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

Last change on this file since 13229 was 13229, checked in by francesca, 3 months ago

dev_r12558_HPC-08_epico_Extra_Halo: merge with trunk@13218, see #2366

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