diff --git a/packages/belos/src/BelosBiCGStabSolMgr.hpp b/packages/belos/src/BelosBiCGStabSolMgr.hpp index d74b2fc8c24b..16f95e573954 100644 --- a/packages/belos/src/BelosBiCGStabSolMgr.hpp +++ b/packages/belos/src/BelosBiCGStabSolMgr.hpp @@ -782,6 +782,15 @@ ReturnType BiCGStabSolMgr::solve () "Belos::BiCGStabSolMgr::solve(): Invalid return from BiCGStabIter::iterate()."); } } + catch (const StatusTestNaNError& e) { + // A NaN was detected in the solver. Set the solution to zero and return unconverged. + achievedTol_ = MT::one(); + Teuchos::RCP X = problem_->getLHS(); + MVT::MvInit( *X, SCT::zero() ); + printer_->stream(Warnings) << "Belos::BiCGStabSolMgr::solve(): Warning! NaN has been detected!" + << std::endl; + return Unconverged; + } catch (const std::exception &e) { printer_->stream(Errors) << "Error! Caught std::exception in BiCGStabIter::iterate() at iteration " << bicgstab_iter->getNumIters() << std::endl diff --git a/packages/belos/src/BelosBlockCGSolMgr.hpp b/packages/belos/src/BelosBlockCGSolMgr.hpp index 33aab3d90679..bd252fe10f10 100644 --- a/packages/belos/src/BelosBlockCGSolMgr.hpp +++ b/packages/belos/src/BelosBlockCGSolMgr.hpp @@ -1015,6 +1015,15 @@ ReturnType BlockCGSolMgr::solve() { "to the Belos developers."); } } + catch (const StatusTestNaNError& e) { + // A NaN was detected in the solver. Set the solution to zero and return unconverged. + achievedTol_ = MT::one(); + Teuchos::RCP X = problem_->getLHS(); + MVT::MvInit( *X, SCT::zero() ); + printer_->stream(Warnings) << "Belos::BlockCGSolMgr::solve(): Warning! NaN has been detected!" + << std::endl; + return Unconverged; + } catch (const std::exception &e) { std::ostream& err = printer_->stream (Errors); err << "Error! Caught std::exception in CGIteration::iterate() at " diff --git a/packages/belos/src/BelosBlockGCRODRSolMgr.hpp b/packages/belos/src/BelosBlockGCRODRSolMgr.hpp index e990f300d1fd..8447ad079b23 100644 --- a/packages/belos/src/BelosBlockGCRODRSolMgr.hpp +++ b/packages/belos/src/BelosBlockGCRODRSolMgr.hpp @@ -1964,7 +1964,7 @@ ReturnType BlockGCRODRSolMgr::solve() { tmpNumBlocks = dim / blockSize_; // Allow for a good breakdown. else{ tmpNumBlocks = ( dim - blockSize_) / blockSize_; // Allow for restarting. - printer_->stream(Warnings) << "Belos::BlockGmresSolMgr::solve(): Warning! Requested Krylov subspace dimension is larger than operator dimension!" + printer_->stream(Warnings) << "Belos::BlockGCRODRSolMgr::solve(): Warning! Requested Krylov subspace dimension is larger than operator dimension!" << std::endl << "The maximum number of blocks allowed for the Krylov subspace will be adjusted to " << tmpNumBlocks << std::endl; primeList.set("Num Blocks",Teuchos::as(tmpNumBlocks)); } @@ -2019,7 +2019,7 @@ ReturnType BlockGCRODRSolMgr::solve() { if ( expConvTest_->getLOADetected() ) { // we don't have convergence loaDetected_ = true; - printer_->stream(Warnings) << "Belos::BlockGmresSolMgr::solve(): Warning! Solver has experienced a loss of accuracy!" << std::endl; + printer_->stream(Warnings) << "Belos::BlockGCRODRSolMgr::solve(): Warning! Solver has experienced a loss of accuracy!" << std::endl; } } // ******************************************* @@ -2054,6 +2054,15 @@ ReturnType BlockGCRODRSolMgr::solve() { isConverged = false; } } // end catch (const GmresIterationOrthoFailure &e) + catch (const StatusTestNaNError& e) { + // A NaN was detected in the solver. Set the solution to zero and return unconverged. + achievedTol_ = MT::one(); + Teuchos::RCP X = problem_->getLHS(); + MVT::MvInit( *X, SCT::zero() ); + printer_->stream(Warnings) << "Belos::BlockGCRODRSolMgr::solve(): Warning! NaN has been detected!" + << std::endl; + return Unconverged; + } catch (const std::exception &e) { printer_->stream(Errors) << "Error! Caught std::exception in BlockGmresIter::iterate() at iteration " << block_gmres_iter->getNumIters() << std::endl diff --git a/packages/belos/src/BelosBlockGmresSolMgr.hpp b/packages/belos/src/BelosBlockGmresSolMgr.hpp index 8204640d1167..6049cff14539 100644 --- a/packages/belos/src/BelosBlockGmresSolMgr.hpp +++ b/packages/belos/src/BelosBlockGmresSolMgr.hpp @@ -1142,6 +1142,15 @@ ReturnType BlockGmresSolMgr::solve() { break; } } + catch (const StatusTestNaNError& e) { + // A NaN was detected in the solver. Set the solution to zero and return unconverged. + achievedTol_ = MT::one(); + Teuchos::RCP X = problem_->getLHS(); + MVT::MvInit( *X, SCT::zero() ); + printer_->stream(Warnings) << "Belos::BlockGmresSolMgr::solve(): Warning! NaN has been detected!" + << std::endl; + return Unconverged; + } catch (const std::exception &e) { printer_->stream(Errors) << "Error! Caught std::exception in BlockGmresIter::iterate() at iteration " << block_gmres_iter->getNumIters() << std::endl diff --git a/packages/belos/src/BelosFixedPointSolMgr.hpp b/packages/belos/src/BelosFixedPointSolMgr.hpp index e3f83f13b63c..92f88d6a1c48 100644 --- a/packages/belos/src/BelosFixedPointSolMgr.hpp +++ b/packages/belos/src/BelosFixedPointSolMgr.hpp @@ -744,6 +744,15 @@ ReturnType FixedPointSolMgr::solve() { "to the Belos developers."); } } + catch (const StatusTestNaNError& e) { + // A NaN was detected in the solver. Set the solution to zero and return unconverged. + achievedTol_ = MT::one(); + Teuchos::RCP X = problem_->getLHS(); + MVT::MvInit( *X, SCT::zero() ); + printer_->stream(Warnings) << "Belos::FixedPointSolMgr::solve(): Warning! NaN has been detected!" + << std::endl; + return Unconverged; + } catch (const std::exception &e) { std::ostream& err = printer_->stream (Errors); err << "Error! Caught std::exception in FixedPointIteration::iterate() at " diff --git a/packages/belos/src/BelosGCRODRSolMgr.hpp b/packages/belos/src/BelosGCRODRSolMgr.hpp index 2344543dbe46..cd1836f72c69 100644 --- a/packages/belos/src/BelosGCRODRSolMgr.hpp +++ b/packages/belos/src/BelosGCRODRSolMgr.hpp @@ -1518,6 +1518,15 @@ ReturnType GCRODRSolMgr::solve() { if (convTest_->getStatus() == Passed) primeConverged = true; } + catch (const StatusTestNaNError& e) { + // A NaN was detected in the solver. Set the solution to zero and return unconverged. + achievedTol_ = MT::one(); + Teuchos::RCP X = problem_->getLHS(); + MVT::MvInit( *X, SCT::zero() ); + printer_->stream(Warnings) << "Belos::GCRODRSolMgr::solve(): Warning! NaN has been detected!" + << std::endl; + return Unconverged; + } catch (const std::exception &e) { printer_->stream(Errors) << "Error! Caught exception in GCRODRIter::iterate() at iteration " << gcrodr_prime_iter->getNumIters() << std::endl diff --git a/packages/belos/src/BelosMinresSolMgr.hpp b/packages/belos/src/BelosMinresSolMgr.hpp index be756a93e7a5..8cec45f7e04a 100644 --- a/packages/belos/src/BelosMinresSolMgr.hpp +++ b/packages/belos/src/BelosMinresSolMgr.hpp @@ -763,7 +763,17 @@ namespace Belos { "nor reached the maximum number of iterations " << maxIters_ << ". That means something went wrong."); } - } catch (const std::exception &e) { + } + catch (const StatusTestNaNError& e) { + // A NaN was detected in the solver. Set the solution to zero and return unconverged. + achievedTol_ = MST::one(); + Teuchos::RCP X = problem_->getLHS(); + MVT::MvInit( *X, SCT::zero() ); + printer_->stream(Warnings) << "Belos::MinresSolMgr::solve(): Warning! NaN has been detected!" + << std::endl; + return Unconverged; + } + catch (const std::exception &e) { printer_->stream (Errors) << "Error! Caught std::exception in MinresIter::iterate() at " << "iteration " << minres_iter->getNumIters() << endl diff --git a/packages/belos/src/BelosPCPGSolMgr.hpp b/packages/belos/src/BelosPCPGSolMgr.hpp index a09a26f21398..bd77822a8cf8 100644 --- a/packages/belos/src/BelosPCPGSolMgr.hpp +++ b/packages/belos/src/BelosPCPGSolMgr.hpp @@ -891,6 +891,15 @@ ReturnType PCPGSolMgr::solve() { "Belos::PCPGSolMgr::solve(): Invalid return from PCPGIter::iterate()."); } // end if } // end try + catch (const StatusTestNaNError& e) { + // A NaN was detected in the solver. Set the solution to zero and return unconverged. + achievedTol_ = MT::one(); + Teuchos::RCP X = problem_->getLHS(); + MVT::MvInit( *X, SCT::zero() ); + printer_->stream(Warnings) << "Belos::PCPG::solve(): Warning! NaN has been detected!" + << std::endl; + return Unconverged; + } catch (const std::exception &e) { printer_->stream(Errors) << "Error! Caught exception in PCPGIter::iterate() at iteration " << pcpg_iter->getNumIters() << std::endl diff --git a/packages/belos/src/BelosPseudoBlockCGSolMgr.hpp b/packages/belos/src/BelosPseudoBlockCGSolMgr.hpp index e6b8c019d7cc..dc21a4d802c3 100644 --- a/packages/belos/src/BelosPseudoBlockCGSolMgr.hpp +++ b/packages/belos/src/BelosPseudoBlockCGSolMgr.hpp @@ -920,6 +920,15 @@ ReturnType PseudoBlockCGSolMgr::solve () "Belos::PseudoBlockCGSolMgr::solve(): Invalid return from PseudoBlockCGIter::iterate()."); } } + catch (const StatusTestNaNError& e) { + // A NaN was detected in the solver. Set the solution to zero and return unconverged. + achievedTol_ = MT::one(); + Teuchos::RCP X = problem_->getLHS(); + MVT::MvInit( *X, SCT::zero() ); + printer_->stream(Warnings) << "Belos::PseudoBlockCGSolMgr::solve(): Warning! NaN has been detected!" + << std::endl; + return Unconverged; + } catch (const std::exception &e) { printer_->stream(Errors) << "Error! Caught std::exception in PseudoBlockCGIter::iterate() at iteration " << block_cg_iter->getNumIters() << std::endl diff --git a/packages/belos/src/BelosPseudoBlockGmresSolMgr.hpp b/packages/belos/src/BelosPseudoBlockGmresSolMgr.hpp index c09c3f55e6be..7149a9504fc7 100644 --- a/packages/belos/src/BelosPseudoBlockGmresSolMgr.hpp +++ b/packages/belos/src/BelosPseudoBlockGmresSolMgr.hpp @@ -1469,6 +1469,15 @@ ReturnType PseudoBlockGmresSolMgr::solve() { isConverged = false; break; } + catch (const StatusTestNaNError& e) { + // A NaN was detected in the solver. Set the solution to zero and return unconverged. + achievedTol_ = MT::one(); + Teuchos::RCP X = problem_->getLHS(); + MVT::MvInit( *X, SCT::zero() ); + printer_->stream(Warnings) << "Belos::PseudoBlockGmresSolMgr::solve(): Warning! NaN has been detected!" + << std::endl; + return Unconverged; + } catch (const std::exception &e) { printer_->stream(Errors) << "Error! Caught std::exception in PseudoBlockGmresIter::iterate() at iteration " << block_gmres_iter->getNumIters() << std::endl diff --git a/packages/belos/src/BelosPseudoBlockTFQMRSolMgr.hpp b/packages/belos/src/BelosPseudoBlockTFQMRSolMgr.hpp index 6c2a3eaffa25..9f84d9a96aee 100644 --- a/packages/belos/src/BelosPseudoBlockTFQMRSolMgr.hpp +++ b/packages/belos/src/BelosPseudoBlockTFQMRSolMgr.hpp @@ -807,6 +807,15 @@ ReturnType PseudoBlockTFQMRSolMgr::solve() { "Belos::PseudoBlockTFQMRSolMgr::solve(): Invalid return from PseudoBlockTFQMRIter::iterate()."); } } + catch (const StatusTestNaNError& e) { + // A NaN was detected in the solver. Set the solution to zero and return unconverged. + achievedTol_ = MT::one(); + Teuchos::RCP X = problem_->getLHS(); + MVT::MvInit( *X, SCT::zero() ); + printer_->stream(Warnings) << "Belos::PseudoBlockTFQMRSolMgr::solve(): Warning! NaN has been detected!" + << std::endl; + return Unconverged; + } catch (const std::exception &e) { printer_->stream(Errors) << "Error! Caught std::exception in PseudoBlockTFQMRIter::iterate() at iteration " << block_tfqmr_iter->getNumIters() << std::endl diff --git a/packages/belos/src/BelosRCGSolMgr.hpp b/packages/belos/src/BelosRCGSolMgr.hpp index fa624d6eecc4..a15bc3f519a3 100644 --- a/packages/belos/src/BelosRCGSolMgr.hpp +++ b/packages/belos/src/BelosRCGSolMgr.hpp @@ -1748,6 +1748,15 @@ ReturnType RCGSolMgr::solve() { "Belos::RCGSolMgr::solve(): Invalid return from RCGIter::iterate()."); } } + catch (const StatusTestNaNError& e) { + // A NaN was detected in the solver. Set the solution to zero and return unconverged. + achievedTol_ = MT::one(); + Teuchos::RCP X = problem_->getLHS(); + MVT::MvInit( *X, SCT::zero() ); + printer_->stream(Warnings) << "Belos::RCGSolMgr::solve(): Warning! NaN has been detected!" + << std::endl; + return Unconverged; + } catch (const std::exception &e) { printer_->stream(Errors) << "Error! Caught std::exception in RCGIter::iterate() at iteration " << rcg_iter->getNumIters() << std::endl diff --git a/packages/belos/src/BelosStatusTest.hpp b/packages/belos/src/BelosStatusTest.hpp index e6e841c7e856..a36d265a4934 100644 --- a/packages/belos/src/BelosStatusTest.hpp +++ b/packages/belos/src/BelosStatusTest.hpp @@ -73,6 +73,9 @@ namespace Belos { class StatusTestError : public BelosError {public: StatusTestError(const std::string& what_arg) : BelosError(what_arg) {}}; + class StatusTestNaNError : public StatusTestError + {public: StatusTestNaNError(const std::string& what_arg) : StatusTestError(what_arg) {}}; + //@} template diff --git a/packages/belos/src/BelosStatusTestGenResNorm.hpp b/packages/belos/src/BelosStatusTestGenResNorm.hpp index d5c7860dc5d2..ce848230d664 100644 --- a/packages/belos/src/BelosStatusTestGenResNorm.hpp +++ b/packages/belos/src/BelosStatusTestGenResNorm.hpp @@ -571,7 +571,7 @@ StatusType StatusTestGenResNorm::checkStatus( Iteration,Thyr } else { // Throw an std::exception if a NaN is found. status_ = Failed; - TEUCHOS_TEST_FOR_EXCEPTION(true,StatusTestError,"StatusTestGenResSubNorm::checkStatus(): NaN has been detected."); + TEUCHOS_TEST_FOR_EXCEPTION(true,StatusTestNaNError,"StatusTestGenResSubNorm::checkStatus(): NaN has been detected."); } } } diff --git a/packages/belos/src/BelosStatusTestImpResNorm.hpp b/packages/belos/src/BelosStatusTestImpResNorm.hpp index f2f7d342561f..39b3e95afb2c 100644 --- a/packages/belos/src/BelosStatusTestImpResNorm.hpp +++ b/packages/belos/src/BelosStatusTestImpResNorm.hpp @@ -630,7 +630,7 @@ checkStatus (Iteration* iSolver) // tolerance is NaN; we assume the former. We also mark the // test as failed, in case you want to catch the exception. status_ = Failed; - TEUCHOS_TEST_FOR_EXCEPTION(true, StatusTestError, "Belos::" + TEUCHOS_TEST_FOR_EXCEPTION(true, StatusTestNaNError, "Belos::" "StatusTestImpResNorm::checkStatus(): One or more of the current " "implicit residual norms is NaN."); } diff --git a/packages/belos/src/BelosTFQMRSolMgr.hpp b/packages/belos/src/BelosTFQMRSolMgr.hpp index 36b2000b3983..db2417005a98 100644 --- a/packages/belos/src/BelosTFQMRSolMgr.hpp +++ b/packages/belos/src/BelosTFQMRSolMgr.hpp @@ -752,6 +752,15 @@ ReturnType TFQMRSolMgr::solve() { "Belos::TFQMRSolMgr::solve(): Invalid return from TFQMRIter::iterate()."); } } + catch (const StatusTestNaNError& e) { + // A NaN was detected in the solver. Set the solution to zero and return unconverged. + achievedTol_ = MT::one(); + Teuchos::RCP X = problem_->getLHS(); + MVT::MvInit( *X, SCT::zero() ); + printer_->stream(Warnings) << "Belos::TFQMRSolMgr::solve(): Warning! NaN has been detected!" + << std::endl; + return Unconverged; + } catch (const std::exception &e) { printer_->stream(Errors) << "Error! Caught std::exception in TFQMRIter::iterate() at iteration " << tfqmr_iter->getNumIters() << std::endl diff --git a/packages/belos/tpetra/test/LinearSolverFactory/CMakeLists.txt b/packages/belos/tpetra/test/LinearSolverFactory/CMakeLists.txt index 810c9f13b611..c392ebfbda73 100644 --- a/packages/belos/tpetra/test/LinearSolverFactory/CMakeLists.txt +++ b/packages/belos/tpetra/test/LinearSolverFactory/CMakeLists.txt @@ -6,6 +6,13 @@ TRIBITS_ADD_EXECUTABLE_AND_TEST( COMM serial mpi ) +TRIBITS_ADD_EXECUTABLE_AND_TEST( + SolverFactoryNaN + SOURCES SolverFactoryNaN.cpp ${TEUCHOS_STD_UNIT_TEST_MAIN} + ARGS + COMM serial mpi + ) + TRIBITS_ADD_EXECUTABLE_AND_TEST( CustomSolverFactory SOURCES CustomSolverFactory.cpp ${TEUCHOS_STD_UNIT_TEST_MAIN} diff --git a/packages/belos/tpetra/test/LinearSolverFactory/SolverFactoryNaN.cpp b/packages/belos/tpetra/test/LinearSolverFactory/SolverFactoryNaN.cpp new file mode 100644 index 000000000000..c4d4d7ba2ebe --- /dev/null +++ b/packages/belos/tpetra/test/LinearSolverFactory/SolverFactoryNaN.cpp @@ -0,0 +1,383 @@ +//@HEADER +// ************************************************************************ +// +// Belos: Block Linear Solvers Package +// Copyright 2004 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. +// +// ************************************************************************ +//@HEADER + +#include "Teuchos_UnitTestHarness.hpp" +#include "Teuchos_TypeNameTraits.hpp" +#include "BelosSolverFactory.hpp" +#include "BelosTpetraAdapter.hpp" +#include "Tpetra_CrsMatrix.hpp" +#include "Tpetra_Core.hpp" +#include "Tpetra_MultiVector.hpp" +#include "TpetraCore_ETIHelperMacros.h" +#include +#include + +// mfh 12 Oct 2017: Test that the Tpetra specialization of +// Belos::SolverFactory builds and runs without throwing. +// See e.g., Trilinos GitHub issue #754. + +namespace { // (anonymous) + +// Create a very simple square test matrix. We use the identity +// matrix here. The point of this test is NOT to exercise the solver; +// it's just to check that Belos' LinearSolverFactory can create +// working solvers. Belos has more rigorous tests for its solvers. +template +Teuchos::RCP > +createTestMatrix (Teuchos::FancyOStream& out, + const Teuchos::RCP >& comm, + const Tpetra::global_size_t gblNumRows) +{ + using Teuchos::Comm; + using Teuchos::RCP; + using Teuchos::rcp; + using std::endl; + typedef Tpetra::CrsMatrix MAT; + typedef Tpetra::Map map_type; + typedef Teuchos::ScalarTraits STS; + typedef Teuchos::ScalarTraits MTS; + + Teuchos::OSTab tab0 (out); + out << "Create test matrix with " << gblNumRows << " row(s)" << endl; + Teuchos::OSTab tab1 (out); + + TEUCHOS_TEST_FOR_EXCEPTION + ( gblNumRows == 0, std::invalid_argument, "gblNumRows = 0"); + + const GO indexBase = 0; + RCP rowMap (new map_type (gblNumRows, indexBase, comm)); + // For this particular row matrix, the row Map and the column Map + // are the same. Giving a column Map to CrsMatrix's constructor + // lets us use local indices. + RCP colMap = rowMap; + const size_t maxNumEntPerRow = 1; + RCP A (new MAT (rowMap, colMap, maxNumEntPerRow)); + + if (rowMap->getLocalNumElements () != 0) { + Teuchos::Array vals (1); + Teuchos::Array inds (1); + for (LO lclRow = rowMap->getMinLocalIndex (); + lclRow <= rowMap->getMaxLocalIndex (); ++lclRow) { + inds[0] = lclRow; + vals[0] = MTS::zero () / MTS::zero (); // Fill matrix with NaNs + A->insertLocalValues (lclRow, inds (), vals ()); + } + } + + RCP domMap = rowMap; + RCP ranMap = rowMap; + A->fillComplete (domMap, ranMap); + return A; +} + +// Create a very simple square test linear system (matrix, right-hand +// side(s), and exact solution(s). We use the identity matrix here. +// The point of this test is NOT to exercise the preconditioner; it's +// just to check that its LinearSolverFactory can create working +// preconditioners. Belos has more rigorous tests for each of its +// preconditioners. +template +void +createTestProblem (Teuchos::FancyOStream& out, + Teuchos::RCP >& X, + Teuchos::RCP >& A, + Teuchos::RCP >& B, + const Teuchos::RCP >& comm, + const Tpetra::global_size_t gblNumRows, + const size_t numVecs) +{ + using Teuchos::Comm; + using Teuchos::RCP; + using Teuchos::rcp; + using std::endl; + typedef Tpetra::MultiVector MV; + typedef Teuchos::ScalarTraits STS; + + A = createTestMatrix (out, comm, gblNumRows); + X = rcp (new MV (A->getDomainMap (), numVecs)); + B = rcp (new MV (A->getRangeMap (), numVecs)); + + B->putScalar (STS::one()); + X->putScalar (STS::zero()); +} + +template +void +testSolver (Teuchos::FancyOStream& out, + bool& success, + const Teuchos::RCP >& X, + const Teuchos::RCP >& A, + const Teuchos::RCP >& B, + Tpetra::MultiVector& /* X_exact */, + const std::string& solverName) +{ + using Teuchos::Comm; + using Teuchos::RCP; + using Teuchos::rcp; + using std::endl; + typedef Tpetra::Operator OP; + typedef Tpetra::MultiVector MV; + typedef Teuchos::ScalarTraits STS; + typedef Teuchos::ScalarTraits MTS; + + Teuchos::OSTab tab0 (out); + out << "Test solver \"" << solverName << "\" from Belos package" << endl; + Teuchos::OSTab tab1 (out); + + RCP > solver; + Belos::SolverFactory factory; + + // Set up Belos solver parameters. + Teuchos::RCP belosList = Teuchos::parameterList (solverName); + belosList->set ("Verbosity", Belos::Errors + Belos::Warnings); + + try { + solver = factory.create (solverName, belosList); + } catch (std::exception& e) { + out << "*** FAILED: Belos::SolverFactory::create threw an exception: " + << e.what () << endl; + success = false; + return; + } + TEST_EQUALITY_CONST( solver.is_null (), false ); + if (solver.is_null ()) { + out << "*** FAILED to create solver \"" << solverName + << "\" from Belos package" << endl; + success = false; + return; + } + + out << "Create the Belos::LinearProblem to solve" << endl; + typedef Belos::LinearProblem linear_problem_type; + X->putScalar (STS::zero ()); + RCP problem (new linear_problem_type (A, X, B)); + problem->setProblem (); + + out << "Set up the solver" << endl; + solver->setProblem (problem); + + out << "Apply solver to \"solve\" AX=B for X. Belos already has solver tests;" + " the point is to check that solve() doesn't throw if it encounters a NaN." << endl; + Belos::ReturnType ret; + + try { + ret = solver->solve (); + } catch (std::exception& e) { + out << "*** FAILED: Belos::SolverFactory::solve threw an unexpected exception: " + << e.what () << endl; + success = false; + return; + } + + // Check that the solution vector is zeros, the achieved tolerance is 1, and the solver return Unconverged. + bool nonZeroX = false; + std::vector normX( X->getNumVectors () ); + Teuchos::ArrayView normAV( normX ); + X->norm2 ( normAV( 0, X->getNumVectors ()) ); + for ( int i=0; i<(int)(X->getNumVectors ()); ++i ) + { + if ( normX[i] != MTS::zero() ) + nonZeroX = true; + } + if ( nonZeroX || (ret != Belos::Unconverged) || ( solver->achievedTol() != MTS::one() ) ) { + success = false; + return; + } + + success = true; +} + +template +void +testCreatingSolver (Teuchos::FancyOStream& out, + bool& success, + const std::string& solverName) +{ + using Teuchos::RCP; + using Teuchos::rcp; + using Teuchos::TypeNameTraits; + using std::endl; + using MV = Tpetra::MultiVector; + using OP = Tpetra::Operator; + + Teuchos::OSTab tab0 (out); + out << "Test Belos solver \"" << solverName + << "\" for Tpetra ::name () + << ", LO=" << TypeNameTraits::name () + << ", GO=" << TypeNameTraits::name () + << ", NT=" << TypeNameTraits::name () + << ">" << endl; + Teuchos::OSTab tab1 (out); + + Belos::SolverFactory factory; + RCP > solver; + TEST_NOTHROW( solver = factory.create (solverName, Teuchos::null) ); + TEST_ASSERT( ! solver.is_null () ); +} + +// +// The actual unit tests start here. +// + +TEUCHOS_UNIT_TEST_TEMPLATE_4_DECL( SolverFactory, CreateSolvers, SC, LO, GO, NT ) +{ + using std::endl; + + out << "Test Belos::SolverFactory with Tpetra for all solvers" << endl; + Teuchos::OSTab tab1 (out); + const std::vector solverNames {{ + "BICGSTAB", + "BLOCK CG", + "BLOCK GMRES", + "TPETRA CG PIPELINE", + "TPETRA CG SINGLE REDUCE", + "FIXED POINT", + "GCRODR", + "TPETRA GMRES PIPELINE", + "HYBRID BLOCK GMRES", // GmresPoly + "TPETRA GMRES SINGLE REDUCE", + "LSQR", + "MINRES", + "PCPG", + "PSEUDOBLOCK CG", + "PSEUDOBLOCK GMRES", + "PSEUDOBLOCK TFQMR", + "TFQMR" + }}; + + for (const std::string& solverName : solverNames) { + const bool isComplex = Teuchos::ScalarTraits::isComplex; + if (isComplex && (solverName == "LSQR" || solverName == "PCPG")) { + continue; // solver not implemented for complex Scalar types + } + testCreatingSolver (out, success, solverName); + if (! success) { + return; + } + } +} + +TEUCHOS_UNIT_TEST_TEMPLATE_4_DECL( SolverFactory, CreateAndSolve, SC, LO, GO, NT ) +{ + using Teuchos::Comm; + using Teuchos::RCP; + using Teuchos::rcp; + using Teuchos::TypeNameTraits; + using std::endl; + typedef Tpetra::CrsMatrix MAT; + typedef Tpetra::MultiVector MV; + typedef Tpetra::Operator OP; + + RCP > comm = Tpetra::getDefaultComm (); + const Tpetra::global_size_t gblNumRows = comm->getSize () * 10; + const size_t numVecs = 3; + + out << "Create test problem" << endl; + RCP X_exact, B; + RCP A; + createTestProblem (out, X_exact, A, B, comm, gblNumRows, numVecs); + + RCP X = rcp (new MV (X_exact->getMap (), numVecs)); + + Belos::SolverFactory factory; + // FIXME (mfh 23 Aug 2015) Not all Belos solvers can handle solves + // with the identity matrix. BiCGSTAB might need a bit of work, for + // example. I'm not so worried about that for now but we should go + // back and revisit this at some point. + // + // Teuchos::Array solverNames = factory.supportedSolverNames (); + // const int numSolvers = static_cast (solverNames.size ()); + const char* solverNames[10] = { + "BICGSTAB", + "BLOCK CG", + "BLOCK GMRES", + "FIXED POINT", + "GCRODR", + "MINRES", + "PSEUDOBLOCK CG", + "PSEUDOBLOCK GMRES", + "PSEUDOBLOCK TFQMR", + "TFQMR"}; + const int numSolvers = 10; + + int numSolversTested = 0; + for (int k = 0; k < numSolvers; ++k) { + const std::string solverName (solverNames[k]); + + // Use Belos' factory to tell us whether the factory supports the + // given combination of template parameters. If create() throws, + // the combination is not supported. + bool skip = false; + try { + (void) factory.create (solverName, Teuchos::null); + } + catch (...) { + skip = true; + } + if (! skip) { + testSolver (out, success, X, A, B, *X_exact, solverName); + ++numSolversTested; + } + } + + out << "Tested " << numSolversTested << " solver(s) of " << numSolvers + << endl; + if (numSolversTested == 0) { + out << "*** ERROR: Tested no solvers for template parameters" + << "SC = " << TypeNameTraits::name () + << ", LO = " << TypeNameTraits::name () + << ", GO = " << TypeNameTraits::name () + << ", NT = " << TypeNameTraits::name () << endl; + success = false; + } +} + +// Define typedefs that make the Tpetra macros work. +TPETRA_ETI_MANGLING_TYPEDEFS() + +// Macro that instantiates the unit test +#define LCLINST( SC, LO, GO, NT ) \ + TEUCHOS_UNIT_TEST_TEMPLATE_4_INSTANT( SolverFactory, CreateSolvers, SC, LO, GO, NT ) \ + TEUCHOS_UNIT_TEST_TEMPLATE_4_INSTANT( SolverFactory, CreateAndSolve, SC, LO, GO, NT ) + +// Tpetra's ETI will instantiate the unit test for all enabled type +// combinations. +TPETRA_INSTANTIATE_SLGN_NO_ORDINAL_SCALAR( LCLINST ) + +} // namespace (anonymous) diff --git a/packages/stratimikos/adapters/belos/src/Thyra_BelosLinearOpWithSolveFactory_def.hpp b/packages/stratimikos/adapters/belos/src/Thyra_BelosLinearOpWithSolveFactory_def.hpp index 322c96349e9a..60f13f254e6f 100644 --- a/packages/stratimikos/adapters/belos/src/Thyra_BelosLinearOpWithSolveFactory_def.hpp +++ b/packages/stratimikos/adapters/belos/src/Thyra_BelosLinearOpWithSolveFactory_def.hpp @@ -63,7 +63,7 @@ #include "BelosFixedPointSolMgr.hpp" #include "BelosThyraAdapter.hpp" -#ifdef HAVE_BELOS_TPETRA +#if defined(HAVE_BELOS_TPETRA) && defined(HAVE_STRATIMIKOS_THYRATPETRAADAPTERS) #include "Thyra_BelosTpetrasSolverAdapter.hpp" #endif @@ -109,7 +109,7 @@ const std::string BelosLinearOpWithSolveFactory::BiCGStab_name = "BiCGSt template const std::string BelosLinearOpWithSolveFactory::FixedPoint_name = "Fixed Point"; -#ifdef HAVE_BELOS_TPETRA +#if defined(HAVE_BELOS_TPETRA) && defined(HAVE_STRATIMIKOS_THYRATPETRAADAPTERS) template const std::string BelosLinearOpWithSolveFactory::TpetraGmres_name = "TPETRA GMRES"; template @@ -419,7 +419,7 @@ Teuchos::ValidatorXMLConverterDB::addConverter( "TFQMR", "BiCGStab", "Fixed Point" -#ifdef HAVE_BELOS_TPETRA +#if defined(HAVE_BELOS_TPETRA) && defined(HAVE_STRATIMIKOS_THYRATPETRAADAPTERS) ,"TPETRA GMRES", "TPETRA GMRES PIPELINE", "TPETRA GMRES SINGLE REDUCE", @@ -471,7 +471,7 @@ Teuchos::ValidatorXMLConverterDB::addConverter( "Fixed point iteration" -#ifdef HAVE_BELOS_TPETRA +#if defined(HAVE_BELOS_TPETRA) && defined(HAVE_STRATIMIKOS_THYRATPETRAADAPTERS) ,"Native Tpetra implementation of GMRES", "Native Tpetra implementation of pipeline GMRES", @@ -493,7 +493,7 @@ Teuchos::ValidatorXMLConverterDB::addConverter( SOLVER_TYPE_TFQMR, SOLVER_TYPE_BICGSTAB, SOLVER_TYPE_FIXEDPOINT -#ifdef HAVE_BELOS_TPETRA +#if defined(HAVE_BELOS_TPETRA) && defined(HAVE_STRATIMIKOS_THYRATPETRAADAPTERS) ,SOLVER_TYPE_TPETRA_GMRES, SOLVER_TYPE_TPETRA_GMRES_PIPELINE, SOLVER_TYPE_TPETRA_GMRES_SINGLE_REDUCE, @@ -583,7 +583,7 @@ Teuchos::ValidatorXMLConverterDB::addConverter( *mgr.getValidParameters() ); } -#ifdef HAVE_BELOS_TPETRA +#if defined(HAVE_BELOS_TPETRA) && defined(HAVE_STRATIMIKOS_THYRATPETRAADAPTERS) { Thyra::BelosTpetraGmres mgr; solverTypesSL.sublist(TpetraGmres_name).setParameters( @@ -1007,7 +1007,7 @@ void BelosLinearOpWithSolveFactory::initializeOpImpl( } break; } -#ifdef HAVE_BELOS_TPETRA +#if defined(HAVE_BELOS_TPETRA) && defined(HAVE_STRATIMIKOS_THYRATPETRAADAPTERS) case SOLVER_TYPE_TPETRA_GMRES: { // Get the PL diff --git a/packages/stratimikos/src/Stratimikos_LinearSolverBuilder_def.hpp b/packages/stratimikos/src/Stratimikos_LinearSolverBuilder_def.hpp index 62a7d205dfb1..b14e5490ecf2 100644 --- a/packages/stratimikos/src/Stratimikos_LinearSolverBuilder_def.hpp +++ b/packages/stratimikos/src/Stratimikos_LinearSolverBuilder_def.hpp @@ -66,7 +66,7 @@ #ifdef HAVE_STRATIMIKOS_IFPACK # include "Thyra_IfpackPreconditionerFactory.hpp" #endif -#ifdef HAVE_STRATIMIKOS_IFPACK2 +#if defined(HAVE_STRATIMIKOS_IFPACK2) && defined(HAVE_STRATIMIKOS_THYRATPETRAADAPTERS) # include "Thyra_Ifpack2PreconditionerFactory.hpp" # include "Tpetra_CrsMatrix.hpp" #endif @@ -541,7 +541,7 @@ void LinearSolverBuilder::initializeDefaults() ); #endif -#ifdef HAVE_STRATIMIKOS_IFPACK2 +#if defined(HAVE_STRATIMIKOS_IFPACK2) && defined(HAVE_STRATIMIKOS_THYRATPETRAADAPTERS) setPreconditioningStrategyFactory( abstractFactoryStd, Thyra::Ifpack2PreconditionerFactory>>(), @@ -624,7 +624,7 @@ void LinearSolverBuilder::initializeDefaults() ); #endif -#ifdef HAVE_STRATIMIKOS_IFPACK2 +#if defined(HAVE_STRATIMIKOS_IFPACK2) && defined(HAVE_STRATIMIKOS_THYRATPETRAADAPTERS) setPreconditioningStrategyFactory( abstractFactoryStd, Thyra::Ifpack2PreconditionerFactory>>(),