New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
crsdom.F90 in branches/2015/dev_r5003_MERCATOR6_CRS/NEMOGCM/NEMO/OPA_SRC/CRS – NEMO

source: branches/2015/dev_r5003_MERCATOR6_CRS/NEMOGCM/NEMO/OPA_SRC/CRS/crsdom.F90 @ 7332

Last change on this file since 7332 was 7332, checked in by cbricaud, 8 years ago

crs branch cleaning

  • Property svn:keywords set to Id
File size: 67.6 KB
Line 
1MODULE crsdom
2   !!===================================================================
3   !!                  ***  crs.F90 ***
4   !!  Purpose: Interface for calculating quantities from a 
5   !!           higher-resolution grid for the coarse grid.
6   !!
7   !!  Method:
8
9   !!  References:  Aumont, O., J.C. Orr, D. Jamous, P. Monfray
10   !!               O. Marti and G. Madec, 1998. A degradation
11   !!               approach to accelerate simulations to steady-state
12   !!               in a 3-D tracer transport model of the global ocean.
13   !!               Climate Dynamics, 14:101-116.
14   !!  History:
15   !!       Original.   May 2012.  (J. Simeon, C. Calone, G. Madec, C. Ethe)
16   !!===================================================================
17
18   USE dom_oce        ! ocean space and time domain and to get jperio
19   USE wrk_nemo       ! work arrays
20   USE crs            ! domain for coarse grid
21   USE in_out_manager 
22   USE par_kind
23   USE crslbclnk
24   USE lib_mpp
25   USE ieee_arithmetic   
26
27   IMPLICIT NONE
28
29   PRIVATE
30
31   PUBLIC crs_dom_ope
32   PUBLIC crs_dom_e3, crs_dom_sfc, crs_dom_msk, crs_dom_hgr, crs_dom_coordinates
33   PUBLIC crs_dom_facvol, crs_dom_def, crs_dom_bat
34
35   INTERFACE crs_dom_ope
36      MODULE PROCEDURE crs_dom_ope_3d, crs_dom_ope_2d
37   END INTERFACE
38
39   REAL(wp),PUBLIC :: r_inf = 1e+7 !cbr 1e+36
40
41   !! Substitutions
42#  include "domzgr_substitute.h90"
43   
44CONTAINS
45
46
47   SUBROUTINE crs_dom_msk
48   !!===================================================================
49   !
50   !
51   !
52   !!===================================================================
53   INTEGER  ::  ji, jj, jk                   ! dummy loop indices
54   INTEGER  ::  ijis,ijie,ijjs,ijje
55   REAL(wp) ::  zmask
56   !!-------------------------------------------------------------------
57     
58   ! Initialize
59   tmask_crs(:,:,:) = 0.0
60   vmask_crs(:,:,:) = 0.0
61   umask_crs(:,:,:) = 0.0
62   fmask_crs(:,:,:) = 0.0
63   !
64   DO jk = 1, jpkm1
65      DO ji = nldi_crs, nlei_crs
66
67         ijis = mis_crs(ji)
68         ijie = mie_crs(ji)
69
70         DO jj = nldj_crs, nlej_crs
71
72            ijjs = mjs_crs(jj)
73            ijje = mje_crs(jj)
74
75            zmask = 0.0
76            zmask = SUM( tmask(ijis:ijie,ijjs:ijje,jk) )
77            IF ( zmask > 0.0 ) tmask_crs(ji,jj,jk) = 1.0
78
79            zmask = 0.0
80            zmask = SUM( vmask(ijis:ijie,ijje     ,jk) )
81            IF ( zmask > 0.0 ) vmask_crs(ji,jj,jk) = 1.0
82
83            zmask = 0.0
84            zmask = SUM( umask(ijie     ,ijjs:ijje,jk) )
85            IF ( zmask > 0.0 ) umask_crs(ji,jj,jk) = 1.0
86
87            fmask_crs(ji,jj,jk) = fmask(ijie,ijje,jk)
88
89         ENDDO
90      ENDDO
91   ENDDO
92
93   CALL crs_lbc_lnk( tmask_crs, 'T', 1.0 )
94   CALL crs_lbc_lnk( vmask_crs, 'V', 1.0 )
95   CALL crs_lbc_lnk( umask_crs, 'U', 1.0 )
96   CALL crs_lbc_lnk( fmask_crs, 'F', 1.0 )
97   !
98   END SUBROUTINE crs_dom_msk
99
100   SUBROUTINE crs_dom_coordinates( p_gphi, p_glam, cd_type, p_gphi_crs, p_glam_crs )
101      !!----------------------------------------------------------------
102      !!               *** SUBROUTINE crs_coordinates ***
103      !! ** Purpose :  Determine the coordinates for the coarse grid
104      !!
105      !! ** Method  :  From the parent grid subset, search for the central
106      !!               point.  For an odd-numbered reduction factor,
107      !!               the coordinate will be that of the central T-cell.
108      !!               For an even-numbered reduction factor, of a non-square
109      !!               coarse grid box, the coordinate will be that of
110      !!               the east or north face or more likely.  For a square
111      !!               coarse grid box, the coordinate will be that of
112      !!               the central f-corner.
113      !!
114      !! ** Input   :  p_gphi = parent grid gphi[t|u|v|f]
115      !!               p_glam = parent grid glam[t|u|v|f]
116      !!               cd_type  = grid type (T,U,V,F)
117      !! ** Output  :  p_gphi_crs = coarse grid gphi[t|u|v|f]
118      !!               p_glam_crs = coarse grid glam[t|u|v|f]
119      !!             
120      !! History. 1 Jun.
121      !!----------------------------------------------------------------
122      !! Arguments
123      REAL(wp), DIMENSION(jpi,jpj)        , INTENT(in)  :: p_gphi  ! Parent grid latitude
124      REAL(wp), DIMENSION(jpi,jpj)        , INTENT(in)  :: p_glam  ! Parent grid longitude
125      CHARACTER(len=1),                     INTENT(in)  :: cd_type   ! grid type (T,U,V,F)
126      REAL(wp), DIMENSION(jpi_crs,jpj_crs), INTENT(out) :: p_gphi_crs  ! Coarse grid latitude
127      REAL(wp), DIMENSION(jpi_crs,jpj_crs), INTENT(out) :: p_glam_crs  ! Coarse grid longitude
128
129      !! Local variables
130      INTEGER :: ji, jj, jk                   ! dummy loop indices
131      INTEGER :: iji, ijj
132      INTEGER  :: ir,jr
133      !!----------------------------------------------------------------
134      p_gphi_crs(:,:)=0._wp
135      p_glam_crs(:,:)=0._wp
136
137 
138      SELECT CASE ( cd_type )
139         CASE ( 'T' )
140            DO jj =  nldj_crs, nlej_crs
141               ijj = mjs_crs(jj) + 1
142               DO ji = nldi_crs, nlei_crs
143                  iji = mis_crs(ji) + 1
144                  IF( ijj .GT. jpj )WRITE(narea+8000-1,*)"BUG ijj ",jj,mjs_crs(jj);CALL FLUSH(narea+8000-1)
145                  p_gphi_crs(ji,jj) = p_gphi(iji,ijj)
146                  p_glam_crs(ji,jj) = p_glam(iji,ijj)
147               ENDDO
148            ENDDO
149         CASE ( 'U' )
150            DO jj =  nldj_crs, nlej_crs
151               ijj = mjs_crs(jj) + 1
152               DO ji = nldi_crs, nlei_crs
153                  iji = mie_crs(ji)
154                  p_gphi_crs(ji,jj) = p_gphi(iji,ijj)
155                  p_glam_crs(ji,jj) = p_glam(iji,ijj)
156 
157               ENDDO
158            ENDDO
159         CASE ( 'V' )
160            DO jj =  nldj_crs, nlej_crs
161               ijj = mje_crs(jj)
162               DO ji = nldi_crs, nlei_crs
163                  iji = mis_crs(ji) + 1
164                  p_gphi_crs(ji,jj) = p_gphi(iji,ijj)
165                  p_glam_crs(ji,jj) = p_glam(iji,ijj)
166               ENDDO
167            ENDDO
168         CASE ( 'F' )
169            DO jj =  nldj_crs, nlej_crs
170               ijj = mje_crs(jj)
171               DO ji = nldi_crs, nlei_crs
172                  iji = mie_crs(ji)
173                  p_gphi_crs(ji,jj) = p_gphi(iji,ijj)
174                  p_glam_crs(ji,jj) = p_glam(iji,ijj)
175               ENDDO
176            ENDDO
177      END SELECT
178
179      ! Retroactively add back the boundary halo cells.
180      CALL crs_lbc_lnk( p_gphi_crs, cd_type, 1.0 )
181      CALL crs_lbc_lnk( p_glam_crs, cd_type, 1.0 )
182      !
183   END SUBROUTINE crs_dom_coordinates
184
185  SUBROUTINE crs_dom_hgr( p_e1, p_e2, cd_type, p_e1_crs, p_e2_crs )
186      !!----------------------------------------------------------------
187      !!               *** SUBROUTINE crs_dom_hgr ***
188      !!
189      !! ** Purpose :  Get coarse grid horizontal scale factors and unmasked fraction
190      !!
191      !! ** Method  :  For grid types T,U,V,Fthe 2D scale factors of
192      !!               the coarse grid are the sum of the east or north faces of the
193      !!               parent grid subset comprising the coarse grid box.     
194      !!               - e1,e2 Scale factors
195      !!                 Valid arguments:
196      !! ** Inputs  : p_e1, p_e2  = parent grid e1 or e2 (t,u,v,f)
197      !!              cd_type     = grid type (T,U,V,F) for scale factors; for velocities (U or V)
198      !! ** Outputs : p_e1_crs, p_e2_crs  = parent grid e1 or e2 (t,u,v,f)
199      !!
200      !! History.     4 Jun.  Write for WGT and scale factors only
201      !!----------------------------------------------------------------
202      !!
203      !!  Arguments
204      REAL(wp), DIMENSION(jpi,jpj)        , INTENT(in)  :: p_e1     ! Parent grid U,V scale factors (e1)
205      REAL(wp), DIMENSION(jpi,jpj)        , INTENT(in)  :: p_e2     ! Parent grid U,V scale factors (e2)
206      CHARACTER(len=1)                    , INTENT(in)  :: cd_type  ! grid type U,V
207
208      REAL(wp), DIMENSION(jpi_crs,jpj_crs), INTENT(out) :: p_e1_crs ! Coarse grid box 2D quantity
209      REAL(wp), DIMENSION(jpi_crs,jpj_crs), INTENT(out) :: p_e2_crs ! Coarse grid box 2D quantity
210
211      !! Local variables
212      INTEGER :: ji, jj, jk     ! dummy loop indices
213      INTEGER :: ijis,ijie,ijjs,ijje
214      INTEGER :: ji1, jj1
215 
216      !!---------------------------------------------------------------- 
217      ! Initialize     
218
219         DO ji = nldi_crs, nlei_crs
220
221            ijis = mis_crs(ji)
222            ijie = mie_crs(ji)
223
224            DO jj = nldj_crs, nlej_crs
225
226               ijjs = mjs_crs(jj)
227               ijje = mje_crs(jj)
228
229               ! Only for a factro 3 coarsening
230               SELECT CASE ( cd_type )
231                   CASE ( 'T' )
232                      p_e1_crs(ji,jj) = REAL(nn_factx,wp)*p_e1(ijis+1,ijjs+1)
233                      p_e2_crs(ji,jj) = REAL(nn_facty,wp)*p_e2(ijis+1,ijjs+1)
234                   CASE ( 'U' )
235                      p_e1_crs(ji,jj) = REAL(nn_factx,wp)*p_e1(ijis+1,ijjs+1       ) 
236                      p_e2_crs(ji,jj) = REAL(nn_facty,wp)*p_e2(ijie  ,ijjs+1       ) 
237
238                   CASE ( 'V' )
239                      p_e1_crs(ji,jj) = REAL(nn_factx,wp)*p_e1(ijis+1,ijje       ) 
240                      p_e2_crs(ji,jj) = REAL(nn_facty,wp)*p_e2(ijis+1,ijjs+1     ) 
241                   CASE ( 'F' )
242                      p_e1_crs(ji,jj) = REAL(nn_factx,wp)*p_e1(ijis+1,ijje       ) 
243                      p_e2_crs(ji,jj) = REAL(nn_facty,wp)*p_e2(ijie  ,ijjs+1     ) 
244               END SELECT
245            ENDDO
246         ENDDO
247
248
249      CALL crs_lbc_lnk( p_e1_crs, cd_type, 1.0 )
250      CALL crs_lbc_lnk( p_e2_crs, cd_type, 1.0 )
251
252   END SUBROUTINE crs_dom_hgr
253
254
255   SUBROUTINE crs_dom_facvol( p_mask, cd_type, p_e1, p_e2, p_e3, p_fld1_crs, p_fld2_crs )
256      !!----------------------------------------------------------------
257      !!               *** SUBROUTINE crsfun_wgt ***
258      !! ** Purpose :  Three applications.
259      !!               1) SUM. Get coarse grid horizontal scale factors and unmasked fraction
260      !!               2) VOL. Get coarse grid box volumes
261      !!               3) WGT. Weighting multiplier for volume-weighted and/or
262      !!                       area-weighted averages.
263      !!                       Weights (i.e. the denominator) calculated here
264      !!                       to avoid IF-tests and division.
265      !! ** Method  :  1) SUM.  For grid types T,U,V,F (and W) the 2D scale factors of
266      !!               the coarse grid are the sum of the east or north faces of the
267      !!               parent grid subset comprising the coarse grid box. 
268      !!               The fractions of masked:total surface (3D) on the east,
269      !!               north and top faces is, optionally, also output.
270      !!               - Top face area sum
271      !!                 Valid arguments: cd_type, cd_op='W', p_pmask, p_e1, p_e2
272      !!               - Top face ocean surface fraction
273      !!                 Valid arguments: cd_type, cd_op='W', p_pmask, p_e1, p_e2       
274      !!               - e1,e2 Scale factors
275      !!                 Valid arguments:
276      !!               2) VOL.  For grid types W and T, the coarse grid box
277      !!               volumes are output. Also optionally, the fraction of 
278      !!               masked:total volume of the parent grid subset is output (i.e. facvol).
279      !!               3) WGT. Based on the grid type, the denominator is pre-determined here to 
280      !!               perform area- or volume- weighted averages,
281      !!               to avoid IF-tests and divisions.
282      !! ** Inputs  : p_e1, p_e2  = parent grid e1 or e2 (t,u,v,f)
283      !!              p_pmask     = parent grid mask (T,U,V,F)
284      !!              cd_type     = grid type (T,U,V,F) for scale factors; for velocities (U or V)
285      !!              cd_op       = applied operation (SUM, VOL, WGT)
286      !!              p_fse3      = (Optional) parent grid vertical level thickness (fse3u or fse3v)
287      !! ** Outputs : p_cfield2d_1 = (Optional) 2D field on coarse grid
288      !!              p_cfield2d_2 = (Optional) 2D field on coarse grid
289      !!              p_cfield3d_1 = (Optional) 3D field on coarse grid
290      !!              p_cfield3d_2 = (Optional) 3D field on coarse grid
291      !!
292      !! History.     4 Jun.  Write for WGT and scale factors only
293      !!----------------------------------------------------------------
294      !!
295      !!  Arguments
296      CHARACTER(len=1),                         INTENT(in)  :: cd_type  ! grid type U,V
297      REAL(wp), DIMENSION(jpi,jpj,jpk)        , INTENT(in)  :: p_mask  ! Parent grid U,V mask
298      REAL(wp), DIMENSION(jpi,jpj)            , INTENT(in)  :: p_e1     ! Parent grid U,V scale factors (e1)
299      REAL(wp), DIMENSION(jpi,jpj)            , INTENT(in)  :: p_e2     ! Parent grid U,V scale factors (e2)
300      REAL(wp), DIMENSION(jpi,jpj,jpk)        , INTENT(in)  :: p_e3     ! Parent grid vertical level thickness (fse3u, fse3v)
301
302      REAL(wp), DIMENSION(jpi_crs,jpj_crs,jpk), INTENT(out) :: p_fld1_crs ! Coarse grid box 3D quantity
303      REAL(wp), DIMENSION(jpi_crs,jpj_crs,jpk), INTENT(out) :: p_fld2_crs ! Coarse grid box 3D quantity
304
305      !! Local variables
306      REAL(wp)                                :: zdAm
307      INTEGER                                 :: ji, jj, jk
308      INTEGER :: ijis,ijie,ijjs,ijje
309
310      REAL(wp), DIMENSION(:,:,:), POINTER     :: zvol, zmask     
311      !!---------------------------------------------------------------- 
312   
313      CALL wrk_alloc( jpi, jpj, jpk, zvol, zmask )
314
315      p_fld1_crs(:,:,:) = 0.0
316      p_fld2_crs(:,:,:) = 0.0
317
318      DO jk = 1, jpk
319         zvol (:,:,jk) = p_e1(:,:) * p_e2(:,:) * p_e3(:,:,jk) 
320         zmask(:,:,jk) = p_mask(:,:,jk) 
321      ENDDO
322
323      DO jk = 1, jpk
324         DO ji = nldi_crs, nlei_crs
325
326            ijis = mis_crs(ji)
327            ijie = mie_crs(ji)
328
329            DO jj = nldj_crs, nlej_crs
330
331               ijjs = mjs_crs(jj)
332               ijje = mje_crs(jj)
333
334               p_fld1_crs(ji,jj,jk) =  SUM( zvol(ijis:ijie,ijjs:ijje,jk) )
335               zdAm                 =  SUM( zvol(ijis:ijie,ijjs:ijje,jk) * zmask(ijis:ijie,ijjs:ijje,jk) )
336               p_fld2_crs(ji,jj,jk) = zdAm / p_fld1_crs(ji,jj,jk) 
337            ENDDO
338         ENDDO
339      ENDDO
340      !                                             !  Retroactively add back the boundary halo cells.
341      CALL crs_lbc_lnk( p_fld1_crs, cd_type, 1.0 ) 
342      CALL crs_lbc_lnk( p_fld2_crs, cd_type, 1.0 ) 
343      !
344      CALL wrk_dealloc( jpi, jpj, jpk, zvol, zmask )
345      !
346   END SUBROUTINE crs_dom_facvol
347
348
349   SUBROUTINE crs_dom_ope_3d( p_fld, cd_op, cd_type, p_mask, p_fld_crs, p_e12, p_e3, p_surf_crs, p_mask_crs, psgn )
350      !!----------------------------------------------------------------
351      !!               *** SUBROUTINE crsfun_UV ***
352      !! ** Purpose :  Average, area-weighted, of U or V on the east and north faces
353      !!
354      !! ** Method  :  The U and V velocities (3D) are determined as the area-weighted averages
355      !!               on the east and north faces, respectively,
356      !!               of the parent grid subset comprising the coarse grid box.
357      !!               In the case of the V and F grid, the last jrow minus 1 is spurious.
358      !! ** Inputs  : p_e1_e2     = parent grid e1 or e2 (t,u,v,f)
359      !!              cd_type     = grid type (T,U,V,F) for scale factors; for velocities (U or V)
360      !!              psgn        = sign change over north fold (See lbclnk.F90)
361      !!              p_pmask     = parent grid mask (T,U,V,F) for scale factors;
362      !!                                       for velocities (U or V)
363      !!              p_fse3      = parent grid vertical level thickness (fse3u or fse3v)
364      !!              p_pfield    = U or V on the parent grid
365      !!              p_surf_crs  = (Optional) Coarse grid weight for averaging
366      !! ** Outputs : p_cfield3d = 3D field on coarse grid
367      !!
368      !! History.  29 May.  completed draft.
369      !!            4 Jun.  Revision for WGT
370      !!            5 Jun.  Streamline for area-weighted average only ; separate scale factor and weights.
371      !!----------------------------------------------------------------
372      !!
373      !!  Arguments
374      REAL(wp), DIMENSION(jpi,jpj,jpk),         INTENT(in)        :: p_fld   ! T, U, V or W on parent grid
375      CHARACTER(len=*),                         INTENT(in)           :: cd_op    ! Operation SUM, MAX or MIN
376      CHARACTER(len=1),                         INTENT(in)           :: cd_type    ! grid type U,V
377      REAL(wp), DIMENSION(jpi,jpj,jpk),         INTENT(in)           :: p_mask    ! Parent grid T,U,V mask
378      REAL(wp), DIMENSION(jpi,jpj),             INTENT(in), OPTIONAL :: p_e12    ! Parent grid T,U,V scale factors (e1 or e2)
379      REAL(wp), DIMENSION(jpi,jpj,jpk),         INTENT(in), OPTIONAL :: p_e3     ! Parent grid vertical level thickness (fse3u, fse3v)
380      REAL(wp), DIMENSION(jpi_crs,jpj_crs,jpk), INTENT(in), OPTIONAL :: p_surf_crs ! Coarse grid area-weighting denominator   
381      REAL(wp), DIMENSION(jpi_crs,jpj_crs,jpk), INTENT(in), OPTIONAL :: p_mask_crs    ! Coarse grid T,U,V maska
382      REAL(wp),                                 INTENT(in)           :: psgn    ! sign
383
384      REAL(wp), DIMENSION(jpi_crs,jpj_crs,jpk), INTENT(out)          :: p_fld_crs ! Coarse grid box 3D quantity
385
386      !! Local variables
387      INTEGER  :: ji, jj, jk 
388      INTEGER  :: ijis, ijie, ijjs, ijje
389      INTEGER  :: ini, inj
390      REAL(wp) :: zflcrs, zsfcrs   
391      REAL(wp), DIMENSION(:,:,:), POINTER :: zsurf, zsurfmsk, zmask,ztabtmp
392      INTEGER  :: ir,jr
393      REAL(wp), DIMENSION(nn_factx,nn_facty):: ztmp
394      REAL(wp), DIMENSION(nn_factx*nn_facty):: ztmp1
395      REAL(wp), DIMENSION(:), ALLOCATABLE   :: ztmp2
396      INTEGER , DIMENSION(1)  :: zdim1
397      REAL(wp) :: zmin,zmax
398      !!---------------------------------------------------------------- 
399   
400      p_fld_crs(:,:,:) = 0.0
401
402      SELECT CASE ( cd_op )
403 
404         CASE ( 'VOL' )
405     
406            CALL wrk_alloc( jpi, jpj, jpk, zsurf, zsurfmsk )
407         
408            SELECT CASE ( cd_type )
409           
410               CASE( 'T', 'W' )
411                  DO jk = 1, jpk
412                     zsurf   (:,:,jk) =  p_e12(:,:) * p_e3(:,:,jk) *  p_mask(:,:,jk) 
413                     zsurfmsk(:,:,jk) =  zsurf(:,:,jk) 
414                  ENDDO
415                  !
416                  DO jk = 1, jpk         
417                     DO jj  = nldj_crs,nlej_crs
418                        ijjs = mjs_crs(jj)
419                        ijje = mje_crs(jj)
420                        DO ji = nldi_crs, nlei_crs
421
422                           ijis = mis_crs(ji)
423                           ijie = mie_crs(ji)
424
425                           zflcrs = SUM( p_fld(ijis:ijie,ijjs:ijje,jk) * zsurfmsk(ijis:ijie,ijjs:ijje,jk) )
426                           zsfcrs = SUM(                                 zsurfmsk(ijis:ijie,ijjs:ijje,jk) )
427
428                           p_fld_crs(ji,jj,jk) = zflcrs
429                           IF( zsfcrs /= 0.0 )  p_fld_crs(ji,jj,jk) = zflcrs / zsfcrs
430                        ENDDO     
431                     ENDDO
432                  ENDDO 
433                  !
434               CASE DEFAULT
435                    STOP
436            END SELECT
437
438            CALL wrk_dealloc( jpi, jpj, jpk, zsurf, zsurfmsk )
439
440         CASE ( 'LOGVOL' )
441
442            CALL wrk_alloc( jpi, jpj, jpk, zsurf, zsurfmsk, ztabtmp )
443
444            ztabtmp(:,:,:)=0._wp
445            WHERE(p_fld* p_mask .NE. 0._wp ) ztabtmp =  LOG10(p_fld * p_mask)*p_mask
446            ztabtmp = ztabtmp * p_mask
447
448            SELECT CASE ( cd_type )
449
450               CASE( 'T', 'W' )
451
452                  DO jk = 1, jpk
453                     zsurf   (:,:,jk) =  p_e12(:,:) * p_e3(:,:,jk) *  p_mask(:,:,jk)
454                     zsurfmsk(:,:,jk) =  zsurf(:,:,jk)
455                  ENDDO
456                  !
457                  DO jk = 1, jpk
458                     DO jj  = nldj_crs,nlej_crs
459                        ijjs = mjs_crs(jj)
460                        ijje = mje_crs(jj)
461                        DO ji = nldi_crs, nlei_crs
462                           ijis = mis_crs(ji)
463                           ijie = mie_crs(ji)
464                           zflcrs = SUM( ztabtmp(ijis:ijie,ijjs:ijje,jk) * zsurfmsk(ijis:ijie,ijjs:ijje,jk) )
465                           zsfcrs = SUM(                                   zsurfmsk(ijis:ijie,ijjs:ijje,jk) )
466                           p_fld_crs(ji,jj,jk) = zflcrs
467                           IF( zsfcrs /= 0.0 )  p_fld_crs(ji,jj,jk) = zflcrs / zsfcrs
468                           p_fld_crs(ji,jj,jk) = 10 ** ( p_fld_crs(ji,jj,jk) *  p_mask_crs(ji,jj,jk) ) * p_mask_crs(ji,jj,jk)
469                        ENDDO
470                     ENDDO
471                  ENDDO
472               CASE DEFAULT
473                    STOP
474            END SELECT
475
476            CALL wrk_dealloc( jpi, jpj, jpk, zsurf, zsurfmsk ,ztabtmp )
477
478         CASE ( 'MED' )
479
480            CALL wrk_alloc( jpi, jpj, jpk, zsurf, zsurfmsk )
481
482            SELECT CASE ( cd_type )
483
484               CASE( 'T', 'W' )
485                  DO jk = 1, jpk
486                     zsurf   (:,:,jk) =  p_e12(:,:) * p_e3(:,:,jk) *  p_mask(:,:,jk)
487                     zsurfmsk(:,:,jk) =  zsurf(:,:,jk)
488                  ENDDO
489                  !
490                  DO jk = 1, jpk
491
492                     DO jj  = nldj_crs,nlej_crs
493
494                        ijjs = mjs_crs(jj)
495                        ijje = mje_crs(jj)
496                        inj  = ijje-ijjs+1
497
498                        DO ji = nldi_crs, nlei_crs
499                           ijis = mis_crs(ji)
500                           ijie = mie_crs(ji)
501                           ini  = ijie-ijis+1
502
503                           ztmp(1:ini,1:inj)= p_fld(ijis:ijie,ijjs:ijje,jk)
504                           zdim1(1) = nn_factx*nn_facty
505                           ztmp1(:) = RESHAPE( ztmp(:,:) , zdim1 )
506                           CALL PIKSRT(nn_factx*nn_facty,ztmp1)
507
508                           ir=0
509                           jr=1
510                           DO WHILE( jr .LE. nn_factx*nn_facty )
511                              IF( ztmp1(jr) == 0. ) THEN
512                                 ir=jr
513                                 jr=jr+1
514                              ELSE
515                                 EXIT
516                              ENDIF
517                           ENDDO
518                           IF( ir .LE. nn_factx*nn_facty-1 )THEN
519                              ALLOCATE( ztmp2(nn_factx*nn_facty-ir) )
520                              ztmp2(1:nn_factx*nn_facty-ir) = ztmp1(1+ir:nn_factx*nn_facty)
521                              jr=INT( 0.5 * REAL(nn_factx*nn_facty-ir,wp) )+1
522                              p_fld_crs(ji,jj,jk) = ztmp2(jr)
523                              DEALLOCATE( ztmp2 )
524                           ELSE
525                              p_fld_crs(ji,jj,jk) = 0._wp
526                           ENDIF
527
528                        ENDDO
529                     ENDDO
530                  ENDDO
531               CASE DEFAULT
532                    STOP
533            END SELECT
534
535           CALL wrk_dealloc( jpi, jpj, jpk, zsurf, zsurfmsk )
536 
537         CASE ( 'SUM' )
538         
539            CALL wrk_alloc( jpi, jpj, jpk, zsurfmsk )
540
541            IF( PRESENT( p_e3 ) ) THEN
542               DO jk = 1, jpk
543                  zsurfmsk(:,:,jk) =  p_e12(:,:) * p_e3(:,:,jk) * p_mask(:,:,jk) 
544               ENDDO
545            ELSE
546               DO jk = 1, jpk
547                  zsurfmsk(:,:,jk) =  p_e12(:,:) * p_mask(:,:,jk) 
548               ENDDO
549            ENDIF
550
551            SELECT CASE ( cd_type )
552           
553               CASE( 'T', 'W' )
554       
555                  DO jk = 1, jpk
556                     DO jj  = nldj_crs,nlej_crs
557                        ijjs = mjs_crs(jj)
558                        ijje = mje_crs(jj)
559                        DO ji = nldi_crs, nlei_crs
560                           ijis = mis_crs(ji)
561                           ijie = mie_crs(ji)
562
563                           p_fld_crs(ji,jj,jk) = SUM( p_fld(ijis:ijie,ijjs:ijje,jk) * zsurfmsk(ijis:ijie,ijjs:ijje,jk) )
564                        ENDDO
565                     ENDDO
566                  ENDDO
567
568               CASE( 'V' )
569
570
571                  DO jk = 1, jpk
572                     DO jj  = nldj_crs,nlej_crs
573                        ijjs = mjs_crs(jj)
574                        ijje = mje_crs(jj)
575                        DO ji = nldi_crs, nlei_crs
576                           ijis = mis_crs(ji)
577                           ijie = mie_crs(ji)
578
579                           p_fld_crs(ji,jj,jk) = SUM( p_fld(ijis:ijie,ijje,jk) * zsurfmsk(ijis:ijie,ijje,jk) )
580                        ENDDO
581                     ENDDO
582                  ENDDO
583
584               CASE( 'U' )
585
586                  DO jk = 1, jpk
587                     DO jj  = nldj_crs,nlej_crs
588                        ijjs = mjs_crs(jj)
589                        ijje = mje_crs(jj)
590                        DO ji = nldi_crs, nlei_crs
591                           ijis = mis_crs(ji)
592                           ijie = mie_crs(ji)
593
594                           p_fld_crs(ji,jj,jk) = SUM( p_fld(ijie,ijjs:ijje,jk) * zsurfmsk(ijie,ijjs:ijje,jk) )
595                        ENDDO
596                     ENDDO
597                  ENDDO
598
599              END SELECT
600
601              IF( PRESENT( p_surf_crs ) ) THEN
602                 WHERE ( p_surf_crs /= 0.0 ) p_fld_crs(:,:,:) = p_fld_crs(:,:,:) / p_surf_crs(:,:,:)
603              ENDIF
604
605              CALL wrk_dealloc( jpi, jpj, jpk, zsurfmsk )
606
607         CASE ( 'MAX' )    !  search the max of unmasked grid cells
608
609            CALL wrk_alloc( jpi, jpj, jpk, zmask )
610
611            DO jk = 1, jpk
612               zmask(:,:,jk) = p_mask(:,:,jk) 
613            ENDDO
614
615            SELECT CASE ( cd_type )
616           
617               CASE( 'T', 'W' )
618       
619                  DO jk = 1, jpk
620                     DO jj  = nldj_crs,nlej_crs
621                        ijjs = mjs_crs(jj)
622                        ijje = mje_crs(jj)
623                        DO ji = nldi_crs, nlei_crs
624                           ijis = mis_crs(ji)
625                           ijie = mie_crs(ji)
626                           p_fld_crs(ji,jj,jk) = MAXVAL( p_fld(ijis:ijie,ijjs:ijje,jk) * zmask(ijis:ijie,ijjs:ijje,jk) - &
627                                                       & ( ( 1._wp - zmask(ijis:ijie,ijjs:ijje,jk))* r_inf )                )
628                        ENDDO
629                     ENDDO
630                  ENDDO
631 
632               CASE( 'V' )
633                  CALL ctl_stop('MAX operator and V case not available')
634           
635               CASE( 'U' )
636                  CALL ctl_stop('MAX operator and U case not available')
637
638            END SELECT
639
640            CALL wrk_dealloc( jpi, jpj, jpk, zmask )
641
642         CASE ( 'MIN' )      !   Search the min of unmasked grid cells
643
644            CALL wrk_alloc( jpi, jpj, jpk, zmask )
645            DO jk = 1, jpk
646               zmask(:,:,jk) = p_mask(:,:,jk)
647            ENDDO
648
649            SELECT CASE ( cd_type )
650
651               CASE( 'T', 'W' )
652
653                  DO jk = 1, jpk
654                     DO jj  = nldj_crs,nlej_crs
655                        ijjs = mjs_crs(jj)
656                        ijje = mje_crs(jj)
657                        DO ji = nldi_crs, nlei_crs
658                           ijis = mis_crs(ji)
659                           ijie = mie_crs(ji)
660
661                           p_fld_crs(ji,jj,jk) = MINVAL( p_fld(ijis:ijie,ijjs:ijje,jk) * zmask(ijis:ijie,ijjs:ijje,jk) + &
662                                                       & ( 1._wp - zmask(ijis:ijie,ijjs:ijje,jk)* r_inf )                )
663                        ENDDO
664                     ENDDO
665                  ENDDO
666
667           
668               CASE( 'V' )
669                  CALL ctl_stop('MIN operator and V case not available')
670           
671               CASE( 'U' )
672                  CALL ctl_stop('MIN operator and U case not available')
673         
674            END SELECT
675            !
676            CALL wrk_dealloc( jpi, jpj, jpk, zmask )
677            !
678         END SELECT
679         !
680         CALL crs_lbc_lnk( p_fld_crs, cd_type, psgn )
681         !
682    END SUBROUTINE crs_dom_ope_3d
683
684    SUBROUTINE crs_dom_ope_2d( p_fld, cd_op, cd_type, p_mask, p_fld_crs, p_e12, p_e3, p_surf_crs, p_mask_crs, psgn )
685      !!----------------------------------------------------------------
686      !!               *** SUBROUTINE crsfun_UV ***
687      !! ** Purpose :  Average, area-weighted, of U or V on the east and north faces
688      !!
689      !! ** Method  :  The U and V velocities (3D) are determined as the area-weighted averages
690      !!               on the east and north faces, respectively,
691      !!               of the parent grid subset comprising the coarse grid box.
692      !!               In the case of the V and F grid, the last jrow minus 1 is spurious.
693      !! ** Inputs  : p_e1_e2     = parent grid e1 or e2 (t,u,v,f)
694      !!              cd_type     = grid type (T,U,V,F) for scale factors; for velocities (U or V)
695      !!              psgn        = sign change over north fold (See lbclnk.F90)
696      !!              p_pmask     = parent grid mask (T,U,V,F) for scale factors;
697      !!                                       for velocities (U or V)
698      !!              p_fse3      = parent grid vertical level thickness (fse3u or fse3v)
699      !!              p_pfield    = U or V on the parent grid
700      !!              p_surf_crs  = (Optional) Coarse grid weight for averaging
701      !! ** Outputs : p_cfield3d = 3D field on coarse grid
702      !!
703      !! History.  29 May.  completed draft.
704      !!            4 Jun.  Revision for WGT
705      !!            5 Jun.  Streamline for area-weighted average only ; separate scale factor and weights.
706      !!----------------------------------------------------------------
707      !!
708      !!  Arguments
709      REAL(wp), DIMENSION(jpi,jpj),             INTENT(in)           :: p_fld   ! T, U, V or W on parent grid
710      CHARACTER(len=3),                         INTENT(in)           :: cd_op    ! Operation SUM, MAX or MIN
711      CHARACTER(len=1),                         INTENT(in)           :: cd_type    ! grid type U,V
712      REAL(wp), DIMENSION(jpi,jpj,jpk),         INTENT(in)           :: p_mask    ! Parent grid T,U,V mask
713      REAL(wp), DIMENSION(jpi,jpj),             INTENT(in), OPTIONAL :: p_e12    ! Parent grid T,U,V scale factors (e1 or e2)
714      REAL(wp), DIMENSION(jpi,jpj,jpk),         INTENT(in), OPTIONAL :: p_e3     ! Parent grid vertical level thickness (fse3u, fse3v)
715      REAL(wp), DIMENSION(jpi_crs,jpj_crs)    , INTENT(in), OPTIONAL :: p_surf_crs ! Coarse grid area-weighting denominator   
716      REAL(wp), DIMENSION(jpi_crs,jpj_crs,jpk), INTENT(in), OPTIONAL :: p_mask_crs    ! Coarse grid T,U,V mask
717      REAL(wp),                                 INTENT(in)           :: psgn   
718
719      REAL(wp), DIMENSION(jpi_crs,jpj_crs)    , INTENT(out)          :: p_fld_crs ! Coarse grid box 3D quantity
720
721      !! Local variables
722      INTEGER  :: ji, jj, jk                 ! dummy loop indices
723      INTEGER ::  ijis, ijie, ijjs, ijje
724      REAL(wp) :: zflcrs, zsfcrs   
725      REAL(wp), DIMENSION(:,:), POINTER :: zsurfmsk   
726
727      !!---------------------------------------------------------------- 
728 
729      p_fld_crs(:,:) = 0.0
730
731      SELECT CASE ( cd_op )
732     
733        CASE ( 'VOL' )
734
735            CALL wrk_alloc( jpi, jpj, zsurfmsk )
736            zsurfmsk(:,:) =  p_e12(:,:) * p_e3(:,:,1) * p_mask(:,:,1)
737
738            DO jj  = nldj_crs,nlej_crs
739               ijjs = mjs_crs(jj)
740               ijje = mje_crs(jj)
741               DO ji = nldi_crs, nlei_crs
742                  ijis = mis_crs(ji)
743                  ijie = mie_crs(ji)
744
745                  zflcrs = SUM( p_fld(ijis:ijie,ijjs:ijje) * zsurfmsk(ijis:ijie,ijjs:ijje) )
746                  zsfcrs = SUM(                              zsurfmsk(ijis:ijie,ijjs:ijje) )
747
748                  p_fld_crs(ji,jj) = zflcrs
749                  IF( zsfcrs /= 0.0 )  p_fld_crs(ji,jj) = zflcrs / zsfcrs
750               ENDDO
751            ENDDO
752            CALL wrk_dealloc( jpi, jpj, zsurfmsk )
753            !
754
755         CASE ( 'SUM' )
756         
757            CALL wrk_alloc( jpi, jpj, zsurfmsk )
758            IF( PRESENT( p_e3 ) ) THEN
759               zsurfmsk(:,:) =  p_e12(:,:) * p_e3(:,:,1) * p_mask(:,:,1)
760            ELSE
761               zsurfmsk(:,:) =  p_e12(:,:) * p_mask(:,:,1)
762            ENDIF
763
764            SELECT CASE ( cd_type )
765
766               CASE( 'T', 'W' )
767
768                  DO jj  = nldj_crs,nlej_crs
769                     ijjs = mjs_crs(jj)
770                     ijje = mje_crs(jj)
771                     DO ji = nldi_crs, nlei_crs
772                        ijis = mis_crs(ji)
773                        ijie = mie_crs(ji)
774                        p_fld_crs(ji,jj) = SUM( p_fld(ijis:ijie,ijjs:ijje) * zsurfmsk(ijis:ijie,ijjs:ijje) )
775                     ENDDO
776                  ENDDO
777           
778               CASE( 'V' )
779
780                  DO jj  = nldj_crs,nlej_crs
781                     ijjs = mjs_crs(jj)
782                     ijje = mje_crs(jj)
783                     DO ji = nldi_crs, nlei_crs
784                        ijis = mis_crs(ji)
785                        ijie = mie_crs(ji)
786                        p_fld_crs(ji,jj) = SUM( p_fld(ijis:ijie,ijje) * zsurfmsk(ijis:ijie,ijje) )
787                     ENDDO
788                  ENDDO
789
790               CASE( 'U' )
791
792                  DO jj  = nldj_crs,nlej_crs
793                     ijjs = mjs_crs(jj)
794                     ijje = mje_crs(jj)
795                     DO ji = nldi_crs, nlei_crs
796                        ijis = mis_crs(ji)
797                        ijie = mie_crs(ji)
798                        p_fld_crs(ji,jj) = SUM( p_fld(ijie,ijjs:ijje) * zsurfmsk(ijie,ijjs:ijje) )
799                     ENDDO
800                  ENDDO
801
802              END SELECT
803
804              IF( PRESENT( p_surf_crs ) ) THEN
805                 WHERE ( p_surf_crs /= 0.0 ) p_fld_crs(:,:) = p_fld_crs(:,:) / p_surf_crs(:,:)
806              ENDIF
807
808              CALL wrk_dealloc( jpi, jpj, zsurfmsk )
809
810         CASE ( 'MAX' )
811
812            SELECT CASE ( cd_type )
813           
814               CASE( 'T', 'W' )
815 
816                  DO jj  = nldj_crs,nlej_crs
817                     ijjs = mjs_crs(jj)
818                     ijje = mje_crs(jj)
819                     DO ji = nldi_crs, nlei_crs
820                        ijis = mis_crs(ji)
821                        ijie = mie_crs(ji)
822                        p_fld_crs(ji,jj) = MAXVAL( p_fld(ijis:ijie,ijjs:ijje) * p_mask(ijis:ijie,ijjs:ijje,1) - &
823                                                 & ( 1._wp - p_mask(ijis:ijie,ijjs:ijje,1) )                    )
824                     ENDDO
825                  ENDDO
826           
827               CASE( 'V' )
828                  CALL ctl_stop('MAX operator and V case not available')
829           
830               CASE( 'U' )
831                  CALL ctl_stop('MAX operator and U case not available')
832
833              END SELECT
834
835         CASE ( 'MIN' )      !   Search the min of unmasked grid cells
836
837           SELECT CASE ( cd_type )
838
839              CASE( 'T', 'W' )
840
841                  DO jj  = nldj_crs,nlej_crs
842                     ijjs = mjs_crs(jj)
843                     ijje = mje_crs(jj)
844                     DO ji = nldi_crs, nlei_crs
845                        ijis = mis_crs(ji)
846                        ijie = mie_crs(ji)
847                        p_fld_crs(ji,jj) = MINVAL( p_fld(ijis:ijie,ijjs:ijje) * p_mask(ijis:ijie,ijjs:ijje,1) + &
848                                                 & ( 1._wp - p_mask(ijis:ijie,ijjs:ijje,1) )                    )
849                     ENDDO
850                  ENDDO
851           
852               CASE( 'V' )
853                  CALL ctl_stop('MIN operator and V case not available')
854           
855               CASE( 'U' )
856                  CALL ctl_stop('MIN operator and U case not available')
857
858              END SELECT
859             !
860         END SELECT
861         !
862         CALL crs_lbc_lnk( p_fld_crs, cd_type, psgn )
863         !
864   END SUBROUTINE crs_dom_ope_2d
865
866   SUBROUTINE crs_dom_e3( p_e1, p_e2, p_e3, p_sfc_2d_crs,  p_sfc_3d_crs, cd_type, p_mask, p_e3_crs, p_e3_max_crs)
867      !!---------------------------------------------------------------- 
868      !!
869      !!
870      !!
871      !!
872      !!----------------------------------------------------------------
873      !!  Arguments
874      CHARACTER(len=1),                         INTENT(in)          :: cd_type           ! grid type T, W ( U, V, F)
875      REAL(wp), DIMENSION(jpi,jpj,jpk),         INTENT(in)          :: p_mask            ! Parent grid T mask
876      REAL(wp), DIMENSION(jpi,jpj)    ,         INTENT(in)          :: p_e1, p_e2        ! 2D tracer T or W on parent grid
877      REAL(wp), DIMENSION(jpi,jpj,jpk),         INTENT(in)          :: p_e3              ! 3D tracer T or W on parent grid
878      REAL(wp), DIMENSION(jpi_crs,jpj_crs)    , INTENT(in),OPTIONAL :: p_sfc_2d_crs      ! Coarse grid box east or north face quantity
879      REAL(wp), DIMENSION(jpi_crs,jpj_crs,jpk), INTENT(in),OPTIONAL :: p_sfc_3d_crs      ! Coarse grid box east or north face quantity
880      REAL(wp), DIMENSION(jpi_crs,jpj_crs,jpk), INTENT(inout)       :: p_e3_crs          ! Coarse grid box east or north face quantity
881      REAL(wp), DIMENSION(jpi_crs,jpj_crs,jpk), INTENT(inout)       :: p_e3_max_crs      ! Coarse grid box east or north face quantity
882
883      !! Local variables
884      INTEGER ::  ji, jj, jk                   ! dummy loop indices
885      INTEGER ::  ijis, ijie, ijjs, ijje 
886      REAL(wp) :: ze3crs 
887
888      !!---------------------------------------------------------------- 
889      p_e3_crs    (:,:,:) = 0._wp
890      p_e3_max_crs(:,:,:) = 0._wp
891   
892
893      SELECT CASE ( cd_type )
894
895         CASE ('T')
896
897            DO jk = 1, jpk
898               DO ji = nldi_crs, nlei_crs
899
900                  ijis = mis_crs(ji)
901                  ijie = mie_crs(ji)
902
903                  DO jj = nldj_crs, nlej_crs
904
905                     ijjs = mjs_crs(jj)
906                     ijje = mje_crs(jj)
907
908                     p_e3_max_crs(ji,jj,jk) = MAXVAL( p_e3(ijis:ijie,ijjs:ijje,jk) * p_mask(ijis:ijie,ijjs:ijje,jk) )
909
910                     ze3crs = SUM( p_e1(ijis:ijie,ijjs:ijje) * p_e2(ijis:ijie,ijjs:ijje) * p_e3(ijis:ijie,ijjs:ijje,jk) * p_mask(ijis:ijie,ijjs:ijje,jk) )
911                     IF( p_sfc_3d_crs(ji,jj,jk) .NE. 0._wp )p_e3_crs(ji,jj,jk) = ze3crs / p_sfc_3d_crs(ji,jj,jk)
912
913                  ENDDO
914               ENDDO
915            ENDDO
916
917         CASE ('U')
918
919            DO jk = 1, jpk
920               DO ji = nldi_crs, nlei_crs
921
922                  ijis = mis_crs(ji)
923                  ijie = mie_crs(ji)
924
925                  DO jj = nldj_crs, nlej_crs
926
927                     ijjs = mjs_crs(jj)
928                     ijje = mje_crs(jj)
929
930                     p_e3_max_crs(ji,jj,jk) = MAXVAL( p_e3(ijie,ijjs:ijje,jk) * p_mask(ijie,ijjs:ijje,jk) )
931
932                     ze3crs = SUM( p_e2(ijie,ijjs:ijje) * p_e3(ijie,ijjs:ijje,jk) * p_mask(ijie,ijjs:ijje,jk) )
933                     IF( p_sfc_2d_crs(ji,jj) .NE. 0._wp )p_e3_crs(ji,jj,jk) = ze3crs / p_sfc_2d_crs(ji,jj)
934                  ENDDO
935               ENDDO
936            ENDDO
937
938         CASE ('V')
939
940            DO jk = 1, jpk
941               DO ji = nldi_crs, nlei_crs
942
943                  ijis = mis_crs(ji)
944                  ijie = mie_crs(ji)
945
946                  DO jj = nldj_crs, nlej_crs
947
948                     ijjs = mjs_crs(jj)
949                     ijje = mje_crs(jj)
950
951                     p_e3_max_crs(ji,jj,jk) = MAXVAL( p_e3(ijis:ijie,ijje,jk) * p_mask(ijis:ijie,ijje,jk) )
952
953                     ze3crs = SUM( p_e1(ijis:ijie,ijje) * p_e3(ijis:ijie,ijje,jk) * p_mask(ijis:ijie,ijje,jk) )
954                     IF( p_sfc_2d_crs(ji,jj) .NE. 0._wp )p_e3_crs(ji,jj,jk) = ze3crs / p_sfc_2d_crs(ji,jj)
955
956                  ENDDO
957               ENDDO
958            ENDDO
959
960         CASE ('W')
961
962            DO jk = 1, jpk
963               DO ji = nldi_crs, nlei_crs
964
965                  ijis = mis_crs(ji)
966                  ijie = mie_crs(ji)
967
968                  DO jj = nldj_crs, nlej_crs
969
970                     ijjs = mjs_crs(jj)
971                     ijje = mje_crs(jj)
972
973                     p_e3_max_crs(ji,jj,jk) = MAXVAL( p_e3(ijis:ijie,ijjs:ijje,jk) * p_mask(ijis:ijie,ijjs:ijje,jk) )
974
975                     ze3crs = SUM( p_e1(ijis:ijie,ijjs:ijje) * p_e2(ijis:ijie,ijjs:ijje) * p_e3(ijis:ijie,ijjs:ijje,jk) * p_mask(ijis:ijie,ijjs:ijje,jk) )
976                     IF( p_sfc_3d_crs(ji,jj,jk) .NE. 0._wp )p_e3_crs(ji,jj,jk) = ze3crs / p_sfc_3d_crs(ji,jj,jk)
977
978                  ENDDO
979               ENDDO
980            ENDDO
981
982      END SELECT
983
984      CALL crs_lbc_lnk( p_e3_max_crs, cd_type, 1.0, pval=1.0 )
985      CALL crs_lbc_lnk( p_e3_crs    , cd_type, 1.0, pval=1.0 )
986
987   END SUBROUTINE crs_dom_e3
988
989   SUBROUTINE crs_dom_sfc(p_mask, cd_type, p_surf_crs, p_surf_crs_msk,  p_e1, p_e2, p_e3 )
990      !!=========================================================================================
991      !!
992      !!
993      !!=========================================================================================
994      !!  Arguments
995      CHARACTER(len=1),                         INTENT(in)           :: cd_type      ! grid type T, W ( U, V, F)
996      REAL(wp), DIMENSION(jpi,jpj,jpk)        , INTENT(in)           :: p_mask       ! Parent grid T mask
997      REAL(wp), DIMENSION(jpi,jpj)            , INTENT(in), OPTIONAL :: p_e1, p_e2         ! 3D tracer T or W on parent grid
998      REAL(wp), DIMENSION(jpi,jpj,jpk)        , INTENT(in), OPTIONAL :: p_e3         ! 3D tracer T or W on parent grid
999      REAL(wp), DIMENSION(jpi_crs,jpj_crs,jpk), INTENT(out)          :: p_surf_crs ! Coarse grid box east or north face quantity
1000      REAL(wp), DIMENSION(jpi_crs,jpj_crs,jpk), INTENT(out)          :: p_surf_crs_msk ! Coarse grid box east or north face quantity
1001
1002      !! Local variables
1003      INTEGER  :: ji, jj, jk                   ! dummy loop indices
1004      INTEGER  :: ijis,ijie,ijjs,ijje
1005      REAL(wp), DIMENSION(:,:,:), POINTER :: zsurf, zsurfmsk   
1006      !!---------------------------------------------------------------- 
1007      ! Initialize
1008      p_surf_crs(:,:,:)=0._wp
1009      p_surf_crs_msk(:,:,:)=0._wp
1010
1011      CALL wrk_alloc( jpi, jpj, jpk, zsurf, zsurfmsk )
1012      !
1013      SELECT CASE ( cd_type )
1014     
1015         CASE ('W')   
1016            DO jk = 1, jpk
1017               zsurf(:,:,jk) = p_e1(:,:) * p_e2(:,:) 
1018            ENDDO
1019
1020         CASE ('V')     
1021            DO jk = 1, jpk
1022               zsurf(:,:,jk) = p_e1(:,:) * p_e3(:,:,jk) 
1023            ENDDO
1024 
1025         CASE ('U')     
1026            DO jk = 1, jpk
1027               zsurf(:,:,jk) = p_e2(:,:) * p_e3(:,:,jk) 
1028            ENDDO
1029
1030         CASE DEFAULT
1031            DO jk = 1, jpk
1032               zsurf(:,:,jk) = p_e1(:,:) * p_e2(:,:) 
1033            ENDDO
1034      END SELECT
1035
1036      DO jk = 1, jpk
1037         zsurfmsk(:,:,jk) = zsurf(:,:,jk) * p_mask(:,:,jk)
1038      ENDDO
1039
1040      SELECT CASE ( cd_type )
1041
1042         CASE ('W')
1043
1044            DO jk = 1, jpk
1045               DO jj = nldj_crs,nlej_crs
1046                  ijjs=mjs_crs(jj)
1047                  ijje=mje_crs(jj)
1048                  DO ji = nldi_crs,nlei_crs
1049                     ijis=mis_crs(ji)
1050                     ijie=mie_crs(ji)
1051                     p_surf_crs    (ji,jj,jk) =  SUM(zsurf   (ijis:ijie,ijjs:ijje,jk) )
1052                     p_surf_crs_msk(ji,jj,jk) =  SUM(zsurfmsk(ijis:ijie,ijjs:ijje,jk) )
1053                  ENDDO     
1054               ENDDO
1055            ENDDO   
1056
1057         CASE ('U')
1058
1059            DO jk = 1, jpk
1060               DO jj = nldj_crs,nlej_crs
1061                  ijjs=mjs_crs(jj)
1062                  ijje=mje_crs(jj)
1063                  DO ji = nldi_crs,nlei_crs
1064                     ijis=mis_crs(ji)
1065                     ijie=mie_crs(ji)
1066                     p_surf_crs    (ji,jj,jk) =  SUM(zsurf   (ijie,ijjs:ijje,jk) )
1067                     p_surf_crs_msk(ji,jj,jk) =  SUM(zsurfmsk(ijie,ijjs:ijje,jk) )
1068                  ENDDO
1069               ENDDO
1070            ENDDO
1071
1072         CASE ('V')
1073
1074            DO jk = 1, jpk
1075               DO jj = nldj_crs,nlej_crs
1076                  ijjs=mjs_crs(jj)
1077                  ijje=mje_crs(jj)
1078                  DO ji = nldi_crs,nlei_crs
1079                     ijis=mis_crs(ji)
1080                     ijie=mie_crs(ji)
1081                     p_surf_crs    (ji,jj,jk) =  SUM(zsurf   (ijis:ijie,ijje,jk) )
1082                     p_surf_crs_msk(ji,jj,jk) =  SUM(zsurfmsk(ijis:ijie,ijje,jk) )
1083                  ENDDO
1084               ENDDO
1085            ENDDO
1086
1087      END SELECT
1088
1089      CALL crs_lbc_lnk( p_surf_crs    , cd_type, 1.0 )
1090      CALL crs_lbc_lnk( p_surf_crs_msk, cd_type, 1.0 )
1091
1092      CALL wrk_dealloc( jpi, jpj, jpk, zsurf, zsurfmsk )
1093
1094   END SUBROUTINE crs_dom_sfc
1095   
1096   SUBROUTINE crs_dom_def
1097      !!----------------------------------------------------------------
1098      !!               *** SUBROUTINE crs_dom_def ***
1099      !! ** Purpose :  Three applications.
1100      !!               1) Define global domain indice of the croasening grid
1101      !!               2) Define local domain indice of the croasening grid
1102      !!               3) Define the processor domain indice for a croasening grid
1103      !!----------------------------------------------------------------
1104      INTEGER :: ji,jj,jk,ijjgloT,ijis,ijie,ijjs,ijje,jn      ! dummy indices
1105      INTEGER :: ierr                                ! allocation error status
1106      INTEGER :: iproci,iprocj,iproc,iprocno,iprocso,iimppt_crs
1107      INTEGER :: ii_start,ii_end,ij_start,ij_end
1108      !!----------------------------------------------------------------
1109 
1110 
1111      !==============================================================================================
1112      ! Define global and local domain sizes
1113      !==============================================================================================
1114      jpiglo_crs   = INT( (jpiglo - 2) / nn_factx ) + 2
1115      jpjglo_crs   = INT( (jpjglo - MOD(jpjglo, nn_facty)) / nn_facty ) + 2
1116      jpiglo_crsm1 = jpiglo_crs - 1
1117      jpjglo_crsm1 = jpjglo_crs - 1 
1118
1119      jpi_crs = ( jpiglo_crs   - 2 * jpreci + (jpni-1) ) / jpni + 2 * jpreci
1120      jpj_crs = ( jpjglo_crsm1 - 2 * jprecj + (jpnj-1) ) / jpnj + 2 * jprecj
1121       
1122      jpi_crsm1   = jpi_crs - 1
1123      jpj_crsm1   = jpj_crs - 1
1124
1125      npolj_crs   = npolj
1126     
1127      ierr = crs_dom_alloc()          ! allocate most coarse grid arrays
1128
1129      !==============================================================================================
1130      ! Define processor domain indices
1131      !==============================================================================================
1132      IF( .NOT. lk_mpp ) THEN
1133
1134         nimpp_crs  = 1
1135         njmpp_crs  = 1
1136         nlci_crs   = jpi_crs
1137         nlcj_crs   = jpj_crs
1138         nldi_crs   = 1
1139         nldj_crs   = 1
1140         nlei_crs   = jpi_crs
1141         nlej_crs   = jpj_crs
1142
1143      ELSE
1144
1145         nimpp_crs  = 1
1146         njmpp_crs  = 1
1147         nlci_crs   = jpi_crs
1148         nlcj_crs   = jpj_crs
1149         nldi_crs   = 1
1150         nldj_crs   = 1
1151         nlei_crs   = jpi_crs
1152         nlej_crs   = jpj_crs
1153
1154         !==============================================================================================
1155         ! mpp_ini2
1156         !==============================================================================================
1157
1158         !order of local domain in i and j directions
1159         DO ji = 1 , jpni
1160            DO jj = 1 ,jpnj
1161               IF( nfipproc(ji,jj)  == narea-1 )THEN
1162                  iproci=ji
1163                  iprocj=jj
1164               ENDIF
1165            ENDDO
1166         ENDDO
1167
1168         !==========================================================================
1169         ! check
1170         !==========================================================================
1171         !WRITE(narea+8000-1,*)"nfipproc(ji,jj),narea :",nfipproc(iproci,iprocj),narea
1172         !WRITE(narea+8000-1,*)"proc i,j ",iproci,iprocj
1173         !WRITE(narea+8000-1,*)"jpni  jpnj jpnij ",jpni,jpnj,jpnij
1174         !WRITE(narea+8000-1,*)"nperio jperio ",nperio,jperio
1175         !WRITE(narea+8000-1,*)"nowe noea",nowe,noea
1176         !WRITE(narea+8000-1,*)"noso nono",noso,nono
1177         !WRITE(narea+8000-1,*)"nbondi nbondj ",nbondi,nbondj
1178         !WRITE(narea+8000-1,*)"jpiglo jpjglo ",jpiglo,jpjglo
1179         !WRITE(narea+8000-1,*)"jpi jpj ",jpi,jpj
1180         !WRITE(narea+8000-1,*)"nbondi nbondj",nbondi,nbondj
1181         !WRITE(narea+8000-1,*)"nimpp njmpp ",nimpp,njmpp
1182         !WRITE(narea+8000-1,*)"loc jpi nldi,nlei,nlci ",jpi, nldi        ,nlei         ,nlci
1183         !WRITE(narea+8000-1,*)"glo jpi nldi,nlei      ",jpi, nldi+nimpp-1,nlei+nimpp-1
1184         !WRITE(narea+8000-1,*)"loc jpj nldj,nlej,nlcj ",jpj, nldj        ,nlej         ,nlcj
1185         !WRITE(narea+8000-1,*)"glo jpj nldj,nlej      ",jpj, nldj+njmpp-1,nlej+njmpp-1
1186         !WRITE(narea+8000-1,*)"jpiglo_crs jpjglo_crs ",jpiglo_crs,jpjglo_crs
1187         !WRITE(narea+8000-1,*)"jpi_crs jpj_crs ",jpi_crs,jpj_crs
1188         !WRITE(narea+8000-1,*)"glamt gphit ",glamt(1,1),gphit(jpi,jpj),glamt(jpi,jpj),gphit(jpi,jpj)
1189         !WRITE(narea+8000-1,*)"min max tmask ",MINVAL(tmask),MAXVAL(tmask)
1190
1191         !==========================================================================
1192         ! coarsened domain: dimensions along I
1193         !==========================================================================
1194
1195         !------------------------------------------------------------------------------------
1196         !I-1 fill mis2_crs and mie2_crs: arrays to switch from physic grid to coarsened grid
1197         !------------------------------------------------------------------------------------
1198
1199         ! !--------!--------!--------!
1200         ! !        !        !        !
1201         ! !        !        !        !
1202         ! !        !        !        !
1203         ! !--------!--------!--------!
1204         ! !        !        !        !
1205         ! !        ! ji,jj  !        !
1206         ! !        !        !        !
1207         ! !--------!--------!--------!
1208         ! !        !        !        !
1209         ! !        !        !        !
1210         ! !        !        !        !
1211         ! !--------!--------!--------!
1212         !  mis2_crs(ji)      mie2_crs(ji)
1213       
1214
1215         SELECT CASE ( jperio )
1216
1217         CASE ( 0, 1 )
1218
1219            DO ji=1,jpiglo_crs
1220               ijis=nn_factx*(ji-1)+1
1221               ijie=nn_factx*(ji-1)+3
1222               mis2_crs(ji)=ijis
1223               mie2_crs(ji)=ijie
1224            ENDDO
1225
1226         CASE ( 3, 4 )    !   3, 4 : T-Pivot at North Fold: make correspondinf the pivot points of the 2 grids
1227
1228            DO ji=1,jpiglo_crs
1229               ijis=nn_factx*(ji-1)-2
1230               ijie=nn_factx*(ji-1)
1231               mis2_crs(ji)=ijis
1232               mie2_crs(ji)=ijie
1233            ENDDO
1234
1235         CASE DEFAULT
1236            WRITE(numout,*) 'crs_init. Only jperio = 0, 1, 3, 4 supported; narea: ',narea
1237         END SELECT
1238
1239         !-------------------------------------------------------------------------------
1240         ! I-2 find the first CRS cell which is inside the physic grid inner domain
1241         !-------------------------------------------------------------------------------
1242         ! ijis           : global indice of the first CRS cell which inside the physic grid inner domain
1243         ! mis2_crs(ijis) : global indice of the bottom-left physic cell corresponding to ijis cell
1244         ! ii_start       : local  ndice of the bottom-left physic cell corresponding to ijis cell
1245
1246         ji=1
1247         DO WHILE( mis2_crs(ji) - nimpp + 1 .LT. 1 ) 
1248            ji=ji+1
1249            IF( ji==jpiglo_crs )EXIT
1250         END DO
1251
1252         ijis=ji
1253         ii_start = mis2_crs(ijis)-nimpp+1
1254
1255         !----------------------------------------------------------------------------------------------
1256         ! I-3 compute nldi_crs and compute mis2_crs and mie2_crs for the first cell of the local domain
1257         !---------------------------------------------------------------------------------------------
1258         nldi_crs = 2
1259         IF( nowe == -1 .AND. ( (jperio==3 .OR. jperio==4 ) .OR. ( (jperio==0 .OR. jperio==1 ) .AND. iproci .NE. 1 )) )THEN
1260
1261            mie2_crs(ijis-1) = mis2_crs(ijis)-1
1262             
1263            SELECT CASE(ii_start)
1264               CASE(1)
1265                  nldi_crs=2
1266                  mie2_crs(ijis-1) = -1
1267                  mis2_crs(ijis-1) = -1
1268               CASE(2)
1269                  nldi_crs=2
1270                  mis2_crs(ijis-1) = mie2_crs(ijis-1)
1271               CASE(3)
1272                  nldi_crs=2
1273                  mis2_crs(ijis-1) = mie2_crs(ijis-1) -1
1274               CASE DEFAULT
1275                  WRITE(narea+8000-1,*)"WRONG VALUE FOR iistart ",ii_start
1276            END SELECT
1277
1278         ENDIF
1279
1280         !----------------------------------------------------------------------------------------------
1281         ! I-4 compute nimpp_crs
1282         !---------------------------------------------------------------------------------------------
1283         nimpp_crs = ijis-1
1284         IF( nimpp==1 )nimpp_crs=1
1285
1286         !-------------------------------------------------------------------------------
1287         ! I-5 find the last CRS cell which is inside the physic grid inner domain
1288         !-------------------------------------------------------------------------------
1289         ! ijie           : global indice of the last CRS cell which inside the physic grid inner domain
1290
1291         ji=jpiglo_crs
1292         DO WHILE( mie2_crs(ji) - nimpp + 1 .GT. jpi )
1293            ji=ji-1
1294            IF( ji==1 )EXIT
1295         END DO
1296         ijie=ji
1297
1298         !-------------------------------------------------------------------------------
1299         ! I-6 compute nlei_crs and nlci_crs
1300         !-------------------------------------------------------------------------------
1301         nlei_crs=ijie-nimpp_crs+1
1302         nlci_crs=nlei_crs+jpreci
1303         IF( iproci == jpni )nlei_crs=nlci_crs
1304
1305         !-------------------------------------------------------------------------------
1306         ! I-7 local to global and global to local indices for CRS grid
1307         !-------------------------------------------------------------------------------
1308         DO ji = 1, jpi_crs
1309            mig_crs(ji) = ji + nimpp_crs - 1
1310         ENDDO
1311         DO ji = 1, jpiglo_crs
1312            mi0_crs(ji) = MAX( 1, MIN( ji - nimpp_crs + 1 , jpi_crs + 1 ) )
1313            mi1_crs(ji) = MAX( 0, MIN( ji - nimpp_crs + 1 , jpi_crs     ) )
1314         ENDDO
1315
1316         !---------------------------------------------------------------------------------------
1317         ! I-8 CRS to physic grid: local indices mis_crs and mie_crs and local coarsening factor
1318         !---------------------------------------------------------------------------------------
1319         DO ji = 1, nlei_crs
1320            mis_crs(ji) = mis2_crs(mig_crs(ji)) - nimpp + 1
1321            mie_crs(ji) = mie2_crs(mig_crs(ji)) - nimpp + 1
1322            IF( iproci == jpni  .AND. ji == nlei_crs )THEN
1323               mie_crs(ji) = nlei
1324               mie2_crs(mig_crs(ji)) = mig(nlei)
1325            ENDIF
1326            nfactx(ji)  = mie_crs(ji)-mis_crs(ji)+1
1327         ENDDO
1328
1329         !==========================================================================
1330         ! coarsened domain: dimensions along J
1331         !==========================================================================
1332
1333         !-----------------------------------------------------------------------------------
1334         !J-1 fill mjs2_crs and mje2_crs: arrays to switch from physic grid to coarsened grid
1335         !-----------------------------------------------------------------------------------
1336
1337         ! !--------!--------!--------!
1338         ! !        !        !        !
1339         ! !        !        !        !
1340         ! !        !        !        ! mje2_crs(jj)
1341         ! !--------!--------!--------!
1342         ! !        !        !        !
1343         ! !        ! ji,jj  !        !
1344         ! !        !        !        !
1345         ! !--------!--------!--------!
1346         ! !        !        !        !
1347         ! !        !        !        ! mjs2_crs(jj)
1348         ! !        !        !        !
1349         ! !--------!--------!--------!
1350
1351         SELECT CASE ( jperio )
1352
1353         CASE ( 0, 1 )    !
1354
1355            DO jj=1,jpjglo_crs
1356               ijjs=nn_facty*(jj-1)+1
1357               ijje=nn_facty*(jj-1)+3
1358               mjs2_crs(jj)=ijjs
1359               mje2_crs(jj)=ijje
1360            ENDDO
1361
1362         CASE ( 3, 4 )    !   3, 4 : T-Pivot at North Fold
1363
1364            DO jj=1,jpjglo_crs
1365               ijjs=nn_facty*(jj)-5
1366               ijje=nn_facty*(jj)-3
1367               mjs2_crs(jj)=ijjs
1368               mje2_crs(jj)=ijje
1369            ENDDO
1370
1371         CASE DEFAULT
1372            WRITE(numout,*) 'crs_init. Only jperio = 0, 1, 3, 4 supported; narea: ',narea
1373         END SELECT
1374
1375         !-------------------------------------------------------------------------------
1376         ! J-2 find the first CRS cell which is inside the physic grid inner domain
1377         !-------------------------------------------------------------------------------
1378         ! ijjs           : global indice of the first CRS cell which inside the physic grid inner domain
1379         ! mis2_crs(ijjs) : global indice of the bottom-left physic cell corresponding to ijis cell
1380         ! ij_start       : local  ndice of the bottom-left physic cell corresponding to ijis cell
1381
1382         jj=1
1383         DO WHILE( mjs2_crs(jj) - njmpp + 1 .LT. 1 )
1384            jj=jj+1
1385            IF( jj==jpjglo_crs )EXIT
1386         END DO
1387
1388         ijjs=jj
1389         ij_start = mjs2_crs(ijjs)-njmpp+1
1390
1391         !----------------------------------------------------------------------------------------------
1392         ! J-3 compute nldj_crs and compute mjs2_crs and mje2_crs for the first cell of the local domain
1393         !---------------------------------------------------------------------------------------------
1394         nldj_crs = 2
1395
1396         IF( jperio==3 .OR. jperio==4 )THEN
1397
1398            IF( noso == -1 )THEN
1399
1400               mje2_crs(ijjs-1) = mjs2_crs(ijjs)-1
1401
1402               SELECT CASE(ij_start)
1403                  CASE(1)
1404                     nldj_crs=2
1405                     mje2_crs(ijjs-1) = -1
1406                     mjs2_crs(ijjs-1) = -1
1407                  CASE(2)
1408                     nldj_crs=2
1409                     mjs2_crs(ijjs-1) = mje2_crs(ijjs-1)
1410                  CASE(3)
1411                     nldj_crs=2
1412                     mjs2_crs(ijjs-1) = mje2_crs(ijjs-1) -1
1413                  CASE DEFAULT
1414                     WRITE(narea+8000-1,*)"WRONG VALUE FOR iistart ",ii_start
1415               END SELECT
1416
1417            ENDIF
1418         ENDIF
1419
1420         !----------------------------------------------------------------------------------------------
1421         ! J-4 compute nimpp_crs
1422         !---------------------------------------------------------------------------------------------
1423         njmpp_crs = ijjs-1
1424         IF( njmpp==1 )njmpp_crs=1
1425
1426         !-------------------------------------------------------------------------------
1427         ! J-5 find the last CRS cell which is inside the physic grid inner domain
1428         !-------------------------------------------------------------------------------
1429         ! ijje           : global indice of the last CRS cell which inside the physic grid inner domain
1430
1431         jj=jpjglo_crs
1432         DO WHILE( mje2_crs(jj) - njmpp + 1 .GT. nlcj )
1433            jj=jj-1
1434            IF( jj==1 )EXIT
1435         END DO
1436         ijje=jj
1437
1438         !-------------------------------------------------------------------------------
1439         ! J-6 compute nlej_crs and nlcj_crs
1440         !-------------------------------------------------------------------------------
1441         nlej_crs=ijje-njmpp_crs+1
1442         nlcj_crs=nlej_crs+jprecj
1443
1444         IF( iprocj == jpnj )THEN
1445            IF( jperio==3 .OR. jperio==4 )THEN
1446               nlej_crs=jpj_crs
1447               nlcj_crs=nlej_crs
1448            ELSE
1449               nlej_crs= nlej_crs+1
1450               nlcj_crs= nlcj_crs+1
1451            ENDIF
1452         ENDIF
1453
1454         !-------------------------------------------------------------------------------
1455         ! J-7 local to global and global to local indices for CRS grid
1456         !-------------------------------------------------------------------------------
1457         DO jj = 1, jpj_crs
1458            mjg_crs(jj) = jj + njmpp_crs - 1
1459         ENDDO
1460         DO jj = 1, jpjglo_crs
1461            mj0_crs(jj) = MAX( 1, MIN( jj - njmpp_crs + 1 , jpj_crs + 1 ) )
1462            mj1_crs(jj) = MAX( 0, MIN( jj - njmpp_crs + 1 , jpj_crs     ) )
1463         ENDDO
1464
1465         !---------------------------------------------------------------------------------------
1466         ! J-8 CRS to physic grid: local indices mis_crs and mie_crs and local coarsening factor
1467         !---------------------------------------------------------------------------------------
1468         DO jj = 1, nlej_crs
1469            mjs_crs(jj) = mjs2_crs(mjg_crs(jj)) - njmpp + 1
1470            mje_crs(jj) = mje2_crs(mjg_crs(jj)) - njmpp + 1
1471            IF( iprocj == jpnj  .AND. jj == nlej_crs )THEN
1472               mjs_crs(jj) = nlej
1473               mjs2_crs(mjg_crs(jj)) = mjg(nlej)
1474               mje_crs(jj) = nlej
1475               mje2_crs(mjg_crs(jj)) = mjg(nlej)
1476            ENDIF
1477            nfacty(jj)  = mje_crs(jj)-mjs_crs(jj)+1
1478         ENDDO
1479
1480         !==========================================================================
1481         ! send local start and end indices to all procs
1482         !==========================================================================
1483
1484         nldit_crs(:)=0 ; nleit_crs(:)=0 ; nlcit_crs(:)=0 ; nimppt_crs(:)=0
1485         nldjt_crs(:)=0 ; nlejt_crs(:)=0 ; nlcjt_crs(:)=0 ; njmppt_crs(:)=0
1486
1487         CALL mppgatheri((/nlci_crs/),0,nlcit_crs) ; CALL mppgatheri((/nlcj_crs/),0,nlcjt_crs) 
1488         CALL mppgatheri((/nldi_crs/),0,nldit_crs) ; CALL mppgatheri((/nldj_crs/),0,nldjt_crs) 
1489         CALL mppgatheri((/nlei_crs/),0,nleit_crs) ; CALL mppgatheri((/nlej_crs/),0,nlejt_crs) 
1490         CALL mppgatheri((/nimpp_crs/),0,nimppt_crs) ; CALL mppgatheri((/njmpp_crs/),0,njmppt_crs) 
1491
1492         DO jj = 1 ,jpnj
1493            DO ji = 1 , jpni
1494               jn=nfipproc(ji,jj)+1
1495               IF( jn .GE. 1 )THEN
1496                  nfiimpp_crs(ji,jj)=nimppt_crs(jn)
1497               ELSE
1498                  nfiimpp_crs(ji,jj) = ANINT( REAL( (nfiimpp(ji,jj) + 1 ) / nn_factx, wp ) ) + 1
1499               ENDIF
1500            ENDDO
1501         ENDDO
1502 
1503         !nogather=T
1504         nfsloop_crs = 1
1505         nfeloop_crs = nlci_crs
1506         DO jn = 2,jpni-1
1507            IF( nfipproc(jn,jpnj) .eq. (narea - 1) )THEN
1508               IF (nfipproc(jn - 1 ,jpnj) .eq. -1 )THEN
1509                  nfsloop_crs = nldi_crs
1510               ENDIF
1511               IF( nfipproc(jn + 1,jpnj) .eq. -1 )THEN
1512                  nfeloop_crs = nlei_crs
1513               ENDIF
1514            ENDIF
1515         END DO
1516
1517         !==========================================================================
1518         ! check
1519         !==========================================================================
1520         !WRITE(narea+8000-1,*)"loc crs jpi nldi,nlei,nlci ",jpi_crs, nldi_crs            ,nlei_crs             ,nlci_crs
1521         !WRITE(narea+8000-1,*)"glo crs jpi nldi,nlei      ",jpi_crs, nldi_crs+nimpp_crs-1,nlei_crs+nimpp_crs-1
1522         !WRITE(narea+8000-1,*)"loc crs jpj nldj,nlej,nlcj ",jpj_crs, nldj_crs            ,nlej_crs             ,nlcj_crs
1523         !WRITE(narea+8000-1,*)"glo crs jpj nldj,nlej      ",jpj_crs, nldj_crs+njmpp_crs-1,nlej_crs+njmpp_crs-1
1524
1525         !==========================================================================
1526         ! Save the parent grid information
1527         !==========================================================================
1528         IF( jpizoom /= 1 .OR. jpjzoom /= 1)    STOP  !cbr mettre un ctlstp et ailleurs ( crsini )
1529         jpi_full    = jpi
1530         jpj_full    = jpj
1531         jpim1_full  = jpim1
1532         jpjm1_full  = jpjm1
1533         npolj_full  = npolj
1534         jpiglo_full = jpiglo
1535         jpjglo_full = jpjglo
1536
1537         nlcj_full   = nlcj
1538         nlci_full   = nlci
1539         nldi_full   = nldi
1540         nldj_full   = nldj
1541         nlei_full   = nlei
1542         nlej_full   = nlej
1543         nimpp_full  = nimpp     
1544         njmpp_full  = njmpp
1545     
1546         nlcit_full(:)  = nlcit(:)
1547         nldit_full(:)  = nldit(:)
1548         nleit_full(:)  = nleit(:)
1549         nimppt_full(:) = nimppt(:)
1550         nlcjt_full(:)  = nlcjt(:)
1551         nldjt_full(:)  = nldjt(:)
1552         nlejt_full(:)  = nlejt(:)
1553         njmppt_full(:) = njmppt(:)
1554     
1555         nfsloop_full = nfsloop
1556         nfeloop_full = nfeloop
1557
1558         nfiimpp_full(:,:) = nfiimpp(:,:) 
1559
1560
1561         !==========================================================================
1562         ! control
1563         !==========================================================================
1564         CALL dom_grid_crs  !swich from mother grid to coarsened grid
1565
1566         IF(lwp) THEN
1567            WRITE(numout,*)
1568            WRITE(numout,*) 'crs_init : coarse grid dimensions'
1569            WRITE(numout,*) '~~~~~~~   coarse domain global j-dimension           jpjglo = ', jpjglo
1570            WRITE(numout,*) '~~~~~~~   coarse domain global i-dimension           jpiglo = ', jpiglo
1571            WRITE(numout,*) '~~~~~~~   coarse domain local  i-dimension              jpi = ', jpi
1572            WRITE(numout,*) '~~~~~~~   coarse domain local  j-dimension              jpj = ', jpj
1573            WRITE(numout,*)
1574            WRITE(numout,*) ' nproc  = '     , nproc
1575            WRITE(numout,*) ' nlci   = '     , nlci
1576            WRITE(numout,*) ' nlcj   = '     , nlcj
1577            WRITE(numout,*) ' nldi   = '     , nldi
1578            WRITE(numout,*) ' nldj   = '     , nldj
1579            WRITE(numout,*) ' nlei   = '     , nlei
1580            WRITE(numout,*) ' nlej   = '     , nlej
1581            WRITE(numout,*) ' nlei_full='    , nlei_full
1582            WRITE(numout,*) ' nldi_full='    , nldi_full
1583            WRITE(numout,*) ' nimpp  = '     , nimpp
1584            WRITE(numout,*) ' njmpp  = '     , njmpp
1585            WRITE(numout,*) ' njmpp_full  = ', njmpp_full
1586            WRITE(numout,*)
1587         ENDIF
1588     
1589         CALL dom_grid_glo ! switch from coarsened grid to mother grid
1590     
1591         nrestx = MOD( nn_factx, 2 )   ! check if even- or odd- numbered reduction factor
1592         nresty = MOD( nn_facty, 2 )
1593
1594         IF( nresty == 0 )THEN
1595            IF( npolj == 3 ) npolj_crs = 5
1596            IF( npolj == 5 ) npolj_crs = 3
1597         ENDIF     
1598     
1599         rfactxy = nn_factx * nn_facty
1600     
1601      ENDIF ! lk_mpp
1602      !
1603      nistr = mis_crs(2)  ;   niend = mis_crs(nlci_crs - 1)
1604      njstr = mjs_crs(3)  ;   njend = mjs_crs(nlcj_crs - 1)
1605      !
1606      !
1607   END SUBROUTINE crs_dom_def
1608   
1609   SUBROUTINE crs_dom_bat
1610      !!----------------------------------------------------------------
1611      !!               *** SUBROUTINE crs_dom_bat ***
1612      !! ** Purpose :  coarsenig bathy
1613      !!----------------------------------------------------------------
1614      !!
1615      !!  local variables
1616      INTEGER  :: ji,jj,jk      ! dummy indices
1617      REAL(wp), DIMENSION(:,:)  , POINTER :: zmbk
1618      !!----------------------------------------------------------------
1619   
1620      CALL wrk_alloc( jpi_crs, jpj_crs, zmbk )
1621   
1622      mbathy_crs(:,:) = jpkm1
1623      mbkt_crs(:,:) = 1
1624      mbku_crs(:,:) = 1
1625      mbkv_crs(:,:) = 1
1626
1627
1628      DO jj = 1, jpj_crs
1629         DO ji = 1, jpi_crs
1630            jk = 0
1631            DO WHILE( tmask_crs(ji,jj,jk+1) > 0.)
1632               jk = jk + 1
1633            ENDDO
1634            mbathy_crs(ji,jj) = float( jk )
1635         ENDDO
1636      ENDDO
1637     
1638      CALL wrk_alloc( jpi_crs, jpj_crs, zmbk )
1639
1640      zmbk(:,:) = 0.0
1641      zmbk(:,:) = REAL( mbathy_crs(:,:), wp ) ;   CALL crs_lbc_lnk(zmbk,'T',1.0)   ;   mbathy_crs(:,:) = INT( zmbk(:,:) )
1642
1643
1644      !
1645      IF(lwp) WRITE(numout,*)
1646      IF(lwp) WRITE(numout,*) '    crsini : mbkt is ocean bottom k-index of T-, U-, V- and W-levels '
1647      IF(lwp) WRITE(numout,*) '    ~~~~~~~~~~~~~'
1648      !
1649      mbkt_crs(:,:) = MAX( mbathy_crs(:,:) , 1 )    ! bottom k-index of T-level (=1 over land)
1650      !                                     ! bottom k-index of W-level = mbkt+1
1651
1652      DO jj = 1, jpj_crsm1                      ! bottom k-index of u- (v-) level
1653         DO ji = 1, jpi_crsm1
1654            mbku_crs(ji,jj) = MIN(  mbkt_crs(ji+1,jj  ) , mbkt_crs(ji,jj)  )
1655            mbkv_crs(ji,jj) = MIN(  mbkt_crs(ji  ,jj+1) , mbkt_crs(ji,jj)  )
1656         END DO
1657      END DO
1658
1659      ! convert into REAL to use lbc_lnk ; impose a min value of 1 as a zero can be set in lbclnk
1660      zmbk(:,:) = 1.e0;   
1661      zmbk(:,:) = REAL( mbku_crs(:,:), wp )   ;   CALL crs_lbc_lnk(zmbk,'U',1.0) ; mbku_crs  (:,:) = MAX( INT( zmbk(:,:) ), 1 ) 
1662      zmbk(:,:) = REAL( mbkv_crs(:,:), wp )   ;   CALL crs_lbc_lnk(zmbk,'V',1.0) ; mbkv_crs  (:,:) = MAX( INT( zmbk(:,:) ), 1 ) 
1663      !
1664      CALL wrk_dealloc( jpi_crs, jpj_crs, zmbk )
1665      !
1666   END SUBROUTINE crs_dom_bat
1667
1668   SUBROUTINE PIKSRT(N,ARR)
1669     INTEGER                  ,INTENT(IN) :: N
1670     REAL(kind=8),DIMENSION(N),INTENT(INOUT) :: ARR
1671
1672     INTEGER      :: i,j
1673     REAL(kind=8) :: a
1674     !!----------------------------------------------------------------
1675
1676     DO j=2, N
1677       a=ARR(j)
1678       DO i=j-1,1,-1
1679          IF(ARR(i)<=a) goto 10
1680          ARR(i+1)=ARR(i)
1681       ENDDO
1682       i=0
168310     ARR(i+1)=a
1684     ENDDO
1685     RETURN
1686
1687   END SUBROUTINE PIKSRT
1688
1689
1690END MODULE crsdom
Note: See TracBrowser for help on using the repository browser.