74 ZOLTAN_ID_PTR globalIDs,
75 ZOLTAN_ID_PTR localIDs,
106 globalIDs[i] = globalMap.
toGlobal(i);
108 for(
int j=0; j<wgt_dim; j++)
110 obj_wgts[wgt_dim*i + j] = weights[wgt_dim*i + j];
130 ZOLTAN_ID_PTR globalIDs,
131 ZOLTAN_ID_PTR localIDs,
148 *ierr = ZOLTAN_FATAL;
174 ZOLTAN_ID_PTR globalIDs,
175 ZOLTAN_ID_PTR localIDs,
185 if ((nGID != 1) || (nLID != 1) || (
nPoints != sizes.
size()))
187 *ierr = ZOLTAN_FATAL;
193 numEdges[i] = sizes[i];
206 ZOLTAN_ID_PTR globalIDs,
207 ZOLTAN_ID_PTR localIDs,
209 ZOLTAN_ID_PTR nborGID,
229 *ierr = ZOLTAN_FATAL;
245 const CompactListList<label>& adjacency,
248 List<label>& finalDecomp
254 <<
" Zoltan cannot redistribute without points"
255 "on all processors." <<
endl;
260 args[0] =
"zoltanDecomp";
264 for (
label i = 0; i < argc; i++)
266 argv[i] = strdup(
args[i].c_str());
270 int rc = Zoltan_Initialize(argc, argv, &ver);
280 const int nWeights = pWeights.size()/
points.size();
283 Zoltan_Set_Param(zz,
"return_lists",
"export");
284 Zoltan_Set_Param(zz,
"obj_weight_dim",
name(nWeights).c_str());
285 Zoltan_Set_Param(zz,
"edge_weight_dim",
"0");
288 Zoltan_Set_Param(zz,
"debug_level",
"0");
289 Zoltan_Set_Param(zz,
"imbalance_tol",
"1.05");
292 Zoltan_Set_Param(zz,
"lb_method",
"graph");
293 Zoltan_Set_Param(zz,
"lb_approach",
"repartition");
297 const dictionary& coeffsDict_ =
302 if (!iter().isDict())
304 const word& key = iter().keyword();
305 const word value(iter().stream());
309 Info<< typeName <<
" : setting parameter " << key
310 <<
" to " << value <<
endl;
313 Zoltan_Set_Param(zz, key.c_str(), value.c_str());
321 Tuple2<const pointField&, const scalarField&> vertexData(
points, pWeights);
329 void* adjacencyPtr = &
const_cast<CompactListList<label>&
>(adjacency);
331 Zoltan_Set_Edge_List_Multi_Fn(zz,
get_edge_list, adjacencyPtr);
333 int changed, num_imp, num_exp, *imp_procs, *exp_procs;
334 int *imp_to_part, *exp_to_part;
335 int num_gid_entries, num_lid_entries;
336 ZOLTAN_ID_PTR imp_global_ids, exp_global_ids;
337 ZOLTAN_ID_PTR imp_local_ids, exp_local_ids;
342 int oldExcepts = fedisableexcept
371 feenableexcept(oldExcepts);
374 if (debug && changed)
376 Pout <<
"Zoltan number to move " << changed <<
" " << num_exp <<
endl;
379 finalDecomp.setSize(
points.size());
384 for(
int i=0; i<num_exp; i++)
386 finalDecomp[exp_local_ids[i]] = exp_procs[i];
434 <<
"Can use this decomposition method only for the whole mesh"
436 <<
"and supply one coordinate (cellCentre) for every cell." <<
endl
437 <<
"The number of coordinates " << cellCentres.
size() <<
endl
438 <<
"The number of cells in the mesh " << mesh.
nCells()
468 decomp[i] = finalDecomp[i];
485 <<
"Size of cell-to-coarse map " << agglom.
size()
486 <<
" differs from number of cells in mesh " << mesh.
nCells()
514 forAll(fineDistribution, i)
516 fineDistribution[i] = finalDecomp[agglom[i]];
519 return fineDistribution;
533 <<
"Inconsistent number of cells (" << globalPointPoints.
size()
534 <<
") and number of cell centres (" <<
points.
size()
#define forAll(list, i)
Loop across all elements in list.
#define forAllConstIter(Container, container, iter)
Iterate across all elements in the container object of type.
Macros for easy insertion into run-time selection tables.
A packed storage unstructured matrix of objects of type <T> using an offset table for access.
void size(const label)
Override size to be inconsistent with allocated storage.
A 2-tuple for storing two objects of different types.
const Type2 & second() const
Return second.
const Type1 & first() const
Return first.
const UList< T > & m() const
Return the packed matrix of data.
label size() const
Return the primary size, i.e. the number of rows.
labelList sizes() const
Return sizes (to be used e.g. for construction)
static bool & parRun()
Is this a parallel run?
static int myProcNo(const label communicator=0)
Number of this process (starting from masterNo() = 0)
static const direction nComponents
Number of components in this vector space.
label size() const
Return the number of arguments.
Abstract base class for decomposition.
dictionary decompositionDict_
A list of keyword definitions, which are a keyword followed by any number of values (e....
const dictionary & subDict(const word &) const
Find and return a sub-dictionary.
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
Calculates a unique integer (label so might not have enough room - 2G max) for processor + local inde...
label whichProcID(const label i) const
Which processor does global come from? Binary search.
label toGlobal(const label i) const
From local to global.
Mesh consisting of general polyhedral cells.
Zoltan redistribution in parallel.
virtual labelList decompose(const pointField &points, const scalarField &pointWeights)
Inherit decompose from decompositionMethod.
zoltanDecomp(const dictionary &decompositionDict)
Construct given the decomposition dictionary and mesh.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
volScalarField scalarField(fieldObject, mesh)
errorManipArg< error, int > exit(error &err, const int errNo=1)
pointField vertices(const blockVertexList &bvl)
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Ostream & endl(Ostream &os)
Add newline and flush stream.
vectorField pointField
pointField is a vectorField.
defineTypeNameAndDebug(combustionModel, 0)
prefixOSstream Pout(cout, "Pout")
List< string > stringList
A List of strings.
labelList identityMap(const label len)
Create identity map (map[i] == i) of given length.
word name(const complex &)
Return a string representation of a complex.
addToRunTimeSelectionTable(ensightPart, ensightPartCells, istream)
Foam::argList args(argc, argv)
static int get_number_of_vertices(void *data, int *ierr)
static void get_num_edges_list(void *data, int nGID, int nLID, int nPoints, ZOLTAN_ID_PTR globalIDs, ZOLTAN_ID_PTR localIDs, int *numEdges, int *ierr)
static void get_geom_list(void *data, int nGID, int nLID, int nPoints, ZOLTAN_ID_PTR globalIDs, ZOLTAN_ID_PTR localIDs, int nDim, double *vertices, int *ierr)
static int get_mesh_dim(void *data, int *ierr)
static void get_vertex_list(void *data, int nGID, int nLID, ZOLTAN_ID_PTR globalIDs, ZOLTAN_ID_PTR localIDs, int wgt_dim, float *obj_wgts, int *ierr)
static void get_edge_list(void *data, int nGID, int nLID, int nPoints, ZOLTAN_ID_PTR globalIDs, ZOLTAN_ID_PTR localIDs, int *num_edges, ZOLTAN_ID_PTR nborGID, int *nborProc, int wgt_dim, float *ewgts, int *ierr)