46#ifndef MUELU_BELOSSMOOTHER_DEF_HPP
47#define MUELU_BELOSSMOOTHER_DEF_HPP
51#if defined(HAVE_MUELU_BELOS)
53#include <Teuchos_ParameterList.hpp>
55#include <Xpetra_CrsMatrix.hpp>
56#include <Xpetra_Matrix.hpp>
57#include <Xpetra_MultiVectorFactory.hpp>
58#ifdef HAVE_XPETRA_TPETRA
59#include <Xpetra_TpetraMultiVector.hpp>
65#include "MueLu_Utilities.hpp"
68#include <BelosLinearProblem.hpp>
69#include <BelosSolverFactory.hpp>
75 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
79 bool solverSupported =
false;
80#ifdef HAVE_MUELU_TPETRA
81 Belos::SolverFactory<Scalar,tMV,tOP> tFactory;
82 solverSupported = solverSupported || tFactory.isSupported(type);
104 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
109 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
112 this->Input(currentLevel,
"A");
116 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
120 A_ = Factory::Get< RCP<Matrix> >(currentLevel,
"A");
121 SetupBelos(currentLevel);
123 this->GetOStream(
Statistics1) << description() << std::endl;
127 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
130 bool useTpetra = A_->getRowMap()->lib() == Xpetra::UseTpetra;
133#ifdef HAVE_MUELU_TPETRA
134 tBelosProblem_ = rcp(
new Belos::LinearProblem<Scalar, tMV, tOP>());
136 tBelosProblem_->setOperator(tA);
138 Belos::SolverFactory<SC, tMV, tOP> solverFactory;
139 tSolver_ = solverFactory.create(type_, rcpFromRef(
const_cast<ParameterList&
>(this->GetParameterList())));
140 tSolver_->setProblem(tBelosProblem_);
143 TEUCHOS_ASSERT(
false);
147 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
151 if (A_->getRowMap()->lib() == Xpetra::UseTpetra) {
152#ifdef HAVE_MUELU_TPETRA
153 if (InitialGuessIsZero) {
159 tBelosProblem_->setInitResVec(tpB);
160 tBelosProblem_->setProblem(tpX, tpB);
164 typedef Teuchos::ScalarTraits<Scalar> TST;
166 RCP<MultiVector> Correction = MultiVectorFactory::Build(A_->getDomainMap(), X.getNumVectors());
171 tBelosProblem_->setInitResVec(tpB);
172 tBelosProblem_->setProblem(tpX, tpB);
175 X.update(TST::one(), *Correction, TST::one());
179 TEUCHOS_ASSERT(
false);
184 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
186 RCP<BelosSmoother> smoother = rcp(
new BelosSmoother(*
this) );
187 smoother->SetParameterList(this->GetParameterList());
191 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
193 std::ostringstream out;
195 if (A_->getRowMap()->lib() == Xpetra::UseTpetra) {
196#ifdef MUELU_HAVE_TPETRA
197 out << tSolver_->description();
201 out <<
"BELOS {type = " << type_ <<
"}";
206 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
211 out0 <<
"Parameter list: " << std::endl;
212 Teuchos::OSTab tab2(out);
213 out << this->GetParameterList();
217#ifdef MUELU_HAVE_TPETRA
218 if (tSolver_ != Teuchos::null) {
219 Teuchos::OSTab tab2(out);
220 out << *tSolver_ << std::endl;
225 if (verbLevel &
Debug) {
226 if (A_->getRowMap()->lib() == Xpetra::UseTpetra) {
227#ifdef MUELU_HAVE_TPETRA
230 <<
"RCP<solver_>: " << tSolver_ << std::endl;
236 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
238 return Teuchos::OrdinalTraits<size_t>::invalid();
#define MUELU_DESCRIBE
Helper macro for implementing Describable::describe() for BaseClass objects.
Class that encapsulates Belos smoothers.
void Apply(MultiVector &X, const MultiVector &B, bool InitialGuessIsZero=false) const
Apply the preconditioner.
void print(Teuchos::FancyOStream &out, const VerbLevel verbLevel=Default) const
Print the object with some verbosity level to an FancyOStream object.
friend class BelosSmoother
Constructor.
void SetParameterList(const Teuchos::ParameterList ¶mList)
Set parameters from a parameter list and return with default values.
void DeclareInput(Level ¤tLevel) const
Input.
void SetupBelos(Level ¤tLevel)
std::string description() const
Return a simple one-line description of this object.
size_t getNodeSmootherComplexity() const
For diagnostic purposes.
void Setup(Level ¤tLevel)
Set up the smoother.
RCP< SmootherPrototype > Copy() const
Exception throws to report errors in the internal logical of the program.
Timer to be used in factories. Similar to Monitor but with additional timers.
Class that holds all level-specific information.
virtual void SetParameterList(const Teuchos::ParameterList ¶mList)=0
Set parameters from a parameter list and return with default values.
void declareConstructionOutcome(bool fail, std::string msg)
bool IsSetup() const
Get the state of a smoother prototype.
static RCP< Tpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op2NonConstTpetraCrs(RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Op)
static RCP< const Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > MV2TpetraMV(RCP< Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > const vec)
Helper utility to pull out the underlying Tpetra objects from an Xpetra object.
static RCP< Tpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > MV2NonConstTpetraMV(RCP< Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > vec)
static RCP< Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Residual(const Xpetra::Operator< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Op, const Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &X, const Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &RHS)
Namespace for MueLu classes and methods.
@ Debug
Print additional debugging information.
@ Statistics1
Print more statistics.
@ External
Print external lib objects.
@ Parameters1
Print class parameters (more parameters, more verbose)