source: branches/2016/dev_r6409_SIMPLIF_2_usrdef/NEMOGCM/NEMO/OPA_SRC/DOM/domain.F90 @ 7188

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

#1692 - branch SIMPLIF_2_usrdef: e3.=dk[dep.] (discret derivative)

  • Property svn:keywords set to Id
File size: 31.6 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   !!   cfg_write     : create the "domain_cfg.nc" file containing all required configuration information   
25   !!----------------------------------------------------------------------
26   USE oce            ! ocean variables
27   USE dom_oce        ! domain: ocean
28   USE sbc_oce        ! surface boundary condition: ocean
29   USE phycst         ! physical constants
30   USE usrdef_closea  ! closed seas
31   USE domhgr         ! domain: set the horizontal mesh
32   USE domzgr         ! domain: set the vertical mesh
33   USE dommsk         ! domain: set the mask system
34   USE domwri         ! domain: write the meshmask file
35   USE domvvl         ! variable volume
36   USE c1d            ! 1D configuration
37   USE domc1d         ! 1D configuration: column location
38   USE dyncor_c1d     ! 1D configuration: Coriolis term    (cor_c1d routine)
39   !
40   USE in_out_manager ! I/O manager
41   USE iom            ! I/O library
42   USE lbclnk         ! ocean lateral boundary condition (or mpp link)
43   USE lib_mpp        ! distributed memory computing library
44   USE wrk_nemo       ! Memory Allocation
45   USE timing         ! Timing
46
47   IMPLICIT NONE
48   PRIVATE
49
50   PUBLIC   dom_init   ! called by opa.F90
51
52   !!-------------------------------------------------------------------------
53   !! NEMO/OPA 3.3 , NEMO Consortium (2010)
54   !! $Id$
55   !! Software governed by the CeCILL licence        (NEMOGCM/NEMO_CeCILL.txt)
56   !!-------------------------------------------------------------------------
57CONTAINS
58
59   SUBROUTINE dom_init
60      !!----------------------------------------------------------------------
61      !!                  ***  ROUTINE dom_init  ***
62      !!                   
63      !! ** Purpose :   Domain initialization. Call the routines that are
64      !!              required to create the arrays which define the space
65      !!              and time domain of the ocean model.
66      !!
67      !! ** Method  : - dom_msk: compute the masks from the bathymetry file
68      !!              - dom_hgr: compute or read the horizontal grid-point position
69      !!                         and scale factors, and the coriolis factor
70      !!              - dom_zgr: define the vertical coordinate and the bathymetry
71      !!              - dom_wri: create the meshmask file if nn_msh=1
72      !!              - 1D configuration, move Coriolis, u and v at T-point
73      !!----------------------------------------------------------------------
74      INTEGER ::   ji, jj, jk, ik   ! dummy loop indices
75      INTEGER ::   iconf = 0    ! local integers
76      CHARACTER (len=64) ::   cform = "(A12, 3(A13, I7))" 
77      INTEGER , DIMENSION(jpi,jpj) ::   ik_top , ik_bot       ! top and bottom ocean level
78      REAL(wp), DIMENSION(jpi,jpj) ::   z1_hu_0, z1_hv_0
79      !!----------------------------------------------------------------------
80      !
81      IF( nn_timing == 1 )   CALL timing_start('dom_init')
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, '   jpreci  : ', jpreci
97         WRITE(numout,*)     '              jpnj    : ', jpnj, '   jprecj  : ', jprecj
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 DEFAULT
109            CALL ctl_stop( 'jperio is out of range' )
110         END SELECT
111         WRITE(numout,*)     '      Ocean model configuration used:'
112         WRITE(numout,*)     '              cp_cfg = ', cp_cfg
113         WRITE(numout,*)     '              jp_cfg = ', jp_cfg
114      ENDIF
115      !
116      !     
117!!gm  This should be removed with the new configuration interface
118      IF( lk_c1d .AND. ln_c1d_locpt )  CALL dom_c1d( rn_lat1d, rn_lon1d )
119!!gm end
120      !
121      !           !==  Reference coordinate system  ==!
122      !
123      CALL dom_glo                     ! global domain versus local domain
124      CALL dom_nam                     ! read namelist ( namrun, namdom )
125      CALL dom_clo( cp_cfg, jp_cfg )   ! Closed seas and lake
126      CALL dom_hgr                     ! Horizontal mesh
127      CALL dom_zgr( ik_top, ik_bot )   ! Vertical mesh and bathymetry
128      IF( nn_closea == 0 )   CALL clo_bat( ik_top, ik_bot )    !==  remove closed seas or lakes  ==!
129      CALL dom_msk( ik_top, ik_bot )   ! Masks
130      !
131      DO jj = 1, jpj                   ! depth of the iceshelves
132         DO ji = 1, jpi
133            ik = mikt(ji,jj)
134            risfdep(ji,jj) = gdepw_0(ji,jj,ik)
135         END DO
136      END DO
137      !
138      ht_0(:,:) = e3t_0(:,:,1) * tmask(:,:,1)   ! Reference ocean thickness
139      hu_0(:,:) = e3u_0(:,:,1) * umask(:,:,1)
140      hv_0(:,:) = e3v_0(:,:,1) * vmask(:,:,1)
141      DO jk = 2, jpk
142         ht_0(:,:) = ht_0(:,:) + e3t_0(:,:,jk) * tmask(:,:,jk)
143         hu_0(:,:) = hu_0(:,:) + e3u_0(:,:,jk) * umask(:,:,jk)
144         hv_0(:,:) = hv_0(:,:) + e3v_0(:,:,jk) * vmask(:,:,jk)
145      END DO
146      !
147      !           !==  time varying part of coordinate system  ==!
148      !
149      IF( ln_linssh ) THEN       != Fix in time : set to the reference one for all
150      !
151         !       before        !          now          !       after         !
152         ;  gdept_b = gdept_0  ;   gdept_n = gdept_0   !        ---          ! depth of grid-points
153         ;  gdepw_b = gdepw_0  ;   gdepw_n = gdepw_0   !        ---          !
154         ;                     ;   gde3w_n = gde3w_0   !        ---          !
155         !                                                                 
156         ;    e3t_b =   e3t_0  ;     e3t_n =   e3t_0   ;   e3t_a =  e3t_0    ! scale factors
157         ;    e3u_b =   e3u_0  ;     e3u_n =   e3u_0   ;   e3u_a =  e3u_0    !
158         ;    e3v_b =   e3v_0  ;     e3v_n =   e3v_0   ;   e3v_a =  e3v_0    !
159         ;                     ;     e3f_n =   e3f_0   !        ---          !
160         ;    e3w_b =   e3w_0  ;     e3w_n =   e3w_0   !        ---          !
161         ;   e3uw_b =  e3uw_0  ;    e3uw_n =  e3uw_0   !        ---          !
162         ;   e3vw_b =  e3vw_0  ;    e3vw_n =  e3vw_0   !        ---          !
163         !
164         z1_hu_0(:,:) = ssumask(:,:) / ( hu_0(:,:) + 1._wp - ssumask(:,:) )     ! _i mask due to ISF
165         z1_hv_0(:,:) = ssvmask(:,:) / ( hv_0(:,:) + 1._wp - ssvmask(:,:) )
166         !
167         !        before       !          now          !       after         !
168         ;                     ;      ht_n =    ht_0   !                     ! water column thickness
169         ;     hu_b =    hu_0  ;      hu_n =    hu_0   ;    hu_a =    hu_0   !
170         ;     hv_b =    hv_0  ;      hv_n =    hv_0   ;    hv_a =    hv_0   !
171         ;  r1_hu_b = z1_hu_0  ;   r1_hu_n = z1_hu_0   ; r1_hu_a = z1_hu_0   ! inverse of water column thickness
172         ;  r1_hv_b = z1_hv_0  ;   r1_hv_n = z1_hv_0   ; r1_hv_a = z1_hv_0   !
173         !
174         !
175      ELSE                       != time varying : initialize before/now/after variables
176         !
177         CALL dom_vvl_init 
178         !
179      ENDIF
180      !
181      IF( lk_c1d         )   CALL cor_c1d       ! 1D configuration: Coriolis set at T-point
182      !
183      IF( nn_msh > 0 .AND. .NOT. ln_iscpl )                         CALL dom_wri      ! Create a domain file
184      IF( nn_msh > 0 .AND.       ln_iscpl .AND. .NOT. ln_rstart )   CALL dom_wri      ! Create a domain file
185      IF( .NOT.ln_rstart )   CALL dom_ctl       ! Domain control
186      !
187     
188      IF(lwp) THEN
189         WRITE(numout,*)
190         WRITE(numout,*) 'dom_init : end of domain initialization nn_msh=', nn_msh
191         WRITE(numout,*) 
192      ENDIF
193      !
194      IF( ln_write_cfg )   CALL cfg_write         ! create the configuration file
195      !
196      IF( nn_timing == 1 )   CALL timing_stop('dom_init')
197      !
198   END SUBROUTINE dom_init
199
200
201   SUBROUTINE dom_glo
202      !!----------------------------------------------------------------------
203      !!                     ***  ROUTINE dom_glo  ***
204      !!
205      !! ** Purpose :   initialization of global domain <--> local domain indices
206      !!
207      !! ** Method  :   
208      !!
209      !! ** Action  : - mig , mjg : local  domain indices ==> global domain indices
210      !!              - mi0 , mi1 : global domain indices ==> local  domain indices
211      !!              - mj0,, mj1   (global point not in the local domain ==> mi0>mi1 and/or mj0>mj1)
212      !!----------------------------------------------------------------------
213      INTEGER ::   ji, jj   ! dummy loop argument
214      !!----------------------------------------------------------------------
215      !
216      DO ji = 1, jpi                 ! local domain indices ==> global domain indices
217        mig(ji) = ji + nimpp - 1
218      END DO
219      DO jj = 1, jpj
220        mjg(jj) = jj + njmpp - 1
221      END DO
222      !                              ! global domain indices ==> local domain indices
223      !                                   ! (return (m.0,m.1)=(1,0) if data domain gridpoint is to the west/south of the
224      !                                   ! local domain, or (m.0,m.1)=(jp.+1,jp.) to the east/north of local domain.
225      DO ji = 1, jpiglo
226        mi0(ji) = MAX( 1 , MIN( ji - nimpp + 1, jpi+1 ) )
227        mi1(ji) = MAX( 0 , MIN( ji - nimpp + 1, jpi   ) )
228      END DO
229      DO jj = 1, jpjglo
230        mj0(jj) = MAX( 1 , MIN( jj - njmpp + 1, jpj+1 ) )
231        mj1(jj) = MAX( 0 , MIN( jj - njmpp + 1, jpj   ) )
232      END DO
233      IF(lwp) THEN                   ! control print
234         WRITE(numout,*)
235         WRITE(numout,*) 'dom_glo : domain: global <<==>> local '
236         WRITE(numout,*) '~~~~~~~ '
237         WRITE(numout,*) '   global domain:   jpiglo = ', jpiglo, ' jpjglo = ', jpjglo, ' jpkglo = ', jpkglo
238         WRITE(numout,*) '   local  domain:   jpi    = ', jpi   , ' jpj    = ', jpj   , ' jpk    = ', jpk
239         WRITE(numout,*)
240         WRITE(numout,*) '   conversion from local to global domain indices (and vise versa) done'
241         IF( nn_print >= 1 ) THEN
242            WRITE(numout,*)
243            WRITE(numout,*) '          conversion local  ==> global i-index domain'
244            WRITE(numout,25)              (mig(ji),ji = 1,jpi)
245            WRITE(numout,*)
246            WRITE(numout,*) '          conversion global ==> local  i-index domain'
247            WRITE(numout,*) '             starting index'
248            WRITE(numout,25)              (mi0(ji),ji = 1,jpiglo)
249            WRITE(numout,*) '             ending index'
250            WRITE(numout,25)              (mi1(ji),ji = 1,jpiglo)
251            WRITE(numout,*)
252            WRITE(numout,*) '          conversion local  ==> global j-index domain'
253            WRITE(numout,25)              (mjg(jj),jj = 1,jpj)
254            WRITE(numout,*)
255            WRITE(numout,*) '          conversion global ==> local  j-index domain'
256            WRITE(numout,*) '             starting index'
257            WRITE(numout,25)              (mj0(jj),jj = 1,jpjglo)
258            WRITE(numout,*) '             ending index'
259            WRITE(numout,25)              (mj1(jj),jj = 1,jpjglo)
260         ENDIF
261      ENDIF
262 25   FORMAT( 100(10x,19i4,/) )
263      !
264   END SUBROUTINE dom_glo
265
266
267   SUBROUTINE dom_nam
268      !!----------------------------------------------------------------------
269      !!                     ***  ROUTINE dom_nam  ***
270      !!                   
271      !! ** Purpose :   read domaine namelists and print the variables.
272      !!
273      !! ** input   : - namrun namelist
274      !!              - namdom namelist
275      !!              - namnc4 namelist   ! "key_netcdf4" only
276      !!----------------------------------------------------------------------
277      USE ioipsl
278      NAMELIST/namrun/ cn_ocerst_indir, cn_ocerst_outdir, nn_stocklist, ln_rst_list,                 &
279         &             nn_no   , cn_exp   , cn_ocerst_in, cn_ocerst_out, ln_rstart , nn_rstctl ,     &
280         &             nn_it000, nn_itend , nn_date0    , nn_time0     , nn_leapy  , nn_istate ,     &
281         &             nn_stock, nn_write , ln_mskland  , ln_clobber   , nn_chunksz, nn_euler  ,     &
282         &             ln_cfmeta, ln_iscpl
283      NAMELIST/namdom/ ln_linssh, nn_closea, nn_msh, rn_isfhmin, rn_rdt, rn_atfp, ln_crs
284#if defined key_netcdf4
285      NAMELIST/namnc4/ nn_nchunks_i, nn_nchunks_j, nn_nchunks_k, ln_nc4zip
286#endif
287      INTEGER  ::   ios                 ! Local integer output status for namelist read
288      !!----------------------------------------------------------------------
289
290      REWIND( numnam_ref )              ! Namelist namrun in reference namelist : Parameters of the run
291      READ  ( numnam_ref, namrun, IOSTAT = ios, ERR = 901)
292901   IF( ios /= 0 )   CALL ctl_nam ( ios , 'namrun in reference namelist', lwp )
293
294      REWIND( numnam_cfg )              ! Namelist namrun in configuration namelist : Parameters of the run
295      READ  ( numnam_cfg, namrun, IOSTAT = ios, ERR = 902 )
296902   IF( ios /= 0 )   CALL ctl_nam ( ios , 'namrun in configuration namelist', lwp )
297      IF(lwm) WRITE ( numond, namrun )
298      !
299      IF(lwp) THEN                  ! control print
300         WRITE(numout,*)
301         WRITE(numout,*) 'dom_nam  : domain initialization through namelist read'
302         WRITE(numout,*) '~~~~~~~ '
303         WRITE(numout,*) '   Namelist namrun'
304         WRITE(numout,*) '      job number                      nn_no      = ', nn_no
305         WRITE(numout,*) '      experiment name for output      cn_exp     = ', cn_exp
306         WRITE(numout,*) '      file prefix restart input       cn_ocerst_in= ', cn_ocerst_in
307         WRITE(numout,*) '      restart input directory         cn_ocerst_indir= ', cn_ocerst_indir
308         WRITE(numout,*) '      file prefix restart output      cn_ocerst_out= ', cn_ocerst_out
309         WRITE(numout,*) '      restart output directory        cn_ocerst_outdir= ', cn_ocerst_outdir
310         WRITE(numout,*) '      restart logical                 ln_rstart  = ', ln_rstart
311         WRITE(numout,*) '      start with forward time step    nn_euler   = ', nn_euler
312         WRITE(numout,*) '      control of time step            nn_rstctl  = ', nn_rstctl
313         WRITE(numout,*) '      number of the first time step   nn_it000   = ', nn_it000
314         WRITE(numout,*) '      number of the last time step    nn_itend   = ', nn_itend
315         WRITE(numout,*) '      initial calendar date aammjj    nn_date0   = ', nn_date0
316         WRITE(numout,*) '      initial time of day in hhmm     nn_time0   = ', nn_time0
317         WRITE(numout,*) '      leap year calendar (0/1)        nn_leapy   = ', nn_leapy
318         WRITE(numout,*) '      initial state output            nn_istate  = ', nn_istate
319         IF( ln_rst_list ) THEN
320            WRITE(numout,*) '      list of restart dump times      nn_stocklist   =', nn_stocklist
321         ELSE
322            WRITE(numout,*) '      frequency of restart file       nn_stock   = ', nn_stock
323         ENDIF
324         WRITE(numout,*) '      frequency of output file        nn_write   = ', nn_write
325         WRITE(numout,*) '      mask land points                ln_mskland = ', ln_mskland
326         WRITE(numout,*) '      additional CF standard metadata ln_cfmeta  = ', ln_cfmeta
327         WRITE(numout,*) '      overwrite an existing file      ln_clobber = ', ln_clobber
328         WRITE(numout,*) '      NetCDF chunksize (bytes)        nn_chunksz = ', nn_chunksz
329         WRITE(numout,*) '      IS coupling at the restart step ln_iscpl   = ', ln_iscpl
330      ENDIF
331
332      no = nn_no                    ! conversion DOCTOR names into model names (this should disappear soon)
333      cexper = cn_exp
334      nrstdt = nn_rstctl
335      nit000 = nn_it000
336      nitend = nn_itend
337      ndate0 = nn_date0
338      nleapy = nn_leapy
339      ninist = nn_istate
340      nstock = nn_stock
341      nstocklist = nn_stocklist
342      nwrite = nn_write
343      neuler = nn_euler
344      IF ( neuler == 1 .AND. .NOT. ln_rstart ) THEN
345         WRITE(ctmp1,*) 'ln_rstart =.FALSE., nn_euler is forced to 0 '
346         CALL ctl_warn( ctmp1 )
347         neuler = 0
348      ENDIF
349      !                             ! control of output frequency
350      IF ( nstock == 0 .OR. nstock > nitend ) THEN
351         WRITE(ctmp1,*) 'nstock = ', nstock, ' it is forced to ', nitend
352         CALL ctl_warn( ctmp1 )
353         nstock = nitend
354      ENDIF
355      IF ( nwrite == 0 ) THEN
356         WRITE(ctmp1,*) 'nwrite = ', nwrite, ' it is forced to ', nitend
357         CALL ctl_warn( ctmp1 )
358         nwrite = nitend
359      ENDIF
360
361#if defined key_agrif
362      IF( Agrif_Root() ) THEN
363#endif
364      SELECT CASE ( nleapy )        ! Choose calendar for IOIPSL
365      CASE (  1 ) 
366         CALL ioconf_calendar('gregorian')
367         IF(lwp) WRITE(numout,*) '   The IOIPSL calendar is "gregorian", i.e. leap year'
368      CASE (  0 )
369         CALL ioconf_calendar('noleap')
370         IF(lwp) WRITE(numout,*) '   The IOIPSL calendar is "noleap", i.e. no leap year'
371      CASE ( 30 )
372         CALL ioconf_calendar('360d')
373         IF(lwp) WRITE(numout,*) '   The IOIPSL calendar is "360d", i.e. 360 days in a year'
374      END SELECT
375#if defined key_agrif
376      ENDIF
377#endif
378
379      REWIND( numnam_ref )              ! Namelist namdom in reference namelist : space & time domain (bathymetry, mesh, timestep)
380      READ  ( numnam_ref, namdom, IOSTAT = ios, ERR = 903)
381903   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namdom in reference namelist', lwp )
382 
383      !
384      REWIND( numnam_cfg )              ! Namelist namdom in configuration namelist : space & time domain (bathymetry, mesh, timestep)
385      READ  ( numnam_cfg, namdom, IOSTAT = ios, ERR = 904 )
386904   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namdom in configuration namelist', lwp )
387      IF(lwm) WRITE ( numond, namdom )
388      !
389      IF(lwp) THEN
390         WRITE(numout,*)
391         WRITE(numout,*) '   Namelist namdom : space & time domain'
392         WRITE(numout,*) '      linear free surface (=T)              ln_linssh  = ', ln_linssh
393         WRITE(numout,*) '      suppression of closed seas (=0)       nn_closea  = ', nn_closea
394         WRITE(numout,*) '      create mesh/mask file(s)              nn_msh     = ', nn_msh
395         WRITE(numout,*) '           = 0   no file created           '
396         WRITE(numout,*) '           = 1   mesh_mask                 '
397         WRITE(numout,*) '           = 2   mesh and mask             '
398         WRITE(numout,*) '           = 3   mesh_hgr, msh_zgr and mask'
399         WRITE(numout,*) '      treshold to open the isf cavity       rn_isfhmin = ', rn_isfhmin, ' (m)'
400         WRITE(numout,*) '      ocean time step                       rn_rdt     = ', rn_rdt
401         WRITE(numout,*) '      asselin time filter parameter         rn_atfp    = ', rn_atfp
402         WRITE(numout,*) '      online coarsening of dynamical fields ln_crs     = ', ln_crs
403      ENDIF
404     
405      call flush( numout )
406      !
407!     !          ! conversion DOCTOR names into model names (this should disappear soon)
408      atfp      = rn_atfp
409      rdt       = rn_rdt
410
411#if defined key_netcdf4
412      !                             ! NetCDF 4 case   ("key_netcdf4" defined)
413      REWIND( numnam_ref )              ! Namelist namnc4 in reference namelist : NETCDF
414      READ  ( numnam_ref, namnc4, IOSTAT = ios, ERR = 907)
415907   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namnc4 in reference namelist', lwp )
416
417      REWIND( numnam_cfg )              ! Namelist namnc4 in configuration namelist : NETCDF
418      READ  ( numnam_cfg, namnc4, IOSTAT = ios, ERR = 908 )
419908   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namnc4 in configuration namelist', lwp )
420      IF(lwm) WRITE( numond, namnc4 )
421
422      IF(lwp) THEN                        ! control print
423         WRITE(numout,*)
424         WRITE(numout,*) '   Namelist namnc4 - Netcdf4 chunking parameters'
425         WRITE(numout,*) '      number of chunks in i-dimension      nn_nchunks_i   = ', nn_nchunks_i
426         WRITE(numout,*) '      number of chunks in j-dimension      nn_nchunks_j   = ', nn_nchunks_j
427         WRITE(numout,*) '      number of chunks in k-dimension      nn_nchunks_k   = ', nn_nchunks_k
428         WRITE(numout,*) '      apply netcdf4/hdf5 chunking & compression ln_nc4zip = ', ln_nc4zip
429      ENDIF
430
431      ! Put the netcdf4 settings into a simple structure (snc4set, defined in in_out_manager module)
432      ! Note the chunk size in the unlimited (time) dimension will be fixed at 1
433      snc4set%ni   = nn_nchunks_i
434      snc4set%nj   = nn_nchunks_j
435      snc4set%nk   = nn_nchunks_k
436      snc4set%luse = ln_nc4zip
437#else
438      snc4set%luse = .FALSE.        ! No NetCDF 4 case
439#endif
440      !
441   END SUBROUTINE dom_nam
442
443
444   SUBROUTINE dom_ctl
445      !!----------------------------------------------------------------------
446      !!                     ***  ROUTINE dom_ctl  ***
447      !!
448      !! ** Purpose :   Domain control.
449      !!
450      !! ** Method  :   compute and print extrema of masked scale factors
451      !!----------------------------------------------------------------------
452      INTEGER ::   iimi1, ijmi1, iimi2, ijmi2, iima1, ijma1, iima2, ijma2
453      INTEGER, DIMENSION(2) ::   iloc   !
454      REAL(wp) ::   ze1min, ze1max, ze2min, ze2max
455      !!----------------------------------------------------------------------
456      !
457      IF(lk_mpp) THEN
458         CALL mpp_minloc( e1t(:,:), tmask_i(:,:), ze1min, iimi1,ijmi1 )
459         CALL mpp_minloc( e2t(:,:), tmask_i(:,:), ze2min, iimi2,ijmi2 )
460         CALL mpp_maxloc( e1t(:,:), tmask_i(:,:), ze1max, iima1,ijma1 )
461         CALL mpp_maxloc( e2t(:,:), tmask_i(:,:), ze2max, iima2,ijma2 )
462      ELSE
463         ze1min = MINVAL( e1t(:,:), mask = tmask_i(:,:) == 1._wp )   
464         ze2min = MINVAL( e2t(:,:), mask = tmask_i(:,:) == 1._wp )   
465         ze1max = MAXVAL( e1t(:,:), mask = tmask_i(:,:) == 1._wp )   
466         ze2max = MAXVAL( e2t(:,:), mask = tmask_i(:,:) == 1._wp )   
467         !
468         iloc  = MINLOC( e1t(:,:), mask = tmask_i(:,:) == 1._wp )
469         iimi1 = iloc(1) + nimpp - 1
470         ijmi1 = iloc(2) + njmpp - 1
471         iloc  = MINLOC( e2t(:,:), mask = tmask_i(:,:) == 1._wp )
472         iimi2 = iloc(1) + nimpp - 1
473         ijmi2 = iloc(2) + njmpp - 1
474         iloc  = MAXLOC( e1t(:,:), mask = tmask_i(:,:) == 1._wp )
475         iima1 = iloc(1) + nimpp - 1
476         ijma1 = iloc(2) + njmpp - 1
477         iloc  = MAXLOC( e2t(:,:), mask = tmask_i(:,:) == 1._wp )
478         iima2 = iloc(1) + nimpp - 1
479         ijma2 = iloc(2) + njmpp - 1
480      ENDIF
481      IF(lwp) THEN
482         WRITE(numout,*)
483         WRITE(numout,*) 'dom_ctl : extrema of the masked scale factors'
484         WRITE(numout,*) '~~~~~~~'
485         WRITE(numout,"(14x,'e1t maxi: ',1f10.2,' at i = ',i5,' j= ',i5)") ze1max, iima1, ijma1
486         WRITE(numout,"(14x,'e1t mini: ',1f10.2,' at i = ',i5,' j= ',i5)") ze1min, iimi1, ijmi1
487         WRITE(numout,"(14x,'e2t maxi: ',1f10.2,' at i = ',i5,' j= ',i5)") ze2max, iima2, ijma2
488         WRITE(numout,"(14x,'e2t mini: ',1f10.2,' at i = ',i5,' j= ',i5)") ze2min, iimi2, ijmi2
489      ENDIF
490      !
491   END SUBROUTINE dom_ctl
492
493
494   SUBROUTINE cfg_write
495      !!----------------------------------------------------------------------
496      !!                  ***  ROUTINE cfg_write  ***
497      !!                   
498      !! ** Purpose :   Create the "domain_cfg" file, a NetCDF file which
499      !!              contains all the ocean domain informations required to
500      !!              define an ocean configuration.
501      !!
502      !! ** Method  :   Write in a file all the arrays required to set up an
503      !!              ocean configuration.
504      !!
505      !! ** output file :   domain_cfg.nc : domain size, characteristics, horizontal mesh,
506      !!                              Coriolis parameter, depth and vertical scale factors
507      !!----------------------------------------------------------------------
508      INTEGER           ::   ji, jj, jk   ! dummy loop indices
509      INTEGER           ::   izco, izps, isco, icav
510      INTEGER           ::   inum     ! temprary units for 'domain_cfg.nc' file
511      CHARACTER(len=21) ::   clnam    ! filename (mesh and mask informations)
512      REAL(wp), DIMENSION(jpi,jpj) ::   z2d   ! workspace
513      !!----------------------------------------------------------------------
514      !
515      IF(lwp) WRITE(numout,*)
516      IF(lwp) WRITE(numout,*) 'cfg_write : create the "domain_cfg.nc" file containing all required configuration information'
517      IF(lwp) WRITE(numout,*) '~~~~~~~~~'
518      !
519      !                       ! ============================= !
520      !                       !  create 'domain_cfg.nc' file  !
521      !                       ! ============================= !
522      !         
523      clnam = 'domain_cfg'  ! filename (configuration information)
524      CALL iom_open( TRIM(clnam), inum, ldwrt = .TRUE., kiolib = jprstlib )
525     
526      !                             !==  global domain size  ==!
527      !
528      CALL iom_rstput( 0, 0, inum, 'jpiglo', REAL( jpiglo, wp), ktype = jp_i4 )
529      CALL iom_rstput( 0, 0, inum, 'jpjglo', REAL( jpjglo, wp), ktype = jp_i4 )
530      CALL iom_rstput( 0, 0, inum, 'jpkglo', REAL( jpk   , wp), ktype = jp_i4 )
531      !
532      !                             !==  domain characteristics  ==!
533      !
534      !                                   ! lateral boundary of the global domain
535      CALL iom_rstput( 0, 0, inum, 'jperio', REAL( jperio, wp), ktype = jp_i4 )
536      !
537      !                                   ! type of vertical coordinate
538      IF( ln_zco    ) THEN   ;   izco = 1   ;   ELSE   ;   izco = 0   ;   ENDIF
539      IF( ln_zps    ) THEN   ;   izps = 1   ;   ELSE   ;   izps = 0   ;   ENDIF
540      IF( ln_sco    ) THEN   ;   isco = 1   ;   ELSE   ;   isco = 0   ;   ENDIF
541      CALL iom_rstput( 0, 0, inum, 'ln_zco'   , REAL( izco, wp), ktype = jp_i4 )
542      CALL iom_rstput( 0, 0, inum, 'ln_zps'   , REAL( izps, wp), ktype = jp_i4 )
543      CALL iom_rstput( 0, 0, inum, 'ln_sco'   , REAL( isco, wp), ktype = jp_i4 )
544      !
545      !                                   ! ocean cavities under iceshelves
546      IF( ln_isfcav ) THEN   ;   icav = 1   ;   ELSE   ;   icav = 0   ;   ENDIF
547      CALL iom_rstput( 0, 0, inum, 'ln_isfcav', REAL( icav, wp), ktype = jp_i4 )
548      !
549      !                             !==  horizontal mesh  !
550      !
551      CALL iom_rstput( 0, 0, inum, 'glamt', glamt, ktype = jp_r8 )   ! latitude
552      CALL iom_rstput( 0, 0, inum, 'glamu', glamu, ktype = jp_r8 )
553      CALL iom_rstput( 0, 0, inum, 'glamv', glamv, ktype = jp_r8 )
554      CALL iom_rstput( 0, 0, inum, 'glamf', glamf, ktype = jp_r8 )
555      !                               
556      CALL iom_rstput( 0, 0, inum, 'gphit', gphit, ktype = jp_r8 )   ! longitude
557      CALL iom_rstput( 0, 0, inum, 'gphiu', gphiu, ktype = jp_r8 )
558      CALL iom_rstput( 0, 0, inum, 'gphiv', gphiv, ktype = jp_r8 )
559      CALL iom_rstput( 0, 0, inum, 'gphif', gphif, ktype = jp_r8 )
560      !                               
561      CALL iom_rstput( 0, 0, inum, 'e1t'  , e1t  , ktype = jp_r8 )   ! i-scale factors (e1.)
562      CALL iom_rstput( 0, 0, inum, 'e1u'  , e1u  , ktype = jp_r8 )
563      CALL iom_rstput( 0, 0, inum, 'e1v'  , e1v  , ktype = jp_r8 )
564      CALL iom_rstput( 0, 0, inum, 'e1f'  , e1f  , ktype = jp_r8 )
565      !
566      CALL iom_rstput( 0, 0, inum, 'e2t'  , e2t  , ktype = jp_r8 )   ! j-scale factors (e2.)
567      CALL iom_rstput( 0, 0, inum, 'e2u'  , e2u  , ktype = jp_r8 )
568      CALL iom_rstput( 0, 0, inum, 'e2v'  , e2v  , ktype = jp_r8 )
569      CALL iom_rstput( 0, 0, inum, 'e2f'  , e2f  , ktype = jp_r8 )
570      !
571      CALL iom_rstput( 0, 0, inum, 'ff_f' , ff_f , ktype = jp_r8 )   ! coriolis factor
572      CALL iom_rstput( 0, 0, inum, 'ff_t' , ff_t , ktype = jp_r8 )
573      !
574      !                             !==  vertical mesh  ==!
575      !                                                     
576      CALL iom_rstput( 0, 0, inum, 'e3t_1d'  , e3t_1d  , ktype = jp_r8 )   ! reference 1D-coordinate
577      CALL iom_rstput( 0, 0, inum, 'e3w_1d'  , e3w_1d  , ktype = jp_r8 )
578      !
579      CALL iom_rstput( 0, 0, inum, 'e3t_0'   , e3t_0   , ktype = jp_r8 )   ! vertical scale factors
580      CALL iom_rstput( 0, 0, inum, 'e3u_0'   , e3u_0   , ktype = jp_r8 )
581      CALL iom_rstput( 0, 0, inum, 'e3v_0'   , e3v_0   , ktype = jp_r8 )
582      CALL iom_rstput( 0, 0, inum, 'e3f_0'   , e3f_0   , ktype = jp_r8 )
583      CALL iom_rstput( 0, 0, inum, 'e3w_0'   , e3w_0   , ktype = jp_r8 )
584      CALL iom_rstput( 0, 0, inum, 'e3uw_0'  , e3uw_0  , ktype = jp_r8 )
585      CALL iom_rstput( 0, 0, inum, 'e3vw_0'  , e3vw_0  , ktype = jp_r8 )
586      !                                         
587      !                             !==  wet top and bottom level  ==!   (caution: multiplied by ssmask)
588      !
589      CALL iom_rstput( 0, 0, inum, 'top_level'    , REAL( mikt, wp )*ssmask , ktype = jp_i4 )   ! nb of ocean T-points (ISF)
590      CALL iom_rstput( 0, 0, inum, 'bottom_level' , REAL( mbkt, wp )*ssmask , ktype = jp_i4 )   ! nb of ocean T-points
591      !
592      IF( ln_sco ) THEN             ! s-coordinate: store grid stiffness ratio  (Not required anyway)
593         CALL dom_stiff( z2d )
594         CALL iom_rstput( 0, 0, inum, 'stiffness', z2d )        !    ! Max. grid stiffness ratio
595      ENDIF
596      !
597      !                                ! ============================
598      !                                !        close the files
599      !                                ! ============================
600      CALL iom_close( inum )
601      !
602   END SUBROUTINE cfg_write
603
604   !!======================================================================
605END MODULE domain
Note: See TracBrowser for help on using the repository browser.