source: NEMO/trunk/src/OCE/DOM/domain.F90

Last change on this file was 15270, checked in by smasson, 5 weeks ago

trunk: forget some cleaning (remove dom_glo), #2724

  • Property svn:keywords set to Id
File size: 37.7 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   !!            4.1  !  2020-02  (G. Madec, S. Techene)  introduce ssh to h0 ratio
18   !!----------------------------------------------------------------------
19
20   !!----------------------------------------------------------------------
21   !!   dom_init      : initialize the space and time domain
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 domtile        ! tiling utilities
30#if defined key_qco
31   USE domqco         ! quasi-eulerian coord.
32#elif defined key_linssh
33   !                  ! fix in time coord.
34#else
35   USE domvvl         ! variable volume coord.
36#endif
37#if defined key_agrif
38   USE agrif_oce_interp, ONLY : Agrif_istate_ssh ! ssh interpolated from parent
39#endif
40   USE sbc_oce        ! surface boundary condition: ocean
41   USE trc_oce        ! shared ocean & passive tracers variab
42   USE phycst         ! physical constants
43   USE domhgr         ! domain: set the horizontal mesh
44   USE domzgr         ! domain: set the vertical mesh
45   USE dommsk         ! domain: set the mask system
46   USE domwri         ! domain: write the meshmask file
47   USE wet_dry , ONLY : ll_wd     ! wet & drying flag
48   USE closea  , ONLY : dom_clo   ! closed seas routine
49   USE c1d
50   !
51   USE in_out_manager ! I/O manager
52   USE iom            ! I/O library
53   USE lbclnk         ! ocean lateral boundary condition (or mpp link)
54   USE lib_mpp        ! distributed memory computing library
55   USE restart        ! only for lrst_oce and rst_read_ssh
56
57   IMPLICIT NONE
58   PRIVATE
59
60   PUBLIC   dom_init     ! called by nemogcm.F90
61   PUBLIC   domain_cfg   ! called by nemogcm.F90
62
63   !! * Substitutions
64#  include "do_loop_substitute.h90"
65   !!-------------------------------------------------------------------------
66   !! NEMO/OCE 4.0 , NEMO Consortium (2018)
67   !! $Id$
68   !! Software governed by the CeCILL license (see ./LICENSE)
69   !!-------------------------------------------------------------------------
70CONTAINS
71
72   SUBROUTINE dom_init( Kbb, Kmm, Kaa )
73      !!----------------------------------------------------------------------
74      !!                  ***  ROUTINE dom_init  ***
75      !!
76      !! ** Purpose :   Domain initialization. Call the routines that are
77      !!              required to create the arrays which define the space
78      !!              and time domain of the ocean model.
79      !!
80      !! ** Method  : - dom_msk: compute the masks from the bathymetry file
81      !!              - dom_hgr: compute or read the horizontal grid-point position
82      !!                         and scale factors, and the coriolis factor
83      !!              - dom_zgr: define the vertical coordinate and the bathymetry
84      !!              - dom_wri: create the meshmask file (ln_meshmask=T)
85      !!              - 1D configuration, move Coriolis, u and v at T-point
86      !!----------------------------------------------------------------------
87      INTEGER          , INTENT(in) :: Kbb, Kmm, Kaa          ! ocean time level indices
88      !
89      INTEGER ::   ji, jj, jk, jt   ! dummy loop indices
90      INTEGER ::   iconf = 0    ! local integers
91      REAL(wp)::   zrdt
92      CHARACTER (len=64) ::   cform = "(A12, 3(A13, I7))"
93      INTEGER , DIMENSION(jpi,jpj) ::   ik_top , ik_bot       ! top and bottom ocean level
94      REAL(wp), DIMENSION(jpi,jpj) ::   z1_hu_0, z1_hv_0
95      !!----------------------------------------------------------------------
96      !
97      IF(lwp) THEN         ! Ocean domain Parameters (control print)
98         WRITE(numout,*)
99         WRITE(numout,*) 'dom_init : domain initialization'
100         WRITE(numout,*) '~~~~~~~~'
101         !
102         WRITE(numout,*)     '   Domain info'
103         WRITE(numout,*)     '      dimension of model:'
104         WRITE(numout,*)     '             Local domain      Global domain       Data domain '
105         WRITE(numout,cform) '        ','   jpi     : ', jpi, '   jpiglo  : ', jpiglo
106         WRITE(numout,cform) '        ','   jpj     : ', jpj, '   jpjglo  : ', jpjglo
107         WRITE(numout,cform) '        ','   jpk     : ', jpk, '   jpkglo  : ', jpkglo
108         WRITE(numout,cform) '       ' ,'   jpij    : ', jpij
109         WRITE(numout,*)     '      mpp local domain info (mpp):'
110         WRITE(numout,*)     '              jpni    : ', jpni, '   nn_hls  : ', nn_hls
111         WRITE(numout,*)     '              jpnj    : ', jpnj, '   nn_hls  : ', nn_hls
112         WRITE(numout,*)     '              jpnij   : ', jpnij
113         WRITE(numout,*)     '      lateral boundary of the Global domain:'
114         WRITE(numout,*)     '              cyclic east-west             :', l_Iperio
115         WRITE(numout,*)     '              cyclic north-south           :', l_Jperio
116         WRITE(numout,*)     '              North Pole folding           :', l_NFold
117         WRITE(numout,*)     '                 type of North pole Folding:', c_NFtype
118         WRITE(numout,*)     '      Ocean model configuration used:'
119         WRITE(numout,*)     '              cn_cfg = ', TRIM( cn_cfg ), '   nn_cfg = ', nn_cfg
120      ENDIF
121
122      !
123      !           !==  Reference coordinate system  ==!
124      !
125      CALL dom_nam                      ! read namelist ( namrun, namdom )
126      CALL dom_tile_init                ! Tile domain
127
128      IF( ln_c1d )   CALL     c1d_init                    ! 1D column configuration
129      !
130      CALL dom_hgr                      ! Horizontal mesh
131
132      IF( ln_closea ) CALL dom_clo      ! Read in masks to define closed seas and lakes
133
134      CALL dom_zgr( ik_top, ik_bot )    ! Vertical mesh and bathymetry (return top and bottom ocean t-level indices)
135
136      CALL dom_msk( ik_top, ik_bot )    ! Masks
137      !
138      ht_0(:,:) = 0._wp  ! Reference ocean thickness
139      hu_0(:,:) = 0._wp
140      hv_0(:,:) = 0._wp
141      hf_0(:,:) = 0._wp
142      DO jk = 1, jpkm1
143         ht_0(:,:) = ht_0(:,:) + e3t_0(:,:,jk) * tmask(:,:,jk)
144         hu_0(:,:) = hu_0(:,:) + e3u_0(:,:,jk) * umask(:,:,jk)
145         hv_0(:,:) = hv_0(:,:) + e3v_0(:,:,jk) * vmask(:,:,jk)
146      END DO
147      !
148      DO_3D( 0, 0, 0, 0, 1, jpkm1 )
149         hf_0(ji,jj) = hf_0(ji,jj) + e3f_0(ji,jj,jk)*vmask(ji,jj,jk)*vmask(ji+1,jj,jk)
150      END_3D
151      CALL lbc_lnk('domain', hf_0, 'F', 1._wp)
152      !
153      IF( lk_SWE ) THEN      ! SWE case redefine hf_0
154         hf_0(:,:) = hf_0(:,:) + e3f_0(:,:,1) * ssfmask(:,:)
155      ENDIF
156      !
157      r1_ht_0(:,:) = ssmask (:,:) / ( ht_0(:,:) + 1._wp -  ssmask (:,:) )
158      r1_hu_0(:,:) = ssumask(:,:) / ( hu_0(:,:) + 1._wp -  ssumask(:,:) )
159      r1_hv_0(:,:) = ssvmask(:,:) / ( hv_0(:,:) + 1._wp -  ssvmask(:,:) )
160      r1_hf_0(:,:) = ssfmask(:,:) / ( hf_0(:,:) + 1._wp -  ssfmask(:,:) )
161      !
162      IF( ll_wd ) THEN       ! wet and drying (check ht_0 >= 0)
163         DO_2D( 1, 1, 1, 1 )
164            IF( ht_0(ji,jj) < 0._wp .AND. ssmask(ji,jj) == 1._wp ) THEN
165               CALL ctl_stop( 'dom_init : ht_0 must be positive at potentially wet points' )
166            ENDIF
167         END_2D
168      ENDIF
169      !
170      !           !==  initialisation of time varying coordinate  ==!
171      !
172      !                                 != ssh initialization
173      !
174      IF( l_SAS ) THEN        !* No ocean dynamics calculation : set to 0
175         ssh(:,:,:) = 0._wp
176#if defined key_agrif
177      ELSEIF( .NOT.Agrif_root() .AND.    &
178         &     ln_init_chfrpar ) THEN        !* Interpolate initial ssh from parent
179         CALL Agrif_istate_ssh( Kbb, Kmm, Kaa )
180#endif
181      ELSE                                   !* Read in restart file or set by user
182         CALL rst_read_ssh( Kbb, Kmm, Kaa )
183      ENDIF
184      !     
185#if defined key_qco
186      !                                 != Quasi-Euerian coordinate case
187      !
188      IF( .NOT.l_offline )   CALL dom_qco_init( Kbb, Kmm, Kaa )
189#elif defined key_linssh
190      !                                 != Fix in time : key_linssh case, set through domzgr_substitute.h90
191#else
192      !
193      IF( ln_linssh ) THEN              != Fix in time : set to the reference one for all
194         !
195         DO jt = 1, jpt                         ! depth of t- and w-grid-points
196            gdept(:,:,:,jt) = gdept_0(:,:,:)
197            gdepw(:,:,:,jt) = gdepw_0(:,:,:)
198         END DO
199            gde3w(:,:,:)    = gde3w_0(:,:,:)    ! = gdept as the sum of e3t
200         !
201         DO jt = 1, jpt                         ! vertical scale factors
202            e3t (:,:,:,jt) =  e3t_0(:,:,:)
203            e3u (:,:,:,jt) =  e3u_0(:,:,:)
204            e3v (:,:,:,jt) =  e3v_0(:,:,:)
205            e3w (:,:,:,jt) =  e3w_0(:,:,:)
206            e3uw(:,:,:,jt) = e3uw_0(:,:,:)
207            e3vw(:,:,:,jt) = e3vw_0(:,:,:)
208         END DO
209            e3f (:,:,:)    =  e3f_0(:,:,:)
210         !
211         DO jt = 1, jpt                         ! water column thickness and its inverse
212               hu(:,:,jt) =    hu_0(:,:)
213               hv(:,:,jt) =    hv_0(:,:)
214            r1_hu(:,:,jt) = r1_hu_0(:,:)
215            r1_hv(:,:,jt) = r1_hv_0(:,:)
216         END DO
217               ht   (:,:) =    ht_0(:,:)
218         !
219      ELSE                              != Time varying : initialize before/now/after variables
220         !
221         IF( .NOT.l_offline )   CALL dom_vvl_init( Kbb, Kmm, Kaa )
222         !
223      ENDIF
224#endif
225
226      !
227
228#if defined key_agrif
229      IF( .NOT. Agrif_Root() ) CALL Agrif_Init_Domain( Kbb, Kmm, Kaa )
230#endif
231      IF( ln_meshmask    )   CALL dom_wri       ! Create a domain file
232      IF( .NOT.ln_rstart )   CALL dom_ctl       ! Domain control
233      !
234      IF( ln_write_cfg   )   CALL cfg_write     ! create the configuration file
235      !
236      IF(lwp) THEN
237         WRITE(numout,*)
238         WRITE(numout,*) 'dom_init :   ==>>>   END of domain initialization'
239         WRITE(numout,*) '~~~~~~~~'
240         WRITE(numout,*)
241      ENDIF
242      !
243   END SUBROUTINE dom_init
244
245
246   SUBROUTINE dom_nam
247      !!----------------------------------------------------------------------
248      !!                     ***  ROUTINE dom_nam  ***
249      !!
250      !! ** Purpose :   read domaine namelists and print the variables.
251      !!
252      !! ** input   : - namrun namelist
253      !!              - namdom namelist
254      !!              - namtile namelist
255      !!              - namnc4 namelist   ! "key_netcdf4" only
256      !!----------------------------------------------------------------------
257      USE ioipsl
258      !!
259      INTEGER ::   ios   ! Local integer
260      REAL(wp)::   zrdt
261      !!----------------------------------------------------------------------
262      !
263      NAMELIST/namrun/ cn_ocerst_indir, cn_ocerst_outdir, nn_stocklist, ln_rst_list,                 &
264         &             nn_no   , cn_exp   , cn_ocerst_in, cn_ocerst_out, ln_rstart , nn_rstctl ,     &
265         &             nn_it000, nn_itend , nn_date0    , nn_time0     , nn_leapy  , nn_istate ,     &
266         &             nn_stock, nn_write , ln_mskland  , ln_clobber   , nn_chunksz, ln_1st_euler  , &
267         &             ln_cfmeta, ln_xios_read, nn_wxios
268      NAMELIST/namdom/ ln_linssh, rn_Dt, rn_atfp, ln_crs, ln_c1d, ln_meshmask
269      NAMELIST/namtile/ ln_tile, nn_ltile_i, nn_ltile_j
270#if defined key_netcdf4
271      NAMELIST/namnc4/ nn_nchunks_i, nn_nchunks_j, nn_nchunks_k, ln_nc4zip
272#endif
273      !!----------------------------------------------------------------------
274      !
275      IF(lwp) THEN
276         WRITE(numout,*)
277         WRITE(numout,*) 'dom_nam : domain initialization through namelist read'
278         WRITE(numout,*) '~~~~~~~ '
279      ENDIF
280      !
281      !                       !=======================!
282      !                       !==  namelist namdom  ==!
283      !                       !=======================!
284      !
285      READ  ( numnam_ref, namdom, IOSTAT = ios, ERR = 903)
286903   IF( ios /= 0 )   CALL ctl_nam ( ios , 'namdom in reference namelist' )
287      READ  ( numnam_cfg, namdom, IOSTAT = ios, ERR = 904 )
288904   IF( ios >  0 )   CALL ctl_nam ( ios , 'namdom in configuration namelist' )
289      IF(lwm) WRITE( numond, namdom )
290      !
291#if defined key_linssh
292      ln_linssh = lk_linssh      ! overwrite ln_linssh with the logical associated with key_linssh
293#endif
294      !
295#if defined key_agrif
296      IF( .NOT. Agrif_Root() ) THEN    ! AGRIF child, subdivide the Parent timestep
297         rn_Dt = Agrif_Parent (rn_Dt ) / Agrif_Rhot()
298      ENDIF
299#endif
300      !
301      IF(lwp) THEN
302         WRITE(numout,*)
303         WRITE(numout,*) '   Namelist : namdom   ---   space & time domain'
304         WRITE(numout,*) '      linear free surface (=T)                ln_linssh   = ', ln_linssh
305         WRITE(numout,*) '      create mesh/mask file                   ln_meshmask = ', ln_meshmask
306         WRITE(numout,*) '      ocean time step                         rn_Dt       = ', rn_Dt
307         WRITE(numout,*) '      asselin time filter parameter           rn_atfp     = ', rn_atfp
308         WRITE(numout,*) '      online coarsening of dynamical fields   ln_crs      = ', ln_crs
309         WRITE(numout,*) '      single column domain (1x1pt)            ln_c1d      = ', ln_c1d
310      ENDIF
311      !
312      ! set current model timestep rDt = 2*rn_Dt if MLF or rDt = rn_Dt if RK3
313      rDt   = 2._wp * rn_Dt
314      r1_Dt = 1._wp / rDt
315      !
316      IF( l_SAS .AND. .NOT.ln_linssh ) THEN
317         CALL ctl_warn( 'SAS requires linear ssh : force ln_linssh = T' )
318         ln_linssh = .TRUE.
319      ENDIF
320      !
321#if defined key_qco
322      IF( ln_linssh )   CALL ctl_stop( 'STOP','domain: key_qco and ln_linssh=T or key_linssh are incompatible' )
323#endif
324      !
325      !                       !=======================!
326      !                       !==  namelist namrun  ==!
327      !                       !=======================!
328      !
329      READ  ( numnam_ref, namrun, IOSTAT = ios, ERR = 901)
330901   IF( ios /= 0 )   CALL ctl_nam ( ios , 'namrun in reference namelist' )
331      READ  ( numnam_cfg, namrun, IOSTAT = ios, ERR = 902 )
332902   IF( ios >  0 )   CALL ctl_nam ( ios , 'namrun in configuration namelist' )
333      IF(lwm) WRITE ( numond, namrun )
334
335#if defined key_agrif
336      IF( .NOT. Agrif_Root() ) THEN
337            nn_it000 = (Agrif_Parent(nn_it000)-1)*Agrif_IRhot() + 1
338            nn_itend =  Agrif_Parent(nn_itend)   *Agrif_IRhot()
339      ENDIF
340#endif
341      !
342      IF(lwp) THEN                  ! control print
343         WRITE(numout,*) '   Namelist : namrun   ---   run parameters'
344         WRITE(numout,*) '      Assimilation cycle              nn_no           = ', nn_no
345         WRITE(numout,*) '      experiment name for output      cn_exp          = ', TRIM( cn_exp           )
346         WRITE(numout,*) '      file prefix restart input       cn_ocerst_in    = ', TRIM( cn_ocerst_in     )
347         WRITE(numout,*) '      restart input directory         cn_ocerst_indir = ', TRIM( cn_ocerst_indir  )
348         WRITE(numout,*) '      file prefix restart output      cn_ocerst_out   = ', TRIM( cn_ocerst_out    )
349         WRITE(numout,*) '      restart output directory        cn_ocerst_outdir= ', TRIM( cn_ocerst_outdir )
350         WRITE(numout,*) '      restart logical                 ln_rstart       = ', ln_rstart
351         WRITE(numout,*) '      start with forward time step    ln_1st_euler    = ', ln_1st_euler
352         WRITE(numout,*) '      control of time step            nn_rstctl       = ', nn_rstctl
353         WRITE(numout,*) '      number of the first time step   nn_it000        = ', nn_it000
354         WRITE(numout,*) '      number of the last time step    nn_itend        = ', nn_itend
355         WRITE(numout,*) '      initial calendar date aammjj    nn_date0        = ', nn_date0
356         WRITE(numout,*) '      initial time of day in hhmm     nn_time0        = ', nn_time0
357         WRITE(numout,*) '      leap year calendar (0/1)        nn_leapy        = ', nn_leapy
358         WRITE(numout,*) '      initial state output            nn_istate       = ', nn_istate
359         IF( ln_rst_list ) THEN
360            WRITE(numout,*) '      list of restart dump times      nn_stocklist    =', nn_stocklist
361         ELSE
362            WRITE(numout,*) '      frequency of restart file       nn_stock        = ', nn_stock
363         ENDIF
364#if ! defined key_xios
365         WRITE(numout,*) '      frequency of output file        nn_write        = ', nn_write
366#endif
367         WRITE(numout,*) '      mask land points                ln_mskland      = ', ln_mskland
368         WRITE(numout,*) '      additional CF standard metadata ln_cfmeta       = ', ln_cfmeta
369         WRITE(numout,*) '      overwrite an existing file      ln_clobber      = ', ln_clobber
370         WRITE(numout,*) '      NetCDF chunksize (bytes)        nn_chunksz      = ', nn_chunksz
371         IF( TRIM(Agrif_CFixed()) == '0' ) THEN
372            WRITE(numout,*) '      READ restart for a single file using XIOS ln_xios_read =', ln_xios_read
373            WRITE(numout,*) '      Write restart using XIOS        nn_wxios   = ', nn_wxios
374         ELSE
375            WRITE(numout,*) "      AGRIF: nn_wxios will be ingored. See setting for parent"
376            WRITE(numout,*) "      AGRIF: ln_xios_read will be ingored. See setting for parent"
377         ENDIF
378      ENDIF
379
380      cexper = cn_exp         ! conversion DOCTOR names into model names (this should disappear soon)
381      nrstdt = nn_rstctl
382      nit000 = nn_it000
383      nitend = nn_itend
384      ndate0 = nn_date0
385      nleapy = nn_leapy
386      ninist = nn_istate
387      !
388      !                                        !==  Set parameters for restart reading using xIOS  ==!
389      !
390      IF( TRIM(Agrif_CFixed()) == '0' ) THEN
391         lrxios = ln_xios_read .AND. ln_rstart
392         IF( nn_wxios > 0 )   lwxios = .TRUE.           !* set output file type for XIOS based on NEMO namelist
393         nxioso = nn_wxios
394      ENDIF
395      !                                        !==  Check consistency between ln_rstart and ln_1st_euler  ==!   (i.e. set l_1st_euler)
396      l_1st_euler = ln_1st_euler
397      !
398      IF( ln_rstart ) THEN                              !*  Restart case
399         !
400         IF(lwp) WRITE(numout,*)
401         IF(lwp) WRITE(numout,*) '   open the restart file'
402         CALL rst_read_open                                              !- Open the restart file
403         !
404         IF( iom_varid( numror, 'rdt', ldstop = .FALSE. ) > 0 ) THEN     !- Check time-step consistency and force Euler restart if changed
405            CALL iom_get( numror, 'rdt', zrdt )
406            IF( zrdt /= rn_Dt ) THEN
407               IF(lwp) WRITE( numout,*)
408               IF(lwp) WRITE( numout,*) '   rn_Dt = ', rn_Dt,' not equal to the READ one rdt = ', zrdt
409               IF(lwp) WRITE( numout,*)
410               IF(lwp) WRITE( numout,*) '      ==>>>   forced euler first time-step'
411               l_1st_euler =  .TRUE.
412            ENDIF
413         ENDIF
414         !
415         IF( .NOT.l_SAS .AND. iom_varid( numror, 'sshb', ldstop = .FALSE. ) <= 0 ) THEN   !- Check absence of one of the Kbb field (here sshb)
416            !                                                                             !  (any Kbb field is missing ==> all Kbb fields are missing)
417            IF( .NOT.l_1st_euler ) THEN
418               CALL ctl_warn('dom_nam : ssh at Kbb not found in restart files ',   &
419                  &                        'l_1st_euler forced to .true. and ' ,   &
420                  &                        'ssh(Kbb) = ssh(Kmm) '                  )
421               l_1st_euler = .TRUE.
422            ENDIF
423         ENDIF
424      ELSEIF( .NOT.l_1st_euler ) THEN                   !*  Initialization case
425         IF(lwp) WRITE(numout,*)
426         IF(lwp) WRITE(numout,*)'   ==>>>   Start from rest (ln_rstart=F)'
427         IF(lwp) WRITE(numout,*)'           an Euler initial time step is used : l_1st_euler is forced to .true. '
428         l_1st_euler = .TRUE.
429      ENDIF
430      !
431      !                                        !==  control of output frequency  ==!
432      !
433      IF( .NOT. ln_rst_list ) THEN   ! we use nn_stock
434         IF( nn_stock == -1 )   CALL ctl_warn( 'nn_stock = -1 --> no restart will be done' )
435         IF( nn_stock == 0 .OR. nn_stock > nitend ) THEN
436            WRITE(ctmp1,*) 'nn_stock = ', nn_stock, ' it is forced to ', nitend
437            CALL ctl_warn( ctmp1 )
438            nn_stock = nitend
439         ENDIF
440      ENDIF
441#if ! defined key_xios
442      IF( nn_write == -1 )   CALL ctl_warn( 'nn_write = -1 --> no output files will be done' )
443      IF ( nn_write == 0 ) THEN
444         WRITE(ctmp1,*) 'nn_write = ', nn_write, ' it is forced to ', nitend
445         CALL ctl_warn( ctmp1 )
446         nn_write = nitend
447      ENDIF
448#endif
449
450      IF( Agrif_Root() ) THEN
451         IF(lwp) WRITE(numout,*)
452         SELECT CASE ( nleapy )                !==  Choose calendar for IOIPSL  ==!
453         CASE (  1 )
454            CALL ioconf_calendar('gregorian')
455            IF(lwp) WRITE(numout,*) '   ==>>>   The IOIPSL calendar is "gregorian", i.e. leap year'
456         CASE (  0 )
457            CALL ioconf_calendar('noleap')
458            IF(lwp) WRITE(numout,*) '   ==>>>   The IOIPSL calendar is "noleap", i.e. no leap year'
459         CASE ( 30 )
460            CALL ioconf_calendar('360d')
461            IF(lwp) WRITE(numout,*) '   ==>>>   The IOIPSL calendar is "360d", i.e. 360 days in a year'
462         END SELECT
463      ENDIF
464      !
465      !                       !========================!
466      !                       !==  namelist namtile  ==!
467      !                       !========================!
468      !
469      READ  ( numnam_ref, namtile, IOSTAT = ios, ERR = 905 )
470905   IF( ios /= 0 )   CALL ctl_nam ( ios , 'namtile in reference namelist' )
471      READ  ( numnam_cfg, namtile, IOSTAT = ios, ERR = 906 )
472906   IF( ios >  0 )   CALL ctl_nam ( ios , 'namtile in configuration namelist' )
473      IF(lwm) WRITE( numond, namtile )
474
475      IF(lwp) THEN
476         WRITE(numout,*)
477         WRITE(numout,*)    '   Namelist : namtile   ---   Domain tiling decomposition'
478         WRITE(numout,*)    '      Tiling (T) or not (F)                ln_tile    = ', ln_tile
479         WRITE(numout,*)    '      Length of tile in i                  nn_ltile_i = ', nn_ltile_i
480         WRITE(numout,*)    '      Length of tile in j                  nn_ltile_j = ', nn_ltile_j
481         WRITE(numout,*)
482         IF( ln_tile ) THEN
483            WRITE(numout,*) '      The domain will be decomposed into tiles of size', nn_ltile_i, 'x', nn_ltile_j
484         ELSE
485            WRITE(numout,*) '      Domain tiling will NOT be used'
486         ENDIF
487      ENDIF
488      !
489#if defined key_netcdf4
490      !                       !=======================!
491      !                       !==  namelist namnc4  ==!   NetCDF 4 case   ("key_netcdf4" defined)
492      !                       !=======================!
493      !
494      READ  ( numnam_ref, namnc4, IOSTAT = ios, ERR = 907)
495907   IF( ios /= 0 )   CALL ctl_nam ( ios , 'namnc4 in reference namelist' )
496      READ  ( numnam_cfg, namnc4, IOSTAT = ios, ERR = 908 )
497908   IF( ios >  0 )   CALL ctl_nam ( ios , 'namnc4 in configuration namelist' )
498      IF(lwm) WRITE( numond, namnc4 )
499
500      IF(lwp) THEN                        ! control print
501         WRITE(numout,*)
502         WRITE(numout,*) '   Namelist namnc4 - Netcdf4 chunking parameters ("key_netcdf4" defined)'
503         WRITE(numout,*) '      number of chunks in i-dimension             nn_nchunks_i = ', nn_nchunks_i
504         WRITE(numout,*) '      number of chunks in j-dimension             nn_nchunks_j = ', nn_nchunks_j
505         WRITE(numout,*) '      number of chunks in k-dimension             nn_nchunks_k = ', nn_nchunks_k
506         WRITE(numout,*) '      apply netcdf4/hdf5 chunking & compression   ln_nc4zip    = ', ln_nc4zip
507      ENDIF
508
509      ! Put the netcdf4 settings into a simple structure (snc4set, defined in in_out_manager module)
510      ! Note the chunk size in the unlimited (time) dimension will be fixed at 1
511      snc4set%ni   = nn_nchunks_i
512      snc4set%nj   = nn_nchunks_j
513      snc4set%nk   = nn_nchunks_k
514      snc4set%luse = ln_nc4zip
515#else
516      snc4set%luse = .FALSE.        ! No NetCDF 4 case
517#endif
518      !
519   END SUBROUTINE dom_nam
520
521
522   SUBROUTINE dom_ctl
523      !!----------------------------------------------------------------------
524      !!                     ***  ROUTINE dom_ctl  ***
525      !!
526      !! ** Purpose :   Domain control.
527      !!
528      !! ** Method  :   compute and print extrema of masked scale factors
529      !!----------------------------------------------------------------------
530      LOGICAL, DIMENSION(jpi,jpj) ::   llmsk
531      INTEGER, DIMENSION(2)       ::   imil, imip, imi1, imi2, imal, imap, ima1, ima2
532      REAL(wp)                    ::   zglmin, zglmax, zgpmin, zgpmax, ze1min, ze1max, ze2min, ze2max
533      !!----------------------------------------------------------------------
534      !
535      llmsk = tmask_h(:,:) == 1._wp
536      !
537      CALL mpp_minloc( 'domain', glamt(:,:), llmsk, zglmin, imil )
538      CALL mpp_minloc( 'domain', gphit(:,:), llmsk, zgpmin, imip )
539      CALL mpp_minloc( 'domain',   e1t(:,:), llmsk, ze1min, imi1 )
540      CALL mpp_minloc( 'domain',   e2t(:,:), llmsk, ze2min, imi2 )
541      CALL mpp_maxloc( 'domain', glamt(:,:), llmsk, zglmax, imal )
542      CALL mpp_maxloc( 'domain', gphit(:,:), llmsk, zgpmax, imap )
543      CALL mpp_maxloc( 'domain',   e1t(:,:), llmsk, ze1max, ima1 )
544      CALL mpp_maxloc( 'domain',   e2t(:,:), llmsk, ze2max, ima2 )
545      !
546      IF(lwp) THEN
547         WRITE(numout,*)
548         WRITE(numout,*) 'dom_ctl : extrema of the masked scale factors'
549         WRITE(numout,*) '~~~~~~~'
550         WRITE(numout,"(14x,'glamt mini: ',1f10.2,' at i = ',i5,' j= ',i5)") zglmin, imil(1), imil(2)
551         WRITE(numout,"(14x,'glamt maxi: ',1f10.2,' at i = ',i5,' j= ',i5)") zglmax, imal(1), imal(2)
552         WRITE(numout,"(14x,'gphit mini: ',1f10.2,' at i = ',i5,' j= ',i5)") zgpmin, imip(1), imip(2)
553         WRITE(numout,"(14x,'gphit maxi: ',1f10.2,' at i = ',i5,' j= ',i5)") zgpmax, imap(1), imap(2)
554         WRITE(numout,"(14x,'  e1t mini: ',1f10.2,' at i = ',i5,' j= ',i5)") ze1min, imi1(1), imi1(2)
555         WRITE(numout,"(14x,'  e1t maxi: ',1f10.2,' at i = ',i5,' j= ',i5)") ze1max, ima1(1), ima1(2)
556         WRITE(numout,"(14x,'  e2t mini: ',1f10.2,' at i = ',i5,' j= ',i5)") ze2min, imi2(1), imi2(2)
557         WRITE(numout,"(14x,'  e2t maxi: ',1f10.2,' at i = ',i5,' j= ',i5)") ze2max, ima2(1), ima2(2)
558      ENDIF
559      !
560   END SUBROUTINE dom_ctl
561
562
563   SUBROUTINE domain_cfg( cd_cfg, kk_cfg, kpi, kpj, kpk, ldIperio, ldJperio, ldNFold, cdNFtype )
564      !!----------------------------------------------------------------------
565      !!                     ***  ROUTINE domain_cfg  ***
566      !!
567      !! ** Purpose :   read the domain size in domain configuration file
568      !!
569      !! ** Method  :   read the cn_domcfg NetCDF file
570      !!----------------------------------------------------------------------
571      CHARACTER(len=*), INTENT(out) ::   cd_cfg               ! configuration name
572      INTEGER         , INTENT(out) ::   kk_cfg               ! configuration resolution
573      INTEGER         , INTENT(out) ::   kpi, kpj, kpk        ! global domain sizes
574      LOGICAL         , INTENT(out) ::   ldIperio, ldJperio   ! i- and j- periodicity
575      LOGICAL         , INTENT(out) ::   ldNFold              ! North pole folding
576      CHARACTER(len=1), INTENT(out) ::   cdNFtype             ! Folding type: T or F
577      !
578      CHARACTER(len=7) ::   catt                  ! 'T', 'F', '-' or 'UNKNOWN'
579      INTEGER ::   inum, iperio, iatt             ! local integer
580      REAL(wp) ::   zorca_res                     ! local scalars
581      REAL(wp) ::   zperio                        !   -      -
582      INTEGER, DIMENSION(4) ::   idvar, idimsz    ! size   of dimensions
583      !!----------------------------------------------------------------------
584      !
585      IF(lwp) THEN
586         WRITE(numout,*) '           '
587         WRITE(numout,*) 'domain_cfg : domain size read in ', TRIM( cn_domcfg ), ' file'
588         WRITE(numout,*) '~~~~~~~~~~ '
589      ENDIF
590      !
591      CALL iom_open( cn_domcfg, inum )
592      !
593      CALL iom_getatt( inum,  'CfgName', cd_cfg )   ! returns 'UNKNOWN' if not found
594      CALL iom_getatt( inum, 'CfgIndex', kk_cfg )   ! returns      -999 if not found
595      !
596      ! ------- keep compatibility with OLD VERSION... start -------
597      IF( cd_cfg == 'UNKNOWN' .AND. kk_cfg == -999 ) THEN
598         IF(  iom_varid( inum, 'ORCA'       , ldstop = .FALSE. ) > 0  .AND.  &
599            & iom_varid( inum, 'ORCA_index' , ldstop = .FALSE. ) > 0    ) THEN
600            !
601            cd_cfg = 'ORCA'
602            CALL iom_get( inum, 'ORCA_index', zorca_res )   ;   kk_cfg = NINT( zorca_res )
603            !
604         ELSE
605            CALL iom_getatt( inum, 'cn_cfg', cd_cfg )  ! returns 'UNKNOWN' if not found
606            CALL iom_getatt( inum, 'nn_cfg', kk_cfg )  ! returns      -999 if not found
607         ENDIF
608      ENDIF
609      ! ------- keep compatibility with OLD VERSION... end -------
610      !
611      idvar = iom_varid( inum, 'e3t_0', kdimsz = idimsz )   ! use e3t_0, that must exist, to get jp(ijk)glo
612      kpi = idimsz(1)
613      kpj = idimsz(2)
614      kpk = idimsz(3)
615      !
616      CALL iom_getatt( inum, 'Iperio', iatt )   ;   ldIperio = iatt == 1   ! returns      -999 if not found -> default = .false.
617      CALL iom_getatt( inum, 'Jperio', iatt )   ;   ldJperio = iatt == 1   ! returns      -999 if not found -> default = .false.
618      CALL iom_getatt( inum,  'NFold', iatt )   ;   ldNFold  = iatt == 1   ! returns      -999 if not found -> default = .false.
619      CALL iom_getatt( inum, 'NFtype', catt )                              ! returns 'UNKNOWN' if not found
620      IF( LEN_TRIM(catt) == 1 ) THEN   ;   cdNFtype = TRIM(catt)
621      ELSE                             ;   cdNFtype = '-'
622      ENDIF
623      !
624      ! ------- keep compatibility with OLD VERSION... start -------
625      IF( iatt == -999 .AND. catt == 'UNKNOWN' .AND. iom_varid( inum, 'jperio', ldstop = .FALSE. ) > 0 ) THEN                   
626         CALL iom_get( inum, 'jperio', zperio )   ;   iperio = NINT( zperio )
627         ldIperio = iperio == 1  .OR. iperio == 4 .OR. iperio == 6 .OR. iperio == 7   ! i-periodicity
628         ldJperio = iperio == 2  .OR. iperio == 7                                     ! j-periodicity
629         ldNFold  = iperio >= 3 .AND. iperio <= 6                                     ! North pole folding
630         IF(     iperio == 3 .OR. iperio == 4 ) THEN   ;   cdNFtype = 'T'             !    folding at T point
631         ELSEIF( iperio == 5 .OR. iperio == 6 ) THEN   ;   cdNFtype = 'F'             !    folding at F point
632         ELSE                                          ;   cdNFtype = '-'             !    default value
633         ENDIF
634      ENDIF
635      ! ------- keep compatibility with OLD VERSION... end -------
636      !
637      CALL iom_close( inum )
638      !
639      IF(lwp) THEN
640         WRITE(numout,*) '   .'
641         WRITE(numout,*) '   ==>>>   ', TRIM(cn_cfg), ' configuration '
642         WRITE(numout,*) '   .'
643         WRITE(numout,*) '      nn_cfg = ', kk_cfg
644         WRITE(numout,*) '      Ni0glo = ', kpi
645         WRITE(numout,*) '      Nj0glo = ', kpj
646         WRITE(numout,*) '      jpkglo = ', kpk
647      ENDIF
648      !
649   END SUBROUTINE domain_cfg
650
651
652   SUBROUTINE cfg_write
653      !!----------------------------------------------------------------------
654      !!                  ***  ROUTINE cfg_write  ***
655      !!
656      !! ** Purpose :   Create the "cn_domcfg_out" file, a NetCDF file which
657      !!              contains all the ocean domain informations required to
658      !!              define an ocean configuration.
659      !!
660      !! ** Method  :   Write in a file all the arrays required to set up an
661      !!              ocean configuration.
662      !!
663      !! ** output file :   domcfg_out.nc : domain size, characteristics, horizontal
664      !!                       mesh, Coriolis parameter, and vertical scale factors
665      !!                    NB: also contain ORCA family information
666      !!----------------------------------------------------------------------
667      INTEGER           ::   ji, jj, jk   ! dummy loop indices
668      INTEGER           ::   inum     ! local units
669      CHARACTER(len=21) ::   clnam    ! filename (mesh and mask informations)
670      REAL(wp), DIMENSION(jpi,jpj) ::   z2d   ! workspace
671      !!----------------------------------------------------------------------
672      !
673      IF(lwp) WRITE(numout,*)
674      IF(lwp) WRITE(numout,*) 'cfg_write : create the domain configuration file (', TRIM(cn_domcfg_out),'.nc)'
675      IF(lwp) WRITE(numout,*) '~~~~~~~~~'
676      !
677      !                       ! ============================= !
678      !                       !  create 'domcfg_out.nc' file  !
679      !                       ! ============================= !
680      !
681      clnam = cn_domcfg_out  ! filename (configuration information)
682      CALL iom_open( TRIM(clnam), inum, ldwrt = .TRUE. )
683      !
684      !                             !==  Configuration specificities  ==!
685      !
686      CALL iom_putatt( inum,  'CfgName', TRIM(cn_cfg) )
687      CALL iom_putatt( inum, 'CfgIndex',      nn_cfg  )
688      !
689      !                             !==  domain characteristics  ==!
690      !
691      !                                   ! lateral boundary of the global domain
692      CALL iom_putatt( inum, 'Iperio', COUNT( (/l_Iperio/) ) )
693      CALL iom_putatt( inum, 'Jperio', COUNT( (/l_Jperio/) ) )
694      CALL iom_putatt( inum,  'NFold', COUNT( (/l_NFold /) ) )
695      CALL iom_putatt( inum, 'NFtype',          c_NFtype     )
696
697      !                                   ! type of vertical coordinate
698      IF(ln_zco)   CALL iom_putatt( inum, 'VertCoord', 'zco' )
699      IF(ln_zps)   CALL iom_putatt( inum, 'VertCoord', 'zps' )
700      IF(ln_sco)   CALL iom_putatt( inum, 'VertCoord', 'sco' )
701     
702      !                                   ! ocean cavities under iceshelves
703      CALL iom_putatt( inum, 'IsfCav', COUNT( (/ln_isfcav/) ) )
704      !
705      !                             !==  horizontal mesh  !
706      !
707      CALL iom_rstput( 0, 0, inum, 'glamt', glamt, ktype = jp_r8 )   ! latitude
708      CALL iom_rstput( 0, 0, inum, 'glamu', glamu, ktype = jp_r8 )
709      CALL iom_rstput( 0, 0, inum, 'glamv', glamv, ktype = jp_r8 )
710      CALL iom_rstput( 0, 0, inum, 'glamf', glamf, ktype = jp_r8 )
711      !
712      CALL iom_rstput( 0, 0, inum, 'gphit', gphit, ktype = jp_r8 )   ! longitude
713      CALL iom_rstput( 0, 0, inum, 'gphiu', gphiu, ktype = jp_r8 )
714      CALL iom_rstput( 0, 0, inum, 'gphiv', gphiv, ktype = jp_r8 )
715      CALL iom_rstput( 0, 0, inum, 'gphif', gphif, ktype = jp_r8 )
716      !
717      CALL iom_rstput( 0, 0, inum, 'e1t'  , e1t  , ktype = jp_r8 )   ! i-scale factors (e1.)
718      CALL iom_rstput( 0, 0, inum, 'e1u'  , e1u  , ktype = jp_r8 )
719      CALL iom_rstput( 0, 0, inum, 'e1v'  , e1v  , ktype = jp_r8 )
720      CALL iom_rstput( 0, 0, inum, 'e1f'  , e1f  , ktype = jp_r8 )
721      !
722      CALL iom_rstput( 0, 0, inum, 'e2t'  , e2t  , ktype = jp_r8 )   ! j-scale factors (e2.)
723      CALL iom_rstput( 0, 0, inum, 'e2u'  , e2u  , ktype = jp_r8 )
724      CALL iom_rstput( 0, 0, inum, 'e2v'  , e2v  , ktype = jp_r8 )
725      CALL iom_rstput( 0, 0, inum, 'e2f'  , e2f  , ktype = jp_r8 )
726      !
727      CALL iom_rstput( 0, 0, inum, 'ff_f' , ff_f , ktype = jp_r8 )   ! coriolis factor
728      CALL iom_rstput( 0, 0, inum, 'ff_t' , ff_t , ktype = jp_r8 )
729      !
730      !                             !==  vertical mesh  ==!
731      !
732      CALL iom_rstput( 0, 0, inum, 'e3t_1d'  , e3t_1d , ktype = jp_r8 )   ! reference 1D-coordinate
733      CALL iom_rstput( 0, 0, inum, 'e3w_1d'  , e3w_1d , ktype = jp_r8 )
734      !
735      CALL iom_rstput( 0, 0, inum, 'e3t_0'   , e3t_0  , ktype = jp_r8 )   ! vertical scale factors
736      CALL iom_rstput( 0, 0, inum, 'e3u_0'   , e3u_0  , ktype = jp_r8 )
737      CALL iom_rstput( 0, 0, inum, 'e3v_0'   , e3v_0  , ktype = jp_r8 )
738      CALL iom_rstput( 0, 0, inum, 'e3f_0'   , e3f_0  , ktype = jp_r8 )
739      CALL iom_rstput( 0, 0, inum, 'e3w_0'   , e3w_0  , ktype = jp_r8 )
740      CALL iom_rstput( 0, 0, inum, 'e3uw_0'  , e3uw_0 , ktype = jp_r8 )
741      CALL iom_rstput( 0, 0, inum, 'e3vw_0'  , e3vw_0 , ktype = jp_r8 )
742      !
743      !                             !==  wet top and bottom level  ==!   (caution: multiplied by ssmask)
744      !
745      CALL iom_rstput( 0, 0, inum, 'top_level'    , REAL( mikt, wp )*ssmask , ktype = jp_i4 )   ! nb of ocean T-points (ISF)
746      CALL iom_rstput( 0, 0, inum, 'bottom_level' , REAL( mbkt, wp )*ssmask , ktype = jp_i4 )   ! nb of ocean T-points
747      !
748      IF( ln_sco ) THEN             ! s-coordinate: store grid stiffness ratio  (Not required anyway)
749         CALL dom_stiff( z2d )
750         CALL iom_rstput( 0, 0, inum, 'stiffness', z2d )        !    ! Max. grid stiffness ratio
751      ENDIF
752      !
753      IF( ll_wd ) THEN              ! wetting and drying domain
754         CALL iom_rstput( 0, 0, inum, 'ht_0'   , ht_0   , ktype = jp_r8 )
755      ENDIF
756      !                       ! ============================ !
757      !                       !        close the files
758      !                       ! ============================ !
759      CALL iom_close( inum )
760      !
761   END SUBROUTINE cfg_write
762
763   !!======================================================================
764END MODULE domain
Note: See TracBrowser for help on using the repository browser.