53 #ifndef MUELU_INDEFBLOCKEDDIAGONALSMOOTHER_DEF_HPP_
54 #define MUELU_INDEFBLOCKEDDIAGONALSMOOTHER_DEF_HPP_
56 #include "Teuchos_ArrayViewDecl.hpp"
62 #include <Xpetra_CrsMatrixWrap.hpp>
64 #include <Xpetra_MultiVectorFactory.hpp>
71 #include "MueLu_HierarchyUtils.hpp"
75 #include "MueLu_SchurComplementFactory.hpp"
76 #include "MueLu_FactoryManager.hpp"
80 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
82 : type_(
"IndefiniteBlockDiagonalSmoother")
86 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
89 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
93 validParamList->
set<
RCP<const FactoryBase> >(
"A", Teuchos::null,
"Generating factory of the matrix A (must be a 2x2 block matrix)");
94 validParamList->
set<
Scalar>(
"Damping factor", 1.0,
"Damping/Scaling factor");
95 validParamList->
set<
LocalOrdinal>(
"Sweeps", 1,
"Number of SIMPLE sweeps (default = 1)");
98 return validParamList;
102 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
106 size_t myPos = Teuchos::as<size_t>(pos);
108 if (myPos < FactManager_.size()) {
110 FactManager_.at(myPos) = FactManager;
111 }
else if (myPos == FactManager_.size()) {
113 FactManager_.push_back(FactManager);
116 *out <<
"Warning: cannot add new FactoryManager at proper position " << pos <<
". The FactoryManager is just appended to the end. Check this!" << std::endl;
119 FactManager_.push_back(FactManager);
123 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
125 currentLevel.
DeclareInput(
"A", this->GetFactory(
"A").
get());
127 TEUCHOS_TEST_FOR_EXCEPTION(FactManager_.size() != 2,
Exceptions::RuntimeError,
"MueLu::IndefBlockedDiagonalSmoother::DeclareInput: You have to declare two FactoryManagers with a \"Smoother\" object: One for predicting the primary variable and one for the SchurComplement system. The smoother for the SchurComplement system needs a SchurComplementFactory as input for variable \"A\"!");
130 std::vector<Teuchos::RCP<const FactoryManagerBase> >::const_iterator it;
131 for (it = FactManager_.begin(); it != FactManager_.end(); ++it) {
135 currentLevel.
DeclareInput(
"PreSmoother", (*it)->GetFactory(
"Smoother").get());
139 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
141 FactoryMonitor m(*
this,
"Setup for indefinite blocked diagonal smoother", currentLevel);
144 this->GetOStream(
Warnings0) <<
"MueLu::IndefBlockedDiagonalSmoother::Setup(): Setup() has already been called";
147 A_ = Factory::Get<RCP<Matrix> >(currentLevel,
"A");
153 rangeMapExtractor_ = bA->getRangeMapExtractor();
154 domainMapExtractor_ = bA->getDomainMapExtractor();
206 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
219 bool bRangeThyraMode = rangeMapExtractor_->getThyraMode();
220 bool bDomainThyraMode = domainMapExtractor_->getThyraMode();
227 bool bCopyResultX =
false;
228 bool bReorderX =
false;
229 bool bReorderB =
false;
239 if (rbA.is_null() ==
false) {
244 if (bX.is_null() ==
true) {
246 rcpX.
swap(vectorWithBlockedMap);
252 if (bB.is_null() ==
true) {
254 rcpB.
swap(vectorWithBlockedMap);
259 if (bX.is_null() ==
true) {
261 rcpX.
swap(vectorWithBlockedMap);
265 if (bB.is_null() ==
true) {
267 rcpB.
swap(vectorWithBlockedMap);
272 bX = Teuchos::rcp_dynamic_cast<BlockedMultiVector>(rcpX);
273 bB = Teuchos::rcp_dynamic_cast<
const BlockedMultiVector>(rcpB);
276 if (rbA.is_null() ==
false) {
280 if (bX->getBlockedMap()->getNumMaps() != bA->getDomainMapExtractor()->NumMaps() || bReorderX) {
283 rcpX.
swap(reorderedBlockedVector);
286 if (bB->getBlockedMap()->getNumMaps() != bA->getRangeMapExtractor()->NumMaps() || bReorderB) {
289 rcpB.
swap(reorderedBlockedVector);
303 RCP<MultiVector> xtilde = MultiVectorFactory::Build(rcpX->getMap(), rcpX->getNumVectors());
311 residual->update(one, *rcpB, zero);
322 bxtilde->putScalar(zero);
323 velPredictSmoo_->Apply(*xtilde1, *r1);
327 schurCompSmoo_->Apply(*xtilde2, *r2);
330 rcpX->update(omega, *bxtilde, one);
337 x1->update(omega, *xtilde1, one);
338 x2->update(omega, *xtilde2, one);
341 domainMapExtractor_->InsertVector(x1, 0, rcpX, bDomainThyraMode);
342 domainMapExtractor_->InsertVector(x2, 1, rcpX, bDomainThyraMode);
346 if (bCopyResultX ==
true) {
348 X.update(one, *Xmerged, zero);
352 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
358 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
360 std::ostringstream out;
362 out <<
"{type = " << type_ <<
"}";
366 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
371 out0 <<
"Prec. type: " << type_ << std::endl;
374 if (verbLevel &
Debug) {
378 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
Important warning messages (one line)
virtual const Teuchos::ParameterList & GetParameterList() const
Exception indicating invalid cast attempted.
std::string description() const
Return a simple one-line description of this object.
MueLu::DefaultLocalOrdinal LocalOrdinal
void Apply(MultiVector &X, const MultiVector &B, bool InitialGuessIsZero=false) const
Apply the Braess Sarazin smoother.
IndefBlockedDiagonalSmoother()
Constructor.
T & get(const std::string &name, T def_value)
ParameterList & set(std::string const &name, T const &value, std::string const &docString="", RCP< const ParameterEntryValidator > const &validator=null)
void Setup(Level ¤tLevel)
Setup routine.
Timer to be used in factories. Similar to Monitor but with additional timers.
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
void swap(RCP< T > &r_ptr)
Print additional debugging information.
void print(Teuchos::FancyOStream &out, const VerbLevel verbLevel=Default) const
Print the object with some verbosity level to an FancyOStream object.
void DeclareInput(Level ¤tLevel) const
Input.
virtual const RCP< const FactoryBase > GetFactory(const std::string &varName) const =0
Get.
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
MueLu::DefaultScalar Scalar
Class that holds all level-specific information.
Cheap Blocked diagonal smoother for indefinite 2x2 block matrices.
RCP< const ParameterList > GetValidParameterList() const
Input.
bool IsSetup() const
Get the state of a smoother prototype.
#define MUELU_DESCRIBE
Helper macro for implementing Describable::describe() for BaseClass objects.
size_t getNodeSmootherComplexity() const
Get a rough estimate of cost per iteration.
void AddFactoryManager(RCP< const FactoryManagerBase > FactManager, int pos=0)
Add a factory manager for BraessSarazin internal SchurComplement handling.
void residual(const Operator< SC, LO, GO, NO > &Aop, const MultiVector< SC, LO, GO, NO > &X_in, const MultiVector< SC, LO, GO, NO > &B_in, MultiVector< SC, LO, GO, NO > &R_in)
Exception throws to report errors in the internal logical of the program.
#define MUELU_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
An exception safe way to call the method 'Level::SetFactoryManager()'.
Teuchos::RCP< const Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > buildReorderedBlockedMultiVector(Teuchos::RCP< const Xpetra::BlockReorderManager > brm, Teuchos::RCP< const Xpetra::BlockedMultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > bvec)
virtual ~IndefBlockedDiagonalSmoother()
Destructor.
RCP< SmootherPrototype > Copy() const
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.
std::string toString(const T &t)