From 6d4f7b9165c99daf70be71aef4db7c657cf67db1 Mon Sep 17 00:00:00 2001 From: Chris Siefert Date: Wed, 20 Feb 2019 11:08:08 -0700 Subject: [PATCH 1/8] MueLu: Refactoring of MueLu Driver (Part 1) --- packages/muelu/test/scaling/Driver.cpp | 102 +------ packages/muelu/test/scaling/DriverCore.hpp | 337 +++++++++++++++++++++ 2 files changed, 352 insertions(+), 87 deletions(-) create mode 100644 packages/muelu/test/scaling/DriverCore.hpp diff --git a/packages/muelu/test/scaling/Driver.cpp b/packages/muelu/test/scaling/Driver.cpp index 5737d38df382..f597cff7c07d 100644 --- a/packages/muelu/test/scaling/Driver.cpp +++ b/packages/muelu/test/scaling/Driver.cpp @@ -77,6 +77,7 @@ #include #include #include +#include #ifdef HAVE_MUELU_BELOS #include @@ -117,44 +118,8 @@ #include "Xpetra_EpetraMultiVector.hpp" #endif -#include -/*********************************************************************/ -// Support for ML interface -#if defined(HAVE_MUELU_ML) and defined(HAVE_MUELU_EPETRA) -#include -#include - -// Helper functions for compilation purposes -template -struct ML_Wrapper{ - static void Generate_ML_MultiLevelPreconditioner(Teuchos::RCP > & A, Teuchos::ParameterList & mueluList, - Teuchos::RCP > & mlopX) { - throw std::runtime_error("Template parameter mismatch"); - } -}; - - -template -struct ML_Wrapper { - static void Generate_ML_MultiLevelPreconditioner(Teuchos::RCP >& A,Teuchos::ParameterList & mueluList, - Teuchos::RCP >& mlopX) { - typedef double SC; - typedef int LO; - typedef GlobalOrdinal GO; - typedef Kokkos::Compat::KokkosSerialWrapperNode NO; - Teuchos::RCP Aep = Xpetra::Helpers::Op2EpetraCrs(A); - Teuchos::RCP mlop = Teuchos::rcp(new ML_Epetra::MultiLevelPreconditioner(*Aep,mueluList)); -#if defined(HAVE_MUELU_BELOS) - // NOTE: Belos needs the Apply() and AppleInverse() routines of ML swapped. So... - mlop = Teuchos::rcp(new Belos::EpetraPrecOp(mlop)); -#endif - - mlopX = Teuchos::rcp(new Xpetra::EpetraOperator(mlop)); - } -}; -#endif /*********************************************************************/ #ifdef HAVE_MUELU_TPETRA @@ -240,7 +205,7 @@ int main_(Teuchos::CommandLineProcessor &clp, Xpetra::UnderlyingLib& lib, int ar // ========================================================================= typedef Teuchos::ScalarTraits STS; SC zero = STS::zero();//, one = STS::one(); - typedef typename STS::magnitudeType real_type; + typedef typename STS::coordinateType real_type; typedef Xpetra::MultiVector RealValuedMultiVector; // ========================================================================= @@ -280,6 +245,9 @@ int main_(Teuchos::CommandLineProcessor &clp, Xpetra::UnderlyingLib& lib, int ar #ifdef HAVE_MUELU_CUDA bool profileSetup = false; clp.setOption("cuda-profile-setup", "no-cuda-profile-setup", &profileSetup, "enable CUDA profiling for setup"); bool profileSolve = false; clp.setOption("cuda-profile-solve", "no-cuda-profile-solve", &profileSolve, "enable CUDA profiling for solve"); +#else + bool profileSetup = false; + bool profileSolve = false; #endif int cacheSize = 0; clp.setOption("cachesize", &cacheSize, "cache size (in KB)"); #ifdef HAVE_MPI @@ -456,60 +424,20 @@ int main_(Teuchos::CommandLineProcessor &clp, Xpetra::UnderlyingLib& lib, int ar // ========================================================================= // Preconditioner construction // ========================================================================= - comm->barrier(); -#ifdef HAVE_MUELU_CUDA - if(profileSetup) cudaProfilerStart(); -#endif - - tm = rcp(new TimeMonitor(*TimeMonitor::getNewTimer("Driver: 2 - MueLu Setup"))); bool useAMGX = mueluList.isParameter("use external multigrid package") && (mueluList.get("use external multigrid package") == "amgx"); - bool useML = mueluList.isParameter("use external multigrid package") && (mueluList.get("use external multigrid package") == "ml"); - if(useML && lib != Xpetra::UseEpetra) throw std::runtime_error("Error: Cannot use ML on non-epetra matrices"); - - RCP H; - // RCP > AMGXprec; - RCP Prec; - for (int i = 0; i <= numRebuilds; i++) { - A->SetMaxEigenvalueEstimate(-Teuchos::ScalarTraits::one()); - if(useAMGX) { -#if defined(HAVE_MUELU_AMGX) and defined(HAVE_MUELU_TPETRA) - Teuchos::ParameterList dummyList; - Teuchos::RCP > Ac = Utilities::Op2NonConstTpetraCrs(A); - Teuchos::RCP > At = Teuchos::rcp_dynamic_cast >(Ac); - Teuchos::RCP > Top = MueLu::CreateTpetraPreconditioner(At, mueluList, dummyList); - Prec = Teuchos::rcp(new Xpetra::TpetraOperator(Top)); -#endif - } else if(useML) { -#if defined(HAVE_MUELU_ML) and defined(HAVE_MUELU_EPETRA) - mueluList.remove("use external multigrid package"); - if(!coordinates.is_null()) { - RCP epetraCoord = MueLu::Utilities::MV2EpetraMV(coordinates); - if(epetraCoord->NumVectors() > 0) mueluList.set("x-coordinates",(*epetraCoord)[0]); - if(epetraCoord->NumVectors() > 1) mueluList.set("y-coordinates",(*epetraCoord)[1]); - if(epetraCoord->NumVectors() > 2) mueluList.set("z-coordinates",(*epetraCoord)[2]); - } - ML_Wrapper::Generate_ML_MultiLevelPreconditioner(A,mueluList,Prec); -#endif - } - else { - const std::string userName = "user data"; - Teuchos::Array lNodesPerDim(3, 10); - Teuchos::ParameterList& userParamList = mueluList.sublist(userName); - userParamList.set >("Coordinates", coordinates); - userParamList.set> >("Nullspace", nullspace); - userParamList.set >("Array lNodesPerDim", lNodesPerDim); + bool useML = mueluList.isParameter("use external multigrid package") && (mueluList.get("use external multigrid package") == "ml"); #ifdef HAVE_MPI - if(provideNodeComm) { - userParamList.set("Node Comm",nodeComm); - } -#endif - - H = MueLu::CreateXpetraPreconditioner(A, mueluList, mueluList); - } + if(provideNodeComm && !useAMGX && !useML) { + Teuchos::ParameterList& userParamList = mueluList.sublist("user data"); + userParamList.set("Node Comm",nodeComm); } -#ifdef HAVE_MUELU_CUDA - if(profileSetup) cudaProfilerStop(); #endif + RCP H; + RCP Prec; + comm->barrier(); + MUELU_SWITCH_TIME_MONITOR(tm,"Driver: 2 - MueLu Setup"); + PreconditionerSetup(A,coordinates,nullspace,mueluList,profileSetup,useAMGX,useML,numRebuilds,H,Prec); + // Get the raw matrices for matvec testing #if defined(HAVE_MUELU_TPETRA) diff --git a/packages/muelu/test/scaling/DriverCore.hpp b/packages/muelu/test/scaling/DriverCore.hpp new file mode 100644 index 000000000000..d62e1c5695b8 --- /dev/null +++ b/packages/muelu/test/scaling/DriverCore.hpp @@ -0,0 +1,337 @@ +// @HEADER +// +// *********************************************************************** +// +// MueLu: A package for multigrid based preconditioning +// Copyright 2012 Sandia Corporation +// +// Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation, +// the U.S. Government retains certain rights in this software. +// +// Redistribution and use in source and binary forms, with or without +// modification, are permitted provided that the following conditions are +// met: +// +// 1. Redistributions of source code must retain the above copyright +// notice, this list of conditions and the following disclaimer. +// +// 2. Redistributions in binary form must reproduce the above copyright +// notice, this list of conditions and the following disclaimer in the +// documentation and/or other materials provided with the distribution. +// +// 3. Neither the name of the Corporation nor the names of the +// contributors may be used to endorse or promote products derived from +// this software without specific prior written permission. +// +// THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY +// EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE +// IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR +// PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE +// CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, +// EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, +// PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR +// PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF +// LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING +// NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS +// SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. +// +// Questions? Contact +// Jonathan Hu (jhu@sandia.gov) +// Andrey Prokopenko (aprokop@sandia.gov) +// Ray Tuminaro (rstumin@sandia.gov) +// +// *********************************************************************** +// +// @HEADER +#ifndef DRIVERCORE_HPP +#include + +// Teuchos +#include + +// Xpetra +#include +#include +#include + + +// Belos +#ifdef HAVE_MUELU_BELOS +#include +#include +#include +#include +#include +#include +#include // => This header defines Belos::XpetraOp +#include // => This header defines Belos::MueLuOp +#ifdef HAVE_MUELU_TPETRA +#include // => This header defines Belos::TpetraOp +#endif +#ifdef HAVE_MUELU_EPETRA +#include // => This header defines Belos::EpetraPrecOp +#endif +#endif + +// Cuda +#ifdef HAVE_MUELU_CUDA +#include "cuda_profiler_api.h" +#endif + +// AMGX +#ifdef HAVE_MUELU_AMGX +#include +#endif + + +#include + + + +//************************************************************************************* +//************************************************************************************* +//************************************************************************************* +// A handy macro to switch time monitors in a StackedTimer-compatible way +#define MUELU_SWITCH_TIME_MONITOR(tm,timername) \ + {tm=Teuchos::null; tm = rcp(new TimeMonitor(*TimeMonitor::getNewTimer("timername")));} + + +//************************************************************************************* +//************************************************************************************* +//************************************************************************************* +// Support for ML interface +#if defined(HAVE_MUELU_ML) and defined(HAVE_MUELU_EPETRA) +#include +#include + +// Helper functions for compilation purposes +template +struct ML_Wrapper{ + static void Generate_ML_MultiLevelPreconditioner(Teuchos::RCP > & A, Teuchos::ParameterList & mueluList, + Teuchos::RCP > & mlopX) { + throw std::runtime_error("Template parameter mismatch"); + } +}; + + +template +struct ML_Wrapper { + static void Generate_ML_MultiLevelPreconditioner(Teuchos::RCP >& A,Teuchos::ParameterList & mueluList, + Teuchos::RCP >& mlopX) { + typedef double SC; + typedef int LO; + typedef GlobalOrdinal GO; + typedef Kokkos::Compat::KokkosSerialWrapperNode NO; + Teuchos::RCP Aep = Xpetra::Helpers::Op2EpetraCrs(A); + Teuchos::RCP mlop = Teuchos::rcp(new ML_Epetra::MultiLevelPreconditioner(*Aep,mueluList)); +#if defined(HAVE_MUELU_BELOS) + // NOTE: Belos needs the Apply() and AppleInverse() routines of ML swapped. So... + mlop = Teuchos::rcp(new Belos::EpetraPrecOp(mlop)); +#endif + + mlopX = Teuchos::rcp(new Xpetra::EpetraOperator(mlop)); + } +}; +#endif + + +//************************************************************************************* +//************************************************************************************* +//************************************************************************************* +// This is a standard setup routine +template +void PreconditionerSetup(Teuchos::RCP > & A, + Teuchos::RCP::coordinateType,LocalOrdinal,GlobalOrdinal,Node> > & coordinates, + Teuchos::RCP > & nullspace, + Teuchos::ParameterList & mueluList, + bool profileSetup, + bool useAMGX, + bool useML, + int numRebuilds, + Teuchos::RCP > & H, + Teuchos::RCP > & Prec) { +#include + using Teuchos::RCP; + Xpetra::UnderlyingLib lib = A->getRowMap()->lib(); + typedef typename Teuchos::ScalarTraits::coordinateType coordinate_type; + typedef Xpetra::MultiVector CoordinateMultiVector; + // ========================================================================= + // Preconditioner construction + // ========================================================================= +#ifdef HAVE_MUELU_CUDA + if(profileSetup) cudaProfilerStart(); +#endif + + if(useML && lib != Xpetra::UseEpetra) throw std::runtime_error("Error: Cannot use ML on non-epetra matrices"); + + for (int i = 0; i <= numRebuilds; i++) { + A->SetMaxEigenvalueEstimate(-Teuchos::ScalarTraits::one()); + if(useAMGX) { +#if defined(HAVE_MUELU_AMGX) and defined(HAVE_MUELU_TPETRA) + Teuchos::ParameterList dummyList; + RCP > Ac = Utilities::Op2NonConstTpetraCrs(A); + RCP > At = Teuchos::rcp_dynamic_cast >(Ac); + RCP > Top = MueLu::CreateTpetraPreconditioner(At, mueluList, dummyList); + Prec = Teuchos::rcp(new Xpetra::TpetraOperator(Top)); +#endif + } else if(useML) { +#if defined(HAVE_MUELU_ML) and defined(HAVE_MUELU_EPETRA) + mueluList.remove("use external multigrid package"); + if(!coordinates.is_null()) { + RCP epetraCoord = MueLu::Utilities::MV2EpetraMV(coordinates); + if(epetraCoord->NumVectors() > 0) mueluList.set("x-coordinates",(*epetraCoord)[0]); + if(epetraCoord->NumVectors() > 1) mueluList.set("y-coordinates",(*epetraCoord)[1]); + if(epetraCoord->NumVectors() > 2) mueluList.set("z-coordinates",(*epetraCoord)[2]); + } + ML_Wrapper::Generate_ML_MultiLevelPreconditioner(A,mueluList,Prec); +#endif + } + else { + Teuchos::Array lNodesPerDim(3, 10); + Teuchos::ParameterList& userParamList = mueluList.sublist("user data"); + userParamList.set >("Coordinates", coordinates); + userParamList.set> >("Nullspace", nullspace); + userParamList.set >("Array lNodesPerDim", lNodesPerDim); + H = MueLu::CreateXpetraPreconditioner(A, mueluList, mueluList); + } + } +#ifdef HAVE_MUELU_CUDA + if(profileSetup) cudaProfilerStop(); +#endif +} + + + +//************************************************************************************* +#if 0 +template +void SystemSolve(Teuchos::RCP > & A, + Teuchos::RCP > & X, + Teuchos::RCP > & B, + Teuchos::RCP & H, + Teuchos::RCP & Prec, + const Teuchos::ParameterList & belosList, + Teuchos::FancyOStream & out, + std::string solveType, + bool profileSolve, + bool useAMGX, + bool useML, + int numResolves, + int maxIts, + double tol) { +#include + using Teuchos::RCP; + Xpetra::UnderlyingLib lib = A->getRowMap()->lib(); + + for(int solveno = 0; solveno<=numResolves; solveno++) { + RCP tm = rcp(new TimeMonitor(*TimeMonitor::getNewTimer("Driver: 3 - LHS and RHS initialization"))); + X->putScalar(zero); + tm = Teuchos::null; + + if (solveType == "none") { + // Do not perform a solve + } else if (solveType == "matvec") { + // Just do matvecs + tm = rcp (new TimeMonitor(*TimeMonitor::getNewTimer("Driver: 6 - Matvec"))); +#if defined(HAVE_MUELU_TPETRA) + if(lib==Xpetra::UseTpetra) Atpetra->apply(*Btpetra,*Xtpetra); +#endif +#if defined(HAVE_MUELU_EPETRA) && !defined(HAVE_MUELU_INST_COMPLEX_INT_INT) + if(lib==Xpetra::UseEpetra) Aepetra->Apply(*Bepetra,*Xepetra); +#endif + //clear the cache (and don't time it) + tm = Teuchos::null; + int ttt = rand(); + for (int i=0; i(max - min + 1))); + } else if (solveType == "standalone") { + tm = rcp (new TimeMonitor(*TimeMonitor::getNewTimer("Driver: 4 - Fixed Point Solve"))); +#ifdef HAVE_MUELU_CUDA + if(profileSolve) cudaProfilerStart(); +#endif + H->IsPreconditioner(false); + H->Iterate(*B, *X, maxIts); +#ifdef HAVE_MUELU_CUDA + if(profileSolve) cudaProfilerStop(); +#endif + } else if (solveType == "belos") { +#ifdef HAVE_MUELU_BELOS + tm = rcp(new TimeMonitor(*TimeMonitor::getNewTimer("Driver: 5 - Belos Solve"))); +#ifdef HAVE_MUELU_CUDA + if(profileSolve) cudaProfilerStart(); +#endif + // Operator and Multivector type that will be used with Belos + typedef MultiVector MV; + typedef Belos::OperatorT OP; + + // Define Operator and Preconditioner + Teuchos::RCP belosOp = Teuchos::rcp(new Belos::XpetraOp(A)); // Turns a Xpetra::Matrix object into a Belos operato + Teuchos::RCP belosPrec; // Turns a MueLu::Hierarchy object into a Belos operator + if(useAMGX) { +#if defined(HAVE_MUELU_AMGX) and defined(HAVE_MUELU_TPETRA) + belosPrec = Teuchos::rcp(new Belos::XpetraOp (Prec)); // Turns an Xpetra::Operator object into a Belos operator +#endif + } + else if(useML) { +#if defined(HAVE_MUELU_ML) and defined(HAVE_MUELU_EPETRA) + belosPrec = Teuchos::rcp(new Belos::XpetraOp (Prec)); // Turns an Xpetra::Operator object into a Belos operator +#endif + } + else { + H->IsPreconditioner(true); + belosPrec = Teuchos::rcp(new Belos::MueLuOp (H)); // Turns a MueLu::Hierarchy object into a Belos operator + } + + // Construct a Belos LinearProblem object + RCP > belosProblem = rcp(new Belos::LinearProblem(belosOp, X, B)); + if(solvePreconditioned) belosProblem->setRightPrec(belosPrec); + + bool set = belosProblem->setProblem(); + if (set == false) { + throw MueLu::Exceptions::RuntimeError("ERROR: Belos::LinearProblem failed to set up correctly!"); + } + + // Belos parameter list + RCP belosList = Teuchos::parameterList(); + belosList->set("Maximum Iterations", maxIts); // Maximum number of iterations allowed + belosList->set("Convergence Tolerance", tol); // Relative convergence tolerance requested + belosList->set("Verbosity", Belos::Errors + Belos::Warnings + Belos::StatusTestDetails); + belosList->set("Output Frequency", 1); + belosList->set("Output Style", Belos::Brief); + if (!scaleResidualHist) + belosList->set("Implicit Residual Scaling", "None"); + + // Create an iterative solver manager + Belos::SolverFactory solverFactory; + RCP< Belos::SolverManager > solver = solverFactory.create(belosType, belosList); + solver->setProblem(belosProblem); + + // Perform solve + Belos::ReturnType ret = Belos::Unconverged; + ret = solver->solve(); + + // Get the number of iterations for this solve. + out << "Number of iterations performed for this solve: " << solver->getNumIters() << std::endl; + // Check convergence + if (ret != Belos::Converged) + out << std::endl << "ERROR: Belos did not converge! " << std::endl; + else + out << std::endl << "SUCCESS: Belos converged!" << std::endl; +#ifdef HAVE_MUELU_CUDA + if(profileSolve) cudaProfilerStop(); +#endif +#endif //ifdef HAVE_MUELU_BELOS + } else { + throw MueLu::Exceptions::RuntimeError("Unknown solver type: \"" + solveType + "\""); + } + }// end resolves + +} +#endif + + + + + + +#endif From 1b1f31526263aea492258b74ecb9167fb73e34a8 Mon Sep 17 00:00:00 2001 From: Chris Siefert Date: Wed, 20 Feb 2019 11:20:04 -0700 Subject: [PATCH 2/8] MueLu: More mods to the driver to extract code for common reuse --- packages/muelu/test/scaling/Driver.cpp | 143 +-------------------- packages/muelu/test/scaling/DriverCore.hpp | 60 ++++++--- 2 files changed, 49 insertions(+), 154 deletions(-) diff --git a/packages/muelu/test/scaling/Driver.cpp b/packages/muelu/test/scaling/Driver.cpp index f597cff7c07d..793095fcbda0 100644 --- a/packages/muelu/test/scaling/Driver.cpp +++ b/packages/muelu/test/scaling/Driver.cpp @@ -436,28 +436,9 @@ int main_(Teuchos::CommandLineProcessor &clp, Xpetra::UnderlyingLib& lib, int ar RCP Prec; comm->barrier(); MUELU_SWITCH_TIME_MONITOR(tm,"Driver: 2 - MueLu Setup"); - PreconditionerSetup(A,coordinates,nullspace,mueluList,profileSetup,useAMGX,useML,numRebuilds,H,Prec); - - // Get the raw matrices for matvec testing -#if defined(HAVE_MUELU_TPETRA) - Teuchos::RCP > Atpetra; - Teuchos::RCP > Xtpetra,Btpetra; - if(lib==Xpetra::UseTpetra) { - Atpetra = Utilities::Op2NonConstTpetraCrs(A); - Xtpetra = Teuchos::rcp(& Xpetra::toTpetra(*X),false); - Btpetra = Teuchos::rcp(& Xpetra::toTpetra(*B),false); - } -#endif -#if defined(HAVE_MUELU_EPETRA) && !defined(HAVE_MUELU_INST_COMPLEX_INT_INT) - Teuchos::RCP Aepetra; - Teuchos::RCP Xepetra,Bepetra; - if(lib==Xpetra::UseEpetra) { - Aepetra = Xpetra::Helpers::Op2EpetraCrs(A); - Xepetra = Teuchos::rcp(& Xpetra::toEpetra(*X),false); - Bepetra = Teuchos::rcp(& Xpetra::toEpetra(*B),false); - } -#endif + // Build the preconditioner numRebuilds+1 times + PreconditionerSetup(A,coordinates,nullspace,mueluList,profileSetup,useAMGX,useML,numRebuilds,H,Prec); comm->barrier(); tm = Teuchos::null; @@ -466,130 +447,14 @@ int main_(Teuchos::CommandLineProcessor &clp, Xpetra::UnderlyingLib& lib, int ar // System solution (Ax = b) // ========================================================================= comm->barrier(); - if (writeMatricesOPT > -2) { tm = rcp(new TimeMonitor(*TimeMonitor::getNewTimer("Driver: 3.5 - Matrix output"))); H->Write(writeMatricesOPT, writeMatricesOPT); tm = Teuchos::null; } - std::vector tempVector; - int min = 0, max = 10; - int numInts = 0; - if (cacheSize > 0) { - cacheSize *= 1024; //convert to bytes - numInts = cacheSize/sizeof(int) + 1; - tempVector.resize(numInts); - } - - for(int solveno = 0; solveno<=numResolves; solveno++) { - tm = rcp(new TimeMonitor(*TimeMonitor::getNewTimer("Driver: 3 - LHS and RHS initialization"))); - X->putScalar(zero); - tm = Teuchos::null; - - if (solveType == "none") { - // Do not perform a solve - } else if (solveType == "matvec") { - // Just do matvecs - tm = rcp (new TimeMonitor(*TimeMonitor::getNewTimer("Driver: 6 - Matvec"))); -#if defined(HAVE_MUELU_TPETRA) - if(lib==Xpetra::UseTpetra) Atpetra->apply(*Btpetra,*Xtpetra); -#endif -#if defined(HAVE_MUELU_EPETRA) && !defined(HAVE_MUELU_INST_COMPLEX_INT_INT) - if(lib==Xpetra::UseEpetra) Aepetra->Apply(*Bepetra,*Xepetra); -#endif - //clear the cache (and don't time it) - tm = Teuchos::null; - int ttt = rand(); - for (int i=0; i(max - min + 1))); - } else if (solveType == "standalone") { - tm = rcp (new TimeMonitor(*TimeMonitor::getNewTimer("Driver: 4 - Fixed Point Solve"))); -#ifdef HAVE_MUELU_CUDA - if(profileSolve) cudaProfilerStart(); -#endif - H->IsPreconditioner(false); - H->Iterate(*B, *X, maxIts); -#ifdef HAVE_MUELU_CUDA - if(profileSolve) cudaProfilerStop(); -#endif - } else if (solveType == "belos") { -#ifdef HAVE_MUELU_BELOS - tm = rcp(new TimeMonitor(*TimeMonitor::getNewTimer("Driver: 5 - Belos Solve"))); -#ifdef HAVE_MUELU_CUDA - if(profileSolve) cudaProfilerStart(); -#endif - // Operator and Multivector type that will be used with Belos - typedef MultiVector MV; - typedef Belos::OperatorT OP; - - // Define Operator and Preconditioner - Teuchos::RCP belosOp = Teuchos::rcp(new Belos::XpetraOp(A)); // Turns a Xpetra::Matrix object into a Belos operato - Teuchos::RCP belosPrec; // Turns a MueLu::Hierarchy object into a Belos operator - if(useAMGX) { -#if defined(HAVE_MUELU_AMGX) and defined(HAVE_MUELU_TPETRA) - belosPrec = Teuchos::rcp(new Belos::XpetraOp (Prec)); // Turns an Xpetra::Operator object into a Belos operator -#endif - } - else if(useML) { -#if defined(HAVE_MUELU_ML) and defined(HAVE_MUELU_EPETRA) - belosPrec = Teuchos::rcp(new Belos::XpetraOp (Prec)); // Turns an Xpetra::Operator object into a Belos operator -#endif - } - else { - H->IsPreconditioner(true); - belosPrec = Teuchos::rcp(new Belos::MueLuOp (H)); // Turns a MueLu::Hierarchy object into a Belos operator - } - - // Construct a Belos LinearProblem object - RCP > belosProblem = rcp(new Belos::LinearProblem(belosOp, X, B)); - if(solvePreconditioned) belosProblem->setRightPrec(belosPrec); - - bool set = belosProblem->setProblem(); - if (set == false) { - out << "\nERROR: Belos::LinearProblem failed to set up correctly!" << std::endl; - // this fixes the resource leak detected by coverity (CID134984) - if (openedOut != NULL) { - fclose(openedOut); - openedOut = NULL; - } - return EXIT_FAILURE; - } - - // Belos parameter list - RCP belosList = Teuchos::parameterList(); - belosList->set("Maximum Iterations", maxIts); // Maximum number of iterations allowed - belosList->set("Convergence Tolerance", tol); // Relative convergence tolerance requested - belosList->set("Verbosity", Belos::Errors + Belos::Warnings + Belos::StatusTestDetails); - belosList->set("Output Frequency", 1); - belosList->set("Output Style", Belos::Brief); - if (!scaleResidualHist) - belosList->set("Implicit Residual Scaling", "None"); - - // Create an iterative solver manager - Belos::SolverFactory solverFactory; - RCP< Belos::SolverManager > solver = solverFactory.create(belosType, belosList); - solver->setProblem(belosProblem); - - // Perform solve - Belos::ReturnType ret = Belos::Unconverged; - ret = solver->solve(); - - // Get the number of iterations for this solve. - out << "Number of iterations performed for this solve: " << solver->getNumIters() << std::endl; - // Check convergence - if (ret != Belos::Converged) - out << std::endl << "ERROR: Belos did not converge! " << std::endl; - else - out << std::endl << "SUCCESS: Belos converged!" << std::endl; -#ifdef HAVE_MUELU_CUDA - if(profileSolve) cudaProfilerStop(); -#endif -#endif //ifdef HAVE_MUELU_BELOS - } else { - throw MueLu::Exceptions::RuntimeError("Unknown solver type: \"" + solveType + "\""); - } - }// end resolves + // Solve the system numResolves+1 times + SystemSolve(A,X,B,H,Prec,out,solveType,belosType,profileSolve,useAMGX,useML,cacheSize,numResolves,scaleResidualHist,solvePreconditioned,maxIts,tol); comm->barrier(); tm = Teuchos::null; diff --git a/packages/muelu/test/scaling/DriverCore.hpp b/packages/muelu/test/scaling/DriverCore.hpp index d62e1c5695b8..4935065bf4f6 100644 --- a/packages/muelu/test/scaling/DriverCore.hpp +++ b/packages/muelu/test/scaling/DriverCore.hpp @@ -203,26 +203,62 @@ void PreconditionerSetup(Teuchos::RCP void SystemSolve(Teuchos::RCP > & A, Teuchos::RCP > & X, Teuchos::RCP > & B, - Teuchos::RCP & H, - Teuchos::RCP & Prec, - const Teuchos::ParameterList & belosList, + Teuchos::RCP > & H, + Teuchos::RCP > & Prec, Teuchos::FancyOStream & out, std::string solveType, + std::string belosType, bool profileSolve, bool useAMGX, bool useML, + int cacheSize, int numResolves, + bool scaleResidualHist, + bool solvePreconditioned, int maxIts, double tol) { #include using Teuchos::RCP; + using Teuchos::rcp; + using Teuchos::TimeMonitor; Xpetra::UnderlyingLib lib = A->getRowMap()->lib(); - + typedef Teuchos::ScalarTraits STS; + SC zero = STS::zero(); + + // Cache clearing + std::vector tempVector; + int min = 0, max = 10; + int numInts = 0; + if (cacheSize > 0) { + cacheSize *= 1024; //convert to bytes + numInts = cacheSize/sizeof(int) + 1; + tempVector.resize(numInts); + } + + // Get the raw matrices for matvec testing +#if defined(HAVE_MUELU_TPETRA) + Teuchos::RCP > Atpetra; + Teuchos::RCP > Xtpetra,Btpetra; + if(lib==Xpetra::UseTpetra) { + Atpetra = Utilities::Op2NonConstTpetraCrs(A); + Xtpetra = rcp(& Xpetra::toTpetra(*X),false); + Btpetra = rcp(& Xpetra::toTpetra(*B),false); + } +#endif +#if defined(HAVE_MUELU_EPETRA) && !defined(HAVE_MUELU_INST_COMPLEX_INT_INT) + Teuchos::RCP Aepetra; + Teuchos::RCP Xepetra,Bepetra; + if(lib==Xpetra::UseEpetra) { + Aepetra = Xpetra::Helpers::Op2EpetraCrs(A); + Xepetra = rcp(& Xpetra::toEpetra(*X),false); + Bepetra = rcp(& Xpetra::toEpetra(*B),false); + } +#endif + for(int solveno = 0; solveno<=numResolves; solveno++) { RCP tm = rcp(new TimeMonitor(*TimeMonitor::getNewTimer("Driver: 3 - LHS and RHS initialization"))); X->putScalar(zero); @@ -265,21 +301,21 @@ void SystemSolve(Teuchos::RCP OP; // Define Operator and Preconditioner - Teuchos::RCP belosOp = Teuchos::rcp(new Belos::XpetraOp(A)); // Turns a Xpetra::Matrix object into a Belos operato + Teuchos::RCP belosOp = rcp(new Belos::XpetraOp(A)); // Turns a Xpetra::Matrix object into a Belos operato Teuchos::RCP belosPrec; // Turns a MueLu::Hierarchy object into a Belos operator if(useAMGX) { #if defined(HAVE_MUELU_AMGX) and defined(HAVE_MUELU_TPETRA) - belosPrec = Teuchos::rcp(new Belos::XpetraOp (Prec)); // Turns an Xpetra::Operator object into a Belos operator + belosPrec = rcp(new Belos::XpetraOp (Prec)); // Turns an Xpetra::Operator object into a Belos operator #endif } else if(useML) { #if defined(HAVE_MUELU_ML) and defined(HAVE_MUELU_EPETRA) - belosPrec = Teuchos::rcp(new Belos::XpetraOp (Prec)); // Turns an Xpetra::Operator object into a Belos operator + belosPrec = rcp(new Belos::XpetraOp (Prec)); // Turns an Xpetra::Operator object into a Belos operator #endif } else { H->IsPreconditioner(true); - belosPrec = Teuchos::rcp(new Belos::MueLuOp (H)); // Turns a MueLu::Hierarchy object into a Belos operator + belosPrec = rcp(new Belos::MueLuOp (H)); // Turns a MueLu::Hierarchy object into a Belos operator } // Construct a Belos LinearProblem object @@ -327,11 +363,5 @@ void SystemSolve(Teuchos::RCP Date: Wed, 20 Feb 2019 11:52:31 -0700 Subject: [PATCH 3/8] MueLu: Adding an example for multiple runs in a single deck --- packages/muelu/test/scaling/CMakeLists.txt | 18 +++- .../muelu/test/scaling/scaling-with-rerun.xml | 84 +++++++++++++++++++ 2 files changed, 101 insertions(+), 1 deletion(-) create mode 100644 packages/muelu/test/scaling/scaling-with-rerun.xml diff --git a/packages/muelu/test/scaling/CMakeLists.txt b/packages/muelu/test/scaling/CMakeLists.txt index 675b6b9ff2aa..ffa47f411c6c 100644 --- a/packages/muelu/test/scaling/CMakeLists.txt +++ b/packages/muelu/test/scaling/CMakeLists.txt @@ -39,7 +39,7 @@ IF (${PACKAGE_NAME}_HAVE_TPETRA_SOLVER_STACK OR ${PACKAGE_NAME}_HAVE_EPETRA_SOLV ) TRIBITS_COPY_FILES_TO_BINARY_DIR(Driver_cp - SOURCE_FILES scaling.xml scaling.yaml scaling-complex.xml scaling-withglobalconstants.xml scaling-complex-withglobalconstants.xml circ_nsp_dependency.xml isorropia.xml iso_poisson.xml conchas_milestone_zoltan.xml conchas_milestone_zoltan2.xml conchas_milestone_zoltan2_complex.xml sa_with_ilu.xml sa_with_Ifpack2_line_detection.xml rap.xml smoother.xml smoother_complex.xml tripleMatrixProduct.xml scaling-ml.xml elasticity3D.xml amgx.json amgx.xml + SOURCE_FILES scaling.xml scaling.yaml scaling-complex.xml scaling-withglobalconstants.xml scaling-complex-withglobalconstants.xml circ_nsp_dependency.xml isorropia.xml iso_poisson.xml conchas_milestone_zoltan.xml conchas_milestone_zoltan2.xml conchas_milestone_zoltan2_complex.xml sa_with_ilu.xml sa_with_Ifpack2_line_detection.xml rap.xml smoother.xml smoother_complex.xml tripleMatrixProduct.xml scaling-ml.xml elasticity3D.xml amgx.json amgx.xml scaling-with-rerun.xml ) TRIBITS_ADD_EXECUTABLE( @@ -125,6 +125,14 @@ IF (${PACKAGE_NAME}_HAVE_EPETRA_SOLVER_STACK AND (NOT Xpetra_INT_LONG_LONG)) COMM mpi # HAVE_MPI required ) + TRIBITS_ADD_TEST( + Driver + NAME "DriverEpetra_Rerun" + ARGS "--linAlgebra=Epetra --xml=scaling-with-rerun.xml" + NUM_MPI_PROCS 4 + COMM mpi # HAVE_MPI required + ) + TRIBITS_ADD_TEST( Driver NAME "DriverEpetra_isotropic_poisson" @@ -225,6 +233,14 @@ IF (${PACKAGE_NAME}_HAVE_TPETRA_SOLVER_STACK) COMM mpi # HAVE_MPI required ) + TRIBITS_ADD_TEST( + Driver + NAME "DriverTpetra_Rerun" + ARGS "--linAlgebra=Tpetra --xml=scaling-with-rerun.xml" + NUM_MPI_PROCS 4 + COMM mpi # HAVE_MPI required + ) + TRIBITS_ADD_TEST( Driver NAME "DriverTpetra_WithGlobalConstants" diff --git a/packages/muelu/test/scaling/scaling-with-rerun.xml b/packages/muelu/test/scaling/scaling-with-rerun.xml new file mode 100644 index 000000000000..b6212ac34c72 --- /dev/null +++ b/packages/muelu/test/scaling/scaling-with-rerun.xml @@ -0,0 +1,84 @@ + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + From 9dd85a65db2b908a42ee1fc9b4ed51e3d21b55a2 Mon Sep 17 00:00:00 2001 From: Chris Siefert Date: Wed, 20 Feb 2019 12:17:49 -0700 Subject: [PATCH 4/8] MueLu: Adding error handling to MueLu Driver --- packages/muelu/test/scaling/Driver.cpp | 98 ++++++++++++++++---------- 1 file changed, 62 insertions(+), 36 deletions(-) diff --git a/packages/muelu/test/scaling/Driver.cpp b/packages/muelu/test/scaling/Driver.cpp index 793095fcbda0..6ad17bfc2103 100644 --- a/packages/muelu/test/scaling/Driver.cpp +++ b/packages/muelu/test/scaling/Driver.cpp @@ -352,12 +352,11 @@ int main_(Teuchos::CommandLineProcessor &clp, Xpetra::UnderlyingLib& lib, int ar int numReruns = 1; if (paramList.isParameter("number of reruns")) numReruns = paramList.get("number of reruns"); - - const bool mustAlreadyExist = true; + for (int rerunCount = 1; rerunCount <= numReruns; rerunCount++) { - ParameterList mueluList, runList; - bool stop = false; + ParameterList mueluList, runList; + const bool mustAlreadyExist = true; if (isDriver) { runList = paramList.sublist("Run1", mustAlreadyExist); mueluList = runList .sublist("MueLu", mustAlreadyExist); @@ -394,13 +393,14 @@ int main_(Teuchos::CommandLineProcessor &clp, Xpetra::UnderlyingLib& lib, int ar } } + int runCount = 1; + int savedOut = -1; + FILE* openedOut = NULL; do { solveType = dsolveType; tol = dtol; - int savedOut = -1; - FILE* openedOut = NULL; if (isDriver) { if (runList.isParameter("filename")) { // Redirect all output into a filename We have to redirect all output, @@ -419,7 +419,13 @@ int main_(Teuchos::CommandLineProcessor &clp, Xpetra::UnderlyingLib& lib, int ar if (runList.isParameter("tol")) tol = runList.get ("tol"); } - out << galeriStream.str(); + RCP fancy2 = Teuchos::fancyOStream(Teuchos::rcpFromRef(std::cout)); + Teuchos::FancyOStream& out2 = *fancy2; + out2.setOutputToRootOnly(0); + + + + out2 << galeriStream.str(); // ========================================================================= // Preconditioner construction @@ -432,31 +438,46 @@ int main_(Teuchos::CommandLineProcessor &clp, Xpetra::UnderlyingLib& lib, int ar userParamList.set("Node Comm",nodeComm); } #endif + out2<<"*********** ParameterList ***********"< H; RCP Prec; - comm->barrier(); - MUELU_SWITCH_TIME_MONITOR(tm,"Driver: 2 - MueLu Setup"); - - // Build the preconditioner numRebuilds+1 times - PreconditionerSetup(A,coordinates,nullspace,mueluList,profileSetup,useAMGX,useML,numRebuilds,H,Prec); - - comm->barrier(); - tm = Teuchos::null; - + try { + comm->barrier(); + // Build the preconditioner numRebuilds+1 times + MUELU_SWITCH_TIME_MONITOR(tm,"Driver: 2 - MueLu Setup"); + PreconditionerSetup(A,coordinates,nullspace,mueluList,profileSetup,useAMGX,useML,numRebuilds,H,Prec); + + comm->barrier(); + tm = Teuchos::null; + } + catch(...) { + out2<<"MueLu_Driver: preconditioner setup crashed!"<barrier(); - if (writeMatricesOPT > -2) { - tm = rcp(new TimeMonitor(*TimeMonitor::getNewTimer("Driver: 3.5 - Matrix output"))); - H->Write(writeMatricesOPT, writeMatricesOPT); - tm = Teuchos::null; + try { + comm->barrier(); + if (writeMatricesOPT > -2) { + tm = rcp(new TimeMonitor(*TimeMonitor::getNewTimer("Driver: 3.5 - Matrix output"))); + H->Write(writeMatricesOPT, writeMatricesOPT); + tm = Teuchos::null; + } + + // Solve the system numResolves+1 times + SystemSolve(A,X,B,H,Prec,out2,solveType,belosType,profileSolve,useAMGX,useML,cacheSize,numResolves,scaleResidualHist,solvePreconditioned,maxIts,tol); + + comm->barrier(); + } + catch(...) { + out2<<"MueLu_Driver: solver crashed!"<barrier(); tm = Teuchos::null; globalTimeMonitor = Teuchos::null; @@ -474,22 +495,17 @@ int main_(Teuchos::CommandLineProcessor &clp, Xpetra::UnderlyingLib& lib, int ar const std::string filter = ""; - std::ios_base::fmtflags ff(out.flags()); - if (timingsFormat == "table-fixed") out << std::fixed; - else out << std::scientific; + std::ios_base::fmtflags ff(out2.flags()); + if (timingsFormat == "table-fixed") out2 << std::fixed; + else out2 << std::scientific; TimeMonitor::report(comm.ptr(), out, filter, reportParams); - out << std::setiosflags(ff); + out2 << std::setiosflags(ff); } TimeMonitor::clearCounters(); + out2 << std::endl; if (isDriver) { - if (openedOut != NULL) { - TEUCHOS_ASSERT(savedOut >= 0); - dup2(savedOut, STDOUT_FILENO); - fclose(openedOut); - openedOut = NULL; - } try { runList = paramList.sublist("Run" + MueLu::toString(++runCount), mustAlreadyExist); mueluList = runList .sublist("MueLu", mustAlreadyExist); @@ -497,8 +513,18 @@ int main_(Teuchos::CommandLineProcessor &clp, Xpetra::UnderlyingLib& lib, int ar stop = true; } } - + comm->barrier(); } while (!stop); + + // Cleanup Output + if (openedOut != NULL) { + TEUCHOS_ASSERT(savedOut >= 0); + dup2(savedOut, STDOUT_FILENO); + fclose(openedOut); + openedOut = NULL; + } + + }//end reruns return EXIT_SUCCESS; From b0c31550c0e29c75af79bd5a9594ab88dbd63322 Mon Sep 17 00:00:00 2001 From: Chris Siefert Date: Wed, 20 Feb 2019 12:17:58 -0700 Subject: [PATCH 5/8] MueLu: File --- packages/muelu/test/scaling/scaling-with-rerun.xml | 1 - 1 file changed, 1 deletion(-) diff --git a/packages/muelu/test/scaling/scaling-with-rerun.xml b/packages/muelu/test/scaling/scaling-with-rerun.xml index b6212ac34c72..d9ec7f288283 100644 --- a/packages/muelu/test/scaling/scaling-with-rerun.xml +++ b/packages/muelu/test/scaling/scaling-with-rerun.xml @@ -1,5 +1,4 @@ - From d8f2528eb2487d1b1b994ba400526e0de0755d14 Mon Sep 17 00:00:00 2001 From: Chris Siefert Date: Wed, 20 Feb 2019 14:24:35 -0700 Subject: [PATCH 6/8] MueLu: Prettify --- packages/muelu/test/scaling/Driver.cpp | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/packages/muelu/test/scaling/Driver.cpp b/packages/muelu/test/scaling/Driver.cpp index 6ad17bfc2103..ba99ad986ff3 100644 --- a/packages/muelu/test/scaling/Driver.cpp +++ b/packages/muelu/test/scaling/Driver.cpp @@ -438,9 +438,9 @@ int main_(Teuchos::CommandLineProcessor &clp, Xpetra::UnderlyingLib& lib, int ar userParamList.set("Node Comm",nodeComm); } #endif - out2<<"*********** ParameterList ***********"< H; RCP Prec; From b8265eadbadc824435a770cd507c79b135e5a9b3 Mon Sep 17 00:00:00 2001 From: Jonathan Hu Date: Wed, 20 Feb 2019 14:31:18 -0700 Subject: [PATCH 7/8] MueLu: Driver: flush all output streams --- packages/muelu/test/scaling/Driver.cpp | 1 + 1 file changed, 1 insertion(+) diff --git a/packages/muelu/test/scaling/Driver.cpp b/packages/muelu/test/scaling/Driver.cpp index ba99ad986ff3..fc054f4af21f 100644 --- a/packages/muelu/test/scaling/Driver.cpp +++ b/packages/muelu/test/scaling/Driver.cpp @@ -513,6 +513,7 @@ int main_(Teuchos::CommandLineProcessor &clp, Xpetra::UnderlyingLib& lib, int ar stop = true; } } + fflush(NULL); comm->barrier(); } while (!stop); From 1ad41f814212cee4f78cdcd4f8aae29489c93a36 Mon Sep 17 00:00:00 2001 From: Chris Siefert Date: Wed, 20 Feb 2019 20:57:54 -0700 Subject: [PATCH 8/8] MueLu: Fixing output --- packages/muelu/test/scaling/Driver.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/packages/muelu/test/scaling/Driver.cpp b/packages/muelu/test/scaling/Driver.cpp index fc054f4af21f..022ee39e11f2 100644 --- a/packages/muelu/test/scaling/Driver.cpp +++ b/packages/muelu/test/scaling/Driver.cpp @@ -440,7 +440,7 @@ int main_(Teuchos::CommandLineProcessor &clp, Xpetra::UnderlyingLib& lib, int ar #endif out2<<"*********** MueLu ParameterList ***********"< H; RCP Prec;