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.
mppini.F90 in NEMO/branches/2020/dev_r12558_HPC-08_epico_Extra_Halo/src/OCE/LBC – NEMO

source: NEMO/branches/2020/dev_r12558_HPC-08_epico_Extra_Halo/src/OCE/LBC/mppini.F90 @ 13174

Last change on this file since 13174 was 13174, checked in by smasson, 4 years ago

Extra_Halo: works if jpni = 1, allows nn_hls >2, remove island in BENCH, see #2366

  • Property svn:keywords set to Id
File size: 63.4 KB
Line 
1MODULE mppini
2   !!======================================================================
3   !!                       ***  MODULE mppini   ***
4   !! Ocean initialization : distributed memory computing initialization
5   !!======================================================================
6   !! History :  6.0  !  1994-11  (M. Guyon)  Original code
7   !!  OPA       7.0  !  1995-04  (J. Escobar, M. Imbard)
8   !!            8.0  !  1998-05  (M. Imbard, J. Escobar, L. Colombet )  SHMEM and MPI versions
9   !!  NEMO      1.0  !  2004-01  (G. Madec, J.M Molines)  F90 : free form , north fold jpni > 1
10   !!            3.4  !  2011-10  (A. C. Coward, NOCS & J. Donners, PRACE)  add init_nfdcom
11   !!            3.   !  2013-06  (I. Epicoco, S. Mocavero, CMCC)  init_nfdcom: setup avoiding MPI communication
12   !!            4.0  !  2016-06  (G. Madec)  use domain configuration file instead of bathymetry file
13   !!            4.0  !  2017-06  (J.M. Molines, T. Lovato) merge of mppini and mppini_2
14   !!----------------------------------------------------------------------
15
16   !!----------------------------------------------------------------------
17   !!  mpp_init       : Lay out the global domain over processors with/without land processor elimination
18   !!      init_ioipsl: IOIPSL initialization in mpp
19   !!      init_nfdcom: Setup for north fold exchanges with explicit point-to-point messaging
20   !!      init_doloop: set the starting/ending indices of DO-loop used in do_loop_substitute
21   !!----------------------------------------------------------------------
22   USE dom_oce        ! ocean space and time domain
23   USE bdy_oce        ! open BounDarY 
24   !
25   USE lbcnfd  , ONLY : isendto, nsndto ! Setup of north fold exchanges
26   USE lib_mpp        ! distribued memory computing library
27   USE iom            ! nemo I/O library
28   USE ioipsl         ! I/O IPSL library
29   USE in_out_manager ! I/O Manager
30
31   IMPLICIT NONE
32   PRIVATE
33
34   PUBLIC   mpp_init   ! called by opa.F90
35
36   INTEGER ::   numbot = -1   ! 'bottom_level' local logical unit
37   INTEGER ::   numbdy = -1   ! 'bdy_msk'      local logical unit
38   
39   !!----------------------------------------------------------------------
40   !! NEMO/OCE 4.0 , NEMO Consortium (2018)
41   !! $Id$
42   !! Software governed by the CeCILL license (see ./LICENSE)
43   !!----------------------------------------------------------------------
44CONTAINS
45
46#if ! defined key_mpp_mpi
47   !!----------------------------------------------------------------------
48   !!   Default option :                            shared memory computing
49   !!----------------------------------------------------------------------
50
51   SUBROUTINE mpp_init
52      !!----------------------------------------------------------------------
53      !!                  ***  ROUTINE mpp_init  ***
54      !!
55      !! ** Purpose :   Lay out the global domain over processors.
56      !!
57      !! ** Method  :   Shared memory computing, set the local processor
58      !!              variables to the value of the global domain
59      !!----------------------------------------------------------------------
60      !
61      jpiglo = Ni0glo
62      jpjglo = Nj0glo
63      jpimax = jpiglo
64      jpjmax = jpjglo
65      jpi    = jpiglo
66      jpj    = jpjglo
67      jpk    = jpkglo
68      jpim1  = jpi-1                         ! inner domain indices
69      jpjm1  = jpj-1                         !   "           "
70      jpkm1  = MAX( 1, jpk-1 )               !   "           "
71      !
72      CALL init_doloop                       ! set start/end indices or do-loop depending on the halo width value (nn_hls)
73      !
74      jpij   = jpi*jpj
75      jpni   = 1
76      jpnj   = 1
77      jpnij  = jpni*jpnj
78      nimpp  = 1           !
79      njmpp  = 1
80      nbondi = 2
81      nbondj = 2
82      nidom  = FLIO_DOM_NONE
83      npolj = 0
84      IF( jperio == 3 .OR. jperio == 4 )   npolj = 3
85      IF( jperio == 5 .OR. jperio == 6 )   npolj = 5
86      l_Iperio = jpni == 1 .AND. (jperio == 1 .OR. jperio == 4 .OR. jperio == 6 .OR. jperio == 7)
87      l_Jperio = jpnj == 1 .AND. (jperio == 2 .OR. jperio == 7)
88      !
89      IF(lwp) THEN
90         WRITE(numout,*)
91         WRITE(numout,*) 'mpp_init : NO massively parallel processing'
92         WRITE(numout,*) '~~~~~~~~ '
93         WRITE(numout,*) '   l_Iperio = ', l_Iperio, '    l_Jperio = ', l_Jperio 
94         WRITE(numout,*) '     npolj  = ',   npolj , '      njmpp  = ', njmpp
95      ENDIF
96      !
97      IF(  jpni /= 1 .OR. jpnj /= 1 .OR. jpnij /= 1 )                                     &
98         CALL ctl_stop( 'mpp_init: equality  jpni = jpnj = jpnij = 1 is not satisfied',   &
99            &           'the domain is lay out for distributed memory computing!' )
100         !
101   END SUBROUTINE mpp_init
102
103#else
104   !!----------------------------------------------------------------------
105   !!   'key_mpp_mpi'                     MPI massively parallel processing
106   !!----------------------------------------------------------------------
107
108
109   SUBROUTINE mpp_init
110      !!----------------------------------------------------------------------
111      !!                  ***  ROUTINE mpp_init  ***
112      !!                   
113      !! ** Purpose :   Lay out the global domain over processors.
114      !!      If land processors are to be eliminated, this program requires the
115      !!      presence of the domain configuration file. Land processors elimination
116      !!      is performed if jpni x jpnj /= jpnij. In this case, using the MPP_PREP
117      !!      preprocessing tool, help for defining the best cutting out.
118      !!
119      !! ** Method  :   Global domain is distributed in smaller local domains.
120      !!      Periodic condition is a function of the local domain position
121      !!      (global boundary or neighbouring domain) and of the global
122      !!      periodic
123      !!      Type :         jperio global periodic condition
124      !!
125      !! ** Action : - set domain parameters
126      !!                    nimpp     : longitudinal index
127      !!                    njmpp     : latitudinal  index
128      !!                    narea     : number for local area
129      !!                    nbondi    : mark for "east-west local boundary"
130      !!                    nbondj    : mark for "north-south local boundary"
131      !!                    nproc     : number for local processor
132      !!                    noea      : number for local neighboring processor
133      !!                    nowe      : number for local neighboring processor
134      !!                    noso      : number for local neighboring processor
135      !!                    nono      : number for local neighboring processor
136      !!----------------------------------------------------------------------
137      INTEGER ::   ji, jj, jn, jproc, jarea   ! dummy loop indices
138      INTEGER ::   inijmin
139      INTEGER ::   i2add
140      INTEGER ::   inum                       ! local logical unit
141      INTEGER ::   idir, ifreq, icont         ! local integers
142      INTEGER ::   ii, il1, ili, imil         !   -       -
143      INTEGER ::   ij, il2, ilj, ijm1         !   -       -
144      INTEGER ::   iino, ijno, iiso, ijso     !   -       -
145      INTEGER ::   iiea, ijea, iiwe, ijwe     !   -       -
146      INTEGER ::   iarea0                     !   -       -
147      INTEGER ::   ierr, ios                  !
148      INTEGER ::   inbi, inbj, iimax,  ijmax, icnt1, icnt2
149      LOGICAL ::   llbest, llauto
150      LOGICAL ::   llwrtlay
151      LOGICAL ::   ln_listonly
152      INTEGER, ALLOCATABLE, DIMENSION(:)     ::   iin, ii_nono, ii_noea          ! 1D workspace
153      INTEGER, ALLOCATABLE, DIMENSION(:)     ::   ijn, ii_noso, ii_nowe          !  -     -
154      INTEGER, ALLOCATABLE, DIMENSION(:,:) ::   iimppt, ijpi, ibondi, ipproc   ! 2D workspace
155      INTEGER, ALLOCATABLE, DIMENSION(:,:) ::   ijmppt, ijpj, ibondj, ipolj    !  -     -
156      INTEGER, ALLOCATABLE, DIMENSION(:,:) ::   iie0, iis0, iono, ioea         !  -     -
157      INTEGER, ALLOCATABLE, DIMENSION(:,:) ::   ije0, ijs0, ioso, iowe         !  -     -
158      LOGICAL, ALLOCATABLE, DIMENSION(:,:) ::   llisoce                        !  -     -
159      NAMELIST/nambdy/ ln_bdy, nb_bdy, ln_coords_file, cn_coords_file,           &
160           &             ln_mask_file, cn_mask_file, cn_dyn2d, nn_dyn2d_dta,     &
161           &             cn_dyn3d, nn_dyn3d_dta, cn_tra, nn_tra_dta,             & 
162           &             ln_tra_dmp, ln_dyn3d_dmp, rn_time_dmp, rn_time_dmp_out, &
163           &             cn_ice, nn_ice_dta,                                     &
164           &             ln_vol, nn_volctl, nn_rimwidth
165      NAMELIST/nammpp/ jpni, jpnj, nn_hls, ln_nnogather, ln_listonly
166      !!----------------------------------------------------------------------
167      !
168      llwrtlay = lwm .OR. sn_cfctl%l_layout
169      !
170      !  0. read namelists parameters
171      ! -----------------------------------
172      !
173      READ  ( numnam_ref, nammpp, IOSTAT = ios, ERR = 901 )
174901   IF( ios /= 0 )   CALL ctl_nam ( ios , 'nammpp in reference namelist' )
175      READ  ( numnam_cfg, nammpp, IOSTAT = ios, ERR = 902 )
176902   IF( ios >  0 )   CALL ctl_nam ( ios , 'nammpp in configuration namelist' )   
177      !
178      nn_hls = MAX(1, nn_hls)   ! nn_hls must be > 0
179      IF(lwp) THEN
180            WRITE(numout,*) '   Namelist nammpp'
181         IF( jpni < 1 .OR. jpnj < 1  ) THEN
182            WRITE(numout,*) '      jpni and jpnj will be calculated automatically'
183         ELSE
184            WRITE(numout,*) '      processor grid extent in i                            jpni = ', jpni
185            WRITE(numout,*) '      processor grid extent in j                            jpnj = ', jpnj
186         ENDIF
187            WRITE(numout,*) '      avoid use of mpi_allgather at the north fold  ln_nnogather = ', ln_nnogather
188      ENDIF
189      !
190      IF(lwm)   WRITE( numond, nammpp )
191      !
192!!!------------------------------------
193!!!  nn_hls shloud be read in nammpp
194!!!------------------------------------
195      jpiglo = Ni0glo + 2 * nn_hls
196      jpjglo = Nj0glo + 2 * nn_hls
197      !
198      ! do we need to take into account bdy_msk?
199      READ  ( numnam_ref, nambdy, IOSTAT = ios, ERR = 903)
200903   IF( ios /= 0 )   CALL ctl_nam ( ios , 'nambdy in reference namelist (mppini)' )
201      READ  ( numnam_cfg, nambdy, IOSTAT = ios, ERR = 904 )
202904   IF( ios >  0 )   CALL ctl_nam ( ios , 'nambdy in configuration namelist (mppini)' )
203      !
204      IF(               ln_read_cfg ) CALL iom_open( cn_domcfg,    numbot )
205      IF( ln_bdy .AND. ln_mask_file ) CALL iom_open( cn_mask_file, numbdy )
206      !
207      IF( ln_listonly )   CALL bestpartition( MAX(mppsize,jpni*jpnj), ldlist = .TRUE. )   ! must be done by all core
208      !
209      !  1. Dimension arrays for subdomains
210      ! -----------------------------------
211      !
212      ! If dimensions of processors grid weren't specified in the namelist file
213      ! then we calculate them here now that we have our communicator size
214      IF(lwp) THEN
215         WRITE(numout,*) 'mpp_init:'
216         WRITE(numout,*) '~~~~~~~~ '
217         WRITE(numout,*)
218      ENDIF
219      IF( jpni < 1 .OR. jpnj < 1 ) THEN
220         CALL bestpartition( mppsize, jpni, jpnj )           ! best mpi decomposition for mppsize mpi processes
221         llauto = .TRUE.
222         llbest = .TRUE.
223      ELSE
224         llauto = .FALSE.
225         CALL bestpartition( mppsize, inbi, inbj, icnt2 )    ! best mpi decomposition for mppsize mpi processes
226         ! largest subdomain size for mpi decoposition jpni*jpnj given in the namelist
227         CALL basic_decomposition( jpiglo, jpjglo, nn_hls, jpni, jpnj, jpimax, jpjmax )
228         ! largest subdomain size for mpi decoposition inbi*inbj given by bestpartition
229         CALL basic_decomposition( jpiglo, jpjglo, nn_hls, inbi, inbj,  iimax,  ijmax )
230         icnt1 = jpni*jpnj - mppsize   ! number of land subdomains that should be removed to use mppsize mpi processes
231         IF(lwp) THEN
232            WRITE(numout,9000) '   The chosen domain decomposition ', jpni, ' x ', jpnj, ' with ', icnt1, ' land subdomains'
233            WRITE(numout,9002) '      - uses a total of ',  mppsize,' mpi process'
234            WRITE(numout,9000) '      - has mpi subdomains with a maximum size of (jpi = ', jpimax, ', jpj = ', jpjmax,   &
235               &                                                                ', jpi*jpj = ', jpimax*jpjmax, ')'
236            WRITE(numout,9000) '   The best domain decompostion ', inbi, ' x ', inbj, ' with ', icnt2, ' land subdomains'
237            WRITE(numout,9002) '      - uses a total of ',  inbi*inbj-icnt2,' mpi process'
238            WRITE(numout,9000) '      - has mpi subdomains with a maximum size of (jpi = ',  iimax, ', jpj = ',  ijmax,   &
239               &                                                             ', jpi*jpj = ',  iimax* ijmax, ')'
240         ENDIF
241         IF( iimax*ijmax < jpimax*jpjmax ) THEN   ! chosen subdomain size is larger that the best subdomain size
242            llbest = .FALSE.
243            IF ( inbi*inbj-icnt2 < mppsize ) THEN
244               WRITE(ctmp1,*) '   ==> You could therefore have smaller mpi subdomains with less mpi processes'
245            ELSE
246               WRITE(ctmp1,*) '   ==> You could therefore have smaller mpi subdomains with the same number of mpi processes'
247            ENDIF
248            CALL ctl_warn( ' ', ctmp1, ' ', '    ---   YOU ARE WASTING CPU...   ---', ' ' )
249         ELSE IF ( iimax*ijmax == jpimax*jpjmax .AND. (inbi*inbj-icnt2) <  mppsize) THEN
250            llbest = .FALSE.
251            WRITE(ctmp1,*) '   ==> You could therefore have the same mpi subdomains size with less mpi processes'
252            CALL ctl_warn( ' ', ctmp1, ' ', '    ---   YOU ARE WASTING CPU...   ---', ' ' )
253         ELSE
254            llbest = .TRUE.
255         ENDIF
256      ENDIF
257     
258      ! look for land mpi subdomains...
259      ALLOCATE( llisoce(jpni,jpnj) )
260      CALL is_ocean( jpni, jpnj, llisoce )
261      inijmin = COUNT( llisoce )   ! number of oce subdomains
262
263      IF( mppsize < inijmin ) THEN   ! too many oce subdomains: can happen only if jpni and jpnj are prescribed...
264         WRITE(ctmp1,9001) '   With this specified domain decomposition: jpni = ', jpni, ' jpnj = ', jpnj
265         WRITE(ctmp2,9002) '   we can eliminate only ', jpni*jpnj - inijmin, ' land mpi subdomains therefore '
266         WRITE(ctmp3,9001) '   the number of ocean mpi subdomains (', inijmin,') exceed the number of MPI processes:', mppsize
267         WRITE(ctmp4,*) '   ==>>> There is the list of best domain decompositions you should use: '
268         CALL ctl_stop( ctmp1, ctmp2, ctmp3, ' ', ctmp4, ' ' )
269         CALL bestpartition( mppsize, ldlist = .TRUE. )   ! must be done by all core
270      ENDIF
271
272      IF( mppsize > jpni*jpnj ) THEN   ! not enough mpi subdomains for the total number of mpi processes
273         IF(lwp) THEN
274            WRITE(numout,9003) '   The number of mpi processes: ', mppsize
275            WRITE(numout,9003) '   exceeds the maximum number of subdomains (ocean+land) = ', jpni*jpnj
276            WRITE(numout,9001) '   defined by the following domain decomposition: jpni = ', jpni, ' jpnj = ', jpnj
277            WRITE(numout,   *) '   You should: '
278           IF( llauto ) THEN
279               WRITE(numout,*) '     - either prescribe your domain decomposition with the namelist variables'
280               WRITE(numout,*) '       jpni and jpnj to match the number of mpi process you want to use, '
281               WRITE(numout,*) '       even IF it not the best choice...'
282               WRITE(numout,*) '     - or keep the automatic and optimal domain decomposition by picking up one'
283               WRITE(numout,*) '       of the number of mpi process proposed in the list bellow'
284            ELSE
285               WRITE(numout,*) '     - either properly prescribe your domain decomposition with jpni and jpnj'
286               WRITE(numout,*) '       in order to be consistent with the number of mpi process you want to use'
287               WRITE(numout,*) '       even IF it not the best choice...'
288               WRITE(numout,*) '     - or use the automatic and optimal domain decomposition and pick up one of'
289               WRITE(numout,*) '       the domain decomposition proposed in the list bellow'
290            ENDIF
291            WRITE(numout,*)
292         ENDIF
293         CALL bestpartition( mppsize, ldlist = .TRUE. )   ! must be done by all core
294      ENDIF
295
296      jpnij = mppsize   ! force jpnij definition <-- remove as much land subdomains as needed to reach this condition
297      IF( mppsize > inijmin ) THEN
298         WRITE(ctmp1,9003) '   The number of mpi processes: ', mppsize
299         WRITE(ctmp2,9003) '   exceeds the maximum number of ocean subdomains = ', inijmin
300         WRITE(ctmp3,9002) '   we suppressed ', jpni*jpnj - mppsize, ' land subdomains '
301         WRITE(ctmp4,9002) '   BUT we had to keep ', mppsize - inijmin, ' land subdomains that are useless...'
302         CALL ctl_warn( ctmp1, ctmp2, ctmp3, ctmp4, ' ', '    --- YOU ARE WASTING CPU... ---', ' ' )
303      ELSE   ! mppsize = inijmin
304         IF(lwp) THEN
305            IF(llbest) WRITE(numout,*) '   ==> you use the best mpi decomposition'
306            WRITE(numout,*)
307            WRITE(numout,9003) '   Number of mpi processes: ', mppsize
308            WRITE(numout,9003) '   Number of ocean subdomains = ', inijmin
309            WRITE(numout,9003) '   Number of suppressed land subdomains = ', jpni*jpnj - inijmin
310            WRITE(numout,*)
311         ENDIF
312      ENDIF
3139000  FORMAT (a, i4, a, i4, a, i7, a)
3149001  FORMAT (a, i4, a, i4)
3159002  FORMAT (a, i4, a)
3169003  FORMAT (a, i5)
317
318      IF( numbot /= -1 )   CALL iom_close( numbot )
319      IF( numbdy /= -1 )   CALL iom_close( numbdy )
320     
321      ALLOCATE(  nfimpp(jpni ) , nfproc(jpni ) ,   nfjpi(jpni ) ,                     &
322         &       nimppt(jpnij) , ibonit(jpnij) ,  jpiall(jpnij) ,  jpjall(jpnij) ,    &
323         &       njmppt(jpnij) , ibonjt(jpnij) , nis0all(jpnij) , njs0all(jpnij) ,    &
324         &                                       nie0all(jpnij) , nje0all(jpnij) ,    &
325         &       iin(jpnij), ii_nono(jpnij), ii_noea(jpnij),   &
326         &       ijn(jpnij), ii_noso(jpnij), ii_nowe(jpnij),   &
327         &       iimppt(jpni,jpnj), ijpi(jpni,jpnj), ibondi(jpni,jpnj), ipproc(jpni,jpnj),   &
328         &       ijmppt(jpni,jpnj), ijpj(jpni,jpnj), ibondj(jpni,jpnj),  ipolj(jpni,jpnj),   &
329         &         iie0(jpni,jpnj), iis0(jpni,jpnj),   iono(jpni,jpnj),   ioea(jpni,jpnj),   &
330         &         ije0(jpni,jpnj), ijs0(jpni,jpnj),   ioso(jpni,jpnj),   iowe(jpni,jpnj),   &
331         &       STAT=ierr )
332      CALL mpp_sum( 'mppini', ierr )
333      IF( ierr /= 0 )   CALL ctl_stop( 'STOP', 'mpp_init: unable to allocate standard ocean arrays' )
334     
335#if defined key_agrif
336      IF( .NOT. Agrif_Root() ) THEN       ! AGRIF children: specific setting (cf. agrif_user.F90)
337         IF( jpiglo /= nbcellsx + 2 * ( nbghostcells + 1 + nn_hls ) )  &
338            CALL ctl_stop( 'STOP', 'mpp_init: Agrif children requires jpiglo == nbcellsx + 2*( nbghostcells + 1 + nn_hls )' )
339         IF( jpjglo /= nbcellsy + 2 * ( nbghostcells + 1 + nn_hls ) )  &
340            CALL ctl_stop( 'STOP', 'mpp_init: Agrif children requires jpjglo == nbcellsy + 2*( nbghostcells + 1 + nn_hls )' )
341         IF( ln_use_jattr )   CALL ctl_stop( 'STOP', 'mpp_init:Agrif children requires ln_use_jattr = .false. ' )
342      ENDIF
343#endif
344      !
345      !  2. Index arrays for subdomains
346      ! -----------------------------------
347      !
348      CALL basic_decomposition( jpiglo, jpjglo, nn_hls, jpni, jpnj, jpimax, jpjmax, iimppt, ijmppt, ijpi, ijpj )
349      nfimpp(:) = iimppt(:,jpnj)
350      nfjpi (:) =   ijpi(:,jpnj)
351      !
352      IF(lwp) THEN
353         WRITE(numout,*)
354         WRITE(numout,*) 'MPI Message Passing MPI - domain lay out over processors'
355         WRITE(numout,*)
356         WRITE(numout,*) '   defines mpp subdomains'
357         WRITE(numout,*) '      jpni = ', jpni 
358         WRITE(numout,*) '      jpnj = ', jpnj
359         WRITE(numout,*)
360         WRITE(numout,*) '      sum ijpi(i,1) = ', sum(ijpi(:,1)), ' jpiglo = ', jpiglo
361         WRITE(numout,*) '      sum ijpj(1,j) = ', sum(ijpj(1,:)), ' jpjglo = ', jpjglo
362      ENDIF
363     
364      ! 3. Subdomain description in the Regular Case
365      ! --------------------------------------------
366      ! specific cases where there is no communication -> must do the periodicity by itself
367      ! Warning: because of potential land-area suppression, do not use nbond[ij] == 2 
368      l_Iperio = jpni == 1 .AND. (jperio == 1 .OR. jperio == 4 .OR. jperio == 6 .OR. jperio == 7)
369      l_Jperio = jpnj == 1 .AND. (jperio == 2 .OR. jperio == 7)
370     
371      DO jarea = 1, jpni*jpnj
372         !
373         iarea0 = jarea - 1
374         ii = 1 + MOD(iarea0,jpni)
375         ij = 1 +     iarea0/jpni
376         ili = ijpi(ii,ij)
377         ilj = ijpj(ii,ij)
378         ibondi(ii,ij) = 0                         ! default: has e-w neighbours
379         IF( ii   ==    1 )   ibondi(ii,ij) = -1   ! first column, has only e neighbour
380         IF( ii   == jpni )   ibondi(ii,ij) =  1   ! last column,  has only w neighbour
381         IF( jpni ==    1 )   ibondi(ii,ij) =  2   ! has no e-w neighbour
382         ibondj(ii,ij) = 0                         ! default: has n-s neighbours
383         IF( ij   ==    1 )   ibondj(ii,ij) = -1   ! first row, has only n neighbour
384         IF( ij   == jpnj )   ibondj(ii,ij) =  1   ! last row,  has only s neighbour
385         IF( jpnj ==    1 )   ibondj(ii,ij) =  2   ! has no n-s neighbour
386
387         ! Subdomain neighbors (get their zone number): default definition
388         ioso(ii,ij) = iarea0 - jpni
389         iowe(ii,ij) = iarea0 - 1
390         ioea(ii,ij) = iarea0 + 1
391         iono(ii,ij) = iarea0 + jpni
392         iis0(ii,ij) =  1  + nn_hls
393         iie0(ii,ij) = ili - nn_hls
394         ijs0(ii,ij) =  1  + nn_hls
395         ije0(ii,ij) = ilj - nn_hls
396
397         ! East-West periodicity: change ibondi, ioea, iowe
398         IF( jperio == 1 .OR. jperio == 4 .OR. jperio == 6 .OR. jperio == 7 ) THEN
399            IF( jpni  /= 1 )   ibondi(ii,ij) = 0                        ! redefine: all have e-w neighbours
400            IF( ii ==    1 )   iowe(ii,ij) = iarea0 +        (jpni-1)   ! redefine: first column, address of w neighbour
401            IF( ii == jpni )   ioea(ii,ij) = iarea0 -        (jpni-1)   ! redefine: last column,  address of e neighbour
402         ENDIF
403
404         ! Simple North-South periodicity: change ibondj, ioso, iono
405         IF( jperio == 2 .OR. jperio == 7 ) THEN
406            IF( jpnj  /= 1 )   ibondj(ii,ij) = 0                        ! redefine: all have n-s neighbours
407            IF( ij ==    1 )   ioso(ii,ij) = iarea0 + jpni * (jpnj-1)   ! redefine: first row, address of s neighbour
408            IF( ij == jpnj )   iono(ii,ij) = iarea0 - jpni * (jpnj-1)   ! redefine: last row,  address of n neighbour
409         ENDIF
410
411         ! North fold: define ipolj, change iono. Warning: we do not change ibondj...
412         ipolj(ii,ij) = 0
413         IF( jperio == 3 .OR. jperio == 4 ) THEN
414            ijm1 = jpni*(jpnj-1)
415            imil = ijm1+(jpni+1)/2
416            IF( jarea > ijm1 ) ipolj(ii,ij) = 3
417            IF( MOD(jpni,2) == 1 .AND. jarea == imil ) ipolj(ii,ij) = 4
418            IF( ipolj(ii,ij) == 3 ) iono(ii,ij) = jpni*jpnj-jarea+ijm1   ! MPI rank of northern neighbour
419         ENDIF
420         IF( jperio == 5 .OR. jperio == 6 ) THEN
421            ijm1 = jpni*(jpnj-1)
422            imil = ijm1+(jpni+1)/2
423            IF( jarea > ijm1) ipolj(ii,ij) = 5
424            IF( MOD(jpni,2) == 1 .AND. jarea == imil ) ipolj(ii,ij) = 6
425            IF( ipolj(ii,ij) == 5) iono(ii,ij) = jpni*jpnj-jarea+ijm1    ! MPI rank of northern neighbour
426         ENDIF
427         !
428      END DO
429
430      ! 4. deal with land subdomains
431      ! ----------------------------
432      !
433      ! specify which subdomains are oce subdomains; other are land subdomains
434      ipproc(:,:) = -1
435      icont = -1
436      DO jarea = 1, jpni*jpnj
437         iarea0 = jarea - 1
438         ii = 1 + MOD(iarea0,jpni)
439         ij = 1 +     iarea0/jpni
440         IF( llisoce(ii,ij) ) THEN
441            icont = icont + 1
442            ipproc(ii,ij) = icont
443            iin(icont+1) = ii
444            ijn(icont+1) = ij
445         ENDIF
446      END DO
447      ! if needed add some land subdomains to reach jpnij active subdomains
448      i2add = jpnij - inijmin
449      DO jarea = 1, jpni*jpnj
450         iarea0 = jarea - 1
451         ii = 1 + MOD(iarea0,jpni)
452         ij = 1 +     iarea0/jpni
453         IF( .NOT. llisoce(ii,ij) .AND. i2add > 0 ) THEN
454            icont = icont + 1
455            ipproc(ii,ij) = icont
456            iin(icont+1) = ii
457            ijn(icont+1) = ij
458            i2add = i2add - 1
459         ENDIF
460      END DO
461      nfproc(:) = ipproc(:,jpnj)
462
463      ! neighbour treatment: change ibondi, ibondj if next to a land zone
464      DO jarea = 1, jpni*jpnj
465         ii = 1 + MOD( jarea-1  , jpni )
466         ij = 1 +     (jarea-1) / jpni
467         ! land-only area with an active n neigbour
468         IF ( ipproc(ii,ij) == -1 .AND. 0 <= iono(ii,ij) .AND. iono(ii,ij) <= jpni*jpnj-1 ) THEN
469            iino = 1 + MOD( iono(ii,ij) , jpni )                    ! ii index of this n neigbour
470            ijno = 1 +      iono(ii,ij) / jpni                      ! ij index of this n neigbour
471            ! In case of north fold exchange: I am the n neigbour of my n neigbour!! (#1057)
472            ! --> for northern neighbours of northern row processors (in case of north-fold)
473            !     need to reverse the LOGICAL direction of communication
474            idir = 1                                           ! we are indeed the s neigbour of this n neigbour
475            IF( ij == jpnj .AND. ijno == jpnj )   idir = -1    ! both are on the last row, we are in fact the n neigbour
476            IF( ibondj(iino,ijno) == idir     )   ibondj(iino,ijno) =   2     ! this n neigbour had only a s/n neigbour -> no more
477            IF( ibondj(iino,ijno) == 0        )   ibondj(iino,ijno) = -idir   ! this n neigbour had both, n-s neighbours -> keep 1
478         ENDIF
479         ! land-only area with an active s neigbour
480         IF( ipproc(ii,ij) == -1 .AND. 0 <= ioso(ii,ij) .AND. ioso(ii,ij) <= jpni*jpnj-1 ) THEN
481            iiso = 1 + MOD( ioso(ii,ij) , jpni )                    ! ii index of this s neigbour
482            ijso = 1 +      ioso(ii,ij) / jpni                      ! ij index of this s neigbour
483            IF( ibondj(iiso,ijso) == -1 )   ibondj(iiso,ijso) = 2   ! this s neigbour had only a n neigbour    -> no more neigbour
484            IF( ibondj(iiso,ijso) ==  0 )   ibondj(iiso,ijso) = 1   ! this s neigbour had both, n-s neighbours -> keep s neigbour
485         ENDIF
486         ! land-only area with an active e neigbour
487         IF( ipproc(ii,ij) == -1 .AND. 0 <= ioea(ii,ij) .AND. ioea(ii,ij) <= jpni*jpnj-1 ) THEN
488            iiea = 1 + MOD( ioea(ii,ij) , jpni )                    ! ii index of this e neigbour
489            ijea = 1 +      ioea(ii,ij) / jpni                      ! ij index of this e neigbour
490            IF( ibondi(iiea,ijea) == 1 )   ibondi(iiea,ijea) =  2   ! this e neigbour had only a w neigbour    -> no more neigbour
491            IF( ibondi(iiea,ijea) == 0 )   ibondi(iiea,ijea) = -1   ! this e neigbour had both, e-w neighbours -> keep e neigbour
492         ENDIF
493         ! land-only area with an active w neigbour
494         IF( ipproc(ii,ij) == -1 .AND. 0 <= iowe(ii,ij) .AND. iowe(ii,ij) <= jpni*jpnj-1) THEN
495            iiwe = 1 + MOD( iowe(ii,ij) , jpni )                    ! ii index of this w neigbour
496            ijwe = 1 +      iowe(ii,ij) / jpni                      ! ij index of this w neigbour
497            IF( ibondi(iiwe,ijwe) == -1 )   ibondi(iiwe,ijwe) = 2   ! this w neigbour had only a e neigbour    -> no more neigbour
498            IF( ibondi(iiwe,ijwe) ==  0 )   ibondi(iiwe,ijwe) = 1   ! this w neigbour had both, e-w neighbours -> keep w neigbour
499         ENDIF
500      END DO
501     
502      ! 5. Subdomain print
503      ! ------------------
504      IF(lwp) THEN
505         ifreq = 4
506         il1 = 1
507         DO jn = 1, (jpni-1)/ifreq+1
508            il2 = MIN(jpni,il1+ifreq-1)
509            WRITE(numout,*)
510            WRITE(numout,9400) ('***',ji=il1,il2-1)
511            DO jj = jpnj, 1, -1
512               WRITE(numout,9403) ('   ',ji=il1,il2-1)
513               WRITE(numout,9402) jj, (ijpi(ji,jj),ijpj(ji,jj),ji=il1,il2)
514               WRITE(numout,9404) (ipproc(ji,jj),ji=il1,il2)
515               WRITE(numout,9403) ('   ',ji=il1,il2-1)
516               WRITE(numout,9400) ('***',ji=il1,il2-1)
517            END DO
518            WRITE(numout,9401) (ji,ji=il1,il2)
519            il1 = il1+ifreq
520         END DO
521 9400    FORMAT('           ***'   ,20('*************',a3)    )
522 9403    FORMAT('           *     ',20('         *   ',a3)    )
523 9401    FORMAT('              '   ,20('   ',i3,'          ') )
524 9402    FORMAT('       ',i3,' *  ',20(i3,'  x',i3,'   *   ') )
525 9404    FORMAT('           *  '   ,20('     ' ,i4,'   *   ') )
526      ENDIF
527         
528      ! just to save nono etc for all proc
529      ! warning ii*ij (zone) /= nproc (processors)!
530      ! ioso = zone number, ii_noso = proc number
531      ii_noso(:) = -1
532      ii_nono(:) = -1
533      ii_noea(:) = -1
534      ii_nowe(:) = -1 
535      DO jproc = 1, jpnij
536         ii = iin(jproc)
537         ij = ijn(jproc)
538         IF( 0 <= ioso(ii,ij) .AND. ioso(ii,ij) <= (jpni*jpnj-1) ) THEN
539            iiso = 1 + MOD( ioso(ii,ij) , jpni )
540            ijso = 1 +      ioso(ii,ij) / jpni
541            ii_noso(jproc) = ipproc(iiso,ijso)
542         ENDIF
543         IF( 0 <= iowe(ii,ij) .AND. iowe(ii,ij) <= (jpni*jpnj-1) ) THEN
544          iiwe = 1 + MOD( iowe(ii,ij) , jpni )
545          ijwe = 1 +      iowe(ii,ij) / jpni
546          ii_nowe(jproc) = ipproc(iiwe,ijwe)
547         ENDIF
548         IF( 0 <= ioea(ii,ij) .AND. ioea(ii,ij) <= (jpni*jpnj-1) ) THEN
549            iiea = 1 + MOD( ioea(ii,ij) , jpni )
550            ijea = 1 +      ioea(ii,ij) / jpni
551            ii_noea(jproc)= ipproc(iiea,ijea)
552         ENDIF
553         IF( 0 <= iono(ii,ij) .AND. iono(ii,ij) <= (jpni*jpnj-1) ) THEN
554            iino = 1 + MOD( iono(ii,ij) , jpni )
555            ijno = 1 +      iono(ii,ij) / jpni
556            ii_nono(jproc)= ipproc(iino,ijno)
557         ENDIF
558      END DO
559   
560      ! 6. Change processor name
561      ! ------------------------
562      ii = iin(narea)
563      ij = ijn(narea)
564      !
565      ! set default neighbours
566      noso = ii_noso(narea)
567      nowe = ii_nowe(narea)
568      noea = ii_noea(narea)
569      nono = ii_nono(narea)
570      jpi    = ijpi(ii,ij) 
571!!$      Nis0  = iis0(ii,ij)
572!!$      Nie0  = iie0(ii,ij)
573      jpj    = ijpj(ii,ij) 
574!!$      Njs0  = ijs0(ii,ij)
575!!$      Nje0  = ije0(ii,ij)
576      nbondi = ibondi(ii,ij)
577      nbondj = ibondj(ii,ij)
578      nimpp = iimppt(ii,ij) 
579      njmpp = ijmppt(ii,ij)
580      jpk = jpkglo                              ! third dim
581#if defined key_agrif
582      ! simple trick to use same vertical grid as parent but different number of levels:
583      ! Save maximum number of levels in jpkglo, then define all vertical grids with this number.
584      ! Suppress once vertical online interpolation is ok
585!!$      IF(.NOT.Agrif_Root())   jpkglo = Agrif_Parent( jpkglo )
586#endif
587      !
588      CALL init_doloop                          ! set start/end indices of do-loop, depending on the halo width value (nn_hls)
589      !
590      jpim1 = jpi-1                             ! inner domain indices
591      jpjm1 = jpj-1                             !   "           "
592      jpkm1 = MAX( 1, jpk-1 )                   !   "           "
593      jpij  = jpi*jpj                           !  jpi x j
594      DO jproc = 1, jpnij
595         ii = iin(jproc)
596         ij = ijn(jproc)
597         jpiall (jproc) = ijpi(ii,ij)
598         nis0all(jproc) = iis0(ii,ij)
599         nie0all(jproc) = iie0(ii,ij)
600         jpjall (jproc) = ijpj(ii,ij)
601         njs0all(jproc) = ijs0(ii,ij)
602         nje0all(jproc) = ije0(ii,ij)
603         ibonit(jproc) = ibondi(ii,ij)
604         ibonjt(jproc) = ibondj(ii,ij)
605         nimppt(jproc) = iimppt(ii,ij) 
606         njmppt(jproc) = ijmppt(ii,ij) 
607      END DO
608
609      ! Save processor layout in ascii file
610      IF (llwrtlay) THEN
611         CALL ctl_opn( inum, 'layout.dat', 'REPLACE', 'FORMATTED', 'SEQUENTIAL', -1, numout, .FALSE., narea )
612         WRITE(inum,'(a)') '   jpnij   jpimax  jpjmax    jpk  jpiglo  jpjglo'//&
613   &           ' ( local:    narea     jpi     jpj )'
614         WRITE(inum,'(6i8,a,3i8,a)') jpnij,jpimax,jpjmax,jpk,jpiglo,jpjglo,&
615   &           ' ( local: ',narea,jpi,jpj,' )'
616         WRITE(inum,'(a)') 'nproc jpi jpj Nis0 Njs0 Nie0 Nje0 nimp njmp nono noso nowe noea nbondi nbondj '
617
618         DO jproc = 1, jpnij
619            WRITE(inum,'(13i5,2i7)')   jproc-1,  jpiall(jproc),  jpjall(jproc),   &
620               &                                nis0all(jproc), njs0all(jproc),   &
621               &                                nie0all(jproc), nje0all(jproc),   &
622               &                                nimppt (jproc), njmppt (jproc),   & 
623               &                                ii_nono(jproc), ii_noso(jproc),   &
624               &                                ii_nowe(jproc), ii_noea(jproc),   &
625               &                                ibonit (jproc), ibonjt (jproc) 
626         END DO
627      END IF
628
629      !                          ! north fold parameter
630      ! Defined npolj, either 0, 3 , 4 , 5 , 6
631      ! In this case the important thing is that npolj /= 0
632      ! Because if we go through these line it is because jpni >1 and thus
633      ! we must use lbcnorthmpp, which tests only npolj =0 or npolj /= 0
634      npolj = 0
635      ij = ijn(narea)
636      IF( jperio == 3 .OR. jperio == 4 ) THEN
637         IF( ij == jpnj )   npolj = 3
638      ENDIF
639      IF( jperio == 5 .OR. jperio == 6 ) THEN
640         IF( ij == jpnj )   npolj = 5
641      ENDIF
642      !
643      nproc = narea-1
644      IF(lwp) THEN
645         WRITE(numout,*)
646         WRITE(numout,*) '   resulting internal parameters : '
647         WRITE(numout,*) '      nproc  = ', nproc
648         WRITE(numout,*) '      nowe   = ', nowe  , '   noea  =  ', noea
649         WRITE(numout,*) '      nono   = ', nono  , '   noso  =  ', noso
650         WRITE(numout,*) '      nbondi = ', nbondi
651         WRITE(numout,*) '      nbondj = ', nbondj
652         WRITE(numout,*) '      npolj  = ', npolj
653         WRITE(numout,*) '    l_Iperio = ', l_Iperio
654         WRITE(numout,*) '    l_Jperio = ', l_Jperio
655         WRITE(numout,*) '      nimpp  = ', nimpp
656         WRITE(numout,*) '      njmpp  = ', njmpp
657         WRITE(numout,*) '      nn_hls = ', nn_hls 
658      ENDIF
659
660      !                          ! Prepare mpp north fold
661      IF( jperio >= 3 .AND. jperio <= 6 .AND. jpni > 1 ) THEN
662         CALL mpp_ini_north
663         IF (lwp) THEN
664            WRITE(numout,*)
665            WRITE(numout,*) '   ==>>>   North fold boundary prepared for jpni >1'
666            ! additional prints in layout.dat
667         ENDIF
668         IF (llwrtlay) THEN
669            WRITE(inum,*)
670            WRITE(inum,*)
671            WRITE(inum,*) 'number of subdomains located along the north fold : ', ndim_rank_north
672            WRITE(inum,*) 'Rank of the subdomains located along the north fold : ', ndim_rank_north
673            DO jproc = 1, ndim_rank_north, 5
674               WRITE(inum,*) nrank_north( jproc:MINVAL( (/jproc+4,ndim_rank_north/) ) )
675            END DO
676         ENDIF
677      ENDIF
678      !
679      CALL init_ioipsl       ! Prepare NetCDF output file (if necessary)
680      !     
681      IF (( jperio >= 3 .AND. jperio <= 6 .AND. jpni > 1 ).AND.( ln_nnogather )) THEN
682         CALL init_nfdcom     ! northfold neighbour lists
683         IF (llwrtlay) THEN
684            WRITE(inum,*)
685            WRITE(inum,*)
686            WRITE(inum,*) 'north fold exchanges with explicit point-to-point messaging :'
687            WRITE(inum,*) 'nsndto : ', nsndto
688            WRITE(inum,*) 'isendto : ', isendto
689         ENDIF
690      ENDIF
691      !
692      IF (llwrtlay) CLOSE(inum)   
693      !
694      DEALLOCATE(iin, ijn, ii_nono, ii_noea, ii_noso, ii_nowe,    &
695         &       iimppt, ijmppt, ibondi, ibondj, ipproc, ipolj,   &
696         &       ijpi, ijpj, iie0, ije0, iis0, ijs0,              &
697         &       iono, ioea, ioso, iowe, llisoce)
698      !
699    END SUBROUTINE mpp_init
700
701
702    SUBROUTINE basic_decomposition( kiglo, kjglo, khls, knbi, knbj, kimax, kjmax, kimppt, kjmppt, klci, klcj)
703      !!----------------------------------------------------------------------
704      !!                  ***  ROUTINE basic_decomposition  ***
705      !!                   
706      !! ** Purpose :   Lay out the global domain over processors.
707      !!
708      !! ** Method  :   Global domain is distributed in smaller local domains.
709      !!
710      !! ** Action : - set for all knbi*knbj domains:
711      !!                    kimppt     : longitudinal index
712      !!                    kjmppt     : latitudinal  index
713      !!                    klci       : first dimension
714      !!                    klcj       : second dimension
715      !!----------------------------------------------------------------------
716      INTEGER,                                 INTENT(in   ) ::   kiglo, kjglo
717      INTEGER,                                 INTENT(in   ) ::   khls
718      INTEGER,                                 INTENT(in   ) ::   knbi, knbj
719      INTEGER,                                 INTENT(  out) ::   kimax, kjmax
720      INTEGER, DIMENSION(knbi,knbj), OPTIONAL, INTENT(  out) ::   kimppt, kjmppt
721      INTEGER, DIMENSION(knbi,knbj), OPTIONAL, INTENT(  out) ::   klci, klcj
722      !
723      INTEGER ::   ji, jj
724      INTEGER ::   i2hls 
725      INTEGER ::   iresti, irestj, irm, ijpjmin
726      !!----------------------------------------------------------------------
727      i2hls = 2*khls
728      !
729#if defined key_nemocice_decomp
730      kimax = ( nx_global+2-i2hls + (knbi-1) ) / knbi + i2hls    ! first  dim.
731      kjmax = ( ny_global+2-i2hls + (knbj-1) ) / knbj + i2hls    ! second dim.
732#else
733      kimax = ( kiglo - i2hls + (knbi-1) ) / knbi + i2hls    ! first  dim.
734      kjmax = ( kjglo - i2hls + (knbj-1) ) / knbj + i2hls    ! second dim.
735#endif
736      IF( .NOT. PRESENT(kimppt) ) RETURN
737      !
738      !  1. Dimension arrays for subdomains
739      ! -----------------------------------
740      !  Computation of local domain sizes klci() klcj()
741      !  These dimensions depend on global sizes knbi,knbj and kiglo,kjglo
742      !  The subdomains are squares lesser than or equal to the global
743      !  dimensions divided by the number of processors minus the overlap array.
744      !
745      iresti = 1 + MOD( kiglo - i2hls - 1 , knbi )
746      irestj = 1 + MOD( kjglo - i2hls - 1 , knbj )
747      !
748      !  Need to use kimax and kjmax here since jpi and jpj not yet defined
749#if defined key_nemocice_decomp
750      ! Change padding to be consistent with CICE
751      klci(1:knbi-1,:       ) = kimax
752      klci(  knbi  ,:       ) = kiglo - (knbi - 1) * (kimax - i2hls)
753      klcj(:       ,1:knbj-1) = kjmax
754      klcj(:       ,  knbj  ) = kjglo - (knbj - 1) * (kjmax - i2hls)
755#else
756      klci(1:iresti      ,:) = kimax
757      klci(iresti+1:knbi ,:) = kimax-1
758      IF( MINVAL(klci) < 2*i2hls ) THEN
759         WRITE(ctmp1,*) '   basic_decomposition: minimum value of jpi must be >= ', 2*i2hls
760         WRITE(ctmp2,*) '   We have ', MINVAL(klci)
761        CALL ctl_stop( 'STOP', ctmp1, ctmp2 )
762      ENDIF
763      IF( jperio == 3 .OR. jperio == 4 .OR. jperio == 5 .OR. jperio == 6 ) THEN
764         ! minimize the size of the last row to compensate for the north pole folding coast
765         IF( jperio == 3 .OR. jperio == 4 )   ijpjmin = 2+3*khls   ! V and F folding must be outside of southern halos
766         IF( jperio == 5 .OR. jperio == 6 )   ijpjmin = 1+3*khls   ! V and F folding must be outside of southern halos
767         irm = knbj - irestj                                       ! total number of lines to be removed
768         klcj(:,knbj) = MAX( ijpjmin, kjmax-irm )                  ! we must have jpj >= ijpjmin in the last row
769         irm = irm - ( kjmax - klcj(1,knbj) )                      ! remaining number of lines to remove
770         irestj = knbj - 1 - irm
771         klcj(:, irestj+1:knbj-1) = kjmax-1
772      ELSE
773         klcj(:, irestj+1:knbj  ) = kjmax-1
774      ENDIF
775      klcj(:,1:irestj) = kjmax
776      IF( MINVAL(klcj) < 2*i2hls ) THEN
777         WRITE(ctmp1,*) '   basic_decomposition: minimum value of jpj must be >= ', 2*i2hls
778         WRITE(ctmp2,*) '   We have ', MINVAL(klcj)
779         CALL ctl_stop( 'STOP', ctmp1, ctmp2 )
780      ENDIF
781#endif
782
783      !  2. Index arrays for subdomains
784      ! -------------------------------
785      kimppt(:,:) = 1
786      kjmppt(:,:) = 1
787      !
788      IF( knbi > 1 ) THEN
789         DO jj = 1, knbj
790            DO ji = 2, knbi
791               kimppt(ji,jj) = kimppt(ji-1,jj) + klci(ji-1,jj) - 2 * nn_hls
792            END DO
793         END DO
794      ENDIF
795      !
796      IF( knbj > 1 )THEN
797         DO jj = 2, knbj
798            DO ji = 1, knbi
799               kjmppt(ji,jj) = kjmppt(ji,jj-1) + klcj(ji,jj-1) - 2 * nn_hls
800            END DO
801         END DO
802      ENDIF
803     
804   END SUBROUTINE basic_decomposition
805
806
807   SUBROUTINE bestpartition( knbij, knbi, knbj, knbcnt, ldlist )
808      !!----------------------------------------------------------------------
809      !!                 ***  ROUTINE bestpartition  ***
810      !!
811      !! ** Purpose :
812      !!
813      !! ** Method  :
814      !!----------------------------------------------------------------------
815      INTEGER,           INTENT(in   ) ::   knbij         ! total number if subdomains               (knbi*knbj)
816      INTEGER, OPTIONAL, INTENT(  out) ::   knbi, knbj    ! number if subdomains along i and j (knbi and knbj)
817      INTEGER, OPTIONAL, INTENT(  out) ::   knbcnt        ! number of land subdomains
818      LOGICAL, OPTIONAL, INTENT(in   ) ::   ldlist        ! .true.: print the list the best domain decompositions (with land)
819      !
820      INTEGER :: ji, jj, ii, iitarget
821      INTEGER :: iszitst, iszjtst
822      INTEGER :: isziref, iszjref
823      INTEGER :: inbij, iszij
824      INTEGER :: inbimax, inbjmax, inbijmax, inbijold
825      INTEGER :: isz0, isz1
826      INTEGER, DIMENSION(  :), ALLOCATABLE :: indexok
827      INTEGER, DIMENSION(  :), ALLOCATABLE :: inbi0, inbj0, inbij0   ! number of subdomains along i,j
828      INTEGER, DIMENSION(  :), ALLOCATABLE :: iszi0, iszj0, iszij0   ! max size of the subdomains along i,j
829      INTEGER, DIMENSION(  :), ALLOCATABLE :: inbi1, inbj1, inbij1   ! number of subdomains along i,j
830      INTEGER, DIMENSION(  :), ALLOCATABLE :: iszi1, iszj1, iszij1   ! max size of the subdomains along i,j
831      LOGICAL :: llist
832      LOGICAL, DIMENSION(:,:), ALLOCATABLE :: llmsk2d                 ! max size of the subdomains along i,j
833      LOGICAL, DIMENSION(:,:), ALLOCATABLE :: llisoce              !  -     -
834      REAL(wp)::   zpropland
835      !!----------------------------------------------------------------------
836      !
837      llist = .FALSE.
838      IF( PRESENT(ldlist) ) llist = ldlist
839
840      CALL mpp_init_landprop( zpropland )                      ! get the proportion of land point over the gloal domain
841      inbij = NINT( REAL(knbij, wp) / ( 1.0 - zpropland ) )    ! define the largest possible value for jpni*jpnj
842      !
843      IF( llist ) THEN   ;   inbijmax = inbij*2
844      ELSE               ;   inbijmax = inbij
845      ENDIF
846      !
847      ALLOCATE(inbi0(inbijmax),inbj0(inbijmax),iszi0(inbijmax),iszj0(inbijmax))
848      !
849      inbimax = 0
850      inbjmax = 0
851      isziref = jpiglo*jpjglo+1
852      iszjref = jpiglo*jpjglo+1
853      !
854      ! get the list of knbi that gives a smaller jpimax than knbi-1
855      ! get the list of knbj that gives a smaller jpjmax than knbj-1
856      DO ji = 1, inbijmax     
857#if defined key_nemocice_decomp
858         iszitst = ( nx_global+2-2*nn_hls + (ji-1) ) / ji + 2*nn_hls    ! first  dim.
859#else
860         iszitst = ( Ni0glo + (ji-1) ) / ji + 2*nn_hls
861#endif
862         IF( iszitst < isziref ) THEN
863            isziref = iszitst
864            inbimax = inbimax + 1
865            inbi0(inbimax) = ji
866            iszi0(inbimax) = isziref
867         ENDIF
868#if defined key_nemocice_decomp
869         iszjtst = ( ny_global+2-2*nn_hls + (ji-1) ) / ji + 2*nn_hls    ! first  dim.
870#else
871         iszjtst = ( Nj0glo + (ji-1) ) / ji + 2*nn_hls
872#endif
873         IF( iszjtst < iszjref ) THEN
874            iszjref = iszjtst
875            inbjmax = inbjmax + 1
876            inbj0(inbjmax) = ji
877            iszj0(inbjmax) = iszjref
878         ENDIF
879      END DO
880
881      ! combine these 2 lists to get all possible knbi*knbj <  inbijmax
882      ALLOCATE( llmsk2d(inbimax,inbjmax) )
883      DO jj = 1, inbjmax
884         DO ji = 1, inbimax
885            IF ( inbi0(ji) * inbj0(jj) <= inbijmax ) THEN   ;   llmsk2d(ji,jj) = .TRUE.
886            ELSE                                            ;   llmsk2d(ji,jj) = .FALSE.
887            ENDIF
888         END DO
889      END DO
890      isz1 = COUNT(llmsk2d)
891      ALLOCATE( inbi1(isz1), inbj1(isz1), iszi1(isz1), iszj1(isz1) )
892      ii = 0
893      DO jj = 1, inbjmax
894         DO ji = 1, inbimax
895            IF( llmsk2d(ji,jj) .EQV. .TRUE. ) THEN
896               ii = ii + 1
897               inbi1(ii) = inbi0(ji)
898               inbj1(ii) = inbj0(jj)
899               iszi1(ii) = iszi0(ji)
900               iszj1(ii) = iszj0(jj)
901            END IF
902         END DO
903      END DO
904      DEALLOCATE( inbi0, inbj0, iszi0, iszj0 )
905      DEALLOCATE( llmsk2d )
906
907      ALLOCATE( inbij1(isz1), iszij1(isz1) )
908      inbij1(:) = inbi1(:) * inbj1(:)
909      iszij1(:) = iszi1(:) * iszj1(:)
910
911      ! if therr is no land and no print
912      IF( .NOT. llist .AND. numbot == -1 .AND. numbdy == -1 ) THEN
913         ! get the smaller partition which gives the smallest subdomain size
914         ii = MINLOC(inbij1, mask = iszij1 == MINVAL(iszij1), dim = 1)
915         knbi = inbi1(ii)
916         knbj = inbj1(ii)
917         IF(PRESENT(knbcnt))   knbcnt = 0
918         DEALLOCATE( inbi1, inbj1, inbij1, iszi1, iszj1, iszij1 )
919         RETURN
920      ENDIF
921
922      ! extract only the partitions which reduce the subdomain size in comparison with smaller partitions
923      ALLOCATE( indexok(isz1) )                                 ! to store indices of the best partitions
924      isz0 = 0                                                  ! number of best partitions     
925      inbij = 1                                                 ! start with the min value of inbij1 => 1
926      iszij = jpiglo*jpjglo+1                                   ! default: larger than global domain
927      DO WHILE( inbij <= inbijmax )                             ! if we did not reach the max of inbij1
928         ii = MINLOC(iszij1, mask = inbij1 == inbij, dim = 1)   ! warning: send back the first occurence if multiple results
929         IF ( iszij1(ii) < iszij ) THEN
930            isz0 = isz0 + 1
931            indexok(isz0) = ii
932            iszij = iszij1(ii)
933         ENDIF
934         inbij = MINVAL(inbij1, mask = inbij1 > inbij)   ! warning: return largest integer value if mask = .false. everywhere
935      END DO
936      DEALLOCATE( inbij1, iszij1 )
937
938      ! keep only the best partitions (sorted by increasing order of subdomains number and decreassing subdomain size)
939      ALLOCATE( inbi0(isz0), inbj0(isz0), iszi0(isz0), iszj0(isz0) )
940      DO ji = 1, isz0
941         ii = indexok(ji)
942         inbi0(ji) = inbi1(ii)
943         inbj0(ji) = inbj1(ii)
944         iszi0(ji) = iszi1(ii)
945         iszj0(ji) = iszj1(ii)
946      END DO
947      DEALLOCATE( indexok, inbi1, inbj1, iszi1, iszj1 )
948
949      IF( llist ) THEN
950         IF(lwp) THEN
951            WRITE(numout,*)
952            WRITE(numout,*) '                  For your information:'
953            WRITE(numout,*) '  list of the best partitions including land supression'
954            WRITE(numout,*) '  -----------------------------------------------------'
955            WRITE(numout,*)
956         END IF
957         ji = isz0   ! initialization with the largest value
958         ALLOCATE( llisoce(inbi0(ji), inbj0(ji)) )
959         CALL is_ocean( inbi0(ji), inbj0(ji), llisoce ) ! Warning: must be call by all cores (call mpp_sum)
960         inbijold = COUNT(llisoce)
961         DEALLOCATE( llisoce )
962         DO ji =isz0-1,1,-1
963            ALLOCATE( llisoce(inbi0(ji), inbj0(ji)) )
964            CALL is_ocean( inbi0(ji), inbj0(ji), llisoce ) ! Warning: must be call by all cores (call mpp_sum)
965            inbij = COUNT(llisoce)
966            DEALLOCATE( llisoce )
967            IF(lwp .AND. inbij < inbijold) THEN
968               WRITE(numout,'(a, i6, a, i6, a, f4.1, a, i9, a, i6, a, i6, a)')                                 &
969                  &   'nb_cores oce: ', inbij, ', land domains excluded: ', inbi0(ji)*inbj0(ji) - inbij,       &
970                  &   ' (', REAL(inbi0(ji)*inbj0(ji) - inbij,wp) / REAL(inbi0(ji)*inbj0(ji),wp) *100.,         &
971                  &   '%), largest oce domain: ', iszi0(ji)*iszj0(ji), ' ( ', iszi0(ji),' x ', iszj0(ji), ' )'
972               inbijold = inbij
973            END IF
974         END DO
975         DEALLOCATE( inbi0, inbj0, iszi0, iszj0 )
976         IF(lwp) THEN
977            WRITE(numout,*)
978            WRITE(numout,*)  '  -----------------------------------------------------------'
979         ENDIF
980         CALL mppsync
981         CALL mppstop( ld_abort = .TRUE. )
982      ENDIF
983     
984      DEALLOCATE( iszi0, iszj0 )
985      inbij = inbijmax + 1        ! default: larger than possible
986      ii = isz0+1                 ! start from the end of the list (smaller subdomains)
987      DO WHILE( inbij > knbij )   ! while the number of ocean subdomains exceed the number of procs
988         ii = ii -1 
989         ALLOCATE( llisoce(inbi0(ii), inbj0(ii)) )
990         CALL is_ocean( inbi0(ii), inbj0(ii), llisoce )            ! must be done by all core
991         inbij = COUNT(llisoce)
992         DEALLOCATE( llisoce )
993      END DO
994      knbi = inbi0(ii)
995      knbj = inbj0(ii)
996      IF(PRESENT(knbcnt))   knbcnt = knbi * knbj - inbij
997      DEALLOCATE( inbi0, inbj0 )
998      !
999   END SUBROUTINE bestpartition
1000   
1001   
1002   SUBROUTINE mpp_init_landprop( propland )
1003      !!----------------------------------------------------------------------
1004      !!                  ***  ROUTINE mpp_init_landprop  ***
1005      !!
1006      !! ** Purpose : the the proportion of land points in the surface land-sea mask
1007      !!
1008      !! ** Method  : read iproc strips (of length Ni0glo) of the land-sea mask
1009      !!----------------------------------------------------------------------
1010      REAL(wp), INTENT(  out) :: propland    ! proportion of land points in the global domain (between 0 and 1)
1011      !
1012      INTEGER, DIMENSION(jpni*jpnj) ::   kusedom_1d
1013      INTEGER :: inboce, iarea
1014      INTEGER :: iproc, idiv, ijsz
1015      INTEGER :: ijstr
1016      LOGICAL, ALLOCATABLE, DIMENSION(:,:) ::   lloce
1017      !!----------------------------------------------------------------------
1018      ! do nothing if there is no land-sea mask
1019      IF( numbot == -1 .and. numbdy == -1 ) THEN
1020         propland = 0.
1021         RETURN
1022      ENDIF
1023
1024      ! number of processes reading the bathymetry file
1025      iproc = MINVAL( (/mppsize, Nj0glo/2, 100/) )  ! read a least 2 lines, no more that 100 processes reading at the same time
1026     
1027      ! we want to read iproc strips of the land-sea mask. -> pick up iproc processes every idiv processes starting at 1
1028      IF( iproc == 1 ) THEN   ;   idiv = mppsize
1029      ELSE                    ;   idiv = ( mppsize - 1 ) / ( iproc - 1 )
1030      ENDIF
1031
1032      iarea = (narea-1)/idiv   ! involed process number (starting counting at 0)
1033      IF( MOD( narea-1, idiv ) == 0 .AND. iarea < iproc ) THEN   ! beware idiv can be = to 1
1034         !
1035         ijsz = Nj0glo / iproc                                               ! width of the stripe to read
1036         IF( iarea < MOD(Nj0glo,iproc) ) ijsz = ijsz + 1
1037         ijstr = iarea*(Nj0glo/iproc) + MIN(iarea, MOD(Nj0glo,iproc)) + 1    ! starting j position of the reading
1038         !
1039         ALLOCATE( lloce(Ni0glo, ijsz) )                                     ! allocate the strip
1040         CALL readbot_strip( ijstr, ijsz, lloce )
1041         inboce = COUNT(lloce)                                               ! number of ocean point in the stripe
1042         DEALLOCATE(lloce)
1043         !
1044      ELSE
1045         inboce = 0
1046      ENDIF
1047      CALL mpp_sum( 'mppini', inboce )   ! total number of ocean points over the global domain
1048      !
1049      propland = REAL( Ni0glo*Nj0glo - inboce, wp ) / REAL( Ni0glo*Nj0glo, wp ) 
1050      !
1051   END SUBROUTINE mpp_init_landprop
1052   
1053   
1054   SUBROUTINE is_ocean( knbi, knbj, ldisoce )
1055      !!----------------------------------------------------------------------
1056      !!                  ***  ROUTINE mpp_init_nboce  ***
1057      !!
1058      !! ** Purpose : Check for a mpi domain decomposition knbi x knbj which
1059      !!              subdomains, including 1 halo (even if nn_hls>1), contain
1060      !!              at least 1 ocean point.
1061      !!              We must indeed ensure that each subdomain that is a neighbour
1062      !!              of a land subdomain as only land points on its boundary
1063      !!              (inside the inner subdomain) with the land subdomain.
1064      !!              This is needed to get the proper bondary conditions on
1065      !!              a subdomain with a closed boundary.
1066      !!
1067      !! ** Method  : read knbj strips (of length Ni0glo) of the land-sea mask
1068      !!----------------------------------------------------------------------
1069      INTEGER,                       INTENT(in   ) ::   knbi, knbj     ! domain decomposition
1070      LOGICAL, DIMENSION(knbi,knbj), INTENT(  out) ::   ldisoce        ! .true. if a sub domain constains 1 ocean point
1071      !
1072      INTEGER, DIMENSION(knbi,knbj) ::   inboce                        ! number oce oce pint in each mpi subdomain
1073      INTEGER, DIMENSION(knbi*knbj) ::   inboce_1d
1074      INTEGER :: idiv, iimax, ijmax, iarea
1075      INTEGER :: inx, iny, inry, isty
1076      INTEGER :: ji, jn
1077      LOGICAL, ALLOCATABLE, DIMENSION(:,:) ::   lloce                  ! lloce(i,j) = .true. if the point (i,j) is ocean
1078      INTEGER, ALLOCATABLE, DIMENSION(:,:) ::   iimppt, ijpi
1079      INTEGER, ALLOCATABLE, DIMENSION(:,:) ::   ijmppt, ijpj
1080      !!----------------------------------------------------------------------
1081      ! do nothing if there is no land-sea mask
1082      IF( numbot == -1 .AND. numbdy == -1 ) THEN
1083         ldisoce(:,:) = .TRUE.
1084         RETURN
1085      ENDIF
1086
1087      ! we want to read knbj strips of the land-sea mask. -> pick up knbj processes every idiv processes starting at 1
1088      IF           ( knbj == 1 ) THEN   ;   idiv = mppsize
1089      ELSE IF ( mppsize < knbj ) THEN   ;   idiv = 1
1090      ELSE                              ;   idiv = ( mppsize - 1 ) / ( knbj - 1 )
1091      ENDIF
1092      inboce(:,:) = 0          ! default no ocean point found
1093
1094      DO jn = 0, (knbj-1)/mppsize   ! if mppsize < knbj : more strips than mpi processes (because of potential land domains)
1095         !
1096         iarea = (narea-1)/idiv + jn * mppsize + 1                     ! involed process number (starting counting at 1)
1097         IF( MOD( narea-1, idiv ) == 0 .AND. iarea <= knbj ) THEN      ! beware idiv can be = to 1
1098            !
1099            ALLOCATE( iimppt(knbi,knbj), ijmppt(knbi,knbj), ijpi(knbi,knbj), ijpj(knbi,knbj) )
1100            CALL basic_decomposition( Ni0glo, Nj0glo, 0, knbi, knbj, iimax, ijmax, iimppt, ijmppt, ijpi, ijpj )
1101            !
1102            inx = Ni0glo + 2   ;   iny = ijpj(1,iarea) + 2             ! strip size + 1 halo on each direction (even if nn_hls>1)
1103            ALLOCATE( lloce(inx, iny) )                                ! allocate the strip
1104            inry = iny - COUNT( (/ iarea == 1, iarea == knbj /) )      ! number of point to read in y-direction
1105            isty = 1 + COUNT( (/ iarea == 1 /) )                       ! read from the first or the second line?
1106            CALL readbot_strip( ijmppt(1,iarea) - 2 + isty, inry, lloce(2:inx-1, isty:inry+isty-1) )   ! read the strip
1107            !
1108            IF( iarea == 1    ) THEN                                   ! the first line was not read
1109               IF( jperio == 2 .OR. jperio == 7 ) THEN                 !   north-south periodocity
1110                  CALL readbot_strip( Nj0glo, 1, lloce(2:inx-1, 1) )   !   read the last line -> first line of lloce
1111               ELSE
1112                  lloce(2:inx-1,  1) = .FALSE.                         !   closed boundary
1113               ENDIF
1114            ENDIF
1115            IF( iarea == knbj ) THEN                                   ! the last line was not read
1116               IF( jperio == 2 .OR. jperio == 7 ) THEN                 !   north-south periodocity
1117                  CALL readbot_strip( 1, 1, lloce(2:inx-1,iny) )       !   read the first line -> last line of lloce
1118               ELSE
1119                  lloce(2:inx-1,iny) = .FALSE.                         !   closed boundary
1120               ENDIF
1121            ENDIF
1122            !                                                          ! first and last column were not read
1123            IF( jperio == 1 .OR. jperio == 4 .OR. jperio == 6 .OR. jperio == 7 ) THEN
1124               lloce(1,:) = lloce(inx-1,:)   ;   lloce(inx,:) = lloce(2,:)   ! east-west periodocity
1125            ELSE
1126               lloce(1,:) = .FALSE.          ;   lloce(inx,:) = .FALSE.      ! closed boundary
1127            ENDIF
1128            !
1129            DO  ji = 1, knbi
1130               inboce(ji,iarea) = COUNT( lloce(iimppt(ji,1):iimppt(ji,1)+ijpi(ji,1)+1,:) )   ! lloce as 2 points more than Ni0glo
1131            END DO
1132            !
1133            DEALLOCATE(lloce)
1134            DEALLOCATE(iimppt, ijmppt, ijpi, ijpj)
1135            !
1136         ENDIF
1137      END DO
1138   
1139      inboce_1d = RESHAPE(inboce, (/ knbi*knbj /))
1140      CALL mpp_sum( 'mppini', inboce_1d )
1141      inboce = RESHAPE(inboce_1d, (/knbi, knbj/))
1142      ldisoce(:,:) = inboce(:,:) /= 0
1143      !
1144   END SUBROUTINE is_ocean
1145   
1146   
1147   SUBROUTINE readbot_strip( kjstr, kjcnt, ldoce )
1148      !!----------------------------------------------------------------------
1149      !!                  ***  ROUTINE readbot_strip  ***
1150      !!
1151      !! ** Purpose : Read relevant bathymetric information in order to
1152      !!              provide a land/sea mask used for the elimination
1153      !!              of land domains, in an mpp computation.
1154      !!
1155      !! ** Method  : read stipe of size (Ni0glo,...)
1156      !!----------------------------------------------------------------------
1157      INTEGER                         , INTENT(in   ) :: kjstr       ! starting j position of the reading
1158      INTEGER                         , INTENT(in   ) :: kjcnt       ! number of lines to read
1159      LOGICAL, DIMENSION(Ni0glo,kjcnt), INTENT(  out) :: ldoce       ! ldoce(i,j) = .true. if the point (i,j) is ocean
1160      !
1161      INTEGER                           ::   inumsave                ! local logical unit
1162      REAL(wp), DIMENSION(Ni0glo,kjcnt) ::   zbot, zbdy 
1163      !!----------------------------------------------------------------------
1164      !
1165      inumsave = numout   ;   numout = numnul   !   redirect all print to /dev/null
1166      !
1167      IF( numbot /= -1 ) THEN   
1168         CALL iom_get( numbot, jpdom_unknown, 'bottom_level', zbot, kstart = (/1,kjstr/), kcount = (/Ni0glo, kjcnt/) )
1169      ELSE
1170         zbot(:,:) = 1._wp                      ! put a non-null value
1171      ENDIF
1172      !
1173      IF( numbdy /= -1 ) THEN                   ! Adjust with bdy_msk if it exists   
1174         CALL iom_get ( numbdy, jpdom_unknown,     'bdy_msk', zbdy, kstart = (/1,kjstr/), kcount = (/Ni0glo, kjcnt/) )
1175         zbot(:,:) = zbot(:,:) * zbdy(:,:)
1176      ENDIF
1177      !
1178      ldoce(:,:) = zbot(:,:) > 0._wp
1179      numout = inumsave
1180      !
1181   END SUBROUTINE readbot_strip
1182
1183
1184   SUBROUTINE init_ioipsl
1185      !!----------------------------------------------------------------------
1186      !!                  ***  ROUTINE init_ioipsl  ***
1187      !!
1188      !! ** Purpose :   
1189      !!
1190      !! ** Method  :   
1191      !!
1192      !! History :
1193      !!   9.0  !  04-03  (G. Madec )  MPP-IOIPSL
1194      !!   " "  !  08-12  (A. Coward)  addition in case of jpni*jpnj < jpnij
1195      !!----------------------------------------------------------------------
1196      INTEGER, DIMENSION(2) ::   iglo, iloc, iabsf, iabsl, ihals, ihale, idid
1197      !!----------------------------------------------------------------------
1198
1199      ! The domain is split only horizontally along i- or/and j- direction
1200      ! So we need at the most only 1D arrays with 2 elements.
1201      ! Set idompar values equivalent to the jpdom_local_noextra definition
1202      ! used in IOM. This works even if jpnij .ne. jpni*jpnj.
1203      iglo( :) = (/ Ni0glo, Nj0glo /)
1204      iloc( :) = (/ Ni_0  , Nj_0   /)
1205      iabsf(:) = (/ Nis0  , Njs0   /) + (/ nimpp, njmpp /) - 1 - nn_hls   ! corresponds to mig0(Nis0) but mig0 is not yet defined!
1206      iabsl(:) = iabsf(:) + iloc(:) - 1
1207      ihals(:) = (/ 0     , 0      /)
1208      ihale(:) = (/ 0     , 0      /)
1209      idid( :) = (/ 1     , 2      /)
1210
1211      IF(lwp) THEN
1212          WRITE(numout,*)
1213          WRITE(numout,*) 'mpp init_ioipsl :   iloc  = ', iloc
1214          WRITE(numout,*) '~~~~~~~~~~~~~~~     iabsf = ', iabsf
1215          WRITE(numout,*) '                    ihals = ', ihals
1216          WRITE(numout,*) '                    ihale = ', ihale
1217      ENDIF
1218      !
1219      CALL flio_dom_set ( jpnij, nproc, idid, iglo, iloc, iabsf, iabsl, ihals, ihale, 'BOX', nidom)
1220      !
1221   END SUBROUTINE init_ioipsl 
1222
1223
1224   SUBROUTINE init_nfdcom
1225      !!----------------------------------------------------------------------
1226      !!                     ***  ROUTINE  init_nfdcom  ***
1227      !! ** Purpose :   Setup for north fold exchanges with explicit
1228      !!                point-to-point messaging
1229      !!
1230      !! ** Method :   Initialization of the northern neighbours lists.
1231      !!----------------------------------------------------------------------
1232      !!    1.0  ! 2011-10  (A. C. Coward, NOCS & J. Donners, PRACE)
1233      !!    2.0  ! 2013-06 Setup avoiding MPI communication (I. Epicoco, S. Mocavero, CMCC)
1234      !!----------------------------------------------------------------------
1235      INTEGER  ::   sxM, dxM, sxT, dxT, jn
1236      !!----------------------------------------------------------------------
1237      !
1238      !initializes the north-fold communication variables
1239      isendto(:) = 0
1240      nsndto     = 0
1241      !
1242      IF ( njmpp == MAXVAL( njmppt ) ) THEN      ! if I am a process in the north
1243         !
1244         !sxM is the first point (in the global domain) needed to compute the north-fold for the current process
1245         sxM = jpiglo - nimppt(narea) - jpiall(narea) + 1
1246         !dxM is the last point (in the global domain) needed to compute the north-fold for the current process
1247         dxM = jpiglo - nimppt(narea) + 2
1248         !
1249         ! loop over the other north-fold processes to find the processes
1250         ! managing the points belonging to the sxT-dxT range
1251         !
1252         DO jn = 1, jpni
1253            !
1254            sxT = nfimpp(jn)                    ! sxT = 1st  point (in the global domain) of the jn process
1255            dxT = nfimpp(jn) + nfjpi(jn) - 1    ! dxT = last point (in the global domain) of the jn process
1256            !
1257            IF    ( sxT < sxM  .AND.  sxM < dxT ) THEN
1258               nsndto          = nsndto + 1
1259               isendto(nsndto) = jn
1260            ELSEIF( sxM <= sxT  .AND.  dxM >= dxT ) THEN
1261               nsndto          = nsndto + 1
1262               isendto(nsndto) = jn
1263            ELSEIF( dxM <  dxT  .AND.  sxT <  dxM ) THEN
1264               nsndto          = nsndto + 1
1265               isendto(nsndto) = jn
1266            ENDIF
1267            !
1268         END DO
1269         !
1270      ENDIF
1271      l_north_nogather = .TRUE.
1272      !
1273   END SUBROUTINE init_nfdcom
1274
1275#endif
1276
1277   SUBROUTINE init_doloop
1278      !!----------------------------------------------------------------------
1279      !!                  ***  ROUTINE init_doloop  ***
1280      !!
1281      !! ** Purpose :   set the starting/ending indices of DO-loop
1282      !!              These indices are used in do_loop_substitute.h90
1283      !!----------------------------------------------------------------------
1284      !
1285      Nis0 =   1+nn_hls   ;   Nis1 = Nis0-1   ;   Nis2 = MAX(  1, Nis0-2)
1286      Njs0 =   1+nn_hls   ;   Njs1 = Njs0-1   ;   Njs2 = MAX(  1, Njs0-2) 
1287      !                                                 
1288      Nie0 = jpi-nn_hls   ;   Nie1 = Nie0+1   ;   Nie2 = MIN(jpi, Nie0+2)
1289      Nje0 = jpj-nn_hls   ;   Nje1 = Nje0+1   ;   Nje2 = MIN(jpj, Nje0+2)
1290      !
1291      IF( nn_hls == 1 ) THEN          !* halo size of 1
1292         !
1293         Nis1nxt2 = Nis0   ;   Njs1nxt2 = Njs0
1294         Nie1nxt2 = Nie0   ;   Nje1nxt2 = Nje0
1295         !
1296      ELSE                            !* larger halo size...
1297         !
1298         Nis1nxt2 = Nis1   ;   Njs1nxt2 = Njs1
1299         Nie1nxt2 = Nie1   ;   Nje1nxt2 = Nje1
1300         !
1301      ENDIF
1302      !
1303      Ni_0 = Nie0 - Nis0 + 1
1304      Nj_0 = Nje0 - Njs0 + 1
1305      Ni_1 = Nie1 - Nis1 + 1
1306      Nj_1 = Nje1 - Njs1 + 1
1307      Ni_2 = Nie2 - Nis2 + 1
1308      Nj_2 = Nje2 - Njs2 + 1
1309      !
1310   END SUBROUTINE init_doloop
1311   
1312   !!======================================================================
1313END MODULE mppini
Note: See TracBrowser for help on using the repository browser.