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/UKMO/dev_r12745_HPC-02_Daley_Tiling_trial_structure/src/OCE/DOM – NEMO

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

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

tra_ldf_iso trial using structures

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