42#ifndef BELOS_PSEUDO_BLOCK_TFQMR_SOLMGR_HPP
43#define BELOS_PSEUDO_BLOCK_TFQMR_SOLMGR_HPP
61#ifdef BELOS_TEUCHOS_TIME_MONITOR
62#include "Teuchos_TimeMonitor.hpp"
93 template<
class ScalarType,
class MV,
class OP>
99 typedef Teuchos::ScalarTraits<ScalarType> SCT;
100 typedef typename Teuchos::ScalarTraits<ScalarType>::magnitudeType MagnitudeType;
101 typedef Teuchos::ScalarTraits<MagnitudeType> MT;
132 const Teuchos::RCP<Teuchos::ParameterList> &pl );
138 Teuchos::RCP<SolverManager<ScalarType, MV, OP> >
clone ()
const override {
163 Teuchos::Array<Teuchos::RCP<Teuchos::Time> >
getTimers()
const {
164 return Teuchos::tuple(timerSolve_);
198 void setParameters(
const Teuchos::RCP<Teuchos::ParameterList> ¶ms )
override;
244 bool checkStatusTest();
247 Teuchos::RCP<LinearProblem<ScalarType,MV,OP> > problem_;
250 Teuchos::RCP<OutputManager<ScalarType> > printer_;
251 Teuchos::RCP<std::ostream> outputStream_;
254 Teuchos::RCP<StatusTest<ScalarType,MV,OP> > sTest_;
255 Teuchos::RCP<StatusTestMaxIters<ScalarType,MV,OP> > maxIterTest_;
256 Teuchos::RCP<StatusTest<ScalarType,MV,OP> > convTest_;
257 Teuchos::RCP<StatusTestGenResNorm<ScalarType,MV,OP> > expConvTest_, impConvTest_;
258 Teuchos::RCP<StatusTestOutput<ScalarType,MV,OP> > outputTest_;
261 Teuchos::RCP<Teuchos::ParameterList> params_;
264 static constexpr int maxIters_default_ = 1000;
265 static constexpr bool expResTest_default_ =
false;
268 static constexpr int outputFreq_default_ = -1;
269 static constexpr int defQuorum_default_ = 1;
270 static constexpr const char * impResScale_default_ =
"Norm of Preconditioned Initial Residual";
271 static constexpr const char * expResScale_default_ =
"Norm of Initial Residual";
272 static constexpr const char * label_default_ =
"Belos";
275 MagnitudeType convtol_, impTolScale_, achievedTol_;
276 int maxIters_, numIters_;
277 int verbosity_, outputStyle_, outputFreq_, defQuorum_;
279 std::string impResScale_, expResScale_;
283 Teuchos::RCP<Teuchos::Time> timerSolve_;
286 bool isSet_, isSTSet_;
291template<
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_),
302 defQuorum_(defQuorum_default_),
303 expResTest_(expResTest_default_),
304 impResScale_(impResScale_default_),
305 expResScale_(expResScale_default_),
306 label_(label_default_),
313template<
class ScalarType,
class MV,
class OP>
316 const Teuchos::RCP<Teuchos::ParameterList> &pl ) :
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_),
327 defQuorum_(defQuorum_default_),
328 expResTest_(expResTest_default_),
329 impResScale_(impResScale_default_),
330 expResScale_(expResScale_default_),
331 label_(label_default_),
335 TEUCHOS_TEST_FOR_EXCEPTION(problem_ == Teuchos::null, std::invalid_argument,
"Problem not given to solver manager.");
338 if ( !is_null(pl) ) {
343template<
class ScalarType,
class MV,
class OP>
347 if (params_ == Teuchos::null) {
348 params_ = Teuchos::rcp(
new Teuchos::ParameterList(*getValidParameters()) );
351 params->validateParameters(*getValidParameters());
355 if (params->isParameter(
"Maximum Iterations")) {
356 maxIters_ = params->get(
"Maximum Iterations",maxIters_default_);
359 params_->set(
"Maximum Iterations", maxIters_);
360 if (maxIterTest_!=Teuchos::null)
361 maxIterTest_->setMaxIters( maxIters_ );
365 if (params->isParameter(
"Timer Label")) {
366 std::string tempLabel = params->get(
"Timer Label", label_default_);
369 if (tempLabel != label_) {
371 params_->set(
"Timer Label", label_);
372 std::string solveLabel = label_ +
": PseudoBlockTFQMRSolMgr total solve time";
373#ifdef BELOS_TEUCHOS_TIME_MONITOR
374 timerSolve_ = Teuchos::TimeMonitor::getNewCounter(solveLabel);
380 if (params->isParameter(
"Verbosity")) {
381 if (Teuchos::isParameterType<int>(*params,
"Verbosity")) {
382 verbosity_ = params->get(
"Verbosity", verbosity_default_);
384 verbosity_ = (int)Teuchos::getParameter<Belos::MsgType>(*params,
"Verbosity");
388 params_->set(
"Verbosity", verbosity_);
389 if (printer_ != Teuchos::null)
390 printer_->setVerbosity(verbosity_);
394 if (params->isParameter(
"Output Style")) {
395 if (Teuchos::isParameterType<int>(*params,
"Output Style")) {
396 outputStyle_ = params->get(
"Output Style", outputStyle_default_);
398 outputStyle_ = (int)Teuchos::getParameter<Belos::OutputType>(*params,
"Output Style");
402 params_->set(
"Output Style", outputStyle_);
407 if (params->isParameter(
"Output Stream")) {
408 outputStream_ = Teuchos::getParameter<Teuchos::RCP<std::ostream> >(*params,
"Output Stream");
411 params_->set(
"Output Stream", outputStream_);
412 if (printer_ != Teuchos::null)
413 printer_->setOStream( outputStream_ );
418 if (params->isParameter(
"Output Frequency")) {
419 outputFreq_ = params->get(
"Output Frequency", outputFreq_default_);
423 params_->set(
"Output Frequency", outputFreq_);
424 if (outputTest_ != Teuchos::null)
425 outputTest_->setOutputFrequency( outputFreq_ );
429 if (printer_ == Teuchos::null) {
434 if (params->isParameter(
"Convergence Tolerance")) {
435 if (params->isType<MagnitudeType> (
"Convergence Tolerance")) {
436 convtol_ = params->get (
"Convergence Tolerance",
444 params_->set(
"Convergence Tolerance", convtol_);
448 if (params->isParameter(
"Implicit Tolerance Scale Factor")) {
449 if (params->isType<MagnitudeType> (
"Implicit Tolerance Scale Factor")) {
450 impTolScale_ = params->get (
"Implicit Tolerance Scale Factor",
455 impTolScale_ = params->get (
"Implicit Tolerance Scale Factor",
460 params_->set(
"Implicit Tolerance Scale Factor", impTolScale_);
464 if (params->isParameter(
"Implicit Residual Scaling")) {
465 std::string tempImpResScale = Teuchos::getParameter<std::string>( *params,
"Implicit Residual Scaling" );
468 if (impResScale_ != tempImpResScale) {
469 impResScale_ = tempImpResScale;
472 params_->set(
"Implicit Residual Scaling", impResScale_);
477 if (params->isParameter(
"Explicit Residual Scaling")) {
478 std::string tempExpResScale = Teuchos::getParameter<std::string>( *params,
"Explicit Residual Scaling" );
481 if (expResScale_ != tempExpResScale) {
482 expResScale_ = tempExpResScale;
485 params_->set(
"Explicit Residual Scaling", expResScale_);
490 if (params->isParameter(
"Explicit Residual Test")) {
491 expResTest_ = Teuchos::getParameter<bool>( *params,
"Explicit Residual Test" );
494 params_->set(
"Explicit Residual Test", expResTest_);
495 if (expConvTest_ == Teuchos::null) {
501 if (params->isParameter(
"Deflation Quorum")) {
502 defQuorum_ = params->get(
"Deflation Quorum", defQuorum_);
503 params_->set (
"Deflation Quorum", defQuorum_);
504 if (! impConvTest_.is_null ()) {
505 impConvTest_->setQuorum (defQuorum_);
507 if (! expConvTest_.is_null ()) {
508 expConvTest_->setQuorum (defQuorum_);
513 if (timerSolve_ == Teuchos::null) {
514 std::string solveLabel = label_ +
": PseudoBlockTFQMRSolMgr total solve time";
515#ifdef BELOS_TEUCHOS_TIME_MONITOR
516 timerSolve_ = Teuchos::TimeMonitor::getNewCounter(solveLabel);
526template<
class ScalarType,
class MV,
class OP>
538 Teuchos::RCP<StatusTestGenResNorm_t> tmpImpConvTest =
539 Teuchos::rcp(
new StatusTestGenResNorm_t( impTolScale_*convtol_, defQuorum_ ) );
541 impConvTest_ = tmpImpConvTest;
544 Teuchos::RCP<StatusTestGenResNorm_t> tmpExpConvTest =
545 Teuchos::rcp(
new StatusTestGenResNorm_t( convtol_, defQuorum_ ) );
546 tmpExpConvTest->defineResForm( StatusTestGenResNorm_t::Explicit,
Belos::TwoNorm );
548 expConvTest_ = tmpExpConvTest;
551 convTest_ = Teuchos::rcp(
new StatusTestCombo_t( StatusTestCombo_t::SEQ, impConvTest_, expConvTest_ ) );
556 Teuchos::RCP<StatusTestGenResNorm_t> tmpImpConvTest =
557 Teuchos::rcp(
new StatusTestGenResNorm_t( convtol_, defQuorum_ ) );
559 impConvTest_ = tmpImpConvTest;
562 expConvTest_ = impConvTest_;
563 convTest_ = impConvTest_;
565 sTest_ = Teuchos::rcp(
new StatusTestCombo_t( StatusTestCombo_t::OR, maxIterTest_, convTest_ ) );
569 StatusTestOutputFactory<ScalarType,MV,OP> stoFactory( outputStyle_ );
573 std::string solverDesc =
" Pseudo Block TFQMR ";
574 outputTest_->setSolverDesc( solverDesc );
584template<
class ScalarType,
class MV,
class OP>
585Teuchos::RCP<const Teuchos::ParameterList>
588 static Teuchos::RCP<const Teuchos::ParameterList> validPL;
591 if(is_null(validPL)) {
592 Teuchos::RCP<Teuchos::ParameterList> pl = Teuchos::parameterList();
597 "The relative residual tolerance that needs to be achieved by the\n"
598 "iterative solver in order for the linear system to be declared converged.");
600 "The scale factor used by the implicit residual test when explicit residual\n"
601 "testing is used. May enable faster convergence when TFQMR bound is too loose.");
602 pl->set(
"Maximum Iterations",
static_cast<int>(maxIters_default_),
603 "The maximum number of block iterations allowed for each\n"
604 "set of RHS solved.");
605 pl->set(
"Verbosity",
static_cast<int>(verbosity_default_),
606 "What type(s) of solver information should be outputted\n"
607 "to the output stream.");
608 pl->set(
"Output Style",
static_cast<int>(outputStyle_default_),
609 "What style is used for the solver information outputted\n"
610 "to the output stream.");
611 pl->set(
"Output Frequency",
static_cast<int>(outputFreq_default_),
612 "How often convergence information should be outputted\n"
613 "to the output stream.");
614 pl->set(
"Deflation Quorum",
static_cast<int>(defQuorum_default_),
615 "The number of linear systems that need to converge before they are deflated.");
616 pl->set(
"Output Stream", Teuchos::rcpFromRef(std::cout),
617 "A reference-counted pointer to the output stream where all\n"
618 "solver output is sent.");
619 pl->set(
"Explicit Residual Test",
static_cast<bool>(expResTest_default_),
620 "Whether the explicitly computed residual should be used in the convergence test.");
621 pl->set(
"Implicit Residual Scaling",
static_cast<const char *
>(impResScale_default_),
622 "The type of scaling used in the implicit residual convergence test.");
623 pl->set(
"Explicit Residual Scaling",
static_cast<const char *
>(expResScale_default_),
624 "The type of scaling used in the explicit residual convergence test.");
625 pl->set(
"Timer Label",
static_cast<const char *
>(label_default_),
626 "The string to use as a prefix for the timer labels.");
634template<
class ScalarType,
class MV,
class OP>
641 setParameters(Teuchos::parameterList(*getValidParameters()));
645 "Belos::PseudoBlockTFQMRSolMgr::solve(): Linear problem is not a valid object.");
648 "Belos::PseudoBlockTFQMRSolMgr::solve(): Linear problem is not ready, setProblem() has not been called.");
652 "Belos::PseudoBlockTFQMRSolMgr::solve(): Linear problem and requested status tests are incompatible.");
657 int numRHS2Solve = MVT::GetNumberVecs( *(problem_->getRHS()) );
658 int numCurrRHS = numRHS2Solve;
660 std::vector<int> currIdx( numRHS2Solve );
661 for (
int i=0; i<numRHS2Solve; ++i) {
662 currIdx[i] = startPtr+i;
666 problem_->setLSIndex( currIdx );
670 Teuchos::ParameterList plist;
673 outputTest_->reset();
676 bool isConverged =
true;
681 Teuchos::RCP<PseudoBlockTFQMRIter<ScalarType,MV,OP> > block_tfqmr_iter =
686#ifdef BELOS_TEUCHOS_TIME_MONITOR
687 Teuchos::TimeMonitor slvtimer(*timerSolve_);
690 while ( numRHS2Solve > 0 ) {
693 std::vector<int> convRHSIdx;
694 std::vector<int> currRHSIdx( currIdx );
695 currRHSIdx.resize(numCurrRHS);
698 block_tfqmr_iter->resetNumIters();
701 outputTest_->resetNumCalls();
704 Teuchos::RCP<MV> R_0 = MVT::CloneViewNonConst( *(Teuchos::rcp_const_cast<MV>(problem_->getInitPrecResVec())), currIdx );
709 block_tfqmr_iter->initializeTFQMR(newstate);
715 block_tfqmr_iter->iterate();
722 if ( convTest_->getStatus() ==
Passed ) {
725 std::vector<int> convIdx = expConvTest_->convIndices();
729 if (convIdx.size() == currRHSIdx.size())
733 problem_->updateSolution( block_tfqmr_iter->getCurrentUpdate(),
true );
736 problem_->setCurrLS();
740 std::vector<int> unconvIdx( currRHSIdx.size() );
741 for (
unsigned int i=0; i<currRHSIdx.size(); ++i) {
743 for (
unsigned int j=0; j<convIdx.size(); ++j) {
744 if (currRHSIdx[i] == convIdx[j]) {
751 currRHSIdx[have++] = currRHSIdx[i];
754 unconvIdx.resize(have);
755 currRHSIdx.resize(have);
758 problem_->setLSIndex( currRHSIdx );
768 defstate.
Rtilde = MVT::CloneView( *currentState.
Rtilde, unconvIdx);
769 defstate.
U = MVT::CloneView( *currentState.
U, unconvIdx );
770 defstate.
AU = MVT::CloneView( *currentState.
AU, unconvIdx );
771 defstate.
V = MVT::CloneView( *currentState.
V, unconvIdx );
772 defstate.
W = MVT::CloneView( *currentState.
W, unconvIdx );
773 defstate.
D = MVT::CloneView( *currentState.
D, unconvIdx );
776 for (std::vector<int>::iterator uIter = unconvIdx.begin(); uIter != unconvIdx.end(); uIter++)
778 defstate.
alpha.push_back( currentState.
alpha[ *uIter ] );
779 defstate.
eta.push_back( currentState.
eta[ *uIter ] );
780 defstate.
rho.push_back( currentState.
rho[ *uIter ] );
781 defstate.
tau.push_back( currentState.
tau[ *uIter ] );
782 defstate.
theta.push_back( currentState.
theta[ *uIter ] );
785 block_tfqmr_iter->initializeTFQMR(defstate);
792 else if ( maxIterTest_->getStatus() ==
Passed ) {
806 TEUCHOS_TEST_FOR_EXCEPTION(
true,std::logic_error,
807 "Belos::PseudoBlockTFQMRSolMgr::solve(): Invalid return from PseudoBlockTFQMRIter::iterate().");
810 catch (
const std::exception &e) {
811 printer_->stream(
Errors) <<
"Error! Caught std::exception in PseudoBlockTFQMRIter::iterate() at iteration "
812 << block_tfqmr_iter->getNumIters() << std::endl
813 << e.what() << std::endl;
819 problem_->updateSolution( block_tfqmr_iter->getCurrentUpdate(),
true );
822 problem_->setCurrLS();
825 startPtr += numCurrRHS;
826 numRHS2Solve -= numCurrRHS;
827 if ( numRHS2Solve > 0 ) {
828 numCurrRHS = numRHS2Solve;
829 currIdx.resize( numCurrRHS );
830 for (
int i=0; i<numCurrRHS; ++i)
831 { currIdx[i] = startPtr+i; }
834 if (defQuorum_ > numCurrRHS) {
835 if (impConvTest_ != Teuchos::null)
836 impConvTest_->setQuorum( numCurrRHS );
837 if (expConvTest_ != Teuchos::null)
838 expConvTest_->setQuorum( numCurrRHS );
842 problem_->setLSIndex( currIdx );
845 currIdx.resize( numRHS2Solve );
856#ifdef BELOS_TEUCHOS_TIME_MONITOR
861 Teuchos::TimeMonitor::summarize( printer_->stream(
TimingDetails) );
865 numIters_ = maxIterTest_->getNumIters();
878 const std::vector<MagnitudeType>* pTestValues = NULL;
880 pTestValues = expConvTest_->getTestValue();
881 if (pTestValues == NULL || pTestValues->size() < 1) {
882 pTestValues = impConvTest_->getTestValue();
887 pTestValues = impConvTest_->getTestValue();
889 TEUCHOS_TEST_FOR_EXCEPTION(pTestValues == NULL, std::logic_error,
890 "Belos::PseudoBlockTFQMRSolMgr::solve(): The implicit convergence test's "
891 "getTestValue() method returned NULL. Please report this bug to the "
892 "Belos developers.");
893 TEUCHOS_TEST_FOR_EXCEPTION(pTestValues->size() < 1, std::logic_error,
894 "Belos::TMQMRSolMgr::solve(): The implicit convergence test's "
895 "getTestValue() method returned a vector of length zero. Please report "
896 "this bug to the Belos developers.");
901 achievedTol_ = *std::max_element (pTestValues->begin(), pTestValues->end());
911template<
class ScalarType,
class MV,
class OP>
914 std::ostringstream oss;
915 oss <<
"Belos::PseudoBlockTFQMRSolMgr<...,"<<Teuchos::ScalarTraits<ScalarType>::name()<<
">";
Belos header file which uses auto-configuration information to include necessary C++ headers.
Class which describes the linear problem to be solved by the iterative solver.
Class which manages the output and verbosity of the Belos solvers.
Belos concrete class for generating iterations with the preconditioned tranpose-free QMR (TFQMR) meth...
Pure virtual base class which describes the basic interface for a solver manager.
Belos::StatusTest for logically combining several status tests.
Belos::StatusTestResNorm for specifying general residual norm stopping criteria.
Belos::StatusTest class for specifying a maximum number of iterations.
A factory class for generating StatusTestOutput objects.
Collection of types and exceptions used within the Belos solvers.
Parent class to all Belos exceptions.
A linear system to solve, and its associated information.
Traits class which defines basic operations on multivectors.
Class which defines basic traits for the operator type.
Belos's basic output manager for sending information of select verbosity levels to the appropriate ou...
This class implements the preconditioned transpose-free QMR algorithm for solving non-Hermitian linea...
PseudoBlockTFQMRSolMgrLinearProblemFailure is thrown when the linear problem is not setup (i....
PseudoBlockTFQMRSolMgrLinearProblemFailure(const std::string &what_arg)
The Belos::PseudoBlockTFQMRSolMgr provides a powerful and fully-featured solver manager over the pseu...
Teuchos::RCP< const Teuchos::ParameterList > getCurrentParameters() const override
Get a parameter list containing the current parameters for this object.
MagnitudeType achievedTol() const override
Tolerance achieved by the last solve() invocation.
bool isLOADetected() const override
Whether loss of accuracy was detected during the last solve() invocation.
void reset(const ResetType type) override
Performs a reset of the solver manager specified by the ResetType. This informs the solver manager th...
Teuchos::RCP< const Teuchos::ParameterList > getValidParameters() const override
Get a parameter list containing the valid parameters for this object.
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)
Teuchos::Array< Teuchos::RCP< Teuchos::Time > > getTimers() const
Return the timers for this object.
std::string description() const override
Method to return description of the pseudo-block TFQMR solver manager.
PseudoBlockTFQMRSolMgr()
Empty constructor for PseudoBlockTFQMRSolMgr. This constructor takes no arguments and sets the defaul...
void setParameters(const Teuchos::RCP< Teuchos::ParameterList > ¶ms) override
Set the parameters the solver manager should use to solve the linear problem.
int getNumIters() const override
Get the iteration count for the most recent call to solve().
virtual ~PseudoBlockTFQMRSolMgr()
Destructor.
void setProblem(const Teuchos::RCP< LinearProblem< ScalarType, MV, OP > > &problem) override
Set the linear problem that needs to be solved.
ReturnType solve() override
This method performs possibly repeated calls to the underlying linear solver's iterate() routine unti...
The Belos::SolverManager is a templated virtual base class that defines the basic interface that any ...
A class for extending the status testing capabilities of Belos via logical combinations.
An implementation of StatusTestResNorm using a family of residual norms.
A Belos::StatusTest class for specifying a maximum number of iterations.
ScaleType convertStringToScaleType(const std::string &scaleType)
Convert the given string to its ScaleType enum value.
ReturnType
Whether the Belos solve converged for all linear systems.
ResetType
How to reset the solver.
Default parameters common to most Belos solvers.
static const double impTolScale
"Implicit Tolerance Scale Factor"
static const double convTol
Default convergence tolerance.
Structure to contain pointers to PseudoBlockTFQMRIter state variables.
std::vector< MagnitudeType > theta
std::vector< MagnitudeType > tau
Teuchos::RCP< const MV > AU
std::vector< ScalarType > rho
Teuchos::RCP< const MV > D
std::vector< ScalarType > alpha
Teuchos::RCP< const MV > Rtilde
Teuchos::RCP< const MV > W
The current residual basis.
Teuchos::RCP< const MV > U
Teuchos::RCP< const MV > V
std::vector< ScalarType > eta