source: NEMO/branches/UKMO/dev_r12745_HPC-02_Daley_Tiling_trial_public/src/OCE/DOM/domain.F90 @ 12879

Last change on this file since 12879 was 12879, checked in by hadcv, 9 months ago

Changes following feedback

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