45 #ifndef _ZOLTAN2_ALGZOLTAN_HPP_ 46 #define _ZOLTAN2_ALGZOLTAN_HPP_ 57 #include <zoltan_cpp.h> 75 template <
typename Adapter>
81 typedef typename Adapter::lno_t lno_t;
82 typedef typename Adapter::gno_t gno_t;
83 typedef typename Adapter::scalar_t scalar_t;
84 typedef typename Adapter::part_t part_t;
85 typedef typename Adapter::user_t user_t;
86 typedef typename Adapter::userCoord_t userCoord_t;
88 const RCP<const Environment> env;
89 const RCP<const Comm<int> > problemComm;
90 const RCP<const typename Adapter::base_adapter_t> adapter;
91 RCP<const Model<Adapter> > model;
96 void setMPIComm(
const RCP<
const Comm<int> > &problemComm__) {
97 # ifdef HAVE_ZOLTAN2_MPI 98 mpicomm = Teuchos::getRawMpiComm(*problemComm__);
100 mpicomm = MPI_COMM_WORLD;
109 Zoltan_Initialize(argc, argv, &ver);
112 void setCallbacksIDs()
114 zz->Set_Num_Obj_Fn(zoltanNumObj<Adapter>, (
void *) &(*adapter));
115 zz->Set_Obj_List_Fn(zoltanObjList<Adapter>, (
void *) &(*adapter));
117 const part_t *myparts;
118 adapter->getPartsView(myparts);
120 zz->Set_Part_Multi_Fn(zoltanParts<Adapter>, (
void *) &(*adapter));
124 zz->Set_Param(
"NUM_GID_ENTRIES", tmp);
126 zz->Set_Param(
"NUM_LID_ENTRIES", tmp);
129 template <
typename AdapterWithCoords>
130 void setCallbacksGeom(
const AdapterWithCoords *ia)
135 zz->Set_Num_Geom_Fn(zoltanNumGeom<AdapterWithCoords>, (
void *) ia);
136 zz->Set_Geom_Multi_Fn(zoltanGeom<AdapterWithCoords>, (
void *) ia);
139 void setCallbacksGraph(
146 void setCallbacksGraph(
153 void setCallbacksGraph(
160 void setCallbacksHypergraph(
166 zz->Set_HG_Size_CS_Fn(zoltanHGSizeCS_withMatrixAdapter<Adapter>,
168 zz->Set_HG_CS_Fn(zoltanHGCS_withMatrixAdapter<Adapter>,
180 const Teuchos::ParameterList &pl = env->getParameters();
182 const Teuchos::ParameterEntry *pe = pl.getEntryPtr(
"hypergraph_model_type");
183 std::string model_type(
"traditional");
185 model_type = pe->getValue<std::string>(&model_type);
188 if (model_type==
"ghosting" ||
189 !adp->areEntityIDsUnique(adp->getPrimaryEntityType())) {
196 zz->Set_Num_Obj_Fn(zoltanHGNumObj_withModel<Adapter>, (
void *) &(*mdl));
197 zz->Set_Obj_List_Fn(zoltanHGObjList_withModel<Adapter>, (
void *) &(*mdl));
199 zz->Set_HG_Size_CS_Fn(zoltanHGSizeCS_withModel<Adapter>, (
void *) &(*mdl));
200 zz->Set_HG_CS_Fn(zoltanHGCS_withModel<Adapter>, (
void *) &(*mdl));
204 zz->Set_HG_Size_CS_Fn(zoltanHGSizeCS_withMeshAdapter<Adapter>,
206 zz->Set_HG_CS_Fn(zoltanHGCS_withMeshAdapter<Adapter>,
223 const RCP<
const Comm<int> > &problemComm__,
225 env(env__), problemComm(problemComm__), adapter(adapter__)
227 setMPIComm(problemComm__);
229 zz = rcp(
new Zoltan(mpicomm));
234 const RCP<
const Comm<int> > &problemComm__,
236 env(env__), problemComm(problemComm__), adapter(adapter__)
238 setMPIComm(problemComm__);
240 zz = rcp(
new Zoltan(mpicomm));
242 setCallbacksGeom(&(*adapter));
246 const RCP<
const Comm<int> > &problemComm__,
248 env(env__), problemComm(problemComm__), adapter(adapter__)
250 setMPIComm(problemComm__);
252 zz = rcp(
new Zoltan(mpicomm));
254 setCallbacksGraph(adapter);
255 if (adapter->coordinatesAvailable()) {
256 setCallbacksGeom(adapter->getCoordinateInput());
261 const RCP<
const Comm<int> > &problemComm__,
263 env(env__), problemComm(problemComm__), adapter(adapter__)
265 setMPIComm(problemComm__);
267 zz = rcp(
new Zoltan(mpicomm));
269 setCallbacksGraph(adapter);
270 setCallbacksHypergraph(adapter);
271 if (adapter->coordinatesAvailable()) {
272 setCallbacksGeom(adapter->getCoordinateInput());
277 const RCP<
const Comm<int> > &problemComm__,
279 env(env__), problemComm(problemComm__), adapter(adapter__)
281 setMPIComm(problemComm__);
283 zz = rcp(
new Zoltan(mpicomm));
285 setCallbacksGraph(adapter);
289 setCallbacksHypergraph(adapter);
290 setCallbacksGeom(&(*adapter));
299 template <
typename Adapter>
307 size_t numGlobalParts = solution->getTargetGlobalNumberOfParts();
309 sprintf(paramstr,
"%lu", numGlobalParts);
310 zz->Set_Param(
"NUM_GLOBAL_PARTS", paramstr);
312 int wdim = adapter->getNumWeightsPerID();
313 sprintf(paramstr,
"%d", wdim);
314 zz->Set_Param(
"OBJ_WEIGHT_DIM", paramstr);
316 const Teuchos::ParameterList &pl = env->getParameters();
319 const Teuchos::ParameterEntry *pe = pl.getEntryPtr(
"imbalance_tolerance");
322 tolerance = pe->getValue<
double>(&tolerance);
323 sprintf(str,
"%f", tolerance);
324 zz->Set_Param(
"IMBALANCE_TOL", str);
327 pe = pl.getEntryPtr(
"partitioning_approach");
329 std::string approach;
330 approach = pe->getValue<std::string>(&approach);
331 if (approach ==
"partition")
332 zz->Set_Param(
"LB_APPROACH",
"PARTITION");
334 zz->Set_Param(
"LB_APPROACH",
"REPARTITION");
337 pe = pl.getEntryPtr(
"partitioning_objective");
339 std::string strChoice = pe->getValue<std::string>(&strChoice);
340 if (strChoice == std::string(
"multicriteria_minimize_total_weight"))
341 zz->Set_Param(
"RCB_MULTICRITERIA_NORM",
"1");
342 else if (strChoice == std::string(
"multicriteria_balance_total_maximum"))
343 zz->Set_Param(
"RCB_MULTICRITERIA_NORM",
"2");
344 else if (strChoice == std::string(
"multicriteria_minimize_maximum_weight"))
345 zz->Set_Param(
"RCB_MULTICRITERIA_NORM",
"3");
348 pe = pl.getEntryPtr(
"rectilinear");
350 bool val = pe->getValue(&val);
352 zz->Set_Param(
"RCB_RECTILINEAR_BLOCKS",
"1");
357 const Teuchos::ParameterList &zpl = pl.sublist(
"zoltan_parameters");
358 for (ParameterList::ConstIterator iter = zpl.begin();
359 iter != zpl.end(); iter++) {
360 const std::string &zname = pl.name(iter);
362 std::string zval = pl.entry(iter).getValue(&zval);
363 zz->Set_Param(zname.c_str(), zval.c_str());
366 catch (std::exception &e) {
371 int pdim = (wdim > 1 ? wdim : 1);
372 for (
int d = 0; d < pdim; d++) {
373 if (!solution->criteriaHasUniformPartSizes(d)) {
374 float *partsizes =
new float[numGlobalParts];
375 int *partidx =
new int[numGlobalParts];
376 int *wgtidx =
new int[numGlobalParts];
377 for (
size_t i=0; i<numGlobalParts; i++) partidx[i] = i;
378 for (
size_t i=0; i<numGlobalParts; i++) wgtidx[i] = d;
379 for (
size_t i=0; i<numGlobalParts; i++)
380 partsizes[i] = solution->getCriteriaPartSize(0, i);
381 zz->LB_Set_Part_Sizes(1, numGlobalParts, partidx, wgtidx, partsizes);
394 ZOLTAN_ID_PTR dGids = NULL, dLids = NULL;
395 int *dProcs = NULL, *dParts = NULL;
397 ZOLTAN_ID_PTR oGids = NULL, oLids = NULL;
398 int *oProcs = NULL, *oParts = NULL;
400 zz->Set_Param(
"RETURN_LISTS",
"PARTS");
403 int ierr = zz->LB_Partition(changed, nGidEnt, nLidEnt,
404 nDummy, dGids, dLids, dProcs, dParts,
405 nObj, oGids, oLids, oProcs, oParts);
407 env->globalInputAssertion(__FILE__, __LINE__,
"Zoltan LB_Partition",
416 numObjects=model->getLocalNumObjects();
420 ArrayRCP<part_t> partList(
new part_t[numObjects], 0, numObjects,
true);
421 for (
int i = 0; i < nObj; i++) {
424 partList[tmp] = oParts[i];
436 Teuchos::RCP<const map_t> mapWithCopies;
437 Teuchos::RCP<const map_t> oneToOneMap;
440 typedef Tpetra::Vector<scalar_t, lno_t, gno_t> vector_t;
441 vector_t vecWithCopies(mapWithCopies);
442 vector_t oneToOneVec(oneToOneMap);
445 assert(nObj == lno_t(oneToOneMap->getNodeNumElements()));
446 for (lno_t i = 0; i < nObj; i++)
447 oneToOneVec.replaceLocalValue(i, oParts[i]);
450 Teuchos::RCP<const Tpetra::Import<lno_t, gno_t> > importer =
451 Tpetra::createImport<lno_t, gno_t>(oneToOneMap, mapWithCopies);
452 vecWithCopies.doImport(oneToOneVec, *importer, Tpetra::REPLACE);
455 lno_t nlocal = lno_t(mapWithCopies->getNodeNumElements());
456 for (lno_t i = 0; i < nlocal; i++)
457 partList[i] = vecWithCopies.getData()[i];
460 solution->setParts(partList);
463 zz->LB_Free_Part(&oGids, &oLids, &oProcs, &oParts);
AlgZoltan(const RCP< const Environment > &env__, const RCP< const Comm< int > > &problemComm__, const RCP< const VectorAdapter< user_t > > &adapter__)
IdentifierAdapter defines the interface for identifiers.
fast typical checks for valid arguments
MatrixAdapter defines the adapter interface for matrices.
Defines the Model interface.
AlgZoltan(const RCP< const Environment > &env__, const RCP< const Comm< int > > &problemComm__, const RCP< const MatrixAdapter< user_t, userCoord_t > > &adapter__)
GraphAdapter defines the interface for graph-based user data.
MeshAdapter defines the interface for mesh input.
std::bitset< NUM_MODEL_FLAGS > modelFlag_t
Defines the PartitioningSolution class.
void partition(const RCP< PartitioningSolution< Adapter > > &solution)
Partitioning method.
callback functions for the Zoltan package (templated on Adapter)
A PartitioningSolution is a solution to a partitioning problem.
VectorAdapter defines the interface for vector input.
AlgZoltan(const RCP< const Environment > &env__, const RCP< const Comm< int > > &problemComm__, const RCP< const MeshAdapter< user_t > > &adapter__)
Algorithm defines the base class for all algorithms.
static void ASSIGN(first_t &a, second_t b)
Traits class to handle conversions between gno_t/lno_t and TPL data types (e.g., ParMETIS's idx_t...
Gathering definitions used in software development.
void getVertexMaps(Teuchos::RCP< const map_t > &copiesMap, Teuchos::RCP< const map_t > &onetooneMap) const
Sets pointers to the vertex map with copies and the vertex map without copies Note: the pointers will...
The base class for all model classes.
AlgZoltan(const RCP< const Environment > &env__, const RCP< const Comm< int > > &problemComm__, const RCP< const GraphAdapter< user_t, userCoord_t > > &adapter__)
A gathering of useful namespace methods.
HyperGraphModel defines the interface required for hyper graph models.
AlgZoltan(const RCP< const Environment > &env__, const RCP< const Comm< int > > &problemComm__, const RCP< const IdentifierAdapter< user_t > > &adapter__)