source: branches/UKMO/dev_r5518_GO6_package/NEMOGCM/NEMO/OPA_SRC/nemogcm.F90

Last change on this file was 11101, checked in by frrh, 16 months ago

Merge changes from Met Office GMED ticket 450 to reduce unnecessary
text output from NEMO.
This output, which is typically not switchable, is rarely of interest
in normal (non-debugging) runs and simply redunantley consumes extra
file space.
Further, the presence of this text output has been shown to
significantly degrade performance of models which are run during
Met Office HPC RAID (disk) checks.
The new code introduces switches which are configurable via the
changes made in the associated Met Office MOCI ticket 399.

File size: 44.4 KB
Line 
1MODULE nemogcm
2   !!======================================================================
3   !!                       ***  MODULE nemogcm   ***
4   !! Ocean system   : NEMO GCM (ocean dynamics, on-line tracers, biochemistry and sea-ice)
5   !!======================================================================
6   !! History :  OPA  ! 1990-10  (C. Levy, G. Madec)  Original code
7   !!            7.0  ! 1991-11  (M. Imbard, C. Levy, G. Madec)
8   !!            7.1  ! 1993-03  (M. Imbard, C. Levy, G. Madec, O. Marti, M. Guyon, A. Lazar,
9   !!                             P. Delecluse, C. Perigaud, G. Caniaux, B. Colot, C. Maes) release 7.1
10   !!             -   ! 1992-06  (L.Terray)  coupling implementation
11   !!             -   ! 1993-11  (M.A. Filiberti) IGLOO sea-ice
12   !!            8.0  ! 1996-03  (M. Imbard, C. Levy, G. Madec, O. Marti, M. Guyon, A. Lazar,
13   !!                             P. Delecluse, L.Terray, M.A. Filiberti, J. Vialar, A.M. Treguier, M. Levy) release 8.0
14   !!            8.1  ! 1997-06  (M. Imbard, G. Madec)
15   !!            8.2  ! 1999-11  (M. Imbard, H. Goosse)  LIM sea-ice model
16   !!                 ! 1999-12  (V. Thierry, A-M. Treguier, M. Imbard, M-A. Foujols)  OPEN-MP
17   !!                 ! 2000-07  (J-M Molines, M. Imbard)  Open Boundary Conditions  (CLIPPER)
18   !!   NEMO     1.0  ! 2002-08  (G. Madec)  F90: Free form and modules
19   !!             -   ! 2004-06  (R. Redler, NEC CCRLE, Germany) add OASIS[3/4] coupled interfaces
20   !!             -   ! 2004-08  (C. Talandier) New trends organization
21   !!             -   ! 2005-06  (C. Ethe) Add the 1D configuration possibility
22   !!             -   ! 2005-11  (V. Garnier) Surface pressure gradient organization
23   !!             -   ! 2006-03  (L. Debreu, C. Mazauric)  Agrif implementation
24   !!             -   ! 2006-04  (G. Madec, R. Benshila)  Step reorganization
25   !!             -   ! 2007-07  (J. Chanut, A. Sellar) Unstructured open boundaries (BDY)
26   !!            3.2  ! 2009-08  (S. Masson)  open/write in the listing file in mpp
27   !!            3.3  ! 2010-05  (K. Mogensen, A. Weaver, M. Martin, D. Lea) Assimilation interface
28   !!             -   ! 2010-10  (C. Ethe, G. Madec) reorganisation of initialisation phase
29   !!            3.3.1! 2011-01  (A. R. Porter, STFC Daresbury) dynamical allocation
30   !!            3.4  ! 2011-11  (C. Harris) decomposition changes for running with CICE
31   !!                 ! 2012-05  (C. Calone, J. Simeon, G. Madec, C. Ethe) Add grid coarsening
32   !!----------------------------------------------------------------------
33
34   !!----------------------------------------------------------------------
35   !!   nemo_gcm       : solve ocean dynamics, tracer, biogeochemistry and/or sea-ice
36   !!   nemo_init      : initialization of the NEMO system
37   !!   nemo_ctl       : initialisation of the contol print
38   !!   nemo_closefile : close remaining open files
39   !!   nemo_alloc     : dynamical allocation
40   !!   nemo_partition : calculate MPP domain decomposition
41   !!   factorise      : calculate the factors of the no. of MPI processes
42   !!----------------------------------------------------------------------
43   USE step_oce        ! module used in the ocean time stepping module
44   USE cla             ! cross land advection               (tra_cla routine)
45   USE domcfg          ! domain configuration               (dom_cfg routine)
46   USE mppini          ! shared/distributed memory setting (mpp_init routine)
47   USE domain          ! domain initialization             (dom_init routine)
48#if defined key_nemocice_decomp
49   USE ice_domain_size, only: nx_global, ny_global
50#endif
51   USE tideini         ! tidal components initialization   (tide_ini routine)
52   USE bdyini          ! open boundary cond. setting       (bdy_init routine)
53   USE bdydta          ! open boundary cond. setting   (bdy_dta_init routine)
54   USE bdytides        ! open boundary cond. setting   (bdytide_init routine)
55   USE istate          ! initial state setting          (istate_init routine)
56   USE ldfdyn          ! lateral viscosity setting      (ldfdyn_init routine)
57   USE ldftra          ! lateral diffusivity setting    (ldftra_init routine)
58   USE zdfini          ! vertical physics setting          (zdf_init routine)
59   USE phycst          ! physical constant                  (par_cst routine)
60   USE trdini          ! dyn/tra trends initialization     (trd_init routine)
61   USE asminc          ! assimilation increments     
62   USE asmbkg          ! writing out state trajectory
63   USE asmbgc          ! biogeochemical assimilation increments
64   USE diaptr          ! poleward transports           (dia_ptr_init routine)
65   USE diadct          ! sections transports           (dia_dct_init routine)
66   USE diaobs          ! Observation diagnostics       (dia_obs_init routine)
67   USE lib_fortran     ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined)
68   USE step            ! NEMO time-stepping                 (stp     routine)
69   USE icbini          ! handle bergs, initialisation
70   USE icbstp          ! handle bergs, calving, themodynamics and transport
71   USE sbccpl 
72   USE cpl_oasis3      ! OASIS3 coupling
73   USE c1d             ! 1D configuration
74   USE step_c1d        ! Time stepping loop for the 1D configuration
75   USE dyndmp          ! Momentum damping
76#if defined key_top
77   USE trcini          ! passive tracer initialisation
78   USE trc, ONLY: numstr  ! tracer stats unit number
79#endif
80   USE lib_mpp         ! distributed memory computing
81#if defined key_iomput
82   USE xios
83#endif
84   USE sbctide, ONLY: lk_tide
85   USE crsini          ! initialise grid coarsening utility
86   USE lbcnfd, ONLY: isendto, nsndto, nfsloop, nfeloop ! Setup of north fold exchanges
87   USE sbc_oce, ONLY: lk_oasis
88   USE stopar
89   USE stopts
90
91   IMPLICIT NONE
92   PRIVATE
93
94   PUBLIC   nemo_gcm    ! called by model.F90
95   PUBLIC   nemo_init   ! needed by AGRIF
96   PUBLIC   nemo_alloc  ! needed by TAM
97
98   CHARACTER(lc) ::   cform_aaa="( /, 'AAAAAAAA', / ) "     ! flag for output listing
99
100   !!----------------------------------------------------------------------
101   !! NEMO/OPA 4.0 , NEMO Consortium (2011)
102   !! $Id$
103   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
104   !!----------------------------------------------------------------------
105CONTAINS
106
107   SUBROUTINE nemo_gcm
108      !!----------------------------------------------------------------------
109      !!                     ***  ROUTINE nemo_gcm  ***
110      !!
111      !! ** Purpose :   NEMO solves the primitive equations on an orthogonal
112      !!              curvilinear mesh on the sphere.
113      !!
114      !! ** Method  : - model general initialization
115      !!              - launch the time-stepping (stp routine)
116      !!              - finalize the run by closing files and communications
117      !!
118      !! References : Madec, Delecluse, Imbard, and Levy, 1997:  internal report, IPSL.
119      !!              Madec, 2008, internal report, IPSL.
120      !!----------------------------------------------------------------------
121      INTEGER ::   istp       ! time step index
122      !!----------------------------------------------------------------------
123      !
124#if defined key_agrif
125      CALL Agrif_Init_Grids()      ! AGRIF: set the meshes
126#endif
127
128      !                            !-----------------------!
129      CALL nemo_init               !==  Initialisations  ==!
130      !                            !-----------------------!
131#if defined key_agrif
132      CALL Agrif_Declare_Var_dom   ! AGRIF: set the meshes for DOM
133      CALL Agrif_Declare_Var       !  "      "   "   "      "  DYN/TRA
134# if defined key_top
135      CALL Agrif_Declare_Var_top   !  "      "   "   "      "  TOP
136# endif
137# if defined key_lim2
138      CALL Agrif_Declare_Var_lim2  !  "      "   "   "      "  LIM
139# endif
140#endif
141      ! check that all process are still there... If some process have an error,
142      ! they will never enter in step and other processes will wait until the end of the cpu time!
143      IF( lk_mpp )   CALL mpp_max( nstop )
144
145      IF(lwp) WRITE(numout,cform_aaa)   ! Flag AAAAAAA
146
147      !                            !-----------------------!
148      !                            !==   time stepping   ==!
149      !                            !-----------------------!
150      istp = nit000
151#if defined key_c1d
152         DO WHILE ( istp <= nitend .AND. nstop == 0 )
153            CALL stp_c1d( istp )
154            istp = istp + 1
155         END DO
156#else
157          IF( lk_asminc ) THEN
158             IF( ln_bkgwri ) CALL asm_bkg_wri( nit000 - 1 )    ! Output background fields
159             IF( ln_asmdin ) THEN                        ! Direct initialization
160                IF( ln_trainc ) CALL tra_asm_inc( nit000 - 1 )    ! Tracers
161                IF( ln_dyninc ) CALL dyn_asm_inc( nit000 - 1 )    ! Dynamics
162                IF( ln_sshinc ) CALL ssh_asm_inc( nit000 - 1 )    ! SSH
163                IF( lk_bgcinc ) CALL bgc_asm_inc( nit000 - 1 )    ! BGC
164             ENDIF
165          ENDIF
166
167#if defined key_agrif
168          CALL Agrif_Regrid()
169#endif
170
171         DO WHILE ( istp <= nitend .AND. nstop == 0 )
172#if defined key_agrif
173            CALL stp                         ! AGRIF: time stepping
174#else
175            IF (lk_oasis) CALL sbc_cpl_snd( istp )  ! Coupling to atmos
176       CALL stp( istp )
177            ! We don't couple on the final timestep because
178            ! our restart file has already been written
179            ! and contains all the necessary data for a
180            ! restart. sbc_cpl_snd could be called here
181            ! but it would require
182            ! a) A test to ensure it was not performed
183            !    on the very last time-step
184            ! b) the presence of another call to
185            !    sbc_cpl_snd call prior to the main DO loop
186            ! This solution produces identical results
187            ! with fewer lines of code.
188#endif
189            istp = istp + 1
190            IF( lk_mpp )   CALL mpp_max( nstop )
191         END DO
192#endif
193
194      IF( lk_diaobs   )   CALL dia_obs_wri
195      !
196      IF( ( lk_asminc ).AND.( ln_balwri ) ) CALL asm_bgc_bal_wri( nitend )  ! Output balancing increments
197      !
198      IF( ln_icebergs )   CALL icb_end( nitend )
199
200      !                            !------------------------!
201      !                            !==  finalize the run  ==!
202      !                            !------------------------!
203      IF(lwp) WRITE(numout,cform_aaa)   ! Flag AAAAAAA
204      !
205      IF( nstop /= 0 .AND. lwp ) THEN   ! error print
206         WRITE(numout,cform_err)
207         WRITE(numout,*) nstop, ' error have been found'
208      ENDIF
209      !
210#if defined key_agrif
211      IF( .NOT. Agrif_Root() ) THEN
212         CALL Agrif_ParentGrid_To_ChildGrid()
213         IF( lk_diaobs ) CALL dia_obs_wri
214         IF( nn_timing > 0 )   CALL timing_finalize
215         CALL Agrif_ChildGrid_To_ParentGrid()
216      ENDIF
217#endif
218      IF( nn_timing > 0 )   CALL timing_finalize
219      !
220      CALL nemo_closefile
221      !
222#if defined key_iomput
223      CALL xios_finalize                ! end mpp communications with xios
224      IF( lk_oasis ) CALL cpl_finalize    ! end coupling and mpp communications with OASIS
225#else
226      IF( lk_oasis ) THEN
227         CALL cpl_finalize              ! end coupling and mpp communications with OASIS
228      ELSE
229         IF( lk_mpp )   CALL mppstop    ! end mpp communications
230      ENDIF
231#endif
232      !
233      ! Met Office addition: if failed, return non-zero exit code
234      IF( nstop /= 0 )  CALL exit( 9 ) 
235      !
236   END SUBROUTINE nemo_gcm
237
238
239   SUBROUTINE nemo_init
240      !!----------------------------------------------------------------------
241      !!                     ***  ROUTINE nemo_init  ***
242      !!
243      !! ** Purpose :   initialization of the NEMO GCM
244      !!----------------------------------------------------------------------
245      INTEGER ::   ji            ! dummy loop indices
246      INTEGER ::   ilocal_comm   ! local integer
247      INTEGER ::   ios
248      CHARACTER(len=80), DIMENSION(16) ::   cltxt
249      !
250      NAMELIST/namctl/ ln_ctl  ,sn_cfctl, nn_print, nn_ictls, nn_ictle,   &
251         &             nn_isplt, nn_jsplt, nn_jctls, nn_jctle,   &
252         &             nn_bench, nn_timing, ln_flush
253      NAMELIST/namcfg/ cp_cfg, cp_cfz, jp_cfg, jpidta, jpjdta, jpkdta, jpiglo, jpjglo, &
254         &             jpizoom, jpjzoom, jperio, ln_use_jattr
255      !!----------------------------------------------------------------------
256      !
257      cltxt = ''
258      cxios_context = 'nemo'
259      !
260      !                             ! Open reference namelist and configuration namelist files
261      CALL ctl_opn( numnam_ref, 'namelist_ref', 'OLD', 'FORMATTED', 'SEQUENTIAL', -1, 6, .FALSE. )
262      CALL ctl_opn( numnam_cfg, 'namelist_cfg', 'OLD', 'FORMATTED', 'SEQUENTIAL', -1, 6, .FALSE. )
263      !
264      REWIND( numnam_ref )              ! Namelist namctl in reference namelist : Control prints & Benchmark
265      READ  ( numnam_ref, namctl, IOSTAT = ios, ERR = 901 )
266901   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namctl in reference namelist', .TRUE. )
267
268      REWIND( numnam_cfg )              ! Namelist namctl in confguration namelist : Control prints & Benchmark
269      READ  ( numnam_cfg, namctl, IOSTAT = ios, ERR = 902 )
270902   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namctl in configuration namelist', .TRUE. )
271
272      !
273      REWIND( numnam_ref )              ! Namelist namcfg in reference namelist : Control prints & Benchmark
274      READ  ( numnam_ref, namcfg, IOSTAT = ios, ERR = 903 )
275903   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namcfg in reference namelist', .TRUE. )
276
277      REWIND( numnam_cfg )              ! Namelist namcfg in confguration namelist : Control prints & Benchmark
278      READ  ( numnam_cfg, namcfg, IOSTAT = ios, ERR = 904 )
279904   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namcfg in configuration namelist', .TRUE. )   
280
281! Force values for AGRIF zoom (cf. agrif_user.F90)
282#if defined key_agrif
283   IF( .NOT. Agrif_Root() ) THEN
284      jpiglo  = nbcellsx + 2 + 2*nbghostcells
285      jpjglo  = nbcellsy + 2 + 2*nbghostcells
286      jpi     = ( jpiglo-2*jpreci + (jpni-1+0) ) / jpni + 2*jpreci
287      jpj     = ( jpjglo-2*jprecj + (jpnj-1+0) ) / jpnj + 2*jprecj
288      jpidta  = jpiglo
289      jpjdta  = jpjglo
290      jpizoom = 1
291      jpjzoom = 1
292      nperio  = 0
293      jperio  = 0
294      ln_use_jattr = .false.
295   ENDIF
296#endif
297      nprint    = nn_print          ! convert DOCTOR namelist names into OLD names
298      lflush    = ln_flush          ! convert namelist variable to model variable
299      !
300      !                             !--------------------------------------------!
301      !                             !  set communicator & select the local node  !
302      !                             !  NB: mynode also opens output.namelist.dyn !
303      !                             !      on unit number numond on first proc   !
304      !                             !--------------------------------------------!
305#if defined key_iomput
306      IF( Agrif_Root() ) THEN
307         IF( lk_oasis ) THEN
308            CALL cpl_init( "toyoce", ilocal_comm )                     ! nemo local communicator given by oasis
309            CALL xios_initialize( "not used",local_comm=ilocal_comm )    ! send nemo communicator to xios
310         ELSE
311            CALL  xios_initialize( "for_xios_mpi_id",return_comm=ilocal_comm )    ! nemo local communicator given by xios
312         ENDIF
313      ENDIF
314      ! Nodes selection (control print return in cltxt)
315      narea = mynode( cltxt, 'output.namelist.dyn', numnam_ref, numnam_cfg, numond , nstop, ilocal_comm )
316#else
317      IF( lk_oasis ) THEN
318         IF( Agrif_Root() ) THEN
319            CALL cpl_init( "toyoce", ilocal_comm )                      ! nemo local communicator given by oasis
320         ENDIF
321         ! Nodes selection (control print return in cltxt)
322         narea = mynode( cltxt, 'output.namelist.dyn', numnam_ref, numnam_cfg, numond , nstop, ilocal_comm )
323      ELSE
324         ilocal_comm = 0
325         ! Nodes selection (control print return in cltxt)
326         narea = mynode( cltxt, 'output.namelist.dyn', numnam_ref, numnam_cfg, numond , nstop )
327      ENDIF
328#endif
329      narea = narea + 1                                     ! mynode return the rank of proc (0 --> jpnij -1 )
330
331      IF( sn_cfctl%l_config ) THEN
332         ! Activate finer control of report outputs
333         ! optionally switch off output from selected areas (note this only
334         ! applies to output which does not involve global communications)
335         IF( ( narea < sn_cfctl%procmin .OR. narea > sn_cfctl%procmax  ) .OR. &
336           & ( MOD( narea - sn_cfctl%procmin, sn_cfctl%procincr ) /= 0 ) )    &
337           &   CALL nemo_set_cfctl( sn_cfctl, .FALSE., .FALSE. )
338      ELSE
339         ! Use ln_ctl to turn on or off all options.
340         CALL nemo_set_cfctl( sn_cfctl, ln_ctl, .TRUE. )
341      ENDIF
342
343      lwm = (narea == 1)                                    ! control of output namelists
344      lwp = (narea == 1) .OR. ln_ctl                        ! control of all listing output print
345
346      IF(lwm .AND. nprint > 2) THEN
347         ! write merged namelists from earlier to output namelist now that the
348         ! file has been opened in call to mynode. nammpp has already been
349         ! written in mynode (if lk_mpp_mpi)
350         WRITE( numond, namctl )
351         WRITE( numond, namcfg )
352      ENDIF
353
354      ! If dimensions of processor grid weren't specified in the namelist file
355      ! then we calculate them here now that we have our communicator size
356      IF( (jpni < 1) .OR. (jpnj < 1) )THEN
357#if   defined key_mpp_mpi
358         IF( Agrif_Root() ) CALL nemo_partition(mppsize)
359#else
360         jpni  = 1
361         jpnj  = 1
362         jpnij = jpni*jpnj
363#endif
364      END IF
365
366      ! Calculate domain dimensions given calculated jpni and jpnj
367      ! This used to be done in par_oce.F90 when they were parameters rather
368      ! than variables
369      IF( Agrif_Root() ) THEN
370#if defined key_nemocice_decomp
371         jpi = ( nx_global+2-2*jpreci + (jpni-1) ) / jpni + 2*jpreci ! first  dim.
372         jpj = ( ny_global+2-2*jprecj + (jpnj-1) ) / jpnj + 2*jprecj ! second dim.
373#else
374         jpi = ( jpiglo-2*jpreci + (jpni-1) ) / jpni + 2*jpreci   ! first  dim.
375         jpj = ( jpjglo-2*jprecj + (jpnj-1) ) / jpnj + 2*jprecj   ! second dim.
376#endif
377      ENDIF         
378         jpk = jpkdta                                             ! third dim
379#if defined key_agrif
380         ! simple trick to use same vertical grid as parent
381         ! but different number of levels:
382         ! Save maximum number of levels in jpkdta, then define all vertical grids
383         ! with this number.
384         ! Suppress once vertical online interpolation is ok
385         IF(.NOT.Agrif_Root()) jpkdta = Agrif_Parent(jpkdta)
386#endif
387         jpim1 = jpi-1                                            ! inner domain indices
388         jpjm1 = jpj-1                                            !   "           "
389         jpkm1 = jpk-1                                            !   "           "
390         jpij  = jpi*jpj                                          !  jpi x j
391
392      IF(lwp) THEN                            ! open listing units
393         !
394         CALL ctl_opn( numout, 'ocean.output', 'REPLACE', 'FORMATTED', 'SEQUENTIAL', -1, 6, .FALSE., narea )
395         !
396         WRITE(numout,*)
397         WRITE(numout,*) '   CNRS - NERC - Met OFFICE - MERCATOR-ocean - INGV - CMCC'
398         WRITE(numout,*) '                       NEMO team'
399         WRITE(numout,*) '            Ocean General Circulation Model'
400         WRITE(numout,*) '                  version 3.6  (2015) '
401         WRITE(numout,*)
402         WRITE(numout,*)
403         DO ji = 1, SIZE(cltxt)
404            IF( TRIM(cltxt(ji)) /= '' )   WRITE(numout,*) cltxt(ji)      ! control print of mynode
405         END DO
406         WRITE(numout,cform_aaa)                                         ! Flag AAAAAAA
407         IF(lflush) CALL flush(numout)
408         !
409      ENDIF
410
411      ! Now we know the dimensions of the grid and numout has been set we can
412      ! allocate arrays
413      CALL nemo_alloc()
414
415      !                             !-------------------------------!
416      !                             !  NEMO general initialization  !
417      !                             !-------------------------------!
418
419      CALL nemo_ctl                          ! Control prints & Benchmark
420
421      !                                      ! Domain decomposition
422      IF( jpni*jpnj == jpnij ) THEN   ;   CALL mpp_init      ! standard cutting out
423      ELSE                            ;   CALL mpp_init2     ! eliminate land processors
424      ENDIF
425      !
426      IF( nn_timing > 0 )  CALL timing_init
427      !
428      !                                      ! General initialization
429                            CALL     phy_cst    ! Physical constants
430                            CALL     eos_init   ! Equation of state
431      IF( lk_c1d        )   CALL     c1d_init   ! 1D column configuration
432                            CALL     dom_cfg    ! Domain configuration
433                            CALL     dom_init   ! Domain
434
435      IF( ln_nnogather )    CALL nemo_northcomms   ! Initialise the northfold neighbour lists (must be done after the masks are defined)
436
437      IF( ln_ctl        )   CALL prt_ctl_init   ! Print control
438
439                            CALL  istate_init   ! ocean initial state (Dynamics and tracers)
440
441      IF( lk_tide       )   CALL    tide_init( nit000 )    ! Initialisation of the tidal harmonics
442
443                            CALL     sbc_init   ! Forcings : surface module (clem: moved here for bdy purpose)
444
445      IF( lk_bdy        )   CALL     bdy_init   ! Open boundaries initialisation
446      IF( lk_bdy        )   CALL bdy_dta_init   ! Open boundaries initialisation of external data arrays
447      IF( lk_bdy .AND. lk_tide )   &
448         &                  CALL bdytide_init   ! Open boundaries initialisation of tidal harmonic forcing
449
450                            CALL dyn_nept_init  ! simplified form of Neptune effect
451      !     
452      IF( ln_crs        )   CALL     crs_init   ! Domain initialization of coarsened grid
453      !
454                                ! Ocean physics
455      !                                         ! Vertical physics
456                            CALL     zdf_init      ! namelist read
457                            CALL zdf_bfr_init      ! bottom friction
458      IF( lk_zdfric     )   CALL zdf_ric_init      ! Richardson number dependent Kz
459      IF( lk_zdftke     )   CALL zdf_tke_init      ! TKE closure scheme
460      IF( lk_zdfgls     )   CALL zdf_gls_init      ! GLS closure scheme
461      IF( lk_zdfkpp     )   CALL zdf_kpp_init      ! KPP closure scheme
462      IF( lk_zdftmx     )   CALL zdf_tmx_init      ! tidal vertical mixing
463      IF( lk_zdfddm .AND. .NOT. lk_zdfkpp )   &
464         &                  CALL zdf_ddm_init      ! double diffusive mixing
465      !                                         ! Lateral physics
466                            CALL ldf_tra_init      ! Lateral ocean tracer physics
467                            CALL ldf_dyn_init      ! Lateral ocean momentum physics
468      IF( lk_ldfslp     )   CALL ldf_slp_init      ! slope of lateral mixing
469
470      !                                     ! Active tracers
471                            CALL tra_qsr_init   ! penetrative solar radiation qsr
472                            CALL tra_bbc_init   ! bottom heat flux
473      IF( lk_trabbl     )   CALL tra_bbl_init   ! advective (and/or diffusive) bottom boundary layer scheme
474                            CALL tra_dmp_init   ! internal damping trends- tracers
475                            CALL tra_adv_init   ! horizontal & vertical advection
476                            CALL tra_ldf_init   ! lateral mixing
477                            CALL tra_zdf_init   ! vertical mixing and after tracer fields
478
479      !                                     ! Dynamics
480      IF( lk_c1d        )   CALL dyn_dmp_init   ! internal damping trends- momentum
481                            CALL dyn_adv_init   ! advection (vector or flux form)
482                            CALL dyn_vor_init   ! vorticity term including Coriolis
483                            CALL dyn_ldf_init   ! lateral mixing
484                            CALL dyn_hpg_init   ! horizontal gradient of Hydrostatic pressure
485                            CALL dyn_zdf_init   ! vertical diffusion
486                            CALL dyn_spg_init   ! surface pressure gradient
487
488      !                                     ! Misc. options
489      IF( nn_cla == 1 .AND. cp_cfg == 'orca' .AND. jp_cfg == 2 )   CALL cla_init       ! Cross Land Advection
490                            CALL icb_init( rdt, nit000)   ! initialise icebergs instance
491                            CALL sto_par_init   ! Stochastic parametrization
492      IF( ln_sto_eos     )  CALL sto_pts_init   ! RRandom T/S fluctuations
493     
494#if defined key_top
495      !                                     ! Passive tracers
496                            CALL     trc_init
497#endif
498      !                                     ! Diagnostics
499      IF( lk_floats     )   CALL     flo_init   ! drifting Floats
500                            CALL dia_ptr_init   ! Poleward TRansports initialization
501      IF( lk_diadct     )   CALL dia_dct_init   ! Sections tranports
502                            CALL dia_hsb_init   ! heat content, salt content and volume budgets
503                            CALL     trd_init   ! Mixed-layer/Vorticity/Integral constraints trends
504                            CALL     bias_init  ! Pressure correction bias
505      IF( lk_diaobs     ) THEN                  ! Observation & model comparison
506                            CALL dia_obs_init            ! Initialize observational data
507                            CALL dia_obs( nit000 - 1 )   ! Observation operator for restart
508      ENDIF
509
510      !                                     ! Assimilation increments
511      IF( lk_asminc     )   CALL asm_inc_init   ! Initialize assimilation increments
512      IF(lwp) WRITE(numout,*) 'Euler time step switch is ', neuler
513      IF(lwp .AND. lflush) CALL flush(numout)
514     
515      IF (nstop > 0) THEN
516        CALL CTL_STOP('STOP','Critical errors in NEMO initialisation')
517      END IF
518
519      !
520   END SUBROUTINE nemo_init
521
522
523   SUBROUTINE nemo_ctl
524      !!----------------------------------------------------------------------
525      !!                     ***  ROUTINE nemo_ctl  ***
526      !!
527      !! ** Purpose :   control print setting
528      !!
529      !! ** Method  : - print namctl information and check some consistencies
530      !!----------------------------------------------------------------------
531      !
532      IF(lwp) THEN                  ! control print
533         WRITE(numout,*)
534         WRITE(numout,*) 'nemo_ctl: Control prints & Benchmark'
535         WRITE(numout,*) '~~~~~~~ '
536         WRITE(numout,*) '   Namelist namctl'
537         WRITE(numout,*) '      run control (for debugging)     ln_ctl     = ', ln_ctl
538         WRITE(numout,*) '       finer control over o/p sn_cfctl%l_config  = ', sn_cfctl%l_config
539         WRITE(numout,*) '                              sn_cfctl%l_runstat = ', sn_cfctl%l_runstat
540         WRITE(numout,*) '                              sn_cfctl%l_trcstat = ', sn_cfctl%l_trcstat
541         WRITE(numout,*) '                              sn_cfctl%l_oceout  = ', sn_cfctl%l_oceout
542         WRITE(numout,*) '                              sn_cfctl%l_layout  = ', sn_cfctl%l_layout
543         WRITE(numout,*) '                              sn_cfctl%l_mppout  = ', sn_cfctl%l_mppout
544         WRITE(numout,*) '                              sn_cfctl%l_mpptop  = ', sn_cfctl%l_mpptop
545         WRITE(numout,*) '                              sn_cfctl%procmin   = ', sn_cfctl%procmin 
546         WRITE(numout,*) '                              sn_cfctl%procmax   = ', sn_cfctl%procmax 
547         WRITE(numout,*) '                              sn_cfctl%procincr  = ', sn_cfctl%procincr 
548         WRITE(numout,*) '                              sn_cfctl%ptimincr  = ', sn_cfctl%ptimincr 
549         WRITE(numout,*) '      level of print                  nn_print   = ', nn_print
550         WRITE(numout,*) '      flush numout/stat               ln_flush   = ', ln_flush
551         WRITE(numout,*) '      Start i indice for SUM control  nn_ictls   = ', nn_ictls
552         WRITE(numout,*) '      End i indice for SUM control    nn_ictle   = ', nn_ictle
553         WRITE(numout,*) '      Start j indice for SUM control  nn_jctls   = ', nn_jctls
554         WRITE(numout,*) '      End j indice for SUM control    nn_jctle   = ', nn_jctle
555         WRITE(numout,*) '      number of proc. following i     nn_isplt   = ', nn_isplt
556         WRITE(numout,*) '      number of proc. following j     nn_jsplt   = ', nn_jsplt
557         WRITE(numout,*) '      benchmark parameter (0/1)       nn_bench   = ', nn_bench
558         WRITE(numout,*) '      timing activated    (0/1)       nn_timing  = ', nn_timing
559         IF(lflush) CALL flush(numout)
560      ENDIF
561      !
562      nictls    = nn_ictls
563      nictle    = nn_ictle
564      njctls    = nn_jctls
565      njctle    = nn_jctle
566      isplt     = nn_isplt
567      jsplt     = nn_jsplt
568      nbench    = nn_bench
569
570      IF(lwp) THEN                  ! control print
571         WRITE(numout,*)
572         WRITE(numout,*) 'namcfg  : configuration initialization through namelist read'
573         WRITE(numout,*) '~~~~~~~ '
574         WRITE(numout,*) '   Namelist namcfg'
575         WRITE(numout,*) '      configuration name              cp_cfg      = ', TRIM(cp_cfg)
576         WRITE(numout,*) '      configuration zoom name         cp_cfz      = ', TRIM(cp_cfz)
577         WRITE(numout,*) '      configuration resolution        jp_cfg      = ', jp_cfg
578         WRITE(numout,*) '      1st lateral dimension ( >= jpi ) jpidta     = ', jpidta
579         WRITE(numout,*) '      2nd    "         "    ( >= jpj ) jpjdta     = ', jpjdta
580         WRITE(numout,*) '      3nd    "         "               jpkdta     = ', jpkdta
581         WRITE(numout,*) '      1st dimension of global domain in i jpiglo  = ', jpiglo
582         WRITE(numout,*) '      2nd    -                  -    in j jpjglo  = ', jpjglo
583         WRITE(numout,*) '      left bottom i index of the zoom (in data domain) jpizoom = ', jpizoom
584         WRITE(numout,*) '      left bottom j index of the zoom (in data domain) jpizoom = ', jpjzoom
585         WRITE(numout,*) '      lateral cond. type (between 0 and 6) jperio = ', jperio   
586         WRITE(numout,*) '      use file attribute if exists as i/p j-start ln_use_jattr = ', ln_use_jattr
587         IF(lflush) CALL flush(numout)
588      ENDIF
589      !                             ! Parameter control
590      !
591      IF( ln_ctl ) THEN                 ! sub-domain area indices for the control prints
592         IF( lk_mpp .AND. jpnij > 1 ) THEN
593            isplt = jpni   ;   jsplt = jpnj   ;   ijsplt = jpni*jpnj   ! the domain is forced to the real split domain
594         ELSE
595            IF( isplt == 1 .AND. jsplt == 1  ) THEN
596               CALL ctl_warn( ' - isplt & jsplt are equal to 1',   &
597                  &           ' - the print control will be done over the whole domain' )
598            ENDIF
599            ijsplt = isplt * jsplt            ! total number of processors ijsplt
600         ENDIF
601         IF(lwp) WRITE(numout,*)'          - The total number of processors over which the'
602         IF(lwp) WRITE(numout,*)'            print control will be done is ijsplt : ', ijsplt
603         IF(lwp .AND. lflush) CALL flush(numout)
604         !
605         !                              ! indices used for the SUM control
606         IF( nictls+nictle+njctls+njctle == 0 )   THEN    ! print control done over the default area
607            lsp_area = .FALSE.
608         ELSE                                             ! print control done over a specific  area
609            lsp_area = .TRUE.
610            IF( nictls < 1 .OR. nictls > jpiglo )   THEN
611               CALL ctl_warn( '          - nictls must be 1<=nictls>=jpiglo, it is forced to 1' )
612               nictls = 1
613            ENDIF
614            IF( nictle < 1 .OR. nictle > jpiglo )   THEN
615               CALL ctl_warn( '          - nictle must be 1<=nictle>=jpiglo, it is forced to jpiglo' )
616               nictle = jpiglo
617            ENDIF
618            IF( njctls < 1 .OR. njctls > jpjglo )   THEN
619               CALL ctl_warn( '          - njctls must be 1<=njctls>=jpjglo, it is forced to 1' )
620               njctls = 1
621            ENDIF
622            IF( njctle < 1 .OR. njctle > jpjglo )   THEN
623               CALL ctl_warn( '          - njctle must be 1<=njctle>=jpjglo, it is forced to jpjglo' )
624               njctle = jpjglo
625            ENDIF
626         ENDIF
627      ENDIF
628      !
629      IF( nbench == 1 ) THEN              ! Benchmark
630         SELECT CASE ( cp_cfg )
631         CASE ( 'gyre' )   ;   CALL ctl_warn( ' The Benchmark is activated ' )
632         CASE DEFAULT      ;   CALL ctl_stop( ' The Benchmark is based on the GYRE configuration:',   &
633            &                                 ' cp_cfg = "gyre" in namelist &namcfg or set nbench = 0' )
634         END SELECT
635      ENDIF
636      !
637      IF( 1_wp /= SIGN(1._wp,-0._wp)  )   CALL ctl_stop( 'nemo_ctl: The intrinsec SIGN function follows ',  &
638         &                                               'f2003 standard. '                              ,  &
639         &                                               'Compile with key_nosignedzero enabled' )
640      !
641   END SUBROUTINE nemo_ctl
642
643
644   SUBROUTINE nemo_closefile
645      !!----------------------------------------------------------------------
646      !!                     ***  ROUTINE nemo_closefile  ***
647      !!
648      !! ** Purpose :   Close the files
649      !!----------------------------------------------------------------------
650      !
651      IF( lk_mpp )   CALL mppsync
652      !
653      CALL iom_close                                 ! close all input/output files managed by iom_*
654      !
655      IF( numstp          /= -1 )   CLOSE( numstp          )   ! time-step file
656      IF( numsol          /= -1 )   CLOSE( numsol          )   ! solver file
657      IF( numnam_ref      /= -1 )   CLOSE( numnam_ref      )   ! oce reference namelist
658      IF( numnam_cfg      /= -1 )   CLOSE( numnam_cfg      )   ! oce configuration namelist
659      IF( lwm.AND.numond  /= -1 )   CLOSE( numond          )   ! oce output namelist
660      IF( numnam_ice_ref  /= -1 )   CLOSE( numnam_ice_ref  )   ! ice reference namelist
661      IF( numnam_ice_cfg  /= -1 )   CLOSE( numnam_ice_cfg  )   ! ice configuration namelist
662      IF( lwm.AND.numoni  /= -1 )   CLOSE( numoni          )   ! ice output namelist
663      IF( numevo_ice      /= -1 )   CLOSE( numevo_ice      )   ! ice variables (temp. evolution)
664      IF( numout          /=  6 )   CLOSE( numout          )   ! standard model output file
665      IF( numdct_vol      /= -1 )   CLOSE( numdct_vol      )   ! volume transports
666      IF( numdct_heat     /= -1 )   CLOSE( numdct_heat     )   ! heat transports
667      IF( numdct_salt     /= -1 )   CLOSE( numdct_salt     )   ! salt transports
668#if defined key_top
669      IF( numstr          /= -1 )   CLOSE( numstr          )   ! tracer statistics
670#endif
671      !
672      numout = 6                                     ! redefine numout in case it is used after this point...
673      !
674   END SUBROUTINE nemo_closefile
675
676
677   SUBROUTINE nemo_alloc
678      !!----------------------------------------------------------------------
679      !!                     ***  ROUTINE nemo_alloc  ***
680      !!
681      !! ** Purpose :   Allocate all the dynamic arrays of the OPA modules
682      !!
683      !! ** Method  :
684      !!----------------------------------------------------------------------
685      USE diawri    , ONLY: dia_wri_alloc
686      USE insitu_tem, ONLY: insitu_tem_alloc
687      USE dom_oce   , ONLY: dom_oce_alloc
688      USE ldfdyn_oce, ONLY: ldfdyn_oce_alloc
689      USE ldftra_oce, ONLY: ldftra_oce_alloc
690      USE trc_oce   , ONLY: trc_oce_alloc
691#if defined key_diadct 
692      USE diadct    , ONLY: diadct_alloc 
693#endif 
694#if defined key_bdy
695      USE bdy_oce   , ONLY: bdy_oce_alloc
696#endif
697      !
698      INTEGER :: ierr
699      !!----------------------------------------------------------------------
700      !
701      ierr =        oce_alloc       ()          ! ocean
702      ierr = ierr + dia_wri_alloc   ()
703      ierr = ierr + insitu_tem_alloc()
704      ierr = ierr + dom_oce_alloc   ()          ! ocean domain
705      ierr = ierr + ldfdyn_oce_alloc()          ! ocean lateral  physics : dynamics
706      ierr = ierr + ldftra_oce_alloc()          ! ocean lateral  physics : tracers
707      ierr = ierr + zdf_oce_alloc   ()          ! ocean vertical physics
708      !
709      ierr = ierr + trc_oce_alloc   ()          ! shared TRC / TRA arrays
710      !
711#if defined key_diadct 
712      ierr = ierr + diadct_alloc    ()          !
713#endif 
714#if defined key_bdy
715      ierr = ierr + bdy_oce_alloc   ()          ! bdy masks (incl. initialization)
716#endif
717      !
718      IF( lk_mpp    )   CALL mpp_sum( ierr )
719      IF( ierr /= 0 )   CALL ctl_stop( 'STOP', 'nemo_alloc : unable to allocate standard ocean arrays' )
720      !
721   END SUBROUTINE nemo_alloc
722
723   SUBROUTINE nemo_set_cfctl(sn_cfctl, setto, for_all )
724      !!----------------------------------------------------------------------
725      !!                     ***  ROUTINE nemo_set_cfctl  ***
726      !!
727      !! ** Purpose :   Set elements of the output control structure to setto.
728      !!                for_all should be .false. unless all areas are to be
729      !!                treated identically.
730      !!
731      !! ** Method  :   Note this routine can be used to switch on/off some
732      !!                types of output for selected areas but any output types
733      !!                that involve global communications (e.g. mpp_max, glob_sum)
734      !!                should be protected from selective switching by the
735      !!                for_all argument
736      !!----------------------------------------------------------------------
737      LOGICAL :: setto, for_all
738      TYPE(sn_ctl) :: sn_cfctl
739      !!----------------------------------------------------------------------
740      IF( for_all ) THEN
741         sn_cfctl%l_runstat = setto
742         sn_cfctl%l_trcstat = setto
743      ENDIF
744      sn_cfctl%l_oceout  = setto
745      sn_cfctl%l_layout  = setto
746      sn_cfctl%l_mppout  = setto
747      sn_cfctl%l_mpptop  = setto
748   END SUBROUTINE nemo_set_cfctl
749
750   SUBROUTINE nemo_partition( num_pes )
751      !!----------------------------------------------------------------------
752      !!                 ***  ROUTINE nemo_partition  ***
753      !!
754      !! ** Purpose :
755      !!
756      !! ** Method  :
757      !!----------------------------------------------------------------------
758      INTEGER, INTENT(in) ::   num_pes   ! The number of MPI processes we have
759      !
760      INTEGER, PARAMETER :: nfactmax = 20
761      INTEGER :: nfact ! The no. of factors returned
762      INTEGER :: ierr  ! Error flag
763      INTEGER :: ji
764      INTEGER :: idiff, mindiff, imin ! For choosing pair of factors that are closest in value
765      INTEGER, DIMENSION(nfactmax) :: ifact ! Array of factors
766      !!----------------------------------------------------------------------
767      !
768      ierr = 0
769      !
770      CALL factorise( ifact, nfactmax, nfact, num_pes, ierr )
771      !
772      IF( nfact <= 1 ) THEN
773         WRITE (numout, *) 'WARNING: factorisation of number of PEs failed'
774         WRITE (numout, *) '       : using grid of ',num_pes,' x 1'
775         IF(lflush) CALL flush(numout)
776         jpnj = 1
777         jpni = num_pes
778      ELSE
779         ! Search through factors for the pair that are closest in value
780         mindiff = 1000000
781         imin    = 1
782         DO ji = 1, nfact-1, 2
783            idiff = ABS( ifact(ji) - ifact(ji+1) )
784            IF( idiff < mindiff ) THEN
785               mindiff = idiff
786               imin = ji
787            ENDIF
788         END DO
789         jpnj = ifact(imin)
790         jpni = ifact(imin + 1)
791      ENDIF
792      !
793      jpnij = jpni*jpnj
794      !
795   END SUBROUTINE nemo_partition
796
797
798   SUBROUTINE factorise( kfax, kmaxfax, knfax, kn, kerr )
799      !!----------------------------------------------------------------------
800      !!                     ***  ROUTINE factorise  ***
801      !!
802      !! ** Purpose :   return the prime factors of n.
803      !!                knfax factors are returned in array kfax which is of
804      !!                maximum dimension kmaxfax.
805      !! ** Method  :
806      !!----------------------------------------------------------------------
807      INTEGER                    , INTENT(in   ) ::   kn, kmaxfax
808      INTEGER                    , INTENT(  out) ::   kerr, knfax
809      INTEGER, DIMENSION(kmaxfax), INTENT(  out) ::   kfax
810      !
811      INTEGER :: ifac, jl, inu
812      INTEGER, PARAMETER :: ntest = 14
813      INTEGER, DIMENSION(ntest) :: ilfax
814      !
815      ! ilfax contains the set of allowed factors.
816      ilfax(:) = (/(2**jl,jl=ntest,1,-1)/)
817      !!----------------------------------------------------------------------
818      ! ilfax contains the set of allowed factors.
819      ilfax(:) = (/(2**jl,jl=ntest,1,-1)/)
820
821      ! Clear the error flag and initialise output vars
822      kerr = 0
823      kfax = 1
824      knfax = 0
825
826      ! Find the factors of n.
827      IF( kn == 1 )   GOTO 20
828
829      ! nu holds the unfactorised part of the number.
830      ! knfax holds the number of factors found.
831      ! l points to the allowed factor list.
832      ! ifac holds the current factor.
833
834      inu   = kn
835      knfax = 0
836
837      DO jl = ntest, 1, -1
838         !
839         ifac = ilfax(jl)
840         IF( ifac > inu )   CYCLE
841
842         ! Test whether the factor will divide.
843
844         IF( MOD(inu,ifac) == 0 ) THEN
845            !
846            knfax = knfax + 1            ! Add the factor to the list
847            IF( knfax > kmaxfax ) THEN
848               kerr = 6
849               write (*,*) 'FACTOR: insufficient space in factor array ', knfax
850               return
851            ENDIF
852            kfax(knfax) = ifac
853            ! Store the other factor that goes with this one
854            knfax = knfax + 1
855            kfax(knfax) = inu / ifac
856            !WRITE (*,*) 'ARPDBG, factors ',knfax-1,' & ',knfax,' are ', kfax(knfax-1),' and ',kfax(knfax)
857         ENDIF
858         !
859      END DO
860
861   20 CONTINUE      ! Label 20 is the exit point from the factor search loop.
862      !
863   END SUBROUTINE factorise
864
865#if defined key_mpp_mpi
866
867   SUBROUTINE nemo_northcomms
868      !!======================================================================
869      !!                     ***  ROUTINE  nemo_northcomms  ***
870      !! nemo_northcomms    :  Setup for north fold exchanges with explicit
871      !!                       point-to-point messaging
872      !!=====================================================================
873      !!----------------------------------------------------------------------
874      !!
875      !! ** Purpose :   Initialization of the northern neighbours lists.
876      !!----------------------------------------------------------------------
877      !!    1.0  ! 2011-10  (A. C. Coward, NOCS & J. Donners, PRACE)
878      !!    2.0  ! 2013-06 Setup avoiding MPI communication (I. Epicoco, S. Mocavero, CMCC)
879      !!----------------------------------------------------------------------
880
881      INTEGER  ::   sxM, dxM, sxT, dxT, jn
882      INTEGER  ::   njmppmax
883
884      njmppmax = MAXVAL( njmppt )
885   
886      !initializes the north-fold communication variables
887      isendto(:) = 0
888      nsndto = 0
889
890      !if I am a process in the north
891      IF ( njmpp == njmppmax ) THEN
892          !sxM is the first point (in the global domain) needed to compute the
893          !north-fold for the current process
894          sxM = jpiglo - nimppt(narea) - nlcit(narea) + 1
895          !dxM is the last point (in the global domain) needed to compute the
896          !north-fold for the current process
897          dxM = jpiglo - nimppt(narea) + 2
898
899          !loop over the other north-fold processes to find the processes
900          !managing the points belonging to the sxT-dxT range
901 
902          DO jn = 1, jpni
903                !sxT is the first point (in the global domain) of the jn
904                !process
905                sxT = nfiimpp(jn, jpnj)
906                !dxT is the last point (in the global domain) of the jn
907                !process
908                dxT = nfiimpp(jn, jpnj) + nfilcit(jn, jpnj) - 1
909                IF ((sxM .gt. sxT) .AND. (sxM .lt. dxT)) THEN
910                   nsndto = nsndto + 1
911                     isendto(nsndto) = jn
912                ELSEIF ((sxM .le. sxT) .AND. (dxM .ge. dxT)) THEN
913                   nsndto = nsndto + 1
914                     isendto(nsndto) = jn
915                ELSEIF ((dxM .lt. dxT) .AND. (sxT .lt. dxM)) THEN
916                   nsndto = nsndto + 1
917                     isendto(nsndto) = jn
918                END IF
919          END DO
920          nfsloop = 1
921          nfeloop = nlci
922          DO jn = 2,jpni-1
923           IF(nfipproc(jn,jpnj) .eq. (narea - 1)) THEN
924              IF (nfipproc(jn - 1 ,jpnj) .eq. -1) THEN
925                 nfsloop = nldi
926              ENDIF
927              IF (nfipproc(jn + 1,jpnj) .eq. -1) THEN
928                 nfeloop = nlei
929              ENDIF
930           ENDIF
931        END DO
932
933      ENDIF
934      l_north_nogather = .TRUE.
935   END SUBROUTINE nemo_northcomms
936#else
937   SUBROUTINE nemo_northcomms      ! Dummy routine
938      WRITE(*,*) 'nemo_northcomms: You should not have seen this print! error?'
939   END SUBROUTINE nemo_northcomms
940#endif
941
942   !!======================================================================
943END MODULE nemogcm
944
945
Note: See TracBrowser for help on using the repository browser.