45 #ifndef THYRA_BELOS_LINEAR_OP_WITH_SOLVE_FACTORY_HPP
46 #define THYRA_BELOS_LINEAR_OP_WITH_SOLVE_FACTORY_HPP
49 #include "Thyra_BelosLinearOpWithSolveFactory_decl.hpp"
50 #include "Thyra_BelosLinearOpWithSolve.hpp"
51 #include "Thyra_ScaledAdjointLinearOpBase.hpp"
66 #if defined(HAVE_BELOS_TPETRA) && defined(HAVE_STRATIMIKOS_THYRATPETRAADAPTERS)
67 #include "Thyra_BelosTpetrasSolverAdapter.hpp"
70 #include "Teuchos_VerboseObjectParameterListHelpers.hpp"
71 #include "Teuchos_StandardParameterEntryValidators.hpp"
72 #include "Teuchos_ParameterList.hpp"
73 #include "Teuchos_dyn_cast.hpp"
74 #include "Teuchos_ValidatorXMLConverterDB.hpp"
75 #include "Teuchos_StandardValidatorXMLConverters.hpp"
83 template<
class Scalar>
84 const std::string BelosLinearOpWithSolveFactory<Scalar>::SolverType_name =
"Solver Type";
85 template<
class Scalar>
86 const std::string BelosLinearOpWithSolveFactory<Scalar>::SolverType_default =
"Pseudo Block GMRES";
87 template<
class Scalar>
88 const std::string BelosLinearOpWithSolveFactory<Scalar>::SolverTypes_name =
"Solver Types";
89 template<
class Scalar>
90 const std::string BelosLinearOpWithSolveFactory<Scalar>::BlockGMRES_name =
"Block GMRES";
91 template<
class Scalar>
92 const std::string BelosLinearOpWithSolveFactory<Scalar>::PseudoBlockGMRES_name =
"Pseudo Block GMRES";
93 template<
class Scalar>
94 const std::string BelosLinearOpWithSolveFactory<Scalar>::BlockCG_name =
"Block CG";
95 template<
class Scalar>
96 const std::string BelosLinearOpWithSolveFactory<Scalar>::PseudoBlockCG_name =
"Pseudo Block CG";
97 template<
class Scalar>
98 const std::string BelosLinearOpWithSolveFactory<Scalar>::PseudoBlockStochasticCG_name =
"Pseudo Block Stochastic CG";
99 template<
class Scalar>
100 const std::string BelosLinearOpWithSolveFactory<Scalar>::GCRODR_name =
"GCRODR";
101 template<
class Scalar>
102 const std::string BelosLinearOpWithSolveFactory<Scalar>::RCG_name =
"RCG";
103 template<
class Scalar>
104 const std::string BelosLinearOpWithSolveFactory<Scalar>::MINRES_name =
"MINRES";
105 template<
class Scalar>
106 const std::string BelosLinearOpWithSolveFactory<Scalar>::TFQMR_name =
"TFQMR";
107 template<
class Scalar>
108 const std::string BelosLinearOpWithSolveFactory<Scalar>::BiCGStab_name =
"BiCGStab";
109 template<
class Scalar>
110 const std::string BelosLinearOpWithSolveFactory<Scalar>::FixedPoint_name =
"Fixed Point";
112 #if defined(HAVE_BELOS_TPETRA) && defined(HAVE_STRATIMIKOS_THYRATPETRAADAPTERS)
113 template<
class Scalar>
114 const std::string BelosLinearOpWithSolveFactory<Scalar>::TpetraGmres_name =
"TPETRA GMRES";
115 template<
class Scalar>
116 const std::string BelosLinearOpWithSolveFactory<Scalar>::TpetraGmresPipeline_name =
"TPETRA GMRES PIPELINE";
117 template<
class Scalar>
118 const std::string BelosLinearOpWithSolveFactory<Scalar>::TpetraGmresSingleReduce_name =
"TPETRA GMRES SINGLE REDUCE";
119 template<
class Scalar>
120 const std::string BelosLinearOpWithSolveFactory<Scalar>::TpetraGmresSstep_name =
"TPETRA GMRES S-STEP";
123 template<
class Scalar>
124 const std::string BelosLinearOpWithSolveFactory<Scalar>::ConvergenceTestFrequency_name =
"Convergence Test Frequency";
127 const std::string LeftPreconditionerIfUnspecified_name =
"Left Preconditioner If Unspecified";
133 template<
class Scalar>
135 :solverType_(SOLVER_TYPE_PSEUDO_BLOCK_GMRES),
136 convergenceTestFrequency_(1)
138 updateThisValidParamList();
142 template<
class Scalar>
146 :solverType_(SOLVER_TYPE_PSEUDO_BLOCK_GMRES)
155 template<
class Scalar>
162 template<
class Scalar>
165 const std::string &precFactoryName
170 precFactoryValidPL = precFactory->getValidParameters();
171 const std::string _precFactoryName =
172 ( precFactoryName !=
""
174 : ( precFactoryValidPL.get() ? precFactoryValidPL->name() :
"GENERIC PRECONDITIONER FACTORY" )
176 precFactory_ = precFactory;
177 precFactoryName_ = _precFactoryName;
178 updateThisValidParamList();
182 template<
class Scalar>
190 template<
class Scalar>
193 std::string *precFactoryName
196 if(precFactory) *precFactory = precFactory_;
197 if(precFactoryName) *precFactoryName = precFactoryName_;
198 precFactory_ = Teuchos::null;
199 precFactoryName_ =
"";
200 updateThisValidParamList();
204 template<
class Scalar>
209 if(precFactory_.get())
210 return precFactory_->isCompatible(fwdOpSrc);
215 template<
class Scalar>
223 template<
class Scalar>
231 initializeOpImpl(fwdOpSrc,null,null,
false,Op,supportSolveUse);
235 template<
class Scalar>
246 template<
class Scalar>
251 if(precFactory_.get())
257 template<
class Scalar>
266 initializeOpImpl(fwdOpSrc,null,prec,
false,Op,supportSolveUse);
270 template<
class Scalar>
279 initializeOpImpl(fwdOpSrc,approxFwdOpSrc,null,
false,Op,supportSolveUse);
283 template<
class Scalar>
308 if(fwdOpSrc) *fwdOpSrc = _fwdOpSrc;
309 if(prec) *prec = _prec;
310 if(approxFwdOpSrc) *approxFwdOpSrc = _approxFwdOpSrc;
311 if(supportSolveUse) *supportSolveUse = _supportSolveUse;
318 template<
class Scalar>
325 paramList_ = paramList;
327 Teuchos::getIntegralValue<EBelosSolverType>(*paramList_, SolverType_name);
328 convergenceTestFrequency_ =
329 Teuchos::getParameter<int>(*paramList_, ConvergenceTestFrequency_name);
330 Teuchos::readVerboseObjectSublist(&*paramList_,
this);
334 template<
class Scalar>
342 template<
class Scalar>
347 paramList_ = Teuchos::null;
352 template<
class Scalar>
360 template<
class Scalar>
364 return thisValidParamList_;
371 template<
class Scalar>
374 std::ostringstream oss;
378 oss << precFactory_->description();
388 template<
class Scalar>
393 using Teuchos::tuple;
394 using Teuchos::setStringToIntegralParameter;
400 EBelosSolverType> >::getDummyObject());
405 if(validParamList.
get()==NULL) {
407 setStringToIntegralParameter<EBelosSolverType>(
408 SolverType_name, SolverType_default,
409 "Type of linear solver algorithm to use.",
412 "Pseudo Block GMRES",
415 "Pseudo Block Stochastic CG",
422 #if defined(HAVE_BELOS_TPETRA) && defined(HAVE_STRATIMIKOS_THYRATPETRAADAPTERS)
424 "TPETRA GMRES PIPELINE",
425 "TPETRA GMRES SINGLE REDUCE",
426 "TPETRA GMRES S-STEP"
430 "Block GMRES solver for nonsymmetric linear systems. It can also solve "
431 "single right-hand side systems, and can also perform Flexible GMRES "
432 "(where the preconditioner may change at every iteration, for example "
433 "for inner-outer iterations), by setting options in the \"Block GMRES\" "
436 "GMRES solver for nonsymmetric linear systems, that performs single "
437 "right-hand side solves on multiple right-hand sides at once. It "
438 "exploits operator multivector multiplication in order to amortize "
439 "global communication costs. Individual linear systems are deflated "
440 "out as they are solved.",
442 "Block CG solver for symmetric (Hermitian in complex arithmetic) "
443 "positive definite linear systems. It can also solve single "
444 "right-hand-side systems.",
446 "CG solver that performs single right-hand side CG on multiple right-hand "
447 "sides at once. It exploits operator multivector multiplication in order "
448 "to amortize global communication costs. Individual linear systems are "
449 "deflated out as they are solved.",
451 "stochastic CG solver that performs single right-hand side CG on multiple right-hand "
452 "sides at once. It exploits operator multivector multiplication in order "
453 "to amortize global communication costs. Individual linear systems are "
454 "deflated out as they are solved. [EXPERIMENTAL]",
456 "Variant of GMRES that performs subspace recycling to accelerate "
457 "convergence for sequences of solves with related linear systems. "
458 "Individual linear systems are deflated out as they are solved. "
459 "The current implementation only supports real-valued Scalar types.",
461 "CG solver for symmetric (Hermitian in complex arithmetic) positive "
462 "definite linear systems, that performs subspace recycling to "
463 "accelerate convergence for sequences of related linear systems.",
465 "MINRES solver for symmetric indefinite linear systems. It performs "
466 "single-right-hand-side solves on multiple right-hand sides sequentially.",
468 "TFQMR (Transpose-Free QMR) solver for nonsymmetric linear systems.",
470 "BiCGStab solver for nonsymmetric linear systems.",
472 "Fixed point iteration"
474 #
if defined(HAVE_BELOS_TPETRA) && defined(HAVE_STRATIMIKOS_THYRATPETRAADAPTERS)
475 ,
"Native Tpetra implementation of GMRES",
477 "Native Tpetra implementation of pipeline GMRES",
479 "Native Tpetra implementation of single-reduce GMRES",
481 "Native Tpetra implementation of s-step GMRES"
484 tuple<EBelosSolverType>(
485 SOLVER_TYPE_BLOCK_GMRES,
486 SOLVER_TYPE_PSEUDO_BLOCK_GMRES,
487 SOLVER_TYPE_BLOCK_CG,
488 SOLVER_TYPE_PSEUDO_BLOCK_CG,
489 SOLVER_TYPE_PSEUDO_BLOCK_STOCHASTIC_CG,
494 SOLVER_TYPE_BICGSTAB,
495 SOLVER_TYPE_FIXEDPOINT
496 #
if defined(HAVE_BELOS_TPETRA) && defined(HAVE_STRATIMIKOS_THYRATPETRAADAPTERS)
497 ,SOLVER_TYPE_TPETRA_GMRES,
498 SOLVER_TYPE_TPETRA_GMRES_PIPELINE,
499 SOLVER_TYPE_TPETRA_GMRES_SINGLE_REDUCE,
500 SOLVER_TYPE_TPETRA_GMRES_SSTEP
505 validParamList->
set(ConvergenceTestFrequency_name, as<int>(1),
506 "Number of linear solver iterations to skip between applying"
507 " user-defined convergence test.");
509 LeftPreconditionerIfUnspecified_name,
false,
510 "If the preconditioner does not specify if it is left or right, and this\n"
511 "option is set to true, put the preconditioner on the left side.\n"
512 "Historically, preconditioning is on the right. Some solvers may not\n"
513 "support left preconditioning.");
515 &solverTypesSL = validParamList->
sublist(SolverTypes_name);
532 if (lapackSupportsScalar) {
535 *mgr.getValidParameters()
538 if (lapackSupportsScalar) {
541 *mgr.getValidParameters()
550 if (lapackSupportsScalar) {
553 *mgr.getValidParameters()
556 if (lapackSupportsScalar && scalarIsReal) {
559 *mgr.getValidParameters()
586 #if defined(HAVE_BELOS_TPETRA) && defined(HAVE_STRATIMIKOS_THYRATPETRAADAPTERS)
588 Thyra::BelosTpetraGmres<Scalar,MV_t,LO_t> mgr;
590 *mgr.getValidParameters()
594 Thyra::BelosTpetraGmresPipeline<Scalar,MV_t,LO_t> mgr;
596 *mgr.getValidParameters()
600 Thyra::BelosTpetraGmresSingleReduce<Scalar,MV_t,LO_t> mgr;
602 *mgr.getValidParameters()
606 Thyra::BelosTpetraGmresSstep<Scalar,MV_t,LO_t> mgr;
608 *mgr.getValidParameters()
613 return validParamList;
617 template<
class Scalar>
618 void BelosLinearOpWithSolveFactory<Scalar>::updateThisValidParamList()
623 Teuchos::setupVerboseObjectSublist(&*thisValidParamList_);
627 template<
class Scalar>
628 void BelosLinearOpWithSolveFactory<Scalar>::initializeOpImpl(
629 const RCP<
const LinearOpSourceBase<Scalar> > &fwdOpSrc,
630 const RCP<
const LinearOpSourceBase<Scalar> > &approxFwdOpSrc,
631 const RCP<
const PreconditionerBase<Scalar> > &prec_in,
632 const bool reusePrec,
633 LinearOpWithSolveBase<Scalar> *Op,
639 using Teuchos::set_extra_data;
641 typedef MultiVectorBase<Scalar> MV_t;
642 typedef LinearOpBase<Scalar> LO_t;
644 const RCP<Teuchos::FancyOStream> out = this->getOStream();
647 if(out.get() &&
static_cast<int>(verbLevel) > static_cast<int>(
Teuchos::VERB_LOW))
648 *out <<
"\nEntering Thyra::BelosLinearOpWithSolveFactory<"<<ST::name()<<
">::initializeOpImpl(...) ...\n";
658 RCP<const LinearOpBase<Scalar> >
659 fwdOp = fwdOpSrc->getOp(),
660 approxFwdOp = ( approxFwdOpSrc.get() ? approxFwdOpSrc->getOp() : Teuchos::null );
666 BelosLinearOpWithSolve<Scalar>
673 RCP<PreconditionerBase<Scalar> > myPrec = Teuchos::null;
674 RCP<const PreconditionerBase<Scalar> > prec = Teuchos::null;
681 if(precFactory_.get()) {
683 ( !belosOp->isExternalPrec()
684 ? Teuchos::rcp_const_cast<PreconditionerBase<Scalar> >(belosOp->extract_prec())
687 bool hasExistingPrec =
false;
689 hasExistingPrec =
true;
694 hasExistingPrec =
false;
695 myPrec = precFactory_->createPrec();
697 if( hasExistingPrec && reusePrec ) {
702 if(approxFwdOp.get())
703 precFactory_->initializePrec(approxFwdOpSrc,&*myPrec);
705 precFactory_->initializePrec(fwdOpSrc,&*myPrec);
715 bool oldIsExternalPrec =
false;
716 RCP<Belos::LinearProblem<Scalar,MV_t,LO_t> > oldLP = Teuchos::null;
717 RCP<Belos::SolverManager<Scalar,MV_t,LO_t> > oldIterSolver = Teuchos::null;
718 RCP<const LinearOpSourceBase<Scalar> > oldFwdOpSrc = Teuchos::null;
719 RCP<const LinearOpSourceBase<Scalar> > oldApproxFwdOpSrc = Teuchos::null;
722 belosOp->uninitialize( &oldLP, NULL, &oldIterSolver, &oldFwdOpSrc,
723 NULL, &oldIsExternalPrec, &oldApproxFwdOpSrc, &oldSupportSolveUse );
732 if (oldLP != Teuchos::null) {
736 lp =
rcp(
new LP_t());
743 lp->setOperator(fwdOp);
750 RCP<const LinearOpBase<Scalar> > unspecified = prec->getUnspecifiedPrecOp();
751 RCP<const LinearOpBase<Scalar> > left = prec->getLeftPrecOp();
752 RCP<const LinearOpBase<Scalar> > right = prec->getRightPrecOp();
754 !( left.get() || right.get() || unspecified.get() ), std::logic_error
755 ,
"Error, at least one preconditoner linear operator objects must be set!"
758 if (paramList_->get<
bool>(LeftPreconditionerIfUnspecified_name,
false))
759 lp->setLeftPrec(unspecified);
761 lp->setRightPrec(unspecified);
764 lp->setLeftPrec(left);
767 lp->setRightPrec(right);
773 ,
"Error, we can not currently handle split preconditioners!"
778 set_extra_data<RCP<PreconditionerBase<Scalar> > >(myPrec,
"Belos::InternalPrec",
779 Teuchos::inOutArg(lp), Teuchos::POST_DESTROY,
false);
781 else if(prec.get()) {
782 set_extra_data<RCP<const PreconditionerBase<Scalar> > >(prec,
"Belos::ExternalPrec",
783 Teuchos::inOutArg(lp), Teuchos::POST_DESTROY,
false);
791 RCP<IterativeSolver_t> iterativeSolver = Teuchos::null;
794 switch(solverType_) {
795 case SOLVER_TYPE_BLOCK_GMRES:
798 if(paramList_.get()) {
804 if (oldIterSolver != Teuchos::null) {
805 iterativeSolver = oldIterSolver;
806 iterativeSolver->setProblem( lp );
807 iterativeSolver->setParameters( solverPL );
814 case SOLVER_TYPE_PSEUDO_BLOCK_GMRES:
817 if(paramList_.get()) {
825 if (oldIterSolver != Teuchos::null) {
826 iterativeSolver = oldIterSolver;
827 iterativeSolver->setProblem( lp );
828 iterativeSolver->setParameters( solverPL );
835 case SOLVER_TYPE_BLOCK_CG:
838 if(paramList_.get()) {
844 if (oldIterSolver != Teuchos::null) {
845 iterativeSolver = oldIterSolver;
846 iterativeSolver->setProblem( lp );
847 iterativeSolver->setParameters( solverPL );
854 case SOLVER_TYPE_PSEUDO_BLOCK_CG:
857 if(paramList_.get()) {
865 if (oldIterSolver != Teuchos::null) {
866 iterativeSolver = oldIterSolver;
867 iterativeSolver->setProblem( lp );
868 iterativeSolver->setParameters( solverPL );
875 case SOLVER_TYPE_PSEUDO_BLOCK_STOCHASTIC_CG:
878 if(paramList_.get()) {
886 if (oldIterSolver != Teuchos::null) {
887 iterativeSolver = oldIterSolver;
888 iterativeSolver->setProblem( lp );
889 iterativeSolver->setParameters( solverPL );
896 case SOLVER_TYPE_GCRODR:
899 if(paramList_.get()) {
905 if (oldIterSolver != Teuchos::null) {
906 iterativeSolver = oldIterSolver;
907 iterativeSolver->setProblem( lp );
908 iterativeSolver->setParameters( solverPL );
915 case SOLVER_TYPE_RCG:
918 if(paramList_.get()) {
924 if (oldIterSolver != Teuchos::null) {
925 iterativeSolver = oldIterSolver;
926 iterativeSolver->setProblem( lp );
927 iterativeSolver->setParameters( solverPL );
934 case SOLVER_TYPE_MINRES:
937 if(paramList_.get()) {
943 if (oldIterSolver != Teuchos::null) {
944 iterativeSolver = oldIterSolver;
945 iterativeSolver->setProblem( lp );
946 iterativeSolver->setParameters( solverPL );
953 case SOLVER_TYPE_TFQMR:
956 if(paramList_.get()) {
962 if (oldIterSolver != Teuchos::null) {
963 iterativeSolver = oldIterSolver;
964 iterativeSolver->setProblem( lp );
965 iterativeSolver->setParameters( solverPL );
972 case SOLVER_TYPE_BICGSTAB:
975 if(paramList_.get()) {
981 if (oldIterSolver != Teuchos::null) {
982 iterativeSolver = oldIterSolver;
983 iterativeSolver->setProblem( lp );
984 iterativeSolver->setParameters( solverPL );
991 case SOLVER_TYPE_FIXEDPOINT:
994 if(paramList_.get()) {
1000 if (oldIterSolver != Teuchos::null) {
1001 iterativeSolver = oldIterSolver;
1002 iterativeSolver->setProblem( lp );
1003 iterativeSolver->setParameters( solverPL );
1010 #if defined(HAVE_BELOS_TPETRA) && defined(HAVE_STRATIMIKOS_THYRATPETRAADAPTERS)
1011 case SOLVER_TYPE_TPETRA_GMRES:
1014 if(paramList_.get()) {
1020 iterativeSolver =
rcp(
new Thyra::BelosTpetraGmres<Scalar,MV_t,LO_t>());
1022 iterativeSolver->setProblem( lp );
1023 iterativeSolver->setParameters( solverPL );
1026 case SOLVER_TYPE_TPETRA_GMRES_PIPELINE:
1029 if(paramList_.get()) {
1032 solverPL =
Teuchos::rcp( &tpetraGmresPipelinePL,
false );
1035 iterativeSolver =
rcp(
new Thyra::BelosTpetraGmresPipeline<Scalar,MV_t,LO_t>());
1037 iterativeSolver->setProblem( lp );
1038 iterativeSolver->setParameters( solverPL );
1041 case SOLVER_TYPE_TPETRA_GMRES_SINGLE_REDUCE:
1044 if(paramList_.get()) {
1047 solverPL =
Teuchos::rcp( &tpetraGmresSingleReducePL,
false );
1050 iterativeSolver =
rcp(
new Thyra::BelosTpetraGmresSingleReduce<Scalar,MV_t,LO_t>());
1052 iterativeSolver->setProblem( lp );
1053 iterativeSolver->setParameters( solverPL );
1056 case SOLVER_TYPE_TPETRA_GMRES_SSTEP:
1059 if(paramList_.get()) {
1065 iterativeSolver =
rcp(
new Thyra::BelosTpetraGmresSstep<Scalar,MV_t,LO_t>());
1067 iterativeSolver->setProblem( lp );
1068 iterativeSolver->setParameters( solverPL );
1082 belosOp->initialize(
1083 lp, solverPL, iterativeSolver,
1084 fwdOpSrc, prec, myPrec.get()==NULL, approxFwdOpSrc,
1085 supportSolveUse, convergenceTestFrequency_
1087 belosOp->setOStream(out);
1088 belosOp->setVerbLevel(verbLevel);
1089 #ifdef TEUCHOS_DEBUG
1090 if(paramList_.get()) {
1092 paramList_->validateParameters(*this->getValidParameters(),1);
1095 if(out.get() &&
static_cast<int>(verbLevel) > static_cast<int>(
Teuchos::VERB_LOW))
1096 *out <<
"\nLeaving Thyra::BelosLinearOpWithSolveFactory<"<<ST::name()<<
">::initializeOpImpl(...) ...\n";
1104 #endif // THYRA_BELOS_LINEAR_OP_WITH_SOLVE_FACTORY_HPP
Teuchos::RCP< const Teuchos::ParameterList > getValidParameters() const
void uninitializeOp(LinearOpWithSolveBase< Scalar > *Op, Teuchos::RCP< const LinearOpSourceBase< Scalar > > *fwdOpSrc, Teuchos::RCP< const PreconditionerBase< Scalar > > *prec, Teuchos::RCP< const LinearOpSourceBase< Scalar > > *approxFwdOpSrc, ESupportSolveUse *supportSolveUse) const
RCP< const LinearOpSourceBase< Scalar > > extract_approxFwdOpSrc()
static void addConverter(RCP< const ParameterEntryValidator > validator, RCP< ValidatorXMLConverter > converterToAdd)
Concrete LinearOpWithSolveBase subclass in terms of Belos.
ParameterList & set(std::string const &name, T const &value, std::string const &docString="", RCP< const ParameterEntryValidator > const &validator=null)
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
Teuchos::RCP< const Teuchos::ParameterList > getValidParameters() const override
RCP< const PreconditionerBase< Scalar > > extract_prec()
T_To & dyn_cast(T_From &from)
void initializeAndReuseOp(const Teuchos::RCP< const LinearOpSourceBase< Scalar > > &fwdOpSrc, LinearOpWithSolveBase< Scalar > *Op) const
bool isExternalPrec() const
void unsetPreconditionerFactory(Teuchos::RCP< PreconditionerFactoryBase< Scalar > > *precFactory, std::string *precFactoryName)
Thyra specializations of MultiVecTraits and OperatorTraits.
void setParameterList(Teuchos::RCP< Teuchos::ParameterList > const ¶mList)
bool isCompatible(const LinearOpSourceBase< Scalar > &fwdOpSrc) const
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
bool acceptsPreconditionerFactory() const
Returns true .
Teuchos::RCP< const Teuchos::ParameterList > getValidParameters() const override
Teuchos::RCP< const Teuchos::ParameterList > getValidParameters() const override
Teuchos::RCP< Teuchos::ParameterList > getNonconstParameterList()
Teuchos::RCP< const Teuchos::ParameterList > getValidParameters() const override
virtual std::string description() const
void initializeOp(const Teuchos::RCP< const LinearOpSourceBase< Scalar > > &fwdOpSrc, LinearOpWithSolveBase< Scalar > *Op, const ESupportSolveUse supportSolveUse) const
std::string description() const
void validateParametersAndSetDefaults(ParameterList const &validParamList, int const depth=1000)
void setPreconditionerFactory(const Teuchos::RCP< PreconditionerFactoryBase< Scalar > > &precFactory, const std::string &precFactoryName)
Teuchos::RCP< const Teuchos::ParameterList > getValidParameters() const override
LinearOpWithSolveFactoryBase subclass implemented in terms of Belos.
Teuchos::RCP< const Teuchos::ParameterList > getValidParameters() const override
ESupportSolveUse supportSolveUse() const
ParameterList & setParameters(const ParameterList &source)
bool supportsPreconditionerInputType(const EPreconditionerInputType precOpType) const
Teuchos::RCP< const Teuchos::ParameterList > getParameterList() const
RCP< const LinearOpSourceBase< Scalar > > extract_fwdOpSrc()
bool nonnull(const boost::shared_ptr< T > &p)
Teuchos::RCP< PreconditionerFactoryBase< Scalar > > getPreconditionerFactory() const
Teuchos::RCP< const Teuchos::ParameterList > getValidParameters() const override
void initializeApproxPreconditionedOp(const Teuchos::RCP< const LinearOpSourceBase< Scalar > > &fwdOpSrc, const Teuchos::RCP< const LinearOpSourceBase< Scalar > > &approxFwdOpSrc, LinearOpWithSolveBase< Scalar > *Op, const ESupportSolveUse supportSolveUse) const
TypeTo as(const TypeFrom &t)
void initializePreconditionedOp(const Teuchos::RCP< const LinearOpSourceBase< Scalar > > &fwdOpSrc, const Teuchos::RCP< const PreconditionerBase< Scalar > > &prec, LinearOpWithSolveBase< Scalar > *Op, const ESupportSolveUse supportSolveUse) const
Teuchos::RCP< LinearOpWithSolveBase< Scalar > > createOp() const
ParameterList & sublist(const std::string &name, bool mustAlreadyExist=false, const std::string &docString="")
BelosLinearOpWithSolveFactory()
Construct without preconditioner factory.
Teuchos::RCP< Teuchos::ParameterList > unsetParameterList()
#define TEUCHOS_TEST_FOR_EXCEPT(throw_exception_test)