New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
nemogcm.F90 in branches/2012/dev_MERGE_2012/NEMOGCM/NEMO/OFF_SRC – NEMO

source: branches/2012/dev_MERGE_2012/NEMOGCM/NEMO/OFF_SRC/nemogcm.F90 @ 3827

Last change on this file since 3827 was 3827, checked in by cetlod, 11 years ago

v3.5alpha: Add missing IF(lwp) before writing in numout, see ticket #1066

File size: 28.6 KB
Line 
1MODULE nemogcm
2   !!======================================================================
3   !!                       ***  MODULE nemogcm   ***
4   !! Off-line Ocean   : passive tracer evolution, dynamics read in files
5   !!======================================================================
6   !! History :  3.3  ! 2010-05  (C. Ethe)  Full reorganization of the off-line: phasing with the on-line
7   !!            4.0  ! 2011-01  (C. Ethe, A. R. Porter, STFC Daresbury) dynamical allocation
8   !!----------------------------------------------------------------------
9
10   !!----------------------------------------------------------------------
11   !!   nemo_gcm        : off-line: solve ocean tracer only
12   !!   nemo_init       : initialization of the nemo model
13   !!   nemo_ctl        : initialisation of algorithm flag
14   !!   nemo_closefile  : close remaining files
15   !!----------------------------------------------------------------------
16   USE dom_oce         ! ocean space domain variables
17   USE oce             ! dynamics and tracers variables
18   USE c1d             ! 1D configuration
19   USE domcfg          ! domain configuration               (dom_cfg routine)
20   USE domain          ! domain initialization             (dom_init routine)
21   USE istate          ! initial state setting          (istate_init routine)
22   USE eosbn2          ! equation of state            (eos bn2 routine)
23   !              ! ocean physics
24   USE ldftra          ! lateral diffusivity setting    (ldf_tra_init routine)
25   USE ldfslp          ! slopes of neutral surfaces     (ldf_slp_init routine)
26   USE traqsr          ! solar radiation penetration    (tra_qsr_init routine)
27   USE trabbl          ! bottom boundary layer          (tra_bbl_init routine)
28   USE zdfini          ! vertical physics: initialization
29   USE sbcmod          ! surface boundary condition       (sbc_init     routine)
30   USE phycst          ! physical constant                  (par_cst routine)
31   USE dtadyn          ! Lecture and Interpolation of the dynamical fields
32   USE trcini          ! Initilization of the passive tracers
33   USE daymod          ! calendar                         (day     routine)
34   USE trcstp          ! passive tracer time-stepping      (trc_stp routine)
35   USE dtadyn          ! Lecture and interpolation of the dynamical fields
36   USE stpctl          ! time stepping control            (stp_ctl routine)
37   !              ! I/O & MPP
38   USE iom             ! I/O library
39   USE in_out_manager  ! I/O manager
40   USE mppini          ! shared/distributed memory setting (mpp_init routine)
41   USE lib_mpp         ! distributed memory computing
42#if defined key_iomput
43   USE xios
44#endif
45   USE prtctl          ! Print control                    (prt_ctl_init routine)
46   USE timing          ! Timing
47   USE lib_fortran     ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined)
48
49   IMPLICIT NONE
50   PRIVATE
51   
52   PUBLIC   nemo_gcm   ! called by nemo.F90
53
54   CHARACTER (len=64) ::   cform_aaa="( /, 'AAAAAAAA', / ) "   ! flag for output listing
55
56   !!----------------------------------------------------------------------
57   !! NEMO/OFF 3.3 , NEMO Consortium (2010)
58   !! $Id: nemogcm.F90 2528 2010-12-27 17:33:53Z rblod $
59   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
60   !!----------------------------------------------------------------------
61CONTAINS
62
63   SUBROUTINE nemo_gcm
64      !!----------------------------------------------------------------------
65      !!                     ***  ROUTINE nemo_gcm  ***
66      !!
67      !! ** Purpose :   nemo solves the primitive equations on an orthogonal
68      !!      curvilinear mesh on the sphere.
69      !!
70      !! ** Method  : - model general initialization
71      !!              - launch the time-stepping (dta_dyn and trc_stp)
72      !!              - finalize the run by closing files and communications
73      !!
74      !! References : Madec, Delecluse,Imbard, and Levy, 1997:  internal report, IPSL.
75      !!              Madec, 2008, internal report, IPSL.
76      !!----------------------------------------------------------------------
77      INTEGER :: istp, indic       ! time step index
78      !!----------------------------------------------------------------------
79
80      CALL nemo_init  ! Initializations
81
82      ! check that all process are still there... If some process have an error,
83      ! they will never enter in step and other processes will wait until the end of the cpu time!
84      IF( lk_mpp )   CALL mpp_max( nstop )
85
86      !                            !-----------------------!
87      !                            !==   time stepping   ==!
88      !                            !-----------------------!
89      istp = nit000
90      !
91      CALL iom_init            ! iom_put initialization (must be done after nemo_init for AGRIF+XIOS+OASIS)
92      !
93      DO WHILE ( istp <= nitend .AND. nstop == 0 )    ! time stepping
94         !
95         IF( istp /= nit000 )   CALL day      ( istp )         ! Calendar (day was already called at nit000 in day_init)
96                                CALL iom_setkt( istp - nit000 + 1 )         ! say to iom that we are at time step kstp
97                                CALL dta_dyn  ( istp )         ! Interpolation of the dynamical fields
98                                CALL trc_stp  ( istp )         ! time-stepping
99                                CALL stp_ctl  ( istp, indic )  ! Time loop: control and print
100         istp = istp + 1
101         IF( lk_mpp )   CALL mpp_max( nstop )
102      END DO
103#if defined key_iomput
104      CALL xios_context_finalize() ! needed for XIOS+AGRIF
105#endif
106
107      !                            !------------------------!
108      !                            !==  finalize the run  ==!
109      !                            !------------------------!
110      IF(lwp) WRITE(numout,cform_aaa)                 ! Flag AAAAAAA
111
112      IF( nstop /= 0 .AND. lwp ) THEN                 ! error print
113         WRITE(numout,cform_err)
114         WRITE(numout,*) nstop, ' error have been found'
115      ENDIF
116      !
117      IF( nn_timing == 1 )   CALL timing_finalize
118      !
119      CALL nemo_closefile
120      !
121# if defined key_iomput
122      CALL xios_finalize             ! end mpp communications
123# else
124      IF( lk_mpp )   CALL mppstop       ! end mpp communications
125# endif
126      !
127   END SUBROUTINE nemo_gcm
128
129
130   SUBROUTINE nemo_init
131      !!----------------------------------------------------------------------
132      !!                     ***  ROUTINE nemo_init ***
133      !!
134      !! ** Purpose :   initialization of the nemo model in off-line mode
135      !!----------------------------------------------------------------------
136      INTEGER ::   ji            ! dummy loop indices
137      INTEGER ::   ilocal_comm   ! local integer
138      CHARACTER(len=80), DIMENSION(16) ::   cltxt
139      !!
140      NAMELIST/namctl/ ln_ctl  , nn_print, nn_ictls, nn_ictle,   &
141         &             nn_isplt, nn_jsplt, nn_jctls, nn_jctle,   &
142         &             nn_bench, nn_timing
143      !!----------------------------------------------------------------------
144      !
145      cltxt = ''
146      !
147      !                             ! open Namelist file
148      CALL ctl_opn( numnam, 'namelist', 'OLD', 'FORMATTED', 'SEQUENTIAL', -1, 6, .FALSE. )
149      !
150      READ( numnam, namctl )        ! Namelist namctl : Control prints & Benchmark
151      !
152      !                             !--------------------------------------------!
153      !                             !  set communicator & select the local node  !
154      !                             !--------------------------------------------!
155#if defined key_iomput
156         CALL  xios_initialize( "nemo",return_comm=ilocal_comm )
157      narea = mynode( cltxt, numnam, nstop, ilocal_comm )   ! Nodes selection
158#else
159      ilocal_comm = 0
160      narea = mynode( cltxt, numnam, nstop )                 ! Nodes selection (control print return in cltxt)
161#endif
162
163      narea = narea + 1                       ! mynode return the rank of proc (0 --> jpnij -1 )
164
165      lwp = (narea == 1) .OR. ln_ctl          ! control of all listing output print
166
167      ! If dimensions of processor grid weren't specified in the namelist file
168      ! then we calculate them here now that we have our communicator size
169      IF( (jpni < 1) .OR. (jpnj < 1) )THEN
170#if   defined key_mpp_mpi
171         CALL nemo_partition(mppsize)
172#else
173         jpni = 1
174         jpnj = 1
175         jpnij = jpni*jpnj
176#endif
177      END IF
178
179      ! Calculate domain dimensions given calculated jpni and jpnj
180      ! This used to be done in par_oce.F90 when they were parameters rather
181      ! than variables
182      jpi = ( jpiglo-2*jpreci + (jpni-1) ) / jpni + 2*jpreci   ! first  dim.
183      jpj = ( jpjglo-2*jprecj + (jpnj-1) ) / jpnj + 2*jprecj   ! second dim.
184      jpk = jpkdta                                             ! third dim
185      jpim1 = jpi-1                                            ! inner domain indices
186      jpjm1 = jpj-1                                            !   "           "
187      jpkm1 = jpk-1                                            !   "           "
188      jpij  = jpi*jpj                                          !  jpi x j
189
190
191      IF(lwp) THEN                            ! open listing units
192         !
193         CALL ctl_opn( numout, 'ocean.output', 'REPLACE', 'FORMATTED', 'SEQUENTIAL', -1, 6, .FALSE., narea )
194         !
195         WRITE(numout,*)
196         WRITE(numout,*) '   CNRS - NERC - Met OFFICE - MERCATOR-ocean - INGV - CMCC'
197         WRITE(numout,*) '                       NEMO team'
198         WRITE(numout,*) '            Ocean General Circulation Model'
199         WRITE(numout,*) '                  version 3.5  (2012) '
200         WRITE(numout,*)
201         WRITE(numout,*)
202         DO ji = 1, SIZE(cltxt) 
203            IF( TRIM(cltxt(ji)) /= '' )   WRITE(numout,*) cltxt(ji)      ! control print of mynode
204         END DO
205         WRITE(numout,cform_aaa)                                         ! Flag AAAAAAA
206         !
207      ENDIF
208
209      ! Now we know the dimensions of the grid and numout has been set we can
210      ! allocate arrays
211      CALL nemo_alloc()
212
213      !                             !--------------------------------!
214      !                             !  Model general initialization  !
215      !                             !--------------------------------!
216
217      CALL nemo_ctl                           ! Control prints & Benchmark
218
219      !                                      ! Domain decomposition
220      IF( jpni*jpnj == jpnij ) THEN   ;   CALL mpp_init      ! standard cutting out
221      ELSE                            ;   CALL mpp_init2     ! eliminate land processors
222      ENDIF
223      !
224      IF( nn_timing == 1 )  CALL timing_init
225      !
226
227      !                                      ! General initialization
228      IF( nn_timing == 1 )  CALL timing_start( 'nemo_init')
229      !
230                            CALL     phy_cst    ! Physical constants
231                            CALL     eos_init   ! Equation of state
232                            CALL     dom_cfg    ! Domain configuration
233                            CALL     dom_init   ! Domain
234                            CALL  istate_init   ! ocean initial state (Dynamics and tracers)
235
236      IF( ln_nnogather )    CALL nemo_northcomms   ! Initialise the northfold neighbour lists (must be done after the masks are defined)
237
238      IF( ln_ctl        )   CALL prt_ctl_init   ! Print control
239
240      !                                     ! Ocean physics
241                            CALL     sbc_init   ! Forcings : surface module
242#if ! defined key_degrad
243                            CALL ldf_tra_init   ! Lateral ocean tracer physics
244#endif
245      IF( lk_ldfslp )       CALL ldf_slp_init   ! slope of lateral mixing
246
247      !                                     ! Active tracers
248                            CALL tra_qsr_init   ! penetrative solar radiation qsr
249      IF( lk_trabbl     )   CALL tra_bbl_init   ! advective (and/or diffusive) bottom boundary layer scheme
250
251      !                                     ! Passive tracers
252                            CALL     trc_init   ! Passive tracers initialization
253      !                                     ! Dynamics
254                            CALL dta_dyn_init   ! Initialization for the dynamics
255
256      IF(lwp) WRITE(numout,cform_aaa)       ! Flag AAAAAAA
257      !
258      IF( nn_timing == 1 )  CALL timing_stop( 'nemo_init')
259      !
260   END SUBROUTINE nemo_init
261
262
263   SUBROUTINE nemo_ctl
264      !!----------------------------------------------------------------------
265      !!                     ***  ROUTINE nemo_ctl  ***
266      !!
267      !! ** Purpose :   control print setting
268      !!
269      !! ** Method  : - print namctl information and check some consistencies
270      !!----------------------------------------------------------------------
271      !
272      IF(lwp) THEN                  ! Parameter print
273         WRITE(numout,*)
274         WRITE(numout,*) 'nemo_flg: Control prints & Benchmark'
275         WRITE(numout,*) '~~~~~~~ '
276         WRITE(numout,*) '   Namelist namctl'
277         WRITE(numout,*) '      run control (for debugging)     ln_ctl     = ', ln_ctl
278         WRITE(numout,*) '      level of print                  nn_print   = ', nn_print
279         WRITE(numout,*) '      Start i indice for SUM control  nn_ictls   = ', nn_ictls
280         WRITE(numout,*) '      End i indice for SUM control    nn_ictle   = ', nn_ictle
281         WRITE(numout,*) '      Start j indice for SUM control  nn_jctls   = ', nn_jctls
282         WRITE(numout,*) '      End j indice for SUM control    nn_jctle   = ', nn_jctle
283         WRITE(numout,*) '      number of proc. following i     nn_isplt   = ', nn_isplt
284         WRITE(numout,*) '      number of proc. following j     nn_jsplt   = ', nn_jsplt
285         WRITE(numout,*) '      benchmark parameter (0/1)       nn_bench   = ', nn_bench
286      ENDIF
287      !
288      nprint    = nn_print          ! convert DOCTOR namelist names into OLD names
289      nictls    = nn_ictls
290      nictle    = nn_ictle
291      njctls    = nn_jctls
292      njctle    = nn_jctle
293      isplt     = nn_isplt
294      jsplt     = nn_jsplt
295      nbench    = nn_bench
296      !                             ! Parameter control
297      !
298      IF( ln_ctl ) THEN                 ! sub-domain area indices for the control prints
299         IF( lk_mpp ) THEN
300            isplt = jpni   ;   jsplt = jpnj   ;   ijsplt = jpni*jpnj   ! the domain is forced to the real splitted domain
301         ELSE
302            IF( isplt == 1 .AND. jsplt == 1  ) THEN
303               CALL ctl_warn( ' - isplt & jsplt are equal to 1',   &
304                  &           ' - the print control will be done over the whole domain' )
305            ENDIF
306            ijsplt = isplt * jsplt            ! total number of processors ijsplt
307         ENDIF
308         IF(lwp) WRITE(numout,*)'          - The total number of processors over which the'
309         IF(lwp) WRITE(numout,*)'            print control will be done is ijsplt : ', ijsplt
310         !
311         !                              ! indices used for the SUM control
312         IF( nictls+nictle+njctls+njctle == 0 )   THEN    ! print control done over the default area
313            lsp_area = .FALSE.
314         ELSE                                             ! print control done over a specific  area
315            lsp_area = .TRUE.
316            IF( nictls < 1 .OR. nictls > jpiglo )   THEN
317               CALL ctl_warn( '          - nictls must be 1<=nictls>=jpiglo, it is forced to 1' )
318               nictls = 1
319            ENDIF
320            IF( nictle < 1 .OR. nictle > jpiglo )   THEN
321               CALL ctl_warn( '          - nictle must be 1<=nictle>=jpiglo, it is forced to jpiglo' )
322               nictle = jpiglo
323            ENDIF
324            IF( njctls < 1 .OR. njctls > jpjglo )   THEN
325               CALL ctl_warn( '          - njctls must be 1<=njctls>=jpjglo, it is forced to 1' )
326               njctls = 1
327            ENDIF
328            IF( njctle < 1 .OR. njctle > jpjglo )   THEN
329               CALL ctl_warn( '          - njctle must be 1<=njctle>=jpjglo, it is forced to jpjglo' )
330               njctle = jpjglo
331            ENDIF
332         ENDIF
333      ENDIF
334      !
335      IF( nbench == 1 )   THEN            ! Benchmark
336         SELECT CASE ( cp_cfg )
337         CASE ( 'gyre' )   ;   CALL ctl_warn( ' The Benchmark is activated ' )
338         CASE DEFAULT      ;   CALL ctl_stop( ' The Benchmark is based on the GYRE configuration:',   &
339            &                                 ' key_gyre must be used or set nbench = 0' )
340         END SELECT
341      ENDIF
342      !
343      IF( lk_c1d .AND. .NOT.lk_iomput )   CALL ctl_stop( 'nemo_ctl: The 1D configuration must be used ',   &
344         &                                               'with the IOM Input/Output manager. '        ,   &
345         &                                               'Compile with key_iomput enabled' )
346      !
347      IF( 1_wp /= SIGN(1._wp,-0._wp)  )   CALL ctl_stop( 'nemo_ctl: The intrinsec SIGN function follows ',  &
348         &                                               'f2003 standard. '                              ,  &
349         &                                               'Compile with key_nosignedzero enabled' )
350      !
351   END SUBROUTINE nemo_ctl
352
353
354   SUBROUTINE nemo_closefile
355      !!----------------------------------------------------------------------
356      !!                     ***  ROUTINE nemo_closefile  ***
357      !!
358      !! ** Purpose :   Close the files
359      !!----------------------------------------------------------------------
360      !
361      IF ( lk_mpp ) CALL mppsync
362      !
363      CALL iom_close                                 ! close all input/output files managed by iom_*
364      !
365      IF( numstp     /= -1 )   CLOSE( numstp     )   ! time-step file
366      IF( numnam     /= -1 )   CLOSE( numnam     )   ! oce namelist
367      IF( numout     /=  6 )   CLOSE( numout     )   ! standard model output file
368      numout = 6                                     ! redefine numout in case it is used after this point...
369      !
370   END SUBROUTINE nemo_closefile
371
372
373   SUBROUTINE nemo_alloc
374      !!----------------------------------------------------------------------
375      !!                     ***  ROUTINE nemo_alloc  ***
376      !!
377      !! ** Purpose :   Allocate all the dynamic arrays of the OPA modules
378      !!
379      !! ** Method  :
380      !!----------------------------------------------------------------------
381      USE diawri,       ONLY: dia_wri_alloc
382      USE dom_oce,      ONLY: dom_oce_alloc
383      USE zdf_oce,      ONLY: zdf_oce_alloc
384      USE ldftra_oce,   ONLY: ldftra_oce_alloc
385      USE trc_oce,      ONLY: trc_oce_alloc
386      !
387      INTEGER :: ierr
388      !!----------------------------------------------------------------------
389      !
390      ierr =        oce_alloc       ()          ! ocean
391      ierr = ierr + dia_wri_alloc   ()
392      ierr = ierr + dom_oce_alloc   ()          ! ocean domain
393      ierr = ierr + ldftra_oce_alloc()          ! ocean lateral  physics : tracers
394      ierr = ierr + zdf_oce_alloc   ()          ! ocean vertical physics
395      !
396      ierr = ierr + lib_mpp_alloc   (numout)    ! mpp exchanges
397      ierr = ierr + trc_oce_alloc   ()          ! shared TRC / TRA arrays
398      !
399      IF( lk_mpp    )   CALL mpp_sum( ierr )
400      IF( ierr /= 0 )   CALL ctl_stop( 'STOP', 'nemo_alloc: unable to allocate standard ocean arrays' )
401      !
402   END SUBROUTINE nemo_alloc
403
404
405   SUBROUTINE nemo_partition( num_pes )
406      !!----------------------------------------------------------------------
407      !!                 ***  ROUTINE nemo_partition  ***
408      !!
409      !! ** Purpose :   
410      !!
411      !! ** Method  :
412      !!----------------------------------------------------------------------
413      INTEGER, INTENT(in) :: num_pes ! The number of MPI processes we have
414      !
415      INTEGER, PARAMETER :: nfactmax = 20
416      INTEGER :: nfact ! The no. of factors returned
417      INTEGER :: ierr  ! Error flag
418      INTEGER :: ji
419      INTEGER :: idiff, mindiff, imin ! For choosing pair of factors that are closest in value
420      INTEGER, DIMENSION(nfactmax) :: ifact ! Array of factors
421      !!----------------------------------------------------------------------
422
423      ierr = 0
424
425      CALL factorise( ifact, nfactmax, nfact, num_pes, ierr )
426
427      IF( nfact <= 1 ) THEN
428         WRITE (numout, *) 'WARNING: factorisation of number of PEs failed'
429         WRITE (numout, *) '       : using grid of ',num_pes,' x 1'
430         jpnj = 1
431         jpni = num_pes
432      ELSE
433         ! Search through factors for the pair that are closest in value
434         mindiff = 1000000
435         imin    = 1
436         DO ji = 1, nfact-1, 2
437            idiff = ABS( ifact(ji) - ifact(ji+1) )
438            IF( idiff < mindiff ) THEN
439               mindiff = idiff
440               imin = ji
441            ENDIF
442         END DO
443         jpnj = ifact(imin)
444         jpni = ifact(imin + 1)
445      ENDIF
446      !
447      jpnij = jpni*jpnj
448      !
449   END SUBROUTINE nemo_partition
450
451
452   SUBROUTINE factorise( kfax, kmaxfax, knfax, kn, kerr )
453      !!----------------------------------------------------------------------
454      !!                     ***  ROUTINE factorise  ***
455      !!
456      !! ** Purpose :   return the prime factors of n.
457      !!                knfax factors are returned in array kfax which is of
458      !!                maximum dimension kmaxfax.
459      !! ** Method  :
460      !!----------------------------------------------------------------------
461      INTEGER                    , INTENT(in   ) ::   kn, kmaxfax
462      INTEGER                    , INTENT(  out) ::   kerr, knfax
463      INTEGER, DIMENSION(kmaxfax), INTENT(  out) ::   kfax
464      !
465      INTEGER :: ifac, jl, inu
466      INTEGER, PARAMETER :: ntest = 14
467      INTEGER :: ilfax(ntest)
468      !
469      ! lfax contains the set of allowed factors.
470      data (ilfax(jl),jl=1,ntest) / 16384, 8192, 4096, 2048, 1024, 512, 256,  &
471         &                            128,   64,   32,   16,    8,   4,   2  /
472      !!----------------------------------------------------------------------
473
474      ! Clear the error flag and initialise output vars
475      kerr = 0
476      kfax = 1
477      knfax = 0
478
479      ! Find the factors of n.
480      IF( kn == 1 )   GOTO 20
481
482      ! nu holds the unfactorised part of the number.
483      ! knfax holds the number of factors found.
484      ! l points to the allowed factor list.
485      ! ifac holds the current factor.
486
487      inu   = kn
488      knfax = 0
489
490      DO jl = ntest, 1, -1
491         !
492         ifac = ilfax(jl)
493         IF( ifac > inu )   CYCLE
494
495         ! Test whether the factor will divide.
496
497         IF( MOD(inu,ifac) == 0 ) THEN
498            !
499            knfax = knfax + 1            ! Add the factor to the list
500            IF( knfax > kmaxfax ) THEN
501               kerr = 6
502               write (*,*) 'FACTOR: insufficient space in factor array ', knfax
503               return
504            ENDIF
505            kfax(knfax) = ifac
506            ! Store the other factor that goes with this one
507            knfax = knfax + 1
508            kfax(knfax) = inu / ifac
509            !WRITE (*,*) 'ARPDBG, factors ',knfax-1,' & ',knfax,' are ', kfax(knfax-1),' and ',kfax(knfax)
510         ENDIF
511         !
512      END DO
513
514   20 CONTINUE      ! Label 20 is the exit point from the factor search loop.
515      !
516   END SUBROUTINE factorise
517
518#if defined key_mpp_mpi
519   SUBROUTINE nemo_northcomms
520      !!======================================================================
521      !!                     ***  ROUTINE  nemo_northcomms  ***
522      !! nemo_northcomms    :  Setup for north fold exchanges with explicit peer to peer messaging
523      !!=====================================================================
524      !!----------------------------------------------------------------------
525      !!
526      !! ** Purpose :   Initialization of the northern neighbours lists.
527      !!----------------------------------------------------------------------
528      !!    1.0  ! 2011-10  (A. C. Coward, NOCS & J. Donners, PRACE)
529      !!----------------------------------------------------------------------
530
531      INTEGER ::   ji, jj, jk, ij, jtyp    ! dummy loop indices
532      INTEGER ::   ijpj                    ! number of rows involved in north-fold exchange
533      INTEGER ::   northcomms_alloc        ! allocate return status
534      REAL(wp), ALLOCATABLE, DIMENSION ( :,: ) ::   znnbrs     ! workspace
535      LOGICAL,  ALLOCATABLE, DIMENSION ( : )   ::   lrankset   ! workspace
536
537      IF(lwp) WRITE(numout,*)
538      IF(lwp) WRITE(numout,*) 'nemo_northcomms : Initialization of the northern neighbours lists'
539      IF(lwp) WRITE(numout,*) '~~~~~~~~~~'
540
541      !!----------------------------------------------------------------------
542      ALLOCATE( znnbrs(jpi,jpj), stat = northcomms_alloc )
543      ALLOCATE( lrankset(jpnij), stat = northcomms_alloc )
544      IF( northcomms_alloc /= 0 ) THEN
545         WRITE(numout,cform_war)
546         WRITE(numout,*) 'northcomms_alloc : failed to allocate arrays'
547         CALL ctl_stop( 'STOP', 'nemo_northcomms : unable to allocate temporary arrays' )
548      ENDIF
549      nsndto = 0
550      isendto = -1
551      ijpj   = 4
552      !
553      ! This routine has been called because ln_nnogather has been set true ( nammpp )
554      ! However, these first few exchanges have to use the mpi_allgather method to
555      ! establish the neighbour lists to use in subsequent peer to peer exchanges.
556      ! Consequently, set l_north_nogather to be false here and set it true only after
557      ! the lists have been established.
558      !
559      l_north_nogather = .FALSE.
560      !
561      ! Exchange and store ranks on northern rows
562
563      DO jtyp = 1,4
564
565         lrankset = .FALSE.
566         znnbrs = narea
567         SELECT CASE (jtyp)
568            CASE(1)
569               CALL lbc_lnk( znnbrs, 'T', 1. )      ! Type 1: T,W-points
570            CASE(2)
571               CALL lbc_lnk( znnbrs, 'U', 1. )      ! Type 2: U-point
572            CASE(3)
573               CALL lbc_lnk( znnbrs, 'V', 1. )      ! Type 3: V-point
574            CASE(4)
575               CALL lbc_lnk( znnbrs, 'F', 1. )      ! Type 4: F-point
576         END SELECT
577
578         IF ( njmppt(narea) .EQ. MAXVAL( njmppt ) ) THEN
579            DO jj = nlcj-ijpj+1, nlcj
580               ij = jj - nlcj + ijpj
581               DO ji = 1,jpi
582                  IF ( INT(znnbrs(ji,jj)) .NE. 0 .AND. INT(znnbrs(ji,jj)) .NE. narea ) &
583               &     lrankset(INT(znnbrs(ji,jj))) = .true.
584               END DO
585            END DO
586
587            DO jj = 1,jpnij
588               IF ( lrankset(jj) ) THEN
589                  nsndto(jtyp) = nsndto(jtyp) + 1
590                  IF ( nsndto(jtyp) .GT. jpmaxngh ) THEN
591                     CALL ctl_stop( ' Too many neighbours in nemo_northcomms ', &
592                  &                 ' jpmaxngh will need to be increased ')
593                  ENDIF
594                  isendto(nsndto(jtyp),jtyp) = jj-1   ! narea converted to MPI rank
595               ENDIF
596            END DO
597         ENDIF
598
599      END DO
600
601      !
602      ! Type 5: I-point
603      !
604      ! ICE point exchanges may involve some averaging. The neighbours list is
605      ! built up using two exchanges to ensure that the whole stencil is covered.
606      ! lrankset should not be reset between these 'J' and 'K' point exchanges
607
608      jtyp = 5
609      lrankset = .FALSE.
610      znnbrs = narea 
611      CALL lbc_lnk( znnbrs, 'J', 1. ) ! first ice U-V point
612
613      IF ( njmppt(narea) .EQ. MAXVAL( njmppt ) ) THEN
614         DO jj = nlcj-ijpj+1, nlcj
615            ij = jj - nlcj + ijpj
616            DO ji = 1,jpi
617               IF ( INT(znnbrs(ji,jj)) .NE. 0 .AND. INT(znnbrs(ji,jj)) .NE. narea ) &
618            &     lrankset(INT(znnbrs(ji,jj))) = .true.
619         END DO
620        END DO
621      ENDIF
622
623      znnbrs = narea 
624      CALL lbc_lnk( znnbrs, 'K', 1. ) ! second ice U-V point
625
626      IF ( njmppt(narea) .EQ. MAXVAL( njmppt )) THEN
627         DO jj = nlcj-ijpj+1, nlcj
628            ij = jj - nlcj + ijpj
629            DO ji = 1,jpi
630               IF ( INT(znnbrs(ji,jj)) .NE. 0 .AND.  INT(znnbrs(ji,jj)) .NE. narea ) &
631            &       lrankset( INT(znnbrs(ji,jj))) = .true.
632            END DO
633         END DO
634
635         DO jj = 1,jpnij
636            IF ( lrankset(jj) ) THEN
637               nsndto(jtyp) = nsndto(jtyp) + 1
638               IF ( nsndto(jtyp) .GT. jpmaxngh ) THEN
639                  CALL ctl_stop( ' Too many neighbours in nemo_northcomms ', &
640               &                 ' jpmaxngh will need to be increased ')
641               ENDIF
642               isendto(nsndto(jtyp),jtyp) = jj-1   ! narea converted to MPI rank
643            ENDIF
644         END DO
645         !
646         ! For northern row areas, set l_north_nogather so that all subsequent exchanges
647         ! can use peer to peer communications at the north fold
648         !
649         l_north_nogather = .TRUE.
650         !
651      ENDIF
652      DEALLOCATE( znnbrs )
653      DEALLOCATE( lrankset )
654
655   END SUBROUTINE nemo_northcomms
656#else
657   SUBROUTINE nemo_northcomms      ! Dummy routine
658      WRITE(*,*) 'nemo_northcomms: You should not have seen this print! error?'
659   END SUBROUTINE nemo_northcomms
660#endif
661   !!======================================================================
662END MODULE nemogcm
Note: See TracBrowser for help on using the repository browser.