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

Last change on this file since 14223 was 14223, checked in by smasson, 3 years ago

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