46 #ifndef MUELU_HIERARCHY_DECL_HPP 47 #define MUELU_HIERARCHY_DECL_HPP 50 #include <Teuchos_Ptr.hpp> 67 #include "MueLu_FactoryManager.hpp" 99 template <class Scalar = Xpetra::Operator<>::scalar_type,
104 #undef MUELU_HIERARCHY_SHORT 158 template<
class S2,
class LO2,
class GO2,
class N2>
234 void Clear(
int startLevel = 0);
255 ReturnType Iterate(
const MultiVector& B, MultiVector& X, ConvData conv = ConvData(),
256 bool InitialGuessIsZero =
false,
LO startLevel = 0);
267 void Write(
const LO &start=-1,
const LO &end=-1,
const std::string &suffix=
"");
319 template<
class Node2>
379 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
380 template<
typename Node2>
389 new_h->Levels_.resize(this->GetNumLevels());
390 new_h->maxCoarseSize_ = maxCoarseSize_;
391 new_h->implicitTranspose_ = implicitTranspose_;
392 new_h->isPreconditioner_ = isPreconditioner_;
393 new_h->isDumpingEnabled_ = isDumpingEnabled_;
394 new_h->dumpLevel_ = dumpLevel_;
395 new_h->dumpFile_ = dumpFile_;
401 for (
int i = 0; i < GetNumLevels(); i++) {
406 A = level->template Get<RCP<Operator> >(
"A");
407 cloneA = Xpetra::clone<Scalar, LocalOrdinal, GlobalOrdinal, Node, Node2>(*A, node2);
408 clonelevel->template Set<RCP<CloneOperator> >(
"A", cloneA);
411 R = level->template Get<RCP<Operator> >(
"R");
412 cloneR = Xpetra::clone<Scalar, LocalOrdinal, GlobalOrdinal, Node, Node2>(*R, node2);
413 clonelevel->template Set<RCP<CloneOperator> >(
"R", cloneR);
416 P = level->template Get<RCP<Operator> >(
"P");
417 cloneP = Xpetra::clone<Scalar, LocalOrdinal, GlobalOrdinal, Node, Node2>(*P, node2);
418 clonelevel->template Set<RCP<CloneOperator> >(
"P", cloneP);
421 Pre = level->template Get<RCP<SmootherBase> >(
"PreSmoother");
422 clonePre = MueLu::clone<Scalar, LocalOrdinal, GlobalOrdinal, Node, Node2> (Pre, cloneA, node2);
423 clonelevel->template Set<RCP<CloneSmoother> >(
"PreSmoother", clonePre);
426 Post = level->template Get<RCP<SmootherBase> >(
"PostSmoother");
427 clonePost = MueLu::clone<Scalar, LocalOrdinal, GlobalOrdinal, Node, Node2> (Post, cloneA, node2);
428 clonelevel->
template Set<RCP<CloneSmoother> >(
"PostSmoother", clonePost);
430 new_h->Levels_[i] = clonelevel;
438 #define MUELU_HIERARCHY_SHORT 439 #endif // MUELU_HIERARCHY_DECL_HPP
void IsPreconditioner(const bool flag)
This class specifies the default factory that should generate some data on a Level if the data does n...
RCP< Level > & GetLevel(const int levelID=0)
Retrieve a certain level from hierarchy.
Hierarchy()
Default constructor.
void AddLevel(const RCP< Level > &level)
Add a level at the end of the hierarchy.
void SetMaxCoarseSize(Xpetra::global_size_t maxCoarseSize)
void AddNewLevel()
Add a new level at the end of the hierarchy.
static Xpetra::global_size_t GetDefaultMaxCoarseSize()
void SetCycle(CycleType Cycle)
Supports VCYCLE and WCYCLE types.
void Clear(int startLevel=0)
Clear impermanent data from previous setup.
Namespace for MueLu classes and methods.
Teuchos::ScalarTraits< SC > STS
MagnitudeType rate_
Convergece rate.
static const NoFactory * get()
void ReplaceCoordinateMap(Level &level)
void AddKeepFlag(const std::string &ename, const FactoryBase *factory=NoFactory::get(), KeepType keep=MueLu::Keep)
Call Level::AddKeepFlag for each level of the Hierarchy.
void print(std::ostream &out=std::cout, const VerbLevel verbLevel=(MueLu::Parameters|MueLu::Statistics0)) const
Xpetra::UnderlyingLib lib()
GlobalOrdinal global_ordinal_type
Base class for factories (e.g., R, P, and A_coarse).
Print statistics that do not involve significant additional computation.
ConvData(std::pair< LO, MagnitudeType > p)
Xpetra::UnderlyingLib lib_
CycleType GetCycle() const
Returns multigrid cycle type (supports VCYCLE and WCYCLE)
static CycleType GetDefaultCycle()
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
void Write(const LO &start=-1, const LO &end=-1, const std::string &suffix="")
Print matrices in the multigrid hierarchy to file.
Class that holds all level-specific information.
LocalOrdinal local_ordinal_type
Class that provides default factories within Needs class.
void RemoveKeepFlag(const std::string &ename, const FactoryBase *factory, KeepType keep=MueLu::All)
Call Level::RemoveKeepFlag for each level of the Hierarchy.
int GetGlobalNumLevels() const
void CheckLevel(Level &level, int levelID)
Helper function.
const RCP< const FactoryManagerBase > & GetLevelManager(const int levelID) const
Always keep data, even accross run. This flag is set by Level::Keep(). This flag is propagated to coa...
ConvData(MagnitudeType tol)
static bool GetDefaultImplicitTranspose()
void describe(Teuchos::FancyOStream &out, const VerbLevel verbLevel=Default) const
Print the Hierarchy with some verbosity level to a FancyOStream object.
bool GetImplicitTranspose() const
void SetPRrebalance(bool doPRrebalance)
Base class for smoothers.
double GetOperatorComplexity() const
static int GetDefaultMaxLevels()
Base class for MueLu classes.
void SetImplicitTranspose(const bool &implicit)
virtual ~Hierarchy()
Destructor.
Data struct for defining stopping criteria of multigrid iteration.
STS::magnitudeType MagnitudeType
void Delete(const std::string &ename, const FactoryBase *factory=NoFactory::get())
Call Level::Delete(ename, factory) for each level of the Hierarchy.
MagnitudeType GetRate() const
Array< RCP< Level > > Levels_
Container for Level objects.
void DumpCurrentGraph() const
static bool GetDefaultPRrebalance()
bool IsAvailable(const std::string &ename, const FactoryBase *factory=NoFactory::get()) const
Test whether a need's value has been saved.
Xpetra::global_size_t GetMaxCoarseSize() const
ReturnType Iterate(const MultiVector &B, MultiVector &X, ConvData conv=ConvData(), bool InitialGuessIsZero=false, LO startLevel=0)
Apply the multigrid preconditioner.
void Keep(const std::string &ename, const FactoryBase *factory=NoFactory::get())
Call Level::Keep(ename, factory) for each level of the Hierarchy.
bool isDumpingEnabled_
Graph dumping.
void setlib(Xpetra::UnderlyingLib inlib)
Teuchos::RCP< Hierarchy< Scalar, LocalOrdinal, GlobalOrdinal, Node2 > > clone(const RCP< Node2 > &node2) const
std::string description() const
Return a simple one-line description of this object.
bool Setup(int coarseLevelID, const RCP< const FactoryManagerBase > fineLevelManager, const RCP< const FactoryManagerBase > coarseLevelManager, const RCP< const FactoryManagerBase > nextLevelManager=Teuchos::null)
Multi-level setup phase: build a new level of the hierarchy.
void EnableGraphDumping(const std::string &filename, int levelID=1)
Provides methods to build a multigrid hierarchy and apply multigrid cycles.
Array< RCP< const FactoryManagerBase > > levelManagers_
Xpetra::global_size_t maxCoarseSize_