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 branches/2017/dev_merge_2017/NEMOGCM/CONFIG/TEST_CASES/WAD/MY_SRC – NEMO

source: branches/2017/dev_merge_2017/NEMOGCM/CONFIG/TEST_CASES/WAD/MY_SRC/domain.F90 @ 9124

Last change on this file since 9124 was 9124, checked in by gm, 6 years ago

dev_merge_2017: ln_timing instead of nn_timing + restricted timing to nemo_init and routine called by step in OPA_SRC

File size: 36.1 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 usrdef_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 domc1d         ! 1D configuration: column location
40   USE dyncor_c1d     ! 1D configuration: Coriolis term    (cor_c1d routine)
41   USE wet_dry        ! wetting and drying
42   !
43   USE in_out_manager ! I/O manager
44   USE iom            ! I/O library
45   USE lbclnk         ! ocean lateral boundary condition (or mpp link)
46   USE lib_mpp        ! distributed memory computing library
47   USE wrk_nemo       ! Memory Allocation
48
49   IMPLICIT NONE
50   PRIVATE
51
52   PUBLIC   dom_init     ! called by nemogcm.F90
53   PUBLIC   domain_cfg   ! called by nemogcm.F90
54
55   !!-------------------------------------------------------------------------
56   !! NEMO/OPA 3.3 , NEMO Consortium (2010)
57   !! $Id: domain.F90 7587 2017-01-20 15:17:56Z acc $
58   !! Software governed by the CeCILL licence        (NEMOGCM/NEMO_CeCILL.txt)
59   !!-------------------------------------------------------------------------
60CONTAINS
61
62   SUBROUTINE dom_init
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 if nn_msh=1
75      !!              - 1D configuration, move Coriolis, u and v at T-point
76      !!----------------------------------------------------------------------
77      INTEGER ::   ji, jj, jk, ik   ! dummy loop indices
78      INTEGER ::   iconf = 0    ! local integers
79      CHARACTER (len=64) ::   cform = "(A12, 3(A13, I7))" 
80      INTEGER , DIMENSION(jpi,jpj) ::   ik_top , ik_bot       ! top and bottom ocean level
81      REAL(wp), DIMENSION(jpi,jpj) ::   z1_hu_0, z1_hv_0
82      !!----------------------------------------------------------------------
83      !
84      IF(lwp) THEN         ! Ocean domain Parameters (control print)
85         WRITE(numout,*)
86         WRITE(numout,*) 'dom_init : domain initialization'
87         WRITE(numout,*) '~~~~~~~~'
88         !
89         WRITE(numout,*)     '   Domain info'
90         WRITE(numout,*)     '      dimension of model:'
91         WRITE(numout,*)     '             Local domain      Global domain       Data domain '
92         WRITE(numout,cform) '        ','   jpi     : ', jpi, '   jpiglo  : ', jpiglo
93         WRITE(numout,cform) '        ','   jpj     : ', jpj, '   jpjglo  : ', jpjglo
94         WRITE(numout,cform) '        ','   jpk     : ', jpk, '   jpkglo  : ', jpkglo
95         WRITE(numout,cform) '       ' ,'   jpij    : ', jpij
96         WRITE(numout,*)     '      mpp local domain info (mpp):'
97         WRITE(numout,*)     '              jpni    : ', jpni, '   jpreci  : ', jpreci
98         WRITE(numout,*)     '              jpnj    : ', jpnj, '   jprecj  : ', jprecj
99         WRITE(numout,*)     '              jpnij   : ', jpnij
100         WRITE(numout,*)     '      lateral boundary of the Global domain : jperio  = ', jperio
101         SELECT CASE ( jperio )
102         CASE( 0 )   ;   WRITE(numout,*) '         (i.e. closed)'
103         CASE( 1 )   ;   WRITE(numout,*) '         (i.e. cyclic east-west)'
104         CASE( 2 )   ;   WRITE(numout,*) '         (i.e. equatorial symmetric)'
105         CASE( 3 )   ;   WRITE(numout,*) '         (i.e. north fold with T-point pivot)'
106         CASE( 4 )   ;   WRITE(numout,*) '         (i.e. cyclic east-west and north fold with T-point pivot)'
107         CASE( 5 )   ;   WRITE(numout,*) '         (i.e. north fold with F-point pivot)'
108         CASE( 6 )   ;   WRITE(numout,*) '         (i.e. cyclic east-west and north fold with F-point pivot)'
109         CASE( 7 )   ;   WRITE(numout,*) '         (i.e. cyclic north-south (jpnij=1 only))'
110         CASE DEFAULT
111            CALL ctl_stop( 'jperio is out of range' )
112         END SELECT
113         WRITE(numout,*)     '      Ocean model configuration used:'
114         WRITE(numout,*)     '              cn_cfg = ', cn_cfg
115         WRITE(numout,*)     '              nn_cfg = ', nn_cfg
116      ENDIF
117      !
118      !     
119!!gm  This should be removed with the new configuration interface
120      IF( lk_c1d .AND. ln_c1d_locpt )  CALL dom_c1d( rn_lat1d, rn_lon1d )
121!!gm end
122      !
123      !           !==  Reference coordinate system  ==!
124      !
125      CALL dom_glo                     ! global domain versus local domain
126      CALL dom_nam                     ! read namelist ( namrun, namdom )
127      CALL dom_clo( cn_cfg, nn_cfg )   ! Closed seas and lake
128      CALL dom_hgr                     ! Horizontal mesh
129      CALL dom_zgr( ik_top, ik_bot )   ! Vertical mesh and bathymetry
130      IF( nn_closea == 0 )   CALL clo_bat( ik_top, ik_bot )    !==  remove closed seas or lakes  ==!
131      CALL dom_msk( ik_top, ik_bot )   ! Masks
132      !
133      DO jj = 1, jpj                   ! depth of the iceshelves
134         DO ji = 1, jpi
135            ik = mikt(ji,jj)
136            risfdep(ji,jj) = gdepw_0(ji,jj,ik)
137         END DO
138      END DO
139      !
140      ht_0(:,:) = 0._wp  ! Reference ocean thickness
141      hu_0(:,:) = 0._wp
142      hv_0(:,:) = 0._wp
143      DO jk = 1, jpk
144         ht_0(:,:) = ht_0(:,:) + e3t_0(:,:,jk) * tmask(:,:,jk)
145         hu_0(:,:) = hu_0(:,:) + e3u_0(:,:,jk) * umask(:,:,jk)
146         hv_0(:,:) = hv_0(:,:) + e3v_0(:,:,jk) * vmask(:,:,jk)
147      END DO
148      !
149      !           !==  time varying part of coordinate system  ==!
150      !
151      IF( ln_linssh ) THEN       != Fix in time : set to the reference one for all
152      !
153         !       before        !          now          !       after         !
154            gdept_b = gdept_0  ;   gdept_n = gdept_0   !        ---          ! depth of grid-points
155            gdepw_b = gdepw_0  ;   gdepw_n = gdepw_0   !        ---          !
156                                   gde3w_n = gde3w_0   !        ---          !
157         !                                                                 
158              e3t_b =   e3t_0  ;     e3t_n =   e3t_0   ;   e3t_a =  e3t_0    ! scale factors
159              e3u_b =   e3u_0  ;     e3u_n =   e3u_0   ;   e3u_a =  e3u_0    !
160              e3v_b =   e3v_0  ;     e3v_n =   e3v_0   ;   e3v_a =  e3v_0    !
161                                     e3f_n =   e3f_0   !        ---          !
162              e3w_b =   e3w_0  ;     e3w_n =   e3w_0   !        ---          !
163             e3uw_b =  e3uw_0  ;    e3uw_n =  e3uw_0   !        ---          !
164             e3vw_b =  e3vw_0  ;    e3vw_n =  e3vw_0   !        ---          !
165         !
166         z1_hu_0(:,:) = ssumask(:,:) / ( hu_0(:,:) + 1._wp - ssumask(:,:) )     ! _i mask due to ISF
167         z1_hv_0(:,:) = ssvmask(:,:) / ( hv_0(:,:) + 1._wp - ssvmask(:,:) )
168         !
169         !        before       !          now          !       after         !
170                                      ht_n =    ht_0   !                     ! water column thickness
171               hu_b =    hu_0  ;      hu_n =    hu_0   ;    hu_a =    hu_0   !
172               hv_b =    hv_0  ;      hv_n =    hv_0   ;    hv_a =    hv_0   !
173            r1_hu_b = z1_hu_0  ;   r1_hu_n = z1_hu_0   ; r1_hu_a = z1_hu_0   ! inverse of water column thickness
174            r1_hv_b = z1_hv_0  ;   r1_hv_n = z1_hv_0   ; r1_hv_a = z1_hv_0   !
175         !
176         !
177      ELSE                       != time varying : initialize before/now/after variables
178         !
179         IF( .NOT.l_offline )  CALL dom_vvl_init 
180         !
181      ENDIF
182      !
183      IF( lk_c1d         )   CALL cor_c1d       ! 1D configuration: Coriolis set at T-point
184      !
185      IF( nn_msh > 0 .AND. .NOT. ln_iscpl )                         CALL dom_wri      ! Create a domain file
186      IF( nn_msh > 0 .AND.       ln_iscpl .AND. .NOT. ln_rstart )   CALL dom_wri      ! Create a domain file
187      IF( .NOT.ln_rstart )   CALL dom_ctl       ! Domain control
188      !
189     
190      IF(lwp) THEN
191         WRITE(numout,*)
192         WRITE(numout,*) 'dom_init : end of domain initialization nn_msh=', nn_msh
193         WRITE(numout,*) 
194      ENDIF
195      !
196      IF( ln_write_cfg )   CALL cfg_write         ! create the configuration file
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      REWIND( numnam_cfg )              ! Namelist namdom in configuration namelist : space & time domain (bathymetry, mesh, timestep)
384      READ  ( numnam_cfg, namdom, IOSTAT = ios, ERR = 904 )
385904   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namdom in configuration namelist', lwp )
386      IF(lwm) WRITE ( numond, namdom )
387      !
388      IF(lwp) THEN
389         WRITE(numout,*)
390         WRITE(numout,*) '   Namelist namdom : space & time domain'
391         WRITE(numout,*) '      linear free surface (=T)              ln_linssh  = ', ln_linssh
392         WRITE(numout,*) '      suppression of closed seas (=0)       nn_closea  = ', nn_closea
393         WRITE(numout,*) '      create mesh/mask file(s)              nn_msh     = ', nn_msh
394         WRITE(numout,*) '           = 0   no file created           '
395         WRITE(numout,*) '           = 1   mesh_mask                 '
396         WRITE(numout,*) '           = 2   mesh and mask             '
397         WRITE(numout,*) '           = 3   mesh_hgr, msh_zgr and mask'
398         WRITE(numout,*) '      treshold to open the isf cavity       rn_isfhmin = ', rn_isfhmin, ' (m)'
399         WRITE(numout,*) '      ocean time step                       rn_rdt     = ', rn_rdt
400         WRITE(numout,*) '      asselin time filter parameter         rn_atfp    = ', rn_atfp
401         WRITE(numout,*) '      online coarsening of dynamical fields ln_crs     = ', ln_crs
402      ENDIF
403     
404      call flush( numout )
405      !
406!     !          ! conversion DOCTOR names into model names (this should disappear soon)
407      atfp      = rn_atfp
408      rdt       = rn_rdt
409
410#if defined key_netcdf4
411      !                             ! NetCDF 4 case   ("key_netcdf4" defined)
412      REWIND( numnam_ref )              ! Namelist namnc4 in reference namelist : NETCDF
413      READ  ( numnam_ref, namnc4, IOSTAT = ios, ERR = 907)
414907   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namnc4 in reference namelist', lwp )
415      !
416      REWIND( numnam_cfg )              ! Namelist namnc4 in configuration namelist : NETCDF
417      READ  ( numnam_cfg, namnc4, IOSTAT = ios, ERR = 908 )
418908   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namnc4 in configuration namelist', lwp )
419      IF(lwm) WRITE( numond, namnc4 )
420
421      IF(lwp) THEN                        ! control print
422         WRITE(numout,*)
423         WRITE(numout,*) '   Namelist namnc4 - Netcdf4 chunking parameters'
424         WRITE(numout,*) '      number of chunks in i-dimension      nn_nchunks_i   = ', nn_nchunks_i
425         WRITE(numout,*) '      number of chunks in j-dimension      nn_nchunks_j   = ', nn_nchunks_j
426         WRITE(numout,*) '      number of chunks in k-dimension      nn_nchunks_k   = ', nn_nchunks_k
427         WRITE(numout,*) '      apply netcdf4/hdf5 chunking & compression ln_nc4zip = ', ln_nc4zip
428      ENDIF
429
430      ! Put the netcdf4 settings into a simple structure (snc4set, defined in in_out_manager module)
431      ! Note the chunk size in the unlimited (time) dimension will be fixed at 1
432      snc4set%ni   = nn_nchunks_i
433      snc4set%nj   = nn_nchunks_j
434      snc4set%nk   = nn_nchunks_k
435      snc4set%luse = ln_nc4zip
436#else
437      snc4set%luse = .FALSE.        ! No NetCDF 4 case
438#endif
439      !
440   END SUBROUTINE dom_nam
441
442
443   SUBROUTINE dom_ctl
444      !!----------------------------------------------------------------------
445      !!                     ***  ROUTINE dom_ctl  ***
446      !!
447      !! ** Purpose :   Domain control.
448      !!
449      !! ** Method  :   compute and print extrema of masked scale factors
450      !!----------------------------------------------------------------------
451      INTEGER ::   iimi1, ijmi1, iimi2, ijmi2, iima1, ijma1, iima2, ijma2
452      INTEGER, DIMENSION(2) ::   iloc   !
453      REAL(wp) ::   ze1min, ze1max, ze2min, ze2max
454      !!----------------------------------------------------------------------
455      !
456      IF(lk_mpp) THEN
457         CALL mpp_minloc( e1t(:,:), tmask_i(:,:), ze1min, iimi1,ijmi1 )
458         CALL mpp_minloc( e2t(:,:), tmask_i(:,:), ze2min, iimi2,ijmi2 )
459         CALL mpp_maxloc( e1t(:,:), tmask_i(:,:), ze1max, iima1,ijma1 )
460         CALL mpp_maxloc( e2t(:,:), tmask_i(:,:), ze2max, iima2,ijma2 )
461      ELSE
462         ze1min = MINVAL( e1t(:,:), mask = tmask_i(:,:) == 1._wp )   
463         ze2min = MINVAL( e2t(:,:), mask = tmask_i(:,:) == 1._wp )   
464         ze1max = MAXVAL( e1t(:,:), mask = tmask_i(:,:) == 1._wp )   
465         ze2max = MAXVAL( e2t(:,:), mask = tmask_i(:,:) == 1._wp )   
466         !
467         iloc  = MINLOC( e1t(:,:), mask = tmask_i(:,:) == 1._wp )
468         iimi1 = iloc(1) + nimpp - 1
469         ijmi1 = iloc(2) + njmpp - 1
470         iloc  = MINLOC( e2t(:,:), mask = tmask_i(:,:) == 1._wp )
471         iimi2 = iloc(1) + nimpp - 1
472         ijmi2 = iloc(2) + njmpp - 1
473         iloc  = MAXLOC( e1t(:,:), mask = tmask_i(:,:) == 1._wp )
474         iima1 = iloc(1) + nimpp - 1
475         ijma1 = iloc(2) + njmpp - 1
476         iloc  = MAXLOC( e2t(:,:), mask = tmask_i(:,:) == 1._wp )
477         iima2 = iloc(1) + nimpp - 1
478         ijma2 = iloc(2) + njmpp - 1
479      ENDIF
480      IF(lwp) THEN
481         WRITE(numout,*)
482         WRITE(numout,*) 'dom_ctl : extrema of the masked scale factors'
483         WRITE(numout,*) '~~~~~~~'
484         WRITE(numout,"(14x,'e1t maxi: ',1f10.2,' at i = ',i5,' j= ',i5)") ze1max, iima1, ijma1
485         WRITE(numout,"(14x,'e1t mini: ',1f10.2,' at i = ',i5,' j= ',i5)") ze1min, iimi1, ijmi1
486         WRITE(numout,"(14x,'e2t maxi: ',1f10.2,' at i = ',i5,' j= ',i5)") ze2max, iima2, ijma2
487         WRITE(numout,"(14x,'e2t mini: ',1f10.2,' at i = ',i5,' j= ',i5)") ze2min, iimi2, ijmi2
488      ENDIF
489      !
490   END SUBROUTINE dom_ctl
491
492
493   SUBROUTINE domain_cfg( ldtxt, cd_cfg, kk_cfg, kpi, kpj, kpk, kperio )
494      !!----------------------------------------------------------------------
495      !!                     ***  ROUTINE dom_nam  ***
496      !!                   
497      !! ** Purpose :   read the domain size in domain configuration file
498      !!
499      !! ** Method  :   
500      !!
501      !!----------------------------------------------------------------------
502      CHARACTER(len=*), DIMENSION(:), INTENT(out) ::   ldtxt           ! stored print information
503      CHARACTER(len=*)              , INTENT(out) ::   cd_cfg          ! configuration name
504      INTEGER                       , INTENT(out) ::   kk_cfg          ! configuration resolution
505      INTEGER                       , INTENT(out) ::   kpi, kpj, kpk   ! global domain sizes
506      INTEGER                       , INTENT(out) ::   kperio          ! lateral global domain b.c.
507      !
508      INTEGER ::   inum, ii   ! local integer
509      REAL(wp) ::   zorca_res                     ! local scalars
510      REAL(wp) ::   ziglo, zjglo, zkglo, zperio   !   -      -
511      !!----------------------------------------------------------------------
512      !
513      ii = 1
514      WRITE(ldtxt(ii),*) '           '                                                    ;   ii = ii+1
515      WRITE(ldtxt(ii),*) 'domain_cfg : domain size read in ', TRIM( cn_domcfg ), ' file'   ;   ii = ii+1
516      WRITE(ldtxt(ii),*) '~~~~~~~~~~ '                                                    ;   ii = ii+1
517      !
518      CALL iom_open( cn_domcfg, inum )
519      !
520      !                                   !- ORCA family specificity
521      IF(  iom_varid( inum, 'ORCA'       , ldstop = .FALSE. ) > 0  .AND.  &
522         & iom_varid( inum, 'ORCA_index' , ldstop = .FALSE. ) > 0    ) THEN
523         !
524         cd_cfg = 'ORCA'
525         CALL iom_get( inum, 'ORCA_index', zorca_res )   ;   kk_cfg = INT( zorca_res )
526         !
527         WRITE(ldtxt(ii),*) '       '                                                    ;   ii = ii+1
528         WRITE(ldtxt(ii),*) '       ==>>>   ORCA configuration '                         ;   ii = ii+1
529         WRITE(ldtxt(ii),*) '       '                                                    ;   ii = ii+1
530         !
531      ELSE                                !- cd_cfg & k_cfg are not used
532         cd_cfg = 'UNKNOWN'
533         kk_cfg = -9999999
534                                          !- or they may be present as global attributes
535                                          !- (netcdf only) 
536         IF( iom_file(inum)%iolib == jpnf90 ) THEN
537            CALL iom_getatt( inum, 'cn_cfg', cd_cfg )  ! returns   !  if not found
538            CALL iom_getatt( inum, 'nn_cfg', kk_cfg )  ! returns -999 if not found
539            IF( TRIM(cd_cfg) .EQ. '!') cd_cfg = 'UNKNOWN'
540            IF( kk_cfg .EQ. -999     ) kk_cfg = -9999999
541         ENDIF
542         !
543      ENDIF
544      !
545      CALL iom_get( inum, 'jpiglo', ziglo  )   ;   kpi = INT( ziglo )
546      CALL iom_get( inum, 'jpjglo', zjglo  )   ;   kpj = INT( zjglo )
547      CALL iom_get( inum, 'jpkglo', zkglo  )   ;   kpk = INT( zkglo )
548      CALL iom_get( inum, 'jperio', zperio )   ;   kperio = INT( zperio )
549      CALL iom_close( inum )
550      !
551      WRITE(ldtxt(ii),*) '   cn_cfg = ', TRIM(cd_cfg), '   nn_cfg = ', kk_cfg             ;   ii = ii+1
552      WRITE(ldtxt(ii),*) '   jpiglo = ', kpi                                              ;   ii = ii+1
553      WRITE(ldtxt(ii),*) '   jpjglo = ', kpj                                              ;   ii = ii+1
554      WRITE(ldtxt(ii),*) '   jpkglo = ', kpk                                              ;   ii = ii+1
555      WRITE(ldtxt(ii),*) '   type of global domain lateral boundary   jperio = ', kperio  ;   ii = ii+1
556      !       
557   END SUBROUTINE domain_cfg
558   
559   
560   SUBROUTINE cfg_write
561      !!----------------------------------------------------------------------
562      !!                  ***  ROUTINE cfg_write  ***
563      !!                   
564      !! ** Purpose :   Create the "cn_domcfg_out" file, a NetCDF file which
565      !!              contains all the ocean domain informations required to
566      !!              define an ocean configuration.
567      !!
568      !! ** Method  :   Write in a file all the arrays required to set up an
569      !!              ocean configuration.
570      !!
571      !! ** output file :   domcfg_out.nc : domain size, characteristics, horizontal
572      !!                       mesh, Coriolis parameter, and vertical scale factors
573      !!                    NB: also contain ORCA family information
574      !!----------------------------------------------------------------------
575      INTEGER           ::   ji, jj, jk   ! dummy loop indices
576      INTEGER           ::   izco, izps, isco, icav
577      INTEGER           ::   inum     ! local units
578      CHARACTER(len=21) ::   clnam    ! filename (mesh and mask informations)
579      REAL(wp), DIMENSION(jpi,jpj) ::   z2d   ! workspace
580      !!----------------------------------------------------------------------
581      !
582      IF(lwp) WRITE(numout,*)
583      IF(lwp) WRITE(numout,*) 'cfg_write : create the domain configuration file (', TRIM(cn_domcfg_out),'.nc)'
584      IF(lwp) WRITE(numout,*) '~~~~~~~~~'
585      !
586      !                       ! ============================= !
587      !                       !  create 'domcfg_out.nc' file  !
588      !                       ! ============================= !
589      !         
590      clnam = 'domcfg_out'  ! filename (configuration information)
591      CALL iom_open( TRIM(clnam), inum, ldwrt = .TRUE., kiolib = jprstlib )
592     
593      !
594      !                             !==  ORCA family specificities  ==!
595      IF( cn_cfg == "ORCA" ) THEN
596         CALL iom_rstput( 0, 0, inum, 'ORCA'      , 1._wp            , ktype = jp_i4 )
597         CALL iom_rstput( 0, 0, inum, 'ORCA_index', REAL( nn_cfg, wp), ktype = jp_i4 )         
598      ENDIF
599      !
600      !                             !==  global domain size  ==!
601      !
602      CALL iom_rstput( 0, 0, inum, 'jpiglo', REAL( jpiglo, wp), ktype = jp_i4 )
603      CALL iom_rstput( 0, 0, inum, 'jpjglo', REAL( jpjglo, wp), ktype = jp_i4 )
604      CALL iom_rstput( 0, 0, inum, 'jpkglo', REAL( jpk   , wp), ktype = jp_i4 )
605      !
606      !                             !==  domain characteristics  ==!
607      !
608      !                                   ! lateral boundary of the global domain
609      CALL iom_rstput( 0, 0, inum, 'jperio', REAL( jperio, wp), ktype = jp_i4 )
610      !
611      !                                   ! type of vertical coordinate
612      IF( ln_zco    ) THEN   ;   izco = 1   ;   ELSE   ;   izco = 0   ;   ENDIF
613      IF( ln_zps    ) THEN   ;   izps = 1   ;   ELSE   ;   izps = 0   ;   ENDIF
614      IF( ln_sco    ) THEN   ;   isco = 1   ;   ELSE   ;   isco = 0   ;   ENDIF
615      CALL iom_rstput( 0, 0, inum, 'ln_zco'   , REAL( izco, wp), ktype = jp_i4 )
616      CALL iom_rstput( 0, 0, inum, 'ln_zps'   , REAL( izps, wp), ktype = jp_i4 )
617      CALL iom_rstput( 0, 0, inum, 'ln_sco'   , REAL( isco, wp), ktype = jp_i4 )
618      !
619      !                                   ! ocean cavities under iceshelves
620      IF( ln_isfcav ) THEN   ;   icav = 1   ;   ELSE   ;   icav = 0   ;   ENDIF
621      CALL iom_rstput( 0, 0, inum, 'ln_isfcav', REAL( icav, wp), ktype = jp_i4 )
622      !
623      !                             !==  horizontal mesh  !
624      !
625      CALL iom_rstput( 0, 0, inum, 'glamt', glamt, ktype = jp_r8 )   ! latitude
626      CALL iom_rstput( 0, 0, inum, 'glamu', glamu, ktype = jp_r8 )
627      CALL iom_rstput( 0, 0, inum, 'glamv', glamv, ktype = jp_r8 )
628      CALL iom_rstput( 0, 0, inum, 'glamf', glamf, ktype = jp_r8 )
629      !                               
630      CALL iom_rstput( 0, 0, inum, 'gphit', gphit, ktype = jp_r8 )   ! longitude
631      CALL iom_rstput( 0, 0, inum, 'gphiu', gphiu, ktype = jp_r8 )
632      CALL iom_rstput( 0, 0, inum, 'gphiv', gphiv, ktype = jp_r8 )
633      CALL iom_rstput( 0, 0, inum, 'gphif', gphif, ktype = jp_r8 )
634      !                               
635      CALL iom_rstput( 0, 0, inum, 'e1t'  , e1t  , ktype = jp_r8 )   ! i-scale factors (e1.)
636      CALL iom_rstput( 0, 0, inum, 'e1u'  , e1u  , ktype = jp_r8 )
637      CALL iom_rstput( 0, 0, inum, 'e1v'  , e1v  , ktype = jp_r8 )
638      CALL iom_rstput( 0, 0, inum, 'e1f'  , e1f  , ktype = jp_r8 )
639      !
640      CALL iom_rstput( 0, 0, inum, 'e2t'  , e2t  , ktype = jp_r8 )   ! j-scale factors (e2.)
641      CALL iom_rstput( 0, 0, inum, 'e2u'  , e2u  , ktype = jp_r8 )
642      CALL iom_rstput( 0, 0, inum, 'e2v'  , e2v  , ktype = jp_r8 )
643      CALL iom_rstput( 0, 0, inum, 'e2f'  , e2f  , ktype = jp_r8 )
644      !
645      CALL iom_rstput( 0, 0, inum, 'ff_f' , ff_f , ktype = jp_r8 )   ! coriolis factor
646      CALL iom_rstput( 0, 0, inum, 'ff_t' , ff_t , ktype = jp_r8 )
647      !
648      !                             !==  vertical mesh  ==!
649      !                                                     
650      CALL iom_rstput( 0, 0, inum, 'e3t_1d'  , e3t_1d , ktype = jp_r8 )   ! reference 1D-coordinate
651      CALL iom_rstput( 0, 0, inum, 'e3w_1d'  , e3w_1d , ktype = jp_r8 )
652      !
653      CALL iom_rstput( 0, 0, inum, 'e3t_0'   , e3t_0  , ktype = jp_r8 )   ! vertical scale factors
654      CALL iom_rstput( 0, 0, inum, 'e3u_0'   , e3u_0  , ktype = jp_r8 )
655      CALL iom_rstput( 0, 0, inum, 'e3v_0'   , e3v_0  , ktype = jp_r8 )
656      CALL iom_rstput( 0, 0, inum, 'e3f_0'   , e3f_0  , ktype = jp_r8 )
657      CALL iom_rstput( 0, 0, inum, 'e3w_0'   , e3w_0  , ktype = jp_r8 )
658      CALL iom_rstput( 0, 0, inum, 'e3uw_0'  , e3uw_0 , ktype = jp_r8 )
659      CALL iom_rstput( 0, 0, inum, 'e3vw_0'  , e3vw_0 , ktype = jp_r8 )
660      !                                         
661      !                             !==  wet top and bottom level  ==!   (caution: multiplied by ssmask)
662      !
663      CALL iom_rstput( 0, 0, inum, 'top_level'    , REAL( mikt, wp )*ssmask , ktype = jp_i4 )   ! nb of ocean T-points (ISF)
664      CALL iom_rstput( 0, 0, inum, 'bottom_level' , REAL( mbkt, wp )*ssmask , ktype = jp_i4 )   ! nb of ocean T-points
665      !
666      IF( ln_sco ) THEN             ! s-coordinate: store grid stiffness ratio  (Not required anyway)
667         CALL dom_stiff( z2d )
668         CALL iom_rstput( 0, 0, inum, 'stiffness', z2d )        !    ! Max. grid stiffness ratio
669      ENDIF
670      !
671      IF( ll_wd ) THEN              ! wetting and drying domain
672         CALL iom_rstput( 0, 0, inum, 'ht_0'   , ht_0   , ktype = jp_r8 )
673      ENDIF
674      !
675      ! Add some global attributes ( netcdf only )
676      IF( iom_file(inum)%iolib == jpnf90 ) THEN
677         CALL iom_putatt( inum, 'nn_cfg', nn_cfg )
678         CALL iom_putatt( inum, 'cn_cfg', TRIM(cn_cfg) )
679      ENDIF
680      !
681      !                                ! ============================
682      !                                !        close the files
683      !                                ! ============================
684      CALL iom_close( inum )
685      !
686   END SUBROUTINE cfg_write
687
688   !!======================================================================
689END MODULE domain
Note: See TracBrowser for help on using the repository browser.