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/2021/ticket2680_C1D_PAPA/src/OCE/LBC – NEMO

source: NEMO/branches/2021/ticket2680_C1D_PAPA/src/OCE/LBC/mppini.F90 @ 14920

Last change on this file since 14920 was 14920, checked in by gsamson, 3 years ago

first steps to get C1D_PAPA working without MPI (key_mpi_off)

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