54 #include <Tpetra_Map.hpp>
62 int main(
int argc,
char *argv[])
64 Tpetra::ScopeGuard tscope(&argc, &argv);
65 Teuchos::RCP<const Teuchos::Comm<int> > comm = Tpetra::getDefaultComm();
67 int rank = comm->getRank();
68 int nprocs = comm->getSize();
72 typedef Tpetra::Map<> Map_t;
73 typedef Map_t::local_ordinal_type localId_t;
74 typedef Map_t::global_ordinal_type globalId_t;
76 typedef Tpetra::Details::DefaultTypes::scalar_type scalar_t;
87 size_t localCount = 40;
90 scalar_t *coords =
new scalar_t [dim * localCount];
93 scalar_t *y = x + localCount;
94 scalar_t *z = y + localCount;
99 scalar_t scalingFactor = 10.0 / RAND_MAX;
101 for (
size_t i=0; i < localCount*dim; i++){
102 coords[i] = scalar_t(rand()) * scalingFactor;
107 globalId_t *globalIds =
new globalId_t [localCount];
108 globalId_t offset = rank * localCount;
110 for (
size_t i=0; i < localCount; i++)
111 globalIds[i] = offset++;
116 double tolerance = 1.1;
119 std::cout <<
"Imbalance tolerance is " << tolerance << std::endl;
121 Teuchos::ParameterList params(
"test params");
122 params.set(
"debug_level",
"basic_status");
123 params.set(
"debug_procs",
"0");
124 params.set(
"error_check_level",
"debug_mode_assertions");
126 params.set(
"algorithm",
"rcb");
127 params.set(
"imbalance_tolerance", tolerance );
128 params.set(
"num_global_parts", nprocs);
138 inputAdapter_t *ia1 =
new inputAdapter_t(localCount,globalIds,x,y,z,1,1,1);
151 quality_t *metricObject1 =
new quality_t(ia1, ¶ms,
161 if (imb <= tolerance)
162 std::cout <<
"pass: " << imb << std::endl;
164 std::cout <<
"fail: " << imb << std::endl;
165 std::cout << std::endl;
167 delete metricObject1;
175 scalar_t *
weights =
new scalar_t [localCount];
176 for (
size_t i=0; i < localCount; i++){
177 weights[i] = 1.0 + scalar_t(rank) / scalar_t(nprocs);
182 std::vector<const scalar_t *>coordVec(2);
183 std::vector<int> coordStrides(2);
185 coordVec[0] = x; coordStrides[0] = 1;
186 coordVec[1] = y; coordStrides[1] = 1;
188 std::vector<const scalar_t *>weightVec(1);
189 std::vector<int> weightStrides(1);
191 weightVec[0] =
weights; weightStrides[0] = 1;
193 inputAdapter_t *ia2=
new inputAdapter_t(localCount, globalIds, coordVec,
194 coordStrides,weightVec,weightStrides);
207 #ifdef HAVE_ZOLTAN2_MPI
208 quality_t *metricObject2 =
new quality_t(ia2, ¶ms,
218 metricObject2->printMetrics(std::cout);
222 scalar_t imb = metricObject2->getWeightImbalance(0);
223 if (imb <= tolerance)
224 std::cout <<
"pass: " << imb << std::endl;
226 std::cout <<
"fail: " << imb << std::endl;
227 std::cout << std::endl;
229 delete metricObject2;
244 params.set(
"partitioning_objective",
"multicriteria_minimize_total_weight");
248 weights =
new scalar_t [localCount*3];
251 for (
size_t i=0; i < localCount*3; i+=3){
252 weights[i] = 1.0 + rank / nprocs;
253 weights[i+1] = rank<nprocs/2 ? 1 : 2;
254 weights[i+2] = rand()/RAND_MAX +.5;
260 weightStrides.resize(3);
262 weightVec[0] =
weights; weightStrides[0] = 3;
263 weightVec[1] = weights+1; weightStrides[1] = 3;
264 weightVec[2] = weights+2; weightStrides[2] = 3;
266 inputAdapter_t *ia3=
new inputAdapter_t(localCount, globalIds, coordVec,
267 coordStrides,weightVec,weightStrides);
290 if (imb <= tolerance)
291 std::cout <<
"pass: " << imb << std::endl;
293 std::cout <<
"fail: " << imb << std::endl;
294 std::cout << std::endl;
296 delete metricObject3;
301 bool dataHasChanged =
false;
303 params.set(
"partitioning_objective",
"multicriteria_minimize_maximum_weight");
305 problem3->
solve(dataHasChanged);
314 if (imb <= tolerance)
315 std::cout <<
"pass: " << imb << std::endl;
317 std::cout <<
"fail: " << imb << std::endl;
318 std::cout << std::endl;
320 delete metricObject3;
322 params.set(
"partitioning_objective",
"multicriteria_balance_total_maximum");
324 problem3->
solve(dataHasChanged);
333 if (imb <= tolerance)
334 std::cout <<
"pass: " << imb << std::endl;
336 std::cout <<
"fail: " << imb << std::endl;
337 std::cout << std::endl;
339 delete metricObject3;
355 params.set(
"num_global_parts", nprocs*2);
363 scalar_t partSizes[2];
365 partIds[0] = rank*2; partSizes[0] = 0;
366 partIds[1] = rank*2+1; partSizes[1] = 1;
374 dataHasChanged =
false;
376 problem1->
solve(dataHasChanged);
387 int numInEmptyParts = 0;
388 for (
size_t i=0; i < localCount; i++){
389 if (partAssignments[i] % 2 == 0)
394 std::cout <<
"Request that " << nprocs <<
" parts be empty." <<std::endl;
398 metricObject1 =
new quality_t(ia1, ¶ms,
408 if (imb <= tolerance)
409 std::cout <<
"pass: " << imb << std::endl;
411 std::cout <<
"fail: " << imb << std::endl;
412 std::cout << std::endl;
414 delete metricObject1;
430 std::cout <<
"PASS" << std::endl;
void printMetrics(std::ostream &os) const
Print all metrics.
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.
A simple class that can be the User template argument for an InputAdapter.
Zoltan2::EvaluatePartition< matrixAdapter_t > quality_t
int main(int narg, char **arg)
void resetParameters(ParameterList *params)
Reset the list of parameters.
Defines the PartitioningSolution class.
scalar_t getWeightImbalance(int weightIndex) const
Return the imbalance for the requested weight.
SparseMatrixAdapter_t::part_t part_t
A PartitioningSolution is a solution to a partitioning problem.
const part_t * getPartListView() const
Returns the part list corresponding to the global ID list.
RCP< const Comm< int > > getComm()
Return the communicator used by the problem.
BasicVectorAdapter represents a vector (plus optional weights) supplied by the user as pointers to st...
const PartitioningSolution< Adapter > & getSolution()
Get the solution to the problem.
PartitioningProblem sets up partitioning problems for the user.
Defines the PartitioningProblem class.
A class that computes and returns quality metrics.
Defines the BasicVectorAdapter class.
void solve(bool updateInputData=true)
Direct the problem to create a solution.
scalar_t getObjectCountImbalance() const
Return the object count imbalance.