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 utils/tools/DOMAINcfg/src – NEMO

source: utils/tools/DOMAINcfg/src/mppini.F90

Last change on this file was 15279, checked in by jchanut, 3 years ago

#2222 and #2638: Enable creating agrif meshes with different vertical grids (geopotential only as a start)

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