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

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

mppini.F90 in NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/OCE/LBC – NEMO

Context Navigation

source: NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/OCE/LBC/mppini.F90 @ 10345

Last change on this file since 10345 was 10345, checked in by smasson, 5 years ago
dev_r10164_HPC09_ESIWACE_PREP_MERGE: merge with trunk@10344, see #2133
Property svn:keywords set to `Id`
File size: 53.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
[9436]	10	!! 3.4 ! 2011-10 (A. C. Coward, NOCS & J. Donners, PRACE) add mpp_init_nfdcom
	11	!! 3. ! 2013-06 (I. Epicoco, S. Mocavero, CMCC) mpp_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	!!----------------------------------------------------------------------
[9436]	17	!! mpp_init : Lay out the global domain over processors with/without land processor elimination
	18	!! mpp_init_mask : Read global bathymetric information to facilitate land suppression
	19	!! mpp_init_ioipsl : IOIPSL initialization in mpp
	20	!! mpp_init_partition: Calculate MPP domain decomposition
	21	!! factorise : Calculate the factors of the no. of MPI processes
	22	!! mpp_init_nfdcom : Setup for north fold exchanges with explicit point-to-point messaging
[3]	23	!!----------------------------------------------------------------------
[9019]	24	USE dom_oce ! ocean space and time domain
	25	USE bdy_oce ! open BounDarY
	26	!
[9436]	27	USE lbcnfd , ONLY : isendto, nsndto, nfsloop, nfeloop ! Setup of north fold exchanges
[9019]	28	USE lib_mpp ! distribued memory computing library
	29	USE iom ! nemo I/O library
	30	USE ioipsl ! I/O IPSL library
	31	USE in_out_manager ! I/O Manager
[3]	32
	33	IMPLICIT NONE
	34	PRIVATE
	35
	36	PUBLIC mpp_init ! called by opa.F90
	37
[10330]	38	INTEGER :: numbot = -1 ! 'bottom_level' local logical unit
	39	INTEGER :: numbdy = -1 ! 'bdy_msk' local logical unit
	40
[3]	41	!!----------------------------------------------------------------------
[9598]	42	!! NEMO/OCE 4.0 , NEMO Consortium (2018)
[1152]	43	!! $Id$
[10068]	44	!! Software governed by the CeCILL license (see ./LICENSE)
[3]	45	!!----------------------------------------------------------------------
	46	CONTAINS
	47
[1556]	48	#if ! defined key_mpp_mpi
[3]	49	!!----------------------------------------------------------------------
	50	!! Default option : shared memory computing
	51	!!----------------------------------------------------------------------
	52
	53	SUBROUTINE mpp_init
	54	!!----------------------------------------------------------------------
	55	!! * ROUTINE mpp_init *
	56	!!
	57	!! ** Purpose : Lay out the global domain over processors.
	58	!!
	59	!! ** Method : Shared memory computing, set the local processor
[9019]	60	!! variables to the value of the global domain
[3]	61	!!----------------------------------------------------------------------
[9019]	62	!
[9436]	63	jpimax = jpiglo
	64	jpjmax = jpjglo
	65	jpi = jpiglo
	66	jpj = jpjglo
[9799]	67	jpk = jpkglo
[9436]	68	jpim1 = jpi-1 ! inner domain indices
	69	jpjm1 = jpj-1 ! " "
	70	jpkm1 = MAX( 1, jpk-1 ) ! " "
	71	jpij = jpi*jpj
	72	jpni = 1
	73	jpnj = 1
	74	jpnij = jpni*jpnj
[9019]	75	nimpp = 1 !
[3]	76	njmpp = 1
	77	nlci = jpi
	78	nlcj = jpj
	79	nldi = 1
	80	nldj = 1
	81	nlei = jpi
	82	nlej = jpj
	83	nbondi = 2
	84	nbondj = 2
[352]	85	nidom = FLIO_DOM_NONE
[1344]	86	npolj = jperio
[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	!
[3]	90	IF(lwp) THEN
	91	WRITE(numout,*)
[9019]	92	WRITE(numout,*) 'mpp_init : NO massively parallel processing'
	93	WRITE(numout,*) '~~~~~~~~ '
[9667]	94	WRITE(numout,*) ' l_Iperio = ', l_Iperio, ' l_Jperio = ', l_Jperio
	95	WRITE(numout,*) ' npolj = ', npolj , ' njmpp = ', njmpp
[3]	96	ENDIF
[9019]	97	!
	98	IF( jpni /= 1 .OR. jpnj /= 1 .OR. jpnij /= 1 ) &
	99	CALL ctl_stop( 'mpp_init: equality jpni = jpnj = jpnij = 1 is not satisfied', &
	100	& 'the domain is lay out for distributed memory computing!' )
	101	!
[3]	102	END SUBROUTINE mpp_init
	103
	104	#else
	105	!!----------------------------------------------------------------------
[9019]	106	!! 'key_mpp_mpi' MPI massively parallel processing
[3]	107	!!----------------------------------------------------------------------
	108
[9449]	109
[3]	110	SUBROUTINE mpp_init
	111	!!----------------------------------------------------------------------
	112	!! * ROUTINE mpp_init *
	113	!!
	114	!! ** Purpose : Lay out the global domain over processors.
[9019]	115	!! If land processors are to be eliminated, this program requires the
	116	!! presence of the domain configuration file. Land processors elimination
	117	!! is performed if jpni x jpnj /= jpnij. In this case, using the MPP_PREP
	118	!! preprocessing tool, help for defining the best cutting out.
[3]	119	!!
	120	!! ** Method : Global domain is distributed in smaller local domains.
	121	!! Periodic condition is a function of the local domain position
	122	!! (global boundary or neighbouring domain) and of the global
	123	!! periodic
	124	!! Type : jperio global periodic condition
	125	!!
[9019]	126	!! ** Action : - set domain parameters
[3]	127	!! nimpp : longitudinal index
	128	!! njmpp : latitudinal index
	129	!! narea : number for local area
	130	!! nlci : first dimension
	131	!! nlcj : second dimension
	132	!! nbondi : mark for "east-west local boundary"
	133	!! nbondj : mark for "north-south local boundary"
	134	!! nproc : number for local processor
	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
[10345]	141	INTEGER :: inijmin
[10330]	142	INTEGER :: i2add
[9019]	143	INTEGER :: inum ! local logical unit
[10330]	144	INTEGER :: idir, ifreq, icont ! 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 ! - -
[10330]	149	INTEGER :: iarea0 ! - -
	150	INTEGER :: ierr, ios !
[9444]	151	INTEGER, ALLOCATABLE, DIMENSION(:) :: iin, ii_nono, ii_noea ! 1D workspace
	152	INTEGER, ALLOCATABLE, DIMENSION(:) :: ijn, ii_noso, ii_nowe ! - -
	153	INTEGER, ALLOCATABLE, DIMENSION(:,:) :: iimppt, ilci, ibondi, ipproc ! 2D workspace
	154	INTEGER, ALLOCATABLE, DIMENSION(:,:) :: ijmppt, ilcj, ibondj, ipolj ! - -
	155	INTEGER, ALLOCATABLE, DIMENSION(:,:) :: ilei, ildi, iono, ioea ! - -
	156	INTEGER, ALLOCATABLE, DIMENSION(:,:) :: ilej, ildj, ioso, iowe ! - -
[10330]	157	LOGICAL, ALLOCATABLE, DIMENSION(:,:) :: llisoce ! - -
	158	NAMELIST/nambdy/ ln_bdy, nb_bdy, ln_coords_file, cn_coords_file, &
	159	& ln_mask_file, cn_mask_file, cn_dyn2d, nn_dyn2d_dta, &
	160	& cn_dyn3d, nn_dyn3d_dta, cn_tra, nn_tra_dta, &
	161	& ln_tra_dmp, ln_dyn3d_dmp, rn_time_dmp, rn_time_dmp_out, &
	162	& cn_ice, nn_ice_dta, &
	163	& rn_ice_tem, rn_ice_sal, rn_ice_age, &
	164	& ln_vol, nn_volctl, nn_rimwidth, nb_jpk_bdy
[3]	165	!!----------------------------------------------------------------------
[9436]	166
[10330]	167	! do we need to take into account bdy_msk?
	168	REWIND( numnam_ref ) ! Namelist nambdy in reference namelist : BDY
	169	READ ( numnam_ref, nambdy, IOSTAT = ios, ERR = 903)
	170	903 IF( ios /= 0 ) CALL ctl_nam ( ios , 'nambdy in reference namelist (mppini)', lwp )
	171	REWIND( numnam_cfg ) ! Namelist nambdy in configuration namelist : BDY
	172	READ ( numnam_cfg, nambdy, IOSTAT = ios, ERR = 904 )
	173	904 IF( ios > 0 ) CALL ctl_nam ( ios , 'nambdy in configuration namelist (mppini)', lwp )
	174	!
	175	IF( ln_read_cfg ) CALL iom_open( cn_domcfg, numbot )
	176	IF( ln_bdy .AND. ln_mask_file ) CALL iom_open( cn_mask_file, numbdy )
	177	!
	178	! 1. Dimension arrays for subdomains
	179	! -----------------------------------
	180	!
[9436]	181	! If dimensions of processor grid weren't specified in the namelist file
	182	! then we calculate them here now that we have our communicator size
[10330]	183	IF( jpni < 1 .OR. jpnj < 1 ) CALL mpp_init_bestpartition( mppsize, jpni, jpnj )
	184
	185	! look for land mpi subdomains...
	186	ALLOCATE( llisoce(jpni,jpnj) )
	187	CALL mpp_init_isoce( jpni, jpnj, llisoce )
	188	inijmin = COUNT( llisoce ) ! number of oce subdomains
	189
	190	IF( mppsize < inijmin ) THEN
	191	WRITE(ctmp1,*) ' With this specified domain decomposition: jpni =', jpni, ' jpnj =', jpnj
	192	WRITE(ctmp2,) ' we can eliminate only ', jpnijpnj - inijmin, ' land mpi subdomains therefore '
	193	WRITE(ctmp3,*) ' the number of ocean mpi subdomains (', inijmin,') exceed the number of MPI processes:', mppsize
	194	WRITE(ctmp4,*) ' ==>>> There is the list of best domain decompositions you should use: '
	195	CALL ctl_stop( 'mpp_init:', '~~~~~~~~ ', ctmp1, ctmp2, ctmp3, ctmp4 )
	196	CALL mpp_init_bestpartition( mppsize, ldlist = .TRUE. ) ! must be done by all core
	197	CALL ctl_stop( 'STOP' )
	198	ENDIF
	199
	200	IF( mppsize > jpni*jpnj ) THEN
	201	WRITE(ctmp1,*) ' The number of mpi processes: ', mppsize
	202	WRITE(ctmp2,) ' exceeds the maximum number of subdomains (ocean+land) = ', jpnijpnj
	203	WRITE(ctmp3,*) ' defined by the following domain decomposition: jpni =', jpni, ' jpnj =', jpnj
	204	WRITE(ctmp4,*) ' ==>>> There is the list of best domain decompositions you should use: '
	205	CALL ctl_stop( 'mpp_init:', '~~~~~~~~ ', ctmp1, ctmp2, ctmp3, ctmp4 )
	206	CALL mpp_init_bestpartition( mppsize, ldlist = .TRUE. ) ! must be done by all core
	207	CALL ctl_stop( 'STOP' )
	208	ENDIF
	209
	210	jpnij = mppsize ! force jpnij definition <-- remove as much land subdomains as needed to reach this condition
	211	IF( mppsize > inijmin ) THEN
	212	WRITE(ctmp1,*) ' The number of mpi processes: ', mppsize
	213	WRITE(ctmp2,*) ' exceeds the maximum number of ocean subdomains = ', inijmin
	214	WRITE(ctmp3,) ' we suppressed ', jpnijpnj - mppsize, ' land subdomains '
	215	WRITE(ctmp4,*) ' BUT we had to keep ', mppsize - inijmin, ' land subdomains that are useless...'
	216	CALL ctl_warn( 'mpp_init:', '~~~~~~~~ ', ctmp1, ctmp2, ctmp3, ctmp4 )
	217	ELSE ! mppsize = inijmin
	218	IF(lwp) THEN
	219	WRITE(numout,*) 'mpp_init: You use an optimal domain decomposition'
	220	WRITE(numout,*) '~~~~~~~~ '
	221	WRITE(numout,*) ' Number of mpi processes: ', mppsize
	222	WRITE(numout,*) ' Number of ocean subdomains = ', inijmin
	223	WRITE(numout,) ' Number of suppressed land subdomains = ', jpnijpnj - inijmin
	224	ENDIF
	225	ENDIF
	226
	227	IF( numbot /= -1 ) CALL iom_close( numbot )
	228	IF( numbdy /= -1 ) CALL iom_close( numbdy )
	229
[9436]	230	ALLOCATE( nfiimpp(jpni,jpnj), nfipproc(jpni,jpnj), nfilcit(jpni,jpnj) , &
	231	& nimppt(jpnij) , ibonit(jpnij) , nlcit(jpnij) , nlcjt(jpnij) , &
	232	& njmppt(jpnij) , ibonjt(jpnij) , nldit(jpnij) , nldjt(jpnij) , &
[9444]	233	& nleit(jpnij) , nlejt(jpnij) , &
	234	& iin(jpnij), ii_nono(jpnij), ii_noea(jpnij), &
	235	& ijn(jpnij), ii_noso(jpnij), ii_nowe(jpnij), &
	236	& iimppt(jpni,jpnj), ilci(jpni,jpnj), ibondi(jpni,jpnj), ipproc(jpni,jpnj), &
	237	& ijmppt(jpni,jpnj), ilcj(jpni,jpnj), ibondj(jpni,jpnj), ipolj(jpni,jpnj), &
	238	& ilei(jpni,jpnj), ildi(jpni,jpnj), iono(jpni,jpnj), ioea(jpni,jpnj), &
	239	& ilej(jpni,jpnj), ildj(jpni,jpnj), ioso(jpni,jpnj), iowe(jpni,jpnj), &
	240	& STAT=ierr )
[10297]	241	CALL mpp_sum( 'mppini', ierr )
[9436]	242	IF( ierr /= 0 ) CALL ctl_stop( 'STOP', 'mpp_init: unable to allocate standard ocean arrays' )
	243
[9508]	244	#if defined key_agrif
[9436]	245	IF( .NOT. Agrif_Root() ) THEN ! AGRIF children: specific setting (cf. agrif_user.F90)
[9449]	246	IF( jpiglo /= nbcellsx + 2 + 2*nbghostcells ) &
	247	CALL ctl_stop( 'STOP', 'mpp_init: Agrif children requires jpiglo == nbcellsx + 2 + 2*nbghostcells' )
	248	IF( jpjglo /= nbcellsy + 2 + 2*nbghostcells ) &
	249	CALL ctl_stop( 'STOP', 'mpp_init: Agrif children requires jpjglo == nbcellsy + 2 + 2*nbghostcells' )
	250	IF( ln_use_jattr ) CALL ctl_stop( 'STOP', 'mpp_init:Agrif children requires ln_use_jattr = .false. ' )
[9436]	251	ENDIF
[9508]	252	#endif
[9436]	253	!
[10330]	254	! 2. Index arrays for subdomains
[3]	255	! -----------------------------------
[9019]	256	!
	257	nreci = 2 * nn_hls
	258	nrecj = 2 * nn_hls
[10330]	259	CALL mpp_basic_decomposition( jpni, jpnj, jpimax, jpjmax, iimppt, ijmppt, ilci, ilcj )
	260	nfiimpp(:,:) = iimppt(:,:)
	261	nfilcit(:,:) = ilci(:,:)
[9019]	262	!
	263	IF(lwp) THEN
	264	WRITE(numout,*)
[10330]	265	WRITE(numout,*) 'MPI Message Passing MPI - domain lay out over processors'
	266	WRITE(numout,*)
[9019]	267	WRITE(numout,*) ' defines mpp subdomains'
[10330]	268	WRITE(numout,*) ' jpni = ', jpni
	269	WRITE(numout,*) ' jpnj = ', jpnj
[9019]	270	WRITE(numout,*)
[10330]	271	WRITE(numout,*) ' sum ilci(i,1) = ', sum(ilci(:,1)), ' jpiglo = ', jpiglo
	272	WRITE(numout,*) ' sum ilcj(1,j) = ', sum(ilcj(1,:)), ' jpjglo = ', jpjglo
[9019]	273	ENDIF
[10330]	274
[9019]	275	! 3. Subdomain description in the Regular Case
	276	! --------------------------------------------
[9667]	277	! specific cases where there is no communication -> must do the periodicity by itself
	278	! Warning: because of potential land-area suppression, do not use nbond[ij] == 2
	279	l_Iperio = jpni == 1 .AND. (jperio == 1 .OR. jperio == 4 .OR. jperio == 6 .OR. jperio == 7)
	280	l_Jperio = jpnj == 1 .AND. (jperio == 2 .OR. jperio == 7)
	281
[9019]	282	DO jarea = 1, jpni*jpnj
[10330]	283	!
[9564]	284	iarea0 = jarea - 1
	285	ii = 1 + MOD(iarea0,jpni)
	286	ij = 1 + iarea0/jpni
[9019]	287	ili = ilci(ii,ij)
	288	ilj = ilcj(ii,ij)
[9667]	289	ibondi(ii,ij) = 0 ! default: has e-w neighbours
	290	IF( ii == 1 ) ibondi(ii,ij) = -1 ! first column, has only e neighbour
	291	IF( ii == jpni ) ibondi(ii,ij) = 1 ! last column, has only w neighbour
	292	IF( jpni == 1 ) ibondi(ii,ij) = 2 ! has no e-w neighbour
	293	ibondj(ii,ij) = 0 ! default: has n-s neighbours
	294	IF( ij == 1 ) ibondj(ii,ij) = -1 ! first row, has only n neighbour
	295	IF( ij == jpnj ) ibondj(ii,ij) = 1 ! last row, has only s neighbour
	296	IF( jpnj == 1 ) ibondj(ii,ij) = 2 ! has no n-s neighbour
[9019]	297
[9436]	298	! Subdomain neighbors (get their zone number): default definition
[9564]	299	ioso(ii,ij) = iarea0 - jpni
	300	iowe(ii,ij) = iarea0 - 1
	301	ioea(ii,ij) = iarea0 + 1
	302	iono(ii,ij) = iarea0 + jpni
[9019]	303	ildi(ii,ij) = 1 + nn_hls
	304	ilei(ii,ij) = ili - nn_hls
	305	ildj(ii,ij) = 1 + nn_hls
	306	ilej(ii,ij) = ilj - nn_hls
	307
[9436]	308	! East-West periodicity: change ibondi, ioea, iowe
[9069]	309	IF( jperio == 1 .OR. jperio == 4 .OR. jperio == 6 .OR. jperio == 7 ) THEN
[9667]	310	IF( jpni /= 1 ) ibondi(ii,ij) = 0 ! redefine: all have e-w neighbours
	311	IF( ii == 1 ) iowe(ii,ij) = iarea0 + (jpni-1) ! redefine: first column, address of w neighbour
	312	IF( ii == jpni ) ioea(ii,ij) = iarea0 - (jpni-1) ! redefine: last column, address of e neighbour
[9019]	313	ENDIF
[9436]	314
[9667]	315	! Simple North-South periodicity: change ibondj, ioso, iono
	316	IF( jperio == 2 .OR. jperio == 7 ) THEN
	317	IF( jpnj /= 1 ) ibondj(ii,ij) = 0 ! redefine: all have n-s neighbours
	318	IF( ij == 1 ) ioso(ii,ij) = iarea0 + jpni * (jpnj-1) ! redefine: first row, address of s neighbour
	319	IF( ij == jpnj ) iono(ii,ij) = iarea0 - jpni * (jpnj-1) ! redefine: last row, address of n neighbour
	320	ENDIF
	321
[9436]	322	! North fold: define ipolj, change iono. Warning: we do not change ibondj...
[9019]	323	ipolj(ii,ij) = 0
	324	IF( jperio == 3 .OR. jperio == 4 ) THEN
	325	ijm1 = jpni*(jpnj-1)
	326	imil = ijm1+(jpni+1)/2
	327	IF( jarea > ijm1 ) ipolj(ii,ij) = 3
	328	IF( MOD(jpni,2) == 1 .AND. jarea == imil ) ipolj(ii,ij) = 4
	329	IF( ipolj(ii,ij) == 3 ) iono(ii,ij) = jpni*jpnj-jarea+ijm1 ! MPI rank of northern neighbour
	330	ENDIF
	331	IF( jperio == 5 .OR. jperio == 6 ) THEN
	332	ijm1 = jpni*(jpnj-1)
	333	imil = ijm1+(jpni+1)/2
	334	IF( jarea > ijm1) ipolj(ii,ij) = 5
	335	IF( MOD(jpni,2) == 1 .AND. jarea == imil ) ipolj(ii,ij) = 6
	336	IF( ipolj(ii,ij) == 5) iono(ii,ij) = jpni*jpnj-jarea+ijm1 ! MPI rank of northern neighbour
	337	ENDIF
	338	!
[10330]	339	END DO
	340
	341	! 4. deal with land subdomains
	342	! ----------------------------
	343	!
	344	! specify which subdomains are oce subdomains; other are land subdomains
	345	ipproc(:,:) = -1
	346	icont = -1
	347	DO jarea = 1, jpni*jpnj
	348	iarea0 = jarea - 1
	349	ii = 1 + MOD(iarea0,jpni)
	350	ij = 1 + iarea0/jpni
	351	IF( llisoce(ii,ij) ) THEN
[9019]	352	icont = icont + 1
	353	ipproc(ii,ij) = icont
	354	iin(icont+1) = ii
	355	ijn(icont+1) = ij
	356	ENDIF
[3]	357	END DO
[10330]	358	! if needed add some land subdomains to reach jpnij active subdomains
	359	i2add = jpnij - inijmin
	360	DO jarea = 1, jpni*jpnj
	361	iarea0 = jarea - 1
	362	ii = 1 + MOD(iarea0,jpni)
	363	ij = 1 + iarea0/jpni
	364	IF( .NOT. llisoce(ii,ij) .AND. i2add > 0 ) THEN
	365	icont = icont + 1
	366	ipproc(ii,ij) = icont
	367	iin(icont+1) = ii
	368	ijn(icont+1) = ij
	369	i2add = i2add - 1
	370	ENDIF
	371	END DO
[9019]	372	nfipproc(:,:) = ipproc(:,:)
[6412]	373
[10330]	374	! neighbour treatment: change ibondi, ibondj if next to a land zone
[9019]	375	DO jarea = 1, jpni*jpnj
	376	ii = 1 + MOD( jarea-1 , jpni )
	377	ij = 1 + (jarea-1) / jpni
[9667]	378	! land-only area with an active n neigbour
[9019]	379	IF ( ipproc(ii,ij) == -1 .AND. 0 <= iono(ii,ij) .AND. iono(ii,ij) <= jpni*jpnj-1 ) THEN
[9667]	380	iino = 1 + MOD( iono(ii,ij) , jpni ) ! ii index of this n neigbour
	381	ijno = 1 + iono(ii,ij) / jpni ! ij index of this n neigbour
	382	! In case of north fold exchange: I am the n neigbour of my n neigbour!! (#1057)
	383	! --> for northern neighbours of northern row processors (in case of north-fold)
	384	! need to reverse the LOGICAL direction of communication
	385	idir = 1 ! we are indeed the s neigbour of this n neigbour
	386	IF( ij == jpnj .AND. ijno == jpnj ) idir = -1 ! both are on the last row, we are in fact the n neigbour
	387	IF( ibondj(iino,ijno) == idir ) ibondj(iino,ijno) = 2 ! this n neigbour had only a s/n neigbour -> no more
	388	IF( ibondj(iino,ijno) == 0 ) ibondj(iino,ijno) = -idir ! this n neigbour had both, n-s neighbours -> keep 1
[9019]	389	ENDIF
[9667]	390	! land-only area with an active s neigbour
[9019]	391	IF( ipproc(ii,ij) == -1 .AND. 0 <= ioso(ii,ij) .AND. ioso(ii,ij) <= jpni*jpnj-1 ) THEN
[9667]	392	iiso = 1 + MOD( ioso(ii,ij) , jpni ) ! ii index of this s neigbour
	393	ijso = 1 + ioso(ii,ij) / jpni ! ij index of this s neigbour
	394	IF( ibondj(iiso,ijso) == -1 ) ibondj(iiso,ijso) = 2 ! this s neigbour had only a n neigbour -> no more neigbour
	395	IF( ibondj(iiso,ijso) == 0 ) ibondj(iiso,ijso) = 1 ! this s neigbour had both, n-s neighbours -> keep s neigbour
[9019]	396	ENDIF
[9667]	397	! land-only area with an active e neigbour
[9019]	398	IF( ipproc(ii,ij) == -1 .AND. 0 <= ioea(ii,ij) .AND. ioea(ii,ij) <= jpni*jpnj-1 ) THEN
[9667]	399	iiea = 1 + MOD( ioea(ii,ij) , jpni ) ! ii index of this e neigbour
	400	ijea = 1 + ioea(ii,ij) / jpni ! ij index of this e neigbour
	401	IF( ibondi(iiea,ijea) == 1 ) ibondi(iiea,ijea) = 2 ! this e neigbour had only a w neigbour -> no more neigbour
	402	IF( ibondi(iiea,ijea) == 0 ) ibondi(iiea,ijea) = -1 ! this e neigbour had both, e-w neighbours -> keep e neigbour
[9019]	403	ENDIF
[9667]	404	! land-only area with an active w neigbour
[9019]	405	IF( ipproc(ii,ij) == -1 .AND. 0 <= iowe(ii,ij) .AND. iowe(ii,ij) <= jpni*jpnj-1) THEN
[9667]	406	iiwe = 1 + MOD( iowe(ii,ij) , jpni ) ! ii index of this w neigbour
	407	ijwe = 1 + iowe(ii,ij) / jpni ! ij index of this w neigbour
	408	IF( ibondi(iiwe,ijwe) == -1 ) ibondi(iiwe,ijwe) = 2 ! this w neigbour had only a e neigbour -> no more neigbour
	409	IF( ibondi(iiwe,ijwe) == 0 ) ibondi(iiwe,ijwe) = 1 ! this w neigbour had both, e-w neighbours -> keep w neigbour
[9019]	410	ENDIF
	411	END DO
[3]	412
[9436]	413	! Update il[de][ij] according to modified ibond[ij]
	414	! ----------------------
[9564]	415	DO jproc = 1, jpnij
	416	ii = iin(jproc)
	417	ij = ijn(jproc)
[9436]	418	IF( ibondi(ii,ij) == -1 .OR. ibondi(ii,ij) == 2 ) ildi(ii,ij) = 1
	419	IF( ibondi(ii,ij) == 1 .OR. ibondi(ii,ij) == 2 ) ilei(ii,ij) = ilci(ii,ij)
	420	IF( ibondj(ii,ij) == -1 .OR. ibondj(ii,ij) == 2 ) ildj(ii,ij) = 1
	421	IF( ibondj(ii,ij) == 1 .OR. ibondj(ii,ij) == 2 ) ilej(ii,ij) = ilcj(ii,ij)
	422	END DO
[10330]	423
	424	! 5. Subdomain print
	425	! ------------------
	426	IF(lwp) THEN
	427	ifreq = 4
	428	il1 = 1
	429	DO jn = 1, (jpni-1)/ifreq+1
	430	il2 = MIN(jpni,il1+ifreq-1)
	431	WRITE(numout,*)
	432	WRITE(numout,9400) ('***',ji=il1,il2-1)
	433	DO jj = jpnj, 1, -1
	434	WRITE(numout,9403) (' ',ji=il1,il2-1)
	435	WRITE(numout,9402) jj, (ilci(ji,jj),ilcj(ji,jj),ji=il1,il2)
	436	WRITE(numout,9404) (ipproc(ji,jj),ji=il1,il2)
	437	WRITE(numout,9403) (' ',ji=il1,il2-1)
	438	WRITE(numout,9400) ('***',ji=il1,il2-1)
	439	END DO
	440	WRITE(numout,9401) (ji,ji=il1,il2)
	441	il1 = il1+ifreq
	442	END DO
	443	9400 FORMAT(' *' ,20('***********',a3) )
	444	9403 FORMAT(' * ',20(' * ',a3) )
	445	9401 FORMAT(' ' ,20(' ',i3,' ') )
	446	9402 FORMAT(' ',i3,' * ',20(i3,' x',i3,' * ') )
	447	9404 FORMAT(' * ' ,20(' ',i3,' * ') )
	448	ENDIF
[9436]	449
[9019]	450	! just to save nono etc for all proc
[9436]	451	! warning ii*ij (zone) /= nproc (processors)!
	452	! ioso = zone number, ii_noso = proc number
[9019]	453	ii_noso(:) = -1
	454	ii_nono(:) = -1
	455	ii_noea(:) = -1
	456	ii_nowe(:) = -1
[9449]	457	DO jproc = 1, jpnij
	458	ii = iin(jproc)
	459	ij = ijn(jproc)
[9019]	460	IF( 0 <= ioso(ii,ij) .AND. ioso(ii,ij) <= (jpni*jpnj-1) ) THEN
	461	iiso = 1 + MOD( ioso(ii,ij) , jpni )
	462	ijso = 1 + ioso(ii,ij) / jpni
[9449]	463	ii_noso(jproc) = ipproc(iiso,ijso)
[9019]	464	ENDIF
	465	IF( 0 <= iowe(ii,ij) .AND. iowe(ii,ij) <= (jpni*jpnj-1) ) THEN
	466	iiwe = 1 + MOD( iowe(ii,ij) , jpni )
	467	ijwe = 1 + iowe(ii,ij) / jpni
[9449]	468	ii_nowe(jproc) = ipproc(iiwe,ijwe)
[9019]	469	ENDIF
	470	IF( 0 <= ioea(ii,ij) .AND. ioea(ii,ij) <= (jpni*jpnj-1) ) THEN
	471	iiea = 1 + MOD( ioea(ii,ij) , jpni )
	472	ijea = 1 + ioea(ii,ij) / jpni
[9449]	473	ii_noea(jproc)= ipproc(iiea,ijea)
[9019]	474	ENDIF
	475	IF( 0 <= iono(ii,ij) .AND. iono(ii,ij) <= (jpni*jpnj-1) ) THEN
	476	iino = 1 + MOD( iono(ii,ij) , jpni )
	477	ijno = 1 + iono(ii,ij) / jpni
[9449]	478	ii_nono(jproc)= ipproc(iino,ijno)
[9019]	479	ENDIF
	480	END DO
	481
	482	! 6. Change processor name
	483	! ------------------------
	484	ii = iin(narea)
	485	ij = ijn(narea)
	486	!
	487	! set default neighbours
	488	noso = ii_noso(narea)
	489	nowe = ii_nowe(narea)
	490	noea = ii_noea(narea)
	491	nono = ii_nono(narea)
	492	nlci = ilci(ii,ij)
	493	nldi = ildi(ii,ij)
	494	nlei = ilei(ii,ij)
[9436]	495	nlcj = ilcj(ii,ij)
[9019]	496	nldj = ildj(ii,ij)
	497	nlej = ilej(ii,ij)
	498	nbondi = ibondi(ii,ij)
	499	nbondj = ibondj(ii,ij)
	500	nimpp = iimppt(ii,ij)
[9436]	501	njmpp = ijmppt(ii,ij)
	502	jpi = nlci
	503	jpj = nlcj
	504	jpk = jpkglo ! third dim
	505	#if defined key_agrif
	506	! simple trick to use same vertical grid as parent but different number of levels:
	507	! Save maximum number of levels in jpkglo, then define all vertical grids with this number.
	508	! Suppress once vertical online interpolation is ok
[9449]	509	!!$ IF(.NOT.Agrif_Root()) jpkglo = Agrif_Parent( jpkglo )
[9436]	510	#endif
	511	jpim1 = jpi-1 ! inner domain indices
	512	jpjm1 = jpj-1 ! " "
	513	jpkm1 = MAX( 1, jpk-1 ) ! " "
	514	jpij = jpi*jpj ! jpi x j
[9019]	515	DO jproc = 1, jpnij
	516	ii = iin(jproc)
	517	ij = ijn(jproc)
	518	nlcit(jproc) = ilci(ii,ij)
	519	nldit(jproc) = ildi(ii,ij)
	520	nleit(jproc) = ilei(ii,ij)
[9436]	521	nlcjt(jproc) = ilcj(ii,ij)
[9019]	522	nldjt(jproc) = ildj(ii,ij)
	523	nlejt(jproc) = ilej(ii,ij)
[9436]	524	ibonit(jproc) = ibondi(ii,ij)
	525	ibonjt(jproc) = ibondj(ii,ij)
	526	nimppt(jproc) = iimppt(ii,ij)
	527	njmppt(jproc) = ijmppt(ii,ij)
[9019]	528	END DO
[3]	529
[9019]	530	! Save processor layout in ascii file
	531	IF (lwp) THEN
[7646]	532	CALL ctl_opn( inum, 'layout.dat', 'REPLACE', 'FORMATTED', 'SEQUENTIAL', -1, numout, .FALSE., narea )
[9019]	533	WRITE(inum,'(a)') ' jpnij jpimax jpjmax jpk jpiglo jpjglo'//&
[9436]	534	& ' ( local: narea jpi jpj )'
[9019]	535	WRITE(inum,'(6i8,a,3i8,a)') jpnij,jpimax,jpjmax,jpk,jpiglo,jpjglo,&
	536	& ' ( local: ',narea,jpi,jpj,' )'
[9436]	537	WRITE(inum,'(a)') 'nproc nlci nlcj nldi nldj nlei nlej nimp njmp nono noso nowe noea nbondi nbondj '
[9019]	538
	539	DO jproc = 1, jpnij
[9436]	540	WRITE(inum,'(13i5,2i7)') jproc-1, nlcit (jproc), nlcjt (jproc), &
	541	& nldit (jproc), nldjt (jproc), &
	542	& nleit (jproc), nlejt (jproc), &
	543	& nimppt (jproc), njmppt (jproc), &
	544	& ii_nono(jproc), ii_noso(jproc), &
	545	& ii_nowe(jproc), ii_noea(jproc), &
[9446]	546	& ibonit (jproc), ibonjt (jproc)
[7646]	547	END DO
[3]	548	END IF
	549
[9019]	550	! ! north fold parameter
	551	! Defined npolj, either 0, 3 , 4 , 5 , 6
	552	! In this case the important thing is that npolj /= 0
	553	! Because if we go through these line it is because jpni >1 and thus
	554	! we must use lbcnorthmpp, which tests only npolj =0 or npolj /= 0
[3]	555	npolj = 0
[9019]	556	ij = ijn(narea)
[3]	557	IF( jperio == 3 .OR. jperio == 4 ) THEN
[9019]	558	IF( ij == jpnj ) npolj = 3
[3]	559	ENDIF
	560	IF( jperio == 5 .OR. jperio == 6 ) THEN
[9019]	561	IF( ij == jpnj ) npolj = 5
[3]	562	ENDIF
[9019]	563	!
[9436]	564	nproc = narea-1
[3]	565	IF(lwp) THEN
[6412]	566	WRITE(numout,*)
[9169]	567	WRITE(numout,*) ' resulting internal parameters : '
	568	WRITE(numout,*) ' nproc = ', nproc
	569	WRITE(numout,*) ' nowe = ', nowe , ' noea = ', noea
	570	WRITE(numout,*) ' nono = ', nono , ' noso = ', noso
	571	WRITE(numout,*) ' nbondi = ', nbondi
	572	WRITE(numout,*) ' nbondj = ', nbondj
	573	WRITE(numout,*) ' npolj = ', npolj
[9667]	574	WRITE(numout,*) ' l_Iperio = ', l_Iperio
	575	WRITE(numout,*) ' l_Jperio = ', l_Jperio
[9169]	576	WRITE(numout,*) ' nlci = ', nlci
	577	WRITE(numout,*) ' nlcj = ', nlcj
	578	WRITE(numout,*) ' nimpp = ', nimpp
	579	WRITE(numout,*) ' njmpp = ', njmpp
	580	WRITE(numout,*) ' nreci = ', nreci
	581	WRITE(numout,*) ' nrecj = ', nrecj
	582	WRITE(numout,*) ' nn_hls = ', nn_hls
[3]	583	ENDIF
	584
[9019]	585	! ! Prepare mpp north fold
[6412]	586	IF( jperio >= 3 .AND. jperio <= 6 .AND. jpni > 1 ) THEN
[3]	587	CALL mpp_ini_north
[9917]	588	IF (lwp) THEN
	589	WRITE(numout,*)
	590	WRITE(numout,*) ' ==>>> North fold boundary prepared for jpni >1'
	591	! additional prints in layout.dat
	592	WRITE(inum,*)
	593	WRITE(inum,*)
	594	WRITE(inum,*) 'number of subdomains located along the north fold : ', ndim_rank_north
	595	WRITE(inum,*) 'Rank of the subdomains located along the north fold : ', ndim_rank_north
	596	DO jproc = 1, ndim_rank_north, 5
	597	WRITE(inum,*) nrank_north( jproc:MINVAL( (/jproc+4,ndim_rank_north/) ) )
	598	END DO
	599	ENDIF
[6412]	600	ENDIF
[9019]	601	!
	602	CALL mpp_init_ioipsl ! Prepare NetCDF output file (if necessary)
	603	!
[9917]	604	IF( ln_nnogather ) THEN
	605	CALL mpp_init_nfdcom ! northfold neighbour lists
	606	IF (lwp) THEN
	607	WRITE(inum,*)
	608	WRITE(inum,*)
	609	WRITE(inum,*) 'north fold exchanges with explicit point-to-point messaging :'
	610	WRITE(inum,*) 'nfsloop : ', nfsloop
	611	WRITE(inum,*) 'nfeloop : ', nfeloop
	612	WRITE(inum,*) 'nsndto : ', nsndto
	613	WRITE(inum,*) 'isendto : ', isendto
	614	ENDIF
	615	ENDIF
[9436]	616	!
[9917]	617	IF (lwp) CLOSE(inum)
	618	!
[9444]	619	DEALLOCATE(iin, ijn, ii_nono, ii_noea, ii_noso, ii_nowe, &
	620	& iimppt, ijmppt, ibondi, ibondj, ipproc, ipolj, &
	621	& ilci, ilcj, ilei, ilej, ildi, ildj, &
[10330]	622	& iono, ioea, ioso, iowe, llisoce)
[9444]	623	!
[9019]	624	END SUBROUTINE mpp_init
[3]	625
	626
[10330]	627	SUBROUTINE mpp_basic_decomposition( knbi, knbj, kimax, kjmax, kimppt, kjmppt, klci, klcj)
[9019]	628	!!----------------------------------------------------------------------
[10330]	629	!! * ROUTINE mpp_basic_decomposition *
	630	!!
	631	!! ** Purpose : Lay out the global domain over processors.
[9019]	632	!!
[10330]	633	!! ** Method : Global domain is distributed in smaller local domains.
[9019]	634	!!
[10330]	635	!! ** Action : - set for all knbi*knbj domains:
	636	!! kimppt : longitudinal index
	637	!! kjmppt : latitudinal index
	638	!! klci : first dimension
	639	!! klcj : second dimension
[9019]	640	!!----------------------------------------------------------------------
[10330]	641	INTEGER, INTENT(in ) :: knbi, knbj
	642	INTEGER, INTENT( out) :: kimax, kjmax
	643	INTEGER, DIMENSION(knbi,knbj), INTENT( out) :: kimppt, kjmppt
	644	INTEGER, DIMENSION(knbi,knbj), INTENT( out) :: klci, klcj
	645	!
	646	INTEGER :: ji, jj
	647	INTEGER :: iresti, irestj
	648	INTEGER :: ireci, irecj
[9019]	649	!!----------------------------------------------------------------------
	650	!
[10330]	651	#if defined key_nemocice_decomp
	652	jpimax = ( nx_global+2-2nn_hls + (knbi-1) ) / knbi + 2nn_hls ! first dim.
	653	jpjmax = ( ny_global+2-2nn_hls + (knbj-1) ) / knbj + 2nn_hls ! second dim.
	654	#else
	655	kimax = ( jpiglo - 2nn_hls + (knbi-1) ) / knbi + 2nn_hls ! first dim.
	656	kjmax = ( jpjglo - 2nn_hls + (knbj-1) ) / knbj + 2nn_hls ! second dim.
	657	#endif
[9019]	658	!
[10330]	659	! 1. Dimension arrays for subdomains
	660	! -----------------------------------
	661	! Computation of local domain sizes klci() klcj()
	662	! These dimensions depend on global sizes knbi,knbj and jpiglo,jpjglo
	663	! The subdomains are squares lesser than or equal to the global
	664	! dimensions divided by the number of processors minus the overlap array.
[9019]	665	!
[10330]	666	ireci = 2 * nn_hls
	667	irecj = 2 * nn_hls
	668	iresti = 1 + MOD( jpiglo - ireci -1 , knbi )
	669	irestj = 1 + MOD( jpjglo - irecj -1 , knbj )
[9168]	670	!
[10330]	671	! Need to use kimax and kjmax here since jpi and jpj not yet defined
	672	#if defined key_nemocice_decomp
	673	! Change padding to be consistent with CICE
	674	klci(1:knbi-1 ,:) = kimax
	675	klci(knbi ,:) = jpiglo - (knbi - 1) * (kimax - nreci)
	676	klcj(:, 1:knbj-1) = kjmax
	677	klcj(:, knbj) = jpjglo - (knbj - 1) * (kjmax - nrecj)
	678	#else
	679	klci(1:iresti ,:) = kimax
	680	klci(iresti+1:knbi ,:) = kimax-1
	681	klcj(:, 1:irestj) = kjmax
	682	klcj(:, irestj+1:knbj) = kjmax-1
	683	#endif
[3]	684
[10330]	685	! 2. Index arrays for subdomains
	686	! -------------------------------
	687	kimppt(:,:) = 1
	688	kjmppt(:,:) = 1
	689	!
	690	IF( knbi > 1 ) THEN
	691	DO jj = 1, knbj
	692	DO ji = 2, knbi
	693	kimppt(ji,jj) = kimppt(ji-1,jj) + klci(ji-1,jj) - ireci
	694	END DO
	695	END DO
[9019]	696	ENDIF
	697	!
[10330]	698	IF( knbj > 1 )THEN
	699	DO jj = 2, knbj
	700	DO ji = 1, knbi
	701	kjmppt(ji,jj) = kjmppt(ji,jj-1) + klcj(ji,jj-1) - irecj
	702	END DO
	703	END DO
	704	ENDIF
	705
	706	END SUBROUTINE mpp_basic_decomposition
[9019]	707
	708
[10330]	709	SUBROUTINE mpp_init_bestpartition( knbij, knbi, knbj, ldlist )
	710	!!----------------------------------------------------------------------
	711	!! * ROUTINE mpp_init_bestpartition *
	712	!!
	713	!! ** Purpose :
	714	!!
	715	!! ** Method :
	716	!!----------------------------------------------------------------------
	717	INTEGER, INTENT(in ) :: knbij ! total number if subdomains (knbi*knbj)
	718	INTEGER, OPTIONAL, INTENT( out) :: knbi, knbj ! number if subdomains along i and j (knbi and knbj)
	719	LOGICAL, OPTIONAL, INTENT(in ) :: ldlist ! .true.: print the list the best domain decompositions (with land)
	720	!
	721	INTEGER :: ji, jj, ii, iitarget
	722	INTEGER :: iszitst, iszjtst
	723	INTEGER :: isziref, iszjref
	724	INTEGER :: inbij, iszij
	725	INTEGER :: inbimax, inbjmax, inbijmax
	726	INTEGER :: isz0, isz1
	727	INTEGER, DIMENSION( :), ALLOCATABLE :: indexok
	728	INTEGER, DIMENSION( :), ALLOCATABLE :: inbi0, inbj0, inbij0 ! number of subdomains along i,j
	729	INTEGER, DIMENSION( :), ALLOCATABLE :: iszi0, iszj0, iszij0 ! max size of the subdomains along i,j
	730	INTEGER, DIMENSION( :), ALLOCATABLE :: inbi1, inbj1, inbij1 ! number of subdomains along i,j
	731	INTEGER, DIMENSION( :), ALLOCATABLE :: iszi1, iszj1, iszij1 ! max size of the subdomains along i,j
	732	LOGICAL :: llist
	733	LOGICAL, DIMENSION(:,:), ALLOCATABLE :: llmsk2d ! max size of the subdomains along i,j
	734	LOGICAL, DIMENSION(:,:), ALLOCATABLE :: llisoce ! - -
	735	REAL(wp):: zpropland
	736	!!----------------------------------------------------------------------
	737	!
	738	llist = .FALSE.
	739	IF( PRESENT(ldlist) ) llist = ldlist
	740
	741	CALL mpp_init_landprop( zpropland ) ! get the proportion of land point over the gloal domain
	742	inbij = NINT( REAL(knbij, wp) / ( 1.0 - zpropland ) ) ! define the largest possible value for jpni*jpnj
	743	!
	744	IF( llist ) THEN ; inbijmax = inbij*2
	745	ELSE ; inbijmax = inbij
	746	ENDIF
	747	!
	748	ALLOCATE(inbi0(inbijmax),inbj0(inbijmax),iszi0(inbijmax),iszj0(inbijmax))
	749	!
	750	inbimax = 0
	751	inbjmax = 0
	752	isziref = jpiglo*jpjglo+1
	753	iszjref = jpiglo*jpjglo+1
	754	!
	755	! get the list of knbi that gives a smaller jpimax than knbi-1
	756	! get the list of knbj that gives a smaller jpjmax than knbj-1
	757	DO ji = 1, inbijmax
	758	#if defined key_nemocice_decomp
	759	iszitst = ( nx_global+2-2nn_hls + (ji-1) ) / ji + 2nn_hls ! first dim.
	760	#else
	761	iszitst = ( jpiglo - 2nn_hls + (ji-1) ) / ji + 2nn_hls
	762	#endif
	763	IF( iszitst < isziref ) THEN
	764	isziref = iszitst
	765	inbimax = inbimax + 1
	766	inbi0(inbimax) = ji
	767	iszi0(inbimax) = isziref
	768	ENDIF
	769	#if defined key_nemocice_decomp
	770	iszjtst = ( ny_global+2-2nn_hls + (ji-1) ) / ji + 2nn_hls ! first dim.
	771	#else
	772	iszjtst = ( jpjglo - 2nn_hls + (ji-1) ) / ji + 2nn_hls
	773	#endif
	774	IF( iszjtst < iszjref ) THEN
	775	iszjref = iszjtst
	776	inbjmax = inbjmax + 1
	777	inbj0(inbjmax) = ji
	778	iszj0(inbjmax) = iszjref
	779	ENDIF
	780	END DO
	781
	782	! combine these 2 lists to get all possible knbi*knbj < inbijmax
	783	ALLOCATE( llmsk2d(inbimax,inbjmax) )
	784	DO jj = 1, inbjmax
	785	DO ji = 1, inbimax
	786	IF ( inbi0(ji) * inbj0(jj) <= inbijmax ) THEN ; llmsk2d(ji,jj) = .TRUE.
	787	ELSE ; llmsk2d(ji,jj) = .FALSE.
	788	ENDIF
	789	END DO
	790	END DO
	791	isz1 = COUNT(llmsk2d)
	792	ALLOCATE( inbi1(isz1), inbj1(isz1), iszi1(isz1), iszj1(isz1) )
	793	ii = 0
	794	DO jj = 1, inbjmax
	795	DO ji = 1, inbimax
	796	IF( llmsk2d(ji,jj) == .TRUE. ) THEN
	797	ii = ii + 1
	798	inbi1(ii) = inbi0(ji)
	799	inbj1(ii) = inbj0(jj)
	800	iszi1(ii) = iszi0(ji)
	801	iszj1(ii) = iszj0(jj)
	802	END IF
	803	END DO
	804	END DO
	805	DEALLOCATE( inbi0, inbj0, iszi0, iszj0 )
	806	DEALLOCATE( llmsk2d )
	807
	808	ALLOCATE( inbij1(isz1), iszij1(isz1) )
	809	inbij1(:) = inbi1(:) * inbj1(:)
	810	iszij1(:) = iszi1(:) * iszj1(:)
	811
	812	! if therr is no land and no print
	813	IF( .NOT. llist .AND. numbot == -1 .AND. numbdy == -1 ) THEN
	814	! get the smaller partition which gives the smallest subdomain size
	815	ii = MINLOC(inbij1, mask = iszij1 == MINVAL(iszij1), dim = 1)
	816	knbi = inbi1(ii)
	817	knbj = inbj1(ii)
	818	DEALLOCATE( inbi1, inbj1, inbij1, iszi1, iszj1, iszij1 )
	819	RETURN
	820	ENDIF
	821
	822	! extract only the partitions which reduce the subdomain size in comparison with smaller partitions
	823	ALLOCATE( indexok(isz1) ) ! to store indices of the best partitions
	824	isz0 = 0 ! number of best partitions
	825	inbij = 1 ! start with the min value of inbij1 => 1
	826	iszij = jpiglo*jpjglo+1 ! default: larger than global domain
	827	DO WHILE( inbij <= inbijmax ) ! if we did not reach the max of inbij1
	828	ii = MINLOC(iszij1, mask = inbij1 == inbij, dim = 1) ! warning: send back the first occurence if multiple results
	829	IF ( iszij1(ii) < iszij ) THEN
	830	isz0 = isz0 + 1
	831	indexok(isz0) = ii
	832	iszij = iszij1(ii)
	833	ENDIF
	834	inbij = MINVAL(inbij1, mask = inbij1 > inbij) ! warning: return largest integer value if mask = .false. everywhere
	835	END DO
	836	DEALLOCATE( inbij1, iszij1 )
	837
	838	! keep only the best partitions (sorted by increasing order of subdomains number and decreassing subdomain size)
	839	ALLOCATE( inbi0(isz0), inbj0(isz0), iszi0(isz0), iszj0(isz0) )
	840	DO ji = 1, isz0
	841	ii = indexok(ji)
	842	inbi0(ji) = inbi1(ii)
	843	inbj0(ji) = inbj1(ii)
	844	iszi0(ji) = iszi1(ii)
	845	iszj0(ji) = iszj1(ii)
	846	END DO
	847	DEALLOCATE( indexok, inbi1, inbj1, iszi1, iszj1 )
	848
	849	IF( llist ) THEN ! we print about 21 best partitions
	850	IF(lwp) THEN
	851	WRITE(numout,*)
	852	WRITE(numout, *) ' For your information:'
	853	WRITE(numout,'(a,i5,a)') ' list of the best partitions around ', knbij, ' mpi processes'
	854	WRITE(numout, *) ' --------------------------------------', '-----', '--------------'
	855	WRITE(numout,*)
	856	END IF
	857	iitarget = MINLOC( inbi0(:)inbj0(:), mask = inbi0(:)inbj0(:) >= knbij, dim = 1 )
	858	DO ji = MAX(1,iitarget-10), MIN(isz0,iitarget+10)
	859	ALLOCATE( llisoce(inbi0(ji), inbj0(ji)) )
	860	CALL mpp_init_isoce( inbi0(ji), inbj0(ji), llisoce ) ! Warning: must be call by all cores (call mpp_sum)
	861	inbij = COUNT(llisoce)
	862	DEALLOCATE( llisoce )
	863	IF(lwp) WRITE(numout,'(a, i5, a, i5, a, i4, a, i4, a, i9, a, i5, a, i5, a)') &
	864	& 'nb_cores ' , inbij,' oce + ', inbi0(ji)*inbj0(ji) - inbij &
	865	& , ' land ( ', inbi0(ji),' x ', inbj0(ji), &
	866	& ' ), nb_points ', iszi0(ji)*iszj0(ji),' ( ', iszi0(ji),' x ', iszj0(ji),' )'
	867	END DO
	868	DEALLOCATE( inbi0, inbj0, iszi0, iszj0 )
	869	RETURN
	870	ENDIF
	871
	872	DEALLOCATE( iszi0, iszj0 )
	873	inbij = inbijmax + 1 ! default: larger than possible
	874	ii = isz0+1 ! start from the end of the list (smaller subdomains)
	875	DO WHILE( inbij > knbij ) ! while the number of ocean subdomains exceed the number of procs
	876	ii = ii -1
	877	ALLOCATE( llisoce(inbi0(ii), inbj0(ii)) )
	878	CALL mpp_init_isoce( inbi0(ii), inbj0(ii), llisoce ) ! must be done by all core
	879	inbij = COUNT(llisoce)
	880	DEALLOCATE( llisoce )
	881	END DO
	882	knbi = inbi0(ii)
	883	knbj = inbj0(ii)
	884	DEALLOCATE( inbi0, inbj0 )
	885	!
	886	END SUBROUTINE mpp_init_bestpartition
	887
	888
	889	SUBROUTINE mpp_init_landprop( propland )
	890	!!----------------------------------------------------------------------
	891	!! * ROUTINE mpp_init_landprop *
	892	!!
	893	!! ** Purpose : the the proportion of land points in the surface land-sea mask
	894	!!
	895	!! ** Method : read iproc strips (of length jpiglo) of the land-sea mask
	896	!!----------------------------------------------------------------------
	897	REAL(wp), INTENT( out) :: propland ! proportion of land points (between 0 and 1)
	898	!
	899	INTEGER, DIMENSION(jpni*jpnj) :: kusedom_1d
	900	INTEGER :: inboce
	901	INTEGER :: iproc, idiv, ijsz
	902	INTEGER :: ijstr, ijend, ijcnt
	903	LOGICAL, ALLOCATABLE, DIMENSION(:,:) :: lloce
	904	!!----------------------------------------------------------------------
	905	! do nothing if there is no land-sea mask
	906	IF( numbot == -1 .and. numbdy == -1 ) THEN
	907	propland = 0.
	908	RETURN
	909	ENDIF
	910
	911	! number of processes reading the bathymetry file
	912	iproc = MINVAL( (/mppsize, jpjglo/2, 100/) ) ! read a least 2 lines, no more that 100 processes reading at the same time
	913
	914	! we want to read iproc strips of the land-sea mask. -> pick up iproc processes among mppsize processes
	915	IF( iproc == 1 ) THEN ; idiv = mppsize
	916	ELSE ; idiv = ( mppsize - 1 ) / ( iproc - 1 )
	917	ENDIF
	918	ijsz = jpjglo / iproc
	919	IF( narea <= MOD(jpjglo,iproc) ) ijsz = ijsz + 1
	920
	921	IF( MOD( narea-1, idiv ) == 0 .AND. (idiv /= 1 .OR. narea <= iproc ) ) THEN
	922	!
	923	ijstr = (narea-1)*(jpjglo/iproc) + MIN(narea-1, MOD(jpjglo,iproc)) + 1
	924	ijend = ijstr + ijsz - 1
	925	ijcnt = ijend - ijstr + 1
	926	!
	927	ALLOCATE( lloce(jpiglo, ijcnt) ) ! allocate the strip
	928	CALL mpp_init_readbot_strip( ijstr, ijcnt, lloce )
	929	inboce = COUNT(lloce)
	930	DEALLOCATE(lloce)
	931	!
	932	ELSE
	933	inboce = 0
	934	ENDIF
	935	CALL mpp_sum( 'mppini', inboce )
	936	!
	937	propland = REAL( jpiglojpjglo - inboce, wp ) / REAL( jpiglojpjglo, wp )
	938	!
	939	END SUBROUTINE mpp_init_landprop
	940
	941
	942	SUBROUTINE mpp_init_isoce( knbi, knbj, ldisoce )
	943	!!----------------------------------------------------------------------
	944	!! * ROUTINE mpp_init_nboce *
	945	!!
	946	!! ** Purpose : check for a mpi domain decomposition knbi x knbj which
	947	!! subdomains contain at least 1 ocean point
	948	!!
	949	!! ** Method : read knbj strips (of length jpiglo) of the land-sea mask
	950	!!----------------------------------------------------------------------
	951	INTEGER, INTENT(in ) :: knbi, knbj
	952	LOGICAL, DIMENSION(knbi,knbj), INTENT( out) :: ldisoce !
	953	!
	954	INTEGER, DIMENSION(knbi,knbj) :: inboce
	955	INTEGER, DIMENSION(knbi*knbj) :: inboce_1d
	956	INTEGER :: idiv, i2read, inj
	957	INTEGER :: iimax, ijmax
	958	INTEGER :: ji,jj
	959	LOGICAL, ALLOCATABLE, DIMENSION(:,:) :: lloce
	960	INTEGER, ALLOCATABLE, DIMENSION(:,:) :: iimppt, ilci
	961	INTEGER, ALLOCATABLE, DIMENSION(:,:) :: ijmppt, ilcj
	962	!!----------------------------------------------------------------------
	963	! do nothing if there is no land-sea mask
	964	IF( numbot == -1 .AND. numbdy == -1 ) THEN
	965	ldisoce(:,:) = .TRUE.
	966	RETURN
	967	ENDIF
	968
	969	! we want to read knbj strips of the land-sea mask. -> pick up knbj processes among mppsize processes
	970	IF( knbj == 1 ) THEN ; idiv = mppsize
	971	ELSE ; idiv = ( mppsize - 1 ) / ( knbj - 1 )
	972	ENDIF
	973	inboce(:,:) = 0
	974	IF( MOD( narea-1, idiv ) == 0 .AND. (idiv /= 1 .OR. narea <= knbj ) ) THEN
	975	!
	976	ALLOCATE( iimppt(knbi,knbj), ijmppt(knbi,knbj), ilci(knbi,knbj), ilcj(knbi,knbj) )
	977	CALL mpp_basic_decomposition( knbi, knbj, iimax, ijmax, iimppt, ijmppt, ilci, ilcj )
	978	!
	979	i2read = knbj / mppsize ! strip number to be read by this process
	980	IF( ( narea - 1 ) / idiv < MOD(knbj,mppsize) ) i2read = i2read + 1
	981	DO jj = 1, i2read
	982	! strip number to be read (from 1 to knbj)
	983	inj = ( narea - 1 ) * ( knbj / mppsize ) + MIN( MOD(knbj,mppsize), ( narea - 1 ) / idiv ) + jj
	984	ALLOCATE( lloce(jpiglo, ilcj(1,inj)) ) ! allocate the strip
	985	CALL mpp_init_readbot_strip( ijmppt(1,inj), ilcj(1,inj), lloce ) ! read the strip
	986	DO ji = 1, knbi
	987	inboce(ji,inj) = COUNT( lloce(iimppt(ji,1):iimppt(ji,1)+ilci(ji,1)-1,:) )
	988	END DO
	989	DEALLOCATE(lloce)
	990	END DO
	991	!
	992	DEALLOCATE(iimppt, ijmppt, ilci, ilcj)
	993	!
	994	ENDIF
	995
	996	inboce_1d = RESHAPE(inboce, (/ knbi*knbj /))
	997	CALL mpp_sum( 'mppini', inboce_1d )
	998	inboce = RESHAPE(inboce_1d, (/knbi, knbj/))
	999	ldisoce = inboce /= 0
	1000	!
	1001	END SUBROUTINE mpp_init_isoce
	1002
	1003
	1004	SUBROUTINE mpp_init_readbot_strip( kjstr, kjcnt, ldoce )
	1005	!!----------------------------------------------------------------------
	1006	!! * ROUTINE mpp_init_readbot_strip *
	1007	!!
	1008	!! ** Purpose : Read relevant bathymetric information in order to
	1009	!! provide a land/sea mask used for the elimination
	1010	!! of land domains, in an mpp computation.
	1011	!!
	1012	!! ** Method : read stipe of size (jpiglo,...)
	1013	!!----------------------------------------------------------------------
	1014	INTEGER , INTENT(in ) :: kjstr !
	1015	INTEGER , INTENT(in ) :: kjcnt !
	1016	LOGICAL, DIMENSION(jpiglo,kjcnt), INTENT( out) :: ldoce !
	1017	!
	1018	INTEGER :: inumsave ! local logical unit
	1019	REAL(wp), DIMENSION(jpiglo,kjcnt) :: zbot, zbdy
	1020	!!----------------------------------------------------------------------
	1021	!
	1022	inumsave = numout ; numout = numnul ! redirect all print to /dev/null
	1023	!
	1024	IF( numbot /= -1 ) THEN
	1025	CALL iom_get( numbot, jpdom_unknown, 'bottom_level', zbot, kstart = (/1,kjstr/), kcount = (/jpiglo, kjcnt/) )
	1026	ELSE
	1027	zbot(:,:) = 1. ! put a non-null value
	1028	ENDIF
	1029
	1030	IF( numbdy /= -1 ) THEN ! Adjust with bdy_msk if it exists
	1031	CALL iom_get ( numbdy, jpdom_unknown, 'bdy_msk', zbdy, kstart = (/1,kjstr/), kcount = (/jpiglo, kjcnt/) )
	1032	zbot(:,:) = zbot(:,:) * zbdy(:,:)
	1033	ENDIF
	1034	!
	1035	ldoce = zbot > 0.
	1036	numout = inumsave
	1037	!
	1038	END SUBROUTINE mpp_init_readbot_strip
	1039
	1040
[88]	1041	SUBROUTINE mpp_init_ioipsl
	1042	!!----------------------------------------------------------------------
	1043	!! * ROUTINE mpp_init_ioipsl *
	1044	!!
	1045	!! ** Purpose :
	1046	!!
	1047	!! ** Method :
	1048	!!
	1049	!! History :
[1238]	1050	!! 9.0 ! 04-03 (G. Madec ) MPP-IOIPSL
	1051	!! " " ! 08-12 (A. Coward) addition in case of jpni*jpnj < jpnij
[88]	1052	!!----------------------------------------------------------------------
[2715]	1053	INTEGER, DIMENSION(2) :: iglo, iloc, iabsf, iabsl, ihals, ihale, idid
[88]	1054	!!----------------------------------------------------------------------
[352]	1055
[1238]	1056	! The domain is split only horizontally along i- or/and j- direction
	1057	! So we need at the most only 1D arrays with 2 elements.
	1058	! Set idompar values equivalent to the jpdom_local_noextra definition
	1059	! used in IOM. This works even if jpnij .ne. jpni*jpnj.
[88]	1060	iglo(1) = jpiglo
	1061	iglo(2) = jpjglo
	1062	iloc(1) = nlci
	1063	iloc(2) = nlcj
	1064	iabsf(1) = nimppt(narea)
	1065	iabsf(2) = njmppt(narea)
	1066	iabsl(:) = iabsf(:) + iloc(:) - 1
[1238]	1067	ihals(1) = nldi - 1
	1068	ihals(2) = nldj - 1
	1069	ihale(1) = nlci - nlei
	1070	ihale(2) = nlcj - nlej
[352]	1071	idid(1) = 1
	1072	idid(2) = 2
	1073
[88]	1074	IF(lwp) THEN
[516]	1075	WRITE(numout,*)
[352]	1076	WRITE(numout,*) 'mpp_init_ioipsl : iloc = ', iloc (1), iloc (2)
	1077	WRITE(numout,*) '~~~~~~~~~~~~~~~ iabsf = ', iabsf(1), iabsf(2)
	1078	WRITE(numout,*) ' ihals = ', ihals(1), ihals(2)
	1079	WRITE(numout,*) ' ihale = ', ihale(1), ihale(2)
[88]	1080	ENDIF
[2715]	1081	!
[352]	1082	CALL flio_dom_set ( jpnij, nproc, idid, iglo, iloc, iabsf, iabsl, ihals, ihale, 'BOX', nidom)
[2715]	1083	!
[88]	1084	END SUBROUTINE mpp_init_ioipsl
	1085
[9436]	1086
	1087	SUBROUTINE mpp_init_nfdcom
	1088	!!----------------------------------------------------------------------
	1089	!! * ROUTINE mpp_init_nfdcom *
	1090	!! ** Purpose : Setup for north fold exchanges with explicit
	1091	!! point-to-point messaging
	1092	!!
	1093	!! ** Method : Initialization of the northern neighbours lists.
	1094	!!----------------------------------------------------------------------
	1095	!! 1.0 ! 2011-10 (A. C. Coward, NOCS & J. Donners, PRACE)
	1096	!! 2.0 ! 2013-06 Setup avoiding MPI communication (I. Epicoco, S. Mocavero, CMCC)
	1097	!!----------------------------------------------------------------------
	1098	INTEGER :: sxM, dxM, sxT, dxT, jn
	1099	INTEGER :: njmppmax
	1100	!!----------------------------------------------------------------------
	1101	!
	1102	njmppmax = MAXVAL( njmppt )
	1103	!
	1104	!initializes the north-fold communication variables
	1105	isendto(:) = 0
	1106	nsndto = 0
	1107	!
	1108	IF ( njmpp == njmppmax ) THEN ! if I am a process in the north
	1109	!
	1110	!sxM is the first point (in the global domain) needed to compute the north-fold for the current process
	1111	sxM = jpiglo - nimppt(narea) - nlcit(narea) + 1
	1112	!dxM is the last point (in the global domain) needed to compute the north-fold for the current process
	1113	dxM = jpiglo - nimppt(narea) + 2
	1114	!
	1115	! loop over the other north-fold processes to find the processes
	1116	! managing the points belonging to the sxT-dxT range
	1117	!
	1118	DO jn = 1, jpni
	1119	!
	1120	sxT = nfiimpp(jn, jpnj) ! sxT = 1st point (in the global domain) of the jn process
	1121	dxT = nfiimpp(jn, jpnj) + nfilcit(jn, jpnj) - 1 ! dxT = last point (in the global domain) of the jn process
	1122	!
	1123	IF ( sxT < sxM .AND. sxM < dxT ) THEN
	1124	nsndto = nsndto + 1
	1125	isendto(nsndto) = jn
	1126	ELSEIF( sxM <= sxT .AND. dxM >= dxT ) THEN
	1127	nsndto = nsndto + 1
	1128	isendto(nsndto) = jn
	1129	ELSEIF( dxM < dxT .AND. sxT < dxM ) THEN
	1130	nsndto = nsndto + 1
	1131	isendto(nsndto) = jn
	1132	ENDIF
	1133	!
	1134	END DO
	1135	nfsloop = 1
	1136	nfeloop = nlci
	1137	DO jn = 2,jpni-1
	1138	IF( nfipproc(jn,jpnj) == (narea - 1) ) THEN
	1139	IF( nfipproc(jn-1,jpnj) == -1 ) nfsloop = nldi
	1140	IF( nfipproc(jn+1,jpnj) == -1 ) nfeloop = nlei
	1141	ENDIF
	1142	END DO
	1143	!
	1144	ENDIF
	1145	l_north_nogather = .TRUE.
	1146	!
	1147	END SUBROUTINE mpp_init_nfdcom
	1148
[10330]	1149
[3]	1150	#endif
[88]	1151
[3]	1152	!!======================================================================
	1153	END MODULE mppini

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

Download in other formats: