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 @ 14325

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

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

Download in other formats: