42 #ifndef __Thyra_BelosTpetraSolverAdapter_hpp
43 #define __Thyra_BelosTpetraSolverAdapter_hpp
50 #include "Thyra_MultiVectorBase.hpp"
51 #include "Thyra_TpetraThyraWrappers.hpp"
53 #include "Teuchos_ParameterList.hpp"
54 #include "Teuchos_ParameterListAcceptorDefaultBase.hpp"
56 #include "Tpetra_CrsMatrix.hpp"
57 #include "BelosTpetraAdapter.hpp"
58 #include "solvers/Belos_Tpetra_Krylov_parameters.hpp"
59 #include "solvers/Belos_Tpetra_Krylov.hpp"
60 #include "solvers/Belos_Tpetra_Gmres.hpp"
61 #include "solvers/Belos_Tpetra_GmresPipeline.hpp"
62 #include "solvers/Belos_Tpetra_GmresSingleReduce.hpp"
63 #include "solvers/Belos_Tpetra_GmresSstep.hpp"
69 template<
class SC,
class MV,
class OP>
72 using tpetra_base_solver_type = BelosTpetra::Impl::Krylov<SC>;
76 BelosTpetraKrylov() =
default;
104 bool defaultValues =
true;
105 tpetra_solver->getParameters (*params, defaultValues);
111 int getNumIters()
const override {
return solver_output.numIters; }
112 bool isLOADetected()
const override {
return false; }
115 problem_->setProblem();
119 return tpetra_solver->achievedTol();
124 auto A = problem_->getOperator();
126 tpetra_solver->setMatrix(tpetraA);
128 auto lM = problem_->getLeftPrec();
129 if (lM != Teuchos::null) {
131 tpetra_solver->setLeftPrec(tpetraLM);
133 auto rM = problem_->getRightPrec();
134 if (rM != Teuchos::null) {
136 tpetra_solver->setRightPrec(tpetraRM);
139 auto B = this->problem_->getRHS();
140 auto X = this->problem_->getLHS();
144 solver_output = tpetra_solver->solve(*tpetraX, *tpetraB);
154 BelosTpetra::Impl::SolverOutput<SC> solver_output;
165 template<
class SC,
class MV,
class OP>
166 class BelosTpetraGmres :
public BelosTpetraKrylov<SC, MV, OP> {
168 using krylov_base_solver_type = BelosTpetraKrylov<SC, MV, OP>;
169 using tpetra_solver_type = BelosTpetra::Impl::Gmres<SC>;
173 krylov_base_solver_type ()
175 this->tpetra_solver =
Teuchos::rcp(
new tpetra_solver_type ());
186 template<
class SC,
class MV,
class OP>
187 class BelosTpetraGmresPipeline :
public BelosTpetraKrylov<SC, MV, OP> {
189 using krylov_base_solver_type = BelosTpetraKrylov<SC, MV, OP>;
190 using tpetra_solver_type = BelosTpetra::Impl::GmresPipeline<SC>;
193 BelosTpetraGmresPipeline() :
194 krylov_base_solver_type ()
196 this->tpetra_solver =
Teuchos::rcp(
new tpetra_solver_type ());
201 return Teuchos::rcp(
new BelosTpetraGmresPipeline<SC, MV, OP>);
207 template<
class SC,
class MV,
class OP>
208 class BelosTpetraGmresSingleReduce :
public BelosTpetraKrylov<SC, MV, OP> {
210 using krylov_base_solver_type = BelosTpetraKrylov<SC, MV, OP>;
211 using tpetra_solver_type = BelosTpetra::Impl::GmresSingleReduce<SC>;
214 BelosTpetraGmresSingleReduce() :
215 krylov_base_solver_type ()
217 this->tpetra_solver =
Teuchos::rcp(
new tpetra_solver_type ());
222 return Teuchos::rcp(
new BelosTpetraGmresSingleReduce<SC, MV, OP>);
228 template<
class SC,
class MV,
class OP>
229 class BelosTpetraGmresSstep :
public BelosTpetraKrylov<SC, MV, OP> {
231 using krylov_base_solver_type = BelosTpetraKrylov<SC, MV, OP>;
232 using tpetra_solver_type = BelosTpetra::Impl::GmresSstep<SC>;
235 BelosTpetraGmresSstep() :
236 krylov_base_solver_type ()
238 this->tpetra_solver =
Teuchos::rcp(
new tpetra_solver_type ());
243 return Teuchos::rcp(
new BelosTpetraGmresSstep<SC, MV, OP>);
249 #endif // __Thyra_TpetraSolverAdapter_hpp
virtual void setProblem(const Teuchos::RCP< LinearProblem< ScalarType, MV, OP > > &problem)=0
bool is_null(const std::shared_ptr< T > &p)
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
ParameterList & setParameters(const ParameterList &source)