8 #ifndef PACKAGES_MUELU_ADAPTERS_XPETRA_MUELU_CREATEXPETRAPRECONDITIONER_HPP_ 9 #define PACKAGES_MUELU_ADAPTERS_XPETRA_MUELU_CREATEXPETRAPRECONDITIONER_HPP_ 11 #include <Teuchos_XMLParameterListHelpers.hpp> 12 #include <Xpetra_CrsMatrix.hpp> 13 #include <Xpetra_MultiVector.hpp> 14 #include <Xpetra_MultiVectorFactory.hpp> 19 #include <MueLu_Hierarchy.hpp> 21 #include <MueLu_MLParameterListInterpreter.hpp> 22 #include <MueLu_ParameterListInterpreter.hpp> 23 #include <MueLu_Utilities.hpp> 24 #include <MueLu_HierarchyUtils.hpp> 28 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
29 Teuchos::RCP<MueLu::Hierarchy<Scalar,LocalOrdinal,GlobalOrdinal,Node> >
31 const Teuchos::ParameterList& inParamList,
32 Teuchos::RCP<Xpetra::MultiVector<double, LocalOrdinal, GlobalOrdinal, Node> > coords = Teuchos::null,
33 Teuchos::RCP<Xpetra::MultiVector<Scalar, LocalOrdinal, GlobalOrdinal, Node> > nullspace = Teuchos::null) {
39 typedef Xpetra::MultiVectorFactory<Scalar,LocalOrdinal,GlobalOrdinal,Node> MultiVectorFactory;
41 std::string timerName =
"MueLu setup time";
42 RCP<Teuchos::Time> tm = Teuchos::TimeMonitor::getNewTimer(timerName);
45 ParameterList paramList = inParamList;
47 bool hasParamList = paramList.numParams();
49 RCP<HierarchyManager> mueLuFactory;
51 std::string syntaxStr =
"parameterlist: syntax";
52 if (hasParamList && paramList.isParameter(syntaxStr) && paramList.get<std::string>(syntaxStr) ==
"ml") {
53 paramList.remove(syntaxStr);
59 RCP<Hierarchy> H = mueLuFactory->CreateHierarchy();
60 H->setlib(op->getDomainMap()->lib());
64 H->GetLevel(0)->Set(
"A", op);
67 if (coords != Teuchos::null) {
68 H->GetLevel(0)->Set(
"Coordinates", coords);
72 if (nullspace == Teuchos::null) {
73 int nPDE = MueLu::MasterList::getDefault<int>(
"number of equations");
74 if (paramList.isSublist(
"Matrix")) {
76 const Teuchos::ParameterList& operatorList = paramList.sublist(
"Matrix");
77 if (operatorList.isParameter(
"PDE equations"))
78 nPDE = operatorList.get<
int>(
"PDE equations");
80 }
else if (paramList.isParameter(
"number of equations")) {
82 nPDE = paramList.get<
int>(
"number of equations");
85 nullspace = MultiVectorFactory::Build(op->getDomainMap(), nPDE);
87 nullspace->putScalar(Teuchos::ScalarTraits<Scalar>::one());
90 for (
int i = 0; i < nPDE; i++) {
91 Teuchos::ArrayRCP<Scalar> nsData = nullspace->getDataNonConst(i);
92 for (
int j = 0; j < nsData.size(); j++) {
93 GlobalOrdinal GID = op->getDomainMap()->getGlobalElement(j) - op->getDomainMap()->getIndexBase();
95 if ((GID-i) % nPDE == 0)
96 nsData[j] = Teuchos::ScalarTraits<Scalar>::one();
101 H->GetLevel(0)->Set(
"Nullspace", nullspace);
104 Teuchos::ParameterList nonSerialList,dummyList;
108 mueLuFactory->SetupHierarchy(*H);
111 tm->incrementNumCalls();
114 const bool alwaysWriteLocal =
true;
115 const bool writeGlobalStats =
true;
116 const bool writeZeroTimers =
false;
117 const bool ignoreZeroTimers =
true;
118 const std::string filter = timerName;
119 Teuchos::TimeMonitor::summarize(op->getRowMap()->getComm().ptr(), std::cout, alwaysWriteLocal, writeGlobalStats,
120 writeZeroTimers, Teuchos::Union, filter, ignoreZeroTimers);
135 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
138 std::string timerName =
"MueLu setup time";
139 RCP<Teuchos::Time> tm = Teuchos::TimeMonitor::getNewTimer(timerName);
143 typedef LocalOrdinal LO;
144 typedef GlobalOrdinal GO;
147 typedef Xpetra::Matrix<SC,LO,GO,NO> Matrix;
148 typedef Xpetra::Operator<SC,LO,GO,NO> Operator;
151 "MueLu::ReuseXpetraPreconditioner: Hierarchy has no levels in it");
153 "MueLu::ReuseXpetraPreconditioner: Hierarchy has no fine level operator");
154 RCP<Level> level0 = H->GetLevel(0);
156 RCP<Operator> O0 = level0->Get<RCP<Operator> >(
"A");
157 RCP<Matrix> A0 = Teuchos::rcp_dynamic_cast<Matrix>(O0);
163 A->SetFixedBlockSize(A0->GetFixedBlockSize());
170 tm->incrementNumCalls();
173 const bool alwaysWriteLocal =
true;
174 const bool writeGlobalStats =
true;
175 const bool writeZeroTimers =
false;
176 const bool ignoreZeroTimers =
true;
177 const std::string filter = timerName;
178 Teuchos::TimeMonitor::summarize(A->getRowMap()->getComm().ptr(), std::cout, alwaysWriteLocal, writeGlobalStats,
179 writeZeroTimers, Teuchos::Union, filter, ignoreZeroTimers);
void ReuseXpetraPreconditioner(const Teuchos::RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &A, Teuchos::RCP< MueLu::Hierarchy< Scalar, LocalOrdinal, GlobalOrdinal, Node >> &H)
Helper function to reuse an existing MueLu preconditioner.
static void AddNonSerializableDataToHierarchy(HierarchyManager &HM, Hierarchy &H, const ParameterList &nonSerialList)
Teuchos::RCP< MueLu::Hierarchy< Scalar, LocalOrdinal, GlobalOrdinal, Node > > CreateXpetraPreconditioner(Teuchos::RCP< Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > op, const Teuchos::ParameterList &inParamList, Teuchos::RCP< Xpetra::MultiVector< double, LocalOrdinal, GlobalOrdinal, Node > > coords=Teuchos::null, Teuchos::RCP< Xpetra::MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > > nullspace=Teuchos::null)
Namespace for MueLu classes and methods.
Print statistics that do not involve significant additional computation.
Class that accepts ML-style parameters and builds a MueLu preconditioner. This interpreter uses the s...
Exception throws to report errors in the internal logical of the program.
Provides methods to build a multigrid hierarchy and apply multigrid cycles.
long ExtractNonSerializableData(const Teuchos::ParameterList &inList, Teuchos::ParameterList &serialList, Teuchos::ParameterList &nonSerialList)