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

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

phaze NEMO routines in CRS branch with nemo_v3_6_STABLE branch at rev 7213 (09-09-2016) (merge -r 5519:7213 )

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