New URL for NEMO forge! http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.

mppini.F90 in NEMO/trunk/src/OCE/LBC – NEMO

Context Navigation

source: NEMO/trunk/src/OCE/LBC/mppini.F90 @ 14072

Last change on this file since 14072 was 14072, checked in by laurent, 3 years ago
Merging branch "2020/dev_r13648_ASINTER-04_laurent_bulk_ice", ticket #2369
Property svn:keywords set to `Id`
File size: 68.1 KB

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

Note: See TracBrowser for help on using the repository browser.

Download in other formats: