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/branches/2021/ticket2680_C1D_PAPA/src/OCE/DOM – NEMO

source: NEMO/branches/2021/ticket2680_C1D_PAPA/src/OCE/DOM/domain.F90 @ 15008

Last change on this file since 15008 was 15008, checked in by gsamson, 3 years ago

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