diff --git a/packages/amesos2/cmake/Amesos2_config.h.in b/packages/amesos2/cmake/Amesos2_config.h.in index 4a6e76ce8f67..273f9f461c57 100644 --- a/packages/amesos2/cmake/Amesos2_config.h.in +++ b/packages/amesos2/cmake/Amesos2_config.h.in @@ -116,3 +116,4 @@ #ifdef HAVE_AMESOS2_ZOLTAN2CORE # define HAVE_AMESOS2_ZOLTAN2 #endif +#cmakedefine HAVE_AMESOS2_GALERI diff --git a/packages/amesos2/cmake/Dependencies.cmake b/packages/amesos2/cmake/Dependencies.cmake index 9c5d753ee8ee..b8bc60738faf 100644 --- a/packages/amesos2/cmake/Dependencies.cmake +++ b/packages/amesos2/cmake/Dependencies.cmake @@ -5,7 +5,7 @@ SET(LIB_REQUIRED_DEP_PACKAGES Teuchos Tpetra TrilinosSS Kokkos) SET(LIB_OPTIONAL_DEP_PACKAGES Epetra EpetraExt ShyLU_NodeBasker ShyLU_NodeTacho) SET(TEST_REQUIRED_DEP_PACKAGES) -SET(TEST_OPTIONAL_DEP_PACKAGES ShyLU_NodeBasker ShyLU_NodeTacho Kokkos TrilinosSS Xpetra Zoltan2Core) +SET(TEST_OPTIONAL_DEP_PACKAGES ShyLU_NodeBasker ShyLU_NodeTacho Kokkos TrilinosSS Xpetra Zoltan2Core Galeri) # SET(LIB_REQUIRED_DEP_TPLS SuperLU) SET(LIB_REQUIRED_DEP_TPLS ) SET(LIB_OPTIONAL_DEP_TPLS MPI SuperLU SuperLUMT SuperLUDist LAPACK UMFPACK PARDISO_MKL CSS_MKL ParMETIS METIS Cholmod MUMPS STRUMPACK CUSPARSE CUSOLVER) diff --git a/packages/amesos2/example/SimpleSolve.cpp b/packages/amesos2/example/SimpleSolve.cpp index a1b3472ad7d8..7a722431aa66 100644 --- a/packages/amesos2/example/SimpleSolve.cpp +++ b/packages/amesos2/example/SimpleSolve.cpp @@ -23,6 +23,9 @@ #include #include #include +#include +#include +#include #include #include @@ -31,7 +34,10 @@ #include "Amesos2.hpp" #include "Amesos2_Version.hpp" - +#if defined(HAVE_AMESOS2_XPETRA) && defined(HAVE_AMESOS2_GALERI) + #include "Galeri_XpetraMaps.hpp" + #include "Galeri_XpetraProblemFactory.hpp" +#endif int main(int argc, char *argv[]) { Tpetra::ScopeGuard tpetraScope(&argc,&argv); @@ -40,6 +46,7 @@ int main(int argc, char *argv[]) { typedef Tpetra::Map<>::local_ordinal_type LO; typedef Tpetra::Map<>::global_ordinal_type GO; + typedef Tpetra::Map MAP; typedef Tpetra::CrsMatrix MAT; typedef Tpetra::MultiVector MV; @@ -48,9 +55,24 @@ int main(int argc, char *argv[]) { using Teuchos::RCP; using Teuchos::rcp; + bool verbose = false; + GO nx = 1; + std::string GaleriName3D {"Laplace3D"}; + std::string solvername("Superlu"); + std::string xml_filename(""); + Teuchos::CommandLineProcessor cmdp(false,true); + cmdp.setOption("verbose","quiet",&verbose,"Print messages and results."); + cmdp.setOption("solvername",&solvername,"Name of solver."); + cmdp.setOption("xml_filename",&xml_filename,"XML Filename for Solver parameters."); + cmdp.setOption("nx",&nx,"Dimension of 3D problem."); + cmdp.setOption ("galeriMatrixName", &GaleriName3D, "Name of 3D Galeri Matrix"); + if (cmdp.parse(argc,argv) != Teuchos::CommandLineProcessor::PARSE_SUCCESSFUL) { + return -1; + } + // Before we do anything, check that SuperLU is enabled - if( !Amesos2::query("SuperLU") ){ - std::cerr << "SuperLU not enabled. Exiting..." << std::endl; + if( !Amesos2::query(solvername) ){ + std::cerr << solvername << " not enabled. Exiting..." << std::endl; return EXIT_SUCCESS; // Otherwise CTest will pick it up as // failure, which it isn't really } @@ -66,82 +88,115 @@ int main(int argc, char *argv[]) { const size_t numVectors = 1; - // create a Map - global_size_t nrows = 6; - RCP > map - = rcp( new Tpetra::Map(nrows,0,comm) ); - - RCP A = rcp( new MAT(map,3) ); // max of three entries in a row - - /* - * We will solve a system with a known solution, for which we will be using - * the following matrix: - * - * [ [ 7, 0, -3, 0, -1, 0 ] - * [ 2, 8, 0, 0, 0, 0 ] - * [ 0, 0, 1, 0, 0, 0 ] - * [ -3, 0, 0, 5, 0, 0 ] - * [ 0, -1, 0, 0, 4, 0 ] - * [ 0, 0, 0, -2, 0, 6 ] ] - * - */ - // Construct matrix - if( myRank == 0 ){ - A->insertGlobalValues(0,tuple(0,2,4),tuple(7,-3,-1)); - A->insertGlobalValues(1,tuple(0,1),tuple(2,8)); - A->insertGlobalValues(2,tuple(2),tuple(1)); - A->insertGlobalValues(3,tuple(0,3),tuple(-3,5)); - A->insertGlobalValues(4,tuple(1,4),tuple(-1,4)); - A->insertGlobalValues(5,tuple(3,5),tuple(-2,6)); + RCP A; + RCP map; + if (nx > 0) { + #if defined(HAVE_AMESOS2_XPETRA) && defined(HAVE_AMESOS2_GALERI) + typedef Galeri::Xpetra::Problem Galeri_t; + + Teuchos::ParameterList galeriList; + Tpetra::global_size_t nGlobalElements = nx * nx * nx; + galeriList.set("nx", nx); + galeriList.set("ny", nx); + galeriList.set("nz", nx); + if (GaleriName3D == "Elasticity3D") { + GO mx = 1; + galeriList.set("mx", mx); + galeriList.set("my", mx); + galeriList.set("mz", mx); + nGlobalElements *= 3; + } + map = rcp(new MAP(nGlobalElements, 0, comm)); + RCP galeriProblem = + Galeri::Xpetra::BuildProblem + (GaleriName3D, map, galeriList); + A = galeriProblem->BuildMatrix(); + #else + std::cerr << "Galeri or Xpetra not enabled. Exiting..." << std::endl; + return EXIT_SUCCESS; // Otherwise CTest will pick it up as + #endif + } else { + // create a Map + global_size_t nrows = 6; + map = rcp( new MAP(nrows,0,comm) ); + + RCP A = rcp( new MAT(map,3) ); // max of three entries in a row + + /* + * We will solve a system with a known solution, for which we will be using + * the following matrix: + * + * [ [ 7, 0, -3, 0, -1, 0 ] + * [ 2, 8, 0, 0, 0, 0 ] + * [ 0, 0, 1, 0, 0, 0 ] + * [ -3, 0, 0, 5, 0, 0 ] + * [ 0, -1, 0, 0, 4, 0 ] + * [ 0, 0, 0, -2, 0, 6 ] ] + * + */ + // Construct matrix + if( myRank == 0 ){ + A->insertGlobalValues(0,tuple(0,2,4),tuple(7,-3,-1)); + A->insertGlobalValues(1,tuple(0,1),tuple(2,8)); + A->insertGlobalValues(2,tuple(2),tuple(1)); + A->insertGlobalValues(3,tuple(0,3),tuple(-3,5)); + A->insertGlobalValues(4,tuple(1,4),tuple(-1,4)); + A->insertGlobalValues(5,tuple(3,5),tuple(-2,6)); + } + A->fillComplete(); } - A->fillComplete(); - // Create random X + // Create X RCP X = rcp(new MV(map,numVectors)); - X->randomize(); + X->putScalar(1); - /* Create B - * - * Use RHS: - * - * [[-7] - * [18] - * [ 3] - * [17] - * [18] - * [28]] - */ + /* Create B */ RCP B = rcp(new MV(map,numVectors)); - int data[6] = {-7,18,3,17,18,28}; - for( int i = 0; i < 6; ++i ){ - if( B->getMap()->isNodeGlobalElement(i) ){ - B->replaceGlobalValue(i,0,data[i]); - } - } - - // Create solver interface to Superlu with Amesos2 factory method - RCP > solver = Amesos2::create("Superlu", A, X, B); - - solver->symbolicFactorization().numericFactorization().solve(); - - - /* Print the solution - * - * Should be: - * - * [[1] - * [2] - * [3] - * [4] - * [5] - * [6]] - */ - RCP fos = Teuchos::fancyOStream(Teuchos::rcpFromRef(out)); + A->apply(*X, *B); + X->randomize(); - *fos << "Solution :" << std::endl; - X->describe(*fos,Teuchos::VERB_EXTREME); - *fos << std::endl; + // Create solver interface with Amesos2 factory method + RCP > solver = Amesos2::create(solvername, A, X, B); + if (xml_filename != "") { + Teuchos::ParameterList test_params = + Teuchos::ParameterXMLFileReader(xml_filename).getParameters(); + Teuchos::ParameterList& amesos2_params = test_params.sublist("Amesos2"); + solver->setParameters( Teuchos::rcpFromRef(amesos2_params) ); + } + RCP stackedTimer; + stackedTimer = rcp(new Teuchos::StackedTimer("Amesos2 SimpleSolve-File")); + Teuchos::TimeMonitor::setStackedTimer(stackedTimer); + { + solver->symbolicFactorization().numericFactorization().solve(); + } + stackedTimer->stopBaseTimer(); + { + Teuchos::StackedTimer::OutputOptions options; + options.num_histogram=3; + options.print_warnings = false; + options.output_histogram = true; + options.output_fraction=true; + options.output_minmax = true; + stackedTimer->report(std::cout, comm, options); + } + if (verbose) { + /* Print the solution + * + * Should be: + * + * [[1] + * [1] + * [1] + * [1] + * [1] + * [1]] + */ + RCP fos = Teuchos::fancyOStream(Teuchos::rcpFromRef(out)); + *fos << "Solution :" << std::endl; + X->describe(*fos,Teuchos::VERB_EXTREME); + *fos << std::endl; + } // We are done. return 0; }