46 #ifndef MUELU_SAPFACTORY_DEF_HPP
47 #define MUELU_SAPFACTORY_DEF_HPP
59 #include "MueLu_PerfUtils.hpp"
61 #include "MueLu_TentativePFactory.hpp"
62 #include "MueLu_Utilities.hpp"
66 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
70 #define SET_VALID_ENTRY(name) validParamList->setEntry(name, MasterList::getEntry(name))
74 #undef SET_VALID_ENTRY
76 validParamList->
set<
RCP<const FactoryBase> >(
"A", Teuchos::null,
"Generating factory of the matrix A used during the prolongator smoothing process");
81 norecurse.disableRecursiveValidation();
82 validParamList->
set<
ParameterList> (
"matrixmatrix: kernel params", norecurse,
"MatrixMatrix kernel parameters");
85 return validParamList;
88 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
90 Input(fineLevel,
"A");
95 if (initialPFact == Teuchos::null) { initialPFact = coarseLevel.
GetFactoryManager()->GetFactory(
"Ptent"); }
99 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
101 return BuildP(fineLevel, coarseLevel);
104 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
110 const std::string prefix =
"MueLu::SaPFactory(" + levelIDs +
"): ";
118 if (initialPFact == Teuchos::null) { initialPFact = coarseLevel.
GetFactoryManager()->GetFactory(
"Ptent"); }
122 RCP<Matrix> A = Get< RCP<Matrix> >(fineLevel,
"A");
125 if (restrictionMode_) {
135 if(pL.
isSublist(
"matrixmatrix: kernel params"))
136 APparams->
sublist(
"matrixmatrix: kernel params") = pL.
sublist(
"matrixmatrix: kernel params");
138 if (coarseLevel.
IsAvailable(
"AP reuse data",
this)) {
139 GetOStream(static_cast<MsgType>(
Runtime0 |
Test)) <<
"Reusing previous AP data" << std::endl;
147 APparams->
set(
"compute global constants: temporaries",APparams->
get(
"compute global constants: temporaries",
false));
148 APparams->
set(
"compute global constants",APparams->
get(
"compute global constants",
false));
150 SC dampingFactor = as<SC>(pL.
get<
double>(
"sa: damping factor"));
151 LO maxEigenIterations = as<LO>(pL.
get<
int> (
"sa: eigenvalue estimate num iterations"));
152 bool estimateMaxEigen = pL.
get<
bool> (
"sa: calculate eigenvalue estimate");
158 lambdaMax = A->GetMaxEigenvalueEstimate();
160 GetOStream(
Statistics1) <<
"Calculating max eigenvalue estimate now (max iters = "<< maxEigenIterations <<
")" << std::endl;
161 Coordinate stopTol = 1e-4;
163 A->SetMaxEigenvalueEstimate(lambdaMax);
165 GetOStream(
Statistics1) <<
"Using cached max eigenvalue estimate" << std::endl;
167 GetOStream(
Statistics1) <<
"Prolongator damping factor = " << dampingFactor/lambdaMax <<
" (" << dampingFactor <<
" / " << lambdaMax <<
")" << std::endl;
174 SC omega = dampingFactor / lambdaMax;
179 GetOStream(
Statistics2), std::string(
"MueLu::SaP-")+levelIDs, APparams);
187 if (!restrictionMode_) {
189 if(!finalP.
is_null()) {std::ostringstream oss; oss <<
"P_" << coarseLevel.
GetLevelID(); finalP->setObjectLabel(oss.str());}
190 Set(coarseLevel,
"P", finalP);
192 APparams->
set(
"graph", finalP);
193 Set(coarseLevel,
"AP reuse data", APparams);
196 if (Ptent->IsView(
"stridedMaps"))
197 finalP->CreateView(
"stridedMaps", Ptent);
201 params->
set(
"printLoadBalancingInfo",
true);
202 params->
set(
"printCommInfo",
true);
212 if(!R.
is_null()) {std::ostringstream oss; oss <<
"R_" << coarseLevel.
GetLevelID(); R->setObjectLabel(oss.str());}
215 Set(coarseLevel,
"R", R);
218 if (Ptent->IsView(
"stridedMaps"))
219 R->CreateView(
"stridedMaps", Ptent,
true);
223 params->
set(
"printLoadBalancingInfo",
true);
224 params->
set(
"printCommInfo",
true);
233 #endif // MUELU_SAPFACTORY_DEF_HPP
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.
std::string toString(const T &what)
Little helper function to convert non-string types to strings.
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)
Timer to be used in factories. Similar to Monitor but with additional timers.
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
One-liner description of what is happening.
static RCP< Xpetra::Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > GetMatrixDiagonalInverse(const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, Magnitude tol=Teuchos::ScalarTraits< Scalar >::eps()*100)
Print even more statistics.
bool isParameter(const std::string &name) const
static RCP< Matrix > Jacobi(SC omega, const Vector &Dinv, const Matrix &A, const Matrix &B, RCP< Matrix > C_in, Teuchos::FancyOStream &fos, const std::string &label, RCP< ParameterList > ¶ms)
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
Class that holds all level-specific information.
bool isSublist(const std::string &name) const
Timer to be used in factories. Similar to SubMonitor but adds a timer level by level.
void Build(Level &fineLevel, Level &coarseLevel) const
Build method.
#define SET_VALID_ENTRY(name)
static RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > Transpose(Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Op, bool optimizeTranspose=false, const std::string &label=std::string(), const Teuchos::RCP< Teuchos::ParameterList > ¶ms=Teuchos::null)
static std::string PrintMatrixInfo(const Matrix &A, const std::string &msgTag, RCP< const Teuchos::ParameterList > params=Teuchos::null)
void BuildP(Level &fineLevel, Level &coarseLevel) const
Abstract Build method.
const RCP< const FactoryManagerBase > GetFactoryManager()
returns the current factory manager
ParameterList & sublist(const std::string &name, bool mustAlreadyExist=false, const std::string &docString="")
int GetLevelID() const
Return level number.
Exception throws to report errors in the internal logical of the program.
void DeclareInput(Level &fineLevel, Level &coarseLevel) const
Input.
void DeclareInput(const std::string &ename, const FactoryBase *factory, const FactoryBase *requestedBy=NoFactory::get())
Callback from FactoryBase::CallDeclareInput() and FactoryBase::DeclareInput()
RCP< const ParameterList > GetValidParameterList() const
Return a const parameter list of valid parameters that setParameterList() will accept.
static Scalar PowerMethod(const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, bool scaleByDiag=true, LocalOrdinal niters=10, Magnitude tolerance=1e-2, bool verbose=false, unsigned int seed=123)
bool IsAvailable(const std::string &ename, const FactoryBase *factory=NoFactory::get()) const
Test whether a need's value has been saved.