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

Last change on this file since 12822 was 12822, checked in by gm, 2 years ago

symmetric sterss tensor and half cell modifications (wet point only, ghost cells)

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