diff --git a/Partition.py b/Partition.py
index b7f24afd04da2aa97e5e78e76439b6df1ffcd45f..f33e9752851f9d84ef620350b6729d2200cd0c95 100644
--- a/Partition.py
+++ b/Partition.py
@@ -62,8 +62,6 @@ def adjustpart(partin, land, nbcore) :
addnbland(partin, land)
nbptproc=np.array([dom['nbi']*dom['nbj'] for dom in partin])
- INFO("Nb of land points per proc : Mean = "+str(np.mean(nbptproc))+" Min = "+str(np.min(nbptproc))+\
- " Max = "+str(np.max(nbptproc))+" Var :"+str((np.max(nbptproc)-np.min(nbptproc))/np.mean(nbptproc)*100)+"%")
if np.min(nbptproc) < 3 :
ERROR("The smallest domain in the partition is of a size less than 5. This will probably not work.")
sys.exit()
@@ -71,7 +69,8 @@ def adjustpart(partin, land, nbcore) :
#
# Partition domain in two steps : 1) halving, 2) adjusting by nb land points
#
-def partitiondom(nig, njg, land, nbcore) :
+def partitiondom(land, nbcore, rank) :
+ njg,nig = land.shape
partout=[{"nbi":nig,"nbj":njg,"istart":0,"jstart":0, "nbland":int(np.sum(land))}]
while len(partout) <= nbcore/2 :
partout = halfpartition(partout)
@@ -79,14 +78,96 @@ def partitiondom(nig, njg, land, nbcore) :
#
adjustpart(partout, land, nbcore)
#
+ gindland=np.zeros((njg,nig), dtype=np.int32)
+ gindland[:,:]=-1
+ n=0
+ for i in range(nig) :
+ for j in range(njg) :
+ if (land[j,i] > 0 ) :
+ gindland[j,i] = n
+ n += 1
+ #
+ for ic in range(len(partout)) :
+ istart = partout[ic]["istart"]
+ nbi = partout[ic]["nbi"]
+ jstart = partout[ic]["jstart"]
+ nbj = partout[ic]["nbj"]
+ locgind=np.copy(gindland[jstart:jstart+nbj,istart:istart+nbi])
+ loccore=np.copy(gindland[jstart:jstart+nbj,istart:istart+nbi])
+ loccore[loccore >= 0] = ic
+ partout[ic]["glandindex"] = locgind
+ partout[ic]["corenb"] = loccore
+ partout[ic]["nblocland"] = np.count_nonzero(locgind >= 0)
+ #
+ return partout
+#
+# Build global map of processors
+#
+def buildprocmap(nig, njg, part) :
procmap=np.zeros((njg,nig), dtype=np.int8)
procmap[:,:] = -1
- for i in range(len(partout)) :
- for ij in range(partout[i]["nbj"]) :
- for ii in range(partout[i]["nbi"]) :
- procmap[ij+partout[i]["jstart"],ii+partout[i]["istart"]] = i
-
- return partout, procmap
+ for i in range(len(part)) :
+ for ij in range(part[i]["nbj"]) :
+ for ii in range(part[i]["nbi"]) :
+ procmap[ij+part[i]["jstart"],ii+part[i]["istart"]] = i
+ #
+ return procmap
+#
+# Partion of only land points in the case where we do not need halos
+#
+def landpartition(land, nbcore, rank) :
+ nbland = np.count_nonzero(land > 0)
+ optland = int(nbland/nbcore)+(nbland%nbcore > 0)
+ print("RANK",rank,"Nbland = ", np.count_nonzero(land > 0), " nbcore = ", nbcore, " optland = ",optland)
+ #
+ distppc=np.zeros(nbcore)
+ distppc[:]=0
+ for i in range(nbcore-1) :
+ ptsleft = int(nbland-np.sum(distppc[0:i]))
+ crsleft = nbcore-i
+ if ptsleft > optland and ptsleft/2 > crsleft :
+ distppc[i] = optland
+ else :
+ print("On proc ", i," Still have ", ptsleft, " points to distribute and ", crsleft, " cores")
+ distppc[i] = int(ptsleft/crsleft)
+ #
+ distppc[-1]=nbland-np.sum(distppc[0:-1])
+ print("Distpcc :", distppc)
+ #
+ njg,nig=land.shape
+ gindland=np.zeros((njg,nig), dtype=np.int32)
+ gindland[:,:]=-1
+ glandcore=np.zeros((njg,nig), dtype=np.int32)
+ glandcore[:,:]=-1
+ n=0
+ for i in range(nig) :
+ for j in range(njg) :
+ if (land[j,i] > 0 ) :
+ gindland[j,i] = n
+ cc = np.argmax(distppc > 0)
+ glandcore[j,i] = cc
+ distppc[cc] = distppc[cc]-1
+ n += 1
+ #
+ # Decompose domain according to land point positions
+ #
+ partout=[]
+ for ic in range(nbcore) :
+ ff = np.nonzero(glandcore == ic)
+ # Region is at least 3x3
+ nbi=max(max(ff[1])-min(ff[1])+1, 3)
+ nbj=max(max(ff[0])-min(ff[0])+1, 3)
+ istart=min(min(ff[1]),nig-nbi)
+ jstart=min(min(ff[0]),njg-nbj)
+ #
+ loccore=np.copy(glandcore[jstart:jstart+nbj,istart:istart+nbi])
+ loccore[loccore != ic] = -1
+ locgind=np.copy(gindland[jstart:jstart+nbj,istart:istart+nbi])
+ locgind[loccore < 0] = -1
+ partout.append({"nbi":nbi, "nbj":nbj, "istart":istart, "jstart":jstart, "nbland":ff[0].shape[0], \
+ "corenb":loccore, "glandindex":locgind, "nblocland":np.count_nonzero(locgind >= 0)})
+ #
+ return partout
#
# Add halo to the partition
#
@@ -217,9 +298,6 @@ def toland_index(x,landmap) :
#
class partition :
def __init__ (self, nig, njg, land, mpicomm, nbcore, halosz, rank, wunit="None") :
- #
- part, procmap = partitiondom(nig, njg, land, nbcore)
- halosource_map, innersend_map = addhalo(nig, njg, part, procmap, halosz)
#
# Self varibales
#
@@ -231,12 +309,27 @@ class partition :
self.nig = nig
self.njg = njg
self.nblandg = np.count_nonzero(land > 0)
+ #
+ part = partitiondom(land, nbcore, rank)
+ procmap = buildprocmap(nig, njg, part)
+ # To be activated when halosz = 0 and we can deal with single point domain grids.
+ #
+ # part = landpartition(land, nbcore, rank)
+ # procmap = buildprocmap(nig, njg, part)
+ #
+ halosource_map, innersend_map = addhalo(nig, njg, part, procmap, halosz)
+ #
+ INFO(" Nb of land points on proc "+str(rank)+" = "+str(part[rank]["nblocland"])+" from "+str(self.nblandg)+\
+ " -- "+"%.2f"%(part[rank]["nblocland"]/self.nblandg*100)+"% of total.")
+ #
self.allistart = []
self.alljstart = []
for i in range(self.size) :
self.allistart.append(part[i]["istart"])
self.alljstart.append(part[i]["jstart"])
- # Local info
+ #
+ # Information of domain on local processor
+ #
self.ni = part[rank]["nbi"]
self.nj = part[rank]["nbj"]
self.istart = part[rank]["istart"]
diff --git a/RoutingPreProc.py b/RoutingPreProc.py
index 52cf319b3dbb9c85a185551c4ba2ca0f60a95389..281f099054e681c91cff7e62a075b8491ae5ec2c 100644
--- a/RoutingPreProc.py
+++ b/RoutingPreProc.py
@@ -27,7 +27,6 @@ import ModelGrid as MG
import Partition as PA
import RPPtools as RPP
import HydroGrid as HG
-import diagplots as DP
import Interface as IF
from Truncate import trunc as TR
#
diff --git a/WeightsOnly.py b/WeightsOnly.py
new file mode 100644
index 0000000000000000000000000000000000000000..a3c80495bf8cec4774a92cad708d99aede506121
--- /dev/null
+++ b/WeightsOnly.py
@@ -0,0 +1,62 @@
+#
+#
+#
+import numpy as np
+import os, sys
+from mpi4py import MPI
+import gc
+#
+import matplotlib as mpl
+mpl.use('Agg')
+from matplotlib.backends.backend_pdf import PdfPages
+import matplotlib.pyplot as plt
+import time
+#
+# Gert the information from the configuration file.
+#
+import configparser
+config=configparser.ConfigParser()
+config.read("run.def")
+nbasmax=config.getint("OverAll", "nbasmax")
+numop=config.getint("OverAll", "numop", fallback=100)
+OutGraphFile=config.get("OverAll","GraphFile")
+DumpHydroSuper=config.getboolean("OverAll","DumpHydroSuper",fallback=False)
+wfile=config.get("OverAll","WeightFile",fallback="Weights.nc")
+#
+import ModelGrid as MG
+import Partition as PA
+import RPPtools as RPP
+import HydroGrid as HG
+#
+# Logging in MPI : https://groups.google.com/forum/#!topic/mpi4py/SaNzc8bdj6U
+#
+import getargs
+log_master, log_world = getargs.getLogger()
+INFO, DEBUG, ERROR = log_master.info, log_master.debug, log_world.error
+INFO_ALL, DEBUG_ALL = log_world.info, log_world.debug
+#
+# Read full grid and partition the domain.
+#
+gg = MG.GlobalGrid()
+#
+comm = MPI.COMM_WORLD
+szhalo=0
+nbcore=comm.Get_size()
+rank=comm.Get_rank()
+#
+# Region of grid to be treated
+#
+part = PA.partition(gg.ni, gg.nj, gg.land, comm, nbcore, szhalo, rank)
+#
+modelgrid=MG.ModelGrid(part.ihstart+gg.igstart,part.nih,part.jhstart+gg.jgstart,part.njh)
+INFO("Longitude interval on proc "+str(rank)+" = "+str(modelgrid.box_land[0][0])+" : "+str(modelgrid.box_land[0][1]))
+INFO("Latitude interval on proc "+str(rank)+" = "+str(modelgrid.box_land[1][0])+" : "+str(modelgrid.box_land[1][1]))
+#
+hydrogrid=HG.HydroGrid(modelgrid.box_land)
+#
+# Computes weights of overlap of modelgrid and hydrogrid
+#
+w = RPP.compweights(wfile, part, modelgrid, hydrogrid)
+#
+#
+