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" 53 #include "BelosBlockGmresSolMgr.hpp" 54 #include "BelosPseudoBlockGmresSolMgr.hpp" 55 #include "BelosBlockCGSolMgr.hpp" 56 #include "BelosPseudoBlockCGSolMgr.hpp" 57 #include "BelosPseudoBlockStochasticCGSolMgr.hpp" 58 #include "BelosGCRODRSolMgr.hpp" 59 #include "BelosRCGSolMgr.hpp" 60 #include "BelosMinresSolMgr.hpp" 61 #include "BelosTFQMRSolMgr.hpp" 64 #include "Teuchos_VerboseObjectParameterListHelpers.hpp" 65 #include "Teuchos_StandardParameterEntryValidators.hpp" 66 #include "Teuchos_ParameterList.hpp" 67 #include "Teuchos_dyn_cast.hpp" 68 #include "Teuchos_ValidatorXMLConverterDB.hpp" 69 #include "Teuchos_StandardValidatorXMLConverters.hpp" 77 template<
class Scalar>
79 template<
class Scalar>
81 template<
class Scalar>
83 template<
class Scalar>
85 template<
class Scalar>
87 template<
class Scalar>
89 template<
class Scalar>
91 template<
class Scalar>
93 template<
class Scalar>
95 template<
class Scalar>
97 template<
class Scalar>
99 template<
class Scalar>
101 template<
class Scalar>
105 const std::string LeftPreconditionerIfUnspecified_name =
"Left Preconditioner If Unspecified";
111 template<
class Scalar>
113 :solverType_(SOLVER_TYPE_PSEUDO_BLOCK_GMRES),
114 convergenceTestFrequency_(1)
116 updateThisValidParamList();
120 template<
class Scalar>
122 const RCP<PreconditionerFactoryBase<Scalar> > &precFactory
124 :solverType_(SOLVER_TYPE_PSEUDO_BLOCK_GMRES)
133 template<
class Scalar>
140 template<
class Scalar>
142 const RCP<PreconditionerFactoryBase<Scalar> > &precFactory,
143 const std::string &precFactoryName
146 TEUCHOS_TEST_FOR_EXCEPT(!precFactory.get());
147 RCP<const Teuchos::ParameterList>
148 precFactoryValidPL = precFactory->getValidParameters();
149 const std::string _precFactoryName =
150 ( precFactoryName !=
"" 152 : ( precFactoryValidPL.get() ? precFactoryValidPL->name() :
"GENERIC PRECONDITIONER FACTORY" )
154 precFactory_ = precFactory;
155 precFactoryName_ = _precFactoryName;
156 updateThisValidParamList();
160 template<
class Scalar>
161 RCP<PreconditionerFactoryBase<Scalar> >
168 template<
class Scalar>
170 RCP<PreconditionerFactoryBase<Scalar> > *precFactory,
171 std::string *precFactoryName
174 if(precFactory) *precFactory = precFactory_;
175 if(precFactoryName) *precFactoryName = precFactoryName_;
176 precFactory_ = Teuchos::null;
177 precFactoryName_ =
"";
178 updateThisValidParamList();
182 template<
class Scalar>
184 const LinearOpSourceBase<Scalar> &fwdOpSrc
187 if(precFactory_.get())
188 return precFactory_->isCompatible(fwdOpSrc);
193 template<
class Scalar>
194 RCP<LinearOpWithSolveBase<Scalar> >
201 template<
class Scalar>
203 const RCP<
const LinearOpSourceBase<Scalar> > &fwdOpSrc,
204 LinearOpWithSolveBase<Scalar> *Op,
205 const ESupportSolveUse supportSolveUse
209 initializeOpImpl(fwdOpSrc,null,null,
false,Op,supportSolveUse);
213 template<
class Scalar>
215 const RCP<
const LinearOpSourceBase<Scalar> > &fwdOpSrc,
216 LinearOpWithSolveBase<Scalar> *Op
220 initializeOpImpl(fwdOpSrc,null,null,
true,Op,SUPPORT_SOLVE_UNSPECIFIED);
224 template<
class Scalar>
226 const EPreconditionerInputType precOpType
229 if(precFactory_.get())
231 return (precOpType==PRECONDITIONER_INPUT_TYPE_AS_OPERATOR);
235 template<
class Scalar>
237 const RCP<
const LinearOpSourceBase<Scalar> > &fwdOpSrc,
238 const RCP<
const PreconditionerBase<Scalar> > &prec,
239 LinearOpWithSolveBase<Scalar> *Op,
240 const ESupportSolveUse supportSolveUse
244 initializeOpImpl(fwdOpSrc,null,prec,
false,Op,supportSolveUse);
248 template<
class Scalar>
250 const RCP<
const LinearOpSourceBase<Scalar> > &fwdOpSrc,
251 const RCP<
const LinearOpSourceBase<Scalar> > &approxFwdOpSrc,
252 LinearOpWithSolveBase<Scalar> *Op,
253 const ESupportSolveUse supportSolveUse
257 initializeOpImpl(fwdOpSrc,approxFwdOpSrc,null,
false,Op,supportSolveUse);
261 template<
class Scalar>
263 LinearOpWithSolveBase<Scalar> *Op,
264 RCP<
const LinearOpSourceBase<Scalar> > *fwdOpSrc,
265 RCP<
const PreconditionerBase<Scalar> > *prec,
266 RCP<
const LinearOpSourceBase<Scalar> > *approxFwdOpSrc,
267 ESupportSolveUse *supportSolveUse
271 TEUCHOS_TEST_FOR_EXCEPT(Op==NULL);
275 RCP<const LinearOpSourceBase<Scalar> >
277 RCP<const PreconditionerBase<Scalar> >
282 RCP<const LinearOpSourceBase<Scalar> >
286 if(fwdOpSrc) *fwdOpSrc = _fwdOpSrc;
287 if(prec) *prec = _prec;
288 if(approxFwdOpSrc) *approxFwdOpSrc = _approxFwdOpSrc;
289 if(supportSolveUse) *supportSolveUse = _supportSolveUse;
296 template<
class Scalar>
298 RCP<Teuchos::ParameterList>
const& paramList
301 TEUCHOS_TEST_FOR_EXCEPT(paramList.get()==NULL);
302 paramList->validateParametersAndSetDefaults(*this->getValidParameters(), 1);
303 paramList_ = paramList;
305 Teuchos::getIntegralValue<EBelosSolverType>(*paramList_, SolverType_name);
306 convergenceTestFrequency_ =
307 Teuchos::getParameter<int>(*paramList_, ConvergenceTestFrequency_name);
308 Teuchos::readVerboseObjectSublist(&*paramList_,
this);
312 template<
class Scalar>
313 RCP<Teuchos::ParameterList>
320 template<
class Scalar>
321 RCP<Teuchos::ParameterList>
324 RCP<Teuchos::ParameterList> _paramList = paramList_;
325 paramList_ = Teuchos::null;
330 template<
class Scalar>
331 RCP<const Teuchos::ParameterList>
338 template<
class Scalar>
339 RCP<const Teuchos::ParameterList>
342 return thisValidParamList_;
349 template<
class Scalar>
352 std::ostringstream oss;
353 oss <<
"Thyra::BelosLinearOpWithSolveFactory";
364 template<
class Scalar>
365 RCP<const Teuchos::ParameterList>
369 using Teuchos::tuple;
370 using Teuchos::setStringToIntegralParameter;
371 Teuchos::ValidatorXMLConverterDB::addConverter(
372 Teuchos::DummyObjectGetter<
373 Teuchos::StringToIntegralParameterEntryValidator<EBelosSolverType>
375 Teuchos::DummyObjectGetter<Teuchos::StringToIntegralValidatorXMLConverter<
376 EBelosSolverType> >::getDummyObject());
378 typedef MultiVectorBase<Scalar> MV_t;
379 typedef LinearOpBase<Scalar> LO_t;
380 static RCP<Teuchos::ParameterList> validParamList;
381 if(validParamList.get()==NULL) {
382 validParamList = Teuchos::rcp(
new Teuchos::ParameterList(
"BelosLinearOpWithSolveFactory"));
383 setStringToIntegralParameter<EBelosSolverType>(
384 SolverType_name, SolverType_default,
385 "Type of linear solver algorithm to use.",
388 "Pseudo Block GMRES",
391 "Pseudo Block Stochastic CG",
398 "Block GMRES solver for nonsymmetric linear systems. It can also solve " 399 "single right-hand side systems, and can also perform Flexible GMRES " 400 "(where the preconditioner may change at every iteration, for example " 401 "for inner-outer iterations), by setting options in the \"Block GMRES\" " 404 "GMRES solver for nonsymmetric linear systems, that performs single " 405 "right-hand side solves on multiple right-hand sides at once. It " 406 "exploits operator multivector multiplication in order to amortize " 407 "global communication costs. Individual linear systems are deflated " 408 "out as they are solved.",
410 "Block CG solver for symmetric (Hermitian in complex arithmetic) " 411 "positive definite linear systems. It can also solve single " 412 "right-hand-side systems.",
414 "CG solver that performs single right-hand side CG on multiple right-hand " 415 "sides at once. It exploits operator multivector multiplication in order " 416 "to amortize global communication costs. Individual linear systems are " 417 "deflated out as they are solved.",
419 "stochastic CG solver that performs single right-hand side CG on multiple right-hand " 420 "sides at once. It exploits operator multivector multiplication in order " 421 "to amortize global communication costs. Individual linear systems are " 422 "deflated out as they are solved. [EXPERIMENTAL]",
424 "Variant of GMRES that performs subspace recycling to accelerate " 425 "convergence for sequences of solves with related linear systems. " 426 "Individual linear systems are deflated out as they are solved. " 427 "The current implementation only supports real-valued Scalar types.",
429 "CG solver for symmetric (Hermitian in complex arithmetic) positive " 430 "definite linear systems, that performs subspace recycling to " 431 "accelerate convergence for sequences of related linear systems.",
433 "MINRES solver for symmetric indefinite linear systems. It performs " 434 "single-right-hand-side solves on multiple right-hand sides sequentially.",
436 "TFQMR (Transpose-Free QMR) solver for nonsymmetric linear systems." 438 tuple<EBelosSolverType>(
439 SOLVER_TYPE_BLOCK_GMRES,
440 SOLVER_TYPE_PSEUDO_BLOCK_GMRES,
441 SOLVER_TYPE_BLOCK_CG,
442 SOLVER_TYPE_PSEUDO_BLOCK_CG,
443 SOLVER_TYPE_PSEUDO_BLOCK_STOCHASTIC_CG,
451 validParamList->set(ConvergenceTestFrequency_name, as<int>(1),
452 "Number of linear solver iterations to skip between applying" 453 " user-defined convergence test.");
455 LeftPreconditionerIfUnspecified_name,
false,
456 "If the preconditioner does not specify if it is left or right, and this\n" 457 "option is set to true, put the preconditioner on the left side.\n" 458 "Historically, preconditioning is on the right. Some solvers may not\n" 459 "support left preconditioning.");
460 Teuchos::ParameterList
461 &solverTypesSL = validParamList->sublist(SolverTypes_name);
463 Belos::BlockGmresSolMgr<Scalar,MV_t,LO_t> mgr;
464 solverTypesSL.sublist(BlockGMRES_name).setParameters(
465 *mgr.getValidParameters()
469 Belos::PseudoBlockGmresSolMgr<Scalar,MV_t,LO_t> mgr;
470 solverTypesSL.sublist(PseudoBlockGMRES_name).setParameters(
471 *mgr.getValidParameters()
475 Belos::BlockCGSolMgr<Scalar,MV_t,LO_t> mgr;
476 solverTypesSL.sublist(BlockCG_name).setParameters(
477 *mgr.getValidParameters()
481 Belos::PseudoBlockCGSolMgr<Scalar,MV_t,LO_t> mgr;
482 solverTypesSL.sublist(PseudoBlockCG_name).setParameters(
483 *mgr.getValidParameters()
487 Belos::PseudoBlockStochasticCGSolMgr<Scalar,MV_t,LO_t> mgr;
488 solverTypesSL.sublist(PseudoBlockStochasticCG_name).setParameters(
489 *mgr.getValidParameters()
493 Belos::GCRODRSolMgr<Scalar,MV_t,LO_t> mgr;
494 solverTypesSL.sublist(GCRODR_name).setParameters(
495 *mgr.getValidParameters()
499 Belos::RCGSolMgr<Scalar,MV_t,LO_t> mgr;
500 solverTypesSL.sublist(RCG_name).setParameters(
501 *mgr.getValidParameters()
505 Belos::MinresSolMgr<Scalar,MV_t,LO_t> mgr;
506 solverTypesSL.sublist(MINRES_name).setParameters(
507 *mgr.getValidParameters()
511 Belos::TFQMRSolMgr<Scalar,MV_t,LO_t> mgr;
512 solverTypesSL.sublist(TFQMR_name).setParameters(
513 *mgr.getValidParameters()
517 return validParamList;
521 template<
class Scalar>
522 void BelosLinearOpWithSolveFactory<Scalar>::updateThisValidParamList()
524 thisValidParamList_ = Teuchos::rcp(
525 new Teuchos::ParameterList(*generateAndGetValidParameters())
527 Teuchos::setupVerboseObjectSublist(&*thisValidParamList_);
531 template<
class Scalar>
532 void BelosLinearOpWithSolveFactory<Scalar>::initializeOpImpl(
533 const RCP<
const LinearOpSourceBase<Scalar> > &fwdOpSrc,
534 const RCP<
const LinearOpSourceBase<Scalar> > &approxFwdOpSrc,
535 const RCP<
const PreconditionerBase<Scalar> > &prec_in,
536 const bool reusePrec,
537 LinearOpWithSolveBase<Scalar> *Op,
538 const ESupportSolveUse supportSolveUse
543 using Teuchos::set_extra_data;
544 typedef Teuchos::ScalarTraits<Scalar> ST;
545 typedef MultiVectorBase<Scalar> MV_t;
546 typedef LinearOpBase<Scalar> LO_t;
548 const RCP<Teuchos::FancyOStream> out = this->getOStream();
549 const Teuchos::EVerbosityLevel verbLevel = this->getVerbLevel();
550 Teuchos::OSTab tab(out);
551 if(out.get() &&
static_cast<int>(verbLevel) > static_cast<int>(Teuchos::VERB_LOW))
552 *out <<
"\nEntering Thyra::BelosLinearOpWithSolveFactory<"<<ST::name()<<
">::initializeOpImpl(...) ...\n";
559 TEUCHOS_TEST_FOR_EXCEPT(Op==NULL);
560 TEUCHOS_TEST_FOR_EXCEPT(fwdOpSrc.get()==NULL);
561 TEUCHOS_TEST_FOR_EXCEPT(fwdOpSrc->getOp().get()==NULL);
562 RCP<const LinearOpBase<Scalar> >
563 fwdOp = fwdOpSrc->getOp(),
564 approxFwdOp = ( approxFwdOpSrc.get() ? approxFwdOpSrc->getOp() : Teuchos::null );
570 BelosLinearOpWithSolve<Scalar>
571 *belosOp = &Teuchos::dyn_cast<BelosLinearOpWithSolve<Scalar> >(*Op);
577 RCP<PreconditionerBase<Scalar> > myPrec = Teuchos::null;
578 RCP<const PreconditionerBase<Scalar> > prec = Teuchos::null;
585 if(precFactory_.get()) {
587 ( !belosOp->isExternalPrec()
588 ? Teuchos::rcp_const_cast<PreconditionerBase<Scalar> >(belosOp->extract_prec())
591 bool hasExistingPrec =
false;
593 hasExistingPrec =
true;
598 hasExistingPrec =
false;
599 myPrec = precFactory_->createPrec();
601 if( hasExistingPrec && reusePrec ) {
606 if(approxFwdOp.get())
607 precFactory_->initializePrec(approxFwdOpSrc,&*myPrec);
609 precFactory_->initializePrec(fwdOpSrc,&*myPrec);
619 bool oldIsExternalPrec =
false;
620 RCP<Belos::LinearProblem<Scalar,MV_t,LO_t> > oldLP = Teuchos::null;
621 RCP<Belos::SolverManager<Scalar,MV_t,LO_t> > oldIterSolver = Teuchos::null;
622 RCP<const LinearOpSourceBase<Scalar> > oldFwdOpSrc = Teuchos::null;
623 RCP<const LinearOpSourceBase<Scalar> > oldApproxFwdOpSrc = Teuchos::null;
624 ESupportSolveUse oldSupportSolveUse = SUPPORT_SOLVE_UNSPECIFIED;
626 belosOp->uninitialize( &oldLP, NULL, &oldIterSolver, &oldFwdOpSrc,
627 NULL, &oldIsExternalPrec, &oldApproxFwdOpSrc, &oldSupportSolveUse );
634 typedef Belos::LinearProblem<Scalar,MV_t,LO_t> LP_t;
636 if (oldLP != Teuchos::null) {
640 lp = rcp(
new LP_t());
647 lp->setOperator(fwdOp);
654 RCP<const LinearOpBase<Scalar> > unspecified = prec->getUnspecifiedPrecOp();
655 RCP<const LinearOpBase<Scalar> > left = prec->getLeftPrecOp();
656 RCP<const LinearOpBase<Scalar> > right = prec->getRightPrecOp();
657 TEUCHOS_TEST_FOR_EXCEPTION(
658 !( left.get() || right.get() || unspecified.get() ), std::logic_error
659 ,
"Error, at least one preconditoner linear operator objects must be set!" 661 if(unspecified.get()) {
662 if (paramList_->get<
bool>(LeftPreconditionerIfUnspecified_name,
false))
663 lp->setLeftPrec(unspecified);
665 lp->setRightPrec(unspecified);
669 TEUCHOS_TEST_FOR_EXCEPTION(
670 left.get(),std::logic_error
671 ,
"Error, we can not currently handle a left preconditioner!" 673 lp->setRightPrec(right);
677 set_extra_data<RCP<PreconditionerBase<Scalar> > >(myPrec,
"Belos::InternalPrec",
678 Teuchos::inOutArg(lp), Teuchos::POST_DESTROY,
false);
680 else if(prec.get()) {
681 set_extra_data<RCP<const PreconditionerBase<Scalar> > >(prec,
"Belos::ExternalPrec",
682 Teuchos::inOutArg(lp), Teuchos::POST_DESTROY,
false);
689 typedef Belos::SolverManager<Scalar,MV_t,LO_t> IterativeSolver_t;
690 RCP<IterativeSolver_t> iterativeSolver = Teuchos::null;
691 RCP<Teuchos::ParameterList> solverPL = Teuchos::rcp(
new Teuchos::ParameterList() );
693 switch(solverType_) {
694 case SOLVER_TYPE_BLOCK_GMRES:
697 if(paramList_.get()) {
698 Teuchos::ParameterList &solverTypesPL = paramList_->sublist(SolverTypes_name);
699 Teuchos::ParameterList &gmresPL = solverTypesPL.sublist(BlockGMRES_name);
700 solverPL = Teuchos::rcp( &gmresPL,
false );
703 if (oldIterSolver != Teuchos::null) {
704 iterativeSolver = oldIterSolver;
705 iterativeSolver->setProblem( lp );
706 iterativeSolver->setParameters( solverPL );
709 iterativeSolver = rcp(
new Belos::BlockGmresSolMgr<Scalar,MV_t,LO_t>(lp,solverPL));
713 case SOLVER_TYPE_PSEUDO_BLOCK_GMRES:
716 if(paramList_.get()) {
717 Teuchos::ParameterList &solverTypesPL = paramList_->sublist(SolverTypes_name);
718 Teuchos::ParameterList &pbgmresPL = solverTypesPL.sublist(PseudoBlockGMRES_name);
719 solverPL = Teuchos::rcp( &pbgmresPL,
false );
724 if (oldIterSolver != Teuchos::null) {
725 iterativeSolver = oldIterSolver;
726 iterativeSolver->setProblem( lp );
727 iterativeSolver->setParameters( solverPL );
730 iterativeSolver = rcp(
new Belos::PseudoBlockGmresSolMgr<Scalar,MV_t,LO_t>(lp,solverPL));
734 case SOLVER_TYPE_BLOCK_CG:
737 if(paramList_.get()) {
738 Teuchos::ParameterList &solverTypesPL = paramList_->sublist(SolverTypes_name);
739 Teuchos::ParameterList &cgPL = solverTypesPL.sublist(BlockCG_name);
740 solverPL = Teuchos::rcp( &cgPL,
false );
743 if (oldIterSolver != Teuchos::null) {
744 iterativeSolver = oldIterSolver;
745 iterativeSolver->setProblem( lp );
746 iterativeSolver->setParameters( solverPL );
749 iterativeSolver = rcp(
new Belos::BlockCGSolMgr<Scalar,MV_t,LO_t>(lp,solverPL));
753 case SOLVER_TYPE_PSEUDO_BLOCK_CG:
756 if(paramList_.get()) {
757 Teuchos::ParameterList &solverTypesPL = paramList_->sublist(SolverTypes_name);
758 Teuchos::ParameterList &pbcgPL = solverTypesPL.sublist(PseudoBlockCG_name);
759 solverPL = Teuchos::rcp( &pbcgPL,
false );
764 if (oldIterSolver != Teuchos::null) {
765 iterativeSolver = oldIterSolver;
766 iterativeSolver->setProblem( lp );
767 iterativeSolver->setParameters( solverPL );
770 iterativeSolver = rcp(
new Belos::PseudoBlockCGSolMgr<Scalar,MV_t,LO_t>(lp,solverPL));
774 case SOLVER_TYPE_PSEUDO_BLOCK_STOCHASTIC_CG:
777 if(paramList_.get()) {
778 Teuchos::ParameterList &solverTypesPL = paramList_->sublist(SolverTypes_name);
779 Teuchos::ParameterList &pbcgPL = solverTypesPL.sublist(PseudoBlockStochasticCG_name);
780 solverPL = Teuchos::rcp( &pbcgPL,
false );
785 if (oldIterSolver != Teuchos::null) {
786 iterativeSolver = oldIterSolver;
787 iterativeSolver->setProblem( lp );
788 iterativeSolver->setParameters( solverPL );
791 iterativeSolver = rcp(
new Belos::PseudoBlockStochasticCGSolMgr<Scalar,MV_t,LO_t>(lp,solverPL));
795 case SOLVER_TYPE_GCRODR:
798 if(paramList_.get()) {
799 Teuchos::ParameterList &solverTypesPL = paramList_->sublist(SolverTypes_name);
800 Teuchos::ParameterList &gcrodrPL = solverTypesPL.sublist(GCRODR_name);
801 solverPL = Teuchos::rcp( &gcrodrPL,
false );
804 if (oldIterSolver != Teuchos::null) {
805 iterativeSolver = oldIterSolver;
806 iterativeSolver->setProblem( lp );
807 iterativeSolver->setParameters( solverPL );
810 iterativeSolver = rcp(
new Belos::GCRODRSolMgr<Scalar,MV_t,LO_t>(lp,solverPL));
814 case SOLVER_TYPE_RCG:
817 if(paramList_.get()) {
818 Teuchos::ParameterList &solverTypesPL = paramList_->sublist(SolverTypes_name);
819 Teuchos::ParameterList &rcgPL = solverTypesPL.sublist(RCG_name);
820 solverPL = Teuchos::rcp( &rcgPL,
false );
823 if (oldIterSolver != Teuchos::null) {
824 iterativeSolver = oldIterSolver;
825 iterativeSolver->setProblem( lp );
826 iterativeSolver->setParameters( solverPL );
829 iterativeSolver = rcp(
new Belos::RCGSolMgr<Scalar,MV_t,LO_t>(lp,solverPL));
833 case SOLVER_TYPE_MINRES:
836 if(paramList_.get()) {
837 Teuchos::ParameterList &solverTypesPL = paramList_->sublist(SolverTypes_name);
838 Teuchos::ParameterList &minresPL = solverTypesPL.sublist(MINRES_name);
839 solverPL = Teuchos::rcp( &minresPL,
false );
842 if (oldIterSolver != Teuchos::null) {
843 iterativeSolver = oldIterSolver;
844 iterativeSolver->setProblem( lp );
845 iterativeSolver->setParameters( solverPL );
848 iterativeSolver = rcp(
new Belos::MinresSolMgr<Scalar,MV_t,LO_t>(lp,solverPL));
852 case SOLVER_TYPE_TFQMR:
855 if(paramList_.get()) {
856 Teuchos::ParameterList &solverTypesPL = paramList_->sublist(SolverTypes_name);
857 Teuchos::ParameterList &minresPL = solverTypesPL.sublist(TFQMR_name);
858 solverPL = Teuchos::rcp( &minresPL,
false );
861 if (oldIterSolver != Teuchos::null) {
862 iterativeSolver = oldIterSolver;
863 iterativeSolver->setProblem( lp );
864 iterativeSolver->setParameters( solverPL );
867 iterativeSolver = rcp(
new Belos::TFQMRSolMgr<Scalar,MV_t,LO_t>(lp,solverPL));
874 TEUCHOS_TEST_FOR_EXCEPT(
true);
883 lp, solverPL, iterativeSolver,
884 fwdOpSrc, prec, myPrec.get()==NULL, approxFwdOpSrc,
885 supportSolveUse, convergenceTestFrequency_
887 belosOp->setOStream(out);
888 belosOp->setVerbLevel(verbLevel);
890 if(paramList_.get()) {
892 paramList_->validateParameters(*this->getValidParameters(),1);
895 if(out.get() &&
static_cast<int>(verbLevel) > static_cast<int>(Teuchos::VERB_LOW))
896 *out <<
"\nLeaving Thyra::BelosLinearOpWithSolveFactory<"<<ST::name()<<
">::initializeOpImpl(...) ...\n";
904 #endif // THYRA_BELOS_LINEAR_OP_WITH_SOLVE_FACTORY_HPP static const std::string BlockGMRES_name
bool isCompatible(const LinearOpSourceBase< Scalar > &fwdOpSrc) const
static const std::string RCG_name
RCP< const LinearOpSourceBase< Scalar > > extract_approxFwdOpSrc()
static const std::string MINRES_name
Concrete LinearOpWithSolveBase subclass in terms of Belos.
void initializeAndReuseOp(const Teuchos::RCP< const LinearOpSourceBase< Scalar > > &fwdOpSrc, LinearOpWithSolveBase< Scalar > *Op) const
RCP< const PreconditionerBase< Scalar > > extract_prec()
bool supportsPreconditionerInputType(const EPreconditionerInputType precOpType) const
static const std::string SolverTypes_name
void unsetPreconditionerFactory(Teuchos::RCP< PreconditionerFactoryBase< Scalar > > *precFactory, std::string *precFactoryName)
Thyra specializations of MultiVecTraits and OperatorTraits.
void setParameterList(Teuchos::RCP< Teuchos::ParameterList > const ¶mList)
static const std::string ConvergenceTestFrequency_name
void initializePreconditionedOp(const Teuchos::RCP< const LinearOpSourceBase< Scalar > > &fwdOpSrc, const Teuchos::RCP< const PreconditionerBase< Scalar > > &prec, LinearOpWithSolveBase< Scalar > *Op, const ESupportSolveUse supportSolveUse) const
static const std::string BlockCG_name
static const std::string PseudoBlockGMRES_name
bool isExternalPrec() const
static const std::string PseudoBlockCG_name
Teuchos::RCP< Teuchos::ParameterList > getNonconstParameterList()
void setPreconditionerFactory(const Teuchos::RCP< PreconditionerFactoryBase< Scalar > > &precFactory, const std::string &precFactoryName)
Teuchos::RCP< const Teuchos::ParameterList > getParameterList() const
bool acceptsPreconditionerFactory() const
Returns true .
Teuchos::RCP< PreconditionerFactoryBase< Scalar > > getPreconditionerFactory() const
LinearOpWithSolveFactoryBase subclass implemented in terms of Belos.
RCP< const LinearOpSourceBase< Scalar > > extract_fwdOpSrc()
static const std::string PseudoBlockStochasticCG_name
static const std::string SolverType_default
Teuchos::RCP< const Teuchos::ParameterList > getValidParameters() const
std::string description() const
BelosLinearOpWithSolveFactory()
Construct without preconditioner factory.
ESupportSolveUse supportSolveUse() const
Teuchos::RCP< LinearOpWithSolveBase< Scalar > > createOp() 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
Teuchos::RCP< Teuchos::ParameterList > unsetParameterList()
static const std::string SolverType_name
static const std::string TFQMR_name
void initializeApproxPreconditionedOp(const Teuchos::RCP< const LinearOpSourceBase< Scalar > > &fwdOpSrc, const Teuchos::RCP< const LinearOpSourceBase< Scalar > > &approxFwdOpSrc, LinearOpWithSolveBase< Scalar > *Op, const ESupportSolveUse supportSolveUse) const
void initializeOp(const Teuchos::RCP< const LinearOpSourceBase< Scalar > > &fwdOpSrc, LinearOpWithSolveBase< Scalar > *Op, const ESupportSolveUse supportSolveUse) const
static const std::string GCRODR_name