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.
domain.F90 in branches/2015/dev_r5056_CMCC4_simplification/NEMOGCM/NEMO/OPA_SRC/DOM – NEMO

source: branches/2015/dev_r5056_CMCC4_simplification/NEMOGCM/NEMO/OPA_SRC/DOM/domain.F90 @ 5282

Last change on this file since 5282 was 5282, checked in by diovino, 9 years ago

Dev. branch CMCC4_simplification ticket #1456

  • Property svn:keywords set to Id
File size: 22.3 KB
Line 
1MODULE domain
2   !!==============================================================================
3   !!                       ***  MODULE domain   ***
4   !! Ocean initialization : domain initialization
5   !!==============================================================================
6   !! History :  OPA  !  1990-10  (C. Levy - G. Madec)  Original code
7   !!                 !  1992-01  (M. Imbard) insert time step initialization
8   !!                 !  1996-06  (G. Madec) generalized vertical coordinate
9   !!                 !  1997-02  (G. Madec) creation of domwri.F
10   !!                 !  2001-05  (E.Durand - G. Madec) insert closed sea
11   !!   NEMO     1.0  !  2002-08  (G. Madec)  F90: Free form and module
12   !!            2.0  !  2005-11  (V. Garnier) Surface pressure gradient organization
13   !!            3.3  !  2010-11  (G. Madec)  initialisation in C1D configuration
14   !!            3.6  !  2013     ( J. Simeon, C. Calone, G. Madec, C. Ethe ) Online coarsening of outputs
15   !!----------------------------------------------------------------------
16   
17   !!----------------------------------------------------------------------
18   !!   dom_init       : initialize the space and time domain
19   !!   dom_nam        : read and contral domain namelists
20   !!   dom_ctl        : control print for the ocean domain
21   !!----------------------------------------------------------------------
22   USE oce             ! ocean variables
23   USE dom_oce         ! domain: ocean
24   USE sbc_oce         ! surface boundary condition: ocean
25   USE phycst          ! physical constants
26   USE closea          ! closed seas
27   USE in_out_manager  ! I/O manager
28   USE lib_mpp         ! distributed memory computing library
29
30   USE domhgr          ! domain: set the horizontal mesh
31   USE domzgr          ! domain: set the vertical mesh
32   USE domstp          ! domain: set the time-step
33   USE dommsk          ! domain: set the mask system
34   USE domwri          ! domain: write the meshmask file
35   USE domvvl          ! variable volume
36   USE c1d             ! 1D vertical configuration
37   USE dyncor_c1d      ! Coriolis term (c1d case)         (cor_c1d routine)
38   USE timing          ! Timing
39   USE lbclnk          ! ocean lateral boundary condition (or mpp link)
40
41   IMPLICIT NONE
42   PRIVATE
43
44   PUBLIC   dom_init   ! called by opa.F90
45
46   !! * Substitutions
47#  include "domzgr_substitute.h90"
48   !!-------------------------------------------------------------------------
49   !! NEMO/OPA 3.3 , NEMO Consortium (2010)
50   !! $Id$
51   !! Software governed by the CeCILL licence        (NEMOGCM/NEMO_CeCILL.txt)
52   !!-------------------------------------------------------------------------
53CONTAINS
54
55   SUBROUTINE dom_init
56      !!----------------------------------------------------------------------
57      !!                  ***  ROUTINE dom_init  ***
58      !!                   
59      !! ** Purpose :   Domain initialization. Call the routines that are
60      !!              required to create the arrays which define the space
61      !!              and time domain of the ocean model.
62      !!
63      !! ** Method  : - dom_msk: compute the masks from the bathymetry file
64      !!              - dom_hgr: compute or read the horizontal grid-point position
65      !!                         and scale factors, and the coriolis factor
66      !!              - dom_zgr: define the vertical coordinate and the bathymetry
67      !!              - dom_stp: defined the model time step
68      !!              - dom_wri: create the meshmask file if nmsh=1
69      !!              - 1D configuration, move Coriolis, u and v at T-point
70      !!----------------------------------------------------------------------
71      INTEGER ::   jk          ! dummy loop argument
72      INTEGER ::   iconf = 0   ! local integers
73      !!----------------------------------------------------------------------
74      !
75      IF( nn_timing == 1 )   CALL timing_start('dom_init')
76      !
77      IF(lwp) THEN
78         WRITE(numout,*)
79         WRITE(numout,*) 'dom_init : domain initialization'
80         WRITE(numout,*) '~~~~~~~~'
81      ENDIF
82      !
83                             CALL dom_nam      ! read namelist ( namrun, namdom )
84                             CALL dom_clo      ! Closed seas and lake
85                             CALL dom_hgr      ! Horizontal mesh
86                             CALL dom_zgr      ! Vertical mesh and bathymetry
87                             CALL dom_msk      ! Masks
88      IF( ln_sco )           CALL dom_stiff    ! Maximum stiffness ratio/hydrostatic consistency
89      !
90      ht_0(:,:) = 0.0_wp                       ! Reference ocean depth at T-points
91      hu_0(:,:) = 0.0_wp                       ! Reference ocean depth at U-points
92      hv_0(:,:) = 0.0_wp                       ! Reference ocean depth at V-points
93      DO jk = 1, jpk
94         ht_0(:,:) = ht_0(:,:) + e3t_0(:,:,jk) * tmask(:,:,jk)
95         hu_0(:,:) = hu_0(:,:) + e3u_0(:,:,jk) * umask(:,:,jk)
96         hv_0(:,:) = hv_0(:,:) + e3v_0(:,:,jk) * vmask(:,:,jk)
97      END DO
98      !
99      IF( lk_vvl )           CALL dom_vvl_init ! Vertical variable mesh
100      !
101      IF( lk_c1d         )   CALL cor_c1d      ! 1D configuration: Coriolis set at T-point
102      !
103      !
104      hu(:,:) = 0._wp                          ! Ocean depth at U-points
105      hv(:,:) = 0._wp                          ! Ocean depth at V-points
106      ht(:,:) = 0._wp                          ! Ocean depth at T-points
107      DO jk = 1, jpkm1
108         hu(:,:) = hu(:,:) + fse3u_n(:,:,jk) * umask(:,:,jk)
109         hv(:,:) = hv(:,:) + fse3v_n(:,:,jk) * vmask(:,:,jk)
110         ht(:,:) = ht(:,:) + fse3t_n(:,:,jk) * tmask(:,:,jk)
111      END DO
112      !                                        ! Inverse of the local depth
113      hur(:,:) = 1._wp / ( hu(:,:) + 1._wp - umask_i(:,:) ) * umask_i(:,:)
114      hvr(:,:) = 1._wp / ( hv(:,:) + 1._wp - vmask_i(:,:) ) * vmask_i(:,:)
115
116                             CALL dom_stp      ! time step
117      IF( nmsh /= 0      )   CALL dom_wri      ! Create a domain file
118      IF( .NOT.ln_rstart )   CALL dom_ctl      ! Domain control
119      !
120      IF( nn_timing == 1 )   CALL timing_stop('dom_init')
121      !
122   END SUBROUTINE dom_init
123
124
125   SUBROUTINE dom_nam
126      !!----------------------------------------------------------------------
127      !!                     ***  ROUTINE dom_nam  ***
128      !!                   
129      !! ** Purpose :   read domaine namelists and print the variables.
130      !!
131      !! ** input   : - namrun namelist
132      !!              - namdom namelist
133      !!              - namnc4 namelist   ! "key_netcdf4" only
134      !!----------------------------------------------------------------------
135      USE ioipsl
136      NAMELIST/namrun/ nn_no   , cn_exp    , cn_ocerst_in, cn_ocerst_out, ln_rstart , nn_rstctl,   &
137         &             nn_it000, nn_itend  , nn_date0    , nn_leapy     , nn_istate , nn_stock ,   &
138         &             nn_write, ln_mskland  , ln_clobber   , nn_chunksz, nn_euler
139      NAMELIST/namdom/ nn_bathy, rn_bathy , rn_e3zps_min, rn_e3zps_rat, nn_msh, rn_hmin,   &
140         &             rn_atfp, rn_rdt ,  nn_closea , ln_crs,  jphgr_msh,       &
141         &             ppglam0, ppgphi0, ppe1_deg, ppe2_deg, ppe1_m, ppe2_m, &
142         &             ppsur, ppa0, ppa1, ppkth, ppacr, ppdzmin, pphmax, ldbletanh, &
143         &             ppa2, ppkth2, ppacr2
144
145#if defined key_netcdf4
146      NAMELIST/namnc4/ nn_nchunks_i, nn_nchunks_j, nn_nchunks_k, ln_nc4zip
147#endif
148      INTEGER  ::   ios                 ! Local integer output status for namelist read
149      !!----------------------------------------------------------------------
150
151      REWIND( numnam_ref )              ! Namelist namrun in reference namelist : Parameters of the run
152      READ  ( numnam_ref, namrun, IOSTAT = ios, ERR = 901)
153901   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namrun in reference namelist', lwp )
154
155      REWIND( numnam_cfg )              ! Namelist namrun in configuration namelist : Parameters of the run
156      READ  ( numnam_cfg, namrun, IOSTAT = ios, ERR = 902 )
157902   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namrun in configuration namelist', lwp )
158      IF(lwm) WRITE ( numond, namrun )
159      !
160      IF(lwp) THEN                  ! control print
161         WRITE(numout,*)
162         WRITE(numout,*) 'dom_nam  : domain initialization through namelist read'
163         WRITE(numout,*) '~~~~~~~ '
164         WRITE(numout,*) '   Namelist namrun'
165         WRITE(numout,*) '      job number                      nn_no      = ', nn_no
166         WRITE(numout,*) '      experiment name for output      cn_exp     = ', cn_exp
167         WRITE(numout,*) '      file prefix restart input       cn_ocerst_in= ', cn_ocerst_in
168         WRITE(numout,*) '      file prefix restart output      cn_ocerst_out= ', cn_ocerst_out
169         WRITE(numout,*) '      restart logical                 ln_rstart  = ', ln_rstart
170         WRITE(numout,*) '      start with forward time step    nn_euler   = ', nn_euler
171         WRITE(numout,*) '      control of time step            nn_rstctl  = ', nn_rstctl
172         WRITE(numout,*) '      number of the first time step   nn_it000   = ', nn_it000
173         WRITE(numout,*) '      number of the last time step    nn_itend   = ', nn_itend
174         WRITE(numout,*) '      initial calendar date aammjj    nn_date0   = ', nn_date0
175         WRITE(numout,*) '      leap year calendar (0/1)        nn_leapy   = ', nn_leapy
176         WRITE(numout,*) '      initial state output            nn_istate  = ', nn_istate
177         WRITE(numout,*) '      frequency of restart file       nn_stock   = ', nn_stock
178         WRITE(numout,*) '      frequency of output file        nn_write   = ', nn_write
179         WRITE(numout,*) '      mask land points                ln_mskland = ', ln_mskland
180         WRITE(numout,*) '      overwrite an existing file      ln_clobber = ', ln_clobber
181         WRITE(numout,*) '      NetCDF chunksize (bytes)        nn_chunksz = ', nn_chunksz
182      ENDIF
183
184      no = nn_no                    ! conversion DOCTOR names into model names (this should disappear soon)
185      cexper = cn_exp
186      nrstdt = nn_rstctl
187      nit000 = nn_it000
188      nitend = nn_itend
189      ndate0 = nn_date0
190      nleapy = nn_leapy
191      ninist = nn_istate
192      nstock = nn_stock
193      nwrite = nn_write
194      neuler = nn_euler
195      IF ( neuler == 1 .AND. .NOT.ln_rstart ) THEN
196         WRITE(ctmp1,*) 'ln_rstart =.FALSE., nn_euler is forced to 0 '
197         CALL ctl_warn( ctmp1 )
198         neuler = 0
199      ENDIF
200
201      !                             ! control of output frequency
202      IF ( nstock == 0 .OR. nstock > nitend ) THEN
203         WRITE(ctmp1,*) 'nstock = ', nstock, ' it is forced to ', nitend
204         CALL ctl_warn( ctmp1 )
205         nstock = nitend
206      ENDIF
207      IF ( nwrite == 0 ) THEN
208         WRITE(ctmp1,*) 'nwrite = ', nwrite, ' it is forced to ', nitend
209         CALL ctl_warn( ctmp1 )
210         nwrite = nitend
211      ENDIF
212
213#if defined key_agrif
214      IF( Agrif_Root() ) THEN
215#endif
216      SELECT CASE ( nleapy )        ! Choose calendar for IOIPSL
217      CASE (  1 ) 
218         CALL ioconf_calendar('gregorian')
219         IF(lwp) WRITE(numout,*) '   The IOIPSL calendar is "gregorian", i.e. leap year'
220      CASE (  0 )
221         CALL ioconf_calendar('noleap')
222         IF(lwp) WRITE(numout,*) '   The IOIPSL calendar is "noleap", i.e. no leap year'
223      CASE ( 30 )
224         CALL ioconf_calendar('360d')
225         IF(lwp) WRITE(numout,*) '   The IOIPSL calendar is "360d", i.e. 360 days in a year'
226      END SELECT
227#if defined key_agrif
228      ENDIF
229#endif
230
231      REWIND( numnam_ref )              ! Namelist namdom in reference namelist : space & time domain (bathymetry, mesh, timestep)
232      READ  ( numnam_ref, namdom, IOSTAT = ios, ERR = 903)
233903   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namdom in reference namelist', lwp )
234 
235      !
236      REWIND( numnam_cfg )              ! Namelist namdom in configuration namelist : space & time domain (bathymetry, mesh, timestep)
237      READ  ( numnam_cfg, namdom, IOSTAT = ios, ERR = 904 )
238904   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namdom in configuration namelist', lwp )
239      IF(lwm) WRITE ( numond, namdom )
240
241      IF(lwp) THEN
242         WRITE(numout,*)
243         WRITE(numout,*) '   Namelist namdom : space & time domain'
244         WRITE(numout,*) '      flag read/compute bathymetry      nn_bathy     = ', nn_bathy
245         WRITE(numout,*) '      Depth (if =0 bathy=jpkm1)         rn_bathy     = ', rn_bathy
246         WRITE(numout,*) '      min depth of the ocean    (>0) or    rn_hmin   = ', rn_hmin
247         WRITE(numout,*) '      min number of ocean level (<0)       '
248         WRITE(numout,*) '      minimum thickness of partial      rn_e3zps_min = ', rn_e3zps_min, ' (m)'
249         WRITE(numout,*) '         step level                     rn_e3zps_rat = ', rn_e3zps_rat
250         WRITE(numout,*) '      create mesh/mask file(s)          nn_msh       = ', nn_msh
251         WRITE(numout,*) '           = 0   no file created           '
252         WRITE(numout,*) '           = 1   mesh_mask                 '
253         WRITE(numout,*) '           = 2   mesh and mask             '
254         WRITE(numout,*) '           = 3   mesh_hgr, msh_zgr and mask'
255         WRITE(numout,*) '      ocean time step                       rn_rdt    = ', rn_rdt
256         WRITE(numout,*) '      asselin time filter parameter         rn_atfp   = ', rn_atfp
257         WRITE(numout,*) '      suppression of closed seas (=0)       nn_closea = ', nn_closea
258         WRITE(numout,*) '      online coarsening of dynamical fields ln_crs    = ', ln_crs
259         WRITE(numout,*) '      type of horizontal mesh jphgr_msh           = ', jphgr_msh
260         WRITE(numout,*) '      longitude of first raw and column T-point ppglam0 = ', ppglam0
261         WRITE(numout,*) '      latitude  of first raw and column T-point ppgphi0 = ', ppgphi0
262         WRITE(numout,*) '      zonal      grid-spacing (degrees) ppe1_deg        = ', ppe1_deg
263         WRITE(numout,*) '      meridional grid-spacing (degrees) ppe2_deg        = ', ppe2_deg
264         WRITE(numout,*) '      zonal      grid-spacing (degrees) ppe1_m          = ', ppe1_m
265         WRITE(numout,*) '      meridional grid-spacing (degrees) ppe2_m          = ', ppe2_m
266         WRITE(numout,*) '      ORCA r4, r2 and r05 coefficients  ppsur           = ', ppsur
267         WRITE(numout,*) '                                        ppa0            = ', ppa0
268         WRITE(numout,*) '                                        ppa1            = ', ppa1
269         WRITE(numout,*) '                                        ppkth           = ', ppkth
270         WRITE(numout,*) '                                        ppacr           = ', ppacr
271         WRITE(numout,*) '      Minimum vertical spacing ppdzmin                  = ', ppdzmin
272         WRITE(numout,*) '      Maximum depth pphmax                              = ', pphmax
273         WRITE(numout,*) '      Use double tanf function for vertical coordinates ldbletanh = ', ldbletanh
274         WRITE(numout,*) '      Double tanh function parameters ppa2              = ', ppa2
275         WRITE(numout,*) '                                      ppkth2            = ', ppkth2
276         WRITE(numout,*) '                                      ppacr2            = ', ppacr2
277      ENDIF
278
279      ntopo     = nn_bathy          ! conversion DOCTOR names into model names (this should disappear soon)
280      e3zps_min = rn_e3zps_min
281      e3zps_rat = rn_e3zps_rat
282      nmsh      = nn_msh
283      atfp      = rn_atfp
284      rdt       = rn_rdt
285
286#if defined key_netcdf4
287      !                             ! NetCDF 4 case   ("key_netcdf4" defined)
288      REWIND( numnam_ref )              ! Namelist namnc4 in reference namelist : NETCDF
289      READ  ( numnam_ref, namnc4, IOSTAT = ios, ERR = 907)
290907   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namnc4 in reference namelist', lwp )
291
292      REWIND( numnam_cfg )              ! Namelist namnc4 in configuration namelist : NETCDF
293      READ  ( numnam_cfg, namnc4, IOSTAT = ios, ERR = 908 )
294908   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namnc4 in configuration namelist', lwp )
295      IF(lwm) WRITE( numond, namnc4 )
296
297      IF(lwp) THEN                        ! control print
298         WRITE(numout,*)
299         WRITE(numout,*) '   Namelist namnc4 - Netcdf4 chunking parameters'
300         WRITE(numout,*) '      number of chunks in i-dimension      nn_nchunks_i   = ', nn_nchunks_i
301         WRITE(numout,*) '      number of chunks in j-dimension      nn_nchunks_j   = ', nn_nchunks_j
302         WRITE(numout,*) '      number of chunks in k-dimension      nn_nchunks_k   = ', nn_nchunks_k
303         WRITE(numout,*) '      apply netcdf4/hdf5 chunking & compression ln_nc4zip = ', ln_nc4zip
304      ENDIF
305
306      ! Put the netcdf4 settings into a simple structure (snc4set, defined in in_out_manager module)
307      ! Note the chunk size in the unlimited (time) dimension will be fixed at 1
308      snc4set%ni   = nn_nchunks_i
309      snc4set%nj   = nn_nchunks_j
310      snc4set%nk   = nn_nchunks_k
311      snc4set%luse = ln_nc4zip
312#else
313      snc4set%luse = .FALSE.        ! No NetCDF 4 case
314#endif
315      !
316   END SUBROUTINE dom_nam
317
318
319   SUBROUTINE dom_ctl
320      !!----------------------------------------------------------------------
321      !!                     ***  ROUTINE dom_ctl  ***
322      !!
323      !! ** Purpose :   Domain control.
324      !!
325      !! ** Method  :   compute and print extrema of masked scale factors
326      !!----------------------------------------------------------------------
327      INTEGER ::   iimi1, ijmi1, iimi2, ijmi2, iima1, ijma1, iima2, ijma2
328      INTEGER, DIMENSION(2) ::   iloc   !
329      REAL(wp) ::   ze1min, ze1max, ze2min, ze2max
330      !!----------------------------------------------------------------------
331      !
332      IF(lk_mpp) THEN
333         CALL mpp_minloc( e1t(:,:), tmask_i(:,:), ze1min, iimi1,ijmi1 )
334         CALL mpp_minloc( e2t(:,:), tmask_i(:,:), ze2min, iimi2,ijmi2 )
335         CALL mpp_maxloc( e1t(:,:), tmask_i(:,:), ze1max, iima1,ijma1 )
336         CALL mpp_maxloc( e2t(:,:), tmask_i(:,:), ze2max, iima2,ijma2 )
337      ELSE
338         ze1min = MINVAL( e1t(:,:), mask = tmask_i(:,:) == 1._wp )   
339         ze2min = MINVAL( e2t(:,:), mask = tmask_i(:,:) == 1._wp )   
340         ze1max = MAXVAL( e1t(:,:), mask = tmask_i(:,:) == 1._wp )   
341         ze2max = MAXVAL( e2t(:,:), mask = tmask_i(:,:) == 1._wp )   
342
343         iloc  = MINLOC( e1t(:,:), mask = tmask_i(:,:) == 1._wp )
344         iimi1 = iloc(1) + nimpp - 1
345         ijmi1 = iloc(2) + njmpp - 1
346         iloc  = MINLOC( e2t(:,:), mask = tmask_i(:,:) == 1._wp )
347         iimi2 = iloc(1) + nimpp - 1
348         ijmi2 = iloc(2) + njmpp - 1
349         iloc  = MAXLOC( e1t(:,:), mask = tmask_i(:,:) == 1._wp )
350         iima1 = iloc(1) + nimpp - 1
351         ijma1 = iloc(2) + njmpp - 1
352         iloc  = MAXLOC( e2t(:,:), mask = tmask_i(:,:) == 1._wp )
353         iima2 = iloc(1) + nimpp - 1
354         ijma2 = iloc(2) + njmpp - 1
355      ENDIF
356      IF(lwp) THEN
357         WRITE(numout,*)
358         WRITE(numout,*) 'dom_ctl : extrema of the masked scale factors'
359         WRITE(numout,*) '~~~~~~~'
360         WRITE(numout,"(14x,'e1t maxi: ',1f10.2,' at i = ',i5,' j= ',i5)") ze1max, iima1, ijma1
361         WRITE(numout,"(14x,'e1t mini: ',1f10.2,' at i = ',i5,' j= ',i5)") ze1min, iimi1, ijmi1
362         WRITE(numout,"(14x,'e2t maxi: ',1f10.2,' at i = ',i5,' j= ',i5)") ze2max, iima2, ijma2
363         WRITE(numout,"(14x,'e2t mini: ',1f10.2,' at i = ',i5,' j= ',i5)") ze2min, iimi2, ijmi2
364      ENDIF
365      !
366   END SUBROUTINE dom_ctl
367
368   SUBROUTINE dom_stiff
369      !!----------------------------------------------------------------------
370      !!                  ***  ROUTINE dom_stiff  ***
371      !!                     
372      !! ** Purpose :   Diagnose maximum grid stiffness/hydrostatic consistency
373      !!
374      !! ** Method  :   Compute Haney (1991) hydrostatic condition ratio
375      !!                Save the maximum in the vertical direction
376      !!                (this number is only relevant in s-coordinates)
377      !!
378      !!                Haney, R. L., 1991: On the pressure gradient force
379      !!                over steep topography in sigma coordinate ocean models.
380      !!                J. Phys. Oceanogr., 21, 610???619.
381      !!----------------------------------------------------------------------
382      INTEGER  ::   ji, jj, jk 
383      REAL(wp) ::   zrxmax
384      REAL(wp), DIMENSION(4) :: zr1
385      !!----------------------------------------------------------------------
386      rx1(:,:) = 0.e0
387      zrxmax   = 0.e0
388      zr1(:)   = 0.e0
389     
390      DO ji = 2, jpim1
391         DO jj = 2, jpjm1
392            DO jk = 1, jpkm1
393               zr1(1) = umask(ji-1,jj  ,jk) *abs( (gdepw_0(ji  ,jj  ,jk  )-gdepw_0(ji-1,jj  ,jk  )  & 
394                    &                         +gdepw_0(ji  ,jj  ,jk+1)-gdepw_0(ji-1,jj  ,jk+1)) &
395                    &                        /(gdepw_0(ji  ,jj  ,jk  )+gdepw_0(ji-1,jj  ,jk  )  &
396                    &                         -gdepw_0(ji  ,jj  ,jk+1)-gdepw_0(ji-1,jj  ,jk+1) + rsmall) )
397               zr1(2) = umask(ji  ,jj  ,jk) *abs( (gdepw_0(ji+1,jj  ,jk  )-gdepw_0(ji  ,jj  ,jk  )  &
398                    &                         +gdepw_0(ji+1,jj  ,jk+1)-gdepw_0(ji  ,jj  ,jk+1)) &
399                    &                        /(gdepw_0(ji+1,jj  ,jk  )+gdepw_0(ji  ,jj  ,jk  )  &
400                    &                         -gdepw_0(ji+1,jj  ,jk+1)-gdepw_0(ji  ,jj  ,jk+1) + rsmall) )
401               zr1(3) = vmask(ji  ,jj  ,jk) *abs( (gdepw_0(ji  ,jj+1,jk  )-gdepw_0(ji  ,jj  ,jk  )  &
402                    &                         +gdepw_0(ji  ,jj+1,jk+1)-gdepw_0(ji  ,jj  ,jk+1)) &
403                    &                        /(gdepw_0(ji  ,jj+1,jk  )+gdepw_0(ji  ,jj  ,jk  )  &
404                    &                         -gdepw_0(ji  ,jj+1,jk+1)-gdepw_0(ji  ,jj  ,jk+1) + rsmall) )
405               zr1(4) = vmask(ji  ,jj-1,jk) *abs( (gdepw_0(ji  ,jj  ,jk  )-gdepw_0(ji  ,jj-1,jk  )  &
406                    &                         +gdepw_0(ji  ,jj  ,jk+1)-gdepw_0(ji  ,jj-1,jk+1)) &
407                    &                        /(gdepw_0(ji  ,jj  ,jk  )+gdepw_0(ji  ,jj-1,jk  )  &
408                    &                         -gdepw_0(ji,  jj  ,jk+1)-gdepw_0(ji  ,jj-1,jk+1) + rsmall) )
409               zrxmax = MAXVAL(zr1(1:4))
410               rx1(ji,jj) = MAX(rx1(ji,jj), zrxmax)
411            END DO
412         END DO
413      END DO
414
415      CALL lbc_lnk( rx1, 'T', 1. )
416
417      zrxmax = MAXVAL(rx1)
418
419      IF( lk_mpp )   CALL mpp_max( zrxmax ) ! max over the global domain
420
421      IF(lwp) THEN
422         WRITE(numout,*)
423         WRITE(numout,*) 'dom_stiff : maximum grid stiffness ratio: ', zrxmax
424         WRITE(numout,*) '~~~~~~~~~'
425      ENDIF
426
427   END SUBROUTINE dom_stiff
428
429
430
431   !!======================================================================
432END MODULE domain
Note: See TracBrowser for help on using the repository browser.