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/UKMO/test_moci_test_suite/NEMOGCM/NEMO/OPA_SRC/DOM – NEMO

source: branches/UKMO/test_moci_test_suite/NEMOGCM/NEMO/OPA_SRC/DOM/domain.F90 @ 8269

Last change on this file since 8269 was 8269, checked in by andmirek, 7 years ago

fix restart read with XIOS whne ln_restart is false

File size: 25.3 KB
RevLine 
[3]1MODULE domain
2   !!==============================================================================
3   !!                       ***  MODULE domain   ***
4   !! Ocean initialization : domain initialization
5   !!==============================================================================
[1438]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
[2528]13   !!            3.3  !  2010-11  (G. Madec)  initialisation in C1D configuration
[4152]14   !!            3.6  !  2013     ( J. Simeon, C. Calone, G. Madec, C. Ethe ) Online coarsening of outputs
[3]15   !!----------------------------------------------------------------------
[1438]16   
17   !!----------------------------------------------------------------------
[3]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   !!----------------------------------------------------------------------
[2528]22   USE oce             ! ocean variables
23   USE dom_oce         ! domain: ocean
[888]24   USE sbc_oce         ! surface boundary condition: ocean
[719]25   USE phycst          ! physical constants
[1601]26   USE closea          ! closed seas
[719]27   USE in_out_manager  ! I/O manager
[3]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
[592]35   USE domvvl          ! variable volume
[2528]36   USE c1d             ! 1D vertical configuration
37   USE dyncor_c1d      ! Coriolis term (c1d case)         (cor_c1d routine)
[3294]38   USE timing          ! Timing
[3680]39   USE lbclnk          ! ocean lateral boundary condition (or mpp link)
[8243]40   USE iom_def, ONLY:lxios_read, lwxios, wxioso
[3]41
42   IMPLICIT NONE
43   PRIVATE
44
[1438]45   PUBLIC   dom_init   ! called by opa.F90
[3]46
47   !! * Substitutions
48#  include "domzgr_substitute.h90"
[1438]49   !!-------------------------------------------------------------------------
[2528]50   !! NEMO/OPA 3.3 , NEMO Consortium (2010)
[888]51   !! $Id$
[2528]52   !! Software governed by the CeCILL licence        (NEMOGCM/NEMO_CeCILL.txt)
[1438]53   !!-------------------------------------------------------------------------
[3]54CONTAINS
55
56   SUBROUTINE dom_init
57      !!----------------------------------------------------------------------
58      !!                  ***  ROUTINE dom_init  ***
59      !!                   
60      !! ** Purpose :   Domain initialization. Call the routines that are
[1601]61      !!              required to create the arrays which define the space
62      !!              and time domain of the ocean model.
[3]63      !!
[1601]64      !! ** Method  : - dom_msk: compute the masks from the bathymetry file
65      !!              - dom_hgr: compute or read the horizontal grid-point position
66      !!                         and scale factors, and the coriolis factor
67      !!              - dom_zgr: define the vertical coordinate and the bathymetry
68      !!              - dom_stp: defined the model time step
69      !!              - dom_wri: create the meshmask file if nmsh=1
[2528]70      !!              - 1D configuration, move Coriolis, u and v at T-point
[3]71      !!----------------------------------------------------------------------
[3764]72      INTEGER ::   jk          ! dummy loop argument
73      INTEGER ::   iconf = 0   ! local integers
[3]74      !!----------------------------------------------------------------------
[1601]75      !
[3764]76      IF( nn_timing == 1 )   CALL timing_start('dom_init')
[3294]77      !
[3]78      IF(lwp) THEN
79         WRITE(numout,*)
80         WRITE(numout,*) 'dom_init : domain initialization'
81         WRITE(numout,*) '~~~~~~~~'
82      ENDIF
[1601]83      !
84                             CALL dom_nam      ! read namelist ( namrun, namdom, namcla )
85                             CALL dom_clo      ! Closed seas and lake
86                             CALL dom_hgr      ! Horizontal mesh
87                             CALL dom_zgr      ! Vertical mesh and bathymetry
88                             CALL dom_msk      ! Masks
[3680]89      IF( ln_sco )           CALL dom_stiff    ! Maximum stiffness ratio/hydrostatic consistency
[4490]90      !
91      ht_0(:,:) = 0.0_wp                       ! Reference ocean depth at T-points
92      hu_0(:,:) = 0.0_wp                       ! Reference ocean depth at U-points
93      hv_0(:,:) = 0.0_wp                       ! Reference ocean depth at V-points
94      DO jk = 1, jpk
95         ht_0(:,:) = ht_0(:,:) + e3t_0(:,:,jk) * tmask(:,:,jk)
96         hu_0(:,:) = hu_0(:,:) + e3u_0(:,:,jk) * umask(:,:,jk)
97         hv_0(:,:) = hv_0(:,:) + e3v_0(:,:,jk) * vmask(:,:,jk)
98      END DO
99      !
[4292]100      IF( lk_vvl )           CALL dom_vvl_init ! Vertical variable mesh
[1601]101      !
[3764]102      IF( lk_c1d         )   CALL cor_c1d      ! 1D configuration: Coriolis set at T-point
[2528]103      !
[4370]104      !
105      hu(:,:) = 0._wp                          ! Ocean depth at U-points
106      hv(:,:) = 0._wp                          ! Ocean depth at V-points
107      ht(:,:) = 0._wp                          ! Ocean depth at T-points
108      DO jk = 1, jpkm1
[4292]109         hu(:,:) = hu(:,:) + fse3u_n(:,:,jk) * umask(:,:,jk)
110         hv(:,:) = hv(:,:) + fse3v_n(:,:,jk) * vmask(:,:,jk)
[4370]111         ht(:,:) = ht(:,:) + fse3t_n(:,:,jk) * tmask(:,:,jk)
[3]112      END DO
[1601]113      !                                        ! Inverse of the local depth
[4990]114      hur(:,:) = 1._wp / ( hu(:,:) + 1._wp - umask_i(:,:) ) * umask_i(:,:)
115      hvr(:,:) = 1._wp / ( hv(:,:) + 1._wp - vmask_i(:,:) ) * vmask_i(:,:)
[216]116
[1601]117                             CALL dom_stp      ! time step
118      IF( nmsh /= 0      )   CALL dom_wri      ! Create a domain file
119      IF( .NOT.ln_rstart )   CALL dom_ctl      ! Domain control
[1438]120      !
[3764]121      IF( nn_timing == 1 )   CALL timing_stop('dom_init')
[3294]122      !
[3]123   END SUBROUTINE dom_init
124
125
126   SUBROUTINE dom_nam
127      !!----------------------------------------------------------------------
128      !!                     ***  ROUTINE dom_nam  ***
129      !!                   
130      !! ** Purpose :   read domaine namelists and print the variables.
131      !!
132      !! ** input   : - namrun namelist
133      !!              - namdom namelist
134      !!              - namcla namelist
[2528]135      !!              - namnc4 namelist   ! "key_netcdf4" only
[3]136      !!----------------------------------------------------------------------
137      USE ioipsl
[5341]138      NAMELIST/namrun/ cn_ocerst_indir, cn_ocerst_outdir, nn_stocklist, ln_rst_list,               &
[6491]139         &             nn_no   , cn_exp    , cn_ocerst_in, cn_ocerst_out, ln_rstart , ln_rstdate, nn_rstctl,   &
[1601]140         &             nn_it000, nn_itend  , nn_date0    , nn_leapy     , nn_istate , nn_stock ,   &
[7924]141         &             nn_write, ln_dimgnnn, ln_mskland  , ln_cfmeta    , ln_clobber, nn_chunksz, nn_euler, &
[8243]142         &             ln_xios_read, nn_wxios
[4245]143      NAMELIST/namdom/ nn_bathy, rn_bathy , rn_e3zps_min, rn_e3zps_rat, nn_msh, rn_hmin,   &
144         &             nn_acc   , rn_atfp     , rn_rdt      , rn_rdtmin ,                  &
[4370]145         &             rn_rdtmax, rn_rdth     , nn_closea , ln_crs,    &
[4147]146         &             jphgr_msh, &
147         &             ppglam0, ppgphi0, ppe1_deg, ppe2_deg, ppe1_m, ppe2_m, &
148         &             ppsur, ppa0, ppa1, ppkth, ppacr, ppdzmin, pphmax, ldbletanh, &
149         &             ppa2, ppkth2, ppacr2
[1601]150      NAMELIST/namcla/ nn_cla
[2528]151#if defined key_netcdf4
152      NAMELIST/namnc4/ nn_nchunks_i, nn_nchunks_j, nn_nchunks_k, ln_nc4zip
153#endif
[4147]154      INTEGER  ::   ios                 ! Local integer output status for namelist read
[3]155      !!----------------------------------------------------------------------
[8001]156      ln_xios_read = .false.            ! set in case ln_xios_read is not in namelist
[8243]157      nn_wxios = 0
[4147]158      REWIND( numnam_ref )              ! Namelist namrun in reference namelist : Parameters of the run
159      READ  ( numnam_ref, namrun, IOSTAT = ios, ERR = 901)
160901   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namrun in reference namelist', lwp )
161
162      REWIND( numnam_cfg )              ! Namelist namrun in configuration namelist : Parameters of the run
163      READ  ( numnam_cfg, namrun, IOSTAT = ios, ERR = 902 )
164902   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namrun in configuration namelist', lwp )
[4624]165      IF(lwm) WRITE ( numond, namrun )
[1601]166      !
167      IF(lwp) THEN                  ! control print
[3]168         WRITE(numout,*)
169         WRITE(numout,*) 'dom_nam  : domain initialization through namelist read'
170         WRITE(numout,*) '~~~~~~~ '
[1601]171         WRITE(numout,*) '   Namelist namrun'
172         WRITE(numout,*) '      job number                      nn_no      = ', nn_no
173         WRITE(numout,*) '      experiment name for output      cn_exp     = ', cn_exp
[4147]174         WRITE(numout,*) '      file prefix restart input       cn_ocerst_in= ', cn_ocerst_in
[5341]175         WRITE(numout,*) '      restart input directory         cn_ocerst_indir= ', cn_ocerst_indir
[4147]176         WRITE(numout,*) '      file prefix restart output      cn_ocerst_out= ', cn_ocerst_out
[5341]177         WRITE(numout,*) '      restart output directory        cn_ocerst_outdir= ', cn_ocerst_outdir
[6491]178         WRITE(numout,*) '      restart logical                 ln_rstart  = ' , ln_rstart
179         WRITE(numout,*) '      datestamping of restarts        ln_rstdate  = ', ln_rstdate
[4370]180         WRITE(numout,*) '      start with forward time step    nn_euler   = ', nn_euler
[1604]181         WRITE(numout,*) '      control of time step            nn_rstctl  = ', nn_rstctl
[1601]182         WRITE(numout,*) '      number of the first time step   nn_it000   = ', nn_it000
183         WRITE(numout,*) '      number of the last time step    nn_itend   = ', nn_itend
184         WRITE(numout,*) '      initial calendar date aammjj    nn_date0   = ', nn_date0
185         WRITE(numout,*) '      leap year calendar (0/1)        nn_leapy   = ', nn_leapy
186         WRITE(numout,*) '      initial state output            nn_istate  = ', nn_istate
[5341]187         IF( ln_rst_list ) THEN
188            WRITE(numout,*) '      list of restart dump times      nn_stocklist   =', nn_stocklist
189         ELSE
190            WRITE(numout,*) '      frequency of restart file       nn_stock   = ', nn_stock
191         ENDIF
[1601]192         WRITE(numout,*) '      frequency of output file        nn_write   = ', nn_write
193         WRITE(numout,*) '      multi file dimgout              ln_dimgnnn = ', ln_dimgnnn
194         WRITE(numout,*) '      mask land points                ln_mskland = ', ln_mskland
[5363]195         WRITE(numout,*) '      additional CF standard metadata ln_cfmeta  = ', ln_cfmeta
[1601]196         WRITE(numout,*) '      overwrite an existing file      ln_clobber = ', ln_clobber
197         WRITE(numout,*) '      NetCDF chunksize (bytes)        nn_chunksz = ', nn_chunksz
[8001]198         WRITE(numout,*) '      READ restart for a single file using XIOS ln_xios_read =', ln_xios_read
[8243]199         WRITE(numout,*) '      Write restart using XIOS        nn_wxios   = ', nn_wxios
[3]200      ENDIF
201
[1601]202      no = nn_no                    ! conversion DOCTOR names into model names (this should disappear soon)
203      cexper = cn_exp
204      nrstdt = nn_rstctl
205      nit000 = nn_it000
206      nitend = nn_itend
207      ndate0 = nn_date0
208      nleapy = nn_leapy
209      ninist = nn_istate
210      nstock = nn_stock
[5341]211      nstocklist = nn_stocklist
[1601]212      nwrite = nn_write
[4370]213      neuler = nn_euler
[8269]214      lxios_read = ln_xios_read.and.ln_rstart
[5341]215      IF ( neuler == 1 .AND. .NOT. ln_rstart ) THEN
[4370]216         WRITE(ctmp1,*) 'ln_rstart =.FALSE., nn_euler is forced to 0 '
217         CALL ctl_warn( ctmp1 )
218         neuler = 0
219      ENDIF
[3]220
[1601]221      !                             ! control of output frequency
[1335]222      IF ( nstock == 0 .OR. nstock > nitend ) THEN
[1601]223         WRITE(ctmp1,*) 'nstock = ', nstock, ' it is forced to ', nitend
[783]224         CALL ctl_warn( ctmp1 )
[1335]225         nstock = nitend
[3]226      ENDIF
227      IF ( nwrite == 0 ) THEN
[1601]228         WRITE(ctmp1,*) 'nwrite = ', nwrite, ' it is forced to ', nitend
[783]229         CALL ctl_warn( ctmp1 )
230         nwrite = nitend
[3]231      ENDIF
232
[2528]233#if defined key_agrif
[1601]234      IF( Agrif_Root() ) THEN
[2528]235#endif
236      SELECT CASE ( nleapy )        ! Choose calendar for IOIPSL
237      CASE (  1 ) 
238         CALL ioconf_calendar('gregorian')
239         IF(lwp) WRITE(numout,*) '   The IOIPSL calendar is "gregorian", i.e. leap year'
240      CASE (  0 )
241         CALL ioconf_calendar('noleap')
242         IF(lwp) WRITE(numout,*) '   The IOIPSL calendar is "noleap", i.e. no leap year'
243      CASE ( 30 )
244         CALL ioconf_calendar('360d')
245         IF(lwp) WRITE(numout,*) '   The IOIPSL calendar is "360d", i.e. 360 days in a year'
246      END SELECT
247#if defined key_agrif
[1601]248      ENDIF
[2528]249#endif
[3]250
[4147]251      REWIND( numnam_ref )              ! Namelist namdom in reference namelist : space & time domain (bathymetry, mesh, timestep)
252      READ  ( numnam_ref, namdom, IOSTAT = ios, ERR = 903)
253903   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namdom in reference namelist', lwp )
[4152]254 
255      !
[4147]256      REWIND( numnam_cfg )              ! Namelist namdom in configuration namelist : space & time domain (bathymetry, mesh, timestep)
257      READ  ( numnam_cfg, namdom, IOSTAT = ios, ERR = 904 )
258904   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namdom in configuration namelist', lwp )
[4624]259      IF(lwm) WRITE ( numond, namdom )
[4147]260
[3]261      IF(lwp) THEN
[72]262         WRITE(numout,*)
[1601]263         WRITE(numout,*) '   Namelist namdom : space & time domain'
264         WRITE(numout,*) '      flag read/compute bathymetry      nn_bathy     = ', nn_bathy
[4245]265         WRITE(numout,*) '      Depth (if =0 bathy=jpkm1)         rn_bathy     = ', rn_bathy
[2528]266         WRITE(numout,*) '      min depth of the ocean    (>0) or    rn_hmin   = ', rn_hmin
267         WRITE(numout,*) '      min number of ocean level (<0)       '
[1601]268         WRITE(numout,*) '      minimum thickness of partial      rn_e3zps_min = ', rn_e3zps_min, ' (m)'
269         WRITE(numout,*) '         step level                     rn_e3zps_rat = ', rn_e3zps_rat
270         WRITE(numout,*) '      create mesh/mask file(s)          nn_msh       = ', nn_msh
[2528]271         WRITE(numout,*) '           = 0   no file created           '
272         WRITE(numout,*) '           = 1   mesh_mask                 '
273         WRITE(numout,*) '           = 2   mesh and mask             '
274         WRITE(numout,*) '           = 3   mesh_hgr, msh_zgr and mask'
[4152]275         WRITE(numout,*) '      ocean time step                       rn_rdt    = ', rn_rdt
276         WRITE(numout,*) '      asselin time filter parameter         rn_atfp   = ', rn_atfp
277         WRITE(numout,*) '      acceleration of converge              nn_acc    = ', nn_acc
278         WRITE(numout,*) '        nn_acc=1: surface tracer rdt        rn_rdtmin = ', rn_rdtmin
279         WRITE(numout,*) '                  bottom  tracer rdt        rdtmax    = ', rn_rdtmax
280         WRITE(numout,*) '                  depth of transition       rn_rdth   = ', rn_rdth
281         WRITE(numout,*) '      suppression of closed seas (=0)       nn_closea = ', nn_closea
282         WRITE(numout,*) '      online coarsening of dynamical fields ln_crs    = ', ln_crs
[4147]283         WRITE(numout,*) '      type of horizontal mesh jphgr_msh           = ', jphgr_msh
284         WRITE(numout,*) '      longitude of first raw and column T-point ppglam0 = ', ppglam0
285         WRITE(numout,*) '      latitude  of first raw and column T-point ppgphi0 = ', ppgphi0
286         WRITE(numout,*) '      zonal      grid-spacing (degrees) ppe1_deg        = ', ppe1_deg
287         WRITE(numout,*) '      meridional grid-spacing (degrees) ppe2_deg        = ', ppe2_deg
288         WRITE(numout,*) '      zonal      grid-spacing (degrees) ppe1_m          = ', ppe1_m
289         WRITE(numout,*) '      meridional grid-spacing (degrees) ppe2_m          = ', ppe2_m
290         WRITE(numout,*) '      ORCA r4, r2 and r05 coefficients  ppsur           = ', ppsur
291         WRITE(numout,*) '                                        ppa0            = ', ppa0
292         WRITE(numout,*) '                                        ppa1            = ', ppa1
293         WRITE(numout,*) '                                        ppkth           = ', ppkth
294         WRITE(numout,*) '                                        ppacr           = ', ppacr
295         WRITE(numout,*) '      Minimum vertical spacing ppdzmin                  = ', ppdzmin
296         WRITE(numout,*) '      Maximum depth pphmax                              = ', pphmax
297         WRITE(numout,*) '      Use double tanf function for vertical coordinates ldbletanh = ', ldbletanh
298         WRITE(numout,*) '      Double tanh function parameters ppa2              = ', ppa2
299         WRITE(numout,*) '                                      ppkth2            = ', ppkth2
300         WRITE(numout,*) '                                      ppacr2            = ', ppacr2
[223]301      ENDIF
302
[1601]303      ntopo     = nn_bathy          ! conversion DOCTOR names into model names (this should disappear soon)
304      e3zps_min = rn_e3zps_min
305      e3zps_rat = rn_e3zps_rat
306      nmsh      = nn_msh
307      nacc      = nn_acc
308      atfp      = rn_atfp
309      rdt       = rn_rdt
310      rdtmin    = rn_rdtmin
311      rdtmax    = rn_rdtmin
312      rdth      = rn_rdth
[8243]313      if (nn_wxios > 0) lwxios = .TRUE. 
314      wxioso = nn_wxios
[1601]315
[4147]316      REWIND( numnam_ref )              ! Namelist namcla in reference namelist : Cross land advection
317      READ  ( numnam_ref, namcla, IOSTAT = ios, ERR = 905)
318905   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namcla in reference namelist', lwp )
319
320      REWIND( numnam_cfg )              ! Namelist namcla in configuration namelist : Cross land advection
321      READ  ( numnam_cfg, namcla, IOSTAT = ios, ERR = 906 )
322906   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namcla in configuration namelist', lwp )
[4624]323      IF(lwm) WRITE( numond, namcla )
[4147]324
[3]325      IF(lwp) THEN
[72]326         WRITE(numout,*)
[1601]327         WRITE(numout,*) '   Namelist namcla'
328         WRITE(numout,*) '      cross land advection                 nn_cla    = ', nn_cla
[3]329      ENDIF
[4147]330      IF ( nn_cla .EQ. 1 ) THEN
331         IF  ( cp_cfg == "orca" .AND. jp_cfg == 2 ) THEN   ! ORCA R2
332            CONTINUE
333         ELSE
334            CALL ctl_stop( 'STOP', 'Cross land advation iplemented only for ORCA2 configuration: cp_cfg = "orca" and jp_cfg = 2 ' )
335         ENDIF
336      ENDIF
[3]337
[2528]338#if defined key_netcdf4
339      !                             ! NetCDF 4 case   ("key_netcdf4" defined)
[4147]340      REWIND( numnam_ref )              ! Namelist namnc4 in reference namelist : NETCDF
341      READ  ( numnam_ref, namnc4, IOSTAT = ios, ERR = 907)
342907   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namnc4 in reference namelist', lwp )
343
344      REWIND( numnam_cfg )              ! Namelist namnc4 in configuration namelist : NETCDF
345      READ  ( numnam_cfg, namnc4, IOSTAT = ios, ERR = 908 )
346908   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namnc4 in configuration namelist', lwp )
[4624]347      IF(lwm) WRITE( numond, namnc4 )
[4147]348
[2528]349      IF(lwp) THEN                        ! control print
350         WRITE(numout,*)
351         WRITE(numout,*) '   Namelist namnc4 - Netcdf4 chunking parameters'
352         WRITE(numout,*) '      number of chunks in i-dimension      nn_nchunks_i   = ', nn_nchunks_i
353         WRITE(numout,*) '      number of chunks in j-dimension      nn_nchunks_j   = ', nn_nchunks_j
354         WRITE(numout,*) '      number of chunks in k-dimension      nn_nchunks_k   = ', nn_nchunks_k
355         WRITE(numout,*) '      apply netcdf4/hdf5 chunking & compression ln_nc4zip = ', ln_nc4zip
356      ENDIF
[1601]357
[2528]358      ! Put the netcdf4 settings into a simple structure (snc4set, defined in in_out_manager module)
359      ! Note the chunk size in the unlimited (time) dimension will be fixed at 1
360      snc4set%ni   = nn_nchunks_i
361      snc4set%nj   = nn_nchunks_j
362      snc4set%nk   = nn_nchunks_k
363      snc4set%luse = ln_nc4zip
364#else
365      snc4set%luse = .FALSE.        ! No NetCDF 4 case
366#endif
[1438]367      !
[3]368   END SUBROUTINE dom_nam
369
370
371   SUBROUTINE dom_ctl
372      !!----------------------------------------------------------------------
373      !!                     ***  ROUTINE dom_ctl  ***
374      !!
375      !! ** Purpose :   Domain control.
376      !!
377      !! ** Method  :   compute and print extrema of masked scale factors
378      !!----------------------------------------------------------------------
379      INTEGER ::   iimi1, ijmi1, iimi2, ijmi2, iima1, ijma1, iima2, ijma2
[1601]380      INTEGER, DIMENSION(2) ::   iloc   !
[3]381      REAL(wp) ::   ze1min, ze1max, ze2min, ze2max
382      !!----------------------------------------------------------------------
[1601]383      !
384      IF(lk_mpp) THEN
[4990]385         CALL mpp_minloc( e1t(:,:), tmask_i(:,:), ze1min, iimi1,ijmi1 )
386         CALL mpp_minloc( e2t(:,:), tmask_i(:,:), ze2min, iimi2,ijmi2 )
387         CALL mpp_maxloc( e1t(:,:), tmask_i(:,:), ze1max, iima1,ijma1 )
388         CALL mpp_maxloc( e2t(:,:), tmask_i(:,:), ze2max, iima2,ijma2 )
[181]389      ELSE
[4990]390         ze1min = MINVAL( e1t(:,:), mask = tmask_i(:,:) == 1._wp )   
391         ze2min = MINVAL( e2t(:,:), mask = tmask_i(:,:) == 1._wp )   
392         ze1max = MAXVAL( e1t(:,:), mask = tmask_i(:,:) == 1._wp )   
393         ze2max = MAXVAL( e2t(:,:), mask = tmask_i(:,:) == 1._wp )   
[32]394
[4990]395         iloc  = MINLOC( e1t(:,:), mask = tmask_i(:,:) == 1._wp )
[181]396         iimi1 = iloc(1) + nimpp - 1
397         ijmi1 = iloc(2) + njmpp - 1
[4990]398         iloc  = MINLOC( e2t(:,:), mask = tmask_i(:,:) == 1._wp )
[181]399         iimi2 = iloc(1) + nimpp - 1
400         ijmi2 = iloc(2) + njmpp - 1
[4990]401         iloc  = MAXLOC( e1t(:,:), mask = tmask_i(:,:) == 1._wp )
[181]402         iima1 = iloc(1) + nimpp - 1
403         ijma1 = iloc(2) + njmpp - 1
[4990]404         iloc  = MAXLOC( e2t(:,:), mask = tmask_i(:,:) == 1._wp )
[181]405         iima2 = iloc(1) + nimpp - 1
406         ijma2 = iloc(2) + njmpp - 1
[32]407      ENDIF
[3]408      IF(lwp) THEN
[1601]409         WRITE(numout,*)
410         WRITE(numout,*) 'dom_ctl : extrema of the masked scale factors'
411         WRITE(numout,*) '~~~~~~~'
[181]412         WRITE(numout,"(14x,'e1t maxi: ',1f10.2,' at i = ',i5,' j= ',i5)") ze1max, iima1, ijma1
413         WRITE(numout,"(14x,'e1t mini: ',1f10.2,' at i = ',i5,' j= ',i5)") ze1min, iimi1, ijmi1
414         WRITE(numout,"(14x,'e2t maxi: ',1f10.2,' at i = ',i5,' j= ',i5)") ze2max, iima2, ijma2
415         WRITE(numout,"(14x,'e2t mini: ',1f10.2,' at i = ',i5,' j= ',i5)") ze2min, iimi2, ijmi2
[3]416      ENDIF
[1438]417      !
[3]418   END SUBROUTINE dom_ctl
419
[3680]420   SUBROUTINE dom_stiff
421      !!----------------------------------------------------------------------
422      !!                  ***  ROUTINE dom_stiff  ***
423      !!                     
424      !! ** Purpose :   Diagnose maximum grid stiffness/hydrostatic consistency
425      !!
426      !! ** Method  :   Compute Haney (1991) hydrostatic condition ratio
427      !!                Save the maximum in the vertical direction
428      !!                (this number is only relevant in s-coordinates)
429      !!
430      !!                Haney, R. L., 1991: On the pressure gradient force
431      !!                over steep topography in sigma coordinate ocean models.
432      !!                J. Phys. Oceanogr., 21, 610???619.
433      !!----------------------------------------------------------------------
434      INTEGER  ::   ji, jj, jk 
435      REAL(wp) ::   zrxmax
436      REAL(wp), DIMENSION(4) :: zr1
437      !!----------------------------------------------------------------------
438      rx1(:,:) = 0.e0
439      zrxmax   = 0.e0
440      zr1(:)   = 0.e0
441     
442      DO ji = 2, jpim1
443         DO jj = 2, jpjm1
444            DO jk = 1, jpkm1
[4292]445               zr1(1) = umask(ji-1,jj  ,jk) *abs( (gdepw_0(ji  ,jj  ,jk  )-gdepw_0(ji-1,jj  ,jk  )  & 
446                    &                         +gdepw_0(ji  ,jj  ,jk+1)-gdepw_0(ji-1,jj  ,jk+1)) &
447                    &                        /(gdepw_0(ji  ,jj  ,jk  )+gdepw_0(ji-1,jj  ,jk  )  &
448                    &                         -gdepw_0(ji  ,jj  ,jk+1)-gdepw_0(ji-1,jj  ,jk+1) + rsmall) )
449               zr1(2) = umask(ji  ,jj  ,jk) *abs( (gdepw_0(ji+1,jj  ,jk  )-gdepw_0(ji  ,jj  ,jk  )  &
450                    &                         +gdepw_0(ji+1,jj  ,jk+1)-gdepw_0(ji  ,jj  ,jk+1)) &
451                    &                        /(gdepw_0(ji+1,jj  ,jk  )+gdepw_0(ji  ,jj  ,jk  )  &
452                    &                         -gdepw_0(ji+1,jj  ,jk+1)-gdepw_0(ji  ,jj  ,jk+1) + rsmall) )
453               zr1(3) = vmask(ji  ,jj  ,jk) *abs( (gdepw_0(ji  ,jj+1,jk  )-gdepw_0(ji  ,jj  ,jk  )  &
454                    &                         +gdepw_0(ji  ,jj+1,jk+1)-gdepw_0(ji  ,jj  ,jk+1)) &
455                    &                        /(gdepw_0(ji  ,jj+1,jk  )+gdepw_0(ji  ,jj  ,jk  )  &
456                    &                         -gdepw_0(ji  ,jj+1,jk+1)-gdepw_0(ji  ,jj  ,jk+1) + rsmall) )
457               zr1(4) = vmask(ji  ,jj-1,jk) *abs( (gdepw_0(ji  ,jj  ,jk  )-gdepw_0(ji  ,jj-1,jk  )  &
458                    &                         +gdepw_0(ji  ,jj  ,jk+1)-gdepw_0(ji  ,jj-1,jk+1)) &
459                    &                        /(gdepw_0(ji  ,jj  ,jk  )+gdepw_0(ji  ,jj-1,jk  )  &
460                    &                         -gdepw_0(ji,  jj  ,jk+1)-gdepw_0(ji  ,jj-1,jk+1) + rsmall) )
[3680]461               zrxmax = MAXVAL(zr1(1:4))
462               rx1(ji,jj) = MAX(rx1(ji,jj), zrxmax)
463            END DO
464         END DO
465      END DO
466
467      CALL lbc_lnk( rx1, 'T', 1. )
468
469      zrxmax = MAXVAL(rx1)
470
471      IF( lk_mpp )   CALL mpp_max( zrxmax ) ! max over the global domain
472
473      IF(lwp) THEN
474         WRITE(numout,*)
475         WRITE(numout,*) 'dom_stiff : maximum grid stiffness ratio: ', zrxmax
476         WRITE(numout,*) '~~~~~~~~~'
477      ENDIF
478
479   END SUBROUTINE dom_stiff
480
481
482
[3]483   !!======================================================================
484END MODULE domain
Note: See TracBrowser for help on using the repository browser.