42 #ifndef BELOS_TFQMR_SOLMGR_HPP
43 #define BELOS_TFQMR_SOLMGR_HPP
61 #ifdef BELOS_TEUCHOS_TIME_MONITOR
93 template<
class ScalarType,
class MV,
class OP>
164 return Teuchos::tuple(timerSolve_);
244 bool checkStatusTest();
264 static constexpr
int maxIters_default_ = 1000;
265 static constexpr
bool expResTest_default_ =
false;
268 static constexpr
int outputFreq_default_ = -1;
269 static constexpr
const char * impResScale_default_ =
"Norm of Preconditioned Initial Residual";
270 static constexpr
const char * expResScale_default_ =
"Norm of Initial Residual";
271 static constexpr
const char * label_default_ =
"Belos";
274 MagnitudeType convtol_, impTolScale_, achievedTol_;
275 int maxIters_, numIters_;
276 int verbosity_, outputStyle_, outputFreq_;
279 std::string impResScale_, expResScale_;
286 bool isSet_, isSTSet_;
291 template<
class ScalarType,
class MV,
class OP>
293 outputStream_(Teuchos::rcpFromRef(std::cout)),
296 achievedTol_(Teuchos::ScalarTraits<typename Teuchos::ScalarTraits<ScalarType>::magnitudeType>::zero()),
297 maxIters_(maxIters_default_),
299 verbosity_(verbosity_default_),
300 outputStyle_(outputStyle_default_),
301 outputFreq_(outputFreq_default_),
303 expResTest_(expResTest_default_),
304 impResScale_(impResScale_default_),
305 expResScale_(expResScale_default_),
306 label_(label_default_),
313 template<
class ScalarType,
class MV,
class OP>
318 outputStream_(Teuchos::rcpFromRef(std::cout)),
321 achievedTol_(Teuchos::ScalarTraits<typename Teuchos::ScalarTraits<ScalarType>::magnitudeType>::zero()),
322 maxIters_(maxIters_default_),
324 verbosity_(verbosity_default_),
325 outputStyle_(outputStyle_default_),
326 outputFreq_(outputFreq_default_),
328 expResTest_(expResTest_default_),
329 impResScale_(impResScale_default_),
330 expResScale_(expResScale_default_),
331 label_(label_default_),
343 template<
class ScalarType,
class MV,
class OP>
347 if (params_ == Teuchos::null) {
356 maxIters_ = params->
get(
"Maximum Iterations",maxIters_default_);
359 params_->set(
"Maximum Iterations", maxIters_);
360 if (maxIterTest_!=Teuchos::null)
361 maxIterTest_->setMaxIters( maxIters_ );
366 blockSize_ = params->
get(
"Block Size",1);
368 "Belos::TFQMRSolMgr: \"Block Size\" must be 1.");
371 params_->set(
"Block Size", blockSize_);
376 std::string tempLabel = params->
get(
"Timer Label", label_default_);
379 if (tempLabel != label_) {
381 params_->set(
"Timer Label", label_);
382 std::string solveLabel = label_ +
": TFQMRSolMgr total solve time";
383 #ifdef BELOS_TEUCHOS_TIME_MONITOR
391 if (Teuchos::isParameterType<int>(*params,
"Verbosity")) {
392 verbosity_ = params->
get(
"Verbosity", verbosity_default_);
394 verbosity_ = (int)Teuchos::getParameter<Belos::MsgType>(*params,
"Verbosity");
398 params_->set(
"Verbosity", verbosity_);
399 if (printer_ != Teuchos::null)
400 printer_->setVerbosity(verbosity_);
405 if (Teuchos::isParameterType<int>(*params,
"Output Style")) {
406 outputStyle_ = params->
get(
"Output Style", outputStyle_default_);
408 outputStyle_ = (int)Teuchos::getParameter<Belos::OutputType>(*params,
"Output Style");
412 params_->set(
"Output Style", outputStyle_);
418 outputStream_ = Teuchos::getParameter<Teuchos::RCP<std::ostream> >(*params,
"Output Stream");
421 params_->set(
"Output Stream", outputStream_);
422 if (printer_ != Teuchos::null)
423 printer_->setOStream( outputStream_ );
429 outputFreq_ = params->
get(
"Output Frequency", outputFreq_default_);
433 params_->set(
"Output Frequency", outputFreq_);
434 if (outputTest_ != Teuchos::null)
435 outputTest_->setOutputFrequency( outputFreq_ );
439 if (printer_ == Teuchos::null) {
444 if (params->
isParameter(
"Convergence Tolerance")) {
445 if (params->
isType<MagnitudeType> (
"Convergence Tolerance")) {
446 convtol_ = params->
get (
"Convergence Tolerance",
454 params_->set(
"Convergence Tolerance", convtol_);
459 if (params->
isParameter(
"Implicit Tolerance Scale Factor")) {
460 if (params->
isType<MagnitudeType> (
"Implicit Tolerance Scale Factor")) {
461 impTolScale_ = params->
get (
"Implicit Tolerance Scale Factor",
466 impTolScale_ = params->
get (
"Implicit Tolerance Scale Factor",
471 params_->set(
"Implicit Tolerance Scale Factor", impTolScale_);
476 if (params->
isParameter(
"Implicit Residual Scaling")) {
477 std::string tempImpResScale = Teuchos::getParameter<std::string>( *params,
"Implicit Residual Scaling" );
480 if (impResScale_ != tempImpResScale) {
481 impResScale_ = tempImpResScale;
484 params_->set(
"Implicit Residual Scaling", impResScale_);
491 if (params->
isParameter(
"Explicit Residual Scaling")) {
492 std::string tempExpResScale = Teuchos::getParameter<std::string>( *params,
"Explicit Residual Scaling" );
495 if (expResScale_ != tempExpResScale) {
496 expResScale_ = tempExpResScale;
499 params_->set(
"Explicit Residual Scaling", expResScale_);
506 if (params->
isParameter(
"Explicit Residual Test")) {
507 expResTest_ = Teuchos::getParameter<bool>( *params,
"Explicit Residual Test" );
510 params_->set(
"Explicit Residual Test", expResTest_);
511 if (expConvTest_ == Teuchos::null) {
517 if (timerSolve_ == Teuchos::null) {
518 std::string solveLabel = label_ +
": TFQMRSolMgr total solve time";
519 #ifdef BELOS_TEUCHOS_TIME_MONITOR
530 template<
class ScalarType,
class MV,
class OP>
543 Teuchos::rcp(
new StatusTestGenResNorm_t( impTolScale_*convtol_ ) );
545 impConvTest_ = tmpImpConvTest;
550 tmpExpConvTest->defineResForm( StatusTestGenResNorm_t::Explicit,
Belos::TwoNorm );
552 expConvTest_ = tmpExpConvTest;
555 convTest_ =
Teuchos::rcp(
new StatusTestCombo_t( StatusTestCombo_t::SEQ, impConvTest_, expConvTest_ ) );
563 impConvTest_ = tmpImpConvTest;
566 expConvTest_ = impConvTest_;
567 convTest_ = impConvTest_;
569 sTest_ =
Teuchos::rcp(
new StatusTestCombo_t( StatusTestCombo_t::OR, maxIterTest_, convTest_ ) );
573 StatusTestOutputFactory<ScalarType,MV,OP> stoFactory( outputStyle_ );
577 std::string solverDesc =
" TFQMR ";
578 outputTest_->setSolverDesc( solverDesc );
588 template<
class ScalarType,
class MV,
class OP>
601 "The relative residual tolerance that needs to be achieved by the\n"
602 "iterative solver in order for the linear system to be declared converged.");
604 "The scale factor used by the implicit residual test when explicit residual\n"
605 "testing is used. May enable faster convergence when TFQMR bound is too loose.");
606 pl->
set(
"Maximum Iterations", static_cast<int>(maxIters_default_),
607 "The maximum number of block iterations allowed for each\n"
608 "set of RHS solved.");
609 pl->
set(
"Verbosity", static_cast<int>(verbosity_default_),
610 "What type(s) of solver information should be outputted\n"
611 "to the output stream.");
612 pl->
set(
"Output Style", static_cast<int>(outputStyle_default_),
613 "What style is used for the solver information outputted\n"
614 "to the output stream.");
615 pl->
set(
"Output Frequency", static_cast<int>(outputFreq_default_),
616 "How often convergence information should be outputted\n"
617 "to the output stream.");
618 pl->
set(
"Output Stream", Teuchos::rcpFromRef(std::cout),
619 "A reference-counted pointer to the output stream where all\n"
620 "solver output is sent.");
621 pl->
set(
"Explicit Residual Test", static_cast<bool>(expResTest_default_),
622 "Whether the explicitly computed residual should be used in the convergence test.");
623 pl->
set(
"Implicit Residual Scaling", static_cast<const char *>(impResScale_default_),
624 "The type of scaling used in the implicit residual convergence test.");
625 pl->
set(
"Explicit Residual Scaling", static_cast<const char *>(expResScale_default_),
626 "The type of scaling used in the explicit residual convergence test.");
627 pl->
set(
"Timer Label", static_cast<const char *>(label_default_),
628 "The string to use as a prefix for the timer labels.");
636 template<
class ScalarType,
class MV,
class OP>
643 setParameters(Teuchos::parameterList(*getValidParameters()));
647 "Belos::TFQMRSolMgr::solve(): Linear problem is not a valid object.");
650 "Belos::TFQMRSolMgr::solve(): Linear problem is not ready, setProblem() has not been called.");
654 "Belos::TFQMRSolMgr::solve(): Linear problem and requested status tests are incompatible.");
659 int numRHS2Solve = MVT::GetNumberVecs( *(problem_->getRHS()) );
660 int numCurrRHS = blockSize_;
662 std::vector<int> currIdx, currIdx2;
665 currIdx.resize( blockSize_ );
666 currIdx2.resize( blockSize_ );
667 for (
int i=0; i<numCurrRHS; ++i)
668 { currIdx[i] = startPtr+i; currIdx2[i]=i; }
671 problem_->setLSIndex( currIdx );
676 plist.
set(
"Block Size",blockSize_);
679 outputTest_->reset();
682 bool isConverged =
true;
692 #ifdef BELOS_TEUCHOS_TIME_MONITOR
696 while ( numRHS2Solve > 0 ) {
699 std::vector<int> convRHSIdx;
700 std::vector<int> currRHSIdx( currIdx );
701 currRHSIdx.resize(numCurrRHS);
704 tfqmr_iter->resetNumIters();
707 outputTest_->resetNumCalls();
710 Teuchos::RCP<MV> R_0 = MVT::CloneViewNonConst( *(Teuchos::rcp_const_cast<MV>(problem_->getInitPrecResVec())), currIdx );
715 tfqmr_iter->initializeTFQMR(newstate);
721 tfqmr_iter->iterate();
728 if ( convTest_->getStatus() ==
Passed ) {
737 else if ( maxIterTest_->getStatus() ==
Passed ) {
752 "Belos::TFQMRSolMgr::solve(): Invalid return from TFQMRIter::iterate().");
757 achievedTol_ = MT::one();
759 MVT::MvInit( *X, SCT::zero() );
760 printer_->stream(
Warnings) <<
"Belos::TFQMRSolMgr::solve(): Warning! NaN has been detected!"
764 catch (
const std::exception &e) {
765 printer_->stream(
Errors) <<
"Error! Caught std::exception in TFQMRIter::iterate() at iteration "
766 << tfqmr_iter->getNumIters() << std::endl
767 << e.what() << std::endl;
773 problem_->updateSolution( tfqmr_iter->getCurrentUpdate(), true );
776 problem_->setCurrLS();
779 startPtr += numCurrRHS;
780 numRHS2Solve -= numCurrRHS;
781 if ( numRHS2Solve > 0 ) {
782 numCurrRHS = blockSize_;
784 currIdx.resize( blockSize_ );
785 currIdx2.resize( blockSize_ );
786 for (
int i=0; i<numCurrRHS; ++i)
787 { currIdx[i] = startPtr+i; currIdx2[i] = i; }
789 problem_->setLSIndex( currIdx );
792 tfqmr_iter->setBlockSize( blockSize_ );
795 currIdx.resize( numRHS2Solve );
806 #ifdef BELOS_TEUCHOS_TIME_MONITOR
815 numIters_ = maxIterTest_->getNumIters();
828 const std::vector<MagnitudeType>* pTestValues = NULL;
830 pTestValues = expConvTest_->getTestValue();
831 if (pTestValues == NULL || pTestValues->size() < 1) {
832 pTestValues = impConvTest_->getTestValue();
837 pTestValues = impConvTest_->getTestValue();
840 "Belos::TFQMRSolMgr::solve(): The implicit convergence test's "
841 "getTestValue() method returned NULL. Please report this bug to the "
842 "Belos developers.");
844 "Belos::TMQMRSolMgr::solve(): The implicit convergence test's "
845 "getTestValue() method returned a vector of length zero. Please report "
846 "this bug to the Belos developers.");
851 achievedTol_ = *std::max_element (pTestValues->begin(), pTestValues->end());
861 template<
class ScalarType,
class MV,
class OP>
864 std::ostringstream oss;
ScaleType convertStringToScaleType(const std::string &scaleType)
Convert the given string to its ScaleType enum value.
Collection of types and exceptions used within the Belos solvers.
This class implements the preconditioned transpose-free QMR algorithm for solving non-Hermitian linea...
int getNumIters() const override
Get the iteration count for the most recent call to solve().
Belos's basic output manager for sending information of select verbosity levels to the appropriate ou...
Class which manages the output and verbosity of the Belos solvers.
void setParameters(const Teuchos::RCP< Teuchos::ParameterList > ¶ms) override
Set the parameters the solver manager should use to solve the linear problem.
bool is_null(const boost::shared_ptr< T > &p)
Teuchos::Array< Teuchos::RCP< Teuchos::Time > > getTimers() const
Return the timers for this object.
bool isLOADetected() const override
Whether loss of accuracy was detected during the last solve() invocation.
T & get(const std::string &name, T def_value)
ParameterList & set(std::string const &name, T const &value, std::string const &docString="", RCP< const ParameterEntryValidator > const &validator=null)
bool is_null(const std::shared_ptr< T > &p)
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
const LinearProblem< ScalarType, MV, OP > & getProblem() const override
Return a reference to the linear problem being solved by this solver manager.
Teuchos::RCP< SolverManager< ScalarType, MV, OP > > clone() const override
clone for Inverted Injection (DII)
An implementation of StatusTestResNorm using a family of residual norms.
Belos concrete class for generating iterations with the preconditioned tranpose-free QMR (TFQMR) meth...
Teuchos::RCP< const MV > R
The current residual basis.
static const double convTol
Default convergence tolerance.
Belos::StatusTest class for specifying a maximum number of iterations.
std::string description() const override
Method to return description of the TFQMR solver manager.
static std::string name()
A factory class for generating StatusTestOutput objects.
ReturnType solve() override
This method performs possibly repeated calls to the underlying linear solver's iterate() routine unti...
Traits class which defines basic operations on multivectors.
Belos::StatusTest for logically combining several status tests.
bool isParameter(const std::string &name) const
MagnitudeType achievedTol() const override
Tolerance achieved by the last solve() invocation.
A Belos::StatusTest class for specifying a maximum number of iterations.
ResetType
How to reset the solver.
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
Pure virtual base class which describes the basic interface for a solver manager. ...
static const double impTolScale
"Implicit Tolerance Scale Factor"
static void summarize(Ptr< const Comm< int > > comm, std::ostream &out=std::cout, const bool alwaysWriteLocal=false, const bool writeGlobalStats=true, const bool writeZeroTimers=true, const ECounterSetOp setOp=Intersection, const std::string &filter="", const bool ignoreZeroTimers=false)
A linear system to solve, and its associated information.
Structure to contain pointers to TFQMRIter state variables.
Class which describes the linear problem to be solved by the iterative solver.
TFQMRSolMgr()
Empty constructor for TFQMRSolMgr. This constructor takes no arguments and sets the default values fo...
virtual ~TFQMRSolMgr()
Destructor.
ReturnType
Whether the Belos solve converged for all linear systems.
void reset(const ResetType type) override
Performs a reset of the solver manager specified by the ResetType. This informs the solver manager th...
The Belos::SolverManager is a templated virtual base class that defines the basic interface that any ...
The Belos::TFQMRSolMgr provides a powerful and fully-featured solver manager over the TFQMR linear so...
void validateParameters(ParameterList const &validParamList, int const depth=1000, EValidateUsed const validateUsed=VALIDATE_USED_ENABLED, EValidateDefaults const validateDefaults=VALIDATE_DEFAULTS_ENABLED) const
Teuchos::RCP< const Teuchos::ParameterList > getValidParameters() const override
Get a parameter list containing the valid parameters for this object.
Belos::StatusTestResNorm for specifying general residual norm stopping criteria.
bool isType(const std::string &name) const
A class for extending the status testing capabilities of Belos via logical combinations.
TFQMRSolMgrLinearProblemFailure(const std::string &what_arg)
Class which defines basic traits for the operator type.
Parent class to all Belos exceptions.
Default parameters common to most Belos solvers.
TFQMRSolMgrLinearProblemFailure is thrown when the linear problem is not setup (i.e.
Belos header file which uses auto-configuration information to include necessary C++ headers...
void setProblem(const Teuchos::RCP< LinearProblem< ScalarType, MV, OP > > &problem) override
Set the linear problem that needs to be solved.
Teuchos::RCP< const Teuchos::ParameterList > getCurrentParameters() const override
Get a parameter list containing the current parameters for this object.