53 #ifndef MUELU_PERMUTINGSMOOTHER_DEF_HPP 54 #define MUELU_PERMUTINGSMOOTHER_DEF_HPP 56 #include <Xpetra_Matrix.hpp> 57 #include <Xpetra_Vector.hpp> 58 #include <Xpetra_VectorFactory.hpp> 59 #include <Xpetra_MultiVector.hpp> 60 #include <Xpetra_MultiVectorFactory.hpp> 66 #include "MueLu_TrilinosSmoother.hpp" 68 #include "MueLu_Utilities.hpp" 69 #include "MueLu_PermutationFactory.hpp" 74 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
76 : type_(type), overlap_(overlap), permQT_(
Teuchos::null), permP_(
Teuchos::null), diagScalingOp_(
Teuchos::null)
83 newPermFact->SetParameter(
"PermutationRowMapName", Teuchos::ParameterEntry(mapName));
84 newPermFact->SetFactory (
"PermutationRowMapFactory", mapFact);
90 #if defined(HAVE_MUELU_EPETRA) && defined(HAVE_MUELU_IFPACK) 104 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
107 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
111 currentLevel.
DeclareInput(
"permScaling", permFact_.get());
113 s_->DeclareInput(currentLevel);
116 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
118 FactoryMonitor monitor(*
this,
"Permuting Smoother", currentLevel);
121 this->GetOStream(
Warnings0) <<
"MueLu::PermutingSmoother::Setup(): Setup() has already been called" << std::endl;
124 permP_ = currentLevel.
Get< RCP<Matrix> > (
"permP", permFact_.get());
125 permQT_ = currentLevel.
Get< RCP<Matrix> > (
"permQT", permFact_.get());
126 diagScalingOp_ = currentLevel.
Get< RCP<Matrix> > (
"permScaling", permFact_.get());
128 s_->Setup(currentLevel);
133 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
137 typedef Teuchos::ScalarTraits<Scalar> STS;
139 Teuchos::RCP<MultiVector> Xtemp = MultiVectorFactory::Build(X.getMap(), 1,
true);
140 Xtemp->update(STS::one(), X, STS::zero());
143 Teuchos::RCP<MultiVector> Btemp = MultiVectorFactory::Build(B.getMap(), 1,
true);
144 Teuchos::RCP<MultiVector> Btemp2 = MultiVectorFactory::Build(B.getMap(), 1,
true);
145 permP_->apply(B, *Btemp, Teuchos::NO_TRANS);
146 diagScalingOp_->apply(*Btemp, *Btemp2, Teuchos::NO_TRANS);
149 s_->Apply(*Xtemp, *Btemp2, InitialGuessIsZero);
152 permQT_->apply(*Xtemp, X, Teuchos::NO_TRANS);
155 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
156 RCP<MueLu::SmootherPrototype<Scalar, LocalOrdinal, GlobalOrdinal, Node> >
161 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
163 std::ostringstream out;
168 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
Important warning messages (one line)
T & Get(const std::string &ename, const FactoryBase *factory=NoFactory::get())
Get data without decrementing associated storage counter (i.e., read-only access). Usage: Level->Get< RCP<Matrix> >("A", factory) if factory == NULL => use default factory.
RCP< SmootherPrototype > Copy() const
virtual ~PermutingSmoother()
Destructor.
bool IsSetup() const
Get the state of a smoother prototype.
Class that encapsulates external library smoothers.
Timer to be used in factories. Similar to Monitor but with additional timers.
void DeclareInput(Level ¤tLevel) const
Input.
factory generates a row- and column permutation operators P and Q such that P*A*Q^T is a (hopefully) ...
Namespace for MueLu classes and methods.
void print(Teuchos::FancyOStream &out, const VerbLevel verbLevel=Default) const
Print the object with some verbosity level to an FancyOStream object.
This class first calculates row- and column permutation operators and applies a smoother to the permu...
LO overlap_
overlap when using the smoother in additive Schwarz mode
void Apply(MultiVector &X, const MultiVector &B, bool InitialGuessIsZero=false) const
Apply the direct solver. Solves the linear system AX=B using the constructed solver.
void Setup(Level ¤tLevel)
Set up the direct solver.
virtual void SetParameterList(const ParameterList ¶mList)
Set parameters from a parameter list and return with default values.
Class that holds all level-specific information.
RCP< FactoryBase > permFact_
Permutation Factory.
std::string type_
ifpack1/2-specific key phrase that denote smoother type
#define MUELU_DESCRIBE
Helper macro for implementing Describable::describe() for BaseClass objects.
std::string description() const
Return a simple one-line description of this object.
virtual const Teuchos::ParameterList & GetParameterList() const
PermutingSmoother(std::string const &mapName, const RCP< const FactoryBase > &mapFact, std::string const &type="", const Teuchos::ParameterList ¶mList=Teuchos::ParameterList(), LO const &overlap=0, RCP< FactoryBase > permFact=Teuchos::null)
Constructor.
Exception throws to report errors in the internal logical of the program.
RCP< SmootherPrototype > s_
Smoother.
void DeclareInput(const std::string &ename, const FactoryBase *factory, const FactoryBase *requestedBy=NoFactory::get())
Callback from FactoryBase::CallDeclareInput() and FactoryBase::DeclareInput()
virtual std::string description() const
Return a simple one-line description of this object.