Zoltan2
Loading...
Searching...
No Matches
Zoltan2_AlgZoltan.hpp
Go to the documentation of this file.
1// @HEADER
2//
3// ***********************************************************************
4//
5// Zoltan2: A package of combinatorial algorithms for scientific computing
6// Copyright 2012 Sandia Corporation
7//
8// Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
9// the U.S. Government retains certain rights in this software.
10//
11// Redistribution and use in source and binary forms, with or without
12// modification, are permitted provided that the following conditions are
13// met:
14//
15// 1. Redistributions of source code must retain the above copyright
16// notice, this list of conditions and the following disclaimer.
17//
18// 2. Redistributions in binary form must reproduce the above copyright
19// notice, this list of conditions and the following disclaimer in the
20// documentation and/or other materials provided with the distribution.
21//
22// 3. Neither the name of the Corporation nor the names of the
23// contributors may be used to endorse or promote products derived from
24// this software without specific prior written permission.
25//
26// THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
27// EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
28// IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
29// PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
30// CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
31// EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
32// PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
33// PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
34// LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
35// NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
36// SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
37//
38// Questions? Contact Karen Devine (kddevin@sandia.gov)
39// Erik Boman (egboman@sandia.gov)
40// Siva Rajamanickam (srajama@sandia.gov)
41//
42// ***********************************************************************
43//
44// @HEADER
45#ifndef _ZOLTAN2_ALGZOLTAN_HPP_
46#define _ZOLTAN2_ALGZOLTAN_HPP_
47
48#include <Zoltan2_Standards.hpp>
49#include <Zoltan2_Algorithm.hpp>
51#include <Zoltan2_Util.hpp>
52#include <Zoltan2_TPLTraits.hpp>
53
54#include <Zoltan2_Model.hpp>
55
57#include <zoltan_cpp.h>
58#include <zoltan_partition_tree.h>
59
63//
64// This first design templates Zoltan's callback functions on the
65// input adapter. This approach has the advantage of simplicity and
66// is most similar to current usage of Zoltan (where the callbacks define
67// the model).
68// A better approach might template them on a model,
69// allowing Zoltan2 greater flexibility in creating models from the input.
70// Alternatively, different callback implementations could be provided to
71// represent different models to Zoltan.
73
74namespace Zoltan2 {
75
76template <typename Adapter>
77class AlgZoltan : public Algorithm<Adapter>
78{
79
80private:
81
82 typedef typename Adapter::lno_t lno_t;
83 typedef typename Adapter::gno_t gno_t;
84 typedef typename Adapter::scalar_t scalar_t;
85 typedef typename Adapter::part_t part_t;
86 typedef typename Adapter::user_t user_t;
87 typedef typename Adapter::userCoord_t userCoord_t;
88
89 const RCP<const Environment> env;
90 const RCP<const Comm<int> > problemComm;
91 const RCP<const typename Adapter::base_adapter_t> adapter;
92 RCP<const Model<Adapter> > model;
93 RCP<Zoltan> zz;
94
95 MPI_Comm mpicomm;
96
97 void setMPIComm(const RCP<const Comm<int> > &problemComm__) {
98# ifdef HAVE_ZOLTAN2_MPI
99 mpicomm = Teuchos::getRawMpiComm(*problemComm__);
100# else
101 mpicomm = MPI_COMM_WORLD; // taken from siMPI
102# endif
103 }
104
105 void zoltanInit() {
106 // call Zoltan_Initialize to make sure MPI_Init is called (in MPI or siMPI).
107 int argc = 0;
108 char **argv = NULL;
109 float ver;
110 Zoltan_Initialize(argc, argv, &ver);
111 }
112
113 void setCallbacksIDs()
114 {
115 zz->Set_Num_Obj_Fn(zoltanNumObj<Adapter>, (void *) &(*adapter));
116 zz->Set_Obj_List_Fn(zoltanObjList<Adapter>, (void *) &(*adapter));
117
118 const part_t *myparts;
119 adapter->getPartsView(myparts);
120 if (myparts != NULL)
121 zz->Set_Part_Multi_Fn(zoltanParts<Adapter>, (void *) &(*adapter));
122
123 char tmp[4];
125 zz->Set_Param("NUM_GID_ENTRIES", tmp);
127 zz->Set_Param("NUM_LID_ENTRIES", tmp);
128 }
129
130 template <typename AdapterWithCoords>
131 void setCallbacksGeom(const AdapterWithCoords *ia)
132 {
133 // Coordinates may be provided by the MeshAdapter or VectorAdapter.
134 // VectorAdapter may be provided directly by user or indirectly through
135 // GraphAdapter or MatrixAdapter. So separate template type is needed.
136 zz->Set_Num_Geom_Fn(zoltanNumGeom<AdapterWithCoords>, (void *) ia);
137 zz->Set_Geom_Multi_Fn(zoltanGeom<AdapterWithCoords>, (void *) ia);
138 }
139
140 void setCallbacksGraph(
141 const RCP<const GraphAdapter<user_t,userCoord_t> > &/* adp */)
142 {
143 // std::cout << "NotReadyForGraphCallbacksYet" << std::endl;
144 // TODO
145 }
146
147 void setCallbacksGraph(
148 const RCP<const MatrixAdapter<user_t,userCoord_t> > &adp)
149 {
150 // std::cout << "NotReadyForGraphCallbacksYet" << std::endl;
151 // TODO
152 }
153
154 void setCallbacksGraph(
155 const RCP<const MeshAdapter<user_t> > &adp)
156 {
157 // std::cout << "NotReadyForGraphCallbacksYet" << std::endl;
158 // TODO
159 }
160
161 void setCallbacksHypergraph(
162 const RCP<const MatrixAdapter<user_t,userCoord_t> > &adp)
163 {
164 // TODO: If add parameter list to this function, can register
165 // TODO: different callbacks depending on the hypergraph model to use
166
167 zz->Set_HG_Size_CS_Fn(zoltanHGSizeCS_withMatrixAdapter<Adapter>,
168 (void *) &(*adp));
169 zz->Set_HG_CS_Fn(zoltanHGCS_withMatrixAdapter<Adapter>,
170 (void *) &(*adp));
171
172 // zz->Set_HG_Size_Edge_Wts_Fn(zoltanHGSizeEdgeWts_withMatrixAdapter<Adapter>,
173 // (void *) &(*adapter));
174 // zz->Set_HG_Edge_Wts_Fn(zoltanHGSizeEdgeWts_withMatrixAdapter<Adapter>,
175 // (void *) &(*adapter));
176 }
177
178 void setCallbacksHypergraph(
179 const RCP<const GraphAdapter<user_t,userCoord_t> > &adp)
180 {
181 zz->Set_HG_Size_CS_Fn(zoltanHGSizeCS_withGraphAdapter<Adapter>,
182 (void *) &(*adp));
183 zz->Set_HG_CS_Fn(zoltanHGCS_withGraphAdapter<Adapter>,
184 (void *) &(*adp));
185
186 if (adp->getNumWeightsPerEdge() != 0) {
187 if (adp->getNumWeightsPerEdge() > 1) {
188 std::cout << "Zoltan2 warning: getNumWeightsPerEdge() returned "
189 << adp->getNumWeightsPerEdge() << " but PHG supports only "
190 << " one weight per edge; only first weight will be used."
191 << std::endl;
192 }
193 zz->Set_HG_Size_Edge_Wts_Fn(zoltanHGSizeEdgeWts_withGraphAdapter<Adapter>,
194 (void *) &(*adapter));
195 zz->Set_HG_Edge_Wts_Fn(zoltanHGEdgeWts_withGraphAdapter<Adapter>,
196 (void *) &(*adapter));
197 }
198 }
199
200 void setCallbacksHypergraph(const RCP<const MeshAdapter<user_t> > &adp)
201 {
202
203 const Teuchos::ParameterList &pl = env->getParameters();
204
205 const Teuchos::ParameterEntry *pe = pl.getEntryPtr("hypergraph_model_type");
206 std::string model_type("traditional");
207 if (pe){
208 model_type = pe->getValue<std::string>(&model_type);
209 }
210
211 if (model_type=="ghosting" ||
212 !adp->areEntityIDsUnique(adp->getPrimaryEntityType())) {
214 HyperGraphModel<Adapter>* mdl = new HyperGraphModel<Adapter>(adp, env,
215 problemComm, flags,
217 model = rcp(static_cast<const Model<Adapter>* >(mdl),true);
218
219 zz->Set_Num_Obj_Fn(zoltanHGNumObj_withModel<Adapter>, (void *) &(*mdl));
220 zz->Set_Obj_List_Fn(zoltanHGObjList_withModel<Adapter>, (void *) &(*mdl));
221
222 zz->Set_HG_Size_CS_Fn(zoltanHGSizeCS_withModel<Adapter>, (void *) &(*mdl));
223 zz->Set_HG_CS_Fn(zoltanHGCS_withModel<Adapter>, (void *) &(*mdl));
224 }
225 else {
226 //If entities are unique we dont need the extra cost of the model
227 zz->Set_HG_Size_CS_Fn(zoltanHGSizeCS_withMeshAdapter<Adapter>,
228 (void *) &(*adp));
229 zz->Set_HG_CS_Fn(zoltanHGCS_withMeshAdapter<Adapter>,
230 (void *) &(*adp));
231 }
232 // zz->Set_HG_Size_Edge_Wts_Fn(zoltanHGSizeEdgeWts_withMeshAdapter<Adapter>,
233 // (void *) &(*adp));
234 // zz->Set_HG_Edge_Wts_Fn(zoltanHGSizeEdgeWts_withMeshAdapter<Adapter>,
235 // (void *) &(*adp));
236 }
237
239 virtual bool isPartitioningTreeBinary() const
240 {
241 return true;
242 }
243
245 void rcb_recursive_partitionTree_calculations(
246 part_t arrayIndex,
247 part_t numParts,
248 std::vector<part_t> & splitRangeBeg,
249 std::vector<part_t> & splitRangeEnd) const
250 {
251 // Note the purpose of the recursive method is make sure the children of a
252 // node have updated their values for splitRangeBeg and splitRangeEnd
253 // Then we can set our own values simply based on the union
254 // first load the rcb data for the node
255 int parent = -1;
256 int left_leaf = -1;
257 int right_leaf = -1;
258 int err = Zoltan_RCB_Partition_Tree(zz->Get_C_Handle(),
259 arrayIndex - numParts + 1, // rcb starts as 1 but does not include terminals
260 &parent, &left_leaf, &right_leaf);
261 if(err != 0) {
262 throw std::logic_error( "Zoltan_RCB_Partition_Tree returned in error." );
263 }
264 // check that children both have their ranges set and if not, do those
265 // range first so we can build them to make our range
266 if(left_leaf > 0) { // neg is terminal and always already built
267 rcb_recursive_partitionTree_calculations(left_leaf+numParts-1, numParts,
268 splitRangeBeg, splitRangeEnd);
269 }
270 if(right_leaf > 0) { // neg is terminal and always already built
271 rcb_recursive_partitionTree_calculations(right_leaf+numParts-1, numParts,
272 splitRangeBeg, splitRangeEnd);
273 }
274 // now we can build our ranges from the children
275 // note this exploits the rcb conventions for right and left so we know
276 // that left_leaf will be our smaller indices
277 int leftIndex = (left_leaf > 0) ? (left_leaf-1+numParts) : (-left_leaf);
278 int rightIndex = (right_leaf > 0) ? (right_leaf-1+numParts) : (-right_leaf);
279 splitRangeBeg[arrayIndex] = splitRangeBeg[leftIndex];
280 splitRangeEnd[arrayIndex] = splitRangeEnd[rightIndex];
281 // for debugging sanity check verify left_leaf is a set of indices which
282 // goes continuously into the right_leaf
283 if(splitRangeBeg[rightIndex] != splitRangeEnd[leftIndex]) { // end is non-inclusive
284 throw std::logic_error( "RCB expected left_leaf indices and right leaf"
285 " indices to be continuous but it was not so." );
286 }
287 }
288
290 void rcb_getPartitionTree(part_t numParts,
291 part_t & numTreeVerts,
292 std::vector<part_t> & permPartNums,
293 std::vector<part_t> & splitRangeBeg,
294 std::vector<part_t> & splitRangeEnd,
295 std::vector<part_t> & treeVertParents) const
296 {
297 // CALCULATE: numTreeVerts
298 // For rcb a tree node always takes 2 nodes and turns them into 1 node
299 // That means it takes numParts - 1 nodes to reduce a tree of numParts to
300 // a single root node - but we do 2 * numParts - 1 because we are currently
301 // treating all of the 'trivial' terminals as tree nodes - something we
302 // discussed we may change later
303 part_t numTreeNodes = 2 * numParts - 1;
304 numTreeVerts = numTreeNodes - 1; // by design convention root not included
305 // CALCULATE: permPartNums
306 permPartNums.resize(numParts);
307 for(part_t n = 0; n < numParts; ++n) {
308 permPartNums[n] = n; // for rcb we can assume this is trivial and in order
309 }
310 // CALCULATE: treeVertParents
311 treeVertParents.resize(numTreeNodes); // allocate space for numTreeNodes array
312 // scan all the non terminal nodes and set all children to have us as parent
313 // that will set all parents except for the root node which we will set to -1
314 // track the children of the root and final node for swapping later. Couple
315 // ways to do this - all seem a bit awkward but at least this is efficient.
316 part_t rootNode = 0; // track index of the root node for swapping
317 // a bit awkward but efficient - save the children of root and final node
318 // for swap at end to satisfy convention that root is highest index node
319 part_t saveRootNodeChildrenA = -1;
320 part_t saveRootNodeChildrenB = -1;
321 part_t saveFinalNodeChildrenA = -1;
322 part_t saveFinalNodeChildrenB = -1;
323 for(part_t n = numParts; n < numTreeNodes; ++n) { // scan and set all parents
324 int parent = -1;
325 int left_leaf = -1;
326 int right_leaf = -1;
327 int err = Zoltan_RCB_Partition_Tree(zz->Get_C_Handle(),
328 static_cast<int>(n - numParts) + 1, // rcb starts as 1 but does not include terminals
329 &parent, &left_leaf, &right_leaf);
330 if(err != 0) {
331 throw std::logic_error("Zoltan_RCB_Partition_Tree returned in error.");
332 }
333 part_t leftIndex = (left_leaf > 0) ? (left_leaf-1+numParts) : (-left_leaf);
334 part_t rightIndex = (right_leaf > 0) ? (right_leaf-1+numParts) : (-right_leaf);
335 treeVertParents[leftIndex] = n;
336 treeVertParents[rightIndex] = n;
337 // save root node for final swap
338 if(parent == 1 || parent == -1) { // is it the root?
339 rootNode = n; // remember I am the root
340 saveRootNodeChildrenA = leftIndex;
341 saveRootNodeChildrenB = rightIndex;
342 }
343 if(n == numTreeNodes-1) {
344 saveFinalNodeChildrenA = leftIndex;
345 saveFinalNodeChildrenB = rightIndex;
346 }
347 }
348 treeVertParents[rootNode] = -1; // convention parent is root -1
349 // splitRangeBeg and splitRangeEnd
350 splitRangeBeg.resize(numTreeNodes);
351 splitRangeEnd.resize(numTreeNodes);
352 // for terminal nodes this is trivial
353 for(part_t n = 0; n < numParts; ++n) {
354 splitRangeBeg[n] = n;
355 splitRangeEnd[n] = n + 1;
356 }
357 if(numParts > 1) { // not relevant for 1 part
358 // build the splitRangeBeg and splitRangeEnd recursively which forces the
359 // children of each node to be calculated before the parent so parent can
360 // just take the union of the two children
361 rcb_recursive_partitionTree_calculations(rootNode, numParts, splitRangeBeg,
362 splitRangeEnd);
363 // now as a final step handle the swap to root is the highest index node
364 // swap the parent of the two nodes
365 std::swap(treeVertParents[rootNode], treeVertParents[numTreeNodes-1]);
366 // get the children of the swapped nodes to have updated parents
367 treeVertParents[saveFinalNodeChildrenA] = rootNode;
368 treeVertParents[saveFinalNodeChildrenB] = rootNode;
369 // handle case where final node is child of the root
370 if(saveRootNodeChildrenA == numTreeNodes - 1) {
371 saveRootNodeChildrenA = rootNode;
372 }
373 if(saveRootNodeChildrenB == numTreeNodes - 1) {
374 saveRootNodeChildrenB = rootNode;
375 }
376 treeVertParents[saveRootNodeChildrenA] = numTreeNodes - 1;
377 treeVertParents[saveRootNodeChildrenB] = numTreeNodes - 1;
378 // update the beg and end indices - simply swap them
379 std::swap(splitRangeBeg[rootNode], splitRangeBeg[numTreeNodes-1]);
380 std::swap(splitRangeEnd[rootNode], splitRangeEnd[numTreeNodes-1]);
381 }
382 }
383
385 void phg_getPartitionTree(part_t numParts,
386 part_t & numTreeVerts,
387 std::vector<part_t> & permPartNums,
388 std::vector<part_t> & splitRangeBeg,
389 std::vector<part_t> & splitRangeEnd,
390 std::vector<part_t> & treeVertParents) const
391 {
392 // First thing is to get the length of the tree from zoltan.
393 // The tree is a list of pairs (lo,hi) for each node.
394 // Here tree_array_size is the number of pairs.
395 // In phg indexing the first pairt (i=0) is empty garbage.
396 // The second pair (index 1) will be the root.
397 // Some nodes will be empty nodes, determined by hi = -1.
398 int tree_array_size = -1; // will be set by Zoltan_PHG_Partition_Tree_Size
399 int err = Zoltan_PHG_Partition_Tree_Size(
400 zz->Get_C_Handle(), &tree_array_size);
401 if(err != 0) {
402 throw std::logic_error("Zoltan_PHG_Partition_Tree_Size returned error.");
403 }
404 // Determine the number of valid nodes (PHG will have empty slots)
405 // We scan the list of pairs and count each node which does not have hi = -1
406 // During the loop we will also construct mapIndex which maps initial n
407 // to final n due to some conversions we apply to meet the design specs.
408 // The conversions implemented by mapIndex are:
409 // Move all terminals to the beginning (terminals have hi = lo)
410 // Resort the terminals in order (simply map to index lo works)
411 // Move non-terminals after the terminals (they start at index numParts)
412 // Map the first pair (root) to the be last to meet the design spec
413 part_t numTreeNodes = 0;
414 std::vector<part_t> mapIndex(tree_array_size); // maps n to final index
415 part_t trackNonTerminalIndex = numParts; // starts after terminals
416 for(part_t n = 0; n < static_cast<part_t>(tree_array_size); ++n) {
417 part_t phgIndex = n + 1; // phg indexing starts at 1
418 int lo_index = -1;
419 int hi_index = -1;
420 err = Zoltan_PHG_Partition_Tree(
421 zz->Get_C_Handle(), phgIndex, &lo_index, &hi_index);
422 if(hi_index != -1) { // hi -1 means it's an unused node
423 ++numTreeNodes; // increase the total count because this is a real node
424 if(n != 0) { // the root is mapped last but we don't know the length yet
425 mapIndex[n] = (lo_index == hi_index) ? // is it a terminal?
426 lo_index : // terminals map in sequence - lo_index is correct
427 (trackNonTerminalIndex++); // set then bump trackNonTerminalIndex +1
428 }
429 }
430 }
431 // now complete the map by setting root to the length-1 for the design specs
432 mapIndex[0] = numTreeNodes - 1;
433 // CALCULATE: numTreeVerts
434 numTreeVerts = numTreeNodes - 1; // this is the design - root not included
435 // CALCULATE: permPartNums
436 permPartNums.resize(numParts);
437 for(part_t n = 0; n < numParts; ++n) {
438 permPartNums[n] = n; // for phg we can assume this is trivial and in order
439 }
440 // CALCULATE: treeVertParents, splitRangeBeg, splitRangeEnd
441 // we will determine all of these in this second loop using mapIndex
442 // First set the arrays to have the proper length
443 treeVertParents.resize(numTreeNodes);
444 splitRangeBeg.resize(numTreeNodes);
445 splitRangeEnd.resize(numTreeNodes);
446 // Now loop a second time
447 for(part_t n = 0; n < tree_array_size; ++n) {
448 part_t phgIndex = n + 1; // phg indexing starts at 1
449 int lo_index = -1; // zoltan Zoltan_PHG_Partition_Tree will set this
450 int hi_index = -1; // zoltan Zoltan_PHG_Partition_Tree will set this
451 err = Zoltan_PHG_Partition_Tree( // access zoltan phg tree data
452 zz->Get_C_Handle(), phgIndex, &lo_index, &hi_index);
453 if(err != 0) {
454 throw std::logic_error("Zoltan_PHG_Partition_Tree returned in error.");
455 }
456 if(hi_index != -1) { // hi -1 means it's an unused node (a gap)
457 // get final index using mapIndex - so convert from phg to design plan
458 part_t finalIndex = mapIndex[n]; // get the index for the final output
459 // now determine the parent
460 // In the original phg indexing, the parent can be directly calculated
461 // from the pair index using the following rules:
462 // if phgIndex even, then parent is phgIndex/2
463 // here we determine even by ((phgIndex%2) == 0)
464 // if phgIndex odd, then parent is (phgIndex-1)/2
465 // but after getting parentPHGIndex we convert back to this array
466 // indexing by subtracting 1
467 part_t parentPHGIndex =
468 ((phgIndex%2) == 0) ? (phgIndex/2) : ((phgIndex-1)/2);
469 // map the parent phg index to the final parent index
470 // however, for the special case of the root (n=1), set the parent to -1
471 treeVertParents[finalIndex] = (n==0) ? -1 : mapIndex[parentPHGIndex-1];
472 // set begin (inclusive) and end (non-inclusive), so end is hi+1
473 splitRangeBeg[finalIndex] = static_cast<part_t>(lo_index);
474 splitRangeEnd[finalIndex] = static_cast<part_t>(hi_index+1);
475 }
476 }
477 }
478
480 void getPartitionTree(part_t numParts,
481 part_t & numTreeVerts,
482 std::vector<part_t> & permPartNums,
483 std::vector<part_t> & splitRangeBeg,
484 std::vector<part_t> & splitRangeEnd,
485 std::vector<part_t> & treeVertParents) const
486 {
487 // first check that our parameters requested we keep the tree
488 const Teuchos::ParameterList &pl = env->getParameters();
489 bool keep_partition_tree = false;
490 const Teuchos::ParameterEntry * pe = pl.getEntryPtr("keep_partition_tree");
491 if(pe) {
492 keep_partition_tree = pe->getValue(&keep_partition_tree);
493 if(!keep_partition_tree) {
494 throw std::logic_error(
495 "Requested tree when param keep_partition_tree is off.");
496 }
497 }
498
499 // now call the appropriate method based on LB_METHOD in the zoltan
500 // parameters list.
501 const Teuchos::ParameterList & zoltan_pl = pl.sublist("zoltan_parameters");
502 std::string lb_method;
503 pe = zoltan_pl.getEntryPtr("LB_METHOD");
504 if(pe) {
505 lb_method = pe->getValue(&lb_method);
506 }
507 if(lb_method == "phg") {
508 phg_getPartitionTree(numParts, numTreeVerts, permPartNums,
509 splitRangeBeg, splitRangeEnd, treeVertParents);
510 }
511 else if(lb_method == "rcb") {
512 rcb_getPartitionTree(numParts, numTreeVerts, permPartNums,
513 splitRangeBeg, splitRangeEnd, treeVertParents);
514 }
515 else {
516 throw std::logic_error("Did not recognize LB_METHOD: " + lb_method);
517 }
518 }
519
520public:
521
527 AlgZoltan(const RCP<const Environment> &env__,
528 const RCP<const Comm<int> > &problemComm__,
529 const RCP<const IdentifierAdapter<user_t> > &adapter__):
530 env(env__), problemComm(problemComm__), adapter(adapter__)
531 {
532 setMPIComm(problemComm__);
533 zoltanInit();
534 zz = rcp(new Zoltan(mpicomm));
535 setCallbacksIDs();
536 }
537
538 AlgZoltan(const RCP<const Environment> &env__,
539 const RCP<const Comm<int> > &problemComm__,
540 const RCP<const VectorAdapter<user_t> > &adapter__) :
541 env(env__), problemComm(problemComm__), adapter(adapter__)
542 {
543 setMPIComm(problemComm__);
544 zoltanInit();
545 zz = rcp(new Zoltan(mpicomm));
546 setCallbacksIDs();
547 setCallbacksGeom(&(*adapter));
548 }
549
550 AlgZoltan(const RCP<const Environment> &env__,
551 const RCP<const Comm<int> > &problemComm__,
552 const RCP<const GraphAdapter<user_t,userCoord_t> > &adapter__) :
553 env(env__), problemComm(problemComm__), adapter(adapter__)
554 {
555 setMPIComm(problemComm__);
556 zoltanInit();
557 zz = rcp(new Zoltan(mpicomm));
558 setCallbacksIDs();
559 setCallbacksGraph(adapter);
560 setCallbacksHypergraph(adapter);
561 if (adapter->coordinatesAvailable()) {
562 setCallbacksGeom(adapter->getCoordinateInput());
563 }
564 }
565
566 AlgZoltan(const RCP<const Environment> &env__,
567 const RCP<const Comm<int> > &problemComm__,
568 const RCP<const MatrixAdapter<user_t,userCoord_t> > &adapter__) :
569 env(env__), problemComm(problemComm__), adapter(adapter__)
570 {
571 setMPIComm(problemComm__);
572 zoltanInit();
573 zz = rcp(new Zoltan(mpicomm));
574 setCallbacksIDs();
575 setCallbacksGraph(adapter);
576 setCallbacksHypergraph(adapter);
577 if (adapter->coordinatesAvailable()) {
578 setCallbacksGeom(adapter->getCoordinateInput());
579 }
580 }
581
582 AlgZoltan(const RCP<const Environment> &env__,
583 const RCP<const Comm<int> > &problemComm__,
584 const RCP<const MeshAdapter<user_t> > &adapter__) :
585 env(env__), problemComm(problemComm__), adapter(adapter__)
586 {
587 setMPIComm(problemComm__);
588 zoltanInit();
589 zz = rcp(new Zoltan(mpicomm));
590 setCallbacksIDs();
591 setCallbacksGraph(adapter);
592 //TODO:: check parameter list to see if hypergraph is needed. We dont want to build the model
593 // if we don't have to and we shouldn't as it can take a decent amount of time if the
594 // primary entity is copied
595 setCallbacksHypergraph(adapter);
596 setCallbacksGeom(&(*adapter));
597 }
598
599 void partition(const RCP<PartitioningSolution<Adapter> > &solution);
600 // void color(const RCP<ColoringSolution<Adapter> > &solution);
601};
602
604template <typename Adapter>
606 const RCP<PartitioningSolution<Adapter> > &solution
607)
608{
609 HELLO;
610 char paramstr[128];
611
612 size_t numGlobalParts = solution->getTargetGlobalNumberOfParts();
613
614 sprintf(paramstr, "%lu", numGlobalParts);
615 zz->Set_Param("NUM_GLOBAL_PARTS", paramstr);
616
617 int wdim = adapter->getNumWeightsPerID();
618 sprintf(paramstr, "%d", wdim);
619 zz->Set_Param("OBJ_WEIGHT_DIM", paramstr);
620
621 const Teuchos::ParameterList &pl = env->getParameters();
622
623 double tolerance;
624 const Teuchos::ParameterEntry *pe = pl.getEntryPtr("imbalance_tolerance");
625 if (pe){
626 char str[30];
627 tolerance = pe->getValue<double>(&tolerance);
628 sprintf(str, "%f", tolerance);
629 zz->Set_Param("IMBALANCE_TOL", str);
630 }
631
632 pe = pl.getEntryPtr("partitioning_approach");
633 if (pe){
634 std::string approach;
635 approach = pe->getValue<std::string>(&approach);
636 if (approach == "partition")
637 zz->Set_Param("LB_APPROACH", "PARTITION");
638 else
639 zz->Set_Param("LB_APPROACH", "REPARTITION");
640 }
641
642 pe = pl.getEntryPtr("partitioning_objective");
643 if (pe){
644 std::string strChoice = pe->getValue<std::string>(&strChoice);
645 if (strChoice == std::string("multicriteria_minimize_total_weight"))
646 zz->Set_Param("RCB_MULTICRITERIA_NORM", "1");
647 else if (strChoice == std::string("multicriteria_balance_total_maximum"))
648 zz->Set_Param("RCB_MULTICRITERIA_NORM", "2");
649 else if (strChoice == std::string("multicriteria_minimize_maximum_weight"))
650 zz->Set_Param("RCB_MULTICRITERIA_NORM", "3");
651 }
652
653 // perhaps make this a bool stored in the AlgZoltan if we want to follow
654 // the pattern of multijagged mj_keep_part_boxes for example. However we can
655 // collect the error straight from Zoltan if we attempt to access the tree
656 // when we never stored it so that may not be necessary
657 bool keep_partition_tree = false;
658 pe = pl.getEntryPtr("keep_partition_tree");
659 if (pe) {
660 keep_partition_tree = pe->getValue(&keep_partition_tree);
661 if (keep_partition_tree) {
662 // need to resolve the organization of this
663 // when do we want to use the zoltan parameters directly versus
664 // using the zoltan2 parameters like this
665 zz->Set_Param("KEEP_CUTS", "1"); // rcb zoltan setting
666 zz->Set_Param("PHG_KEEP_TREE", "1"); // phg zoltan setting
667 }
668 }
669
670 pe = pl.getEntryPtr("rectilinear");
671 if (pe) {
672 bool val = pe->getValue(&val);
673 if (val)
674 zz->Set_Param("RCB_RECTILINEAR_BLOCKS", "1");
675 }
676
677 // Look for zoltan_parameters sublist; pass all zoltan parameters to Zoltan
678 try {
679 const Teuchos::ParameterList &zpl = pl.sublist("zoltan_parameters");
680 for (ParameterList::ConstIterator iter = zpl.begin();
681 iter != zpl.end(); iter++) {
682 const std::string &zname = pl.name(iter);
683 // Convert the value to a string to pass to Zoltan
684 std::string zval = pl.entry(iter).getValue(&zval);
685 zz->Set_Param(zname.c_str(), zval.c_str());
686 }
687 }
688 catch (std::exception &) {
689 // No zoltan_parameters sublist found; no Zoltan parameters to register
690 }
691
692 // Get target part sizes
693 int pdim = (wdim > 1 ? wdim : 1);
694 for (int d = 0; d < pdim; d++) {
695 if (!solution->criteriaHasUniformPartSizes(d)) {
696 float *partsizes = new float[numGlobalParts];
697 int *partidx = new int[numGlobalParts];
698 int *wgtidx = new int[numGlobalParts];
699 for (size_t i=0; i<numGlobalParts; i++) partidx[i] = i;
700 for (size_t i=0; i<numGlobalParts; i++) wgtidx[i] = d;
701 for (size_t i=0; i<numGlobalParts; i++)
702 partsizes[i] = solution->getCriteriaPartSize(0, i);
703 zz->LB_Set_Part_Sizes(1, numGlobalParts, partidx, wgtidx, partsizes);
704 delete [] partsizes;
705 delete [] partidx;
706 delete [] wgtidx;
707 }
708 }
709
710 // Make the call to LB_Partition
711 int changed = 0;
714
715 int nDummy = -1; // Dummy vars to satisfy arglist
716 ZOLTAN_ID_PTR dGids = NULL, dLids = NULL;
717 int *dProcs = NULL, *dParts = NULL;
718 int nObj = -1; // Output vars with new part info
719 ZOLTAN_ID_PTR oGids = NULL, oLids = NULL;
720 int *oProcs = NULL, *oParts = NULL;
721
722 zz->Set_Param("RETURN_LISTS", "PARTS"); // required format for Zoltan2;
723 // results in last five arguments
724
725 int ierr = zz->LB_Partition(changed, nGidEnt, nLidEnt,
726 nDummy, dGids, dLids, dProcs, dParts,
727 nObj, oGids, oLids, oProcs, oParts);
728
729 env->globalInputAssertion(__FILE__, __LINE__, "Zoltan LB_Partition",
730 (ierr==ZOLTAN_OK || ierr==ZOLTAN_WARN), BASIC_ASSERTION, problemComm);
731
732 int numObjects=nObj;
733 // The number of objects may be larger than zoltan knows due to copies that
734 // were removed by the hypergraph model
735 if (model!=RCP<const Model<Adapter> >() &&
736 dynamic_cast<const HyperGraphModel<Adapter>* >(&(*model)) &&
737 !(dynamic_cast<const HyperGraphModel<Adapter>* >(&(*model))->areVertexIDsUnique())) {
738 numObjects=model->getLocalNumObjects();
739 }
740
741 // Load answer into the solution.
742 ArrayRCP<part_t> partList(new part_t[numObjects], 0, numObjects, true);
743 for (int i = 0; i < nObj; i++) {
744 lno_t tmp;
745 TPL_Traits<lno_t, ZOLTAN_ID_PTR>::ASSIGN(tmp, &(oLids[i*nLidEnt]));
746 partList[tmp] = oParts[i];
747 }
748
749 if (model!=RCP<const Model<Adapter> >() &&
750 dynamic_cast<const HyperGraphModel<Adapter>* >(&(*model)) &&
751 !(dynamic_cast<const HyperGraphModel<Adapter>* >(&(*model))->areVertexIDsUnique())) {
752 // Setup the part ids for copied entities removed by ownership in
753 // hypergraph model.
754 const HyperGraphModel<Adapter>* mdl =
755 static_cast<const HyperGraphModel<Adapter>* >(&(*model));
756
757 typedef typename HyperGraphModel<Adapter>::map_t map_t;
758 Teuchos::RCP<const map_t> mapWithCopies;
759 Teuchos::RCP<const map_t> oneToOneMap;
760 mdl->getVertexMaps(mapWithCopies,oneToOneMap);
761
762 typedef Tpetra::Vector<scalar_t, lno_t, gno_t> vector_t;
763 vector_t vecWithCopies(mapWithCopies);
764 vector_t oneToOneVec(oneToOneMap);
765
766 // Set values in oneToOneVec: each entry == rank
767 assert(nObj == lno_t(oneToOneMap->getLocalNumElements()));
768 for (lno_t i = 0; i < nObj; i++)
769 oneToOneVec.replaceLocalValue(i, oParts[i]);
770
771 // Now import oneToOneVec's values back to vecWithCopies
772 Teuchos::RCP<const Tpetra::Import<lno_t, gno_t> > importer =
773 Tpetra::createImport<lno_t, gno_t>(oneToOneMap, mapWithCopies);
774 vecWithCopies.doImport(oneToOneVec, *importer, Tpetra::REPLACE);
775
776 // Should see copied vector values when print VEC WITH COPIES
777 lno_t nlocal = lno_t(mapWithCopies->getLocalNumElements());
778 for (lno_t i = 0; i < nlocal; i++)
779 partList[i] = vecWithCopies.getData()[i];
780 }
781
782 solution->setParts(partList);
783
784 // Clean up
785 zz->LB_Free_Part(&oGids, &oLids, &oProcs, &oParts);
786}
787
788} // namespace Zoltan2
789
790#endif
callback functions for the Zoltan package (templated on Adapter)
Defines the Model interface.
Defines the PartitioningSolution class.
Gathering definitions used in software development.
#define HELLO
Traits class to handle conversions between gno_t/lno_t and TPL data types (e.g., ParMETIS's idx_t,...
A gathering of useful namespace methods.
AlgZoltan(const RCP< const Environment > &env__, const RCP< const Comm< int > > &problemComm__, const RCP< const MeshAdapter< user_t > > &adapter__)
AlgZoltan(const RCP< const Environment > &env__, const RCP< const Comm< int > > &problemComm__, const RCP< const IdentifierAdapter< user_t > > &adapter__)
void partition(const RCP< PartitioningSolution< Adapter > > &solution)
Partitioning method.
AlgZoltan(const RCP< const Environment > &env__, const RCP< const Comm< int > > &problemComm__, const RCP< const MatrixAdapter< user_t, userCoord_t > > &adapter__)
AlgZoltan(const RCP< const Environment > &env__, const RCP< const Comm< int > > &problemComm__, const RCP< const VectorAdapter< user_t > > &adapter__)
AlgZoltan(const RCP< const Environment > &env__, const RCP< const Comm< int > > &problemComm__, const RCP< const GraphAdapter< user_t, userCoord_t > > &adapter__)
Algorithm defines the base class for all algorithms.
HyperGraphModel defines the interface required for hyper graph models.
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...
map_t::local_ordinal_type lno_t
Tpetra::Map map_t
Created by mbenlioglu on Aug 31, 2020.
std::bitset< NUM_MODEL_FLAGS > modelFlag_t
@ BASIC_ASSERTION
fast typical checks for valid arguments
static void ASSIGN(first_t &a, second_t b)