source: trunk/NEMOGCM/NEMO/OOO_SRC/nemogcm.F90 @ 4624

Last change on this file since 4624 was 4624, checked in by acc, 7 years ago

#1305. Fix slow start-up problems on some systems by introducing and using lwm logical to restrict output of merged namelists to the first (or only) processor. lwm is true only on the first processor regardless of ln_ctl. Small changes to all flavours of nemogcm.F90 are also required to write namctl and namcfg after the call to mynode which now opens output.namelist.dyn and writes nammpp.

  • Property svn:keywords set to Id
File size: 30.6 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   !!----------------------------------------------------------------------
32
33   !!----------------------------------------------------------------------
34   !!   nemo_gcm       : solve ocean dynamics, tracer, biogeochemistry and/or sea-ice
35   !!   nemo_init      : initialization of the NEMO system
36   !!   nemo_ctl       : initialisation of the contol print
37   !!   nemo_closefile : close remaining open files
38   !!   nemo_alloc     : dynamical allocation
39   !!   nemo_partition : calculate MPP domain decomposition
40   !!   factorise      : calculate the factors of the no. of MPI processes
41   !!----------------------------------------------------------------------
42   USE step_oce        ! module used in the ocean time stepping module
43   USE domcfg          ! domain configuration               (dom_cfg routine)
44   USE mppini          ! shared/distributed memory setting (mpp_init routine)
45   USE domain          ! domain initialization             (dom_init routine)
46#if defined key_nemocice_decomp
47   USE ice_domain_size, only: nx_global, ny_global
48#endif
49   USE istate          ! initial state setting          (istate_init routine)
50   USE phycst          ! physical constant                  (par_cst routine)
51   USE diaobs          ! Observation diagnostics       (dia_obs_init routine)
52   USE lib_fortran     ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined)
53   USE step            ! NEMO time-stepping                 (stp     routine)
54   USE icbini          ! handle bergs, initialisation
55   USE icbstp          ! handle bergs, calving, themodynamics and transport
56#if defined key_oasis3
57   USE cpl_oasis3      ! OASIS3 coupling
58#elif defined key_oasis4
59   USE cpl_oasis4      ! OASIS4 coupling (not working)
60#endif
61   USE lib_mpp         ! distributed memory computing
62#if defined key_iomput
63   USE xios
64#endif
65   USE ooo_data        ! Offline obs_oper data
66   USE ooo_read        ! Offline obs_oper read routines
67   USE ooo_intp        ! Offline obs_oper interpolation
68
69   IMPLICIT NONE
70   PRIVATE
71
72   PUBLIC   nemo_gcm    ! called by nemo.f90
73   PUBLIC   nemo_init   ! needed by AGRIF
74   PUBLIC   nemo_alloc  ! needed by TAM
75
76   CHARACTER(lc) ::   cform_aaa="( /, 'AAAAAAAA', / ) "     ! flag for output listing
77
78   !!----------------------------------------------------------------------
79   !! NEMO/OPA 4.0 , NEMO Consortium (2011)
80   !! $Id$
81   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
82   !!----------------------------------------------------------------------
83CONTAINS
84
85   SUBROUTINE nemo_gcm
86         !!----------------------------------------------------------------------
87         !!                    ***  SUBROUTINE offline_obs_oper ***
88         !!
89         !! ** Purpose : To use NEMO components to interpolate model fields
90         !!              to observation space.
91         !!
92         !! ** Method : 1. Initialise NEMO
93         !!             2. Initialise offline obs_oper
94         !!             3. Cycle through match ups
95         !!             4. Write results to file
96         !!
97         !!----------------------------------------------------------------------
98         !! Class 4 output stream switch
99         USE obs_fbm, ONLY: ln_cl4
100         !! Initialise NEMO
101         CALL nemo_init
102         !! Initialise Offline obs_oper data
103         CALL ooo_data_init( ln_cl4 )
104         !! Loop over various model counterparts
105         DO jimatch = 1, cl4_match_len
106            IF (jimatch .GT. 1) THEN
107               !! Initialise obs_oper
108               CALL dia_obs_init
109            END IF
110            !! Interpolate to observation space
111            CALL ooo_interp
112            !! Pipe to output files
113            CALL dia_obs_wri
114            !! Reset the obs_oper between
115            CALL dia_obs_dealloc
116         END DO
117         !! Safely stop MPI
118         IF(lk_mpp) CALL mppstop  ! end mpp communications
119   END SUBROUTINE nemo_gcm
120
121   SUBROUTINE nemo_init
122      !!----------------------------------------------------------------------
123      !!                     ***  ROUTINE nemo_init  ***
124      !!
125      !! ** Purpose :   initialization of the NEMO GCM
126      !!----------------------------------------------------------------------
127      INTEGER ::   ji            ! dummy loop indices
128      INTEGER ::   ilocal_comm   ! local integer
129      CHARACTER(len=80), DIMENSION(16) ::   cltxt
130      !!
131      NAMELIST/namctl/ ln_ctl, nn_print, nn_ictls, nn_ictle,   &
132         &             nn_isplt, nn_jsplt, nn_jctls, nn_jctle,   &
133         &             nn_bench, nn_timing
134      NAMELIST/namcfg/ cp_cfg, cp_cfz, jp_cfg, jpidta, jpjdta, jpkdta, jpiglo, jpjglo, &
135         &             jpizoom, jpjzoom, jperio
136      !!----------------------------------------------------------------------
137      !
138      cltxt = ''
139      !
140      !                             ! Open reference namelist and configuration namelist files
141      CALL ctl_opn( numnam_ref, 'namelist_ref', 'OLD', 'FORMATTED', 'SEQUENTIAL', -1, 6, .FALSE. )
142      CALL ctl_opn( numnam_cfg, 'namelist_cfg', 'OLD', 'FORMATTED', 'SEQUENTIAL', -1, 6, .FALSE. )
143      !
144      REWIND( numnam_ref )              ! Namelist namctl in reference namelist : Control prints & Benchmark
145      READ  ( numnam_ref, namctl, IOSTAT = ios, ERR = 901 )
146901   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namctl in reference namelist', .TRUE. )
147
148      REWIND( numnam_cfg )              ! Namelist namctl in confguration namelist : Control prints & Benchmark
149      READ  ( numnam_cfg, namctl, IOSTAT = ios, ERR = 902 )
150902   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namctl in configuration namelist', .TRUE. )
151
152      !
153      REWIND( numnam_ref )              ! Namelist namcfg in reference namelist : Control prints & Benchmark
154      READ  ( numnam_ref, namcfg, IOSTAT = ios, ERR = 903 )
155903   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namcfg in reference namelist', .TRUE. )
156
157      REWIND( numnam_cfg )              ! Namelist namcfg in confguration namelist : Control prints & Benchmark
158      READ  ( numnam_cfg, namcfg, IOSTAT = ios, ERR = 904 )
159904   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namcfg in configuration namelist', .TRUE. )   
160
161      !                             !--------------------------------------------!
162      !                             !  set communicator & select the local node  !
163      !                             !  NB: mynode also opens output.namelist.dyn !
164      !                             !      on unit number numond on first proc   !
165      !                             !--------------------------------------------!
166#if defined key_iomput
167      IF( Agrif_Root() ) THEN
168# if defined key_oasis3 || defined key_oasis4
169         CALL cpl_prism_init( ilocal_comm )      ! nemo local communicator given by oasis
170         CALL xios_initialize( "oceanx",local_comm=ilocal_comm )
171# else
172         CALL  xios_initialize( "nemo",return_comm=ilocal_comm )
173# endif
174      ENDIF
175      narea = mynode( cltxt, numnam_ref, numnam_cfg, numond , nstop, ilocal_comm )   ! Nodes selection
176#else
177# if defined key_oasis3 || defined key_oasis4
178      IF( Agrif_Root() ) THEN
179         CALL cpl_prism_init( ilocal_comm )                 ! nemo local communicator given by oasis
180      ENDIF
181      narea = mynode( cltxt, numnam_ref, numnam_cfg, numond , nstop, ilocal_comm )   ! Nodes selection (control print return in cltxt)
182# else
183      ilocal_comm = 0
184      narea = mynode( cltxt, numnam_ref, numnam_cfg, numond , nstop )                ! Nodes selection (control print return in cltxt)
185# endif
186#endif
187      narea = narea + 1                                     ! mynode return the rank of proc (0 --> jpnij -1 )
188
189      lwm = (narea == 1)                                    ! control of output namelists
190      lwp = (narea == 1) .OR. ln_ctl                        ! control of all listing output print
191
192      IF(lwm) THEN
193         ! write merged namelists from earlier to output namelist now that the
194         ! file has been opened in call to mynode. nammpp has already been
195         ! written in mynode (if lk_mpp_mpi)
196         WRITE( numond, namctl )
197         WRITE( numond, namcfg )
198      ENDIF
199
200      ! If dimensions of processor grid weren't specified in the namelist file
201      ! then we calculate them here now that we have our communicator size
202      IF( (jpni < 1) .OR. (jpnj < 1) )THEN
203#if   defined key_mpp_mpi
204         IF( Agrif_Root() ) CALL nemo_partition(mppsize)
205#else
206         jpni  = 1
207         jpnj  = 1
208         jpnij = jpni*jpnj
209#endif
210      END IF
211
212      ! Calculate domain dimensions given calculated jpni and jpnj
213      ! This used to be done in par_oce.F90 when they were parameters rather
214      ! than variables
215      IF( Agrif_Root() ) THEN
216#if defined key_nemocice_decomp
217         jpi = ( nx_global+2-2*jpreci + (jpni-1) ) / jpni + 2*jpreci ! first  dim.
218         jpj = ( ny_global+2-2*jprecj + (jpnj-1) ) / jpnj + 2*jprecj ! second dim.
219#else
220         jpi = ( jpiglo-2*jpreci + (jpni-1) ) / jpni + 2*jpreci   ! first  dim.
221         jpj = ( jpjglo-2*jprecj + (jpnj-1) ) / jpnj + 2*jprecj   ! second dim.
222#endif
223         jpk = jpkdta                                             ! third dim
224         jpim1 = jpi-1                                            ! inner domain indices
225         jpjm1 = jpj-1                                            !   "           "
226         jpkm1 = jpk-1                                            !   "           "
227         jpij  = jpi*jpj                                          !  jpi x j
228      ENDIF
229
230      IF(lwp) THEN                            ! open listing units
231         !
232         CALL ctl_opn( numout, 'ocean.output', 'REPLACE', 'FORMATTED', 'SEQUENTIAL', -1, 6, .FALSE., narea )
233         !
234         WRITE(numout,*)
235         WRITE(numout,*) '   CNRS - NERC - Met OFFICE - MERCATOR-ocean - INGV - CMCC'
236         WRITE(numout,*) '                       NEMO team'
237         WRITE(numout,*) '            Ocean General Circulation Model'
238         WRITE(numout,*) '                  version 3.4  (2011) '
239         WRITE(numout,*)
240         WRITE(numout,*)
241         DO ji = 1, SIZE(cltxt)
242            IF( TRIM(cltxt(ji)) /= '' )   WRITE(numout,*) cltxt(ji)      ! control print of mynode
243         END DO
244         WRITE(numout,cform_aaa)                                         ! Flag AAAAAAA
245         !
246      ENDIF
247
248      ! Now we know the dimensions of the grid and numout has been set we can
249      ! allocate arrays
250      CALL nemo_alloc()
251
252      !                             !-------------------------------!
253      !                             !  NEMO general initialization  !
254      !                             !-------------------------------!
255
256      CALL nemo_ctl                          ! Control prints & Benchmark
257
258      !                                      ! Domain decomposition
259      IF( jpni*jpnj == jpnij ) THEN   ;   CALL mpp_init      ! standard cutting out
260      ELSE                            ;   CALL mpp_init2     ! eliminate land processors
261      ENDIF
262      !
263      IF( nn_timing == 1 )  CALL timing_init
264      !
265      !                                      ! General initialization
266                            CALL     phy_cst    ! Physical constants
267                            CALL     eos_init   ! Equation of state
268                            CALL     dom_cfg    ! Domain configuration
269                            CALL     dom_init   ! Domain
270
271      IF( ln_nnogather )    CALL nemo_northcomms   ! Initialise the northfold neighbour lists (must be done after the masks are defined)
272
273      IF( ln_ctl        )   CALL prt_ctl_init   ! Print control
274
275                            CALL  istate_init   ! ocean initial state (Dynamics and tracers)
276
277      IF( lk_diaobs     ) THEN                  ! Observation & model comparison
278                            CALL dia_obs_init            ! Initialize observational data
279                            CALL dia_obs( nit000 - 1 )   ! Observation operator for restart
280      ENDIF
281   END SUBROUTINE nemo_init
282
283
284   SUBROUTINE nemo_ctl
285      !!----------------------------------------------------------------------
286      !!                     ***  ROUTINE nemo_ctl  ***
287      !!
288      !! ** Purpose :   control print setting
289      !!
290      !! ** Method  : - print namctl information and check some consistencies
291      !!----------------------------------------------------------------------
292      !
293      IF(lwp) THEN                  ! control print
294         WRITE(numout,*)
295         WRITE(numout,*) 'nemo_ctl: Control prints & Benchmark'
296         WRITE(numout,*) '~~~~~~~ '
297         WRITE(numout,*) '   Namelist namctl'
298         WRITE(numout,*) '      run control (for debugging)     ln_ctl     = ', ln_ctl
299         WRITE(numout,*) '      level of print                  nn_print   = ', nn_print
300         WRITE(numout,*) '      Start i indice for SUM control  nn_ictls   = ', nn_ictls
301         WRITE(numout,*) '      End i indice for SUM control    nn_ictle   = ', nn_ictle
302         WRITE(numout,*) '      Start j indice for SUM control  nn_jctls   = ', nn_jctls
303         WRITE(numout,*) '      End j indice for SUM control    nn_jctle   = ', nn_jctle
304         WRITE(numout,*) '      number of proc. following i     nn_isplt   = ', nn_isplt
305         WRITE(numout,*) '      number of proc. following j     nn_jsplt   = ', nn_jsplt
306         WRITE(numout,*) '      benchmark parameter (0/1)       nn_bench   = ', nn_bench
307         WRITE(numout,*) '      timing activated    (0/1)       nn_timing  = ', nn_timing
308      ENDIF
309      !
310      nprint    = nn_print          ! convert DOCTOR namelist names into OLD names
311      nictls    = nn_ictls
312      nictle    = nn_ictle
313      njctls    = nn_jctls
314      njctle    = nn_jctle
315      isplt     = nn_isplt
316      jsplt     = nn_jsplt
317      nbench    = nn_bench
318      !                             ! Parameter control
319      !
320      IF( ln_ctl ) THEN                 ! sub-domain area indices for the control prints
321         IF( lk_mpp .AND. jpnij > 1 ) THEN
322            isplt = jpni   ;   jsplt = jpnj   ;   ijsplt = jpni*jpnj   ! the domain is forced to the real split domain
323         ELSE
324            IF( isplt == 1 .AND. jsplt == 1  ) THEN
325               CALL ctl_warn( ' - isplt & jsplt are equal to 1',   &
326                  &           ' - the print control will be done over the whole domain' )
327            ENDIF
328            ijsplt = isplt * jsplt            ! total number of processors ijsplt
329         ENDIF
330         IF(lwp) WRITE(numout,*)'          - The total number of processors over which the'
331         IF(lwp) WRITE(numout,*)'            print control will be done is ijsplt : ', ijsplt
332         !
333         !                              ! indices used for the SUM control
334         IF( nictls+nictle+njctls+njctle == 0 )   THEN    ! print control done over the default area
335            lsp_area = .FALSE.
336         ELSE                                             ! print control done over a specific  area
337            lsp_area = .TRUE.
338            IF( nictls < 1 .OR. nictls > jpiglo )   THEN
339               CALL ctl_warn( '          - nictls must be 1<=nictls>=jpiglo, it is forced to 1' )
340               nictls = 1
341            ENDIF
342            IF( nictle < 1 .OR. nictle > jpiglo )   THEN
343               CALL ctl_warn( '          - nictle must be 1<=nictle>=jpiglo, it is forced to jpiglo' )
344               nictle = jpiglo
345            ENDIF
346            IF( njctls < 1 .OR. njctls > jpjglo )   THEN
347               CALL ctl_warn( '          - njctls must be 1<=njctls>=jpjglo, it is forced to 1' )
348               njctls = 1
349            ENDIF
350            IF( njctle < 1 .OR. njctle > jpjglo )   THEN
351               CALL ctl_warn( '          - njctle must be 1<=njctle>=jpjglo, it is forced to jpjglo' )
352               njctle = jpjglo
353            ENDIF
354         ENDIF
355      ENDIF
356      !
357      IF( nbench == 1 ) THEN              ! Benchmark
358         SELECT CASE ( cp_cfg )
359         CASE ( 'gyre' )   ;   CALL ctl_warn( ' The Benchmark is activated ' )
360         CASE DEFAULT      ;   CALL ctl_stop( ' The Benchmark is based on the GYRE configuration:',   &
361            &                                 ' key_gyre must be used or set nbench = 0' )
362         END SELECT
363      ENDIF
364      !
365      IF( lk_c1d .AND. .NOT.lk_iomput )   CALL ctl_stop( 'nemo_ctl: The 1D configuration must be used ',   &
366         &                                               'with the IOM Input/Output manager. '         ,   &
367         &                                               'Compile with key_iomput enabled' )
368      !
369      IF( 1_wp /= SIGN(1._wp,-0._wp)  )   CALL ctl_stop( 'nemo_ctl: The intrinsec SIGN function follows ',  &
370         &                                               'f2003 standard. '                              ,  &
371         &                                               'Compile with key_nosignedzero enabled' )
372      !
373   END SUBROUTINE nemo_ctl
374
375
376   SUBROUTINE nemo_closefile
377      !!----------------------------------------------------------------------
378      !!                     ***  ROUTINE nemo_closefile  ***
379      !!
380      !! ** Purpose :   Close the files
381      !!----------------------------------------------------------------------
382      !
383      IF( lk_mpp )   CALL mppsync
384      !
385      CALL iom_close                                 ! close all input/output files managed by iom_*
386      !
387      IF( numstp      /= -1 )   CLOSE( numstp      )   ! time-step file
388      IF( numsol      /= -1 )   CLOSE( numsol      )   ! solver file
389      IF( numnam      /= -1 )   CLOSE( numnam      )   ! oce namelist
390      IF( numnam_ice  /= -1 )   CLOSE( numnam_ice  )   ! ice namelist
391      IF( numevo_ice  /= -1 )   CLOSE( numevo_ice  )   ! ice variables (temp. evolution)
392      IF( numout      /=  6 )   CLOSE( numout      )   ! standard model output file
393      IF( numdct_vol  /= -1 )   CLOSE( numdct_vol  )   ! volume transports
394      IF( numdct_heat /= -1 )   CLOSE( numdct_heat )   ! heat transports
395      IF( numdct_salt /= -1 )   CLOSE( numdct_salt )   ! salt transports
396
397      !
398      numout = 6                                     ! redefine numout in case it is used after this point...
399      !
400   END SUBROUTINE nemo_closefile
401
402
403   SUBROUTINE nemo_alloc
404      !!----------------------------------------------------------------------
405      !!                     ***  ROUTINE nemo_alloc  ***
406      !!
407      !! ** Purpose :   Allocate all the dynamic arrays of the OPA modules
408      !!
409      !! ** Method  :
410      !!----------------------------------------------------------------------
411      USE diawri    , ONLY: dia_wri_alloc
412      USE dom_oce   , ONLY: dom_oce_alloc
413      !
414      INTEGER :: ierr
415      !!----------------------------------------------------------------------
416      !
417      ierr =        oce_alloc       ()          ! ocean
418      ierr = ierr + dia_wri_alloc   ()
419      ierr = ierr + dom_oce_alloc   ()          ! ocean domain
420      !
421      ierr = ierr + lib_mpp_alloc   (numout)    ! mpp exchanges
422      !
423      IF( lk_mpp    )   CALL mpp_sum( ierr )
424      IF( ierr /= 0 )   CALL ctl_stop( 'STOP', 'nemo_alloc : unable to allocate standard ocean arrays' )
425      !
426   END SUBROUTINE nemo_alloc
427
428
429   SUBROUTINE nemo_partition( num_pes )
430      !!----------------------------------------------------------------------
431      !!                 ***  ROUTINE nemo_partition  ***
432      !!
433      !! ** Purpose :
434      !!
435      !! ** Method  :
436      !!----------------------------------------------------------------------
437      INTEGER, INTENT(in) :: num_pes ! The number of MPI processes we have
438      !
439      INTEGER, PARAMETER :: nfactmax = 20
440      INTEGER :: nfact ! The no. of factors returned
441      INTEGER :: ierr  ! Error flag
442      INTEGER :: ji
443      INTEGER :: idiff, mindiff, imin ! For choosing pair of factors that are closest in value
444      INTEGER, DIMENSION(nfactmax) :: ifact ! Array of factors
445      !!----------------------------------------------------------------------
446
447      ierr = 0
448
449      CALL factorise( ifact, nfactmax, nfact, num_pes, ierr )
450
451      IF( nfact <= 1 ) THEN
452         WRITE (numout, *) 'WARNING: factorisation of number of PEs failed'
453         WRITE (numout, *) '       : using grid of ',num_pes,' x 1'
454         jpnj = 1
455         jpni = num_pes
456      ELSE
457         ! Search through factors for the pair that are closest in value
458         mindiff = 1000000
459         imin    = 1
460         DO ji = 1, nfact-1, 2
461            idiff = ABS( ifact(ji) - ifact(ji+1) )
462            IF( idiff < mindiff ) THEN
463               mindiff = idiff
464               imin = ji
465            ENDIF
466         END DO
467         jpnj = ifact(imin)
468         jpni = ifact(imin + 1)
469      ENDIF
470      !
471      jpnij = jpni*jpnj
472      !
473   END SUBROUTINE nemo_partition
474
475
476   SUBROUTINE factorise( kfax, kmaxfax, knfax, kn, kerr )
477      !!----------------------------------------------------------------------
478      !!                     ***  ROUTINE factorise  ***
479      !!
480      !! ** Purpose :   return the prime factors of n.
481      !!                knfax factors are returned in array kfax which is of
482      !!                maximum dimension kmaxfax.
483      !! ** Method  :
484      !!----------------------------------------------------------------------
485      INTEGER                    , INTENT(in   ) ::   kn, kmaxfax
486      INTEGER                    , INTENT(  out) ::   kerr, knfax
487      INTEGER, DIMENSION(kmaxfax), INTENT(  out) ::   kfax
488      !
489      INTEGER :: ifac, jl, inu
490      INTEGER, PARAMETER :: ntest = 14
491      INTEGER :: ilfax(ntest)
492
493      ! lfax contains the set of allowed factors.
494      data (ilfax(jl),jl=1,ntest) / 16384, 8192, 4096, 2048, 1024, 512, 256,  &
495         &                            128,   64,   32,   16,    8,   4,   2  /
496      !!----------------------------------------------------------------------
497
498      ! Clear the error flag and initialise output vars
499      kerr = 0
500      kfax = 1
501      knfax = 0
502
503      ! Find the factors of n.
504      IF( kn == 1 )   GOTO 20
505
506      ! nu holds the unfactorised part of the number.
507      ! knfax holds the number of factors found.
508      ! l points to the allowed factor list.
509      ! ifac holds the current factor.
510
511      inu   = kn
512      knfax = 0
513
514      DO jl = ntest, 1, -1
515         !
516         ifac = ilfax(jl)
517         IF( ifac > inu )   CYCLE
518
519         ! Test whether the factor will divide.
520
521         IF( MOD(inu,ifac) == 0 ) THEN
522            !
523            knfax = knfax + 1            ! Add the factor to the list
524            IF( knfax > kmaxfax ) THEN
525               kerr = 6
526               write (*,*) 'FACTOR: insufficient space in factor array ', knfax
527               return
528            ENDIF
529            kfax(knfax) = ifac
530            ! Store the other factor that goes with this one
531            knfax = knfax + 1
532            kfax(knfax) = inu / ifac
533            !WRITE (*,*) 'ARPDBG, factors ',knfax-1,' & ',knfax,' are ', kfax(knfax-1),' and ',kfax(knfax)
534         ENDIF
535         !
536      END DO
537
538   20 CONTINUE      ! Label 20 is the exit point from the factor search loop.
539      !
540   END SUBROUTINE factorise
541
542#if defined key_mpp_mpi
543   SUBROUTINE nemo_northcomms
544      !!======================================================================
545      !!                     ***  ROUTINE  nemo_northcomms  ***
546      !! nemo_northcomms    :  Setup for north fold exchanges with explicit peer to peer messaging
547      !!=====================================================================
548      !!----------------------------------------------------------------------
549      !!
550      !! ** Purpose :   Initialization of the northern neighbours lists.
551      !!----------------------------------------------------------------------
552      !!    1.0  ! 2011-10  (A. C. Coward, NOCS & J. Donners, PRACE)
553      !!----------------------------------------------------------------------
554
555      INTEGER ::   ji, jj, jk, ij, jtyp    ! dummy loop indices
556      INTEGER ::   ijpj                    ! number of rows involved in north-fold exchange
557      INTEGER ::   northcomms_alloc        ! allocate return status
558      REAL(wp), ALLOCATABLE, DIMENSION ( :,: ) ::   znnbrs     ! workspace
559      LOGICAL,  ALLOCATABLE, DIMENSION ( : )   ::   lrankset   ! workspace
560
561      IF(lwp) WRITE(numout,*)
562      IF(lwp) WRITE(numout,*) 'nemo_northcomms : Initialization of the northern neighbours lists'
563      IF(lwp) WRITE(numout,*) '~~~~~~~~~~'
564
565      !!----------------------------------------------------------------------
566      ALLOCATE( znnbrs(jpi,jpj), stat = northcomms_alloc )
567      ALLOCATE( lrankset(jpnij), stat = northcomms_alloc )
568      IF( northcomms_alloc /= 0 ) THEN
569         WRITE(numout,cform_war)
570         WRITE(numout,*) 'northcomms_alloc : failed to allocate arrays'
571         CALL ctl_stop( 'STOP', 'nemo_northcomms : unable to allocate temporary arrays' )
572      ENDIF
573      nsndto = 0
574      isendto = -1
575      ijpj   = 4
576      !
577      ! This routine has been called because ln_nnogather has been set true ( nammpp )
578      ! However, these first few exchanges have to use the mpi_allgather method to
579      ! establish the neighbour lists to use in subsequent peer to peer exchanges.
580      ! Consequently, set l_north_nogather to be false here and set it true only after
581      ! the lists have been established.
582      !
583      l_north_nogather = .FALSE.
584      !
585      ! Exchange and store ranks on northern rows
586
587      DO jtyp = 1,4
588
589         lrankset = .FALSE.
590         znnbrs = narea
591         SELECT CASE (jtyp)
592            CASE(1)
593               CALL lbc_lnk( znnbrs, 'T', 1. )      ! Type 1: T,W-points
594            CASE(2)
595               CALL lbc_lnk( znnbrs, 'U', 1. )      ! Type 2: U-point
596            CASE(3)
597               CALL lbc_lnk( znnbrs, 'V', 1. )      ! Type 3: V-point
598            CASE(4)
599               CALL lbc_lnk( znnbrs, 'F', 1. )      ! Type 4: F-point
600         END SELECT
601
602         IF ( njmppt(narea) .EQ. MAXVAL( njmppt ) ) THEN
603            DO jj = nlcj-ijpj+1, nlcj
604               ij = jj - nlcj + ijpj
605               DO ji = 1,jpi
606                  IF ( INT(znnbrs(ji,jj)) .NE. 0 .AND. INT(znnbrs(ji,jj)) .NE. narea ) &
607               &     lrankset(INT(znnbrs(ji,jj))) = .true.
608               END DO
609            END DO
610
611            DO jj = 1,jpnij
612               IF ( lrankset(jj) ) THEN
613                  nsndto(jtyp) = nsndto(jtyp) + 1
614                  IF ( nsndto(jtyp) .GT. jpmaxngh ) THEN
615                     CALL ctl_stop( ' Too many neighbours in nemo_northcomms ', &
616                  &                 ' jpmaxngh will need to be increased ')
617                  ENDIF
618                  isendto(nsndto(jtyp),jtyp) = jj-1   ! narea converted to MPI rank
619               ENDIF
620            END DO
621         ENDIF
622
623      END DO
624
625      !
626      ! Type 5: I-point
627      !
628      ! ICE point exchanges may involve some averaging. The neighbours list is
629      ! built up using two exchanges to ensure that the whole stencil is covered.
630      ! lrankset should not be reset between these 'J' and 'K' point exchanges
631
632      jtyp = 5
633      lrankset = .FALSE.
634      znnbrs = narea
635      CALL lbc_lnk( znnbrs, 'J', 1. ) ! first ice U-V point
636
637      IF ( njmppt(narea) .EQ. MAXVAL( njmppt ) ) THEN
638         DO jj = nlcj-ijpj+1, nlcj
639            ij = jj - nlcj + ijpj
640            DO ji = 1,jpi
641               IF ( INT(znnbrs(ji,jj)) .NE. 0 .AND. INT(znnbrs(ji,jj)) .NE. narea ) &
642            &     lrankset(INT(znnbrs(ji,jj))) = .true.
643         END DO
644        END DO
645      ENDIF
646
647      znnbrs = narea
648      CALL lbc_lnk( znnbrs, 'K', 1. ) ! second ice U-V point
649
650      IF ( njmppt(narea) .EQ. MAXVAL( njmppt )) THEN
651         DO jj = nlcj-ijpj+1, nlcj
652            ij = jj - nlcj + ijpj
653            DO ji = 1,jpi
654               IF ( INT(znnbrs(ji,jj)) .NE. 0 .AND.  INT(znnbrs(ji,jj)) .NE. narea ) &
655            &       lrankset( INT(znnbrs(ji,jj))) = .true.
656            END DO
657         END DO
658
659         DO jj = 1,jpnij
660            IF ( lrankset(jj) ) THEN
661               nsndto(jtyp) = nsndto(jtyp) + 1
662               IF ( nsndto(jtyp) .GT. jpmaxngh ) THEN
663                  CALL ctl_stop( ' Too many neighbours in nemo_northcomms ', &
664               &                 ' jpmaxngh will need to be increased ')
665               ENDIF
666               isendto(nsndto(jtyp),jtyp) = jj-1   ! narea converted to MPI rank
667            ENDIF
668         END DO
669         !
670         ! For northern row areas, set l_north_nogather so that all subsequent exchanges
671         ! can use peer to peer communications at the north fold
672         !
673         l_north_nogather = .TRUE.
674         !
675      ENDIF
676      DEALLOCATE( znnbrs )
677      DEALLOCATE( lrankset )
678
679   END SUBROUTINE nemo_northcomms
680#else
681   SUBROUTINE nemo_northcomms      ! Dummy routine
682      WRITE(*,*) 'nemo_northcomms: You should not have seen this print! error?'
683   END SUBROUTINE nemo_northcomms
684#endif
685   !!======================================================================
686END MODULE nemogcm
Note: See TracBrowser for help on using the repository browser.