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

source: NEMO/branches/2021/dev_r14312_MPI_Interface/src/OCE/LBC/mppini.F90 @ 14363

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

dev_r14312_MPI_Interface: suppress communications involving only land points, #2598

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