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/2020/dev_r12527_Gurvan_ShallowWater/src/SWE – NEMO

source: NEMO/branches/2020/dev_r12527_Gurvan_ShallowWater/src/SWE/domain.F90 @ 12614

Last change on this file since 12614 was 12614, checked in by gm, 4 years ago

first Shallow Water Eq. update

  • Property svn:keywords set to Id
File size: 36.4 KB
Line 
1MODULE domain
2   !!==============================================================================
3   !!                       ***  MODULE domain   ***
4   !! Ocean initialization : domain initialization
5   !!==============================================================================
6   !! History :  OPA  !  1990-10  (C. Levy - G. Madec)  Original code
7   !!                 !  1992-01  (M. Imbard) insert time step initialization
8   !!                 !  1996-06  (G. Madec) generalized vertical coordinate
9   !!                 !  1997-02  (G. Madec) creation of domwri.F
10   !!                 !  2001-05  (E.Durand - G. Madec) insert closed sea
11   !!   NEMO     1.0  !  2002-08  (G. Madec)  F90: Free form and module
12   !!            2.0  !  2005-11  (V. Garnier) Surface pressure gradient organization
13   !!            3.3  !  2010-11  (G. Madec)  initialisation in C1D configuration
14   !!            3.6  !  2013     ( J. Simeon, C. Calone, G. Madec, C. Ethe ) Online coarsening of outputs
15   !!            3.7  !  2015-11  (G. Madec, A. Coward)  time varying zgr by default
16   !!            4.0  !  2016-10  (G. Madec, S. Flavoni)  domain configuration / user defined interface
17   !!----------------------------------------------------------------------
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      !
125      IF( lwxios ) THEN
126!define names for restart write and set core output (restart.F90)
127         CALL iom_set_rst_vars(rst_wfields)
128         CALL iom_set_rstw_core(cdstr)
129      ENDIF
130!reset namelist for SAS
131      IF(cdstr == 'SAS') THEN
132         IF(lrxios) THEN
133               IF(lwp) write(numout,*) 'Disable reading restart file using XIOS for SAS'
134               lrxios = .FALSE.
135         ENDIF
136      ENDIF
137      !
138      CALL dom_hgr                      ! Horizontal mesh
139
140      IF( ln_closea ) CALL dom_clo      ! Read in masks to define closed seas and lakes
141
142      CALL dom_zgr( ik_top, ik_bot )    ! Vertical mesh and bathymetry
143
144      CALL dom_msk( ik_top, ik_bot )    ! Masks
145
146      ht_0(:,:) = 0._wp                 ! Reference ocean thickness
147      hu_0(:,:) = 0._wp
148      hv_0(:,:) = 0._wp
149      hf_0(:,:) = 0._wp
150      DO jk = 1, jpkm1
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         hf_0(:,:) = hf_0(:,:) + e3f_0(:,:,jk) * ssfmask(:,:)     ! CAUTION : only valid in SWE, not with bathymetry
155      END DO
156      !                                 ! Inverse of reference ocean thickness
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      !           !==  time varying part of coordinate system  ==!
163      !
164      IF( ln_linssh ) THEN       != Fix in time : set to the reference one for all
165      !
166         !       before        !          now          !       after         !
167            gdept(:,:,:,Kbb) = gdept_0  ;   gdept(:,:,:,Kmm) = gdept_0   ;   gdept(:,:,:,Kaa) = gdept_0   ! depth of grid-points
168            gdepw(:,:,:,Kbb) = gdepw_0  ;   gdepw(:,:,:,Kmm) = gdepw_0   ;   gdepw(:,:,:,Kaa) = gdepw_0   !
169                                   gde3w = gde3w_0   !        ---          !
170         !                                                                 
171              e3t(:,:,:,Kbb) =   e3t_0  ;     e3t(:,:,:,Kmm) =   e3t_0   ;   e3t(:,:,:,Kaa) =  e3t_0    ! scale factors
172              e3u(:,:,:,Kbb) =   e3u_0  ;     e3u(:,:,:,Kmm) =   e3u_0   ;   e3u(:,:,:,Kaa) =  e3u_0    !
173              e3v(:,:,:,Kbb) =   e3v_0  ;     e3v(:,:,:,Kmm) =   e3v_0   ;   e3v(:,:,:,Kaa) =  e3v_0    !
174                                     e3f =   e3f_0   !        ---          !
175              e3w(:,:,:,Kbb) =   e3w_0  ;     e3w(:,:,:,Kmm) =   e3w_0   ;    e3w(:,:,:,Kaa) =   e3w_0   !
176             e3uw(:,:,:,Kbb) =  e3uw_0  ;    e3uw(:,:,:,Kmm) =  e3uw_0   ;   e3uw(:,:,:,Kaa) =  e3uw_0   
177             e3vw(:,:,:,Kbb) =  e3vw_0  ;    e3vw(:,:,:,Kmm) =  e3vw_0   ;   e3vw(:,:,:,Kaa) =  e3vw_0   !
178         !
179         z1_hu_0(:,:) = ssumask(:,:) / ( hu_0(:,:) + 1._wp - ssumask(:,:) )     ! _i mask due to ISF
180         z1_hv_0(:,:) = ssvmask(:,:) / ( hv_0(:,:) + 1._wp - ssvmask(:,:) )
181         !
182         !        before       !          now          !       after         !
183                                      ht =    ht_0   !                     ! water column thickness
184               hu(:,:,Kbb) =    hu_0  ;      hu(:,:,Kmm) =    hu_0   ;    hu(:,:,Kaa) =    hu_0   !
185               hv(:,:,Kbb) =    hv_0  ;      hv(:,:,Kmm) =    hv_0   ;    hv(:,:,Kaa) =    hv_0   !
186            r1_hu(:,:,Kbb) = z1_hu_0  ;   r1_hu(:,:,Kmm) = z1_hu_0   ; r1_hu(:,:,Kaa) = z1_hu_0   ! inverse of water column thickness
187            r1_hv(:,:,Kbb) = z1_hv_0  ;   r1_hv(:,:,Kmm) = z1_hv_0   ; r1_hv(:,:,Kaa) = z1_hv_0   !
188         !
189         !
190      ELSE                       != time varying : initialize before/now/after variables
191         !
192         IF( .NOT.l_offline )  CALL dom_vvl_init( Kbb, Kmm, Kaa )
193         !
194      ENDIF
195      !
196      IF( lk_c1d         )   CALL cor_c1d       ! 1D configuration: Coriolis set at T-point
197      !
198      IF( ln_meshmask    )   CALL dom_wri       ! Create a domain file
199
200      IF( .NOT.ln_rstart )   CALL dom_ctl       ! Domain control
201      !
202      IF( ln_write_cfg   )   CALL cfg_write     ! create the configuration file
203      !
204      IF(lwp) THEN
205         WRITE(numout,*)
206         WRITE(numout,*) 'dom_init :   ==>>>   END of domain initialization'
207         WRITE(numout,*) '~~~~~~~~'
208         WRITE(numout,*) 
209      ENDIF
210     
211      !
212   END SUBROUTINE dom_init
213
214
215   SUBROUTINE dom_glo
216      !!----------------------------------------------------------------------
217      !!                     ***  ROUTINE dom_glo  ***
218      !!
219      !! ** Purpose :   initialization of global domain <--> local domain indices
220      !!
221      !! ** Method  :   
222      !!
223      !! ** Action  : - mig , mjg : local  domain indices ==> global domain indices
224      !!              - mi0 , mi1 : global domain indices ==> local  domain indices
225      !!              - mj0,, mj1   (global point not in the local domain ==> mi0>mi1 and/or mj0>mj1)
226      !!----------------------------------------------------------------------
227      INTEGER ::   ji, jj   ! dummy loop argument
228      !!----------------------------------------------------------------------
229      !
230      DO ji = 1, jpi                 ! local domain indices ==> global domain indices
231        mig(ji) = ji + nimpp - 1
232      END DO
233      DO jj = 1, jpj
234        mjg(jj) = jj + njmpp - 1
235      END DO
236      !                              ! global domain indices ==> local domain indices
237      !                                   ! (return (m.0,m.1)=(1,0) if data domain gridpoint is to the west/south of the
238      !                                   ! local domain, or (m.0,m.1)=(jp.+1,jp.) to the east/north of local domain.
239      DO ji = 1, jpiglo
240        mi0(ji) = MAX( 1 , MIN( ji - nimpp + 1, jpi+1 ) )
241        mi1(ji) = MAX( 0 , MIN( ji - nimpp + 1, jpi   ) )
242      END DO
243      DO jj = 1, jpjglo
244        mj0(jj) = MAX( 1 , MIN( jj - njmpp + 1, jpj+1 ) )
245        mj1(jj) = MAX( 0 , MIN( jj - njmpp + 1, jpj   ) )
246      END DO
247      IF(lwp) THEN                   ! control print
248         WRITE(numout,*)
249         WRITE(numout,*) 'dom_glo : domain: global <<==>> local '
250         WRITE(numout,*) '~~~~~~~ '
251         WRITE(numout,*) '   global domain:   jpiglo = ', jpiglo, ' jpjglo = ', jpjglo, ' jpkglo = ', jpkglo
252         WRITE(numout,*) '   local  domain:   jpi    = ', jpi   , ' jpj    = ', jpj   , ' jpk    = ', jpk
253         WRITE(numout,*)
254         WRITE(numout,*) '   conversion from local to global domain indices (and vise versa) done'
255         IF( nn_print >= 1 ) THEN
256            WRITE(numout,*)
257            WRITE(numout,*) '          conversion local  ==> global i-index domain (mig)'
258            WRITE(numout,25)              (mig(ji),ji = 1,jpi)
259            WRITE(numout,*)
260            WRITE(numout,*) '          conversion global ==> local  i-index domain'
261            WRITE(numout,*) '             starting index (mi0)'
262            WRITE(numout,25)              (mi0(ji),ji = 1,jpiglo)
263            WRITE(numout,*) '             ending index (mi1)'
264            WRITE(numout,25)              (mi1(ji),ji = 1,jpiglo)
265            WRITE(numout,*)
266            WRITE(numout,*) '          conversion local  ==> global j-index domain (mjg)'
267            WRITE(numout,25)              (mjg(jj),jj = 1,jpj)
268            WRITE(numout,*)
269            WRITE(numout,*) '          conversion global ==> local  j-index domain'
270            WRITE(numout,*) '             starting index (mj0)'
271            WRITE(numout,25)              (mj0(jj),jj = 1,jpjglo)
272            WRITE(numout,*) '             ending index (mj1)'
273            WRITE(numout,25)              (mj1(jj),jj = 1,jpjglo)
274         ENDIF
275      ENDIF
276 25   FORMAT( 100(10x,19i4,/) )
277      !
278   END SUBROUTINE dom_glo
279
280
281   SUBROUTINE dom_nam
282      !!----------------------------------------------------------------------
283      !!                     ***  ROUTINE dom_nam  ***
284      !!                   
285      !! ** Purpose :   read domaine namelists and print the variables.
286      !!
287      !! ** input   : - namrun namelist
288      !!              - namdom namelist
289      !!              - namnc4 namelist   ! "key_netcdf4" only
290      !!----------------------------------------------------------------------
291      USE ioipsl
292      !!
293      INTEGER  ::   ios   ! Local integer
294      !
295      NAMELIST/namrun/ cn_ocerst_indir, cn_ocerst_outdir, nn_stocklist, ln_rst_list,                 &
296         &             nn_no   , cn_exp   , cn_ocerst_in, cn_ocerst_out, ln_rstart , nn_rstctl ,     &
297         &             nn_it000, nn_itend , nn_date0    , nn_time0     , nn_leapy  , nn_istate ,     &
298         &             nn_stock, nn_write , ln_mskland  , ln_clobber   , nn_chunksz, ln_1st_euler  , &
299         &             ln_cfmeta, ln_xios_read, nn_wxios
300      NAMELIST/namdom/ ln_linssh, rn_Dt, rn_atfp, ln_crs, ln_meshmask
301#if defined key_netcdf4
302      NAMELIST/namnc4/ nn_nchunks_i, nn_nchunks_j, nn_nchunks_k, ln_nc4zip
303#endif
304      !!----------------------------------------------------------------------
305      !
306      IF(lwp) THEN
307         WRITE(numout,*)
308         WRITE(numout,*) 'dom_nam : domain initialization through namelist read'
309         WRITE(numout,*) '~~~~~~~ '
310      ENDIF
311      !
312      !
313      READ  ( numnam_ref, namrun, IOSTAT = ios, ERR = 901)
314901   IF( ios /= 0 )   CALL ctl_nam ( ios , 'namrun in reference namelist' )
315      READ  ( numnam_cfg, namrun, IOSTAT = ios, ERR = 902 )
316902   IF( ios >  0 )   CALL ctl_nam ( ios , 'namrun in configuration namelist' )
317      IF(lwm) WRITE ( numond, namrun )
318      !
319      IF(lwp) THEN                  ! control print
320         WRITE(numout,*) '   Namelist : namrun   ---   run parameters'
321         WRITE(numout,*) '      Assimilation cycle              nn_no           = ', nn_no
322         WRITE(numout,*) '      experiment name for output      cn_exp          = ', TRIM( cn_exp           )
323         WRITE(numout,*) '      file prefix restart input       cn_ocerst_in    = ', TRIM( cn_ocerst_in     )
324         WRITE(numout,*) '      restart input directory         cn_ocerst_indir = ', TRIM( cn_ocerst_indir  )
325         WRITE(numout,*) '      file prefix restart output      cn_ocerst_out   = ', TRIM( cn_ocerst_out    )
326         WRITE(numout,*) '      restart output directory        cn_ocerst_outdir= ', TRIM( cn_ocerst_outdir )
327         WRITE(numout,*) '      restart logical                 ln_rstart       = ', ln_rstart
328         WRITE(numout,*) '      start with forward time step    ln_1st_euler    = ', ln_1st_euler
329         WRITE(numout,*) '      control of time step            nn_rstctl       = ', nn_rstctl
330         WRITE(numout,*) '      number of the first time step   nn_it000        = ', nn_it000
331         WRITE(numout,*) '      number of the last time step    nn_itend        = ', nn_itend
332         WRITE(numout,*) '      initial calendar date aammjj    nn_date0        = ', nn_date0
333         WRITE(numout,*) '      initial time of day in hhmm     nn_time0        = ', nn_time0
334         WRITE(numout,*) '      leap year calendar (0/1)        nn_leapy        = ', nn_leapy
335         WRITE(numout,*) '      initial state output            nn_istate       = ', nn_istate
336         IF( ln_rst_list ) THEN
337            WRITE(numout,*) '      list of restart dump times      nn_stocklist    =', nn_stocklist
338         ELSE
339            WRITE(numout,*) '      frequency of restart file       nn_stock        = ', nn_stock
340         ENDIF
341#if ! defined key_iomput
342         WRITE(numout,*) '      frequency of output file        nn_write        = ', nn_write
343#endif
344         WRITE(numout,*) '      mask land points                ln_mskland      = ', ln_mskland
345         WRITE(numout,*) '      additional CF standard metadata ln_cfmeta       = ', ln_cfmeta
346         WRITE(numout,*) '      overwrite an existing file      ln_clobber      = ', ln_clobber
347         WRITE(numout,*) '      NetCDF chunksize (bytes)        nn_chunksz      = ', nn_chunksz
348         IF( TRIM(Agrif_CFixed()) == '0' ) THEN
349            WRITE(numout,*) '      READ restart for a single file using XIOS ln_xios_read =', ln_xios_read
350            WRITE(numout,*) '      Write restart using XIOS        nn_wxios   = ', nn_wxios
351         ELSE
352            WRITE(numout,*) "      AGRIF: nn_wxios will be ingored. See setting for parent"
353            WRITE(numout,*) "      AGRIF: ln_xios_read will be ingored. See setting for parent"
354         ENDIF
355      ENDIF
356
357      cexper = cn_exp         ! conversion DOCTOR names into model names (this should disappear soon)
358      nrstdt = nn_rstctl
359      nit000 = nn_it000
360      nitend = nn_itend
361      ndate0 = nn_date0
362      nleapy = nn_leapy
363      ninist = nn_istate
364      l_1st_euler = ln_1st_euler
365      IF( .NOT. l_1st_euler .AND. .NOT. ln_rstart ) THEN
366         IF(lwp) WRITE(numout,*) 
367         IF(lwp) WRITE(numout,*)'   ==>>>   Start from rest (ln_rstart=F)'
368         IF(lwp) WRITE(numout,*)'           an Euler initial time step is used : l_1st_euler is forced to .true. '   
369         l_1st_euler = .true.
370      ENDIF
371      !                             ! control of output frequency
372      IF( .NOT. ln_rst_list ) THEN     ! we use nn_stock
373         IF( nn_stock == -1 )   CALL ctl_warn( 'nn_stock = -1 --> no restart will be done' )
374         IF( nn_stock == 0 .OR. nn_stock > nitend ) THEN
375            WRITE(ctmp1,*) 'nn_stock = ', nn_stock, ' it is forced to ', nitend
376            CALL ctl_warn( ctmp1 )
377            nn_stock = nitend
378         ENDIF
379      ENDIF
380#if ! defined key_iomput
381      IF( nn_write == -1 )   CALL ctl_warn( 'nn_write = -1 --> no output files will be done' )
382      IF ( nn_write == 0 ) THEN
383         WRITE(ctmp1,*) 'nn_write = ', nn_write, ' it is forced to ', nitend
384         CALL ctl_warn( ctmp1 )
385         nn_write = nitend
386      ENDIF
387#endif
388
389#if defined key_agrif
390      IF( Agrif_Root() ) THEN
391#endif
392      IF(lwp) WRITE(numout,*)
393      SELECT CASE ( nleapy )        ! Choose calendar for IOIPSL
394      CASE (  1 ) 
395         CALL ioconf_calendar('gregorian')
396         IF(lwp) WRITE(numout,*) '   ==>>>   The IOIPSL calendar is "gregorian", i.e. leap year'
397      CASE (  0 )
398         CALL ioconf_calendar('noleap')
399         IF(lwp) WRITE(numout,*) '   ==>>>   The IOIPSL calendar is "noleap", i.e. no leap year'
400      CASE ( 30 )
401         CALL ioconf_calendar('360d')
402         IF(lwp) WRITE(numout,*) '   ==>>>   The IOIPSL calendar is "360d", i.e. 360 days in a year'
403      END SELECT
404#if defined key_agrif
405      ENDIF
406#endif
407
408      READ  ( numnam_ref, namdom, IOSTAT = ios, ERR = 903)
409903   IF( ios /= 0 )   CALL ctl_nam ( ios , 'namdom in reference namelist' )
410      READ  ( numnam_cfg, namdom, IOSTAT = ios, ERR = 904 )
411904   IF( ios >  0 )   CALL ctl_nam ( ios , 'namdom in configuration namelist' )
412      IF(lwm) WRITE( numond, namdom )
413      !
414      IF(lwp) THEN
415         WRITE(numout,*)
416         WRITE(numout,*) '   Namelist : namdom   ---   space & time domain'
417         WRITE(numout,*) '      linear free surface (=T)                ln_linssh   = ', ln_linssh
418         WRITE(numout,*) '      create mesh/mask file                   ln_meshmask = ', ln_meshmask
419         WRITE(numout,*) '      ocean time step                         rn_Dt       = ', rn_Dt
420         WRITE(numout,*) '      asselin time filter parameter           rn_atfp     = ', rn_atfp
421         WRITE(numout,*) '      online coarsening of dynamical fields   ln_crs      = ', ln_crs
422      ENDIF
423      !
424      !! Initialise current model timestep rDt = 2*rn_Dt if MLF or rDt = rn_Dt if RK3
425      rDt  = 2._wp * rn_Dt
426      r1_Dt = 1._wp / rDt
427
428      IF( TRIM(Agrif_CFixed()) == '0' ) THEN
429         lrxios = ln_xios_read.AND.ln_rstart
430!set output file type for XIOS based on NEMO namelist
431         IF (nn_wxios > 0) lwxios = .TRUE. 
432         nxioso = nn_wxios
433      ENDIF
434
435#if defined key_netcdf4
436      !                             ! NetCDF 4 case   ("key_netcdf4" defined)
437      READ  ( numnam_ref, namnc4, IOSTAT = ios, ERR = 907)
438907   IF( ios /= 0 )   CALL ctl_nam ( ios , 'namnc4 in reference namelist' )
439      READ  ( numnam_cfg, namnc4, IOSTAT = ios, ERR = 908 )
440908   IF( ios >  0 )   CALL ctl_nam ( ios , 'namnc4 in configuration namelist' )
441      IF(lwm) WRITE( numond, namnc4 )
442
443      IF(lwp) THEN                        ! control print
444         WRITE(numout,*)
445         WRITE(numout,*) '   Namelist namnc4 - Netcdf4 chunking parameters'
446         WRITE(numout,*) '      number of chunks in i-dimension             nn_nchunks_i = ', nn_nchunks_i
447         WRITE(numout,*) '      number of chunks in j-dimension             nn_nchunks_j = ', nn_nchunks_j
448         WRITE(numout,*) '      number of chunks in k-dimension             nn_nchunks_k = ', nn_nchunks_k
449         WRITE(numout,*) '      apply netcdf4/hdf5 chunking & compression   ln_nc4zip    = ', ln_nc4zip
450      ENDIF
451
452      ! Put the netcdf4 settings into a simple structure (snc4set, defined in in_out_manager module)
453      ! Note the chunk size in the unlimited (time) dimension will be fixed at 1
454      snc4set%ni   = nn_nchunks_i
455      snc4set%nj   = nn_nchunks_j
456      snc4set%nk   = nn_nchunks_k
457      snc4set%luse = ln_nc4zip
458#else
459      snc4set%luse = .FALSE.        ! No NetCDF 4 case
460#endif
461      !
462   END SUBROUTINE dom_nam
463
464
465   SUBROUTINE dom_ctl
466      !!----------------------------------------------------------------------
467      !!                     ***  ROUTINE dom_ctl  ***
468      !!
469      !! ** Purpose :   Domain control.
470      !!
471      !! ** Method  :   compute and print extrema of masked scale factors
472      !!----------------------------------------------------------------------
473      INTEGER, DIMENSION(2) ::   imi1, imi2, ima1, ima2
474      INTEGER, DIMENSION(2) ::   iloc   !
475      REAL(wp) ::   ze1min, ze1max, ze2min, ze2max
476      !!----------------------------------------------------------------------
477      !
478      IF(lk_mpp) THEN
479         CALL mpp_minloc( 'domain', e1t(:,:), tmask_i(:,:), ze1min, imi1 )
480         CALL mpp_minloc( 'domain', e2t(:,:), tmask_i(:,:), ze2min, imi2 )
481         CALL mpp_maxloc( 'domain', e1t(:,:), tmask_i(:,:), ze1max, ima1 )
482         CALL mpp_maxloc( 'domain', e2t(:,:), tmask_i(:,:), ze2max, ima2 )
483      ELSE
484         ze1min = MINVAL( e1t(:,:), mask = tmask_i(:,:) == 1._wp )   
485         ze2min = MINVAL( e2t(:,:), mask = tmask_i(:,:) == 1._wp )   
486         ze1max = MAXVAL( e1t(:,:), mask = tmask_i(:,:) == 1._wp )   
487         ze2max = MAXVAL( e2t(:,:), mask = tmask_i(:,:) == 1._wp )   
488         !
489         iloc  = MINLOC( e1t(:,:), mask = tmask_i(:,:) == 1._wp )
490         imi1(1) = iloc(1) + nimpp - 1
491         imi1(2) = iloc(2) + njmpp - 1
492         iloc  = MINLOC( e2t(:,:), mask = tmask_i(:,:) == 1._wp )
493         imi2(1) = iloc(1) + nimpp - 1
494         imi2(2) = iloc(2) + njmpp - 1
495         iloc  = MAXLOC( e1t(:,:), mask = tmask_i(:,:) == 1._wp )
496         ima1(1) = iloc(1) + nimpp - 1
497         ima1(2) = iloc(2) + njmpp - 1
498         iloc  = MAXLOC( e2t(:,:), mask = tmask_i(:,:) == 1._wp )
499         ima2(1) = iloc(1) + nimpp - 1
500         ima2(2) = iloc(2) + njmpp - 1
501      ENDIF
502      IF(lwp) THEN
503         WRITE(numout,*)
504         WRITE(numout,*) 'dom_ctl : extrema of the masked scale factors'
505         WRITE(numout,*) '~~~~~~~'
506         WRITE(numout,"(14x,'e1t maxi: ',1f10.2,' at i = ',i5,' j= ',i5)") ze1max, ima1(1), ima1(2)
507         WRITE(numout,"(14x,'e1t mini: ',1f10.2,' at i = ',i5,' j= ',i5)") ze1min, imi1(1), imi1(2)
508         WRITE(numout,"(14x,'e2t maxi: ',1f10.2,' at i = ',i5,' j= ',i5)") ze2max, ima2(1), ima2(2)
509         WRITE(numout,"(14x,'e2t mini: ',1f10.2,' at i = ',i5,' j= ',i5)") ze2min, imi2(1), imi2(2)
510      ENDIF
511      !
512   END SUBROUTINE dom_ctl
513
514
515   SUBROUTINE domain_cfg( cd_cfg, kk_cfg, kpi, kpj, kpk, kperio )
516      !!----------------------------------------------------------------------
517      !!                     ***  ROUTINE dom_nam  ***
518      !!                   
519      !! ** Purpose :   read the domain size in domain configuration file
520      !!
521      !! ** Method  :   read the cn_domcfg NetCDF file
522      !!----------------------------------------------------------------------
523      CHARACTER(len=*)              , INTENT(out) ::   cd_cfg          ! configuration name
524      INTEGER                       , INTENT(out) ::   kk_cfg          ! configuration resolution
525      INTEGER                       , INTENT(out) ::   kpi, kpj, kpk   ! global domain sizes
526      INTEGER                       , INTENT(out) ::   kperio          ! lateral global domain b.c.
527      !
528      INTEGER ::   inum   ! local integer
529      REAL(wp) ::   zorca_res                     ! local scalars
530      REAL(wp) ::   zperio                        !   -      -
531      INTEGER, DIMENSION(4) ::   idvar, idimsz    ! size   of dimensions
532      !!----------------------------------------------------------------------
533      !
534      IF(lwp) THEN
535         WRITE(numout,*) '           '
536         WRITE(numout,*) 'domain_cfg : domain size read in ', TRIM( cn_domcfg ), ' file'
537         WRITE(numout,*) '~~~~~~~~~~ '
538      ENDIF
539      !
540      CALL iom_open( cn_domcfg, inum )
541      !
542      !                                   !- ORCA family specificity
543      IF(  iom_varid( inum, 'ORCA'       , ldstop = .FALSE. ) > 0  .AND.  &
544         & iom_varid( inum, 'ORCA_index' , ldstop = .FALSE. ) > 0    ) THEN
545         !
546         cd_cfg = 'ORCA'
547         CALL iom_get( inum, 'ORCA_index', zorca_res )   ;   kk_cfg = NINT( zorca_res )
548         !
549         IF(lwp) THEN
550            WRITE(numout,*) '   .'
551            WRITE(numout,*) '   ==>>>   ORCA configuration '
552            WRITE(numout,*) '   .'
553         ENDIF
554         !
555      ELSE                                !- cd_cfg & k_cfg are not used
556         cd_cfg = 'UNKNOWN'
557         kk_cfg = -9999999
558                                          !- or they may be present as global attributes
559                                          !- (netcdf only) 
560         CALL iom_getatt( inum, 'cn_cfg', cd_cfg )  ! returns   !  if not found
561         CALL iom_getatt( inum, 'nn_cfg', kk_cfg )  ! returns -999 if not found
562         IF( TRIM(cd_cfg) == '!') cd_cfg = 'UNKNOWN'
563         IF( kk_cfg == -999     ) kk_cfg = -9999999
564         !
565      ENDIF
566       !
567      idvar = iom_varid( inum, 'e3t_0', kdimsz = idimsz )   ! use e3t_0, that must exist, to get jp(ijk)glo
568      kpi = idimsz(1)
569      kpj = idimsz(2)
570      kpk = idimsz(3)
571      CALL iom_get( inum, 'jperio', zperio )   ;   kperio = NINT( zperio )
572      CALL iom_close( inum )
573      !
574      IF(lwp) THEN
575         WRITE(numout,*) '      cn_cfg = ', TRIM(cd_cfg), '   nn_cfg = ', kk_cfg
576         WRITE(numout,*) '      jpiglo = ', kpi
577         WRITE(numout,*) '      jpjglo = ', kpj
578         WRITE(numout,*) '      jpkglo = ', kpk
579         WRITE(numout,*) '      type of global domain lateral boundary   jperio = ', kperio
580      ENDIF
581      !       
582   END SUBROUTINE domain_cfg
583   
584   
585   SUBROUTINE cfg_write
586      !!----------------------------------------------------------------------
587      !!                  ***  ROUTINE cfg_write  ***
588      !!                   
589      !! ** Purpose :   Create the "cn_domcfg_out" file, a NetCDF file which
590      !!              contains all the ocean domain informations required to
591      !!              define an ocean configuration.
592      !!
593      !! ** Method  :   Write in a file all the arrays required to set up an
594      !!              ocean configuration.
595      !!
596      !! ** output file :   domcfg_out.nc : domain size, characteristics, horizontal
597      !!                       mesh, Coriolis parameter, and vertical scale factors
598      !!                    NB: also contain ORCA family information
599      !!----------------------------------------------------------------------
600      INTEGER           ::   ji, jj, jk   ! dummy loop indices
601      INTEGER           ::   izco, izps, isco, icav
602      INTEGER           ::   inum     ! local units
603      CHARACTER(len=21) ::   clnam    ! filename (mesh and mask informations)
604      REAL(wp), DIMENSION(jpi,jpj) ::   z2d   ! workspace
605      !!----------------------------------------------------------------------
606      !
607      IF(lwp) WRITE(numout,*)
608      IF(lwp) WRITE(numout,*) 'cfg_write : create the domain configuration file (', TRIM(cn_domcfg_out),'.nc)'
609      IF(lwp) WRITE(numout,*) '~~~~~~~~~'
610      !
611      !                       ! ============================= !
612      !                       !  create 'domcfg_out.nc' file  !
613      !                       ! ============================= !
614      !         
615      clnam = cn_domcfg_out  ! filename (configuration information)
616      CALL iom_open( TRIM(clnam), inum, ldwrt = .TRUE. )
617     
618      !
619      !                             !==  ORCA family specificities  ==!
620      IF( cn_cfg == "ORCA" ) THEN
621         CALL iom_rstput( 0, 0, inum, 'ORCA'      , 1._wp            , ktype = jp_i4 )
622         CALL iom_rstput( 0, 0, inum, 'ORCA_index', REAL( nn_cfg, wp), ktype = jp_i4 )         
623      ENDIF
624      !
625      !                             !==  global domain size  ==!
626      !
627      CALL iom_rstput( 0, 0, inum, 'jpiglo', REAL( jpiglo, wp), ktype = jp_i4 )
628      CALL iom_rstput( 0, 0, inum, 'jpjglo', REAL( jpjglo, wp), ktype = jp_i4 )
629      CALL iom_rstput( 0, 0, inum, 'jpkglo', REAL( jpk   , wp), ktype = jp_i4 )
630      !
631      !                             !==  domain characteristics  ==!
632      !
633      !                                   ! lateral boundary of the global domain
634      CALL iom_rstput( 0, 0, inum, 'jperio', REAL( jperio, wp), ktype = jp_i4 )
635      !
636      !                                   ! type of vertical coordinate
637      IF( ln_zco    ) THEN   ;   izco = 1   ;   ELSE   ;   izco = 0   ;   ENDIF
638      IF( ln_zps    ) THEN   ;   izps = 1   ;   ELSE   ;   izps = 0   ;   ENDIF
639      IF( ln_sco    ) THEN   ;   isco = 1   ;   ELSE   ;   isco = 0   ;   ENDIF
640      CALL iom_rstput( 0, 0, inum, 'ln_zco'   , REAL( izco, wp), ktype = jp_i4 )
641      CALL iom_rstput( 0, 0, inum, 'ln_zps'   , REAL( izps, wp), ktype = jp_i4 )
642      CALL iom_rstput( 0, 0, inum, 'ln_sco'   , REAL( isco, wp), ktype = jp_i4 )
643      !
644      !                                   ! ocean cavities under iceshelves
645      IF( ln_isfcav ) THEN   ;   icav = 1   ;   ELSE   ;   icav = 0   ;   ENDIF
646      CALL iom_rstput( 0, 0, inum, 'ln_isfcav', REAL( icav, wp), ktype = jp_i4 )
647      !
648      !                             !==  horizontal mesh  !
649      !
650      CALL iom_rstput( 0, 0, inum, 'glamt', glamt, ktype = jp_r8 )   ! latitude
651      CALL iom_rstput( 0, 0, inum, 'glamu', glamu, ktype = jp_r8 )
652      CALL iom_rstput( 0, 0, inum, 'glamv', glamv, ktype = jp_r8 )
653      CALL iom_rstput( 0, 0, inum, 'glamf', glamf, ktype = jp_r8 )
654      !                               
655      CALL iom_rstput( 0, 0, inum, 'gphit', gphit, ktype = jp_r8 )   ! longitude
656      CALL iom_rstput( 0, 0, inum, 'gphiu', gphiu, ktype = jp_r8 )
657      CALL iom_rstput( 0, 0, inum, 'gphiv', gphiv, ktype = jp_r8 )
658      CALL iom_rstput( 0, 0, inum, 'gphif', gphif, ktype = jp_r8 )
659      !                               
660      CALL iom_rstput( 0, 0, inum, 'e1t'  , e1t  , ktype = jp_r8 )   ! i-scale factors (e1.)
661      CALL iom_rstput( 0, 0, inum, 'e1u'  , e1u  , ktype = jp_r8 )
662      CALL iom_rstput( 0, 0, inum, 'e1v'  , e1v  , ktype = jp_r8 )
663      CALL iom_rstput( 0, 0, inum, 'e1f'  , e1f  , ktype = jp_r8 )
664      !
665      CALL iom_rstput( 0, 0, inum, 'e2t'  , e2t  , ktype = jp_r8 )   ! j-scale factors (e2.)
666      CALL iom_rstput( 0, 0, inum, 'e2u'  , e2u  , ktype = jp_r8 )
667      CALL iom_rstput( 0, 0, inum, 'e2v'  , e2v  , ktype = jp_r8 )
668      CALL iom_rstput( 0, 0, inum, 'e2f'  , e2f  , ktype = jp_r8 )
669      !
670      CALL iom_rstput( 0, 0, inum, 'ff_f' , ff_f , ktype = jp_r8 )   ! coriolis factor
671      CALL iom_rstput( 0, 0, inum, 'ff_t' , ff_t , ktype = jp_r8 )
672      !
673      !                             !==  vertical mesh  ==!
674      !                                                     
675      CALL iom_rstput( 0, 0, inum, 'e3t_1d'  , e3t_1d , ktype = jp_r8 )   ! reference 1D-coordinate
676      CALL iom_rstput( 0, 0, inum, 'e3w_1d'  , e3w_1d , ktype = jp_r8 )
677      !
678      CALL iom_rstput( 0, 0, inum, 'e3t_0'   , e3t_0  , ktype = jp_r8 )   ! vertical scale factors
679      CALL iom_rstput( 0, 0, inum, 'e3u_0'   , e3u_0  , ktype = jp_r8 )
680      CALL iom_rstput( 0, 0, inum, 'e3v_0'   , e3v_0  , ktype = jp_r8 )
681      CALL iom_rstput( 0, 0, inum, 'e3f_0'   , e3f_0  , ktype = jp_r8 )
682      CALL iom_rstput( 0, 0, inum, 'e3w_0'   , e3w_0  , ktype = jp_r8 )
683      CALL iom_rstput( 0, 0, inum, 'e3uw_0'  , e3uw_0 , ktype = jp_r8 )
684      CALL iom_rstput( 0, 0, inum, 'e3vw_0'  , e3vw_0 , ktype = jp_r8 )
685      !                                         
686      !                             !==  wet top and bottom level  ==!   (caution: multiplied by ssmask)
687      !
688      CALL iom_rstput( 0, 0, inum, 'top_level'    , REAL( mikt, wp )*ssmask , ktype = jp_i4 )   ! nb of ocean T-points (ISF)
689      CALL iom_rstput( 0, 0, inum, 'bottom_level' , REAL( mbkt, wp )*ssmask , ktype = jp_i4 )   ! nb of ocean T-points
690      !
691      IF( ln_sco ) THEN             ! s-coordinate: store grid stiffness ratio  (Not required anyway)
692         CALL dom_stiff( z2d )
693         CALL iom_rstput( 0, 0, inum, 'stiffness', z2d )        !    ! Max. grid stiffness ratio
694      ENDIF
695      !
696      IF( ll_wd ) THEN              ! wetting and drying domain
697         CALL iom_rstput( 0, 0, inum, 'ht_0'   , ht_0   , ktype = jp_r8 )
698      ENDIF
699      !
700      ! Add some global attributes ( netcdf only )
701      CALL iom_putatt( inum, 'nn_cfg', nn_cfg )
702      CALL iom_putatt( inum, 'cn_cfg', TRIM(cn_cfg) )
703      !
704      !                                ! ============================
705      !                                !        close the files
706      !                                ! ============================
707      CALL iom_close( inum )
708      !
709   END SUBROUTINE cfg_write
710
711   !!======================================================================
712END MODULE domain
Note: See TracBrowser for help on using the repository browser.