Zoltan2
Loading...
Searching...
No Matches
MultiJaggedTest.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
45
57
59
60#include "Teuchos_XMLParameterListHelpers.hpp"
61
62#include <Teuchos_LAPACK.hpp>
63#include <fstream>
64#include <string>
65using Teuchos::RCP;
66using Teuchos::rcp;
67
68
69//#define hopper_separate_test
70#ifdef hopper_separate_test
71#include "stdio.h"
72#endif
73#define CATCH_EXCEPTIONS_AND_RETURN(pp) \
74 catch (std::runtime_error &e) { \
75 std::cout << "Runtime exception returned from " << pp << ": " \
76 << e.what() << " FAIL" << std::endl; \
77 return -1; \
78 } \
79 catch (std::logic_error &e) { \
80 std::cout << "Logic exception returned from " << pp << ": " \
81 << e.what() << " FAIL" << std::endl; \
82 return -1; \
83 } \
84 catch (std::bad_alloc &e) { \
85 std::cout << "Bad_alloc exception returned from " << pp << ": " \
86 << e.what() << " FAIL" << std::endl; \
87 return -1; \
88 } \
89 catch (std::exception &e) { \
90 std::cout << "Unknown exception returned from " << pp << ": " \
91 << e.what() << " FAIL" << std::endl; \
92 return -1; \
93 }
94
95#define CATCH_EXCEPTIONS_WITH_COUNT(ierr, pp) \
96 catch (std::runtime_error &e) { \
97 std::cout << "Runtime exception returned from " << pp << ": " \
98 << e.what() << " FAIL" << std::endl; \
99 (ierr)++; \
100 } \
101 catch (std::logic_error &e) { \
102 std::cout << "Logic exception returned from " << pp << ": " \
103 << e.what() << " FAIL" << std::endl; \
104 (ierr)++; \
105 } \
106 catch (std::bad_alloc &e) { \
107 std::cout << "Bad_alloc exception returned from " << pp << ": " \
108 << e.what() << " FAIL" << std::endl; \
109 (ierr)++; \
110 } \
111 catch (std::exception &e) { \
112 std::cout << "Unknown exception returned from " << pp << ": " \
113 << e.what() << " FAIL" << std::endl; \
114 (ierr)++; \
115 }
116
122const char param_comment = '#';
123
125 const string& s,
126 const string& delimiters = " \f\n\r\t\v" )
127{
128 return s.substr( 0, s.find_last_not_of( delimiters ) + 1 );
129}
130
132 const string& s,
133 const string& delimiters = " \f\n\r\t\v" )
134{
135 return s.substr( s.find_first_not_of( delimiters ) );
136}
137
139 const string& s,
140 const string& delimiters = " \f\n\r\t\v" )
141{
142 return trim_left_copy( trim_right_copy( s, delimiters ), delimiters );
143}
144
145template <typename Adapter>
147 const char *str,
148 int dim,
149 typename Adapter::scalar_t *lower,
150 typename Adapter::scalar_t *upper,
151 size_t nparts,
152 typename Adapter::part_t *parts
153)
154{
155 std::cout << "boxAssign test " << str << ": Box (";
156 for (int j = 0; j < dim; j++) std::cout << lower[j] << " ";
157 std::cout << ") x (";
158 for (int j = 0; j < dim; j++) std::cout << upper[j] << " ";
159
160 if (nparts == 0)
161 std::cout << ") does not overlap any parts" << std::endl;
162 else {
163 std::cout << ") overlaps parts ";
164 for (size_t k = 0; k < nparts; k++) std::cout << parts[k] << " ";
165 std::cout << std::endl;
166 }
167}
168
169template <typename Adapter>
172 RCP<tMVector_t> &coords,
173 bool print_details)
174{
175 int ierr = 0;
176
177 // pointAssign tests
178 int coordDim = coords->getNumVectors();
179 zscalar_t *pointDrop = new zscalar_t[coordDim];
180 typename Adapter::part_t part = -1;
181
182 char mechar[10];
183 sprintf(mechar, "%d", problem->getComm()->getRank());
184 string me(mechar);
185
186 // test correctness of pointAssign for owned points
187 {
188 const typename Adapter::part_t *solnPartView =
189 problem->getSolution().getPartListView();
190
191 size_t numPoints = coords->getLocalLength();
192 for (size_t localID = 0; localID < numPoints; localID++) {
193
194 typename Adapter::part_t solnPart = solnPartView[localID];
195
196 for (int i = 0; i < coordDim; i++)
197 pointDrop[i] = coords->getData(i)[localID];
198
199 try {
200 part = problem->getSolution().pointAssign(coordDim, pointDrop);
201 }
202 CATCH_EXCEPTIONS_WITH_COUNT(ierr, me + ": pointAssign -- OwnedPoints");
203
204 if(print_details) {
205 std::cout << me << " Point " << localID
206 << " gid " << coords->getMap()->getGlobalElement(localID)
207 << " (" << pointDrop[0];
208 if (coordDim > 1) std::cout << " " << pointDrop[1];
209 if (coordDim > 2) std::cout << " " << pointDrop[2];
210 std::cout << ") in boxPart " << part
211 << " in solnPart " << solnPart
212 << std::endl;
213 }
214
215// this error test does not work for points that fall on the cuts.
216// like Zoltan's RCB, pointAssign arbitrarily picks a part along the cut.
217// the arbitrarily chosen part will not necessarily be the one to which
218// the coordinate was assigned in partitioning.
219//
220// if (part != solnPart) {
221// std::cout << me << " pointAssign: incorrect part " << part
222// << " found; should be " << solnPart
223// << " for point " << j << std::endl;
224// ierr++;
225// }
226 }
227 }
228
229 {
230 const std::vector<Zoltan2::coordinateModelPartBox>
231 pBoxes = problem->getSolution().getPartBoxesView();
232 if(print_details) {
233 for (size_t i = 0; i < pBoxes.size(); i++) {
234 typename Zoltan2::coordinateModelPartBox::coord_t *lmin = pBoxes[i].getlmins();
235 typename Zoltan2::coordinateModelPartBox::coord_t *lmax = pBoxes[i].getlmaxs();;
236 std::cout << me << " pBox " << i << " pid " << pBoxes[i].getpId()
237 << " (" << lmin[0] << "," << lmin[1] << ","
238 << (coordDim > 2 ? lmin[2] : 0) << ") x "
239 << " (" << lmax[0] << "," << lmax[1] << ","
240 << (coordDim > 2 ? lmax[2] : 0) << ")" << std::endl;
241 }
242 }
243 }
244
245 // test the origin
246 {
247 for (int i = 0; i < coordDim; i++) pointDrop[i] = 0.;
248 try {
249 part = problem->getSolution().pointAssign(coordDim, pointDrop);
250 }
251 CATCH_EXCEPTIONS_WITH_COUNT(ierr, me + " pointAssign -- Origin");
252
253 std::cout << me << " OriginPoint (" << pointDrop[0];
254 if (coordDim > 1) std::cout << " " << pointDrop[1];
255 if (coordDim > 2) std::cout << " " << pointDrop[2];
256 std::cout << ") part " << part << std::endl;
257 }
258
259 // test point with negative coordinates
260 {
261 for (int i = 0; i < coordDim; i++) pointDrop[i] = -100.+i;
262 try {
263 part = problem->getSolution().pointAssign(coordDim, pointDrop);
264 }
265 CATCH_EXCEPTIONS_WITH_COUNT(ierr, me + " pointAssign -- Negative Point");
266 std::cout << me << " NegativePoint (" << pointDrop[0];
267 if (coordDim > 1) std::cout << " " << pointDrop[1];
268 if (coordDim > 2) std::cout << " " << pointDrop[2];
269 std::cout << ") part " << part << std::endl;
270 }
271
272 // test a point that's way out there
273 {
274 for (int i = 0; i < coordDim; i++) pointDrop[i] = i*5;
275 try {
276 part = problem->getSolution().pointAssign(coordDim, pointDrop);
277 }
278 CATCH_EXCEPTIONS_WITH_COUNT(ierr, me + " pointAssign -- i*5 Point");
279 std::cout << me << " i*5-Point (" << pointDrop[0];
280 if (coordDim > 1) std::cout << " " << pointDrop[1];
281 if (coordDim > 2) std::cout << " " << pointDrop[2];
282 std::cout << ") part " << part << std::endl;
283 }
284
285 // test a point that's way out there
286 {
287 for (int i = 0; i < coordDim; i++) pointDrop[i] = 10+i*5;
288 try {
289 part = problem->getSolution().pointAssign(coordDim, pointDrop);
290 }
291 CATCH_EXCEPTIONS_WITH_COUNT(ierr, me + " pointAssign -- WoopWoop");
292 std::cout << me << " WoopWoop-Point (" << pointDrop[0];
293 if (coordDim > 1) std::cout << " " << pointDrop[1];
294 if (coordDim > 2) std::cout << " " << pointDrop[2];
295 std::cout << ") part " << part << std::endl;
296 }
297
298 delete [] pointDrop;
299 return ierr;
300}
301
302template <typename Adapter>
305 RCP<tMVector_t> &coords)
306{
307 int ierr = 0;
308
309 // boxAssign tests
310 int coordDim = coords->getNumVectors();
311 zscalar_t *lower = new zscalar_t[coordDim];
312 zscalar_t *upper = new zscalar_t[coordDim];
313
314 char mechar[10];
315 sprintf(mechar, "%d", problem->getComm()->getRank());
316 string me(mechar);
317
318 const std::vector<Zoltan2::coordinateModelPartBox>
319 pBoxes = problem->getSolution().getPartBoxesView();
320 size_t nBoxes = pBoxes.size();
321
322 // test a box that is smaller than a part
323 {
324 size_t nparts;
325 typename Adapter::part_t *parts;
326 size_t pickabox = nBoxes / 2;
327 for (int i = 0; i < coordDim; i++) {
328 zscalar_t dd = 0.2 * (pBoxes[pickabox].getlmaxs()[i] -
329 pBoxes[pickabox].getlmins()[i]);
330 lower[i] = pBoxes[pickabox].getlmins()[i] + dd;
331 upper[i] = pBoxes[pickabox].getlmaxs()[i] - dd;
332 }
333 try {
334 problem->getSolution().boxAssign(coordDim, lower, upper,
335 nparts, &parts);
336 }
337 CATCH_EXCEPTIONS_WITH_COUNT(ierr, me + " boxAssign -- smaller");
338 if (nparts > 1) {
339 std::cout << me << " FAIL boxAssign error: smaller test, nparts > 1"
340 << std::endl;
341 ierr++;
342 }
343 print_boxAssign_result<Adapter>("smallerbox", coordDim,
344 lower, upper, nparts, parts);
345 delete [] parts;
346 }
347
348 // test a box that is larger than a part
349 {
350 size_t nparts;
351 typename Adapter::part_t *parts;
352 size_t pickabox = nBoxes / 2;
353 for (int i = 0; i < coordDim; i++) {
354 zscalar_t dd = 0.2 * (pBoxes[pickabox].getlmaxs()[i] -
355 pBoxes[pickabox].getlmins()[i]);
356 lower[i] = pBoxes[pickabox].getlmins()[i] - dd;
357 upper[i] = pBoxes[pickabox].getlmaxs()[i] + dd;
358 }
359 try {
360 problem->getSolution().boxAssign(coordDim, lower, upper,
361 nparts, &parts);
362 }
363 CATCH_EXCEPTIONS_WITH_COUNT(ierr, me + " boxAssign -- larger");
364
365 // larger box should have at least two parts in it for k > 1.
366 if ((nBoxes > 1) && (nparts < 2)) {
367 std::cout << me << " FAIL boxAssign error: "
368 << "larger test, nparts < 2"
369 << std::endl;
370 ierr++;
371 }
372
373 // parts should include pickabox's part
374 bool found_pickabox = 0;
375 for (size_t i = 0; i < nparts; i++)
376 if (parts[i] == pBoxes[pickabox].getpId()) {
377 found_pickabox = 1;
378 break;
379 }
380 if (!found_pickabox) {
381 std::cout << me << " FAIL boxAssign error: "
382 << "larger test, pickabox not found"
383 << std::endl;
384 ierr++;
385 }
386
387 print_boxAssign_result<Adapter>("largerbox", coordDim,
388 lower, upper, nparts, parts);
389 delete [] parts;
390 }
391
392 // test a box that includes all parts
393 {
394 size_t nparts;
395 typename Adapter::part_t *parts;
396 for (int i = 0; i < coordDim; i++) {
397 lower[i] = std::numeric_limits<zscalar_t>::max();
398 upper[i] = std::numeric_limits<zscalar_t>::min();
399 }
400 for (size_t j = 0; j < nBoxes; j++) {
401 for (int i = 0; i < coordDim; i++) {
402 if (pBoxes[j].getlmins()[i] <= lower[i])
403 lower[i] = pBoxes[j].getlmins()[i];
404 if (pBoxes[j].getlmaxs()[i] >= upper[i])
405 upper[i] = pBoxes[j].getlmaxs()[i];
406 }
407 }
408 try {
409 problem->getSolution().boxAssign(coordDim, lower, upper,
410 nparts, &parts);
411 }
412 CATCH_EXCEPTIONS_WITH_COUNT(ierr, me + " boxAssign -- global");
413
414 // global box should have all parts
415 if (nparts != nBoxes) {
416 std::cout << me << " FAIL boxAssign error: "
417 << "global test, nparts found " << nparts
418 << " != num global parts " << nBoxes
419 << std::endl;
420 ierr++;
421 }
422 print_boxAssign_result<Adapter>("globalbox", coordDim,
423 lower, upper, nparts, parts);
424 delete [] parts;
425 }
426
427 // test a box that is bigger than the entire domain
428 // Assuming lower and upper are still set to the global box boundary
429 // from the previous test
430 {
431 size_t nparts;
432 typename Adapter::part_t *parts;
433 for (int i = 0; i < coordDim; i++) {
434 lower[i] -= 2.;
435 upper[i] += 2.;
436 }
437
438 try {
439 problem->getSolution().boxAssign(coordDim, lower, upper,
440 nparts, &parts);
441 }
442 CATCH_EXCEPTIONS_WITH_COUNT(ierr, me + " boxAssign -- bigdomain");
443
444 // bigdomain box should have all parts
445 if (nparts != nBoxes) {
446 std::cout << me << " FAIL boxAssign error: "
447 << "bigdomain test, nparts found " << nparts
448 << " != num global parts " << nBoxes
449 << std::endl;
450 ierr++;
451 }
452 print_boxAssign_result<Adapter>("bigdomainbox", coordDim,
453 lower, upper, nparts, parts);
454 delete [] parts;
455 }
456
457 // test a box that is way out there
458 // Assuming lower and upper are still set to at least the global box
459 // boundary from the previous test
460 {
461 size_t nparts;
462 typename Adapter::part_t *parts;
463 for (int i = 0; i < coordDim; i++) {
464 lower[i] = upper[i] + 10;
465 upper[i] += 20;
466 }
467
468 try {
469 problem->getSolution().boxAssign(coordDim, lower, upper,
470 nparts, &parts);
471 }
472 CATCH_EXCEPTIONS_WITH_COUNT(ierr, me + " boxAssign -- out there");
473
474 // For now, boxAssign returns zero if there is no box overlap.
475 // TODO: this result should be changed in boxAssign definition
476 if (nparts != 0) {
477 std::cout << me << " FAIL boxAssign error: "
478 << "outthere test, nparts found " << nparts
479 << " != zero"
480 << std::endl;
481 ierr++;
482 }
483 print_boxAssign_result<Adapter>("outthere box", coordDim,
484 lower, upper, nparts, parts);
485 delete [] parts;
486 }
487
488 delete [] lower;
489 delete [] upper;
490 return ierr;
491}
492
493void readGeoGenParams(string paramFileName, Teuchos::ParameterList &geoparams, const RCP<const Teuchos::Comm<int> > & comm){
494 std::string input = "";
495 char inp[25000];
496 for(int i = 0; i < 25000; ++i){
497 inp[i] = 0;
498 }
499
500 bool fail = false;
501 if(comm->getRank() == 0){
502
503 std::fstream inParam(paramFileName.c_str());
504 if (inParam.fail())
505 {
506 fail = true;
507 }
508 if(!fail)
509 {
510 std::string tmp = "";
511 getline (inParam,tmp);
512 while (!inParam.eof()){
513 if(tmp != ""){
514 tmp = trim_copy(tmp);
515 if(tmp != ""){
516 input += tmp + "\n";
517 }
518 }
519 getline (inParam,tmp);
520 }
521 inParam.close();
522 for (size_t i = 0; i < input.size(); ++i){
523 inp[i] = input[i];
524 }
525 }
526 }
527
528
529
530 int size = input.size();
531 if(fail){
532 size = -1;
533 }
534 comm->broadcast(0, sizeof(int), (char*) &size);
535 if(size == -1){
536 throw std::runtime_error("File " + paramFileName + " cannot be opened.");
537 }
538 comm->broadcast(0, size, inp);
539 std::istringstream inParam(inp);
540 string str;
541 getline (inParam,str);
542 while (!inParam.eof()){
543 if(str[0] != param_comment){
544 size_t pos = str.find('=');
545 if(pos == string::npos){
546 throw std::runtime_error("Invalid Line:" + str + " in parameter file");
547 }
548 string paramname = trim_copy(str.substr(0,pos));
549 string paramvalue = trim_copy(str.substr(pos + 1));
550 geoparams.set(paramname, paramvalue);
551 }
552 getline (inParam,str);
553 }
554}
555
556template<class bv_use_node_t>
557int compareWithBasicVectorAdapterTest(RCP<const Teuchos::Comm<int> > &comm,
558 Teuchos::RCP<Teuchos::ParameterList> params,
560 RCP<tMVector_t> coords,
561 Zoltan2::XpetraMultiVectorAdapter<tMVector_t>::scalar_t ** weights = NULL, int numWeightsPerCoord = 0) {
562
564
565 // Run a test with BasicVectorAdapter and xyzxyz format coordinates
566 const int bvme = comm->getRank();
567 const inputAdapter_t::lno_t bvlen =
568 inputAdapter_t::lno_t(coords->getLocalLength());
569 const size_t bvnvecs = coords->getNumVectors();
570 const size_t bvsize = coords->getNumVectors() * coords->getLocalLength();
571
572 ArrayRCP<inputAdapter_t::scalar_t> *bvtpetravectors =
573 new ArrayRCP<inputAdapter_t::scalar_t>[bvnvecs];
574 for (size_t i = 0; i < bvnvecs; i++)
575 bvtpetravectors[i] = coords->getDataNonConst(i);
576 int idx = 0;
577 inputAdapter_t::gno_t *bvgids = new
578 inputAdapter_t::gno_t[coords->getLocalLength()];
579 inputAdapter_t::scalar_t *bvcoordarr = new inputAdapter_t::scalar_t[bvsize];
580 for (inputAdapter_t::lno_t j = 0; j < bvlen; j++) {
581 bvgids[j] = coords->getMap()->getGlobalElement(j);
582 for (size_t i = 0; i < bvnvecs; i++) {
583 bvcoordarr[idx++] = bvtpetravectors[i][j];
584 }
585 }
586
587 // my test node type
588 typedef Zoltan2::BasicUserTypes<inputAdapter_t::scalar_t,
589 inputAdapter_t::lno_t,
590 inputAdapter_t::gno_t,
591 bv_use_node_t> bvtypes_t;
592 typedef Zoltan2::BasicVectorAdapter<bvtypes_t> bvadapter_t;
593 std::vector<const inputAdapter_t::scalar_t *> bvcoords(bvnvecs);
594 std::vector<int> bvstrides(bvnvecs);
595 for (size_t i = 0; i < bvnvecs; i++) {
596 bvcoords[i] = &bvcoordarr[i];
597 bvstrides[i] = bvnvecs;
598 }
599 std::vector<const inputAdapter_t::scalar_t *> bvwgts;
600 std::vector<int> bvwgtstrides;
601
602 if(numWeightsPerCoord > 0) {
603 bvwgts = std::vector<const inputAdapter_t::scalar_t *>(numWeightsPerCoord);
604 bvwgtstrides = std::vector<int>(coords->getLocalLength());
605 for (size_t i = 0; i < coords->getLocalLength(); i++) {
606 bvwgtstrides[i] = numWeightsPerCoord;
607 }
608 for (int i = 0; i < numWeightsPerCoord; i++) {
609 bvwgts[i] = weights[i];
610 }
611 }
612
613 bvadapter_t bvia(bvlen, bvgids, bvcoords, bvstrides,
614 bvwgts, bvwgtstrides);
615
617 try {
618 bvproblem = new Zoltan2::PartitioningProblem<bvadapter_t>(&bvia,
619 params.getRawPtr(),
620 comm);
621 }
622 CATCH_EXCEPTIONS_AND_RETURN("PartitioningProblem()")
623
624 try {
625 bvproblem->solve();
626 }
628
629 int ierr = 0;
630
631 // Compare with MultiVectorAdapter result
632 for (inputAdapter_t::lno_t i = 0; i < bvlen; i++) {
633 if (problem->getSolution().getPartListView()[i] !=
634 bvproblem->getSolution().getPartListView()[i]) {
635 std::cout << bvme << " " << i << " "
636 << coords->getMap()->getGlobalElement(i) << " " << bvgids[i]
637 << ": XMV " << problem->getSolution().getPartListView()[i]
638 << "; BMV " << bvproblem->getSolution().getPartListView()[i]
639 << " : FAIL" << std::endl;
640 ++ierr;
641 }
642 /* For debugging - plot all success as well
643 else {
644 std::cout << bvme << " " << i << " "
645 << coords->getMap()->getGlobalElement(i) << " " << bvgids[i]
646 << ": XMV " << problem->getSolution().getPartListView()[i]
647 << "; BMV " << bvproblem->getSolution().getPartListView()[i]
648 << " : PASS" << std::endl;
649 }
650 */
651 }
652
653 delete [] bvgids;
654 delete [] bvcoordarr;
655 delete [] bvtpetravectors;
656 delete bvproblem;
657
658 if (coords->getGlobalLength() < 40) {
659 int len = coords->getLocalLength();
660 const inputAdapter_t::part_t *zparts =
661 problem->getSolution().getPartListView();
662 for (int i = 0; i < len; i++)
663 std::cout << comm->getRank()
664 << " lid " << i
665 << " gid " << coords->getMap()->getGlobalElement(i)
666 << " part " << zparts[i] << std::endl;
667 }
668
669 return ierr;
670}
671
672template<class bv_use_node_t>
673int GeometricGenInterface(RCP<const Teuchos::Comm<int> > &comm,
674 int numTeams, int numParts, float imbalance,
675 std::string paramFile, std::string pqParts,
676 std::string pfname,
677 int k,
678 int migration_check_option,
679 int migration_all_to_all_type,
680 zscalar_t migration_imbalance_cut_off,
681 int migration_processor_assignment_type,
682 int migration_doMigration_type,
683 bool uvm,
684 bool print_details,
685 bool test_boxes,
686 bool rectilinear,
687 int mj_premigration_option,
688 int mj_premigration_coordinate_cutoff
689)
690{
691 int ierr = 0;
692 Teuchos::ParameterList geoparams("geo params");
693 readGeoGenParams(paramFile, geoparams, comm);
694
697 comm);
698
699 int coord_dim = gg->getCoordinateDimension();
700 int numWeightsPerCoord = gg->getNumWeights();
701 zlno_t numLocalPoints = gg->getNumLocalCoords();
702 zgno_t numGlobalPoints = gg->getNumGlobalCoords();
703 zscalar_t **coords = new zscalar_t * [coord_dim];
704 for(int i = 0; i < coord_dim; ++i){
705 coords[i] = new zscalar_t[numLocalPoints];
706 }
707 gg->getLocalCoordinatesCopy(coords);
708 zscalar_t **weight = NULL;
709 if (numWeightsPerCoord) {
710 weight= new zscalar_t * [numWeightsPerCoord];
711 for(int i = 0; i < numWeightsPerCoord; ++i){
712 weight[i] = new zscalar_t[numLocalPoints];
713 }
714 gg->getLocalWeightsCopy(weight);
715 }
716
717 delete gg;
718
719 // Run 1st test with MV which always runs UVM on
720 RCP<Tpetra::Map<zlno_t, zgno_t, znode_t> > mp = rcp(
721 new Tpetra::Map<zlno_t, zgno_t, znode_t>(numGlobalPoints,
722 numLocalPoints, 0, comm));
723 Teuchos::Array<Teuchos::ArrayView<const zscalar_t> > coordView(coord_dim);
724 for (int i=0; i < coord_dim; i++){
725 if(numLocalPoints > 0){
726 Teuchos::ArrayView<const zscalar_t> a(coords[i], numLocalPoints);
727 coordView[i] = a;
728 }
729 else {
730 Teuchos::ArrayView<const zscalar_t> a;
731 coordView[i] = a;
732 }
733 }
734 RCP<tMVector_t> tmVector = RCP<tMVector_t>(new
735 tMVector_t(mp, coordView.view(0, coord_dim),
736 coord_dim));
737 std::vector<const zscalar_t *> weights;
738 if(numWeightsPerCoord){
739 for (int i = 0; i < numWeightsPerCoord;++i){
740 weights.push_back(weight[i]);
741 }
742 }
743 std::vector<int> stride;
746 //inputAdapter_t ia(coordsConst);
747 inputAdapter_t *ia = new inputAdapter_t(tmVector, weights, stride);
748
749 Teuchos::RCP<Teuchos::ParameterList> params;
750
751 //Teuchos::ParameterList params("test params");
752 if(pfname != ""){
753 params = Teuchos::getParametersFromXmlFile(pfname);
754 }
755 else {
756 params =RCP<Teuchos::ParameterList>(new Teuchos::ParameterList, true);
757 }
758/*
759 params->set("memory_output_stream" , "std::cout");
760 params->set("memory_procs" , 0);
761 */
762
763 params->set("timer_output_stream" , "std::cout");
764
765 params->set("algorithm", "multijagged");
766 if (test_boxes)
767 params->set("mj_keep_part_boxes", true); // bool parameter
768 if (rectilinear)
769 params->set("rectilinear", true ); // bool parameter
770
771 if(imbalance > 1)
772 params->set("imbalance_tolerance", double(imbalance));
773 params->set("mj_premigration_option", mj_premigration_option);
774 if (mj_premigration_coordinate_cutoff > 0){
775 params->set("mj_premigration_coordinate_count",
776 mj_premigration_coordinate_cutoff);
777 }
778
779 if(pqParts != "")
780 params->set("mj_parts", pqParts);
781 if(numTeams > 0) {
782 params->set("mj_num_teams", numTeams);
783 }
784 if(numParts > 0)
785 params->set("num_global_parts", numParts);
786 if (k > 0)
787 params->set("mj_concurrent_part_count", k);
788 if(migration_check_option >= 0)
789 params->set("mj_migration_option", migration_check_option);
790 if(migration_imbalance_cut_off >= 0)
791 params->set("mj_minimum_migration_imbalance",
792 double(migration_imbalance_cut_off));
793
795 try {
797 params.getRawPtr(),
798 comm);
799 }
800 CATCH_EXCEPTIONS_AND_RETURN("PartitioningProblem()")
801
802 try {
803 problem->solve();
804 }
806 {
807 ierr += compareWithBasicVectorAdapterTest<bv_use_node_t>(
808 comm, params, problem, tmVector,
809 weight, numWeightsPerCoord);
810 }
811
812 // create metric object
813
814 RCP<quality_t> metricObject =
815 rcp(new quality_t(ia,params.getRawPtr(),comm,&problem->getSolution()));
816
817 if (comm->getRank() == 0){
818 metricObject->printMetrics(std::cout);
819 }
820
821 problem->printTimers();
822
823 // run pointAssign tests
824 if (test_boxes) {
825 ierr += run_pointAssign_tests<inputAdapter_t>(problem, tmVector,
826 print_details);
827 ierr += run_boxAssign_tests<inputAdapter_t>(problem, tmVector);
828 }
829
830 if(numWeightsPerCoord){
831 for(int i = 0; i < numWeightsPerCoord; ++i) {
832 delete [] weight[i];
833 }
834 delete [] weight;
835 }
836 if(coord_dim){
837 for(int i = 0; i < coord_dim; ++i) {
838 delete [] coords[i];
839 }
840 delete [] coords;
841 }
842
843 delete problem;
844 delete ia;
845 return ierr;
846}
847
848template<class bv_use_node_t>
850 RCP<const Teuchos::Comm<int> > &comm,
851 int numTeams,
852 int numParts,
853 float imbalance,
854 std::string fname,
855 std::string pqParts,
856 std::string pfname,
857 int k,
858 int migration_check_option,
859 int migration_all_to_all_type,
860 zscalar_t migration_imbalance_cut_off,
861 int migration_processor_assignment_type,
862 int migration_doMigration_type,
863 bool uvm,
864 bool print_details,
865 bool test_boxes,
866 bool rectilinear,
867 int mj_premigration_option,
868 int mj_premigration_coordinate_cutoff
869)
870{
871 int ierr = 0;
872 //std::string fname("simple");
873 //std::cout << "running " << fname << std::endl;
874
875 UserInputForTests uinput(testDataFilePath, fname, comm, true);
876
877 RCP<tMVector_t> coords = uinput.getUICoordinates();
878
881 inputAdapter_t *ia = new inputAdapter_t(coords);
882
883 Teuchos::RCP <Teuchos::ParameterList> params ;
884
885 //Teuchos::ParameterList params("test params");
886 if(pfname != ""){
887 params = Teuchos::getParametersFromXmlFile(pfname);
888 }
889 else {
890 params =RCP <Teuchos::ParameterList> (new Teuchos::ParameterList, true);
891 }
892
893 //params->set("timer_output_stream" , "std::cout");
894 if (test_boxes)
895 params->set("mj_keep_part_boxes", true); // bool parameter
896 if (rectilinear)
897 params->set("rectilinear", true); // bool parameter
898 params->set("algorithm", "multijagged");
899 if(imbalance > 1){
900 params->set("imbalance_tolerance", double(imbalance));
901 }
902
903 if(pqParts != ""){
904 params->set("mj_parts", pqParts);
905 }
906 params->set("mj_premigration_option", mj_premigration_option);
907
908 if(numParts > 0){
909 params->set("num_global_parts", numParts);
910 }
911 if (k > 0){
912 params->set("mj_concurrent_part_count", k);
913 }
914 if(migration_check_option >= 0){
915 params->set("mj_migration_option", migration_check_option);
916 }
917 if(migration_imbalance_cut_off >= 0){
918 params->set("mj_minimum_migration_imbalance",
919 double (migration_imbalance_cut_off));
920 }
921 if (mj_premigration_coordinate_cutoff > 0){
922 params->set("mj_premigration_coordinate_count",
923 mj_premigration_coordinate_cutoff);
924 }
925
927 try {
929 params.getRawPtr(),
930 comm);
931 }
932 CATCH_EXCEPTIONS_AND_RETURN("PartitioningProblem()")
933 try {
934 problem->solve();
935 }
937 {
938 ierr += compareWithBasicVectorAdapterTest<bv_use_node_t>(
939 comm, params, problem, coords);
940 }
941
942 // create metric object
943
944 RCP<quality_t> metricObject =
945 rcp(new quality_t(ia,params.getRawPtr(),comm,&problem->getSolution()));
946
947 if (comm->getRank() == 0){
948 metricObject->printMetrics(std::cout);
949 std::cout << "testFromDataFile is done " << std::endl;
950 }
951
952 problem->printTimers();
953
954 // run pointAssign tests
955 if (test_boxes) {
956 ierr += run_pointAssign_tests<inputAdapter_t>(problem, coords,
957 print_details);
958 ierr += run_boxAssign_tests<inputAdapter_t>(problem, coords);
959 }
960
961 delete problem;
962 delete ia;
963 return ierr;
964}
965
966#ifdef hopper_separate_test
967
968template <typename zscalar_t, typename zlno_t>
969void getCoords(zscalar_t **&coords, zlno_t &numLocal, int &dim, string fileName){
970 FILE *f = fopen(fileName.c_str(), "r");
971 if (f == NULL){
972 std::cout << fileName << " cannot be opened" << std::endl;
973 std::terminate();
974 }
975 fscanf(f, "%d", &numLocal);
976 fscanf(f, "%d", &dim);
977 coords = new zscalar_t *[ dim];
978 for (int i = 0; i < dim; ++i){
979 coords[i] = new zscalar_t[numLocal];
980 }
981 for (int i = 0; i < dim; ++i){
982 for (zlno_t j = 0; j < numLocal; ++j){
983 fscanf(f, "%lf", &(coords[i][j]));
984 }
985 }
986 fclose(f);
987 std::cout << "done reading:" << fileName << std::endl;
988}
989
990int testFromSeparateDataFiles(
991 RCP<const Teuchos::Comm<int> > &comm,
992 int numTeams,
993 int numParts,
994 float imbalance,
995 std::string fname,
996 std::string pqParts,
997 std::string pfname,
998 int k,
999 int migration_check_option,
1000 int migration_all_to_all_type,
1001 zscalar_t migration_imbalance_cut_off,
1002 int migration_processor_assignment_type,
1003 int migration_doMigration_type,
1004 bool uvm,
1005 bool test_boxes,
1006 bool rectilinear,
1007 int mj_premigration_option
1008)
1009{
1010 //std::string fname("simple");
1011 //std::cout << "running " << fname << std::endl;
1012
1013 int ierr = 0;
1014 int mR = comm->getRank();
1015 if (mR == 0) std::cout << "size of zscalar_t:" << sizeof(zscalar_t) << std::endl;
1016 string tFile = fname +"_" + std::to_string(mR) + ".mtx";
1017 zscalar_t **double_coords;
1018 zlno_t numLocal = 0;
1019 int dim = 0;
1020 getCoords<zscalar_t, zlno_t>(double_coords, numLocal, dim, tFile);
1021 //UserInputForTests uinput(testDataFilePath, fname, comm, true);
1022 Teuchos::Array<Teuchos::ArrayView<const zscalar_t> > coordView(dim);
1023 for (int i=0; i < dim; i++){
1024 if(numLocal > 0){
1025 Teuchos::ArrayView<const zscalar_t> a(double_coords[i], numLocal);
1026 coordView[i] = a;
1027 } else{
1028 Teuchos::ArrayView<const zscalar_t> a;
1029 coordView[i] = a;
1030 }
1031 }
1032
1033 zgno_t numGlobal;
1034 zgno_t nL = numLocal;
1035 Teuchos::Comm<int> *tcomm = (Teuchos::Comm<int> *)comm.getRawPtr();
1036
1037 reduceAll<int, zgno_t>(
1038 *tcomm,
1039 Teuchos::REDUCE_SUM,
1040 1,
1041 &nL,
1042 &numGlobal
1043 );
1044
1045
1046 RCP<Tpetra::Map<zlno_t, zgno_t, znode_t> > mp = rcp(
1047 new Tpetra::Map<zlno_t, zgno_t, znode_t> (numGlobal, numLocal, 0, comm));
1048 RCP< Tpetra::MultiVector<zscalar_t, zlno_t, zgno_t, znode_t> >coords = RCP< Tpetra::MultiVector<zscalar_t, zlno_t, zgno_t, znode_t> >(
1049 new Tpetra::MultiVector<zscalar_t, zlno_t, zgno_t, znode_t>( mp, coordView.view(0, dim), dim));
1050
1051
1052 RCP<const tMVector_t> coordsConst = rcp_const_cast<const tMVector_t>(coords);
1053
1054 typedef Zoltan2::XpetraMultiVectorAdapter<tMVector_t> inputAdapter_t;
1056
1057 inputAdapter_t *ia = new inputAdapter_t(coordsConst);
1058
1059
1060
1061 Teuchos::RCP <Teuchos::ParameterList> params ;
1062
1063 //Teuchos::ParameterList params("test params");
1064 if(pfname != ""){
1065 params = Teuchos::getParametersFromXmlFile(pfname);
1066 }
1067 else {
1068 params =RCP <Teuchos::ParameterList> (new Teuchos::ParameterList, true);
1069 }
1070
1071 //params->set("timer_output_stream" , "std::cout");
1072 params->set("algorithm", "multijagged");
1073 if(imbalance > 1){
1074 params->set("imbalance_tolerance", double(imbalance));
1075 }
1076
1077 params->set("mj_premigration_option", mj_premigration_option);
1078 if(pqParts != ""){
1079 params->set("mj_parts", pqParts);
1080 }
1081 if(numTeams > 0){
1082 params->set("mj_num_teams", numTeams);
1083 }
1084 if(numParts > 0){
1085 params->set("num_global_parts", numParts);
1086 }
1087 if (k > 0){
1088 params->set("parallel_part_calculation_count", k);
1089 }
1090 if(migration_processor_assignment_type >= 0){
1091 params->set("migration_processor_assignment_type", migration_processor_assignment_type);
1092 }
1093 if(migration_check_option >= 0){
1094 params->set("migration_check_option", migration_check_option);
1095 }
1096 if(migration_all_to_all_type >= 0){
1097 params->set("migration_all_to_all_type", migration_all_to_all_type);
1098 }
1099 if(migration_imbalance_cut_off >= 0){
1100 params->set("migration_imbalance_cut_off",
1101 double (migration_imbalance_cut_off));
1102 }
1103 if (migration_doMigration_type >= 0){
1104 params->set("migration_doMigration_type", int (migration_doMigration_type));
1105 }
1106 if (test_boxes)
1107 params->set("mj_keep_part_boxes", true); // bool parameter
1108 if (rectilinear)
1109 params->set("rectilinear", true); // bool parameter
1110
1112 try {
1113 problem =
1115 params.getRawPtr(),
1116 comm);
1117 }
1118 CATCH_EXCEPTIONS_AND_RETURN("PartitioningProblem()")
1119
1120 try {
1121 problem->solve();
1122 }
1124 {
1125 ierr += compareWithBasicVectorAdapterTest<bv_use_node_t>(
1126 comm, params, problem, coords);
1127 }
1128
1129 if (coordsConst->getGlobalLength() < 40) {
1130 int len = coordsConst->getLocalLength();
1131 const inputAdapter_t::part_t *zparts =
1132 problem->getSolution().getPartListView();
1133 for (int i = 0; i < len; i++)
1134 std::cout << comm->getRank()
1135 << " gid " << coords->getMap()->getGlobalElement(i)
1136 << " part " << zparts[i] << std::endl;
1137 }
1138
1139 //create metric object
1140
1141 RCP<quality_t> metricObject =
1142 rcp(new quality_t(ia,params.getRawPtr(),comm,&problem->getSolution()));
1143
1144 if (comm->getRank() == 0){
1145 metricObject->printMetrics(std::cout);
1146 std::cout << "testFromDataFile is done " << std::endl;
1147 }
1148
1149 problem->printTimers();
1150
1151 // run pointAssign tests
1152 if (test_boxes) {
1153 ierr += run_pointAssign_tests<inputAdapter_t>(problem, coords,
1154 print_details);
1155 ierr += run_boxAssign_tests<inputAdapter_t>(problem, coords);
1156 }
1157
1158 delete problem;
1159 delete ia;
1160 return ierr;
1161}
1162#endif
1163
1164
1165
1166string convert_to_string(char *args){
1167 string tmp = "";
1168 for(int i = 0; args[i] != 0; i++)
1169 tmp += args[i];
1170 return tmp;
1171}
1172bool getArgumentValue(string &argumentid, double &argumentValue, string argumentline){
1173 std::stringstream stream(std::stringstream::in | std::stringstream::out);
1174 stream << argumentline;
1175 getline(stream, argumentid, '=');
1176 if (stream.eof()){
1177 return false;
1178 }
1179 stream >> argumentValue;
1180 return true;
1181}
1182
1184 int narg,
1185 char **arg,
1186 int &numTeams,
1187 int &numParts,
1188 float &imbalance ,
1189 string &pqParts,
1190 int &opt,
1191 std::string &fname,
1192 std::string &pfname,
1193 int &k,
1194 int &migration_check_option,
1195 int &migration_all_to_all_type,
1196 zscalar_t &migration_imbalance_cut_off,
1197 int &migration_processor_assignment_type,
1198 int &migration_doMigration_type,
1199 bool &uvm,
1200 bool &print_details,
1201 bool &test_boxes,
1202 bool &rectilinear,
1203 int &mj_premigration_option,
1204 int &mj_coordinate_cutoff
1205)
1206{
1207 bool isCset = false;
1208 bool isPset = false;
1209 bool isFset = false;
1210 bool isPFset = false;
1211
1212 for(int i = 0; i < narg; ++i){
1213 string tmp = convert_to_string(arg[i]);
1214 string identifier = "";
1215 long long int value = -1; double fval = -1;
1216 if(!getArgumentValue(identifier, fval, tmp)) continue;
1217 value = (long long int) (fval);
1218
1219 if(identifier == "W"){
1220 if(value == 0 || value == 1){
1221 print_details = (value == 0 ? false : true);
1222 } else {
1223 throw std::runtime_error("Invalid argument at " + tmp);
1224 }
1225 }
1226 else if(identifier == "UVM"){
1227 if(value == 0 || value == 1){
1228 uvm = (value == 0 ? false : true);
1229 } else {
1230 throw std::runtime_error("Invalid argument at " + tmp);
1231 }
1232 }
1233 else if(identifier == "T"){
1234 if(value > 0){
1235 numTeams=value;
1236 } else {
1237 throw std::runtime_error("Invalid argument at " + tmp);
1238 }
1239 } else if(identifier == "C"){
1240 if(value > 0){
1241 numParts=value;
1242 isCset = true;
1243 } else {
1244 throw std::runtime_error("Invalid argument at " + tmp);
1245 }
1246 } else if(identifier == "P"){
1247 std::stringstream stream(std::stringstream::in | std::stringstream::out);
1248 stream << tmp;
1249 string ttmp;
1250 getline(stream, ttmp, '=');
1251 stream >> pqParts;
1252 isPset = true;
1253 }else if(identifier == "I"){
1254 if(fval > 0){
1255 imbalance=fval;
1256 } else {
1257 throw std::runtime_error("Invalid argument at " + tmp);
1258 }
1259 } else if(identifier == "MI"){
1260 if(fval > 0){
1261 migration_imbalance_cut_off=fval;
1262 } else {
1263 throw std::runtime_error("Invalid argument at " + tmp);
1264 }
1265 } else if(identifier == "MO"){
1266 if(value >=0 ){
1267 migration_check_option = value;
1268 } else {
1269 throw std::runtime_error("Invalid argument at " + tmp);
1270 }
1271 } else if(identifier == "AT"){
1272 if(value >=0 ){
1273 migration_processor_assignment_type = value;
1274 } else {
1275 throw std::runtime_error("Invalid argument at " + tmp);
1276 }
1277 }
1278
1279 else if(identifier == "MT"){
1280 if(value >=0 ){
1281 migration_all_to_all_type = value;
1282 } else {
1283 throw std::runtime_error("Invalid argument at " + tmp);
1284 }
1285 }
1286 else if(identifier == "PCC"){
1287 if(value >=0 ){
1288 mj_coordinate_cutoff = value;
1289 } else {
1290 throw std::runtime_error("Invalid argument at " + tmp);
1291 }
1292 }
1293
1294 else if(identifier == "PM"){
1295 if(value >=0 ){
1296 mj_premigration_option = value;
1297 } else {
1298 throw std::runtime_error("Invalid argument at " + tmp);
1299 }
1300 }
1301
1302 else if(identifier == "DM"){
1303 if(value >=0 ){
1304 migration_doMigration_type = value;
1305 } else {
1306 throw std::runtime_error("Invalid argument at " + tmp);
1307 }
1308 }
1309 else if(identifier == "F"){
1310 std::stringstream stream(std::stringstream::in | std::stringstream::out);
1311 stream << tmp;
1312 getline(stream, fname, '=');
1313
1314 stream >> fname;
1315 isFset = true;
1316 }
1317 else if(identifier == "PF"){
1318 std::stringstream stream(std::stringstream::in | std::stringstream::out);
1319 stream << tmp;
1320 getline(stream, pfname, '=');
1321
1322 stream >> pfname;
1323 isPFset = true;
1324 }
1325
1326 else if(identifier == "O"){
1327 if(value >= 0 && value <= 3){
1328 opt = value;
1329 } else {
1330 throw std::runtime_error("Invalid argument at " + tmp);
1331 }
1332 }
1333 else if(identifier == "K"){
1334 if(value >=0 ){
1335 k = value;
1336 } else {
1337 throw std::runtime_error("Invalid argument at " + tmp);
1338 }
1339 }
1340 else if(identifier == "TB"){
1341 if(value >=0 ){
1342 test_boxes = (value == 0 ? false : true);
1343 } else {
1344 throw std::runtime_error("Invalid argument at " + tmp);
1345 }
1346 }
1347 else if(identifier == "R"){
1348 if(value >=0 ){
1349 rectilinear = (value == 0 ? false : true);
1350 } else {
1351 throw std::runtime_error("Invalid argument at " + tmp);
1352 }
1353 }
1354 else {
1355 throw std::runtime_error("Invalid argument at " + tmp);
1356 }
1357
1358 }
1359 if(!( (isCset || isPset || isPFset) && isFset)){
1360 throw std::runtime_error("(C || P || PF) && F are mandatory arguments.");
1361 }
1362
1363}
1364
1365void print_usage(char *executable){
1366 std::cout << "\nUsage:" << std::endl;
1367 std::cout << executable << " arglist" << std::endl;
1368 std::cout << "arglist:" << std::endl;
1369 std::cout << "\tT=numTeams: numTeams > 0" << std::endl;
1370 std::cout << "\tC=numParts: numParts > 0" << std::endl;
1371 std::cout << "\tP=MultiJaggedPart: Example: P=512,512" << std::endl;
1372 std::cout << "\tI=imbalance: Example I=1.03 (ignored for now.)" << std::endl;
1373 std::cout << "\tF=filePath: When O=0 the path of the coordinate input file, for O>1 the path to the geometric generator parameter file." << std::endl;
1374 std::cout << "\tO=input option: O=0 for reading coordinate from file, O>0 for generating coordinate from coordinate generator file. Default will run geometric generator." << std::endl;
1375 std::cout << "\tK=concurrent part calculation input: K>0." << std::endl;
1376 std::cout << "\tMI=migration_imbalance_cut_off: MI=1.35. " << std::endl;
1377 std::cout << "\tMT=migration_all_to_all_type: 0 for alltoallv, 1 for Zoltan_Comm, 2 for Zoltan2 Distributor object(Default 1)." << std::endl;
1378 std::cout << "\tMO=migration_check_option: 0 for decision on imbalance, 1 for forcing migration, >1 for avoiding migration. (Default-0)" << std::endl;
1379 std::cout << "\tAT=migration_processor_assignment_type. 0-for assigning procs with respect to proc ownment, otherwise, assignment with respect to proc closeness." << std::endl;
1380 std::cout << "Example:\n" << executable << " P=2,2,2 C=8 F=simple O=0" << std::endl;
1381}
1382
1383int main(int narg, char *arg[])
1384{
1385 Tpetra::ScopeGuard tscope(&narg, &arg);
1386 Teuchos::RCP<const Teuchos::Comm<int> > tcomm = Tpetra::getDefaultComm();
1387
1388 int rank = tcomm->getRank();
1389
1390
1391 int numTeams = 0; // will use default if not set
1392 int numParts = -10;
1393 float imbalance = -1.03;
1394 int k = -1;
1395
1396 string pqParts = "";
1397 int opt = 1;
1398 std::string fname = "";
1399 std::string paramFile = "";
1400
1401
1402 int migration_check_option = -2;
1403 int migration_all_to_all_type = -1;
1404 zscalar_t migration_imbalance_cut_off = -1.15;
1405 int migration_processor_assignment_type = -1;
1406 int migration_doMigration_type = -1;
1407 int mj_premigration_option = 0;
1408 int mj_premigration_coordinate_cutoff = 0;
1409
1410 bool uvm = true;
1411 bool print_details = true;
1412 bool test_boxes = false;
1413 bool rectilinear = false;
1414
1415 // make a new node type so we can run BasicVectorAdapter with UVM off
1416 // The Tpetra MV will still run with UVM on and we'll compare the results.
1417 // For Serial/OpenMP the 2nd test will be turned off at the CMake level.
1418 // For CUDA we control uvm on/off with parameter uvm set to 0 or 1.
1419 // For HIP uvm is simply always off.
1420 // TODO: Probably this should all change eventually so we don't have a node
1421 // declared like this.
1422#if defined(KOKKOS_ENABLE_CUDA)
1423 using uvm_off_node_t = Kokkos::Compat::KokkosDeviceWrapperNode<
1424 Kokkos::Cuda, Kokkos::CudaSpace>;
1425#elif defined(KOKKOS_ENABLE_HIP)
1426 using uvm_off_node_t = Kokkos::Compat::KokkosDeviceWrapperNode<
1427 Kokkos::Experimental::HIP, Kokkos::Experimental::HIPSpace>;
1428#endif
1429
1430 try{
1431 try {
1432 getArgVals(
1433 narg,
1434 arg,
1435 numTeams,
1436 numParts,
1437 imbalance ,
1438 pqParts,
1439 opt,
1440 fname,
1441 paramFile,
1442 k,
1443 migration_check_option,
1444 migration_all_to_all_type,
1445 migration_imbalance_cut_off,
1446 migration_processor_assignment_type,
1447 migration_doMigration_type,
1448 uvm,
1449 print_details,
1450 test_boxes,
1451 rectilinear,
1452 mj_premigration_option, mj_premigration_coordinate_cutoff);
1453 }
1454 catch(std::exception &s){
1455 if(tcomm->getRank() == 0){
1456 print_usage(arg[0]);
1457 }
1458 throw s;
1459 }
1460
1461 int ierr = 0;
1462
1463 switch (opt){
1464
1465 case 0:
1466 if(uvm == true) { // true by default, if not CUDA it should be unset
1467 ierr = testFromDataFile<znode_t>(tcomm, numTeams, numParts,
1468 imbalance, fname, pqParts, paramFile, k,
1469 migration_check_option,
1470 migration_all_to_all_type,
1471 migration_imbalance_cut_off,
1472 migration_processor_assignment_type,
1473 migration_doMigration_type, uvm, print_details, test_boxes,
1474 rectilinear,
1475 mj_premigration_option, mj_premigration_coordinate_cutoff);
1476 }
1477 else {
1478#if defined(KOKKOS_ENABLE_CUDA) || defined(KOKKOS_ENABLE_HIP)
1479 ierr = testFromDataFile<uvm_off_node_t>(tcomm, numTeams, numParts,
1480 imbalance, fname, pqParts, paramFile, k,
1481 migration_check_option,
1482 migration_all_to_all_type,
1483 migration_imbalance_cut_off,
1484 migration_processor_assignment_type,
1485 migration_doMigration_type, uvm, print_details, test_boxes,
1486 rectilinear,
1487 mj_premigration_option, mj_premigration_coordinate_cutoff);
1488#else
1489 throw std::logic_error("uvm set off but this is not a cuda/hip test.");
1490#endif
1491 }
1492 break;
1493
1494#ifdef hopper_separate_test
1495 case 1:
1496 // TODO: Note made changes here but did not actually run this
1497 // method.
1498 ierr = testFromSeparateDataFiles(tcomm, numTeams, numParts,
1499 imbalance, fname, pqParts, paramFile, k,
1500 migration_check_option,
1501 migration_all_to_all_type,
1502 migration_imbalance_cut_off,
1503 migration_processor_assignment_type,
1504 migration_doMigration_type,uvm, print_details, test_boxes,
1505 rectilinear,
1506 mj_premigration_option, mj_premigration_coordinate_cutoff);
1507 break;
1508#endif
1509 default:
1510 if(uvm == true) { // true by default, if not CUDA it should be unset
1511 ierr = GeometricGenInterface<znode_t>(tcomm, numTeams, numParts,
1512 imbalance, fname, pqParts, paramFile, k,
1513 migration_check_option,
1514 migration_all_to_all_type,
1515 migration_imbalance_cut_off,
1516 migration_processor_assignment_type,
1517 migration_doMigration_type, uvm, print_details, test_boxes,
1518 rectilinear,
1519 mj_premigration_option, mj_premigration_coordinate_cutoff);
1520 }
1521 else {
1522#if defined(KOKKOS_ENABLE_CUDA) || defined(KOKKOS_ENABLE_HIP)
1523 ierr = GeometricGenInterface<uvm_off_node_t>(tcomm, numTeams,
1524 numParts, imbalance, fname, pqParts, paramFile, k,
1525 migration_check_option,
1526 migration_all_to_all_type,
1527 migration_imbalance_cut_off,
1528 migration_processor_assignment_type,
1529 migration_doMigration_type, uvm, print_details, test_boxes,
1530 rectilinear,
1531 mj_premigration_option, mj_premigration_coordinate_cutoff);
1532#else
1533 throw std::logic_error("uvm set off but this is not a cuda test.");
1534#endif
1535 }
1536 break;
1537 }
1538
1539 if (rank == 0) {
1540 if (ierr == 0) std::cout << "PASS" << std::endl;
1541 else std::cout << "FAIL" << std::endl;
1542 }
1543 }
1544
1545
1546 catch(std::string &s){
1547 if (rank == 0)
1548 std::cerr << s << std::endl;
1549 }
1550
1551 catch(char * s){
1552 if (rank == 0)
1553 std::cerr << s << std::endl;
1554 }
1555
1556 return 0;
1557}
#define CATCH_EXCEPTIONS_WITH_COUNT(ierr, pp)
int run_pointAssign_tests(Zoltan2::PartitioningProblem< Adapter > *problem, RCP< tMVector_t > &coords, bool print_details)
void readGeoGenParams(string paramFileName, Teuchos::ParameterList &geoparams, const RCP< const Teuchos::Comm< int > > &comm)
string trim_left_copy(const string &s, const string &delimiters=" \f\n\r\t\v")
string convert_to_string(char *args)
bool getArgumentValue(string &argumentid, double &argumentValue, string argumentline)
void getArgVals(int narg, char **arg, int &numTeams, int &numParts, float &imbalance, string &pqParts, int &opt, std::string &fname, std::string &pfname, int &k, int &migration_check_option, int &migration_all_to_all_type, zscalar_t &migration_imbalance_cut_off, int &migration_processor_assignment_type, int &migration_doMigration_type, bool &uvm, bool &print_details, bool &test_boxes, bool &rectilinear, int &mj_premigration_option, int &mj_coordinate_cutoff)
int testFromDataFile(RCP< const Teuchos::Comm< int > > &comm, int numTeams, int numParts, float imbalance, std::string fname, std::string pqParts, std::string pfname, int k, int migration_check_option, int migration_all_to_all_type, zscalar_t migration_imbalance_cut_off, int migration_processor_assignment_type, int migration_doMigration_type, bool uvm, bool print_details, bool test_boxes, bool rectilinear, int mj_premigration_option, int mj_premigration_coordinate_cutoff)
string trim_right_copy(const string &s, const string &delimiters=" \f\n\r\t\v")
#define CATCH_EXCEPTIONS_AND_RETURN(pp)
string trim_copy(const string &s, const string &delimiters=" \f\n\r\t\v")
const char param_comment
void print_usage(char *executable)
void print_boxAssign_result(const char *str, int dim, typename Adapter::scalar_t *lower, typename Adapter::scalar_t *upper, size_t nparts, typename Adapter::part_t *parts)
int compareWithBasicVectorAdapterTest(RCP< const Teuchos::Comm< int > > &comm, Teuchos::RCP< Teuchos::ParameterList > params, Zoltan2::PartitioningProblem< Zoltan2::XpetraMultiVectorAdapter< tMVector_t > > *problem, RCP< tMVector_t > coords, Zoltan2::XpetraMultiVectorAdapter< tMVector_t >::scalar_t **weights=NULL, int numWeightsPerCoord=0)
int GeometricGenInterface(RCP< const Teuchos::Comm< int > > &comm, int numTeams, int numParts, float imbalance, std::string paramFile, std::string pqParts, std::string pfname, int k, int migration_check_option, int migration_all_to_all_type, zscalar_t migration_imbalance_cut_off, int migration_processor_assignment_type, int migration_doMigration_type, bool uvm, bool print_details, bool test_boxes, bool rectilinear, int mj_premigration_option, int mj_premigration_coordinate_cutoff)
int run_boxAssign_tests(Zoltan2::PartitioningProblem< Adapter > *problem, RCP< tMVector_t > &coords)
Defines the BasicVectorAdapter class.
Defines the EvaluatePartition class.
Defines the PartitioningProblem class.
Defines the PartitioningSolution 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 XpetraMultiVectorAdapter.
int main()
InputTraits< User >::scalar_t scalar_t
A simple class that can be the User template argument for an InputAdapter.
BasicVectorAdapter represents a vector (plus optional weights) supplied by the user as pointers to st...
A class that computes and returns quality metrics.
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.
RCP< const Comm< int > > getComm()
Return the communicator used by the problem.
void printTimers() const
Return the communicator passed to the problem.
static const std::string fail
Tpetra::MultiVector< zscalar_t, zlno_t, zgno_t, znode_t > tMVector_t
fname
Begin.
Definition validXML.py:19
static ArrayRCP< ArrayRCP< zscalar_t > > weights
void getCoords(void *data, int numGid, int numLid, int numObj, zgno_t *gids, zgno_t *lids, int dim, double *coords, int *ierr)
Zoltan2::EvaluatePartition< matrixAdapter_t > quality_t