50 #ifndef _ZOLTAN2_PARTITIONINGPROBLEM_HPP_
51 #define _ZOLTAN2_PARTITIONINGPROBLEM_HPP_
62 #ifdef ZOLTAN2_TASKMAPPING_MOVE
74 #ifdef HAVE_ZOLTAN2_OVIS
103 template<
typename Adapter>
117 const RCP<
const Teuchos::Comm<int> > &comm):
128 #ifdef HAVE_ZOLTAN2_MPI
133 rcp<const Comm<int> >(new Teuchos::MpiComm<int>(
134 Teuchos::opaqueWrapper(mpicomm))))
164 void solve(
bool updateInputData=
true);
253 scalar_t *partSizes,
bool makeCopy=
true) ;
276 Array<std::string> algorithm_names(19);
277 algorithm_names[0] =
"rcb";
278 algorithm_names[1] =
"multijagged";
279 algorithm_names[2] =
"rib";
280 algorithm_names[3] =
"hsfc";
281 algorithm_names[4] =
"patoh";
282 algorithm_names[5] =
"phg";
283 algorithm_names[6] =
"metis";
284 algorithm_names[7] =
"parmetis";
285 algorithm_names[8] =
"quotient";
286 algorithm_names[9] =
"pulp";
287 algorithm_names[10] =
"parma";
288 algorithm_names[11] =
"scotch";
289 algorithm_names[12] =
"ptscotch";
290 algorithm_names[13] =
"block";
291 algorithm_names[14] =
"cyclic";
292 algorithm_names[15] =
"sarma";
293 algorithm_names[16] =
"random";
294 algorithm_names[17] =
"zoltan";
295 algorithm_names[18] =
"forTestingOnly";
296 RCP<Teuchos::StringValidator> algorithm_Validator = Teuchos::rcp(
297 new Teuchos::StringValidator( algorithm_names ));
298 pl.set(
"algorithm",
"random",
"partitioning algorithm",
299 algorithm_Validator);
302 pl.set(
"keep_partition_tree",
false,
"If true, will keep partition tree",
306 pl.set(
"rectilinear",
false,
"If true, then when a cut is made, all of the "
307 "dots located on the cut are moved to the same side of the cut. The "
308 "resulting regions are then rectilinear. The resulting load balance may "
309 "not be as good as when the group of dots is split by the cut. ",
312 RCP<Teuchos::StringValidator> partitioning_objective_Validator =
313 Teuchos::rcp(
new Teuchos::StringValidator(
314 Teuchos::tuple<std::string>(
"balance_object_count",
315 "balance_object_weight",
"multicriteria_minimize_total_weight",
316 "multicriteria_minimize_maximum_weight",
317 "multicriteria_balance_total_maximum",
"minimize_cut_edge_count",
318 "minimize_cut_edge_weight",
"minimize_neighboring_parts",
319 "minimize_boundary_vertices" )));
320 pl.set(
"partitioning_objective",
"balance_object_weight",
321 "objective of partitioning", partitioning_objective_Validator);
323 pl.set(
"imbalance_tolerance", 1.1,
"imbalance tolerance, ratio of "
327 RCP<Teuchos::EnhancedNumberValidator<int>> num_global_parts_Validator =
328 Teuchos::rcp(
new Teuchos::EnhancedNumberValidator<int>(
329 1, Teuchos::EnhancedNumberTraits<int>::max()) );
330 pl.set(
"num_global_parts", 1,
"global number of parts to compute "
331 "(0 means use the number of processes)", num_global_parts_Validator);
334 RCP<Teuchos::EnhancedNumberValidator<int>> num_local_parts_Validator =
335 Teuchos::rcp(
new Teuchos::EnhancedNumberValidator<int>(
336 0, Teuchos::EnhancedNumberTraits<int>::max()) );
337 pl.set(
"num_local_parts", 0,
"number of parts to compute for this "
338 "process (num_global_parts == sum of all num_local_parts)",
339 num_local_parts_Validator);
341 RCP<Teuchos::StringValidator> partitioning_approach_Validator =
342 Teuchos::rcp(
new Teuchos::StringValidator(
343 Teuchos::tuple<std::string>(
"partition",
"repartition",
344 "maximize_overlap" )));
345 pl.set(
"partitioning_approach",
"partition",
"Partition from scratch, "
346 "partition incrementally from current partition, of partition from "
347 "scratch but maximize overlap with the current partition",
348 partitioning_approach_Validator);
350 RCP<Teuchos::StringValidator> objects_to_partition_Validator =
351 Teuchos::rcp(
new Teuchos::StringValidator(
352 Teuchos::tuple<std::string>(
"matrix_rows",
"matrix_columns",
353 "matrix_nonzeros",
"mesh_elements",
"mesh_nodes",
"graph_edges",
354 "graph_vertices",
"coordinates",
"identifiers" )));
355 pl.set(
"objects_to_partition",
"graph_vertices",
"Objects to be partitioned",
356 objects_to_partition_Validator);
358 RCP<Teuchos::StringValidator> model_Validator = Teuchos::rcp(
359 new Teuchos::StringValidator(
360 Teuchos::tuple<std::string>(
"hypergraph",
"graph",
361 "geometry",
"ids" )));
362 pl.set(
"model",
"graph",
"This is a low level parameter. Normally the "
363 "library will choose a computational model based on the algorithm or "
364 "objective specified by the user.", model_Validator);
367 pl.set(
"remap_parts",
false,
"remap part numbers to minimize migration "
370 pl.set(
"mapping_type", -1,
"Mapping of solution to the processors. -1 No"
373 RCP<Teuchos::EnhancedNumberValidator<int>> ghost_layers_Validator =
374 Teuchos::rcp(
new Teuchos::EnhancedNumberValidator<int>(1, 10, 1, 0) );
375 pl.set(
"ghost_layers", 2,
"number of layers for ghosting used in "
376 "hypergraph ghost method", ghost_layers_Validator);
381 virtual void processAlgorithmName(
const std::string& algorithm,
const std::string& defString,
const std::string& model,
382 Environment &env,
bool& removeSelfEdges,
bool &isGraphType,
bool& needConsecutiveGlobalIds);
388 #ifdef ZOLTAN2_TASKMAPPING_MOVE
389 RCP<MachineRepresentation<scalar_t,part_t> > machine_;
429 template <
typename Adapter>
434 this->env_->debug(
DETAILED_STATUS,
"PartitioningProblem::initializeProblem");
436 if (getenv(
"DEBUGME")){
438 std::cout << getpid() << std::endl;
441 std::cout << _getpid() << std::endl;
446 #ifdef HAVE_ZOLTAN2_OVIS
447 ovis_enabled(this->comm_->getRank());
452 #ifdef ZOLTAN2_TASKMAPPING_MOVE
453 machine_ = RCP<MachineRepresentation<scalar_t,part_t> >(
460 numberOfWeights_ = this->inputAdapter_->getNumWeightsPerID();
462 numberOfCriteria_ = (numberOfWeights_ > 1) ? numberOfWeights_ : 1;
464 inputType_ = this->inputAdapter_->adapterType();
469 ArrayRCP<part_t> *noIds =
new ArrayRCP<part_t> [numberOfCriteria_];
470 ArrayRCP<scalar_t> *noSizes =
new ArrayRCP<scalar_t> [numberOfCriteria_];
472 partIds_ = arcp(noIds, 0, numberOfCriteria_,
true);
473 partSizes_ = arcp(noSizes, 0, numberOfCriteria_,
true);
476 std::ostringstream msg;
477 msg << this->comm_->getSize() <<
" procs,"
478 << numberOfWeights_ <<
" user-defined weights\n";
482 this->env_->memory(
"After initializeProblem");
485 template <
typename Adapter>
487 int criteria,
int len,
part_t *partIds,
scalar_t *partSizes,
bool makeCopy)
489 this->env_->localInputAssertion(__FILE__, __LINE__,
"invalid length",
492 this->env_->localInputAssertion(__FILE__, __LINE__,
"invalid criteria",
496 partIds_[criteria] = ArrayRCP<part_t>();
497 partSizes_[criteria] = ArrayRCP<scalar_t>();
501 this->env_->localInputAssertion(__FILE__, __LINE__,
"invalid arrays",
508 part_t *z2_partIds = NULL;
510 bool own_memory =
false;
513 z2_partIds =
new part_t [len];
515 this->env_->localMemoryAssertion(__FILE__, __LINE__, len, z2_partSizes);
516 memcpy(z2_partIds, partIds, len *
sizeof(
part_t));
517 memcpy(z2_partSizes, partSizes, len *
sizeof(
scalar_t));
521 z2_partIds = partIds;
522 z2_partSizes = partSizes;
525 partIds_[criteria] = arcp(z2_partIds, 0, len, own_memory);
526 partSizes_[criteria] = arcp(z2_partSizes, 0, len, own_memory);
529 template <
typename Adapter>
533 if (algName_ == std::string(
"multijagged")) {
536 this->baseInputAdapter_));
538 else if (algName_ == std::string(
"zoltan")) {
541 this->baseInputAdapter_));
543 else if (algName_ == std::string(
"parma")) {
546 this->baseInputAdapter_));
548 else if (algName_ == std::string(
"scotch")) {
551 this->baseInputAdapter_));
553 else if (algName_ == std::string(
"parmetis")) {
557 this->baseInputAdapter_,
560 else if (algName_ == std::string(
"quotient")) {
563 this->baseInputAdapter_,
567 else if (algName_ == std::string(
"pulp")) {
570 this->baseInputAdapter_));
572 else if (algName_ == std::string(
"block")) {
574 this->comm_, this->baseInputAdapter_));
576 else if (algName_ == std::string(
"phg") ||
577 algName_ == std::string(
"patoh")) {
579 Teuchos::ParameterList &pl = this->env_->getParametersNonConst();
580 Teuchos::ParameterList &zparams = pl.sublist(
"zoltan_parameters",
false);
581 if (numberOfWeights_ > 0) {
583 sprintf(strval,
"%d", numberOfWeights_);
584 zparams.set(
"OBJ_WEIGHT_DIM", strval);
586 zparams.set(
"LB_METHOD", algName_.c_str());
587 zparams.set(
"LB_APPROACH",
"PARTITION");
588 algName_ = std::string(
"zoltan");
592 this->baseInputAdapter_));
594 else if (algName_ == std::string(
"sarma")) {
597 this->baseInputAdapter_));
599 else if (algName_ == std::string(
"forTestingOnly")) {
602 this->baseInputAdapter_));
609 throw std::logic_error(
"partitioning algorithm not supported");
613 template <
typename Adapter>
622 createPartitioningProblem(updateInputData);
634 this->createAlgorithm();
638 this->env_->timerStart(
MACRO_TIMERS,
"create solution");
643 this->envConst_, this->comm_, numberOfWeights_,
644 partIds_.view(0, numberOfCriteria_),
645 partSizes_.view(0, numberOfCriteria_), this->algorithm_);
649 solution_ = rcp(soln);
652 this->env_->memory(
"After creating Solution");
657 this->algorithm_->partition(solution_);
662 const Teuchos::ParameterEntry *pe = this->envConst_->getParameters().getEntryPtr(
"mapping_type");
663 int mapping_type = -1;
665 mapping_type = pe->getValue(&mapping_type);
669 #if ZOLTAN2_TASKMAPPING_MOVE
670 if (mapping_type == 0){
676 this->comm_.getRawPtr(),
677 machine_.getRawPtr(),
678 this->baseInputAdapter_.getRawPtr(),
679 solution_.getRawPtr(),
680 this->envConst_.getRawPtr()
693 const part_t *oldParts = solution_->getPartListView();
694 size_t nLocal = ia->getNumLocalIds();
695 for (
size_t i = 0; i < nLocal; i++) {
707 else if (mapping_type == 1){
714 template <
typename Adapter>
716 const std::string &algorithm,
const std::string &defString,
717 const std::string &model,
Environment &env,
bool &removeSelfEdges,
718 bool &isGraphType,
bool &needConsecutiveGlobalIds) {
721 if (algorithm != defString) {
722 if (algorithm == std::string(
"block") ||
723 algorithm == std::string(
"random") ||
724 algorithm == std::string(
"cyclic") ||
725 algorithm == std::string(
"zoltan") ||
726 algorithm == std::string(
"parma") ||
727 algorithm == std::string(
"forTestingOnly") ||
728 algorithm == std::string(
"quotient") ||
729 algorithm == std::string(
"scotch") ||
730 algorithm == std::string(
"ptscotch") ||
731 algorithm == std::string(
"pulp") || algorithm == std::string(
"sarma") ||
732 algorithm == std::string(
"patoh") || algorithm == std::string(
"phg") ||
733 algorithm == std::string(
"multijagged")) {
734 algName_ = algorithm;
735 }
else if (algorithm == std::string(
"rcb") ||
736 algorithm == std::string(
"rib") ||
737 algorithm == std::string(
"hsfc")) {
739 Teuchos::ParameterList &zparams = pl.sublist(
"zoltan_parameters",
false);
740 zparams.set(
"LB_METHOD", algorithm);
741 if (numberOfWeights_ > 0) {
743 sprintf(strval,
"%d", numberOfWeights_);
744 zparams.set(
"OBJ_WEIGHT_DIM", strval);
746 algName_ = std::string(
"zoltan");
747 }
else if (algorithm == std::string(
"metis") ||
748 algorithm == std::string(
"parmetis")) {
749 algName_ = algorithm;
750 removeSelfEdges =
true;
751 needConsecutiveGlobalIds =
true;
755 throw std::logic_error(
"parameter list algorithm is invalid");
757 }
else if (model != defString) {
759 if (model == std::string(
"hypergraph")) {
760 algName_ = std::string(
"phg");
761 }
else if (model == std::string(
"graph")) {
762 #ifdef HAVE_ZOLTAN2_SCOTCH
763 if (this->comm_->getSize() > 1)
764 algName_ = std::string(
"ptscotch");
766 algName_ = std::string(
"scotch");
768 #ifdef HAVE_ZOLTAN2_PARMETIS
769 if (this->comm_->getSize() > 1)
770 algName_ = std::string(
"parmetis");
772 algName_ = std::string(
"metis");
773 removeSelfEdges =
true;
774 needConsecutiveGlobalIds =
true;
777 #ifdef HAVE_ZOLTAN2_PULP
782 algName_ = std::string(
"pulp");
784 algName_ = std::string(
"phg");
788 }
else if (model == std::string(
"geometry")) {
789 algName_ = std::string(
"multijagged");
790 }
else if (model == std::string(
"ids")) {
791 algName_ = std::string(
"block");
795 "parameter list model type is invalid", 1,
804 algName_ = std::string(
"phg");
807 algName_ = std::string(
"phg");
809 algName_ = std::string(
"multijagged");
811 algName_ = std::string(
"block");
814 throw std::logic_error(
"input type is invalid");
819 template <
typename Adapter>
823 "PartitioningProblem::createPartitioningProblem");
825 using Teuchos::ParameterList;
854 std::string defString(
"default");
858 std::string model(defString);
859 const Teuchos::ParameterEntry *pe = pl.getEntryPtr(
"model");
861 model = pe->getValue<std::string>(&model);
865 std::string algorithm(defString);
866 pe = pl.getEntryPtr(
"algorithm");
868 algorithm = pe->getValue<std::string>(&algorithm);
872 bool needConsecutiveGlobalIds =
false;
873 bool removeSelfEdges=
false;
874 bool isGraphModel =
false;
880 this->processAlgorithmName(algorithm, defString, model, env, removeSelfEdges,
881 isGraphModel, needConsecutiveGlobalIds);
885 Array<int> valueList;
886 pe = pl.getEntryPtr(
"topology");
889 valueList = pe->getValue<Array<int> >(&valueList);
891 if (!Zoltan2::noValuesAreInRangeList<int>(valueList)){
892 int *n =
new int [valueList.size() + 1];
893 levelNumberParts_ = arcp(n, 0, valueList.size() + 1,
true);
894 int procsPerNode = 1;
895 for (
int i=0; i < valueList.size(); i++){
896 levelNumberParts_[i+1] = valueList[i];
897 procsPerNode *= valueList[i];
900 levelNumberParts_[0] = env.
numProcs_ / procsPerNode;
903 levelNumberParts_[0]++;
907 levelNumberParts_.clear();
910 hierarchical_ = levelNumberParts_.size() > 0;
914 std::string objectOfInterest(defString);
915 pe = pl.getEntryPtr(
"objects_to_partition");
917 objectOfInterest = pe->getValue<std::string>(&objectOfInterest);
924 if (isGraphModel ==
true)
929 std::string symParameter(defString);
930 pe = pl.getEntryPtr(
"symmetrize_graph");
932 symParameter = pe->getValue<std::string>(&symParameter);
933 if (symParameter == std::string(
"transpose"))
935 else if (symParameter == std::string(
"bipartite"))
939 bool sgParameter =
false;
940 pe = pl.getEntryPtr(
"subset_graph");
942 sgParameter = pe->getValue(&sgParameter);
944 if (sgParameter == 1)
952 if (needConsecutiveGlobalIds)
958 if (objectOfInterest == defString ||
959 objectOfInterest == std::string(
"matrix_rows") )
961 else if (objectOfInterest == std::string(
"matrix_columns"))
963 else if (objectOfInterest == std::string(
"matrix_nonzeros"))
968 if (objectOfInterest == defString ||
969 objectOfInterest == std::string(
"mesh_nodes") )
971 else if (objectOfInterest == std::string(
"mesh_elements"))
use columns as graph vertices
Serial greedy first-fit graph coloring (local graph only)
~PartitioningProblem()
Destructor.
Zoltan2::BaseAdapter< userTypes_t > base_adapter_t
Time an algorithm (or other entity) as a whole.
Adapter::scalar_t scalar_t
BaseAdapterType inputType_
void setPartSizes(int len, part_t *partIds, scalar_t *partSizes, bool makeCopy=true)
Set or reset relative sizes for the parts that Zoltan2 will create.
use mesh nodes as vertices
ArrayRCP< ArrayRCP< part_t > > partIds_
fast typical checks for valid arguments
ArrayRCP< ArrayRCP< scalar_t > > partSizes_
#define Z2_FORWARD_EXCEPTIONS
Forward an exception back through call stack.
static void getValidParameters(ParameterList &pl)
Set up validators specific to this algorithm.
static RCP< Teuchos::BoolParameterEntryValidator > getBoolValidator()
Exists to make setting up validators less cluttered.
ArrayRCP< int > levelNumberParts_
algorithm requires consecutive ids
std::bitset< NUM_MODEL_FLAGS > modelFlag_t
PartitioningProblem(Adapter *A, ParameterList *p)
Constructor where communicator is the Teuchos default.
model must symmetrize input
map_t::global_ordinal_type gno_t
model must symmetrize input
PartitioningProblem(Adapter *A, ParameterList *p, const RCP< const Teuchos::Comm< int > > &comm)
Constructor where Teuchos communicator is specified.
Defines the PartitioningSolution class.
static RCP< Teuchos::AnyNumberParameterEntryValidator > getAnyDoubleValidator()
Exists to make setting up validators less cluttered.
part_t getAssignedProcForTask(part_t taskId)
getAssignedProcForTask function, returns the assigned processor id for the given task ...
sub-steps, each method's entry and exit
static void getValidParameters(ParameterList &pl)
Set up validators specific to this Problem.
void setPartSizesForCriteria(int criteria, int len, part_t *partIds, scalar_t *partSizes, bool makeCopy=true)
Set or reset the relative sizes (per weight) for the parts that Zoltan2 will create.
int numProcs_
number of processes (relative to comm_)
SparseMatrixAdapter_t::part_t part_t
void createPartitioningProblem(bool newData)
use matrix rows as graph vertices
Teuchos::ParameterList & getParametersNonConst()
Returns a reference to a non-const copy of the parameters.
algorithm requires no self edges
Defines the IdentifierModel interface.
virtual void createAlgorithm()
A PartitioningSolution is a solution to a partitioning problem.
use nonzeros as graph vertices
Defines the EvaluatePartition class.
virtual void processAlgorithmName(const std::string &algorithm, const std::string &defString, const std::string &model, Environment &env, bool &removeSelfEdges, bool &isGraphType, bool &needConsecutiveGlobalIds)
BaseAdapterType
An enum to identify general types of adapters.
identifier data, just a list of IDs
void localBugAssertion(const char *file, int lineNum, const char *msg, bool ok, AssertionLevel level) const
Test for valid library behavior on local process only.
Problem base class from which other classes (PartitioningProblem, ColoringProblem, OrderingProblem, MatchingProblem, etc.) derive.
Defines the Problem base class.
map_t::local_ordinal_type lno_t
RCP< PartitioningSolution< Adapter > > solution_
The user parameters, debug, timing and memory profiling output objects, and error checking methods...
static RCP< Teuchos::AnyNumberParameterEntryValidator > getAnyIntValidator()
Exists to make setting up validators less cluttered.
static void getValidParameters(ParameterList &pl)
Set up validators specific to this algorithm.
MachineRepresentation Class Base class for representing machine coordinates, networks, etc.
Adapter::base_adapter_t base_adapter_t
Define IntegerRangeList validator.
const PartitioningSolution< Adapter > & getSolution()
Get the solution to the problem.
PartitioningProblem sets up partitioning problems for the user.
use mesh elements as vertices
GraphModel defines the interface required for graph models.
Defines the GraphModel interface.
static void getValidParameters(ParameterList &pl)
Set up validators specific to this algorithm.
static void getValidParameters(ParameterList &pl)
Set up validators specific to this algorithm.
static void getValidParameters(ParameterList &)
Set up validators specific to this algorithm.
void solve(bool updateInputData=true)
Direct the problem to create a solution.
Multi Jagged coordinate partitioning algorithm.
static void getValidParameters(ParameterList &pl)
Set up validators specific to this algorithm.
Zoltan2::BasicUserTypes< zscalar_t, zlno_t, zgno_t > user_t