46 #ifndef MUELU_USERPFACTORY_DEF_HPP 47 #define MUELU_USERPFACTORY_DEF_HPP 49 #include <Xpetra_MultiVectorFactory.hpp> 50 #include <Xpetra_Matrix.hpp> 51 #include <Xpetra_IO.hpp> 54 #include "MueLu_PerfUtils.hpp" 55 #include "MueLu_UserPFactory.hpp" 56 #include "MueLu_Utilities.hpp" 60 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
62 RCP<ParameterList> validParamList = rcp(
new ParameterList());
64 validParamList->set< RCP<const FactoryBase> >(
"A", Teuchos::null,
"Generating factory of the matrix A");
65 validParamList->set< RCP<const FactoryBase> >(
"Nullspace", Teuchos::null,
"Generating factory of the nullspace");
66 validParamList->set< std::string > (
"matrixFileName",
"",
"File with matrix data");
67 validParamList->set< std::string > (
"mapFileName",
"",
"File with coarse map data");
69 return validParamList;
72 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
74 Input(fineLevel,
"A");
75 Input(fineLevel,
"Nullspace");
78 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
80 return BuildP(fineLevel, coarseLevel);
83 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
87 RCP<Matrix> A = Get< RCP<Matrix> > (fineLevel,
"A");
88 RCP<MultiVector> fineNullspace = Get< RCP<MultiVector> > (fineLevel,
"Nullspace");
90 TEUCHOS_TEST_FOR_EXCEPTION(A->GetFixedBlockSize() != 1,
Exceptions::RuntimeError,
"Block size > 1 has not been implemented");
92 const Teuchos::ParameterList& pL = GetParameterList();
94 std::string mapFile = pL.get<std::string>(
"mapFileName");
95 RCP<const Map> rowMap = A->getRowMap();
96 RCP<const Map> coarseMap = Xpetra::IO<Scalar, LocalOrdinal, GlobalOrdinal, Node>::ReadMap(mapFile, rowMap->lib(), rowMap->getComm());
97 Set(coarseLevel,
"CoarseMap", coarseMap);
99 std::string matrixFile = pL.get<std::string>(
"matrixFileName");
100 RCP<Matrix> P = Xpetra::IO<Scalar, LocalOrdinal, GlobalOrdinal, Node>::Read(matrixFile, rowMap, coarseMap, coarseMap, rowMap);
102 Set(coarseLevel,
"P", P);
105 RCP<Matrix> P1 = Xpetra::MatrixMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>::Multiply(*A,
false, *P,
false);
106 P = Xpetra::IO<Scalar, LocalOrdinal, GlobalOrdinal, Node>::Read(matrixFile, rowMap, P1->getColMap(), coarseMap, rowMap);
107 Set(coarseLevel,
"P", P);
110 RCP<MultiVector> coarseNullspace = MultiVectorFactory::Build(coarseMap, fineNullspace->getNumVectors());
111 P->apply(*fineNullspace, *coarseNullspace, Teuchos::TRANS, Teuchos::ScalarTraits<SC>::one(), Teuchos::ScalarTraits<SC>::zero());
112 Set(coarseLevel,
"Nullspace", coarseNullspace);
115 size_t n = Teuchos::as<size_t>(sqrt(coarseMap->getGlobalNumElements()));
116 TEUCHOS_TEST_FOR_EXCEPTION(n*n != coarseMap->getGlobalNumElements(),
Exceptions::RuntimeError,
"Unfortunately, this is not the case, don't know what to do");
118 RCP<MultiVector> coarseCoords = MultiVectorFactory::Build(coarseMap, 2);
119 ArrayRCP<Scalar> x = coarseCoords->getDataNonConst(0), y = coarseCoords->getDataNonConst(1);
120 for (
size_t LID = 0; LID < coarseMap->getNodeNumElements(); ++LID) {
121 GlobalOrdinal GID = coarseMap->getGlobalElement(LID) - coarseMap->getIndexBase();
122 GlobalOrdinal i = GID % n, j = GID/n;
126 Set(coarseLevel,
"Coordinates", coarseCoords);
129 RCP<ParameterList> params = rcp(
new ParameterList());
130 params->set(
"printLoadBalancingInfo",
true);
138 #define MUELU_USERPFACTORY_SHORT 139 #endif // MUELU_USERPFACTORY_DEF_HPP Timer to be used in factories. Similar to Monitor but with additional timers.
Namespace for MueLu classes and methods.
void DeclareInput(Level &fineLevel, Level &coarseLevel) const
Input.
Class that holds all level-specific information.
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.
void Build(Level &fineLevel, Level &coarseLevel) const
Build an object with this factory.
Exception throws to report errors in the internal logical of the program.
RCP< const ParameterList > GetValidParameterList() const
Return a const parameter list of valid parameters that setParameterList() will accept.