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 @ 7521

Last change on this file since 7521 was 7521, checked in by cbricaud, 7 years ago

correc bugs in crs branch

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