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>;
73 using converter = Thyra::TpetraOperatorVectorExtraction<SC>;
76 BelosTpetraKrylov() =
default;
79 virtual Teuchos::RCP<Belos::SolverManager<SC, MV, OP> > clone ()
const override = 0;
91 void setParameters(
const Teuchos::RCP<Teuchos::ParameterList> ¶ms )
override {
93 tpetra_solver->setParameters(*params);
96 Teuchos::RCP<const Teuchos::ParameterList> getCurrentParameters()
const override {
100 Teuchos::RCP<const Teuchos::ParameterList> getValidParameters()
const override {
102 Teuchos::RCP<Teuchos::ParameterList> params = Teuchos::parameterList();
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();
118 typename Teuchos::ScalarTraits<SC>::magnitudeType achievedTol()
const {
119 return tpetra_solver->achievedTol();
124 auto A = problem_->getOperator();
125 auto tpetraA = converter::getConstTpetraOperator(A);
126 tpetra_solver->setMatrix(tpetraA);
128 auto lM = problem_->getLeftPrec();
129 if (lM != Teuchos::null) {
130 auto tpetraLM = converter::getConstTpetraOperator(lM);
131 tpetra_solver->setLeftPrec(tpetraLM);
133 auto rM = problem_->getRightPrec();
134 if (rM != Teuchos::null) {
135 auto tpetraRM = converter::getConstTpetraOperator(rM);
136 tpetra_solver->setRightPrec(tpetraRM);
139 auto B = this->problem_->getRHS();
140 auto X = this->problem_->getLHS();
141 auto tpetraB = converter::getConstTpetraMultiVector(B);
142 auto tpetraX = converter::getTpetraMultiVector(X);
144 solver_output = tpetra_solver->solve(*tpetraX, *tpetraB);
153 Teuchos::RCP<tpetra_base_solver_type> tpetra_solver;
154 BelosTpetra::Impl::SolverOutput<SC> solver_output;
157 Teuchos::RCP<Teuchos::ParameterList> params_;
160 Teuchos::RCP<Belos::LinearProblem<SC, MV, OP> > problem_;
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 ());
179 Teuchos::RCP<Belos::SolverManager<SC, MV, OP> > clone ()
const override {
180 return Teuchos::rcp(
new BelosTpetraGmres<SC, MV, OP>);
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 ());
200 Teuchos::RCP<Belos::SolverManager<SC, MV, OP> > clone ()
const override {
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 ());
221 Teuchos::RCP<Belos::SolverManager<SC, MV, OP> > clone ()
const override {
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 ());
242 Teuchos::RCP<Belos::SolverManager<SC, MV, OP> > clone ()
const override {
243 return Teuchos::rcp(
new BelosTpetraGmresSstep<SC, MV, OP>);
virtual void setProblem(const Teuchos::RCP< LinearProblem< ScalarType, MV, OP > > &problem)=0