Zoltan2
Loading...
Searching...
No Matches
partitioningTree.cpp
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
49#include <iostream>
50#include <limits>
51#include <Teuchos_ParameterList.hpp>
52#include <Teuchos_RCP.hpp>
53#include <Teuchos_FancyOStream.hpp>
54#include <Teuchos_CommandLineProcessor.hpp>
55#include <Tpetra_CrsMatrix.hpp>
56#include <Tpetra_Vector.hpp>
57#include <MatrixMarket_Tpetra.hpp>
58
59using Teuchos::RCP;
60
64
65typedef Tpetra::CrsMatrix<z2TestScalar, z2TestLO, z2TestGO> SparseMatrix_t;
66typedef Tpetra::Vector<z2TestScalar, z2TestLO, z2TestGO> Vector;
67typedef Vector::node_type Node;
68
69typedef Tpetra::MultiVector<z2TestScalar, z2TestLO, z2TestGO,znode_t> tMVector_t;
70
71
73
75
77
78
79
80int testForRCB(SparseMatrixAdapter_t &matAdapter, int myrank, part_t numparts,
81 RCP<tMVector_t> coords, RCP<const Teuchos::Comm<int> > comm);
82int testForPHG(SparseMatrixAdapter_t &matAdapter, int myrank, part_t numparts,
83 RCP<tMVector_t> coords, RCP<const Teuchos::Comm<int> >);
84int testForMJ(SparseMatrixAdapter_t &matAdapter, int myrank, part_t numparts,
85 RCP<tMVector_t> coords, RCP<const Teuchos::Comm<int> >);
86
87
90int main(int narg, char** arg)
91{
92 Tpetra::ScopeGuard tscope(&narg, &arg);
93 Teuchos::RCP<const Teuchos::Comm<int> > comm = Tpetra::getDefaultComm();
94 int me = comm->getRank();
95
97
98 std::string inputFile = ""; // Matrix Market or Zoltan file to read
99 std::string inputPath = testDataFilePath; // Directory with input file
100 bool distributeInput = true;
101 int success = 0;
102 part_t numParts = 8;
103
104
106 // Read run-time options.
108 Teuchos::CommandLineProcessor cmdp (false, false);
109 cmdp.setOption("inputPath", &inputPath,
110 "Path to the MatrixMarket or Zoltan file to be read; "
111 "if not specified, a default path will be used.");
112 cmdp.setOption("inputFile", &inputFile,
113 "Name of the Matrix Market or Zoltan file to read; "
114 "");
115 cmdp.setOption("distribute", "no-distribute", &distributeInput,
116 "for Zoltan input files only, "
117 "indicate whether or not to distribute "
118 "input across the communicator");
119 cmdp.setOption("numParts", &numParts,
120 "Global number of parts;");
121
122 Teuchos::CommandLineProcessor::EParseCommandLineReturn
123 parseReturn= cmdp.parse( narg, arg );
124
125 if( parseReturn == Teuchos::CommandLineProcessor::PARSE_HELP_PRINTED )
126 {
127 return 0;
128 }
130
132 // Construct matrix from file
134 RCP<UserInputForTests> uinput;
135
136 if (inputFile != "") // Input file specified; read a matrix
137 {
138 uinput = rcp(new UserInputForTests(inputPath, inputFile, comm,
139 true, distributeInput));
140 }
141 else
142 {
143 std::cout << "Input file must be specified." << std::endl;
144 }
145
146 RCP<SparseMatrix_t> origMatrix = uinput->getUITpetraCrsMatrix();
147
148
149 if (me == 0)
150 {
151 std::cout << "NumRows = " << origMatrix->getGlobalNumRows() << std::endl
152 << "NumNonzeros = " << origMatrix->getGlobalNumEntries() << std::endl
153 << "NumProcs = " << comm->getSize() << std::endl
154 << "NumParts = " << numParts << std::endl;
155 }
156
157 if (origMatrix->getGlobalNumRows() < 40)
158 {
159 Teuchos::FancyOStream out(Teuchos::rcp(&std::cout,false));
160 origMatrix->describe(out, Teuchos::VERB_EXTREME);
161 }
163
165 // Get coordinates, corresponding to graph vertices
167 RCP<tMVector_t> coords;
168 try
169 {
170 coords = uinput->getUICoordinates();
171 }
172 catch(...)
173 {
174 if (me == 0)
175 std::cout << "FAIL: get coordinates" << std::endl;
176 return 1;
177 }
179
183 SparseMatrixAdapter_t matAdapter(origMatrix);
184
185 MultiVectorAdapter_t *ca = NULL;
186
187 try
188 {
189 ca = new MultiVectorAdapter_t(coords);
190 }
191 catch(...)
192 {
193 if (me == 0)
194 std::cout << "FAIL: vector adapter" << std::endl;
195 return 1;
196 }
197
198 matAdapter.setCoordinateInput(ca);
200
201 // MDM - MJ disabled - not implemented yet
202 // int errMJ = testForMJ(matAdapter, me, numParts, coords, comm); -check err
203 int errRCB = testForRCB(matAdapter, me, numParts, coords, comm);
204 int errPHG = testForPHG(matAdapter, me, numParts, coords, comm);
205
206 // Currently just for internal development - sweep from 1 - 25 numParts
207 // and validate each to check for any bad results.
208 // #define BRUTE_FORCE_SEARCH
209 #ifdef BRUTE_FORCE_SEARCH
210 for(int setParts = 1; setParts <= 25; ++setParts) {
211 // numParts is a parameter so should not be set like for for production
212 // code - this is just a 2nd development check to run try a bunch of tests
213 numParts = setParts;
214 std::cout << "Brute force testing num parts: " << numParts << std::endl;
215 // MJ not implemented yet
216 // if(testForMJ(matAdapter, me, numParts, coords, comm) != 0) { return 1; }
217 if(testForRCB(matAdapter, me, numParts, coords, comm) != 0) { return 1; }
218 if(testForPHG(matAdapter, me, numParts, coords, comm) != 0) { return 1; }
219 }
220 #endif
221
222 delete ca;
223
224 // if(errMJ==0)
225 if(errRCB==0)
226 if(errPHG==0)
227 {
228 std::cout << "PASS" << std::endl;
229 return success;
230 }
231 return 1;
232}
234
235
237// Runs validation on the results to make sure the arrays are sensible
239bool validate(part_t numTreeVerts,
240 std::vector<part_t> permPartNums,
241 std::vector<part_t> splitRangeBeg,
242 std::vector<part_t> splitRangeEnd,
243 std::vector<part_t> treeVertParents
244)
245{
246 // by design numTreeVerts does not include the root
247 if(numTreeVerts != static_cast<part_t>(splitRangeBeg.size()) - 1) {
248 return false;
249 }
250 if(numTreeVerts != static_cast<part_t>(splitRangeEnd.size()) - 1) {
251 return false;
252 }
253 if(numTreeVerts != static_cast<part_t>(treeVertParents.size()) - 1) {
254 return false;
255 }
256 // now search every node and validate it's properties
257 for(part_t n = 0; n <= numTreeVerts; ++n) {
258 if(n < static_cast<part_t>(permPartNums.size())) {
259 // terminal so children should just be range 1
260 if(splitRangeEnd[n] != splitRangeBeg[n] + 1) {
261 std::cout << "Invalid terminal - range should be 1" << std::endl;
262 return false;
263 }
264 if(splitRangeBeg[n] != n) {
265 std::cout << "Invalid terminal - not pointing to myself!" << std::endl;
266 return false;
267 }
268 }
269 else {
270 part_t beg = splitRangeBeg[n];
271 part_t end = splitRangeEnd[n];
272 part_t count = end - beg;
273 std::vector<bool> findChildren(count, false);
274 for(part_t n2 = 0; n2 <= numTreeVerts; ++n2) {
275 if(treeVertParents[n2] == n) {
276 part_t beg2 = splitRangeBeg[n2];
277 part_t end2 = splitRangeEnd[n2];
278 for(part_t q = beg2; q < end2; ++q) {
279 part_t setIndex = q - beg; // rel to parent so beg, not beg2
280 if(setIndex < 0 || setIndex >= count) {
281 std::cout << "Found bad child index - invalid results!" << std::endl;
282 return false;
283 }
284 if(findChildren[setIndex]) {
285 std::cout << "Found child twice - invalid results!" << std::endl;
286 return false;
287 }
288 findChildren[setIndex] = true;
289 }
290 }
291 }
292 for(part_t q = 0; q < count; ++q) {
293 if(!findChildren[q]) {
294 std::cout << "Did not find all children. Invalid results!" << std::endl;
295 return false;
296 }
297 }
298 }
299 if(n == numTreeVerts) {
300 // this is the root
301 if(splitRangeBeg[n] != 0) {
302 std::cout << "Root must start at 0!" << std::endl;
303 return false;
304 }
305 if(splitRangeEnd[n] != static_cast<part_t>(permPartNums.size())) {
306 std::cout << "Root must contain all parts!" << std::endl;
307 return false;
308 }
309 }
310 }
311 return true;
312}
314
316// General test function to analyze solution and print out some information
319 RCP<const Teuchos::Comm<int> > comm)
320{
321 part_t numTreeVerts = 0;
322 std::vector<part_t> permPartNums; // peritab in Scotch
323 std::vector<part_t> splitRangeBeg;
324 std::vector<part_t> splitRangeEnd;
325 std::vector<part_t> treeVertParents;
326
328 problem.getSolution();
329
330 solution.getPartitionTree(numTreeVerts,permPartNums,splitRangeBeg,
331 splitRangeEnd,treeVertParents);
332
333 comm->barrier(); // for tidy output...
334 if(comm->getRank() == 0) { // for now just plot for rank 0
335
336 // print the acquired information about the tree
337
338 // Header
339 std::cout << std::endl << "Printing partition tree info..." << std::endl << std::endl;
340
341 // numTreeVerts
342 std::cout << " numTreeVerts: " << numTreeVerts << std::endl << std::endl;
343
344 // array index values 0 1 2 3 ...
345 std::cout << " part array index:";
346 for(part_t n = 0; n < static_cast<part_t>(permPartNums.size()); ++n) {
347 std::cout << std::setw(4) << n << " ";
348 }
349 std::cout << std::endl;
350
351 // permParNums
352 std::cout << " permPartNums: ";
353 for(part_t n = 0; n < static_cast<part_t>(permPartNums.size()); ++n) {
354 std::cout << std::setw(4) << permPartNums[n] << " ";
355 }
356 std::cout << std::endl << std::endl;
357
358 // node index values 0 1 2 3 ...
359 std::cout << " node index: ";
360 for(part_t n = 0; n < static_cast<part_t>(splitRangeBeg.size()); ++n) {
361 std::cout << std::setw(4) << n << " ";
362 }
363 std::cout << std::endl;
364
365 // splitRangeBeg
366 std::cout << " splitRangeBeg: ";
367 for(part_t n = 0; n < static_cast<part_t>(splitRangeBeg.size()); ++n) {
368 std::cout << std::setw(4) << splitRangeBeg[n] << " ";
369 }
370 std::cout << std::endl;
371
372 // splitRangeEnd
373 std::cout << " splitRangeEnd: ";
374 for(part_t n = 0; n < static_cast<part_t>(splitRangeEnd.size()); ++n) {
375 std::cout << std::setw(4) << splitRangeEnd[n] << " ";
376 }
377 std::cout << std::endl;
378
379 // treeVertParents
380 std::cout << " treeVertParents: ";
381 for(part_t n = 0; n < static_cast<part_t>(treeVertParents.size()); ++n) {
382 std::cout << std::setw(4) << treeVertParents[n] << " ";
383 }
384 std::cout << std::endl << std::endl;
385 }
386 comm->barrier(); // for tidy output...
387
388 if(!validate(numTreeVerts, permPartNums, splitRangeBeg, splitRangeEnd,
389 treeVertParents)) {
390 return 1;
391 }
392
393 return 0;
394}
395
397// Test partitioning tree access for RCB partitioning. This partitioning tree
398// should be binary..
400int testForRCB(SparseMatrixAdapter_t &matAdapter, int me, part_t numParts,
401 RCP<tMVector_t> coords, RCP<const Teuchos::Comm<int> > comm)
402{
403
407 Teuchos::ParameterList params;
408
409 params.set("num_global_parts", numParts);
410 params.set("partitioning_approach", "partition");
411 params.set("algorithm", "rcb");
412 params.set("keep_partition_tree", true);
413
415
419 Zoltan2::PartitioningProblem<SparseMatrixAdapter_t> problem(&matAdapter, &params);
420
421 try
422 {
423 if (me == 0) std::cout << "Calling solve() " << std::endl;
424 problem.solve();
425 if (me == 0) std::cout << "Done solve() " << std::endl;
426 }
427 catch (std::runtime_error &e)
428 {
429 std::cout << "Runtime exception returned from solve(): " << e.what();
430 if (!strncmp(e.what(), "BUILD ERROR", 11)) {
431 // Catching build errors as exceptions is OK in the tests
432 std::cout << " PASS" << std::endl;
433 return 0;
434 }
435 else {
436 // All other runtime_errors are failures
437 std::cout << " FAIL" << std::endl;
438 return -1;
439 }
440 }
441 catch (std::logic_error &e)
442 {
443 std::cout << "Logic exception returned from solve(): " << e.what()
444 << " FAIL" << std::endl;
445 return -1;
446 }
447 catch (std::bad_alloc &e)
448 {
449 std::cout << "Bad_alloc exception returned from solve(): " << e.what()
450 << " FAIL" << std::endl;
451 return -1;
452 }
453 catch (std::exception &e)
454 {
455 std::cout << "Unknown exception returned from solve(). " << e.what()
456 << " FAIL" << std::endl;
457 return -1;
458 }
460
461
463
464 bool binary = solution.isPartitioningTreeBinary();
465
466 if (binary == false)
467 {
468 std::cout << "RCB should produce a binary partitioning tree. FAIL" << std::endl;
469 return -1;
470 }
471
472 return analyze(problem, comm);
473}
475
477// Test partitioning tree access for PHG partitioning. This partitioning tree
478// should be balanced.
480int testForPHG(SparseMatrixAdapter_t &matAdapter, int me, part_t numParts,
481 RCP<tMVector_t> coords, RCP<const Teuchos::Comm<int> > comm)
482{
483
487 Teuchos::ParameterList params;
488
489 params.set("num_global_parts", numParts);
490 params.set("partitioning_approach", "partition");
491 params.set("algorithm", "zoltan");
492 params.set("keep_partition_tree", true);
493
494 Teuchos::ParameterList &zparams = params.sublist("zoltan_parameters",false);
495 zparams.set("LB_METHOD","phg");
496 zparams.set("FINAL_OUTPUT", "1");
497
499
503 Zoltan2::PartitioningProblem<SparseMatrixAdapter_t> problem(&matAdapter, &params);
504
505 try
506 {
507 if (me == 0) std::cout << "Calling solve() " << std::endl;
508 problem.solve();
509 if (me == 0) std::cout << "Done solve() " << std::endl;
510 }
511 catch (std::runtime_error &e)
512 {
513 std::cout << "Runtime exception returned from solve(): " << e.what();
514 if (!strncmp(e.what(), "BUILD ERROR", 11)) {
515 // Catching build errors as exceptions is OK in the tests
516 std::cout << " PASS" << std::endl;
517 return 0;
518 }
519 else {
520 // All other runtime_errors are failures
521 std::cout << " FAIL" << std::endl;
522 return -1;
523 }
524 }
525 catch (std::logic_error &e)
526 {
527 std::cout << "Logic exception returned from solve(): " << e.what()
528 << " FAIL" << std::endl;
529 return -1;
530 }
531 catch (std::bad_alloc &e)
532 {
533 std::cout << "Bad_alloc exception returned from solve(): " << e.what()
534 << " FAIL" << std::endl;
535 return -1;
536 }
537 catch (std::exception &e)
538 {
539 std::cout << "Unknown exception returned from solve(). " << e.what()
540 << " FAIL" << std::endl;
541 return -1;
542 }
544
546
547 bool binary = solution.isPartitioningTreeBinary();
548
549 if (binary == false)
550 {
551 std::cout << "PHG should produce a binary partitioning tree. FAIL" << std::endl;
552 return -1;
553 }
554
555 return analyze(problem, comm);
556}
558
560// Test partitioning tree access for MJ partitioning. This partitioning tree
561// should be balanced.
563int testForMJ(SparseMatrixAdapter_t &matAdapter, int me, part_t numParts,
564 RCP<tMVector_t> coords, RCP<const Teuchos::Comm<int> > comm)
565{
566
570 Teuchos::ParameterList params;
571
572 params.set("num_global_parts", numParts);
573 params.set("algorithm", "multijagged");
574 params.set("rectilinear", true);
575 // params.set("mj_keep_part_boxes", true); // allows getPartBoxesView on solution
576 params.set("keep_partition_tree", true);
577
579
583 Zoltan2::PartitioningProblem<SparseMatrixAdapter_t> problem(&matAdapter, &params);
584
585 try
586 {
587 if (me == 0) std::cout << "Calling solve() " << std::endl;
588 problem.solve();
589 if (me == 0) std::cout << "Done solve() " << std::endl;
590 }
591 catch (std::runtime_error &e)
592 {
593 std::cout << "Runtime exception returned from solve(): " << e.what();
594 if (!strncmp(e.what(), "BUILD ERROR", 11)) {
595 // Catching build errors as exceptions is OK in the tests
596 std::cout << " PASS" << std::endl;
597 return 0;
598 }
599 else {
600 // All other runtime_errors are failures
601 std::cout << " FAIL" << std::endl;
602 return -1;
603 }
604 }
605 catch (std::logic_error &e)
606 {
607 std::cout << "Logic exception returned from solve(): " << e.what()
608 << " FAIL" << std::endl;
609 return -1;
610 }
611 catch (std::bad_alloc &e)
612 {
613 std::cout << "Bad_alloc exception returned from solve(): " << e.what()
614 << " FAIL" << std::endl;
615 return -1;
616 }
617 catch (std::exception &e)
618 {
619 std::cout << "Unknown exception returned from solve(). " << e.what()
620 << " FAIL" << std::endl;
621 return -1;
622 }
624
626
627 // copied from the MultiJaggedTest.cpp to inspect boxes
628 // To activate set mj_keep_part_boxes true above
629 /*
630 std::vector<Zoltan2::coordinateModelPartBox<scalar_t, part_t> >
631 & pBoxes = solution.getPartBoxesView();
632 int coordDim = coords->getNumVectors();
633 std::cout << std::endl;
634 std::cout << "Plot final boxes..." << std::endl;
635 for (size_t i = 0; i < pBoxes.size(); i++) {
636 zscalar_t *lmin = pBoxes[i].getlmins();
637 zscalar_t *lmax = pBoxes[i].getlmaxs();;
638 int dim = pBoxes[i].getDim();
639 std::set<int> * pNeighbors = pBoxes[i].getNeighbors();
640 std::vector<int> * pGridIndices = pBoxes[i].getGridIndices();
641
642 std::cout << me << " pBox " << i << " pid " << pBoxes[i].getpId()
643 << " dim: " << dim
644 << " (" << lmin[0] << "," << lmin[1] << ") "
645 << "x"
646 << " (" << lmax[0] << "," << lmax[1] << ")" << std::endl;
647 }
648 */
649
650 bool binary = solution.isPartitioningTreeBinary();
651
652 if (binary == true)
653 {
654 std::cout << "MJ should not produce a binary partitioning tree for this problem. FAIL" << std::endl;
655 return -1;
656 }
657
658 return analyze(problem, comm);
659}
661
662
663
664
Defines the PartitioningProblem class.
common code used by tests
float zscalar_t
Tpetra::Map ::local_ordinal_type zlno_t
std::string testDataFilePath(".")
Tpetra::Map ::global_ordinal_type zgno_t
Defines the XpetraCrsMatrixAdapter class.
Defines the XpetraMultiVectorAdapter.
int main()
InputTraits< User >::part_t part_t
void setCoordinateInput(VectorAdapter< UserCoord > *coordData) override
Allow user to provide additional data that contains coordinate info associated with the MatrixAdapter...
PartitioningProblem sets up partitioning problems for the user.
const PartitioningSolution< Adapter > & getSolution()
Get the solution to the problem.
void solve(bool updateInputData=true)
Direct the problem to create a solution.
A PartitioningSolution is a solution to a partitioning problem.
virtual bool isPartitioningTreeBinary() const
calculate if partition tree is binary.
void getPartitionTree(part_t &numTreeVerts, std::vector< part_t > &permPartNums, std::vector< part_t > &splitRangeBeg, std::vector< part_t > &splitRangeEnd, std::vector< part_t > &treeVertParents) const
get the partition tree - fill the relevant arrays
Provides access for Zoltan2 to Xpetra::CrsMatrix data.
Vector::node_type Node
Tpetra::MultiVector< z2TestScalar, z2TestLO, z2TestGO, znode_t > tMVector_t
zlno_t z2TestLO
Tpetra::CrsMatrix< z2TestScalar, z2TestLO, z2TestGO > SparseMatrix_t
int analyze(Zoltan2::PartitioningProblem< SparseMatrixAdapter_t > &problem, RCP< const Teuchos::Comm< int > > comm)
bool validate(part_t numTreeVerts, std::vector< part_t > permPartNums, std::vector< part_t > splitRangeBeg, std::vector< part_t > splitRangeEnd, std::vector< part_t > treeVertParents)
int testForRCB(SparseMatrixAdapter_t &matAdapter, int myrank, part_t numparts, RCP< tMVector_t > coords, RCP< const Teuchos::Comm< int > > comm)
Zoltan2::XpetraMultiVectorAdapter< tMVector_t > MultiVectorAdapter_t
Zoltan2::XpetraCrsMatrixAdapter< SparseMatrix_t, tMVector_t > SparseMatrixAdapter_t
int testForPHG(SparseMatrixAdapter_t &matAdapter, int myrank, part_t numparts, RCP< tMVector_t > coords, RCP< const Teuchos::Comm< int > >)
Tpetra::Vector< z2TestScalar, z2TestLO, z2TestGO > Vector
zgno_t z2TestGO
SparseMatrixAdapter_t::part_t part_t
zscalar_t z2TestScalar
int testForMJ(SparseMatrixAdapter_t &matAdapter, int myrank, part_t numparts, RCP< tMVector_t > coords, RCP< const Teuchos::Comm< int > >)