10 #ifndef MUELU_RAPFACTORY_DEF_HPP
11 #define MUELU_RAPFACTORY_DEF_HPP
25 #include "MueLu_PerfUtils.hpp"
29 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
31 : hasDeclaredInput_(false) {}
33 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
36 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
40 #define SET_VALID_ENTRY(name) validParamList->setEntry(name, MasterList::getEntry(name))
47 #undef SET_VALID_ENTRY
48 validParamList->
set<
RCP<const FactoryBase> >(
"A", null,
"Generating factory of the matrix A used during the prolongator smoothing process");
52 validParamList->
set<
bool>(
"CheckMainDiagonal",
false,
"Check main diagonal for zeros");
53 validParamList->
set<
bool>(
"RepairMainDiagonal",
false,
"Repair zeros on main diagonal");
58 validParamList->
set<
ParameterList>(
"matrixmatrix: kernel params", norecurse,
"MatrixMatrix kernel parameters");
60 return validParamList;
63 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
66 if (pL.
get<
bool>(
"transpose: use implicit") ==
false)
67 Input(coarseLevel,
"R");
69 Input(fineLevel,
"A");
70 Input(coarseLevel,
"P");
73 for (std::vector<
RCP<const FactoryBase> >::const_iterator it = transferFacts_.begin(); it != transferFacts_.end(); ++it)
74 (*it)->CallDeclareInput(coarseLevel);
76 hasDeclaredInput_ =
true;
79 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
81 const bool doTranspose =
true;
82 const bool doFillComplete =
true;
83 const bool doOptimizeStorage =
true;
87 std::ostringstream levelstr;
92 "MueLu::RAPFactory::Build(): CallDeclareInput has not been called before Build!");
96 RCP<Matrix> P = Get<RCP<Matrix> >(coarseLevel,
"P"), AP;
99 if (P == Teuchos::null) {
101 Set(coarseLevel,
"A", Ac);
106 bool isGPU = Node::is_gpu;
108 if (pL.
get<
bool>(
"rap: triple product") ==
false || isEpetra || isGPU) {
109 if (pL.
get<
bool>(
"rap: triple product") && isEpetra)
110 GetOStream(
Warnings1) <<
"Switching from triple product to R x (A x P) since triple product has not been implemented for Epetra.\n";
111 if (pL.
get<
bool>(
"rap: triple product") && isGPU)
112 GetOStream(
Warnings1) <<
"Switching from triple product to R x (A x P) since triple product has not been implemented for "
113 << Node::execution_space::name() << std::endl;
117 if (pL.
isSublist(
"matrixmatrix: kernel params"))
121 APparams->
set(
"compute global constants: temporaries", APparams->
get(
"compute global constants: temporaries",
false));
122 APparams->
set(
"compute global constants", APparams->
get(
"compute global constants",
false));
124 if (coarseLevel.IsAvailable(
"AP reuse data",
this)) {
125 GetOStream(static_cast<MsgType>(
Runtime0 |
Test)) <<
"Reusing previous AP data" << std::endl;
136 AP = MatrixMatrix::Multiply(*A, !doTranspose, *P, !doTranspose, AP, GetOStream(
Statistics2),
137 doFillComplete, doOptimizeStorage, labelstr + std::string(
"MueLu::A*P-") + levelstr.str(), APparams);
142 if (pL.
isSublist(
"matrixmatrix: kernel params"))
145 if (coarseLevel.IsAvailable(
"RAP reuse data",
this)) {
146 GetOStream(static_cast<MsgType>(
Runtime0 |
Test)) <<
"Reusing previous RAP data" << std::endl;
159 RAPparams->
set(
"compute global constants: temporaries", RAPparams->
get(
"compute global constants: temporaries",
false));
160 RAPparams->
set(
"compute global constants",
true);
166 if (pL.
get<
bool>(
"transpose: use implicit") ==
true) {
169 Ac = MatrixMatrix::Multiply(*P, doTranspose, *AP, !doTranspose, Ac, GetOStream(
Statistics2),
170 doFillComplete, doOptimizeStorage, labelstr + std::string(
"MueLu::R*(AP)-implicit-") + levelstr.str(), RAPparams);
173 RCP<Matrix> R = Get<RCP<Matrix> >(coarseLevel,
"R");
177 Ac = MatrixMatrix::Multiply(*R, !doTranspose, *AP, !doTranspose, Ac, GetOStream(
Statistics2),
178 doFillComplete, doOptimizeStorage, labelstr + std::string(
"MueLu::R*(AP)-explicit-") + levelstr.str(), RAPparams);
182 if (relativeFloor.
size() > 0) {
186 bool repairZeroDiagonals = pL.
get<
bool>(
"RepairMainDiagonal") || pL.
get<
bool>(
"rap: fix zero diagonals");
187 bool checkAc = pL.
get<
bool>(
"CheckMainDiagonal") || pL.
get<
bool>(
"rap: fix zero diagonals");
189 if (checkAc || repairZeroDiagonals) {
191 magnitudeType threshold;
192 if (pL.
isType<magnitudeType>(
"rap: fix zero diagonals threshold"))
193 threshold = pL.
get<magnitudeType>(
"rap: fix zero diagonals threshold");
195 threshold = Teuchos::as<magnitudeType>(pL.
get<
double>(
"rap: fix zero diagonals threshold"));
196 Scalar replacement = Teuchos::as<Scalar>(pL.
get<
double>(
"rap: fix zero diagonals replacement"));
203 params->
set(
"printLoadBalancingInfo",
true);
204 params->
set(
"printCommInfo",
true);
209 std::ostringstream oss;
210 oss <<
"A_" << coarseLevel.GetLevelID();
211 Ac->setObjectLabel(oss.str());
213 Set(coarseLevel,
"A", Ac);
216 APparams->
set(
"graph", AP);
217 Set(coarseLevel,
"AP reuse data", APparams);
220 RAPparams->
set(
"graph", Ac);
221 Set(coarseLevel,
"RAP reuse data", RAPparams);
225 if (pL.
isSublist(
"matrixmatrix: kernel params"))
226 RAPparams->
sublist(
"matrixmatrix: kernel params") = pL.
sublist(
"matrixmatrix: kernel params");
228 if (coarseLevel.IsAvailable(
"RAP reuse data",
this)) {
229 GetOStream(static_cast<MsgType>(
Runtime0 |
Test)) <<
"Reusing previous RAP data" << std::endl;
242 RAPparams->
set(
"compute global constants: temporaries", RAPparams->
get(
"compute global constants: temporaries",
false));
243 RAPparams->
set(
"compute global constants",
true);
245 if (pL.
get<
bool>(
"transpose: use implicit") ==
true) {
246 Ac = MatrixFactory::Build(P->getDomainMap(), Teuchos::as<LO>(0));
251 MultiplyRAP(*P, doTranspose, *A, !doTranspose, *P, !doTranspose, *Ac, doFillComplete,
252 doOptimizeStorage, labelstr + std::string(
"MueLu::R*A*P-implicit-") + levelstr.str(),
255 RCP<Matrix> R = Get<RCP<Matrix> >(coarseLevel,
"R");
256 Ac = MatrixFactory::Build(R->getRowMap(), Teuchos::as<LO>(0));
261 MultiplyRAP(*R, !doTranspose, *A, !doTranspose, *P, !doTranspose, *Ac, doFillComplete,
262 doOptimizeStorage, labelstr + std::string(
"MueLu::R*A*P-explicit-") + levelstr.str(),
267 if (relativeFloor.
size() > 0) {
271 bool repairZeroDiagonals = pL.
get<
bool>(
"RepairMainDiagonal") || pL.
get<
bool>(
"rap: fix zero diagonals");
272 bool checkAc = pL.
get<
bool>(
"CheckMainDiagonal") || pL.
get<
bool>(
"rap: fix zero diagonals");
274 if (checkAc || repairZeroDiagonals) {
276 magnitudeType threshold;
277 if (pL.
isType<magnitudeType>(
"rap: fix zero diagonals threshold"))
278 threshold = pL.
get<magnitudeType>(
"rap: fix zero diagonals threshold");
280 threshold = Teuchos::as<magnitudeType>(pL.
get<
double>(
"rap: fix zero diagonals threshold"));
281 Scalar replacement = Teuchos::as<Scalar>(pL.
get<
double>(
"rap: fix zero diagonals replacement"));
288 params->
set(
"printLoadBalancingInfo",
true);
289 params->
set(
"printCommInfo",
true);
294 std::ostringstream oss;
295 oss <<
"A_" << coarseLevel.GetLevelID();
296 Ac->setObjectLabel(oss.str());
298 Set(coarseLevel,
"A", Ac);
301 RAPparams->
set(
"graph", Ac);
302 Set(coarseLevel,
"RAP reuse data", RAPparams);
307 #ifdef HAVE_MUELU_DEBUG
308 MatrixUtils::checkLocalRowMapMatchesColMap(*Ac);
309 #endif // HAVE_MUELU_DEBUG
311 if (transferFacts_.begin() != transferFacts_.end()) {
315 for (std::vector<
RCP<const FactoryBase> >::const_iterator it = transferFacts_.begin(); it != transferFacts_.end(); ++it) {
317 GetOStream(
Runtime0) <<
"RAPFactory: call transfer factory: " << fac->
description() << std::endl;
352 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
356 "MueLu::RAPFactory::AddTransferFactory: Transfer factory is not derived from TwoLevelFactoryBase. "
357 "This is very strange. (Note: you can remove this exception if there's a good reason for)");
359 transferFacts_.push_back(factory);
364 #define MUELU_RAPFACTORY_SHORT
365 #endif // MUELU_RAPFACTORY_DEF_HPP
#define SET_VALID_ENTRY(name)
Exception indicating invalid cast attempted.
virtual std::string getObjectLabel() const
static void RelativeDiagonalBoost(RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &A, const Teuchos::ArrayView< const double > &relativeThreshold, Teuchos::FancyOStream &fos)
virtual void CallBuild(Level &requestedLevel) const =0
ParameterList & disableRecursiveValidation()
void Release(const FactoryBase &factory)
Decrement the storage counter for all the inputs of a factory.
T & get(const std::string &name, T def_value)
static void MultiplyRAP(const Matrix &R, bool transposeR, const Matrix &A, bool transposeA, const Matrix &P, bool transposeP, Matrix &Ac, bool call_FillComplete_on_result=true, bool doOptimizeStorage=true, const std::string &label=std::string(), const RCP< ParameterList > ¶ms=null)
Timer to be used in factories. Similar to Monitor but with additional timers.
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
RCP< const ParameterList > GetValidParameterList() const
Return a const parameter list of valid parameters that setParameterList() will accept.
One-liner description of what is happening.
void AddTransferFactory(const RCP< const FactoryBase > &factory)
Add transfer factory in the end of list of transfer factories in RepartitionAcFactory.
ParameterList & set(std::string const &name, T &&value, std::string const &docString="", RCP< const ParameterEntryValidator > const &validator=null)
static void CheckRepairMainDiagonal(RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >> &Ac, bool const &repairZeroDiagonals, Teuchos::FancyOStream &fos, const typename Teuchos::ScalarTraits< Scalar >::magnitudeType threshold=Teuchos::ScalarTraits< typename Teuchos::ScalarTraits< Scalar >::magnitudeType >::zero())
Print even more statistics.
bool isParameter(const std::string &name) const
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
MueLu::DefaultScalar Scalar
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.
static std::string PrintMatrixInfo(const Matrix &A, const std::string &msgTag, RCP< const Teuchos::ParameterList > params=Teuchos::null)
bool isType(const std::string &name) const
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.
virtual std::string description() const
Return a simple one-line description of this object.
void Build(Level &fineLevel, Level &coarseLevel) const
Build an object with this factory.
void DeclareInput(Level &fineLevel, Level &coarseLevel) const
Input.