33 #include <Teuchos_DefaultComm.hpp>
34 #include <Teuchos_XMLObject.hpp>
35 #include <Teuchos_FileInputSource.hpp>
36 #include <Teuchos_StackedTimer.hpp>
44 using Teuchos::ParameterList;
47 using Teuchos::ArrayRCP;
48 using Teuchos::XMLObject;
49 using namespace Zoltan2_TestingFramework;
55 using std::ostringstream;
58 #define ERRMSG(msg) if (rank == 0){ std::cerr << "FAIL: " << msg << std::endl; }
59 #define EXC_ERRMSG(msg, e) \
60 if (rank==0){ std::cerr << "FAIL: " << msg << std::endl << e.what() << std::endl;}
63 Teuchos::ParameterList & plist)
66 Teuchos::XMLParameterListReader reader;
67 plist = reader.toParameterList(xml);
75 Teuchos::ParameterList zoltan2Parameters;
78 if (plist.isSublist(
"Zoltan2Parameters")) {
80 ParameterList &sub = plist.sublist(
"Zoltan2Parameters");
81 zoltan2Parameters.setParameters(sub);
86 queue<ParameterList> &problems,
87 queue<ParameterList> &comparisons,
88 const RCP<
const Teuchos::Comm<int> > & comm)
90 int rank = comm->getRank();
93 Teuchos::FileInputSource inputSource(inputFileName);
95 std::cout <<
"Input file source: " << inputFileName << std::endl;
102 xmlInput = inputSource.getObject();
111 for(
int i = 0; i < xmlInput.numChildren(); i++)
116 if(plist.name() ==
"Comparison") {
117 comparisons.emplace(plist);
120 problems.emplace(plist);
129 std::ostringstream & msg,
130 const ParameterList &problem_parameters) {
131 #define ANALYZE_METRICS(adapterClass, metricAnalyzerClass) \
132 RCP<EvaluateBaseClass<adapterClass>> pCast = \
133 rcp_dynamic_cast<EvaluateBaseClass<adapterClass>>( \
134 evaluateFactory->getEvaluateClass()); \
135 if(pCast == Teuchos::null) throw std::logic_error( \
136 "Bad evaluate class cast in analyzeMetrics!" ); \
137 metricAnalyzerClass analyzer(pCast); \
138 return analyzer.analyzeMetrics( \
139 problem_parameters.sublist("Metrics"), msg);
141 #define ANALYZE_METRICS_PARTITIONING(adapterClass) \
142 ANALYZE_METRICS(adapterClass, \
143 MetricAnalyzerEvaluatePartition<adapterClass>)
145 #define ANALYZE_METRICS_ORDERING(adapterClass) \
146 ANALYZE_METRICS(adapterClass, \
147 MetricAnalyzerEvaluateOrdering<adapterClass>)
149 if(evaluateFactory->getProblemName() ==
"partitioning") {
152 else if(evaluateFactory->getProblemName() ==
"ordering") {
156 throw std::logic_error(
157 "analyzeMetrics not implemented for this problem type!" );
163 RCP<ProblemFactory> problemFactory) {
164 #define GET_LOCAL_ORDERING(adapterClass) \
165 return (rcp_dynamic_cast<OrderingProblem<adapterClass>>( \
166 problemFactory->getProblem()))->getLocalOrderingSolution();
172 #define GET_PROBLEM_PARTS(adapterClass) \
173 return (rcp_dynamic_cast<PartitioningProblem<adapterClass>>( \
174 problemFactory->getProblem()))->getSolution().getPartListView();
180 #define GET_IDS_VIEW(adapterClass) \
181 return dynamic_cast<adapterClass*>( \
182 adapterFactory->getMainAdapter())->getIDsView(Ids);
184 throw std::logic_error(
"getIDsView() failed to match adapter name" );
188 const ParameterList &problem_parameters,
189 bool bHasComparisons,
190 RCP<ComparisonHelper> & comparison_helper,
191 const RCP<
const Teuchos::Comm<int> > & comm)
200 int rank = comm->getRank();
201 if(!problem_parameters.isParameter(
"kind"))
204 std::cout <<
"Problem kind not provided" << std::endl;
208 if(!problem_parameters.isParameter(
"InputAdapterParameters"))
211 std::cout <<
"Input adapter parameters not provided" << std::endl;
215 if(!problem_parameters.isParameter(
"Zoltan2Parameters"))
218 std::cout <<
"Zoltan2 problem parameters not provided" << std::endl;
224 std::cout <<
"\n\nRunning test: " << problem_parameters.name() << std::endl;
230 auto comparison_source = comparison_helper->AddSource(problem_parameters.name());
232 comparison_source->startBaseTimer();
237 const ParameterList &adapterPlist =
238 problem_parameters.sublist(
"InputAdapterParameters");
239 comparison_source->startTimer(
"adapter construction time");
243 const_cast<UserInputForTests*>(&uinput), adapterPlist, comm));
245 comparison_source->stopTimer(
"adapter construction time");
247 comparison_source->adapterFactory = adapterFactory;
253 string adapter_name = adapterPlist.get<
string>(
"input adapter");
255 ParameterList zoltan2_parameters =
256 const_cast<ParameterList &
>(problem_parameters.sublist(
"Zoltan2Parameters"));
258 std::cout << std::endl;
261 comparison_source->startTimer(
"problem construction time");
262 std::string problem_kind = problem_parameters.get<std::string>(
"kind");
264 std::cout <<
"Creating a new " << problem_kind <<
" problem." << std::endl;
271 #ifdef HAVE_ZOLTAN2_MPI
277 std::cout <<
"Using input adapter type: " + adapter_name << std::endl;
280 comparison_source->stopTimer(
"problem construction time");
282 comparison_source->problemFactory = problemFactory;
287 comparison_source->startTimer(
"solve time");
289 problemFactory->getProblem()->solve();
291 comparison_source->stopTimer(
"solve time");
293 std::cout << problem_kind +
" problem solved." << std::endl;
298 if(problem_kind ==
"partitioning") {
302 i < adapterFactory->getMainAdapter()->getLocalNumIDs(); i++) {
303 std::cout << rank <<
" LID " << i
304 <<
" GID " << kddIDs[i]
310 if (adapter_name ==
"XpetraCrsGraph") {
315 dynamic_cast<const xCG_xCG_t *
>(adapterFactory->getMainAdapter());
319 const gno_t *vertexIds;
321 const offset_t *offsets;
324 for (
int edim = 0; edim < ewgtDim; edim++) {
328 for (lno_t i=0; i < localNumObj; i++)
329 for (offset_t j=offsets[i]; j < offsets[i+1]; j++)
330 std::cout << edim <<
" " << vertexIds[i] <<
" "
331 << adjIds[j] <<
" " << weights[stride*j] << std::endl;
339 bool bSuccess =
true;
340 if(problem_parameters.isSublist(
"Metrics") || bHasComparisons) {
348 std::cout <<
"Create evaluate class for: " + problem_kind << std::endl;
352 comparison_source->evaluateFactory = evaluateFactory;
354 std::ostringstream msgSummary;
356 evaluateFactory->getEvaluateClass()->printMetrics(msgSummary);
358 std::cout << msgSummary.str();
361 std::ostringstream msgResults;
362 if (!
analyzeMetrics(evaluateFactory, msgResults, problem_parameters)) {
364 std::cout <<
"MetricAnalyzer::analyzeMetrics() "
365 <<
"returned false and the test is FAILED." << std::endl;
368 std::cout << msgResults.str();
373 if (problem_kind ==
"ordering") {
374 std::cout <<
"\nLet's examine the solution..." << std::endl;
378 std::cout <<
"Number of column blocks: "
381 if (localOrderingSolution->
havePerm()) {
382 std::cout <<
"permutation: {";
384 std::cout <<
" " << x;
385 std::cout <<
"}" << std::endl;
389 std::cout <<
"inverse permutation: {";
391 std::cout <<
" " << x;
392 std::cout <<
"}" << std::endl;
396 std::cout <<
"separator range: {";
398 std::cout <<
" " << x;
399 std::cout <<
"}" << std::endl;
403 std::cout <<
"separator tree: {";
405 std::cout <<
" " << x;
406 std::cout <<
"}" << std::endl;
413 comparison_source->stopBaseTimer();
425 bool mainExecute(
int narg,
char *arg[], RCP<
const Comm<int> > &comm)
427 Teuchos::RCP<Teuchos::StackedTimer> stacked_timer = Teuchos::rcp(
new Teuchos::StackedTimer(
"Zoltan2_Driver"));;
432 int rank = comm->getRank();
438 string inputFileName(
"");
440 inputFileName = arg[1];
443 std::cout <<
"\nFAILED to specify xml input file!" << std::endl;
444 std::ostringstream msg;
445 msg <<
"\nStandard use of test_driver.cpp:\n";
446 msg <<
"mpiexec -n <procs> ./Zoltan2_test_driver.exe <input_file.xml>\n";
447 std::cout << msg.str() << std::endl;
455 queue<ParameterList> problems, comparisons;
465 const ParameterList inputParameters = problems.front();
466 if(inputParameters.name() !=
"InputParameters")
469 std::cout <<
"InputParameters not defined. Testing FAILED." << std::endl;
489 RCP<ComparisonHelper> comparison_manager = rcp(
new ComparisonHelper(stacked_timer));
490 while (!problems.empty()) {
494 if (!
run(uinput, problems.front(), !comparisons.empty(),
495 comparison_manager, comm)) {
496 std::cout <<
"Problem run returned false" << std::endl;
507 while (!comparisons.empty()) {
508 if (!comparison_manager->Compare(comparisons.front(),comm)) {
510 std::cout <<
"Comparison manager returned false so the "
511 <<
"test should fail." << std::endl;
520 std::cout <<
"\nFAILED to load input data source. Skipping "
521 "all tests." << std::endl;
526 stacked_timer->stopBaseTimer();
527 Teuchos::StackedTimer::OutputOptions options;
528 options.output_fraction = options.output_histogram = options.output_minmax =
true;
529 stacked_timer->report(std::cout, comm, options);
530 std::string watchrProblemName = std::string(
"Zoltan2 test driver ") + std::to_string(comm->getSize()) +
" ranks";
531 auto xmlOut = stacked_timer->reportWatchrXML(watchrProblemName, comm);
533 std::cout <<
"\nAlso created Watchr performance report " << xmlOut <<
'\n';
538 int main(
int narg,
char *arg[])
540 Tpetra::ScopeGuard tscope(&narg, &arg);
541 Teuchos::RCP<const Teuchos::Comm<int> > comm = Tpetra::getDefaultComm();
544 int rank = comm->getRank();
549 catch(std::logic_error &e) {
554 std::cout <<
"Test driver for rank " << rank
555 <<
" caught the following exception: " << e.what() << std::endl;
559 catch(std::runtime_error &e) {
560 std::cout <<
"Test driver for rank " << rank
561 <<
" caught the following exception: " << e.what() << std::endl;
564 catch(std::exception &e) {
565 std::cout <<
"Test driver for rank " << rank
566 <<
" caught the following exception: " << e.what() << std::endl;
576 reduceAll<int,int>(*comm, Teuchos::EReductionType::REDUCE_MAX, 1,
577 &result, &resultReduced);
582 std::cout <<
"Test Driver with " << comm->getSize()
583 <<
" processes has reduced result code " << resultReduced
584 <<
" and is exiting in the "
585 << ((resultReduced == 0 ) ?
"PASSED" :
"FAILED") <<
" state."
keep typedefs that commonly appear in many places localized
#define GET_LOCAL_ORDERING(adapterClass)
const zpart_t * getPartListView(RCP< ProblemFactory > problemFactory)
LocalOrderingSolution< zlno_t > * getLocalOrderingSolution(RCP< ProblemFactory > problemFactory)
#define ANALYZE_METRICS_PARTITIONING(adapterClass)
bool getParameterLists(const string &inputFileName, queue< ParameterList > &problems, queue< ParameterList > &comparisons, const RCP< const Teuchos::Comm< int > > &comm)
typename InputTraits< User >::scalar_t scalar_t
ProblemFactory class contains 1 static factory method.
Defines Parameter related enumerators, declares functions.
ArrayRCP< lno_t > & getPermutationRCPConst(bool inverse=false) const
Get (local) permuted GIDs by const RCP.
Returns a pointer to new test classes. Is not responsible for memory management!
ArrayRCP< lno_t > & getSeparatorRangeRCPConst() const
Get separator range by const RCP.
Store and compare solution sets from different problems.
map_t::global_ordinal_type gno_t
Provides access for Zoltan2 to Xpetra::CrsGraph data.
bool havePerm() const
Do we have the direct permutation?
int main(int narg, char **arg)
size_t getLocalNumVertices() const override
Returns the number of vertices on this process.
void createAllParameters(Teuchos::ParameterList &pList)
Create a list of all Zoltan2 parameters and validators.
bool haveSeparatorRange() const
Do we have the separator range?
ProblemFactory class contains 1 static factory method.
Defines the XpetraMultiVectorAdapter.
bool analyzeMetrics(RCP< EvaluateFactory > evaluateFactory, std::ostringstream &msg, const ParameterList &problem_parameters)
Defines XpetraCrsGraphAdapter class.
Defines the XpetraCrsMatrixAdapter class.
#define Z2_TEST_UPCAST(adptr, TEMPLATE_ACTION)
int getNumWeightsPerEdge() const override
Returns the number (0 or greater) of edge weights.
lno_t getNumSeparatorBlocks() const
Get number of separator column blocks.
typename InputTraits< User >::gno_t gno_t
map_t::local_ordinal_type lno_t
bool haveSeparators() const
Do we have the separators?
A class for comparing solutions, metrics, and timing data of Zoltan2 problems.
Defines the BasicIdentifierAdapter class.
void getVertexIDsView(const gno_t *&ids) const override
void xmlToModelPList(const Teuchos::XMLObject &xml, Teuchos::ParameterList &plist)
void getEdgeWeightsView(const scalar_t *&weights, int &stride, int idx) const override
Provide a pointer to the edge weights, if any.
bool haveSeparatorTree() const
Do we have the separator tree?
void getIDsView(RCP< AdapterFactory > adapterFactory, const zgno_t *&Ids)
int run(const RCP< const Comm< int > > &comm, int numGlobalParts, int testCnt, std::string *thisTest)
void getEdgesView(const offset_t *&offsets, const gno_t *&adjIds) const override
#define ANALYZE_METRICS_ORDERING(adapterClass)
bool mainExecute(int narg, char *arg[], RCP< const Comm< int > > &comm)
#define EXC_ERRMSG(msg, e)
#define GET_IDS_VIEW(adapterClass)
Tpetra::Map::global_ordinal_type zgno_t
typename InputTraits< User >::lno_t lno_t
ArrayRCP< lno_t > & getSeparatorTreeRCPConst() const
Get separator tree by const RCP.
#define GET_PROBLEM_PARTS(adapterClass)
Generate Adapter for testing purposes.
bool haveInverse() const
Do we have the inverse permutation?