46 #ifndef MUELU_HIERARCHY_DECL_HPP
47 #define MUELU_HIERARCHY_DECL_HPP
50 #include <Teuchos_Ptr.hpp>
67 #include "MueLu_FactoryManager.hpp"
100 template <class Scalar = Xpetra::Operator<>::scalar_type,
105 #undef MUELU_HIERARCHY_SHORT
164 template<
class S2,
class LO2,
class GO2,
class N2>
243 void Clear(
int startLevel = 0);
267 ReturnType Iterate(
const MultiVector& B, MultiVector& X, ConvData conv = ConvData(),
268 bool InitialGuessIsZero =
false,
LO startLevel = 0);
279 void Write(
const LO &start=-1,
const LO &end=-1,
const std::string &suffix=
"");
331 #ifdef HAVE_MUELU_DEPRECATED_CODE
332 template<
class Node2>
413 #ifdef HAVE_MUELU_DEPRECATED_CODE
414 template<
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
415 template<
typename Node2>
425 new_h->Levels_.resize(this->GetNumLevels());
426 new_h->maxCoarseSize_ = maxCoarseSize_;
427 new_h->implicitTranspose_ = implicitTranspose_;
428 new_h->isPreconditioner_ = isPreconditioner_;
429 new_h->isDumpingEnabled_ = isDumpingEnabled_;
430 new_h->dumpLevel_ = dumpLevel_;
431 new_h->dumpFile_ = dumpFile_;
437 for (
int i = 0; i < GetNumLevels(); i++) {
442 A = level->template Get<RCP<Operator> >(
"A");
443 cloneA = Xpetra::clone<Scalar, LocalOrdinal, GlobalOrdinal, Node, Node2>(*A, node2);
444 clonelevel->template Set<RCP<CloneOperator> >(
"A", cloneA);
447 R = level->template Get<RCP<Operator> >(
"R");
448 cloneR = Xpetra::clone<Scalar, LocalOrdinal, GlobalOrdinal, Node, Node2>(*R, node2);
449 clonelevel->template Set<RCP<CloneOperator> >(
"R", cloneR);
452 P = level->template Get<RCP<Operator> >(
"P");
453 cloneP = Xpetra::clone<Scalar, LocalOrdinal, GlobalOrdinal, Node, Node2>(*P, node2);
454 clonelevel->template Set<RCP<CloneOperator> >(
"P", cloneP);
457 Pre = level->template Get<RCP<SmootherBase> >(
"PreSmoother");
458 clonePre = MueLu::clone<Scalar, LocalOrdinal, GlobalOrdinal, Node, Node2> (Pre, cloneA, node2);
459 clonelevel->template Set<RCP<CloneSmoother> >(
"PreSmoother", clonePre);
462 Post = level->template Get<RCP<SmootherBase> >(
"PostSmoother");
463 clonePost = MueLu::clone<Scalar, LocalOrdinal, GlobalOrdinal, Node, Node2> (Post, cloneA, node2);
464 clonelevel->
template Set<RCP<CloneSmoother> >(
"PostSmoother", clonePost);
466 new_h->Levels_[i] = clonelevel;
475 #define MUELU_HIERARCHY_SHORT
476 #endif // MUELU_HIERARCHY_DECL_HPP
void IsPreconditioner(const bool flag)
RCP< Map< LocalOrdinal, GlobalOrdinal, Node2 > > clone(const Map< LocalOrdinal, GlobalOrdinal, Node1 > &map, const RCP< Node2 > &node2)
Array< RCP< MultiVector > > coarseExport_
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 DeleteLevelMultiVectors()
Array< RCP< MultiVector > > coarseImport_
void AddLevel(const RCP< Level > &level)
Add a level at the end of the hierarchy.
int sizeOfAllocatedLevelMultiVectors_
void SetMaxCoarseSize(Xpetra::global_size_t maxCoarseSize)
Xpetra::global_size_t GetMaxCoarseSize() const
void AddNewLevel()
Add a new level at the end of the hierarchy.
static Xpetra::global_size_t GetDefaultMaxCoarseSize()
void DumpCurrentGraph() const
void SetCycle(CycleType Cycle)
Supports VCYCLE and WCYCLE types.
void Clear(int startLevel=0)
Clear impermanent data from previous setup.
bool GetImplicitTranspose() const
const RCP< const FactoryManagerBase > & GetLevelManager(const int levelID) const
Teuchos::ScalarTraits< SC > STS
MagnitudeType rate_
Convergece rate.
static const NoFactory * get()
void ReplaceCoordinateMap(Level &level)
void AllocateLevelMultiVectors(int numvecs)
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()
double GetSmootherComplexity() const
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_
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.
std::string description() const
Return a simple one-line description of this object.
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.
void CheckLevel(Level &level, int levelID)
Helper function.
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 SetPRrebalance(bool doPRrebalance)
Base class for smoothers.
int GetGlobalNumLevels() const
static int GetDefaultMaxLevels()
CycleType GetCycle() const
Returns multigrid cycle type (supports VCYCLE and WCYCLE)
Base class for MueLu classes.
void SetImplicitTranspose(const bool &implicit)
virtual ~Hierarchy()
Destructor.
Data struct for defining stopping criteria of multigrid iteration.
Array< RCP< MultiVector > > coarseRhs_
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.
Array< RCP< Level > > Levels_
Container for Level objects.
static bool GetDefaultPRrebalance()
MagnitudeType GetRate() const
void SetProlongatorScalingFactor(double scalingFactor)
Specify damping factor alpha such that x = x + alpha*P*c, where c is the coarse grid correction...
void describe(Teuchos::FancyOStream &out, const VerbLevel verbLevel=Default) const
Print the Hierarchy with some verbosity level to a FancyOStream object.
ReturnType Iterate(const MultiVector &B, MultiVector &X, ConvData conv=ConvData(), bool InitialGuessIsZero=false, LO startLevel=0)
Apply the multigrid preconditioner.
Array< RCP< MultiVector > > residual_
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)
double GetOperatorComplexity() const
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.
Array< RCP< MultiVector > > coarseX_
Array< RCP< MultiVector > > correction_
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_
bool IsAvailable(const std::string &ename, const FactoryBase *factory=NoFactory::get()) const
Test whether a need's value has been saved.