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 NEMO/trunk/src/OCE/DOM – NEMO

source: NEMO/trunk/src/OCE/DOM/domain.F90 @ 14090

Last change on this file since 14090 was 14090, checked in by hadcv, 4 years ago

#2365: Move dom_tile into domtile.F90

  • Property svn:keywords set to Id
File size: 39.1 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   !!            3.7  !  2015-11  (G. Madec, A. Coward)  time varying zgr by default
16   !!            4.0  !  2016-10  (G. Madec, S. Flavoni)  domain configuration / user defined interface
17   !!            4.1  !  2020-02  (G. Madec, S. Techene)  introduce ssh to h0 ratio
18   !!----------------------------------------------------------------------
19
20   !!----------------------------------------------------------------------
21   !!   dom_init      : initialize the space and time domain
22   !!   dom_glo       : initialize global domain <--> local domain indices
23   !!   dom_nam       : read and contral domain namelists
24   !!   dom_ctl       : control print for the ocean domain
25   !!   domain_cfg    : read the global domain size in domain configuration file
26   !!   cfg_write     : create the domain configuration file
27   !!----------------------------------------------------------------------
28   USE oce            ! ocean variables
29   USE dom_oce        ! domain: ocean
30#if defined key_qco
31   USE domqco         ! quasi-eulerian
32#else
33   USE domvvl         ! variable volume
34#endif
35   USE sshwzv  , ONLY : ssh_init_rst   ! set initial ssh
36   USE sbc_oce        ! surface boundary condition: ocean
37   USE trc_oce        ! shared ocean & passive tracers variab
38   USE phycst         ! physical constants
39   USE domhgr         ! domain: set the horizontal mesh
40   USE domzgr         ! domain: set the vertical mesh
41   USE domtile
42   USE dommsk         ! domain: set the mask system
43   USE domwri         ! domain: write the meshmask file
44   USE c1d            ! 1D configuration
45   USE dyncor_c1d     ! 1D configuration: Coriolis term    (cor_c1d routine)
46   USE wet_dry , ONLY : ll_wd     ! wet & drying flag
47   USE closea  , ONLY : dom_clo   ! closed seas routine
48   !
49   USE in_out_manager ! I/O manager
50   USE iom            ! I/O library
51   USE lbclnk         ! ocean lateral boundary condition (or mpp link)
52   USE lib_mpp        ! distributed memory computing library
53   USE restart        ! only for lrst_oce
54
55   IMPLICIT NONE
56   PRIVATE
57
58   PUBLIC   dom_init     ! called by nemogcm.F90
59   PUBLIC   domain_cfg   ! called by nemogcm.F90
60
61   !! * Substitutions
62#  include "do_loop_substitute.h90"
63   !!-------------------------------------------------------------------------
64   !! NEMO/OCE 4.0 , NEMO Consortium (2018)
65   !! $Id$
66   !! Software governed by the CeCILL license (see ./LICENSE)
67   !!-------------------------------------------------------------------------
68CONTAINS
69
70   SUBROUTINE dom_init( Kbb, Kmm, Kaa )
71      !!----------------------------------------------------------------------
72      !!                  ***  ROUTINE dom_init  ***
73      !!
74      !! ** Purpose :   Domain initialization. Call the routines that are
75      !!              required to create the arrays which define the space
76      !!              and time domain of the ocean model.
77      !!
78      !! ** Method  : - dom_msk: compute the masks from the bathymetry file
79      !!              - dom_hgr: compute or read the horizontal grid-point position
80      !!                         and scale factors, and the coriolis factor
81      !!              - dom_zgr: define the vertical coordinate and the bathymetry
82      !!              - dom_wri: create the meshmask file (ln_meshmask=T)
83      !!              - 1D configuration, move Coriolis, u and v at T-point
84      !!----------------------------------------------------------------------
85      INTEGER          , INTENT(in) :: Kbb, Kmm, Kaa          ! ocean time level indices
86      !
87      INTEGER ::   ji, jj, jk, jt   ! dummy loop indices
88      INTEGER ::   iconf = 0    ! local integers
89      REAL(wp)::   zrdt
90      CHARACTER (len=64) ::   cform = "(A12, 3(A13, I7))"
91      INTEGER , DIMENSION(jpi,jpj) ::   ik_top , ik_bot       ! top and bottom ocean level
92      REAL(wp), DIMENSION(jpi,jpj) ::   z1_hu_0, z1_hv_0
93      !!----------------------------------------------------------------------
94      !
95      IF(lwp) THEN         ! Ocean domain Parameters (control print)
96         WRITE(numout,*)
97         WRITE(numout,*) 'dom_init : domain initialization'
98         WRITE(numout,*) '~~~~~~~~'
99         !
100         WRITE(numout,*)     '   Domain info'
101         WRITE(numout,*)     '      dimension of model:'
102         WRITE(numout,*)     '             Local domain      Global domain       Data domain '
103         WRITE(numout,cform) '        ','   jpi     : ', jpi, '   jpiglo  : ', jpiglo
104         WRITE(numout,cform) '        ','   jpj     : ', jpj, '   jpjglo  : ', jpjglo
105         WRITE(numout,cform) '        ','   jpk     : ', jpk, '   jpkglo  : ', jpkglo
106         WRITE(numout,cform) '       ' ,'   jpij    : ', jpij
107         WRITE(numout,*)     '      mpp local domain info (mpp):'
108         WRITE(numout,*)     '              jpni    : ', jpni, '   nn_hls  : ', nn_hls
109         WRITE(numout,*)     '              jpnj    : ', jpnj, '   nn_hls  : ', nn_hls
110         WRITE(numout,*)     '              jpnij   : ', jpnij
111         WRITE(numout,*)     '      lateral boundary of the Global domain : jperio  = ', jperio
112         SELECT CASE ( jperio )
113         CASE( 0 )   ;   WRITE(numout,*) '         (i.e. closed)'
114         CASE( 1 )   ;   WRITE(numout,*) '         (i.e. cyclic east-west)'
115         CASE( 2 )   ;   WRITE(numout,*) '         (i.e. cyclic north-south)'
116         CASE( 3 )   ;   WRITE(numout,*) '         (i.e. north fold with T-point pivot)'
117         CASE( 4 )   ;   WRITE(numout,*) '         (i.e. cyclic east-west and north fold with T-point pivot)'
118         CASE( 5 )   ;   WRITE(numout,*) '         (i.e. north fold with F-point pivot)'
119         CASE( 6 )   ;   WRITE(numout,*) '         (i.e. cyclic east-west and north fold with F-point pivot)'
120         CASE( 7 )   ;   WRITE(numout,*) '         (i.e. cyclic east-west and north-south)'
121         CASE DEFAULT
122            CALL ctl_stop( 'dom_init:   jperio is out of range' )
123         END SELECT
124         WRITE(numout,*)     '      Ocean model configuration used:'
125         WRITE(numout,*)     '         cn_cfg = ', TRIM( cn_cfg ), '   nn_cfg = ', nn_cfg
126      ENDIF
127
128      !
129      !           !==  Reference coordinate system  ==!
130      !
131      CALL dom_glo                            ! global domain versus local domain
132      CALL dom_nam                            ! read namelist ( namrun, namdom )
133      CALL dom_tile( ntsi, ntsj, ntei, ntej ) ! Tile domain
134
135      !
136      CALL dom_hgr                      ! Horizontal mesh
137
138      IF( ln_closea ) CALL dom_clo      ! Read in masks to define closed seas and lakes
139
140      CALL dom_zgr( ik_top, ik_bot )    ! Vertical mesh and bathymetry (return top and bottom ocean t-level indices)
141
142      CALL dom_msk( ik_top, ik_bot )    ! Masks
143      !
144      ht_0(:,:) = 0._wp  ! Reference ocean thickness
145      hu_0(:,:) = 0._wp
146      hv_0(:,:) = 0._wp
147      hf_0(:,:) = 0._wp
148      DO jk = 1, jpkm1
149         ht_0(:,:) = ht_0(:,:) + e3t_0(:,:,jk) * tmask(:,:,jk)
150         hu_0(:,:) = hu_0(:,:) + e3u_0(:,:,jk) * umask(:,:,jk)
151         hv_0(:,:) = hv_0(:,:) + e3v_0(:,:,jk) * vmask(:,:,jk)
152      END DO
153      !
154      DO jk = 1, jpkm1
155         hf_0(1:jpim1,:) = hf_0(1:jpim1,:) + e3f_0(1:jpim1,:,jk)*vmask(1:jpim1,:,jk)*vmask(2:jpi,:,jk)
156      END DO
157      CALL lbc_lnk('domain', hf_0, 'F', 1._wp)
158      !
159      IF( lk_SWE ) THEN      ! SWE case redefine hf_0
160         hf_0(:,:) = hf_0(:,:) + e3f_0(:,:,1) * ssfmask(:,:)
161      ENDIF
162      !
163      r1_ht_0(:,:) = ssmask (:,:) / ( ht_0(:,:) + 1._wp -  ssmask (:,:) )
164      r1_hu_0(:,:) = ssumask(:,:) / ( hu_0(:,:) + 1._wp -  ssumask(:,:) )
165      r1_hv_0(:,:) = ssvmask(:,:) / ( hv_0(:,:) + 1._wp -  ssvmask(:,:) )
166      r1_hf_0(:,:) = ssfmask(:,:) / ( hf_0(:,:) + 1._wp -  ssfmask(:,:) )
167      !
168      IF( ll_wd ) THEN       ! wet and drying (check ht_0 >= 0)
169         DO_2D( 1, 1, 1, 1 )
170            IF( ht_0(ji,jj) < 0._wp .AND. ssmask(ji,jj) == 1._wp ) THEN
171               CALL ctl_stop( 'ssh_init_rst : ht_0 must be positive at potentially wet points' )
172            ENDIF
173         END_2D
174      ENDIF
175      !
176      !           !==  initialisation of time varying coordinate  ==!
177      !
178      !                                 != ssh initialization
179      IF( .NOT.l_offline .AND. .NOT.l_SAS ) THEN
180         CALL ssh_init_rst( Kbb, Kmm, Kaa )
181      ELSE
182         ssh(:,:,:) = 0._wp
183      ENDIF
184      !
185#if defined key_qco
186      !                                 != Quasi-Euerian coordinate case
187      !
188      IF( .NOT.l_offline )   CALL dom_qco_init( Kbb, Kmm, Kaa )
189#else
190      !
191      IF( ln_linssh ) THEN              != Fix in time : set to the reference one for all
192         !
193         DO jt = 1, jpt                         ! depth of t- and w-grid-points
194            gdept(:,:,:,jt) = gdept_0(:,:,:)
195            gdepw(:,:,:,jt) = gdepw_0(:,:,:)
196         END DO
197            gde3w(:,:,:)    = gde3w_0(:,:,:)    ! = gdept as the sum of e3t
198         !
199         DO jt = 1, jpt                         ! vertical scale factors
200            e3t (:,:,:,jt) =  e3t_0(:,:,:)
201            e3u (:,:,:,jt) =  e3u_0(:,:,:)
202            e3v (:,:,:,jt) =  e3v_0(:,:,:)
203            e3w (:,:,:,jt) =  e3w_0(:,:,:)
204            e3uw(:,:,:,jt) = e3uw_0(:,:,:)
205            e3vw(:,:,:,jt) = e3vw_0(:,:,:)
206         END DO
207            e3f (:,:,:)    =  e3f_0(:,:,:)
208         !
209         DO jt = 1, jpt                         ! water column thickness and its inverse
210               hu(:,:,jt) =    hu_0(:,:)
211               hv(:,:,jt) =    hv_0(:,:)
212            r1_hu(:,:,jt) = r1_hu_0(:,:)
213            r1_hv(:,:,jt) = r1_hv_0(:,:)
214         END DO
215               ht   (:,:) =    ht_0(:,:)
216         !
217      ELSE                              != Time varying : initialize before/now/after variables
218         !
219         IF( .NOT.l_offline )   CALL dom_vvl_init( Kbb, Kmm, Kaa )
220         !
221      ENDIF
222#endif
223
224      !
225
226      IF( lk_c1d         )   CALL cor_c1d       ! 1D configuration: Coriolis set at T-point
227      !
228
229#if defined key_agrif
230      IF( .NOT. Agrif_Root() ) CALL Agrif_Init_Domain( Kbb, Kmm, Kaa )
231#endif
232      IF( ln_meshmask    )   CALL dom_wri       ! Create a domain file
233      IF( .NOT.ln_rstart )   CALL dom_ctl       ! Domain control
234      !
235      IF( ln_write_cfg   )   CALL cfg_write     ! create the configuration file
236      !
237      IF(lwp) THEN
238         WRITE(numout,*)
239         WRITE(numout,*) 'dom_init :   ==>>>   END of domain initialization'
240         WRITE(numout,*) '~~~~~~~~'
241         WRITE(numout,*)
242      ENDIF
243      !
244   END SUBROUTINE dom_init
245
246
247   SUBROUTINE dom_glo
248      !!----------------------------------------------------------------------
249      !!                     ***  ROUTINE dom_glo  ***
250      !!
251      !! ** Purpose :   initialization of global domain <--> local domain indices
252      !!
253      !! ** Method  :
254      !!
255      !! ** Action  : - mig , mjg : local  domain indices ==> global domain, including halos, indices
256      !!              - mig0, mjg0: local  domain indices ==> global domain, excluding halos, indices
257      !!              - mi0 , mi1 : global domain indices ==> local  domain indices
258      !!              - mj0 , mj1   (if global point not in the local domain ==> mi0>mi1 and/or mj0>mj1)
259      !!----------------------------------------------------------------------
260      INTEGER ::   ji, jj   ! dummy loop argument
261      !!----------------------------------------------------------------------
262      !
263      DO ji = 1, jpi                 ! local domain indices ==> global domain indices, including halos
264        mig(ji) = ji + nimpp - 1
265      END DO
266      DO jj = 1, jpj
267        mjg(jj) = jj + njmpp - 1
268      END DO
269      !                              ! local domain indices ==> global domain indices, excluding halos
270      !
271      mig0(:) = mig(:) - nn_hls
272      mjg0(:) = mjg(:) - nn_hls
273      ! WARNING: to keep compatibility with the trunk that was including periodocity into the input data,
274      ! we must define mig0 and mjg0 as bellow.
275      ! Once we decide to forget trunk compatibility, we must simply define mig0 and mjg0 as:
276      mig0_oldcmp(:) = mig0(:) + COUNT( (/ jperio == 1 .OR. jperio == 4 .OR. jperio == 6 .OR. jperio == 7 /) )
277      mjg0_oldcmp(:) = mjg0(:) + COUNT( (/ jperio == 2 .OR. jperio == 7 /) )
278      !
279      !                              ! global domain, including halos, indices ==> local domain indices
280      !                                   ! (return (m.0,m.1)=(1,0) if data domain gridpoint is to the west/south of the
281      !                                   ! local domain, or (m.0,m.1)=(jp.+1,jp.) to the east/north of local domain.
282      DO ji = 1, jpiglo
283        mi0(ji) = MAX( 1 , MIN( ji - nimpp + 1, jpi+1 ) )
284        mi1(ji) = MAX( 0 , MIN( ji - nimpp + 1, jpi   ) )
285      END DO
286      DO jj = 1, jpjglo
287        mj0(jj) = MAX( 1 , MIN( jj - njmpp + 1, jpj+1 ) )
288        mj1(jj) = MAX( 0 , MIN( jj - njmpp + 1, jpj   ) )
289      END DO
290      IF(lwp) THEN                   ! control print
291         WRITE(numout,*)
292         WRITE(numout,*) 'dom_glo : domain: global <<==>> local '
293         WRITE(numout,*) '~~~~~~~ '
294         WRITE(numout,*) '   global domain:   jpiglo = ', jpiglo, ' jpjglo = ', jpjglo, ' jpkglo = ', jpkglo
295         WRITE(numout,*) '   local  domain:   jpi    = ', jpi   , ' jpj    = ', jpj   , ' jpk    = ', jpk
296         WRITE(numout,*)
297      ENDIF
298      !
299   END SUBROUTINE dom_glo
300
301
302   SUBROUTINE dom_nam
303      !!----------------------------------------------------------------------
304      !!                     ***  ROUTINE dom_nam  ***
305      !!
306      !! ** Purpose :   read domaine namelists and print the variables.
307      !!
308      !! ** input   : - namrun namelist
309      !!              - namdom namelist
310      !!              - namtile namelist
311      !!              - namnc4 namelist   ! "key_netcdf4" only
312      !!----------------------------------------------------------------------
313      USE ioipsl
314      !!
315      INTEGER ::   ios   ! Local integer
316      REAL(wp)::   zrdt
317      !!----------------------------------------------------------------------
318      !
319      NAMELIST/namrun/ cn_ocerst_indir, cn_ocerst_outdir, nn_stocklist, ln_rst_list,                 &
320         &             nn_no   , cn_exp   , cn_ocerst_in, cn_ocerst_out, ln_rstart , nn_rstctl ,     &
321         &             nn_it000, nn_itend , nn_date0    , nn_time0     , nn_leapy  , nn_istate ,     &
322         &             nn_stock, nn_write , ln_mskland  , ln_clobber   , nn_chunksz, ln_1st_euler  , &
323         &             ln_cfmeta, ln_xios_read, nn_wxios
324      NAMELIST/namdom/ ln_linssh, rn_Dt, rn_atfp, ln_crs, ln_meshmask
325      NAMELIST/namtile/ ln_tile, nn_ltile_i, nn_ltile_j
326#if defined key_netcdf4
327      NAMELIST/namnc4/ nn_nchunks_i, nn_nchunks_j, nn_nchunks_k, ln_nc4zip
328#endif
329      !!----------------------------------------------------------------------
330      !
331      IF(lwp) THEN
332         WRITE(numout,*)
333         WRITE(numout,*) 'dom_nam : domain initialization through namelist read'
334         WRITE(numout,*) '~~~~~~~ '
335      ENDIF
336      !
337      !                       !=======================!
338      !                       !==  namelist namdom  ==!
339      !                       !=======================!
340      !
341      READ  ( numnam_ref, namdom, IOSTAT = ios, ERR = 903)
342903   IF( ios /= 0 )   CALL ctl_nam ( ios , 'namdom in reference namelist' )
343      READ  ( numnam_cfg, namdom, IOSTAT = ios, ERR = 904 )
344904   IF( ios >  0 )   CALL ctl_nam ( ios , 'namdom in configuration namelist' )
345      IF(lwm) WRITE( numond, namdom )
346      !
347#if defined key_agrif
348      IF( .NOT. Agrif_Root() ) THEN    ! AGRIF child, subdivide the Parent timestep
349         rn_Dt = Agrif_Parent (rn_Dt ) / Agrif_Rhot()
350      ENDIF
351#endif
352      !
353      IF(lwp) THEN
354         WRITE(numout,*)
355         WRITE(numout,*) '   Namelist : namdom   ---   space & time domain'
356         WRITE(numout,*) '      linear free surface (=T)                ln_linssh   = ', ln_linssh
357         WRITE(numout,*) '      create mesh/mask file                   ln_meshmask = ', ln_meshmask
358         WRITE(numout,*) '      ocean time step                         rn_Dt       = ', rn_Dt
359         WRITE(numout,*) '      asselin time filter parameter           rn_atfp     = ', rn_atfp
360         WRITE(numout,*) '      online coarsening of dynamical fields   ln_crs      = ', ln_crs
361      ENDIF
362      !
363      ! set current model timestep rDt = 2*rn_Dt if MLF or rDt = rn_Dt if RK3
364      rDt   = 2._wp * rn_Dt
365      r1_Dt = 1._wp / rDt
366      !
367      IF( l_SAS .AND. .NOT.ln_linssh ) THEN
368         CALL ctl_warn( 'SAS requires linear ssh : force ln_linssh = T' )
369         ln_linssh = .TRUE.
370      ENDIF
371      !
372#if defined key_qco
373      IF( ln_linssh )   CALL ctl_stop( 'STOP','domain: key_qco and ln_linssh = T are incompatible' )
374#endif
375      !
376      !                       !=======================!
377      !                       !==  namelist namrun  ==!
378      !                       !=======================!
379      !
380      READ  ( numnam_ref, namrun, IOSTAT = ios, ERR = 901)
381901   IF( ios /= 0 )   CALL ctl_nam ( ios , 'namrun in reference namelist' )
382      READ  ( numnam_cfg, namrun, IOSTAT = ios, ERR = 902 )
383902   IF( ios >  0 )   CALL ctl_nam ( ios , 'namrun in configuration namelist' )
384      IF(lwm) WRITE ( numond, namrun )
385
386#if defined key_agrif
387      IF( .NOT. Agrif_Root() ) THEN
388            nn_it000 = (Agrif_Parent(nn_it000)-1)*Agrif_IRhot() + 1
389            nn_itend =  Agrif_Parent(nn_itend)   *Agrif_IRhot()
390      ENDIF
391#endif
392      !
393      IF(lwp) THEN                  ! control print
394         WRITE(numout,*) '   Namelist : namrun   ---   run parameters'
395         WRITE(numout,*) '      Assimilation cycle              nn_no           = ', nn_no
396         WRITE(numout,*) '      experiment name for output      cn_exp          = ', TRIM( cn_exp           )
397         WRITE(numout,*) '      file prefix restart input       cn_ocerst_in    = ', TRIM( cn_ocerst_in     )
398         WRITE(numout,*) '      restart input directory         cn_ocerst_indir = ', TRIM( cn_ocerst_indir  )
399         WRITE(numout,*) '      file prefix restart output      cn_ocerst_out   = ', TRIM( cn_ocerst_out    )
400         WRITE(numout,*) '      restart output directory        cn_ocerst_outdir= ', TRIM( cn_ocerst_outdir )
401         WRITE(numout,*) '      restart logical                 ln_rstart       = ', ln_rstart
402         WRITE(numout,*) '      start with forward time step    ln_1st_euler    = ', ln_1st_euler
403         WRITE(numout,*) '      control of time step            nn_rstctl       = ', nn_rstctl
404         WRITE(numout,*) '      number of the first time step   nn_it000        = ', nn_it000
405         WRITE(numout,*) '      number of the last time step    nn_itend        = ', nn_itend
406         WRITE(numout,*) '      initial calendar date aammjj    nn_date0        = ', nn_date0
407         WRITE(numout,*) '      initial time of day in hhmm     nn_time0        = ', nn_time0
408         WRITE(numout,*) '      leap year calendar (0/1)        nn_leapy        = ', nn_leapy
409         WRITE(numout,*) '      initial state output            nn_istate       = ', nn_istate
410         IF( ln_rst_list ) THEN
411            WRITE(numout,*) '      list of restart dump times      nn_stocklist    =', nn_stocklist
412         ELSE
413            WRITE(numout,*) '      frequency of restart file       nn_stock        = ', nn_stock
414         ENDIF
415#if ! defined key_iomput
416         WRITE(numout,*) '      frequency of output file        nn_write        = ', nn_write
417#endif
418         WRITE(numout,*) '      mask land points                ln_mskland      = ', ln_mskland
419         WRITE(numout,*) '      additional CF standard metadata ln_cfmeta       = ', ln_cfmeta
420         WRITE(numout,*) '      overwrite an existing file      ln_clobber      = ', ln_clobber
421         WRITE(numout,*) '      NetCDF chunksize (bytes)        nn_chunksz      = ', nn_chunksz
422         IF( TRIM(Agrif_CFixed()) == '0' ) THEN
423            WRITE(numout,*) '      READ restart for a single file using XIOS ln_xios_read =', ln_xios_read
424            WRITE(numout,*) '      Write restart using XIOS        nn_wxios   = ', nn_wxios
425         ELSE
426            WRITE(numout,*) "      AGRIF: nn_wxios will be ingored. See setting for parent"
427            WRITE(numout,*) "      AGRIF: ln_xios_read will be ingored. See setting for parent"
428         ENDIF
429      ENDIF
430
431      cexper = cn_exp         ! conversion DOCTOR names into model names (this should disappear soon)
432      nrstdt = nn_rstctl
433      nit000 = nn_it000
434      nitend = nn_itend
435      ndate0 = nn_date0
436      nleapy = nn_leapy
437      ninist = nn_istate
438      !
439      !                                        !==  Set parameters for restart reading using xIOS  ==!
440      !
441      IF( TRIM(Agrif_CFixed()) == '0' ) THEN
442         lrxios = ln_xios_read .AND. ln_rstart
443         IF( nn_wxios > 0 )   lwxios = .TRUE.           !* set output file type for XIOS based on NEMO namelist
444         nxioso = nn_wxios
445      ENDIF
446      !                                        !==  Check consistency between ln_rstart and ln_1st_euler  ==!   (i.e. set l_1st_euler)
447      l_1st_euler = ln_1st_euler
448      !
449      IF( ln_rstart ) THEN                              !*  Restart case
450         !
451         IF(lwp) WRITE(numout,*)
452         IF(lwp) WRITE(numout,*) '   open the restart file'
453         CALL rst_read_open                                              !- Open the restart file
454         !
455         IF( iom_varid( numror, 'rdt', ldstop = .FALSE. ) > 0 ) THEN     !- Check time-step consistency and force Euler restart if changed
456            CALL iom_get( numror, 'rdt', zrdt )
457            IF( zrdt /= rn_Dt ) THEN
458               IF(lwp) WRITE( numout,*)
459               IF(lwp) WRITE( numout,*) '   rn_Dt = ', rn_Dt,' not equal to the READ one rdt = ', zrdt
460               IF(lwp) WRITE( numout,*)
461               IF(lwp) WRITE( numout,*) '      ==>>>   forced euler first time-step'
462               l_1st_euler =  .TRUE.
463            ENDIF
464         ENDIF
465         !
466         IF( .NOT.l_SAS .AND. iom_varid( numror, 'sshb', ldstop = .FALSE. ) <= 0 ) THEN   !- Check absence of one of the Kbb field (here sshb)
467            !                                                                             !  (any Kbb field is missing ==> all Kbb fields are missing)
468            IF( .NOT.l_1st_euler ) THEN
469               CALL ctl_warn('dom_nam : ssh at Kbb not found in restart files ',   &
470                  &                        'l_1st_euler forced to .true. and ' ,   &
471                  &                        'ssh(Kbb) = ssh(Kmm) '                  )
472               l_1st_euler = .TRUE.
473            ENDIF
474         ENDIF
475      ELSEIF( .NOT.l_1st_euler ) THEN                   !*  Initialization case
476         IF(lwp) WRITE(numout,*)
477         IF(lwp) WRITE(numout,*)'   ==>>>   Start from rest (ln_rstart=F)'
478         IF(lwp) WRITE(numout,*)'           an Euler initial time step is used : l_1st_euler is forced to .true. '
479         l_1st_euler = .TRUE.
480      ENDIF
481      !
482      !                                        !==  control of output frequency  ==!
483      !
484      IF( .NOT. ln_rst_list ) THEN   ! we use nn_stock
485         IF( nn_stock == -1 )   CALL ctl_warn( 'nn_stock = -1 --> no restart will be done' )
486         IF( nn_stock == 0 .OR. nn_stock > nitend ) THEN
487            WRITE(ctmp1,*) 'nn_stock = ', nn_stock, ' it is forced to ', nitend
488            CALL ctl_warn( ctmp1 )
489            nn_stock = nitend
490         ENDIF
491      ENDIF
492#if ! defined key_iomput
493      IF( nn_write == -1 )   CALL ctl_warn( 'nn_write = -1 --> no output files will be done' )
494      IF ( nn_write == 0 ) THEN
495         WRITE(ctmp1,*) 'nn_write = ', nn_write, ' it is forced to ', nitend
496         CALL ctl_warn( ctmp1 )
497         nn_write = nitend
498      ENDIF
499#endif
500
501      IF( Agrif_Root() ) THEN
502         IF(lwp) WRITE(numout,*)
503         SELECT CASE ( nleapy )                !==  Choose calendar for IOIPSL  ==!
504         CASE (  1 )
505            CALL ioconf_calendar('gregorian')
506            IF(lwp) WRITE(numout,*) '   ==>>>   The IOIPSL calendar is "gregorian", i.e. leap year'
507         CASE (  0 )
508            CALL ioconf_calendar('noleap')
509            IF(lwp) WRITE(numout,*) '   ==>>>   The IOIPSL calendar is "noleap", i.e. no leap year'
510         CASE ( 30 )
511            CALL ioconf_calendar('360d')
512            IF(lwp) WRITE(numout,*) '   ==>>>   The IOIPSL calendar is "360d", i.e. 360 days in a year'
513         END SELECT
514      ENDIF
515      !
516      !                       !========================!
517      !                       !==  namelist namtile  ==!
518      !                       !========================!
519      !
520      READ  ( numnam_ref, namtile, IOSTAT = ios, ERR = 905 )
521905   IF( ios /= 0 )   CALL ctl_nam ( ios , 'namtile in reference namelist' )
522      READ  ( numnam_cfg, namtile, IOSTAT = ios, ERR = 906 )
523906   IF( ios >  0 )   CALL ctl_nam ( ios , 'namtile in configuration namelist' )
524      IF(lwm) WRITE( numond, namtile )
525
526      IF(lwp) THEN
527         WRITE(numout,*)
528         WRITE(numout,*)    '   Namelist : namtile   ---   Domain tiling decomposition'
529         WRITE(numout,*)    '      Tiling (T) or not (F)                ln_tile    = ', ln_tile
530         WRITE(numout,*)    '      Length of tile in i                  nn_ltile_i = ', nn_ltile_i
531         WRITE(numout,*)    '      Length of tile in j                  nn_ltile_j = ', nn_ltile_j
532         WRITE(numout,*)
533         IF( ln_tile ) THEN
534            WRITE(numout,*) '      The domain will be decomposed into tiles of size', nn_ltile_i, 'x', nn_ltile_j
535         ELSE
536            WRITE(numout,*) '      Domain tiling will NOT be used'
537         ENDIF
538      ENDIF
539      !
540#if defined key_netcdf4
541      !                       !=======================!
542      !                       !==  namelist namnc4  ==!   NetCDF 4 case   ("key_netcdf4" defined)
543      !                       !=======================!
544      !
545      READ  ( numnam_ref, namnc4, IOSTAT = ios, ERR = 907)
546907   IF( ios /= 0 )   CALL ctl_nam ( ios , 'namnc4 in reference namelist' )
547      READ  ( numnam_cfg, namnc4, IOSTAT = ios, ERR = 908 )
548908   IF( ios >  0 )   CALL ctl_nam ( ios , 'namnc4 in configuration namelist' )
549      IF(lwm) WRITE( numond, namnc4 )
550
551      IF(lwp) THEN                        ! control print
552         WRITE(numout,*)
553         WRITE(numout,*) '   Namelist namnc4 - Netcdf4 chunking parameters ("key_netcdf4" defined)'
554         WRITE(numout,*) '      number of chunks in i-dimension             nn_nchunks_i = ', nn_nchunks_i
555         WRITE(numout,*) '      number of chunks in j-dimension             nn_nchunks_j = ', nn_nchunks_j
556         WRITE(numout,*) '      number of chunks in k-dimension             nn_nchunks_k = ', nn_nchunks_k
557         WRITE(numout,*) '      apply netcdf4/hdf5 chunking & compression   ln_nc4zip    = ', ln_nc4zip
558      ENDIF
559
560      ! Put the netcdf4 settings into a simple structure (snc4set, defined in in_out_manager module)
561      ! Note the chunk size in the unlimited (time) dimension will be fixed at 1
562      snc4set%ni   = nn_nchunks_i
563      snc4set%nj   = nn_nchunks_j
564      snc4set%nk   = nn_nchunks_k
565      snc4set%luse = ln_nc4zip
566#else
567      snc4set%luse = .FALSE.        ! No NetCDF 4 case
568#endif
569      !
570   END SUBROUTINE dom_nam
571
572
573   SUBROUTINE dom_ctl
574      !!----------------------------------------------------------------------
575      !!                     ***  ROUTINE dom_ctl  ***
576      !!
577      !! ** Purpose :   Domain control.
578      !!
579      !! ** Method  :   compute and print extrema of masked scale factors
580      !!----------------------------------------------------------------------
581      LOGICAL, DIMENSION(jpi,jpj) ::   llmsk
582      INTEGER, DIMENSION(2)       ::   imil, imip, imi1, imi2, imal, imap, ima1, ima2
583      REAL(wp)                    ::   zglmin, zglmax, zgpmin, zgpmax, ze1min, ze1max, ze2min, ze2max
584      !!----------------------------------------------------------------------
585      !
586      llmsk = tmask_h(:,:) == 1._wp
587      !
588      CALL mpp_minloc( 'domain', glamt(:,:), llmsk, zglmin, imil )
589      CALL mpp_minloc( 'domain', gphit(:,:), llmsk, zgpmin, imip )
590      CALL mpp_minloc( 'domain',   e1t(:,:), llmsk, ze1min, imi1 )
591      CALL mpp_minloc( 'domain',   e2t(:,:), llmsk, ze2min, imi2 )
592      CALL mpp_maxloc( 'domain', glamt(:,:), llmsk, zglmax, imal )
593      CALL mpp_maxloc( 'domain', gphit(:,:), llmsk, zgpmax, imap )
594      CALL mpp_maxloc( 'domain',   e1t(:,:), llmsk, ze1max, ima1 )
595      CALL mpp_maxloc( 'domain',   e2t(:,:), llmsk, ze2max, ima2 )
596      !
597      IF(lwp) THEN
598         WRITE(numout,*)
599         WRITE(numout,*) 'dom_ctl : extrema of the masked scale factors'
600         WRITE(numout,*) '~~~~~~~'
601         WRITE(numout,"(14x,'glamt mini: ',1f10.2,' at i = ',i5,' j= ',i5)") zglmin, imil(1), imil(2)
602         WRITE(numout,"(14x,'glamt maxi: ',1f10.2,' at i = ',i5,' j= ',i5)") zglmax, imal(1), imal(2)
603         WRITE(numout,"(14x,'gphit mini: ',1f10.2,' at i = ',i5,' j= ',i5)") zgpmin, imip(1), imip(2)
604         WRITE(numout,"(14x,'gphit maxi: ',1f10.2,' at i = ',i5,' j= ',i5)") zgpmax, imap(1), imap(2)
605         WRITE(numout,"(14x,'  e1t mini: ',1f10.2,' at i = ',i5,' j= ',i5)") ze1min, imi1(1), imi1(2)
606         WRITE(numout,"(14x,'  e1t maxi: ',1f10.2,' at i = ',i5,' j= ',i5)") ze1max, ima1(1), ima1(2)
607         WRITE(numout,"(14x,'  e2t mini: ',1f10.2,' at i = ',i5,' j= ',i5)") ze2min, imi2(1), imi2(2)
608         WRITE(numout,"(14x,'  e2t maxi: ',1f10.2,' at i = ',i5,' j= ',i5)") ze2max, ima2(1), ima2(2)
609      ENDIF
610      !
611   END SUBROUTINE dom_ctl
612
613
614   SUBROUTINE domain_cfg( cd_cfg, kk_cfg, kpi, kpj, kpk, kperio )
615      !!----------------------------------------------------------------------
616      !!                     ***  ROUTINE domain_cfg  ***
617      !!
618      !! ** Purpose :   read the domain size in domain configuration file
619      !!
620      !! ** Method  :   read the cn_domcfg NetCDF file
621      !!----------------------------------------------------------------------
622      CHARACTER(len=*)              , INTENT(out) ::   cd_cfg          ! configuration name
623      INTEGER                       , INTENT(out) ::   kk_cfg          ! configuration resolution
624      INTEGER                       , INTENT(out) ::   kpi, kpj, kpk   ! global domain sizes
625      INTEGER                       , INTENT(out) ::   kperio          ! lateral global domain b.c.
626      !
627      INTEGER ::   inum   ! local integer
628      REAL(wp) ::   zorca_res                     ! local scalars
629      REAL(wp) ::   zperio                        !   -      -
630      INTEGER, DIMENSION(4) ::   idvar, idimsz    ! size   of dimensions
631      !!----------------------------------------------------------------------
632      !
633      IF(lwp) THEN
634         WRITE(numout,*) '           '
635         WRITE(numout,*) 'domain_cfg : domain size read in ', TRIM( cn_domcfg ), ' file'
636         WRITE(numout,*) '~~~~~~~~~~ '
637      ENDIF
638      !
639      CALL iom_open( cn_domcfg, inum )
640      !
641      !                                   !- ORCA family specificity
642      IF(  iom_varid( inum, 'ORCA'       , ldstop = .FALSE. ) > 0  .AND.  &
643         & iom_varid( inum, 'ORCA_index' , ldstop = .FALSE. ) > 0    ) THEN
644         !
645         cd_cfg = 'ORCA'
646         CALL iom_get( inum, 'ORCA_index', zorca_res )   ;   kk_cfg = NINT( zorca_res )
647         !
648         IF(lwp) THEN
649            WRITE(numout,*) '   .'
650            WRITE(numout,*) '   ==>>>   ORCA configuration '
651            WRITE(numout,*) '   .'
652         ENDIF
653         !
654      ELSE                                !- cd_cfg & k_cfg are not used
655         cd_cfg = 'UNKNOWN'
656         kk_cfg = -9999999
657                                          !- or they may be present as global attributes
658                                          !- (netcdf only)
659         CALL iom_getatt( inum, 'cn_cfg', cd_cfg )  ! returns   !  if not found
660         CALL iom_getatt( inum, 'nn_cfg', kk_cfg )  ! returns -999 if not found
661         IF( TRIM(cd_cfg) == '!') cd_cfg = 'UNKNOWN'
662         IF( kk_cfg == -999     ) kk_cfg = -9999999
663         !
664      ENDIF
665       !
666      idvar = iom_varid( inum, 'e3t_0', kdimsz = idimsz )   ! use e3t_0, that must exist, to get jp(ijk)glo
667      kpi = idimsz(1)
668      kpj = idimsz(2)
669      kpk = idimsz(3)
670      CALL iom_get( inum, 'jperio', zperio )   ;   kperio = NINT( zperio )
671      CALL iom_close( inum )
672      !
673      IF(lwp) THEN
674         WRITE(numout,*) '      cn_cfg = ', TRIM(cd_cfg), '   nn_cfg = ', kk_cfg
675         WRITE(numout,*) '      Ni0glo = ', kpi
676         WRITE(numout,*) '      Nj0glo = ', kpj
677         WRITE(numout,*) '      jpkglo = ', kpk
678         WRITE(numout,*) '      type of global domain lateral boundary   jperio = ', kperio
679      ENDIF
680      !
681   END SUBROUTINE domain_cfg
682
683
684   SUBROUTINE cfg_write
685      !!----------------------------------------------------------------------
686      !!                  ***  ROUTINE cfg_write  ***
687      !!
688      !! ** Purpose :   Create the "cn_domcfg_out" file, a NetCDF file which
689      !!              contains all the ocean domain informations required to
690      !!              define an ocean configuration.
691      !!
692      !! ** Method  :   Write in a file all the arrays required to set up an
693      !!              ocean configuration.
694      !!
695      !! ** output file :   domcfg_out.nc : domain size, characteristics, horizontal
696      !!                       mesh, Coriolis parameter, and vertical scale factors
697      !!                    NB: also contain ORCA family information
698      !!----------------------------------------------------------------------
699      INTEGER           ::   ji, jj, jk   ! dummy loop indices
700      INTEGER           ::   inum     ! local units
701      CHARACTER(len=21) ::   clnam    ! filename (mesh and mask informations)
702      REAL(wp), DIMENSION(jpi,jpj) ::   z2d   ! workspace
703      !!----------------------------------------------------------------------
704      !
705      IF(lwp) WRITE(numout,*)
706      IF(lwp) WRITE(numout,*) 'cfg_write : create the domain configuration file (', TRIM(cn_domcfg_out),'.nc)'
707      IF(lwp) WRITE(numout,*) '~~~~~~~~~'
708      !
709      !                       ! ============================= !
710      !                       !  create 'domcfg_out.nc' file  !
711      !                       ! ============================= !
712      !
713      clnam = cn_domcfg_out  ! filename (configuration information)
714      CALL iom_open( TRIM(clnam), inum, ldwrt = .TRUE. )
715      !
716      !                             !==  ORCA family specificities  ==!
717      IF( TRIM(cn_cfg) == "orca" .OR. TRIM(cn_cfg) == "ORCA" ) THEN
718         CALL iom_rstput( 0, 0, inum, 'ORCA'      , 1._wp            , ktype = jp_i4 )
719         CALL iom_rstput( 0, 0, inum, 'ORCA_index', REAL( nn_cfg, wp), ktype = jp_i4 )
720      ENDIF
721      !
722      !                             !==  domain characteristics  ==!
723      !
724      !                                   ! lateral boundary of the global domain
725      CALL iom_rstput( 0, 0, inum, 'jperio', REAL( jperio, wp), ktype = jp_i4 )
726      !
727      !                                   ! type of vertical coordinate
728      CALL iom_rstput( 0, 0, inum, 'ln_zco', REAL(COUNT((/ln_zco/)), wp), ktype = jp_i4 )
729      CALL iom_rstput( 0, 0, inum, 'ln_zps', REAL(COUNT((/ln_zps/)), wp), ktype = jp_i4 )
730      CALL iom_rstput( 0, 0, inum, 'ln_sco', REAL(COUNT((/ln_sco/)), wp), ktype = jp_i4 )
731      !
732      !                                   ! ocean cavities under iceshelves
733      CALL iom_rstput( 0, 0, inum, 'ln_isfcav', REAL(COUNT((/ln_isfcav/)), wp), ktype = jp_i4 )
734      !
735      !                             !==  horizontal mesh  !
736      !
737      CALL iom_rstput( 0, 0, inum, 'glamt', glamt, ktype = jp_r8 )   ! latitude
738      CALL iom_rstput( 0, 0, inum, 'glamu', glamu, ktype = jp_r8 )
739      CALL iom_rstput( 0, 0, inum, 'glamv', glamv, ktype = jp_r8 )
740      CALL iom_rstput( 0, 0, inum, 'glamf', glamf, ktype = jp_r8 )
741      !
742      CALL iom_rstput( 0, 0, inum, 'gphit', gphit, ktype = jp_r8 )   ! longitude
743      CALL iom_rstput( 0, 0, inum, 'gphiu', gphiu, ktype = jp_r8 )
744      CALL iom_rstput( 0, 0, inum, 'gphiv', gphiv, ktype = jp_r8 )
745      CALL iom_rstput( 0, 0, inum, 'gphif', gphif, ktype = jp_r8 )
746      !
747      CALL iom_rstput( 0, 0, inum, 'e1t'  , e1t  , ktype = jp_r8 )   ! i-scale factors (e1.)
748      CALL iom_rstput( 0, 0, inum, 'e1u'  , e1u  , ktype = jp_r8 )
749      CALL iom_rstput( 0, 0, inum, 'e1v'  , e1v  , ktype = jp_r8 )
750      CALL iom_rstput( 0, 0, inum, 'e1f'  , e1f  , ktype = jp_r8 )
751      !
752      CALL iom_rstput( 0, 0, inum, 'e2t'  , e2t  , ktype = jp_r8 )   ! j-scale factors (e2.)
753      CALL iom_rstput( 0, 0, inum, 'e2u'  , e2u  , ktype = jp_r8 )
754      CALL iom_rstput( 0, 0, inum, 'e2v'  , e2v  , ktype = jp_r8 )
755      CALL iom_rstput( 0, 0, inum, 'e2f'  , e2f  , ktype = jp_r8 )
756      !
757      CALL iom_rstput( 0, 0, inum, 'ff_f' , ff_f , ktype = jp_r8 )   ! coriolis factor
758      CALL iom_rstput( 0, 0, inum, 'ff_t' , ff_t , ktype = jp_r8 )
759      !
760      !                             !==  vertical mesh  ==!
761      !
762      CALL iom_rstput( 0, 0, inum, 'e3t_1d'  , e3t_1d , ktype = jp_r8 )   ! reference 1D-coordinate
763      CALL iom_rstput( 0, 0, inum, 'e3w_1d'  , e3w_1d , ktype = jp_r8 )
764      !
765      CALL iom_rstput( 0, 0, inum, 'e3t_0'   , e3t_0  , ktype = jp_r8 )   ! vertical scale factors
766      CALL iom_rstput( 0, 0, inum, 'e3u_0'   , e3u_0  , ktype = jp_r8 )
767      CALL iom_rstput( 0, 0, inum, 'e3v_0'   , e3v_0  , ktype = jp_r8 )
768      CALL iom_rstput( 0, 0, inum, 'e3f_0'   , e3f_0  , ktype = jp_r8 )
769      CALL iom_rstput( 0, 0, inum, 'e3w_0'   , e3w_0  , ktype = jp_r8 )
770      CALL iom_rstput( 0, 0, inum, 'e3uw_0'  , e3uw_0 , ktype = jp_r8 )
771      CALL iom_rstput( 0, 0, inum, 'e3vw_0'  , e3vw_0 , ktype = jp_r8 )
772      !
773      !                             !==  wet top and bottom level  ==!   (caution: multiplied by ssmask)
774      !
775      CALL iom_rstput( 0, 0, inum, 'top_level'    , REAL( mikt, wp )*ssmask , ktype = jp_i4 )   ! nb of ocean T-points (ISF)
776      CALL iom_rstput( 0, 0, inum, 'bottom_level' , REAL( mbkt, wp )*ssmask , ktype = jp_i4 )   ! nb of ocean T-points
777      !
778      IF( ln_sco ) THEN             ! s-coordinate: store grid stiffness ratio  (Not required anyway)
779         CALL dom_stiff( z2d )
780         CALL iom_rstput( 0, 0, inum, 'stiffness', z2d )        !    ! Max. grid stiffness ratio
781      ENDIF
782      !
783      IF( ll_wd ) THEN              ! wetting and drying domain
784         CALL iom_rstput( 0, 0, inum, 'ht_0'   , ht_0   , ktype = jp_r8 )
785      ENDIF
786      !
787      ! Add some global attributes ( netcdf only )
788      CALL iom_putatt( inum, 'nn_cfg', nn_cfg )
789      CALL iom_putatt( inum, 'cn_cfg', TRIM(cn_cfg) )
790      !
791      !                                ! ============================
792      !                                !        close the files
793      !                                ! ============================
794      CALL iom_close( inum )
795      !
796   END SUBROUTINE cfg_write
797
798   !!======================================================================
799END MODULE domain
Note: See TracBrowser for help on using the repository browser.