69 #include <Teuchos_DefaultComm.hpp>
70 #include <Teuchos_XMLObject.hpp>
71 #include <Teuchos_FileInputSource.hpp>
79 using Teuchos::ParameterList;
82 using Teuchos::ArrayRCP;
83 using Teuchos::XMLObject;
84 using namespace Zoltan2_TestingFramework;
90 using std::ostringstream;
93 #define ERRMSG(msg) if (rank == 0){ std::cerr << "FAIL: " << msg << std::endl; }
94 #define EXC_ERRMSG(msg, e) \
95 if (rank==0){ std::cerr << "FAIL: " << msg << std::endl << e.what() << std::endl;}
98 Teuchos::ParameterList & plist)
101 Teuchos::XMLParameterListReader reader;
102 plist = reader.toParameterList(xml);
110 Teuchos::ParameterList zoltan2Parameters;
113 if (plist.isSublist(
"Zoltan2Parameters")) {
115 ParameterList &sub = plist.sublist(
"Zoltan2Parameters");
116 zoltan2Parameters.setParameters(sub);
121 queue<ParameterList> &problems,
122 queue<ParameterList> &comparisons,
123 const RCP<
const Teuchos::Comm<int> > & comm)
125 int rank = comm->getRank();
128 Teuchos::FileInputSource inputSource(inputFileName);
130 std::cout <<
"Input file source: " << inputFileName << std::endl;
137 xmlInput = inputSource.getObject();
146 for(
int i = 0; i < xmlInput.numChildren(); i++)
151 if(plist.name() ==
"Comparison") {
152 comparisons.emplace(plist);
155 problems.emplace(plist);
164 std::ostringstream & msg,
165 const ParameterList &problem_parameters) {
166 #define ANALYZE_METRICS(adapterClass, metricAnalyzerClass) \
167 RCP<EvaluateBaseClass<adapterClass>> pCast = \
168 rcp_dynamic_cast<EvaluateBaseClass<adapterClass>>( \
169 evaluateFactory->getEvaluateClass()); \
170 if(pCast == Teuchos::null) throw std::logic_error( \
171 "Bad evaluate class cast in analyzeMetrics!" ); \
172 metricAnalyzerClass analyzer(pCast); \
173 return analyzer.analyzeMetrics( \
174 problem_parameters.sublist("Metrics"), msg);
176 #define ANALYZE_METRICS_PARTITIONING(adapterClass) \
177 ANALYZE_METRICS(adapterClass, \
178 MetricAnalyzerEvaluatePartition<adapterClass>)
180 #define ANALYZE_METRICS_ORDERING(adapterClass) \
181 ANALYZE_METRICS(adapterClass, \
182 MetricAnalyzerEvaluateOrdering<adapterClass>)
184 if(evaluateFactory->getProblemName() ==
"partitioning") {
187 else if(evaluateFactory->getProblemName() ==
"ordering") {
191 throw std::logic_error(
192 "analyzeMetrics not implemented for this problem type!" );
198 RCP<ProblemFactory> problemFactory) {
199 #define GET_LOCAL_ORDERING(adapterClass) \
200 return (rcp_dynamic_cast<OrderingProblem<adapterClass>>( \
201 problemFactory->getProblem()))->getLocalOrderingSolution();
207 #define GET_PROBLEM_PARTS(adapterClass) \
208 return (rcp_dynamic_cast<PartitioningProblem<adapterClass>>( \
209 problemFactory->getProblem()))->getSolution().getPartListView();
215 #define GET_IDS_VIEW(adapterClass) \
216 return dynamic_cast<adapterClass*>( \
217 adapterFactory->getMainAdapter())->getIDsView(Ids);
219 throw std::logic_error(
"getIDsView() failed to match adapter name" );
223 const ParameterList &problem_parameters,
224 bool bHasComparisons,
225 RCP<ComparisonHelper> & comparison_helper,
226 const RCP<
const Teuchos::Comm<int> > & comm)
235 int rank = comm->getRank();
236 if(!problem_parameters.isParameter(
"kind"))
239 std::cout <<
"Problem kind not provided" << std::endl;
243 if(!problem_parameters.isParameter(
"InputAdapterParameters"))
246 std::cout <<
"Input adapter parameters not provided" << std::endl;
250 if(!problem_parameters.isParameter(
"Zoltan2Parameters"))
253 std::cout <<
"Zoltan2 problem parameters not provided" << std::endl;
259 std::cout <<
"\n\nRunning test: " << problem_parameters.name() << std::endl;
267 comparison_helper->AddSource(problem_parameters.name(), comparison_source);
268 comparison_source->addTimer(
"adapter construction time");
269 comparison_source->addTimer(
"problem construction time");
270 comparison_source->addTimer(
"solve time");
275 const ParameterList &adapterPlist =
276 problem_parameters.sublist(
"InputAdapterParameters");
277 comparison_source->timers[
"adapter construction time"]->start();
281 const_cast<UserInputForTests*>(&uinput), adapterPlist, comm));
283 comparison_source->timers[
"adapter construction time"]->stop();
285 comparison_source->adapterFactory = adapterFactory;
291 string adapter_name = adapterPlist.get<
string>(
"input adapter");
293 ParameterList zoltan2_parameters =
294 const_cast<ParameterList &
>(problem_parameters.sublist(
"Zoltan2Parameters"));
296 std::cout << std::endl;
299 comparison_source->timers[
"problem construction time"]->start();
300 std::string problem_kind = problem_parameters.get<std::string>(
"kind");
302 std::cout <<
"Creating a new " << problem_kind <<
" problem." << std::endl;
309 #ifdef HAVE_ZOLTAN2_MPI
315 std::cout <<
"Using input adapter type: " + adapter_name << std::endl;
318 comparison_source->problemFactory = problemFactory;
323 comparison_source->timers[
"solve time"]->start();
325 problemFactory->getProblem()->solve();
327 comparison_source->timers[
"solve time"]->stop();
329 std::cout << problem_kind +
" problem solved." << std::endl;
334 if(problem_kind ==
"partitioning") {
338 i < adapterFactory->getMainAdapter()->getLocalNumIDs(); i++) {
339 std::cout << rank <<
" LID " << i
340 <<
" GID " << kddIDs[i]
346 if (adapter_name ==
"XpetraCrsGraph") {
351 dynamic_cast<const xCG_xCG_t *
>(adapterFactory->getMainAdapter());
355 const gno_t *vertexIds;
357 const offset_t *offsets;
360 for (
int edim = 0; edim < ewgtDim; edim++) {
364 for (lno_t i=0; i < localNumObj; i++)
365 for (offset_t j=offsets[i]; j < offsets[i+1]; j++)
366 std::cout << edim <<
" " << vertexIds[i] <<
" "
367 << adjIds[j] <<
" " << weights[stride*j] << std::endl;
375 bool bSuccess =
true;
376 if(problem_parameters.isSublist(
"Metrics") || bHasComparisons) {
384 std::cout <<
"Create evaluate class for: " + problem_kind << std::endl;
388 comparison_source->evaluateFactory = evaluateFactory;
390 std::ostringstream msgSummary;
392 evaluateFactory->getEvaluateClass()->printMetrics(msgSummary);
394 std::cout << msgSummary.str();
397 std::ostringstream msgResults;
398 if (!
analyzeMetrics(evaluateFactory, msgResults, problem_parameters)) {
400 std::cout <<
"MetricAnalyzer::analyzeMetrics() "
401 <<
"returned false and the test is FAILED." << std::endl;
404 std::cout << msgResults.str();
409 if (problem_kind ==
"ordering") {
410 std::cout <<
"\nLet's examine the solution..." << std::endl;
414 std::cout <<
"Number of column blocks: "
417 if (localOrderingSolution->
havePerm()) {
418 std::cout <<
"permutation: {";
420 std::cout <<
" " << x;
421 std::cout <<
"}" << std::endl;
425 std::cout <<
"inverse permutation: {";
427 std::cout <<
" " << x;
428 std::cout <<
"}" << std::endl;
432 std::cout <<
"separator range: {";
434 std::cout <<
" " << x;
435 std::cout <<
"}" << std::endl;
439 std::cout <<
"separator tree: {";
441 std::cout <<
" " << x;
442 std::cout <<
"}" << std::endl;
449 comparison_source->printTimers();
461 bool mainExecute(
int narg,
char *arg[], RCP<
const Comm<int> > &comm)
466 int rank = comm->getRank();
472 string inputFileName(
"");
474 inputFileName = arg[1];
477 std::cout <<
"\nFAILED to specify xml input file!" << std::endl;
478 std::ostringstream msg;
479 msg <<
"\nStandard use of test_driver.cpp:\n";
480 msg <<
"mpiexec -n <procs> ./Zoltan2_test_driver.exe <input_file.xml>\n";
481 std::cout << msg.str() << std::endl;
489 queue<ParameterList> problems, comparisons;
499 const ParameterList inputParameters = problems.front();
500 if(inputParameters.name() !=
"InputParameters")
503 std::cout <<
"InputParameters not defined. Testing FAILED." << std::endl;
524 while (!problems.empty()) {
528 if (!
run(uinput, problems.front(), !comparisons.empty(),
529 comparison_manager, comm)) {
530 std::cout <<
"Problem run returned false" << std::endl;
541 while (!comparisons.empty()) {
542 if (!comparison_manager->Compare(comparisons.front(),comm)) {
544 std::cout <<
"Comparison manager returned false so the "
545 <<
"test should fail." << std::endl;
554 std::cout <<
"\nFAILED to load input data source. Skipping "
555 "all tests." << std::endl;
563 int main(
int narg,
char *arg[])
565 Tpetra::ScopeGuard tscope(&narg, &arg);
566 Teuchos::RCP<const Teuchos::Comm<int> > comm = Tpetra::getDefaultComm();
569 int rank = comm->getRank();
574 catch(std::logic_error &e) {
579 std::cout <<
"Test driver for rank " << rank
580 <<
" caught the following exception: " << e.what() << std::endl;
584 catch(std::runtime_error &e) {
585 std::cout <<
"Test driver for rank " << rank
586 <<
" caught the following exception: " << e.what() << std::endl;
589 catch(std::exception &e) {
590 std::cout <<
"Test driver for rank " << rank
591 <<
" caught the following exception: " << e.what() << std::endl;
601 reduceAll<int,int>(*comm, Teuchos::EReductionType::REDUCE_MAX, 1,
602 &result, &resultReduced);
607 std::cout <<
"Test Driver with " << comm->getSize()
608 <<
" processes has reduced result code " << resultReduced
609 <<
" and is exiting in the "
610 << ((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
A class used to save problem solutions and timers.
#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?