42 #ifndef BELOS_BLOCK_CG_SOLMGR_HPP
43 #define BELOS_BLOCK_CG_SOLMGR_HPP
65 #ifdef BELOS_TEUCHOS_TIME_MONITOR
68 #if defined(HAVE_TEUCHOSCORE_CXX11)
69 # include <type_traits>
70 #endif // defined(HAVE_TEUCHOSCORE_CXX11)
104 template<
class ScalarType,
class MV,
class OP,
105 const bool lapackSupportsScalarType =
130 template<
class ScalarType,
class MV,
class OP>
242 return Teuchos::tuple(timerSolve_);
313 std::string description()
const override;
351 static constexpr
int maxIters_default_ = 1000;
352 static constexpr
bool adaptiveBlockSize_default_ =
true;
353 static constexpr
bool showMaxResNormOnly_default_ =
false;
354 static constexpr
bool useSingleReduction_default_ =
false;
355 static constexpr
int blockSize_default_ = 1;
358 static constexpr
int outputFreq_default_ = -1;
359 static constexpr
const char * resNorm_default_ =
"TwoNorm";
360 static constexpr
bool foldConvergenceDetectionIntoAllreduce_default_ =
false;
361 static constexpr
const char * resScale_default_ =
"Norm of Initial Residual";
362 static constexpr
const char * label_default_ =
"Belos";
363 static constexpr
const char * orthoType_default_ =
"ICGS";
364 static constexpr
bool assertPositiveDefiniteness_default_ =
true;
408 template<
class ScalarType,
class MV,
class OP>
410 outputStream_(Teuchos::rcpFromRef(std::cout)),
414 maxIters_(maxIters_default_),
416 blockSize_(blockSize_default_),
417 verbosity_(verbosity_default_),
418 outputStyle_(outputStyle_default_),
419 outputFreq_(outputFreq_default_),
420 adaptiveBlockSize_(adaptiveBlockSize_default_),
421 showMaxResNormOnly_(showMaxResNormOnly_default_),
422 useSingleReduction_(useSingleReduction_default_),
423 orthoType_(orthoType_default_),
424 resScale_(resScale_default_),
425 assertPositiveDefiniteness_(assertPositiveDefiniteness_default_),
426 foldConvergenceDetectionIntoAllreduce_(foldConvergenceDetectionIntoAllreduce_default_),
427 label_(label_default_),
433 template<
class ScalarType,
class MV,
class OP>
438 outputStream_(Teuchos::rcpFromRef(std::cout)),
442 maxIters_(maxIters_default_),
444 blockSize_(blockSize_default_),
445 verbosity_(verbosity_default_),
446 outputStyle_(outputStyle_default_),
447 outputFreq_(outputFreq_default_),
448 adaptiveBlockSize_(adaptiveBlockSize_default_),
449 showMaxResNormOnly_(showMaxResNormOnly_default_),
450 useSingleReduction_(useSingleReduction_default_),
451 orthoType_(orthoType_default_),
452 resScale_(resScale_default_),
453 assertPositiveDefiniteness_(assertPositiveDefiniteness_default_),
454 foldConvergenceDetectionIntoAllreduce_(foldConvergenceDetectionIntoAllreduce_default_),
455 label_(label_default_),
459 "BlockCGSolMgr's constructor requires a nonnull LinearProblem instance.");
469 template<
class ScalarType,
class MV,
class OP>
484 maxIters_ = params->
get(
"Maximum Iterations",maxIters_default_);
487 params_->set(
"Maximum Iterations", maxIters_);
489 maxIterTest_->setMaxIters( maxIters_ );
494 blockSize_ = params->
get(
"Block Size",blockSize_default_);
496 "Belos::BlockCGSolMgr: \"Block Size\" must be strictly positive.");
499 params_->set(
"Block Size", blockSize_);
504 adaptiveBlockSize_ = params->
get(
"Adaptive Block Size",adaptiveBlockSize_default_);
507 params_->set(
"Adaptive Block Size", adaptiveBlockSize_);
512 useSingleReduction_ = params->
get(
"Use Single Reduction", useSingleReduction_default_);
515 if (params->
isParameter(
"Fold Convergence Detection Into Allreduce")) {
516 foldConvergenceDetectionIntoAllreduce_ = params->
get(
"Fold Convergence Detection Into Allreduce",
517 foldConvergenceDetectionIntoAllreduce_default_);
522 std::string tempLabel = params->
get(
"Timer Label", label_default_);
525 if (tempLabel != label_) {
527 params_->set(
"Timer Label", label_);
528 std::string solveLabel = label_ +
": BlockCGSolMgr total solve time";
529 #ifdef BELOS_TEUCHOS_TIME_MONITOR
533 ortho_->setLabel( label_ );
540 if (Teuchos::isParameterType<int>(*params,
"Verbosity")) {
541 verbosity_ = params->
get(
"Verbosity", verbosity_default_);
543 verbosity_ = (int)Teuchos::getParameter<Belos::MsgType>(*params,
"Verbosity");
547 params_->set(
"Verbosity", verbosity_);
549 printer_->setVerbosity(verbosity_);
554 if (Teuchos::isParameterType<int>(*params,
"Output Style")) {
555 outputStyle_ = params->
get(
"Output Style", outputStyle_default_);
557 outputStyle_ = (int)Teuchos::getParameter<Belos::OutputType>(*params,
"Output Style");
561 params_->set(
"Output Style", outputStyle_);
567 outputStream_ = Teuchos::getParameter<Teuchos::RCP<std::ostream> >(*params,
"Output Stream");
570 params_->set(
"Output Stream", outputStream_);
572 printer_->setOStream( outputStream_ );
578 outputFreq_ = params->
get(
"Output Frequency", outputFreq_default_);
582 params_->set(
"Output Frequency", outputFreq_);
584 outputTest_->setOutputFrequency( outputFreq_ );
593 bool changedOrthoType =
false;
595 std::string tempOrthoType = params->
get(
"Orthogonalization",orthoType_default_);
596 if (tempOrthoType != orthoType_) {
597 orthoType_ = tempOrthoType;
598 changedOrthoType =
true;
601 params_->set(
"Orthogonalization", orthoType_);
604 if (params->
isParameter(
"Orthogonalization Constant")) {
606 orthoKappa_ = params->
get (
"Orthogonalization Constant",
610 orthoKappa_ = params->
get (
"Orthogonalization Constant",
615 params_->set(
"Orthogonalization Constant",orthoKappa_);
616 if (orthoType_==
"DGKS") {
617 if (orthoKappa_ > 0 && ortho_ !=
Teuchos::null && !changedOrthoType) {
627 if (orthoType_==
"DGKS" && orthoKappa_ > 0) {
629 paramsOrtho->
set (
"depTol", orthoKappa_ );
640 if (params->
isParameter(
"Convergence Tolerance")) {
642 convtol_ = params->
get (
"Convergence Tolerance",
650 params_->set(
"Convergence Tolerance", convtol_);
652 convTest_->setTolerance( convtol_ );
655 if (params->
isParameter(
"Show Maximum Residual Norm Only")) {
656 showMaxResNormOnly_ = Teuchos::getParameter<bool>(*params,
"Show Maximum Residual Norm Only");
659 params_->set(
"Show Maximum Residual Norm Only", showMaxResNormOnly_);
661 convTest_->setShowMaxResNormOnly( showMaxResNormOnly_ );
665 bool newResTest =
false;
667 std::string tempResScale = resScale_;
668 if (params->
isParameter (
"Implicit Residual Scaling")) {
669 tempResScale = params->
get<std::string> (
"Implicit Residual Scaling");
673 if (resScale_ != tempResScale) {
676 resScale_ = tempResScale;
679 params_->set (
"Implicit Residual Scaling", resScale_);
681 if (! convTest_.is_null ()) {
685 if (params->
isType<std::string> (
"Residual Norm")) {
689 convTest_->defineResForm(StatusTestResNorm_t::Implicit, normType);
692 catch (std::exception& e) {
707 if (convTest_.is_null () || newResTest) {
711 if (params->
isType<std::string> (
"Residual Norm")) {
716 convTest_ =
rcp (
new StatusTestResNorm_t (convtol_, 1, showMaxResNormOnly_));
717 convTest_->defineResForm(StatusTestResNorm_t::Implicit, normType);
721 if (sTest_.is_null () || newResTest)
722 sTest_ =
Teuchos::rcp(
new StatusTestCombo_t( StatusTestCombo_t::OR, maxIterTest_, convTest_ ) );
724 if (outputTest_.is_null () || newResTest) {
732 std::string solverDesc =
" Block CG ";
733 outputTest_->setSolverDesc( solverDesc );
738 if (params->
isParameter(
"Assert Positive Definiteness")) {
739 assertPositiveDefiniteness_ = Teuchos::getParameter<bool>(*params,
"Assert Positive Definiteness");
740 params_->set(
"Assert Positive Definiteness", assertPositiveDefiniteness_);
745 std::string solveLabel = label_ +
": BlockCGSolMgr total solve time";
746 #ifdef BELOS_TEUCHOS_TIME_MONITOR
756 template<
class ScalarType,
class MV,
class OP>
766 "The relative residual tolerance that needs to be achieved by the\n"
767 "iterative solver in order for the linear system to be declared converged.");
768 pl->
set(
"Maximum Iterations", static_cast<int>(maxIters_default_),
769 "The maximum number of block iterations allowed for each\n"
770 "set of RHS solved.");
771 pl->
set(
"Block Size", static_cast<int>(blockSize_default_),
772 "The number of vectors in each block.");
773 pl->
set(
"Adaptive Block Size", static_cast<bool>(adaptiveBlockSize_default_),
774 "Whether the solver manager should adapt to the block size\n"
775 "based on the number of RHS to solve.");
776 pl->
set(
"Verbosity", static_cast<int>(verbosity_default_),
777 "What type(s) of solver information should be outputted\n"
778 "to the output stream.");
779 pl->
set(
"Output Style", static_cast<int>(outputStyle_default_),
780 "What style is used for the solver information outputted\n"
781 "to the output stream.");
782 pl->
set(
"Output Frequency", static_cast<int>(outputFreq_default_),
783 "How often convergence information should be outputted\n"
784 "to the output stream.");
785 pl->
set(
"Output Stream", Teuchos::rcpFromRef(std::cout),
786 "A reference-counted pointer to the output stream where all\n"
787 "solver output is sent.");
788 pl->
set(
"Show Maximum Residual Norm Only", static_cast<bool>(showMaxResNormOnly_default_),
789 "When convergence information is printed, only show the maximum\n"
790 "relative residual norm when the block size is greater than one.");
791 pl->
set(
"Use Single Reduction", static_cast<bool>(useSingleReduction_default_),
792 "Use single reduction iteration when the block size is one.");
793 pl->
set(
"Implicit Residual Scaling", resScale_default_,
794 "The type of scaling used in the residual convergence test.");
795 pl->
set(
"Timer Label", static_cast<const char *>(label_default_),
796 "The string to use as a prefix for the timer labels.");
797 pl->
set(
"Orthogonalization", static_cast<const char *>(orthoType_default_),
798 "The type of orthogonalization to use: DGKS, ICGS, or IMGS.");
799 pl->
set(
"Assert Positive Definiteness",static_cast<bool>(assertPositiveDefiniteness_default_),
800 "Assert for positivity of p^H*A*p in CG iteration.");
802 "The constant used by DGKS orthogonalization to determine\n"
803 "whether another step of classical Gram-Schmidt is necessary.");
804 pl->
set(
"Residual Norm",static_cast<const char *>(resNorm_default_),
805 "Norm used for the convergence check on the residual.");
806 pl->
set(
"Fold Convergence Detection Into Allreduce",static_cast<bool>(foldConvergenceDetectionIntoAllreduce_default_),
807 "Merge the allreduce for convergence detection with the one for CG.\n"
808 "This saves one all-reduce, but incurs more computation.");
816 template<
class ScalarType,
class MV,
class OP>
820 using Teuchos::rcp_const_cast;
821 using Teuchos::rcp_dynamic_cast;
828 setParameters(Teuchos::parameterList(*getValidParameters()));
835 "Belos::BlockCGSolMgr::solve(): Linear problem is not ready, setProblem() "
836 "has not been called.");
840 int numRHS2Solve = MVT::GetNumberVecs( *(problem_->getRHS()) );
841 int numCurrRHS = ( numRHS2Solve < blockSize_) ? numRHS2Solve : blockSize_;
843 std::vector<int> currIdx, currIdx2;
848 if ( adaptiveBlockSize_ ) {
849 blockSize_ = numCurrRHS;
850 currIdx.resize( numCurrRHS );
851 currIdx2.resize( numCurrRHS );
852 for (
int i=0; i<numCurrRHS; ++i)
853 { currIdx[i] = startPtr+i; currIdx2[i]=i; }
857 currIdx.resize( blockSize_ );
858 currIdx2.resize( blockSize_ );
859 for (
int i=0; i<numCurrRHS; ++i)
860 { currIdx[i] = startPtr+i; currIdx2[i]=i; }
861 for (
int i=numCurrRHS; i<blockSize_; ++i)
862 { currIdx[i] = -1; currIdx2[i] = i; }
866 problem_->setLSIndex( currIdx );
871 plist.
set(
"Block Size",blockSize_);
874 outputTest_->reset();
878 bool isConverged =
true;
883 plist.
set(
"Assert Positive Definiteness", assertPositiveDefiniteness_);
885 RCP<CGIteration<ScalarType,MV,OP> > block_cg_iter;
886 if (blockSize_ == 1) {
890 plist.
set(
"Fold Convergence Detection Into Allreduce",
891 foldConvergenceDetectionIntoAllreduce_);
892 if (useSingleReduction_) {
895 outputTest_, convTest_, plist));
900 outputTest_, convTest_, plist));
911 #ifdef BELOS_TEUCHOS_TIME_MONITOR
915 while ( numRHS2Solve > 0 ) {
918 std::vector<int> convRHSIdx;
919 std::vector<int> currRHSIdx( currIdx );
920 currRHSIdx.resize(numCurrRHS);
923 block_cg_iter->resetNumIters();
926 outputTest_->resetNumCalls();
929 RCP<MV> R_0 = MVT::CloneViewNonConst( *(rcp_const_cast<MV>(problem_->getInitResVec())), currIdx );
934 block_cg_iter->initializeCG(newstate);
940 block_cg_iter->iterate();
944 if (convTest_->getStatus() ==
Passed) {
949 std::vector<int> convIdx =
950 rcp_dynamic_cast<conv_test_type>(convTest_)->convIndices();
955 if (convIdx.size() == currRHSIdx.size())
960 problem_->setCurrLS();
965 std::vector<int> unconvIdx(currRHSIdx.size());
966 for (
unsigned int i=0; i<currRHSIdx.size(); ++i) {
968 for (
unsigned int j=0; j<convIdx.size(); ++j) {
969 if (currRHSIdx[i] == convIdx[j]) {
975 currIdx2[have] = currIdx2[i];
976 currRHSIdx[have++] = currRHSIdx[i];
981 currRHSIdx.resize(have);
982 currIdx2.resize(have);
985 problem_->setLSIndex( currRHSIdx );
988 std::vector<MagnitudeType> norms;
989 R_0 = MVT::CloneCopy( *(block_cg_iter->getNativeResiduals(&norms)),currIdx2 );
990 for (
int i=0; i<have; ++i) { currIdx2[i] = i; }
993 block_cg_iter->setBlockSize( have );
998 block_cg_iter->initializeCG(defstate);
1004 else if (maxIterTest_->getStatus() ==
Passed) {
1005 isConverged =
false;
1014 "Belos::BlockCGSolMgr::solve(): Neither the convergence test nor "
1015 "the maximum iteration count test passed. Please report this bug "
1016 "to the Belos developers.");
1021 achievedTol_ = MT::one();
1023 MVT::MvInit( *X, SCT::zero() );
1024 printer_->stream(
Warnings) <<
"Belos::BlockCGSolMgr::solve(): Warning! NaN has been detected!"
1028 catch (
const std::exception &e) {
1029 std::ostream& err = printer_->stream (
Errors);
1030 err <<
"Error! Caught std::exception in CGIteration::iterate() at "
1031 <<
"iteration " << block_cg_iter->getNumIters() << std::endl
1032 << e.what() << std::endl;
1039 problem_->setCurrLS();
1042 startPtr += numCurrRHS;
1043 numRHS2Solve -= numCurrRHS;
1044 if ( numRHS2Solve > 0 ) {
1045 numCurrRHS = ( numRHS2Solve < blockSize_) ? numRHS2Solve : blockSize_;
1047 if ( adaptiveBlockSize_ ) {
1048 blockSize_ = numCurrRHS;
1049 currIdx.resize( numCurrRHS );
1050 currIdx2.resize( numCurrRHS );
1051 for (
int i=0; i<numCurrRHS; ++i)
1052 { currIdx[i] = startPtr+i; currIdx2[i] = i; }
1055 currIdx.resize( blockSize_ );
1056 currIdx2.resize( blockSize_ );
1057 for (
int i=0; i<numCurrRHS; ++i)
1058 { currIdx[i] = startPtr+i; currIdx2[i] = i; }
1059 for (
int i=numCurrRHS; i<blockSize_; ++i)
1060 { currIdx[i] = -1; currIdx2[i] = i; }
1063 problem_->setLSIndex( currIdx );
1066 block_cg_iter->setBlockSize( blockSize_ );
1069 currIdx.resize( numRHS2Solve );
1080 #ifdef BELOS_TEUCHOS_TIME_MONITOR
1091 numIters_ = maxIterTest_->getNumIters();
1097 const std::vector<MagnitudeType>* pTestValues =
1098 rcp_dynamic_cast<conv_test_type>(convTest_)->getTestValue();
1101 "Belos::BlockCGSolMgr::solve(): The convergence test's getTestValue() "
1102 "method returned NULL. Please report this bug to the Belos developers.");
1105 "Belos::BlockCGSolMgr::solve(): The convergence test's getTestValue() "
1106 "method returned a vector of length zero. Please report this bug to the "
1107 "Belos developers.");
1112 achievedTol_ = *std::max_element (pTestValues->begin(), pTestValues->end());
1122 template<
class ScalarType,
class MV,
class OP>
1125 std::ostringstream oss;
1128 oss <<
"Ortho Type='"<<orthoType_<<
"\', Block Size=" << blockSize_;
static const double orthoKappa
DGKS orthogonalization constant.
ScaleType convertStringToScaleType(const std::string &scaleType)
Convert the given string to its ScaleType enum value.
Teuchos::RCP< const MV > R
The current residual.
bool isLOADetected() const override
Return whether a loss of accuracy was detected by this solver during the most current solve...
Collection of types and exceptions used within the Belos solvers.
Teuchos::Array< Teuchos::RCP< Teuchos::Time > > getTimers() const
Return the timers for this object.
Belos's basic output manager for sending information of select verbosity levels to the appropriate ou...
int maxIters_
Maximum iteration count (read from parameter list).
Class which manages the output and verbosity of the Belos solvers.
bool is_null(const boost::shared_ptr< T > &p)
bool isSet_
Whether or not the parameters have been set (via setParameters()).
Teuchos::RCP< Teuchos::ParameterList > params_
Current parameter list.
Teuchos::RCP< StatusTestOutput< ScalarType, MV, OP > > outputTest_
Output "status test" that controls all the other status tests.
ScaleType
The type of scaling to use on the residual norm value.
Teuchos::RCP< SolverManager< ScalarType, MV, OP > > clone() const override
clone for Inverted Injection (DII)
Teuchos::ScalarTraits< ScalarType > SCT
This class implements the preconditioned Conjugate Gradient (CG) iteration.
T & get(ParameterList &l, const std::string &name)
Structure to contain pointers to CGIteration state variables.
Belos concrete class for performing the conjugate-gradient (CG) iteration.
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)
A factory class for generating StatusTestOutput objects.
Belos concrete class for performing the block conjugate-gradient (CG) iteration.
An implementation of StatusTestResNorm using a family of residual norms.
void setProblem(const Teuchos::RCP< LinearProblem< ScalarType, MV, OP > > &problem) override
Set the linear problem that needs to be solved.
std::string label_
Prefix label for all the timers.
The Belos::BlockCGSolMgr provides a powerful and fully-featured solver manager over the CG and BlockC...
static const double convTol
Default convergence tolerance.
Belos::StatusTest class for specifying a maximum number of iterations.
static std::string name()
BlockCGSolMgrLinearProblemFailure(const std::string &what_arg)
Teuchos::RCP< const Teuchos::ParameterList > getCurrentParameters() const override
Get a parameter list containing the current parameters for this object.
BlockCGSolMgrLinearProblemFailure is thrown when the linear problem is not setup (i.e.
A factory class for generating StatusTestOutput objects.
Teuchos::RCP< MatOrthoManager< ScalarType, MV, OP > > ortho_
Orthogonalization manager.
BlockCGSolMgr(const Teuchos::RCP< LinearProblem< ScalarType, MV, OP > > &problem, const Teuchos::RCP< Teuchos::ParameterList > &pl)
Teuchos::RCP< StatusTestMaxIters< ScalarType, MV, OP > > maxIterTest_
Maximum iteration count stopping criterion.
Traits class which defines basic operations on multivectors.
Belos::StatusTest for logically combining several status tests.
Teuchos::RCP< StatusTest< ScalarType, MV, OP > > sTest_
Aggregate stopping criterion.
bool assertPositiveDefiniteness_
Details::SolverManagerRequiresLapack< ScalarType, MV, OP, requiresLapack > base_type
bool isParameter(const std::string &name) const
virtual ~BlockCGSolMgr()
Destructor.
OperatorTraits< ScalarType, MV, OP > OPT
MultiVecTraits< ScalarType, MV > MVT
A Belos::StatusTest class for specifying a maximum number of iterations.
ResetType
How to reset the solver.
void reset(const ResetType type) override
Performs a reset of the solver manager specified by the ResetType. This informs the solver manager th...
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. ...
Belos concrete class for performing a single-reduction conjugate-gradient (CG) iteration.
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)
bool is_null(const RCP< T > &p)
NormType convertStringToNormType(const std::string &normType)
Convert the given string to its NormType enum value.
A linear system to solve, and its associated information.
Class which describes the linear problem to be solved by the iterative solver.
MagnitudeType achievedTol_
Tolerance achieved by the last solve() invocation.
MagnitudeType achievedTol() const override
Tolerance achieved by the last solve() invocation.
Type traits class that says whether Teuchos::LAPACK has a valid implementation for the given ScalarTy...
Teuchos::ScalarTraits< ScalarType >::magnitudeType MagnitudeType
ReturnType
Whether the Belos solve converged for all linear systems.
Teuchos::RCP< StatusTestGenResNorm< ScalarType, MV, OP > > convTest_
Convergence stopping criterion.
int getNumIters() const override
Get the iteration count for the most recent call to solve().
Teuchos::RCP< LinearProblem< ScalarType, MV, OP > > problem_
The linear problem to solve.
void validateParameters(ParameterList const &validParamList, int const depth=1000, EValidateUsed const validateUsed=VALIDATE_USED_ENABLED, EValidateDefaults const validateDefaults=VALIDATE_DEFAULTS_ENABLED) const
NormType
The type of vector norm to compute.
Teuchos::RCP< StatusTestOutput< ScalarType, MV, OP > > create(const Teuchos::RCP< OutputManager< ScalarType > > &printer, Teuchos::RCP< StatusTest< ScalarType, MV, OP > > test, int mod, int printStates)
Create the StatusTestOutput object specified by the outputStyle.
Belos::StatusTestResNorm for specifying general residual norm stopping criteria.
MagnitudeType convtol_
Convergence tolerance (read from parameter list).
bool isType(const std::string &name) const
Teuchos::RCP< Teuchos::Time > timerSolve_
Solve timer.
An implementation of the Belos::MatOrthoManager that performs orthogonalization using (potentially) m...
A class for extending the status testing capabilities of Belos via logical combinations.
const LinearProblem< ScalarType, MV, OP > & getProblem() const override
Return a reference to the linear problem being solved by this solver manager.
Class which defines basic traits for the operator type.
bool foldConvergenceDetectionIntoAllreduce_
Parent class to all Belos exceptions.
Default parameters common to most Belos solvers.
This class implements the preconditioned single-reduction Conjugate Gradient (CG) iteration...
Teuchos::RCP< OutputManager< ScalarType > > printer_
Output manager, that handles printing of different kinds of messages.
Base class for Belos::SolverManager subclasses which normally can only compile with ScalarType types ...
Belos header file which uses auto-configuration information to include necessary C++ headers...
int numIters_
Number of iterations taken by the last solve() invocation.
Teuchos::ScalarTraits< MagnitudeType > MT
static const bool requiresLapack
Teuchos::RCP< std::ostream > outputStream_
Output stream to which the output manager prints.
MagnitudeType orthoKappa_
Orthogonalization parameter (read from parameter list).
Stub implementation of BlockCGIter, for ScalarType types for which Teuchos::LAPACK does NOT have a va...
Teuchos::RCP< Belos::MatOrthoManager< Scalar, MV, OP > > makeMatOrthoManager(const std::string &ortho, const Teuchos::RCP< const OP > &M, const Teuchos::RCP< OutputManager< Scalar > > &, const std::string &label, const Teuchos::RCP< Teuchos::ParameterList > ¶ms)
Return an instance of the specified MatOrthoManager subclass.