60 #include "Teuchos_XMLParameterListHelpers.hpp"
62 #include <Teuchos_LAPACK.hpp>
70 #ifdef hopper_separate_test
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; \
79 catch (std::logic_error &e) { \
80 std::cout << "Logic exception returned from " << pp << ": " \
81 << e.what() << " FAIL" << std::endl; \
84 catch (std::bad_alloc &e) { \
85 std::cout << "Bad_alloc exception returned from " << pp << ": " \
86 << e.what() << " FAIL" << std::endl; \
89 catch (std::exception &e) { \
90 std::cout << "Unknown exception returned from " << pp << ": " \
91 << e.what() << " FAIL" << std::endl; \
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; \
101 catch (std::logic_error &e) { \
102 std::cout << "Logic exception returned from " << pp << ": " \
103 << e.what() << " FAIL" << std::endl; \
106 catch (std::bad_alloc &e) { \
107 std::cout << "Bad_alloc exception returned from " << pp << ": " \
108 << e.what() << " FAIL" << std::endl; \
111 catch (std::exception &e) { \
112 std::cout << "Unknown exception returned from " << pp << ": " \
113 << e.what() << " FAIL" << std::endl; \
126 const string& delimiters =
" \f\n\r\t\v" )
128 return s.substr( 0, s.find_last_not_of( delimiters ) + 1 );
133 const string& delimiters =
" \f\n\r\t\v" )
135 return s.substr( s.find_first_not_of( delimiters ) );
140 const string& delimiters =
" \f\n\r\t\v" )
145 template <
typename Adapter>
149 typename Adapter::scalar_t *lower,
150 typename Adapter::scalar_t *upper,
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] <<
" ";
161 std::cout <<
") does not overlap any parts" << std::endl;
163 std::cout <<
") overlaps parts ";
164 for (
size_t k = 0; k < nparts; k++) std::cout << parts[k] <<
" ";
165 std::cout << std::endl;
169 template <
typename Adapter>
172 RCP<tMVector_t> &coords,
178 int coordDim = coords->getNumVectors();
183 sprintf(mechar,
"%d", problem->
getComm()->getRank());
191 size_t numPoints = coords->getLocalLength();
192 for (
size_t localID = 0; localID < numPoints; localID++) {
196 for (
int i = 0; i < coordDim; i++)
197 pointDrop[i] = coords->getData(i)[localID];
200 part = problem->
getSolution().pointAssign(coordDim, pointDrop);
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
230 const std::vector<Zoltan2::coordinateModelPartBox>
231 pBoxes = problem->
getSolution().getPartBoxesView();
233 for (
size_t i = 0; i < pBoxes.size(); i++) {
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;
247 for (
int i = 0; i < coordDim; i++) pointDrop[i] = 0.;
249 part = problem->
getSolution().pointAssign(coordDim, pointDrop);
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;
261 for (
int i = 0; i < coordDim; i++) pointDrop[i] = -100.+i;
263 part = problem->
getSolution().pointAssign(coordDim, pointDrop);
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;
274 for (
int i = 0; i < coordDim; i++) pointDrop[i] = i*5;
276 part = problem->
getSolution().pointAssign(coordDim, pointDrop);
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;
287 for (
int i = 0; i < coordDim; i++) pointDrop[i] = 10+i*5;
289 part = problem->
getSolution().pointAssign(coordDim, pointDrop);
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;
302 template <
typename Adapter>
305 RCP<tMVector_t> &coords)
310 int coordDim = coords->getNumVectors();
315 sprintf(mechar,
"%d", problem->
getComm()->getRank());
318 const std::vector<Zoltan2::coordinateModelPartBox>
319 pBoxes = problem->
getSolution().getPartBoxesView();
320 size_t nBoxes = pBoxes.size();
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;
334 problem->
getSolution().boxAssign(coordDim, lower, upper,
339 std::cout << me <<
" FAIL boxAssign error: smaller test, nparts > 1"
343 print_boxAssign_result<Adapter>(
"smallerbox", coordDim,
344 lower, upper, nparts, 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;
360 problem->
getSolution().boxAssign(coordDim, lower, upper,
366 if ((nBoxes > 1) && (nparts < 2)) {
367 std::cout << me <<
" FAIL boxAssign error: "
368 <<
"larger test, nparts < 2"
374 bool found_pickabox = 0;
375 for (
size_t i = 0; i < nparts; i++)
376 if (parts[i] == pBoxes[pickabox].getpId()) {
380 if (!found_pickabox) {
381 std::cout << me <<
" FAIL boxAssign error: "
382 <<
"larger test, pickabox not found"
387 print_boxAssign_result<Adapter>(
"largerbox", coordDim,
388 lower, upper, nparts, 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();
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];
409 problem->
getSolution().boxAssign(coordDim, lower, upper,
415 if (nparts != nBoxes) {
416 std::cout << me <<
" FAIL boxAssign error: "
417 <<
"global test, nparts found " << nparts
418 <<
" != num global parts " << nBoxes
422 print_boxAssign_result<Adapter>(
"globalbox", coordDim,
423 lower, upper, nparts, parts);
433 for (
int i = 0; i < coordDim; i++) {
439 problem->
getSolution().boxAssign(coordDim, lower, upper,
445 if (nparts != nBoxes) {
446 std::cout << me <<
" FAIL boxAssign error: "
447 <<
"bigdomain test, nparts found " << nparts
448 <<
" != num global parts " << nBoxes
452 print_boxAssign_result<Adapter>(
"bigdomainbox", coordDim,
453 lower, upper, nparts, parts);
463 for (
int i = 0; i < coordDim; i++) {
464 lower[i] = upper[i] + 10;
469 problem->
getSolution().boxAssign(coordDim, lower, upper,
477 std::cout << me <<
" FAIL boxAssign error: "
478 <<
"outthere test, nparts found " << nparts
483 print_boxAssign_result<Adapter>(
"outthere box", coordDim,
484 lower, upper, nparts, parts);
493 void readGeoGenParams(
string paramFileName, Teuchos::ParameterList &geoparams,
const RCP<
const Teuchos::Comm<int> > & comm){
494 std::string input =
"";
496 for(
int i = 0; i < 25000; ++i){
501 if(comm->getRank() == 0){
503 std::fstream inParam(paramFileName.c_str());
510 std::string tmp =
"";
511 getline (inParam,tmp);
512 while (!inParam.eof()){
519 getline (inParam,tmp);
522 for (
size_t i = 0; i < input.size(); ++i){
530 int size = input.size();
534 comm->broadcast(0,
sizeof(
int), (
char*) &size);
536 throw std::runtime_error(
"File " + paramFileName +
" cannot be opened.");
538 comm->broadcast(0, size, inp);
539 std::istringstream inParam(inp);
541 getline (inParam,str);
542 while (!inParam.eof()){
544 size_t pos = str.find(
'=');
545 if(pos == string::npos){
546 throw std::runtime_error(
"Invalid Line:" + str +
" in parameter file");
548 string paramname =
trim_copy(str.substr(0,pos));
549 string paramvalue =
trim_copy(str.substr(pos + 1));
550 geoparams.set(paramname, paramvalue);
552 getline (inParam,str);
556 template<
class bv_use_node_t>
558 Teuchos::RCP<Teuchos::ParameterList> params,
560 RCP<tMVector_t> coords,
566 const int bvme = comm->getRank();
569 const size_t bvnvecs = coords->getNumVectors();
570 const size_t bvsize = coords->getNumVectors() * coords->getLocalLength();
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);
579 inputAdapter_t::scalar_t *bvcoordarr =
new inputAdapter_t::scalar_t[bvsize];
581 bvgids[j] = coords->getMap()->getGlobalElement(j);
582 for (
size_t i = 0; i < bvnvecs; i++) {
583 bvcoordarr[idx++] = bvtpetravectors[i][j];
591 bv_use_node_t> bvtypes_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;
599 std::vector<const inputAdapter_t::scalar_t *> bvwgts;
600 std::vector<int> bvwgtstrides;
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;
608 for (
int i = 0; i < numWeightsPerCoord; i++) {
613 bvadapter_t bvia(bvlen, bvgids, bvcoords, bvstrides,
614 bvwgts, bvwgtstrides);
632 for (inputAdapter_t::lno_t i = 0; i < bvlen; i++) {
633 if (problem->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;
654 delete [] bvcoordarr;
655 delete [] bvtpetravectors;
658 if (coords->getGlobalLength() < 40) {
659 int len = coords->getLocalLength();
661 problem->getSolution().getPartListView();
662 for (
int i = 0; i < len; i++)
663 std::cout << comm->getRank()
665 <<
" gid " << coords->getMap()->getGlobalElement(i)
666 <<
" part " << zparts[i] << std::endl;
672 template<
class bv_use_node_t>
674 int numTeams,
int numParts,
float imbalance,
675 std::string paramFile, std::string pqParts,
678 int migration_check_option,
679 int migration_all_to_all_type,
681 int migration_processor_assignment_type,
682 int migration_doMigration_type,
687 int mj_premigration_option,
688 int mj_premigration_coordinate_cutoff
692 Teuchos::ParameterList geoparams(
"geo params");
704 for(
int i = 0; i < coord_dim; ++i){
705 coords[i] =
new zscalar_t[numLocalPoints];
709 if (numWeightsPerCoord) {
710 weight=
new zscalar_t * [numWeightsPerCoord];
711 for(
int i = 0; i < numWeightsPerCoord; ++i){
712 weight[i] =
new zscalar_t[numLocalPoints];
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);
730 Teuchos::ArrayView<const zscalar_t> a;
734 RCP<tMVector_t> tmVector = RCP<tMVector_t>(
new
737 std::vector<const zscalar_t *>
weights;
738 if(numWeightsPerCoord){
739 for (
int i = 0; i < numWeightsPerCoord;++i){
740 weights.push_back(weight[i]);
743 std::vector<int> stride;
747 inputAdapter_t *ia =
new inputAdapter_t(tmVector, weights, stride);
749 Teuchos::RCP<Teuchos::ParameterList> params;
753 params = Teuchos::getParametersFromXmlFile(pfname);
756 params =RCP<Teuchos::ParameterList>(
new Teuchos::ParameterList,
true);
763 params->set(
"timer_output_stream" ,
"std::cout");
765 params->set(
"algorithm",
"multijagged");
767 params->set(
"mj_keep_part_boxes",
true);
769 params->set(
"rectilinear",
true );
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);
780 params->set(
"mj_parts", pqParts);
782 params->set(
"mj_num_teams", numTeams);
785 params->set(
"num_global_parts", numParts);
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));
807 ierr += compareWithBasicVectorAdapterTest<bv_use_node_t>(
808 comm, params, problem, tmVector,
809 weight, numWeightsPerCoord);
814 RCP<quality_t> metricObject =
817 if (comm->getRank() == 0){
818 metricObject->printMetrics(std::cout);
825 ierr += run_pointAssign_tests<inputAdapter_t>(problem, tmVector,
827 ierr += run_boxAssign_tests<inputAdapter_t>(problem, tmVector);
830 if(numWeightsPerCoord){
831 for(
int i = 0; i < numWeightsPerCoord; ++i) {
837 for(
int i = 0; i < coord_dim; ++i) {
848 template<
class bv_use_node_t>
850 RCP<
const Teuchos::Comm<int> > &comm,
858 int migration_check_option,
859 int migration_all_to_all_type,
861 int migration_processor_assignment_type,
862 int migration_doMigration_type,
867 int mj_premigration_option,
868 int mj_premigration_coordinate_cutoff
881 inputAdapter_t *ia =
new inputAdapter_t(coords);
883 Teuchos::RCP <Teuchos::ParameterList> params ;
887 params = Teuchos::getParametersFromXmlFile(pfname);
890 params =RCP <Teuchos::ParameterList> (
new Teuchos::ParameterList,
true);
895 params->set(
"mj_keep_part_boxes",
true);
897 params->set(
"rectilinear",
true);
898 params->set(
"algorithm",
"multijagged");
900 params->set(
"imbalance_tolerance",
double(imbalance));
904 params->set(
"mj_parts", pqParts);
906 params->set(
"mj_premigration_option", mj_premigration_option);
909 params->set(
"num_global_parts", numParts);
912 params->set(
"mj_concurrent_part_count", k);
914 if(migration_check_option >= 0){
915 params->set(
"mj_migration_option", migration_check_option);
917 if(migration_imbalance_cut_off >= 0){
918 params->set(
"mj_minimum_migration_imbalance",
919 double (migration_imbalance_cut_off));
921 if (mj_premigration_coordinate_cutoff > 0){
922 params->set(
"mj_premigration_coordinate_count",
923 mj_premigration_coordinate_cutoff);
938 ierr += compareWithBasicVectorAdapterTest<bv_use_node_t>(
939 comm, params, problem, coords);
944 RCP<quality_t> metricObject =
947 if (comm->getRank() == 0){
948 metricObject->printMetrics(std::cout);
949 std::cout <<
"testFromDataFile is done " << std::endl;
956 ierr += run_pointAssign_tests<inputAdapter_t>(problem, coords,
958 ierr += run_boxAssign_tests<inputAdapter_t>(problem, coords);
966 #ifdef hopper_separate_test
968 template <
typename zscalar_t,
typename zlno_t>
970 FILE *f = fopen(fileName.c_str(),
"r");
972 std::cout << fileName <<
" cannot be opened" << std::endl;
975 fscanf(f,
"%d", &numLocal);
976 fscanf(f,
"%d", &dim);
978 for (
int i = 0; i < dim; ++i){
981 for (
int i = 0; i < dim; ++i){
982 for (
zlno_t j = 0; j < numLocal; ++j){
983 fscanf(f,
"%lf", &(coords[i][j]));
987 std::cout <<
"done reading:" << fileName << std::endl;
990 int testFromSeparateDataFiles(
991 RCP<
const Teuchos::Comm<int> > &comm,
999 int migration_check_option,
1000 int migration_all_to_all_type,
1002 int migration_processor_assignment_type,
1003 int migration_doMigration_type,
1007 int mj_premigration_option
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";
1020 getCoords<zscalar_t, zlno_t>(double_coords, numLocal, dim, tFile);
1022 Teuchos::Array<Teuchos::ArrayView<const zscalar_t> > coordView(dim);
1023 for (
int i=0; i < dim; i++){
1025 Teuchos::ArrayView<const zscalar_t> a(double_coords[i], numLocal);
1028 Teuchos::ArrayView<const zscalar_t> a;
1035 Teuchos::Comm<int> *tcomm = (Teuchos::Comm<int> *)comm.getRawPtr();
1037 reduceAll<int, zgno_t>(
1039 Teuchos::REDUCE_SUM,
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));
1052 RCP<const tMVector_t> coordsConst = rcp_const_cast<
const tMVector_t>(coords);
1057 inputAdapter_t *ia =
new inputAdapter_t(coordsConst);
1061 Teuchos::RCP <Teuchos::ParameterList> params ;
1065 params = Teuchos::getParametersFromXmlFile(pfname);
1068 params =RCP <Teuchos::ParameterList> (
new Teuchos::ParameterList,
true);
1072 params->set(
"algorithm",
"multijagged");
1074 params->set(
"imbalance_tolerance",
double(imbalance));
1077 params->set(
"mj_premigration_option", mj_premigration_option);
1079 params->set(
"mj_parts", pqParts);
1082 params->set(
"mj_num_teams", numTeams);
1085 params->set(
"num_global_parts", numParts);
1088 params->set(
"parallel_part_calculation_count", k);
1090 if(migration_processor_assignment_type >= 0){
1091 params->set(
"migration_processor_assignment_type", migration_processor_assignment_type);
1093 if(migration_check_option >= 0){
1094 params->set(
"migration_check_option", migration_check_option);
1096 if(migration_all_to_all_type >= 0){
1097 params->set(
"migration_all_to_all_type", migration_all_to_all_type);
1099 if(migration_imbalance_cut_off >= 0){
1100 params->set(
"migration_imbalance_cut_off",
1101 double (migration_imbalance_cut_off));
1103 if (migration_doMigration_type >= 0){
1104 params->set(
"migration_doMigration_type",
int (migration_doMigration_type));
1107 params->set(
"mj_keep_part_boxes",
true);
1109 params->set(
"rectilinear",
true);
1125 ierr += compareWithBasicVectorAdapterTest<bv_use_node_t>(
1126 comm, params, problem, coords);
1129 if (coordsConst->getGlobalLength() < 40) {
1130 int len = coordsConst->getLocalLength();
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;
1141 RCP<quality_t> metricObject =
1144 if (comm->getRank() == 0){
1145 metricObject->printMetrics(std::cout);
1146 std::cout <<
"testFromDataFile is done " << std::endl;
1153 ierr += run_pointAssign_tests<inputAdapter_t>(problem, coords,
1155 ierr += run_boxAssign_tests<inputAdapter_t>(problem, coords);
1168 for(
int i = 0; args[i] != 0; i++)
1173 std::stringstream stream(std::stringstream::in | std::stringstream::out);
1174 stream << argumentline;
1175 getline(stream, argumentid,
'=');
1179 stream >> argumentValue;
1192 std::string &pfname,
1194 int &migration_check_option,
1195 int &migration_all_to_all_type,
1197 int &migration_processor_assignment_type,
1198 int &migration_doMigration_type,
1200 bool &print_details,
1203 int &mj_premigration_option,
1204 int &mj_coordinate_cutoff
1207 bool isCset =
false;
1208 bool isPset =
false;
1209 bool isFset =
false;
1210 bool isPFset =
false;
1212 for(
int i = 0; i < narg; ++i){
1214 string identifier =
"";
1215 long long int value = -1;
double fval = -1;
1217 value = (
long long int) (fval);
1219 if(identifier ==
"W"){
1220 if(value == 0 || value == 1){
1221 print_details = (value == 0 ?
false :
true);
1223 throw std::runtime_error(
"Invalid argument at " + tmp);
1226 else if(identifier ==
"UVM"){
1227 if(value == 0 || value == 1){
1228 uvm = (value == 0 ?
false :
true);
1230 throw std::runtime_error(
"Invalid argument at " + tmp);
1233 else if(identifier ==
"T"){
1237 throw std::runtime_error(
"Invalid argument at " + tmp);
1239 }
else if(identifier ==
"C"){
1244 throw std::runtime_error(
"Invalid argument at " + tmp);
1246 }
else if(identifier ==
"P"){
1247 std::stringstream stream(std::stringstream::in | std::stringstream::out);
1250 getline(stream, ttmp,
'=');
1253 }
else if(identifier ==
"I"){
1257 throw std::runtime_error(
"Invalid argument at " + tmp);
1259 }
else if(identifier ==
"MI"){
1261 migration_imbalance_cut_off=fval;
1263 throw std::runtime_error(
"Invalid argument at " + tmp);
1265 }
else if(identifier ==
"MO"){
1267 migration_check_option = value;
1269 throw std::runtime_error(
"Invalid argument at " + tmp);
1271 }
else if(identifier ==
"AT"){
1273 migration_processor_assignment_type = value;
1275 throw std::runtime_error(
"Invalid argument at " + tmp);
1279 else if(identifier ==
"MT"){
1281 migration_all_to_all_type = value;
1283 throw std::runtime_error(
"Invalid argument at " + tmp);
1286 else if(identifier ==
"PCC"){
1288 mj_coordinate_cutoff = value;
1290 throw std::runtime_error(
"Invalid argument at " + tmp);
1294 else if(identifier ==
"PM"){
1296 mj_premigration_option = value;
1298 throw std::runtime_error(
"Invalid argument at " + tmp);
1302 else if(identifier ==
"DM"){
1304 migration_doMigration_type = value;
1306 throw std::runtime_error(
"Invalid argument at " + tmp);
1309 else if(identifier ==
"F"){
1310 std::stringstream stream(std::stringstream::in | std::stringstream::out);
1312 getline(stream, fname,
'=');
1317 else if(identifier ==
"PF"){
1318 std::stringstream stream(std::stringstream::in | std::stringstream::out);
1320 getline(stream, pfname,
'=');
1326 else if(identifier ==
"O"){
1327 if(value >= 0 && value <= 3){
1330 throw std::runtime_error(
"Invalid argument at " + tmp);
1333 else if(identifier ==
"K"){
1337 throw std::runtime_error(
"Invalid argument at " + tmp);
1340 else if(identifier ==
"TB"){
1342 test_boxes = (value == 0 ?
false :
true);
1344 throw std::runtime_error(
"Invalid argument at " + tmp);
1347 else if(identifier ==
"R"){
1349 rectilinear = (value == 0 ?
false :
true);
1351 throw std::runtime_error(
"Invalid argument at " + tmp);
1355 throw std::runtime_error(
"Invalid argument at " + tmp);
1359 if(!( (isCset || isPset || isPFset) && isFset)){
1360 throw std::runtime_error(
"(C || P || PF) && F are mandatory arguments.");
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;
1385 Tpetra::ScopeGuard tscope(&narg, &arg);
1386 Teuchos::RCP<const Teuchos::Comm<int> > tcomm = Tpetra::getDefaultComm();
1388 int rank = tcomm->getRank();
1393 float imbalance = -1.03;
1396 string pqParts =
"";
1398 std::string fname =
"";
1399 std::string paramFile =
"";
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;
1411 bool print_details =
true;
1412 bool test_boxes =
false;
1413 bool rectilinear =
false;
1422 #if defined(KOKKOS_ENABLE_CUDA)
1423 using uvm_off_node_t = Tpetra::KokkosCompat::KokkosDeviceWrapperNode<
1424 Kokkos::Cuda, Kokkos::CudaSpace>;
1425 #elif defined(KOKKOS_ENABLE_HIP)
1426 using uvm_off_node_t = Tpetra::KokkosCompat::KokkosDeviceWrapperNode<
1427 Kokkos::HIP, Kokkos::HIPSpace>;
1443 migration_check_option,
1444 migration_all_to_all_type,
1445 migration_imbalance_cut_off,
1446 migration_processor_assignment_type,
1447 migration_doMigration_type,
1452 mj_premigration_option, mj_premigration_coordinate_cutoff);
1454 catch(std::exception &s){
1455 if(tcomm->getRank() == 0){
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,
1475 mj_premigration_option, mj_premigration_coordinate_cutoff);
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,
1487 mj_premigration_option, mj_premigration_coordinate_cutoff);
1489 throw std::logic_error(
"uvm set off but this is not a cuda/hip test.");
1494 #ifdef hopper_separate_test
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,
1506 mj_premigration_option, mj_premigration_coordinate_cutoff);
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,
1519 mj_premigration_option, mj_premigration_coordinate_cutoff);
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,
1531 mj_premigration_option, mj_premigration_coordinate_cutoff);
1533 throw std::logic_error(
"uvm set off but this is not a cuda test.");
1540 if (ierr == 0) std::cout <<
"PASS" << std::endl;
1541 else std::cout <<
"FAIL" << std::endl;
1546 catch(std::string &s){
1548 std::cerr << s << std::endl;
1553 std::cerr << s << std::endl;
void getCoords(void *data, int numGid, int numLid, int numObj, ZOLTAN_ID_PTR gids, ZOLTAN_ID_PTR lids, int dim, double *coords_, int *ierr)
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)
#define CATCH_EXCEPTIONS_WITH_COUNT(ierr, pp)
void getLocalCoordinatesCopy(scalar_t **c)
typename InputTraits< User >::scalar_t scalar_t
string trim_right_copy(const string &s, const string &delimiters=" \f\n\r\t\v")
void print_usage(char *executable)
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)
void readGeoGenParams(string paramFileName, Teuchos::ParameterList &geoparams, const RCP< const Teuchos::Comm< int > > &comm)
map_t::global_ordinal_type gno_t
A simple class that can be the User template argument for an InputAdapter.
int getCoordinateDimension()
Zoltan2::EvaluatePartition< matrixAdapter_t > quality_t
int main(int narg, char **arg)
Defines the PartitioningSolution class.
common code used by tests
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)
Defines the XpetraMultiVectorAdapter.
string trim_copy(const string &s, const string &delimiters=" \f\n\r\t\v")
SparseMatrixAdapter_t::part_t part_t
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)
Defines the EvaluatePartition class.
int run_boxAssign_tests(Zoltan2::PartitioningProblem< Adapter > *problem, RCP< tMVector_t > &coords)
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...
map_t::local_ordinal_type lno_t
void printTimers() const
Return the communicator passed to the problem.
An adapter for Xpetra::MultiVector.
#define CATCH_EXCEPTIONS_AND_RETURN(pp)
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)
Tpetra::Map::local_ordinal_type zlno_t
int run_pointAssign_tests(Zoltan2::PartitioningProblem< Adapter > *problem, RCP< tMVector_t > &coords, bool print_details)
static const std::string fail
string convert_to_string(char *args)
const PartitioningSolution< Adapter > & getSolution()
Get the solution to the problem.
PartitioningProblem sets up partitioning problems for the user.
string trim_left_copy(const string &s, const string &delimiters=" \f\n\r\t\v")
int getNumWeights()
##END Predistribution functions######################//
bool getArgumentValue(string &argumentid, double &argumentValue, string argumentline)
lno_t getNumLocalCoords()
Defines the PartitioningProblem class.
A class that computes and returns quality metrics.
gno_t getNumGlobalCoords()
Defines the BasicVectorAdapter class.
Tpetra::MultiVector< zscalar_t, zlno_t, zgno_t, znode_t > tMVector_t
Tpetra::Map::global_ordinal_type zgno_t
void solve(bool updateInputData=true)
Direct the problem to create a solution.
void getLocalWeightsCopy(scalar_t **w)
std::string testDataFilePath(".")