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/2019/fix_vvl_ticket1791/src/OCE/DOM – NEMO

source: NEMO/branches/2019/fix_vvl_ticket1791/src/OCE/DOM/domain.F90 @ 11837

Last change on this file since 11837 was 11002, checked in by jchanut, 5 years ago

vvl interp, #1791: 1) Add resting depth at F-points. 2) Rewrite dom_vvl_interpol with 2 aditionnal methods for horizontal interpolation at U-V-F points (nmet=0,1,2 in subroutine header). Set nmet=1, i.e. interpolate interfaces instead of scale factors anomaly. 3) Changes in the rest of the code: consider input array in dom_vvl_interpol at T point for F scale factors interpolation.

  • Property svn:keywords set to Id
File size: 37.2 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 closea         ! closed seas
33   USE domhgr         ! domain: set the horizontal mesh
34   USE domzgr         ! domain: set the vertical mesh
35   USE dommsk         ! domain: set the mask system
36   USE domwri         ! domain: write the meshmask file
37   USE domvvl         ! variable volume
38   USE c1d            ! 1D configuration
39   USE dyncor_c1d     ! 1D configuration: Coriolis term    (cor_c1d routine)
40   USE wet_dry,  ONLY : ll_wd
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(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 ::   ji, jj, jk, ik   ! dummy loop indices
76      INTEGER ::   iconf = 0    ! local integers
77      CHARACTER (len=64) ::   cform = "(A12, 3(A13, I7))" 
78      CHARACTER (len=*), INTENT(IN) :: cdstr                  ! model: NEMO or SAS. Determines core restart variables
79      INTEGER , DIMENSION(jpi,jpj) ::   ik_top , ik_bot       ! top and bottom ocean level
80      REAL(wp), DIMENSION(jpi,jpj) ::   z1_hu_0, z1_hv_0
81      !!----------------------------------------------------------------------
82      !
83      IF(lwp) THEN         ! Ocean domain Parameters (control print)
84         WRITE(numout,*)
85         WRITE(numout,*) 'dom_init : domain initialization'
86         WRITE(numout,*) '~~~~~~~~'
87         !
88         WRITE(numout,*)     '   Domain info'
89         WRITE(numout,*)     '      dimension of model:'
90         WRITE(numout,*)     '             Local domain      Global domain       Data domain '
91         WRITE(numout,cform) '        ','   jpi     : ', jpi, '   jpiglo  : ', jpiglo
92         WRITE(numout,cform) '        ','   jpj     : ', jpj, '   jpjglo  : ', jpjglo
93         WRITE(numout,cform) '        ','   jpk     : ', jpk, '   jpkglo  : ', jpkglo
94         WRITE(numout,cform) '       ' ,'   jpij    : ', jpij
95         WRITE(numout,*)     '      mpp local domain info (mpp):'
96         WRITE(numout,*)     '              jpni    : ', jpni, '   nn_hls  : ', nn_hls
97         WRITE(numout,*)     '              jpnj    : ', jpnj, '   nn_hls  : ', nn_hls
98         WRITE(numout,*)     '              jpnij   : ', jpnij
99         WRITE(numout,*)     '      lateral boundary of the Global domain : jperio  = ', jperio
100         SELECT CASE ( jperio )
101         CASE( 0 )   ;   WRITE(numout,*) '         (i.e. closed)'
102         CASE( 1 )   ;   WRITE(numout,*) '         (i.e. cyclic east-west)'
103         CASE( 2 )   ;   WRITE(numout,*) '         (i.e. equatorial symmetric)'
104         CASE( 3 )   ;   WRITE(numout,*) '         (i.e. north fold with T-point pivot)'
105         CASE( 4 )   ;   WRITE(numout,*) '         (i.e. cyclic east-west and north fold with T-point pivot)'
106         CASE( 5 )   ;   WRITE(numout,*) '         (i.e. north fold with F-point pivot)'
107         CASE( 6 )   ;   WRITE(numout,*) '         (i.e. cyclic east-west and north fold with F-point pivot)'
108         CASE( 7 )   ;   WRITE(numout,*) '         (i.e. cyclic east-west and north-south)'
109         CASE DEFAULT
110            CALL ctl_stop( 'jperio is out of range' )
111         END SELECT
112         WRITE(numout,*)     '      Ocean model configuration used:'
113         WRITE(numout,*)     '         cn_cfg = ', TRIM( cn_cfg ), '   nn_cfg = ', nn_cfg
114      ENDIF
115      lwxios = .FALSE.
116      ln_xios_read = .FALSE.
117      !
118      !           !==  Reference coordinate system  ==!
119      !
120      CALL dom_glo                     ! global domain versus local domain
121      CALL dom_nam                     ! read namelist ( namrun, namdom )
122      !
123      IF( lwxios ) THEN
124!define names for restart write and set core output (restart.F90)
125         CALL iom_set_rst_vars(rst_wfields)
126         CALL iom_set_rstw_core(cdstr)
127      ENDIF
128!reset namelist for SAS
129      IF(cdstr == 'SAS') THEN
130         IF(lrxios) THEN
131               IF(lwp) write(numout,*) 'Disable reading restart file using XIOS for SAS'
132               lrxios = .FALSE.
133         ENDIF
134      ENDIF
135      !
136      CALL dom_hgr                     ! Horizontal mesh
137      CALL dom_zgr( ik_top, ik_bot )   ! Vertical mesh and bathymetry
138      CALL dom_msk( ik_top, ik_bot )   ! Masks
139      IF( ln_closea )   CALL dom_clo   ! ln_closea=T : closed seas included in the simulation
140                                       ! Read in masks to define closed seas and lakes
141      !
142      DO jj = 1, jpj                   ! depth of the iceshelves
143         DO ji = 1, jpi
144            ik = mikt(ji,jj)
145            risfdep(ji,jj) = gdepw_0(ji,jj,ik)
146         END DO
147      END DO
148      !
149      ht_0(:,:) = 0._wp  ! Reference ocean thickness
150      hu_0(:,:) = 0._wp
151      hv_0(:,:) = 0._wp
152      hf_0(:,:) = 0._wp
153      DO jk = 1, jpk
154         ht_0(:,:) = ht_0(:,:) + e3t_0(:,:,jk) * tmask(:,:,jk)
155         hu_0(:,:) = hu_0(:,:) + e3u_0(:,:,jk) * umask(:,:,jk)
156         hv_0(:,:) = hv_0(:,:) + e3v_0(:,:,jk) * vmask(:,:,jk)
157      END DO
158      !
159      DO jj = 1, jpjm1
160         DO jk = 1, jpk
161            hf_0(:,jj) = hf_0(:,jj) + e3f_0(:,jj,jk) * umask(:,jj,jk) * umask(:,jj+1,jk)
162         END DO
163      END DO
164      CALL lbc_lnk( 'dom_init' , hf_0, 'F', 1._wp )
165      !           !==  time varying part of coordinate system  ==!
166      !
167      IF( ln_linssh ) THEN       != Fix in time : set to the reference one for all
168      !
169         !       before        !          now          !       after         !
170            gdept_b = gdept_0  ;   gdept_n = gdept_0   !        ---          ! depth of grid-points
171            gdepw_b = gdepw_0  ;   gdepw_n = gdepw_0   !        ---          !
172                                   gde3w_n = gde3w_0   !        ---          !
173         !                                                                 
174              e3t_b =   e3t_0  ;     e3t_n =   e3t_0   ;   e3t_a =  e3t_0    ! scale factors
175              e3u_b =   e3u_0  ;     e3u_n =   e3u_0   ;   e3u_a =  e3u_0    !
176              e3v_b =   e3v_0  ;     e3v_n =   e3v_0   ;   e3v_a =  e3v_0    !
177                                     e3f_n =   e3f_0   !        ---          !
178              e3w_b =   e3w_0  ;     e3w_n =   e3w_0   !        ---          !
179             e3uw_b =  e3uw_0  ;    e3uw_n =  e3uw_0   !        ---          !
180             e3vw_b =  e3vw_0  ;    e3vw_n =  e3vw_0   !        ---          !
181         !
182         z1_hu_0(:,:) = ssumask(:,:) / ( hu_0(:,:) + 1._wp - ssumask(:,:) )     ! _i mask due to ISF
183         z1_hv_0(:,:) = ssvmask(:,:) / ( hv_0(:,:) + 1._wp - ssvmask(:,:) )
184         !
185         !        before       !          now          !       after         !
186                                      ht_n =    ht_0   !                     ! water column thickness
187               hu_b =    hu_0  ;      hu_n =    hu_0   ;    hu_a =    hu_0   !
188               hv_b =    hv_0  ;      hv_n =    hv_0   ;    hv_a =    hv_0   !
189            r1_hu_b = z1_hu_0  ;   r1_hu_n = z1_hu_0   ; r1_hu_a = z1_hu_0   ! inverse of water column thickness
190            r1_hv_b = z1_hv_0  ;   r1_hv_n = z1_hv_0   ; r1_hv_a = z1_hv_0   !
191         !
192         !
193      ELSE                       != time varying : initialize before/now/after variables
194         !
195         IF( .NOT.l_offline )  CALL dom_vvl_init 
196         !
197      ENDIF
198      !
199      IF( lk_c1d         )   CALL cor_c1d       ! 1D configuration: Coriolis set at T-point
200      !
201      IF( ln_meshmask .AND. .NOT.ln_iscpl )                        CALL dom_wri     ! Create a domain file
202      IF( ln_meshmask .AND.      ln_iscpl .AND. .NOT.ln_rstart )   CALL dom_wri     ! Create a domain file
203      IF(                                       .NOT.ln_rstart )   CALL dom_ctl     ! Domain control
204      !
205      IF( ln_write_cfg )   CALL cfg_write         ! create the configuration file
206      !
207      IF(lwp) THEN
208         WRITE(numout,*)
209         WRITE(numout,*) 'dom_init :   ==>>>   END of domain initialization'
210         WRITE(numout,*) '~~~~~~~~'
211         WRITE(numout,*) 
212      ENDIF
213      !
214   END SUBROUTINE dom_init
215
216
217   SUBROUTINE dom_glo
218      !!----------------------------------------------------------------------
219      !!                     ***  ROUTINE dom_glo  ***
220      !!
221      !! ** Purpose :   initialization of global domain <--> local domain indices
222      !!
223      !! ** Method  :   
224      !!
225      !! ** Action  : - mig , mjg : local  domain indices ==> global domain indices
226      !!              - mi0 , mi1 : global domain indices ==> local  domain indices
227      !!              - mj0,, mj1   (global point not in the local domain ==> mi0>mi1 and/or mj0>mj1)
228      !!----------------------------------------------------------------------
229      INTEGER ::   ji, jj   ! dummy loop argument
230      !!----------------------------------------------------------------------
231      !
232      DO ji = 1, jpi                 ! local domain indices ==> global domain indices
233        mig(ji) = ji + nimpp - 1
234      END DO
235      DO jj = 1, jpj
236        mjg(jj) = jj + njmpp - 1
237      END DO
238      !                              ! global domain indices ==> local domain indices
239      !                                   ! (return (m.0,m.1)=(1,0) if data domain gridpoint is to the west/south of the
240      !                                   ! local domain, or (m.0,m.1)=(jp.+1,jp.) to the east/north of local domain.
241      DO ji = 1, jpiglo
242        mi0(ji) = MAX( 1 , MIN( ji - nimpp + 1, jpi+1 ) )
243        mi1(ji) = MAX( 0 , MIN( ji - nimpp + 1, jpi   ) )
244      END DO
245      DO jj = 1, jpjglo
246        mj0(jj) = MAX( 1 , MIN( jj - njmpp + 1, jpj+1 ) )
247        mj1(jj) = MAX( 0 , MIN( jj - njmpp + 1, jpj   ) )
248      END DO
249      IF(lwp) THEN                   ! control print
250         WRITE(numout,*)
251         WRITE(numout,*) 'dom_glo : domain: global <<==>> local '
252         WRITE(numout,*) '~~~~~~~ '
253         WRITE(numout,*) '   global domain:   jpiglo = ', jpiglo, ' jpjglo = ', jpjglo, ' jpkglo = ', jpkglo
254         WRITE(numout,*) '   local  domain:   jpi    = ', jpi   , ' jpj    = ', jpj   , ' jpk    = ', jpk
255         WRITE(numout,*)
256         WRITE(numout,*) '   conversion from local to global domain indices (and vise versa) done'
257         IF( nn_print >= 1 ) THEN
258            WRITE(numout,*)
259            WRITE(numout,*) '          conversion local  ==> global i-index domain (mig)'
260            WRITE(numout,25)              (mig(ji),ji = 1,jpi)
261            WRITE(numout,*)
262            WRITE(numout,*) '          conversion global ==> local  i-index domain'
263            WRITE(numout,*) '             starting index (mi0)'
264            WRITE(numout,25)              (mi0(ji),ji = 1,jpiglo)
265            WRITE(numout,*) '             ending index (mi1)'
266            WRITE(numout,25)              (mi1(ji),ji = 1,jpiglo)
267            WRITE(numout,*)
268            WRITE(numout,*) '          conversion local  ==> global j-index domain (mjg)'
269            WRITE(numout,25)              (mjg(jj),jj = 1,jpj)
270            WRITE(numout,*)
271            WRITE(numout,*) '          conversion global ==> local  j-index domain'
272            WRITE(numout,*) '             starting index (mj0)'
273            WRITE(numout,25)              (mj0(jj),jj = 1,jpjglo)
274            WRITE(numout,*) '             ending index (mj1)'
275            WRITE(numout,25)              (mj1(jj),jj = 1,jpjglo)
276         ENDIF
277      ENDIF
278 25   FORMAT( 100(10x,19i4,/) )
279      !
280   END SUBROUTINE dom_glo
281
282
283   SUBROUTINE dom_nam
284      !!----------------------------------------------------------------------
285      !!                     ***  ROUTINE dom_nam  ***
286      !!                   
287      !! ** Purpose :   read domaine namelists and print the variables.
288      !!
289      !! ** input   : - namrun namelist
290      !!              - namdom namelist
291      !!              - namnc4 namelist   ! "key_netcdf4" only
292      !!----------------------------------------------------------------------
293      USE ioipsl
294      !!
295      INTEGER  ::   ios   ! Local integer
296      !
297      NAMELIST/namrun/ cn_ocerst_indir, cn_ocerst_outdir, nn_stocklist, ln_rst_list,                 &
298         &             nn_no   , cn_exp   , cn_ocerst_in, cn_ocerst_out, ln_rstart , nn_rstctl ,     &
299         &             nn_it000, nn_itend , nn_date0    , nn_time0     , nn_leapy  , nn_istate ,     &
300         &             nn_stock, nn_write , ln_mskland  , ln_clobber   , nn_chunksz, nn_euler  ,     &
301         &             ln_cfmeta, ln_iscpl, ln_xios_read, nn_wxios
302      NAMELIST/namdom/ ln_linssh, rn_isfhmin, rn_rdt, rn_atfp, ln_crs, ln_meshmask
303#if defined key_netcdf4
304      NAMELIST/namnc4/ nn_nchunks_i, nn_nchunks_j, nn_nchunks_k, ln_nc4zip
305#endif
306      !!----------------------------------------------------------------------
307      !
308      IF(lwp) THEN
309         WRITE(numout,*)
310         WRITE(numout,*) 'dom_nam : domain initialization through namelist read'
311         WRITE(numout,*) '~~~~~~~ '
312      ENDIF
313      !
314      !
315      REWIND( numnam_ref )              ! Namelist namrun in reference namelist : Parameters of the run
316      READ  ( numnam_ref, namrun, IOSTAT = ios, ERR = 901)
317901   IF( ios /= 0 )   CALL ctl_nam ( ios , 'namrun in reference namelist', lwp )
318      REWIND( numnam_cfg )              ! Namelist namrun in configuration namelist : Parameters of the run
319      READ  ( numnam_cfg, namrun, IOSTAT = ios, ERR = 902 )
320902   IF( ios >  0 )   CALL ctl_nam ( ios , 'namrun in configuration namelist', lwp )
321      IF(lwm) WRITE ( numond, namrun )
322      !
323      IF(lwp) THEN                  ! control print
324         WRITE(numout,*) '   Namelist : namrun   ---   run parameters'
325         WRITE(numout,*) '      Assimilation cycle              nn_no           = ', nn_no
326         WRITE(numout,*) '      experiment name for output      cn_exp          = ', TRIM( cn_exp           )
327         WRITE(numout,*) '      file prefix restart input       cn_ocerst_in    = ', TRIM( cn_ocerst_in     )
328         WRITE(numout,*) '      restart input directory         cn_ocerst_indir = ', TRIM( cn_ocerst_indir  )
329         WRITE(numout,*) '      file prefix restart output      cn_ocerst_out   = ', TRIM( cn_ocerst_out    )
330         WRITE(numout,*) '      restart output directory        cn_ocerst_outdir= ', TRIM( cn_ocerst_outdir )
331         WRITE(numout,*) '      restart logical                 ln_rstart       = ', ln_rstart
332         WRITE(numout,*) '      start with forward time step    nn_euler        = ', nn_euler
333         WRITE(numout,*) '      control of time step            nn_rstctl       = ', nn_rstctl
334         WRITE(numout,*) '      number of the first time step   nn_it000        = ', nn_it000
335         WRITE(numout,*) '      number of the last time step    nn_itend        = ', nn_itend
336         WRITE(numout,*) '      initial calendar date aammjj    nn_date0        = ', nn_date0
337         WRITE(numout,*) '      initial time of day in hhmm     nn_time0        = ', nn_time0
338         WRITE(numout,*) '      leap year calendar (0/1)        nn_leapy        = ', nn_leapy
339         WRITE(numout,*) '      initial state output            nn_istate       = ', nn_istate
340         IF( ln_rst_list ) THEN
341            WRITE(numout,*) '      list of restart dump times      nn_stocklist    =', nn_stocklist
342         ELSE
343            WRITE(numout,*) '      frequency of restart file       nn_stock        = ', nn_stock
344         ENDIF
345         WRITE(numout,*) '      frequency of output file        nn_write        = ', nn_write
346         WRITE(numout,*) '      mask land points                ln_mskland      = ', ln_mskland
347         WRITE(numout,*) '      additional CF standard metadata ln_cfmeta       = ', ln_cfmeta
348         WRITE(numout,*) '      overwrite an existing file      ln_clobber      = ', ln_clobber
349         WRITE(numout,*) '      NetCDF chunksize (bytes)        nn_chunksz      = ', nn_chunksz
350         WRITE(numout,*) '      IS coupling at the restart step ln_iscpl        = ', ln_iscpl
351         IF( TRIM(Agrif_CFixed()) == '0' ) THEN
352            WRITE(numout,*) '      READ restart for a single file using XIOS ln_xios_read =', ln_xios_read
353            WRITE(numout,*) '      Write restart using XIOS        nn_wxios   = ', nn_wxios
354         ELSE
355            WRITE(numout,*) "      AGRIF: nn_wxios will be ingored. See setting for parent"
356            WRITE(numout,*) "      AGRIF: ln_xios_read will be ingored. See setting for parent"
357         ENDIF
358      ENDIF
359
360      cexper = cn_exp         ! conversion DOCTOR names into model names (this should disappear soon)
361      nrstdt = nn_rstctl
362      nit000 = nn_it000
363      nitend = nn_itend
364      ndate0 = nn_date0
365      nleapy = nn_leapy
366      ninist = nn_istate
367      nstock = nn_stock
368      nstocklist = nn_stocklist
369      nwrite = nn_write
370      neuler = nn_euler
371      IF( neuler == 1 .AND. .NOT. ln_rstart ) THEN
372         IF(lwp) WRITE(numout,*) 
373         IF(lwp) WRITE(numout,*)'   ==>>>   Start from rest (ln_rstart=F)'
374         IF(lwp) WRITE(numout,*)'           an Euler initial time step is used : nn_euler is forced to 0 '   
375         neuler = 0
376      ENDIF
377      !                             ! control of output frequency
378      IF( nstock == 0 .OR. nstock > nitend ) THEN
379         WRITE(ctmp1,*) 'nstock = ', nstock, ' it is forced to ', nitend
380         CALL ctl_warn( ctmp1 )
381         nstock = nitend
382      ENDIF
383      IF ( nwrite == 0 ) THEN
384         WRITE(ctmp1,*) 'nwrite = ', nwrite, ' it is forced to ', nitend
385         CALL ctl_warn( ctmp1 )
386         nwrite = nitend
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      REWIND( numnam_ref )              ! Namelist namdom in reference namelist : space & time domain (bathymetry, mesh, timestep)
409      READ  ( numnam_ref, namdom, IOSTAT = ios, ERR = 903)
410903   IF( ios /= 0 )   CALL ctl_nam ( ios , 'namdom in reference namelist', lwp )
411      REWIND( numnam_cfg )              ! Namelist namdom in configuration namelist : space & time domain (bathymetry, mesh, timestep)
412      READ  ( numnam_cfg, namdom, IOSTAT = ios, ERR = 904 )
413904   IF( ios >  0 )   CALL ctl_nam ( ios , 'namdom in configuration namelist', lwp )
414      IF(lwm) WRITE( numond, namdom )
415      !
416      IF(lwp) THEN
417         WRITE(numout,*)
418         WRITE(numout,*) '   Namelist : namdom   ---   space & time domain'
419         WRITE(numout,*) '      linear free surface (=T)                ln_linssh   = ', ln_linssh
420         WRITE(numout,*) '      create mesh/mask file                   ln_meshmask = ', ln_meshmask
421         WRITE(numout,*) '      treshold to open the isf cavity         rn_isfhmin  = ', rn_isfhmin, ' [m]'
422         WRITE(numout,*) '      ocean time step                         rn_rdt      = ', rn_rdt
423         WRITE(numout,*) '      asselin time filter parameter           rn_atfp     = ', rn_atfp
424         WRITE(numout,*) '      online coarsening of dynamical fields   ln_crs      = ', ln_crs
425      ENDIF
426      !
427      !          ! conversion DOCTOR names into model names (this should disappear soon)
428      atfp = rn_atfp
429      rdt  = rn_rdt
430
431      IF( TRIM(Agrif_CFixed()) == '0' ) THEN
432         lrxios = ln_xios_read.AND.ln_rstart
433!set output file type for XIOS based on NEMO namelist
434         IF (nn_wxios > 0) lwxios = .TRUE. 
435         nxioso = nn_wxios
436      ENDIF
437
438#if defined key_netcdf4
439      !                             ! NetCDF 4 case   ("key_netcdf4" defined)
440      REWIND( numnam_ref )              ! Namelist namnc4 in reference namelist : NETCDF
441      READ  ( numnam_ref, namnc4, IOSTAT = ios, ERR = 907)
442907   IF( ios /= 0 )   CALL ctl_nam ( ios , 'namnc4 in reference namelist', lwp )
443      REWIND( numnam_cfg )              ! Namelist namnc4 in configuration namelist : NETCDF
444      READ  ( numnam_cfg, namnc4, IOSTAT = ios, ERR = 908 )
445908   IF( ios >  0 )   CALL ctl_nam ( ios , 'namnc4 in configuration namelist', lwp )
446      IF(lwm) WRITE( numond, namnc4 )
447
448      IF(lwp) THEN                        ! control print
449         WRITE(numout,*)
450         WRITE(numout,*) '   Namelist namnc4 - Netcdf4 chunking parameters'
451         WRITE(numout,*) '      number of chunks in i-dimension             nn_nchunks_i = ', nn_nchunks_i
452         WRITE(numout,*) '      number of chunks in j-dimension             nn_nchunks_j = ', nn_nchunks_j
453         WRITE(numout,*) '      number of chunks in k-dimension             nn_nchunks_k = ', nn_nchunks_k
454         WRITE(numout,*) '      apply netcdf4/hdf5 chunking & compression   ln_nc4zip    = ', ln_nc4zip
455      ENDIF
456
457      ! Put the netcdf4 settings into a simple structure (snc4set, defined in in_out_manager module)
458      ! Note the chunk size in the unlimited (time) dimension will be fixed at 1
459      snc4set%ni   = nn_nchunks_i
460      snc4set%nj   = nn_nchunks_j
461      snc4set%nk   = nn_nchunks_k
462      snc4set%luse = ln_nc4zip
463#else
464      snc4set%luse = .FALSE.        ! No NetCDF 4 case
465#endif
466      !
467   END SUBROUTINE dom_nam
468
469
470   SUBROUTINE dom_ctl
471      !!----------------------------------------------------------------------
472      !!                     ***  ROUTINE dom_ctl  ***
473      !!
474      !! ** Purpose :   Domain control.
475      !!
476      !! ** Method  :   compute and print extrema of masked scale factors
477      !!----------------------------------------------------------------------
478      INTEGER, DIMENSION(2) ::   imi1, imi2, ima1, ima2
479      INTEGER, DIMENSION(2) ::   iloc   !
480      REAL(wp) ::   ze1min, ze1max, ze2min, ze2max
481      !!----------------------------------------------------------------------
482      !
483      IF(lk_mpp) THEN
484         CALL mpp_minloc( 'domain', e1t(:,:), tmask_i(:,:), ze1min, imi1 )
485         CALL mpp_minloc( 'domain', e2t(:,:), tmask_i(:,:), ze2min, imi2 )
486         CALL mpp_maxloc( 'domain', e1t(:,:), tmask_i(:,:), ze1max, ima1 )
487         CALL mpp_maxloc( 'domain', e2t(:,:), tmask_i(:,:), ze2max, ima2 )
488      ELSE
489         ze1min = MINVAL( e1t(:,:), mask = tmask_i(:,:) == 1._wp )   
490         ze2min = MINVAL( e2t(:,:), mask = tmask_i(:,:) == 1._wp )   
491         ze1max = MAXVAL( e1t(:,:), mask = tmask_i(:,:) == 1._wp )   
492         ze2max = MAXVAL( e2t(:,:), mask = tmask_i(:,:) == 1._wp )   
493         !
494         iloc  = MINLOC( e1t(:,:), mask = tmask_i(:,:) == 1._wp )
495         imi1(1) = iloc(1) + nimpp - 1
496         imi1(2) = iloc(2) + njmpp - 1
497         iloc  = MINLOC( e2t(:,:), mask = tmask_i(:,:) == 1._wp )
498         imi2(1) = iloc(1) + nimpp - 1
499         imi2(2) = iloc(2) + njmpp - 1
500         iloc  = MAXLOC( e1t(:,:), mask = tmask_i(:,:) == 1._wp )
501         ima1(1) = iloc(1) + nimpp - 1
502         ima1(2) = iloc(2) + njmpp - 1
503         iloc  = MAXLOC( e2t(:,:), mask = tmask_i(:,:) == 1._wp )
504         ima2(1) = iloc(1) + nimpp - 1
505         ima2(2) = iloc(2) + njmpp - 1
506      ENDIF
507      IF(lwp) THEN
508         WRITE(numout,*)
509         WRITE(numout,*) 'dom_ctl : extrema of the masked scale factors'
510         WRITE(numout,*) '~~~~~~~'
511         WRITE(numout,"(14x,'e1t maxi: ',1f10.2,' at i = ',i5,' j= ',i5)") ze1max, ima1(1), ima1(2)
512         WRITE(numout,"(14x,'e1t mini: ',1f10.2,' at i = ',i5,' j= ',i5)") ze1min, imi1(1), imi1(2)
513         WRITE(numout,"(14x,'e2t maxi: ',1f10.2,' at i = ',i5,' j= ',i5)") ze2max, ima2(1), ima2(2)
514         WRITE(numout,"(14x,'e2t mini: ',1f10.2,' at i = ',i5,' j= ',i5)") ze2min, imi2(1), imi2(2)
515      ENDIF
516      !
517   END SUBROUTINE dom_ctl
518
519
520   SUBROUTINE domain_cfg( ldtxt, cd_cfg, kk_cfg, kpi, kpj, kpk, kperio )
521      !!----------------------------------------------------------------------
522      !!                     ***  ROUTINE dom_nam  ***
523      !!                   
524      !! ** Purpose :   read the domain size in domain configuration file
525      !!
526      !! ** Method  :   read the cn_domcfg NetCDF file
527      !!----------------------------------------------------------------------
528      CHARACTER(len=*), DIMENSION(:), INTENT(out) ::   ldtxt           ! stored print information
529      CHARACTER(len=*)              , INTENT(out) ::   cd_cfg          ! configuration name
530      INTEGER                       , INTENT(out) ::   kk_cfg          ! configuration resolution
531      INTEGER                       , INTENT(out) ::   kpi, kpj, kpk   ! global domain sizes
532      INTEGER                       , INTENT(out) ::   kperio          ! lateral global domain b.c.
533      !
534      INTEGER ::   inum, ii   ! local integer
535      REAL(wp) ::   zorca_res                     ! local scalars
536      REAL(wp) ::   ziglo, zjglo, zkglo, zperio   !   -      -
537      !!----------------------------------------------------------------------
538      !
539      ii = 1
540      WRITE(ldtxt(ii),*) '           '                                                    ;   ii = ii+1
541      WRITE(ldtxt(ii),*) 'domain_cfg : domain size read in ', TRIM( cn_domcfg ), ' file'  ;   ii = ii+1
542      WRITE(ldtxt(ii),*) '~~~~~~~~~~ '                                                    ;   ii = ii+1
543      !
544      CALL iom_open( cn_domcfg, inum )
545      !
546      !                                   !- ORCA family specificity
547      IF(  iom_varid( inum, 'ORCA'       , ldstop = .FALSE. ) > 0  .AND.  &
548         & iom_varid( inum, 'ORCA_index' , ldstop = .FALSE. ) > 0    ) THEN
549         !
550         cd_cfg = 'ORCA'
551         CALL iom_get( inum, 'ORCA_index', zorca_res )   ;   kk_cfg = NINT( zorca_res )
552         !
553         WRITE(ldtxt(ii),*) '   .'                                                     ;   ii = ii+1
554         WRITE(ldtxt(ii),*) '   ==>>>   ORCA configuration '                           ;   ii = ii+1
555         WRITE(ldtxt(ii),*) '   .'                                                     ;   ii = ii+1
556         !
557      ELSE                                !- cd_cfg & k_cfg are not used
558         cd_cfg = 'UNKNOWN'
559         kk_cfg = -9999999
560                                          !- or they may be present as global attributes
561                                          !- (netcdf only) 
562         CALL iom_getatt( inum, 'cn_cfg', cd_cfg )  ! returns   !  if not found
563         CALL iom_getatt( inum, 'nn_cfg', kk_cfg )  ! returns -999 if not found
564         IF( TRIM(cd_cfg) == '!') cd_cfg = 'UNKNOWN'
565         IF( kk_cfg == -999     ) kk_cfg = -9999999
566         !
567      ENDIF
568      !
569      CALL iom_get( inum, 'jpiglo', ziglo  )   ;   kpi = NINT( ziglo )
570      CALL iom_get( inum, 'jpjglo', zjglo  )   ;   kpj = NINT( zjglo )
571      CALL iom_get( inum, 'jpkglo', zkglo  )   ;   kpk = NINT( zkglo )
572      CALL iom_get( inum, 'jperio', zperio )   ;   kperio = NINT( zperio )
573      CALL iom_close( inum )
574      !
575      WRITE(ldtxt(ii),*) '      cn_cfg = ', TRIM(cd_cfg), '   nn_cfg = ', kk_cfg             ;   ii = ii+1
576      WRITE(ldtxt(ii),*) '      jpiglo = ', kpi                                              ;   ii = ii+1
577      WRITE(ldtxt(ii),*) '      jpjglo = ', kpj                                              ;   ii = ii+1
578      WRITE(ldtxt(ii),*) '      jpkglo = ', kpk                                              ;   ii = ii+1
579      WRITE(ldtxt(ii),*) '      type of global domain lateral boundary   jperio = ', kperio  ;   ii = ii+1
580      !       
581   END SUBROUTINE domain_cfg
582   
583   
584   SUBROUTINE cfg_write
585      !!----------------------------------------------------------------------
586      !!                  ***  ROUTINE cfg_write  ***
587      !!                   
588      !! ** Purpose :   Create the "cn_domcfg_out" file, a NetCDF file which
589      !!              contains all the ocean domain informations required to
590      !!              define an ocean configuration.
591      !!
592      !! ** Method  :   Write in a file all the arrays required to set up an
593      !!              ocean configuration.
594      !!
595      !! ** output file :   domcfg_out.nc : domain size, characteristics, horizontal
596      !!                       mesh, Coriolis parameter, and vertical scale factors
597      !!                    NB: also contain ORCA family information
598      !!----------------------------------------------------------------------
599      INTEGER           ::   ji, jj, jk   ! dummy loop indices
600      INTEGER           ::   izco, izps, isco, icav
601      INTEGER           ::   inum     ! local units
602      CHARACTER(len=21) ::   clnam    ! filename (mesh and mask informations)
603      REAL(wp), DIMENSION(jpi,jpj) ::   z2d   ! workspace
604      !!----------------------------------------------------------------------
605      !
606      IF(lwp) WRITE(numout,*)
607      IF(lwp) WRITE(numout,*) 'cfg_write : create the domain configuration file (', TRIM(cn_domcfg_out),'.nc)'
608      IF(lwp) WRITE(numout,*) '~~~~~~~~~'
609      !
610      !                       ! ============================= !
611      !                       !  create 'domcfg_out.nc' file  !
612      !                       ! ============================= !
613      !         
614      clnam = cn_domcfg_out  ! filename (configuration information)
615      CALL iom_open( TRIM(clnam), inum, ldwrt = .TRUE. )
616     
617      !
618      !                             !==  ORCA family specificities  ==!
619      IF( cn_cfg == "ORCA" ) THEN
620         CALL iom_rstput( 0, 0, inum, 'ORCA'      , 1._wp            , ktype = jp_i4 )
621         CALL iom_rstput( 0, 0, inum, 'ORCA_index', REAL( nn_cfg, wp), ktype = jp_i4 )         
622      ENDIF
623      !
624      !                             !==  global domain size  ==!
625      !
626      CALL iom_rstput( 0, 0, inum, 'jpiglo', REAL( jpiglo, wp), ktype = jp_i4 )
627      CALL iom_rstput( 0, 0, inum, 'jpjglo', REAL( jpjglo, wp), ktype = jp_i4 )
628      CALL iom_rstput( 0, 0, inum, 'jpkglo', REAL( jpk   , wp), ktype = jp_i4 )
629      !
630      !                             !==  domain characteristics  ==!
631      !
632      !                                   ! lateral boundary of the global domain
633      CALL iom_rstput( 0, 0, inum, 'jperio', REAL( jperio, wp), ktype = jp_i4 )
634      !
635      !                                   ! type of vertical coordinate
636      IF( ln_zco    ) THEN   ;   izco = 1   ;   ELSE   ;   izco = 0   ;   ENDIF
637      IF( ln_zps    ) THEN   ;   izps = 1   ;   ELSE   ;   izps = 0   ;   ENDIF
638      IF( ln_sco    ) THEN   ;   isco = 1   ;   ELSE   ;   isco = 0   ;   ENDIF
639      CALL iom_rstput( 0, 0, inum, 'ln_zco'   , REAL( izco, wp), ktype = jp_i4 )
640      CALL iom_rstput( 0, 0, inum, 'ln_zps'   , REAL( izps, wp), ktype = jp_i4 )
641      CALL iom_rstput( 0, 0, inum, 'ln_sco'   , REAL( isco, wp), ktype = jp_i4 )
642      !
643      !                                   ! ocean cavities under iceshelves
644      IF( ln_isfcav ) THEN   ;   icav = 1   ;   ELSE   ;   icav = 0   ;   ENDIF
645      CALL iom_rstput( 0, 0, inum, 'ln_isfcav', REAL( icav, wp), ktype = jp_i4 )
646      !
647      !                             !==  horizontal mesh  !
648      !
649      CALL iom_rstput( 0, 0, inum, 'glamt', glamt, ktype = jp_r8 )   ! latitude
650      CALL iom_rstput( 0, 0, inum, 'glamu', glamu, ktype = jp_r8 )
651      CALL iom_rstput( 0, 0, inum, 'glamv', glamv, ktype = jp_r8 )
652      CALL iom_rstput( 0, 0, inum, 'glamf', glamf, ktype = jp_r8 )
653      !                               
654      CALL iom_rstput( 0, 0, inum, 'gphit', gphit, ktype = jp_r8 )   ! longitude
655      CALL iom_rstput( 0, 0, inum, 'gphiu', gphiu, ktype = jp_r8 )
656      CALL iom_rstput( 0, 0, inum, 'gphiv', gphiv, ktype = jp_r8 )
657      CALL iom_rstput( 0, 0, inum, 'gphif', gphif, ktype = jp_r8 )
658      !                               
659      CALL iom_rstput( 0, 0, inum, 'e1t'  , e1t  , ktype = jp_r8 )   ! i-scale factors (e1.)
660      CALL iom_rstput( 0, 0, inum, 'e1u'  , e1u  , ktype = jp_r8 )
661      CALL iom_rstput( 0, 0, inum, 'e1v'  , e1v  , ktype = jp_r8 )
662      CALL iom_rstput( 0, 0, inum, 'e1f'  , e1f  , ktype = jp_r8 )
663      !
664      CALL iom_rstput( 0, 0, inum, 'e2t'  , e2t  , ktype = jp_r8 )   ! j-scale factors (e2.)
665      CALL iom_rstput( 0, 0, inum, 'e2u'  , e2u  , ktype = jp_r8 )
666      CALL iom_rstput( 0, 0, inum, 'e2v'  , e2v  , ktype = jp_r8 )
667      CALL iom_rstput( 0, 0, inum, 'e2f'  , e2f  , ktype = jp_r8 )
668      !
669      CALL iom_rstput( 0, 0, inum, 'ff_f' , ff_f , ktype = jp_r8 )   ! coriolis factor
670      CALL iom_rstput( 0, 0, inum, 'ff_t' , ff_t , ktype = jp_r8 )
671      !
672      !                             !==  vertical mesh  ==!
673      !                                                     
674      CALL iom_rstput( 0, 0, inum, 'e3t_1d'  , e3t_1d , ktype = jp_r8 )   ! reference 1D-coordinate
675      CALL iom_rstput( 0, 0, inum, 'e3w_1d'  , e3w_1d , ktype = jp_r8 )
676      !
677      CALL iom_rstput( 0, 0, inum, 'e3t_0'   , e3t_0  , ktype = jp_r8 )   ! vertical scale factors
678      CALL iom_rstput( 0, 0, inum, 'e3u_0'   , e3u_0  , ktype = jp_r8 )
679      CALL iom_rstput( 0, 0, inum, 'e3v_0'   , e3v_0  , ktype = jp_r8 )
680      CALL iom_rstput( 0, 0, inum, 'e3f_0'   , e3f_0  , ktype = jp_r8 )
681      CALL iom_rstput( 0, 0, inum, 'e3w_0'   , e3w_0  , ktype = jp_r8 )
682      CALL iom_rstput( 0, 0, inum, 'e3uw_0'  , e3uw_0 , ktype = jp_r8 )
683      CALL iom_rstput( 0, 0, inum, 'e3vw_0'  , e3vw_0 , ktype = jp_r8 )
684      !                                         
685      !                             !==  wet top and bottom level  ==!   (caution: multiplied by ssmask)
686      !
687      CALL iom_rstput( 0, 0, inum, 'top_level'    , REAL( mikt, wp )*ssmask , ktype = jp_i4 )   ! nb of ocean T-points (ISF)
688      CALL iom_rstput( 0, 0, inum, 'bottom_level' , REAL( mbkt, wp )*ssmask , ktype = jp_i4 )   ! nb of ocean T-points
689      !
690      IF( ln_sco ) THEN             ! s-coordinate: store grid stiffness ratio  (Not required anyway)
691         CALL dom_stiff( z2d )
692         CALL iom_rstput( 0, 0, inum, 'stiffness', z2d )        !    ! Max. grid stiffness ratio
693      ENDIF
694      !
695      IF( ll_wd ) THEN              ! wetting and drying domain
696         CALL iom_rstput( 0, 0, inum, 'ht_0'   , ht_0   , ktype = jp_r8 )
697      ENDIF
698      !
699      ! Add some global attributes ( netcdf only )
700      CALL iom_putatt( inum, 'nn_cfg', nn_cfg )
701      CALL iom_putatt( inum, 'cn_cfg', TRIM(cn_cfg) )
702      !
703      !                                ! ============================
704      !                                !        close the files
705      !                                ! ============================
706      CALL iom_close( inum )
707      !
708   END SUBROUTINE cfg_write
709
710   !!======================================================================
711END MODULE domain
Note: See TracBrowser for help on using the repository browser.