46#ifndef MUELU_REFMAXWELLSMOOTHER_DEF_HPP
47#define MUELU_REFMAXWELLSMOOTHER_DEF_HPP
51#include <Teuchos_ParameterList.hpp>
53#include <Xpetra_CrsMatrix.hpp>
54#include <Xpetra_Matrix.hpp>
55#include <Xpetra_MultiVectorFactory.hpp>
60#include "MueLu_Utilities.hpp"
67 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
71 const bool solverSupported = (
type_ ==
"RefMaxwell");
77 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
82 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
84 this->Input(currentLevel,
"A");
85 this->Input(currentLevel,
"D0");
86 this->Input(currentLevel,
"M1");
87 this->Input(currentLevel,
"Ms");
88 this->Input(currentLevel,
"M0inv");
89 this->Input(currentLevel,
"Coordinates");
92 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
96 SetupRefMaxwell(currentLevel);
98 this->GetOStream(
Statistics1) << description() << std::endl;
102 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
105 using coordinateType =
typename Teuchos::ScalarTraits<Scalar>::coordinateType;
106 using RealValuedMultiVector = Xpetra::MultiVector<coordinateType,LO,GO,NO>;
108 auto SM = currentLevel.
Get<RCP<Matrix> >(
"A");
109 auto D0 = currentLevel.
Get<RCP<Matrix> >(
"D0");
110 auto M1 = currentLevel.
Get<RCP<Matrix> >(
"M1");
111 auto Ms = currentLevel.
Get<RCP<Matrix> >(
"Ms");
112 auto M0inv = currentLevel.
Get<RCP<Matrix> >(
"M0inv");
113 auto coords = currentLevel.
Get<RCP<RealValuedMultiVector> >(
"Coordinates");
114 auto params = this->GetParameterList();
117 Teuchos::null, coords,
119 std::ostringstream oss;
120 op_->describe(*Teuchos::fancyOStream(Teuchos::rcpFromRef(oss)));
121 cachedDescription_ = oss.str();
124 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
128 if (InitialGuessIsZero) {
134 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
137 smoother->SetParameterList(this->GetParameterList());
141 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
143 std::ostringstream out;
145 out << op_->description();
146 out << std::endl << std::endl;
154 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
159 out0 <<
"Parameter list: " << std::endl;
160 Teuchos::OSTab tab2(out);
161 out << this->GetParameterList();
165 if (op_ != Teuchos::null) {
166 Teuchos::OSTab tab2(out);
167 out << *op_ << std::endl << std::endl;
178 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
180 return Teuchos::OrdinalTraits<size_t>::invalid();
#define MUELU_DESCRIBE
Helper macro for implementing Describable::describe() for BaseClass objects.
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.
T & Get(const std::string &ename, const FactoryBase *factory=NoFactory::get())
Get data without decrementing associated storage counter (i.e., read-only access)....
virtual void SetParameterList(const Teuchos::ParameterList ¶mList)=0
Set parameters from a parameter list and return with default values.
Class that encapsulates Operator smoothers.
void print(Teuchos::FancyOStream &out, const VerbLevel verbLevel=Default) const
Print the object with some verbosity level to an FancyOStream object.
void SetupRefMaxwell(Level ¤tLevel)
void Setup(Level ¤tLevel)
Set up the smoother.
void SetParameterList(const Teuchos::ParameterList ¶mList)
Set parameters from a parameter list and return with default values.
size_t getNodeSmootherComplexity() const
Get a rough estimate of cost per iteration.
void DeclareInput(Level ¤tLevel) const
Input.
RefMaxwellSmoother(const std::string type, const Teuchos::ParameterList ¶mList)
Constructor.
std::string description() const
Return a simple one-line description of this object.
void Apply(MultiVector &X, const MultiVector &B, bool InitialGuessIsZero=false) const
Apply the preconditioner.
RCP< SmootherPrototype > Copy() const
Preconditioner (wrapped as a Xpetra::Operator) for Maxwell's equations in curl-curl form.
void declareConstructionOutcome(bool fail, std::string msg)
bool IsSetup() const
Get the state of a smoother prototype.
Namespace for MueLu classes and methods.
@ Statistics1
Print more statistics.
@ External
Print external lib objects.
@ Parameters1
Print class parameters (more parameters, more verbose)