42#ifndef __Thyra_BelosTpetraSolverAdapter_hpp
43#define __Thyra_BelosTpetraSolverAdapter_hpp
45#include "BelosConfigDefs.hpp"
46#include "BelosLinearProblem.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>
73 using converter = Thyra::TpetraOperatorVectorExtraction<SC>;
79 virtual Teuchos::RCP<Belos::SolverManager<SC, MV, OP> >
clone ()
const override = 0;
82 void setProblem(
const Teuchos::RCP<Belos::LinearProblem<SC, MV, OP> > &problem )
override {
86 const Belos::LinearProblem<SC, MV, OP>&
getProblem()
const override {
91 void setParameters(
const Teuchos::RCP<Teuchos::ParameterList> ¶ms )
override {
102 Teuchos::RCP<Teuchos::ParameterList> params = Teuchos::parameterList();
104 bool defaultValues =
true;
113 void reset(
const Belos::ResetType type )
override {
114 if ((type & Belos::Problem) && !Teuchos::is_null(
problem_))
118 typename Teuchos::ScalarTraits<SC>::magnitudeType
achievedTol()
const {
123 Belos::ReturnType
solve()
override {
125 auto tpetraA = converter::getConstTpetraOperator(A);
129 if (lM != Teuchos::null) {
130 auto tpetraLM = converter::getConstTpetraOperator(lM);
134 if (rM != Teuchos::null) {
135 auto tpetraRM = converter::getConstTpetraOperator(rM);
141 auto tpetraB = converter::getConstTpetraMultiVector(B);
142 auto tpetraX = converter::getTpetraMultiVector(X);
147 Belos::ReturnType belos_output = (
solver_output.converged ? Belos::Converged : Belos::Unconverged);
160 Teuchos::RCP<Belos::LinearProblem<SC, MV, OP> >
problem_;
165 template<
class SC,
class MV,
class OP>
179 Teuchos::RCP<Belos::SolverManager<SC, MV, OP> >
clone ()
const override {
186 template<
class SC,
class MV,
class OP>
200 Teuchos::RCP<Belos::SolverManager<SC, MV, OP> >
clone ()
const override {
207 template<
class SC,
class MV,
class OP>
221 Teuchos::RCP<Belos::SolverManager<SC, MV, OP> >
clone ()
const override {
228 template<
class SC,
class MV,
class OP>
242 Teuchos::RCP<Belos::SolverManager<SC, MV, OP> >
clone ()
const override {
BelosTpetra::Impl::GmresPipeline< SC > tpetra_solver_type
BelosTpetraGmresPipeline()
constructor
Teuchos::RCP< Belos::SolverManager< SC, MV, OP > > clone() const override
clone for Inverted Injection (DII)
BelosTpetraGmresSingleReduce()
constructor
BelosTpetra::Impl::GmresSingleReduce< SC > tpetra_solver_type
Teuchos::RCP< Belos::SolverManager< SC, MV, OP > > clone() const override
clone for Inverted Injection (DII)
BelosTpetra::Impl::GmresSstep< SC > tpetra_solver_type
Teuchos::RCP< Belos::SolverManager< SC, MV, OP > > clone() const override
clone for Inverted Injection (DII)
BelosTpetraGmresSstep()
constructor
BelosTpetra::Impl::Gmres< SC > tpetra_solver_type
BelosTpetraGmres()
constructor
Teuchos::RCP< Belos::SolverManager< SC, MV, OP > > clone() const override
clone for Inverted Injection (DII)
bool isLOADetected() const override
virtual Teuchos::RCP< Belos::SolverManager< SC, MV, OP > > clone() const override=0
clone for Inverted Injection (DII)
Teuchos::RCP< Belos::LinearProblem< SC, MV, OP > > problem_
The linear problem to solve.
void setProblem(const Teuchos::RCP< Belos::LinearProblem< SC, MV, OP > > &problem) override
set/get problem
void setParameters(const Teuchos::RCP< Teuchos::ParameterList > ¶ms) override
set/get parameters
Belos::ReturnType solve() override
solve
BelosTpetra::Impl::Krylov< SC > tpetra_base_solver_type
Teuchos::RCP< const Teuchos::ParameterList > getCurrentParameters() const override
void reset(const Belos::ResetType type) override
Thyra::TpetraOperatorVectorExtraction< SC > converter
BelosTpetra::Impl::SolverOutput< SC > solver_output
Teuchos::RCP< const Teuchos::ParameterList > getValidParameters() const override
Teuchos::RCP< Teuchos::ParameterList > params_
Current parameter list.
const Belos::LinearProblem< SC, MV, OP > & getProblem() const override
Teuchos::ScalarTraits< SC >::magnitudeType achievedTol() const
Teuchos::RCP< tpetra_base_solver_type > tpetra_solver
int getNumIters() const override
Get the iteration count for the most recent call to solve().
BelosTpetraKrylov()=default
constructor