47 #ifndef PACKAGES_MUELU_SRC_GRAPH_MUELU_UNSMOOSHFACTORY_DEF_HPP_
48 #define PACKAGES_MUELU_SRC_GRAPH_MUELU_UNSMOOSHFACTORY_DEF_HPP_
56 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
59 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
62 validParamList->
set<
RCP<const FactoryBase> >(
"A", Teuchos::null,
"Generating factory for unamalgamated matrix. Row map of (unamalgamted) output prolongation operator should match row map of this A.");
63 validParamList->
set<
RCP<const FactoryBase> >(
"P", Teuchos::null,
"Generating factory of the (amalgamated) prolongator P");
64 validParamList->
set<
RCP<const FactoryBase> >(
"DofStatus", Teuchos::null,
"Generating factory for dofStatus array (usually the VariableDofLaplacdianFactory)");
66 validParamList->
set<
int > (
"maxDofPerNode", 1,
"Maximum number of DOFs per node");
67 validParamList->
set<
bool > (
"fineIsPadded" ,
false,
"true if finest level input matrix is padded");
69 return validParamList;
72 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
75 Input(fineLevel,
"A");
76 Input(coarseLevel,
"P");
81 Input(fineLevel,
"DofStatus");
84 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
92 RCP<Matrix> unamalgA = Get< RCP<Matrix> >(fineLevel,
"A");
93 RCP<Matrix> amalgP = Get< RCP<Matrix> >(coarseLevel,
"P");
96 int maxDofPerNode = pL.
get<
int> (
"maxDofPerNode");
97 bool fineIsPadded = pL.
get<
bool>(
"fineIsPadded");
103 if(fineLevel.GetLevelID() == 0) {
104 dofStatus = Get<Teuchos::Array<char> >(fineLevel,
"DofStatus");
109 bool bHasZeroDiagonal =
false;
113 for(decltype(dirOrNot.
size()) i = 0; i < dirOrNot.
size(); ++i) {
114 if(dirOrNot[i] ==
true) dofStatus[i] =
'p';
127 amalgPcrs->getAllValues(amalgRowPtr, amalgCols, amalgVals);
130 size_t paddedNrows = amalgP->getRowMap()->getNodeNumElements() * Teuchos::as<size_t>(maxDofPerNode);
138 if(fineIsPadded ==
true || fineLevel.GetLevelID() > 0) {
147 for (decltype(amalgRowPtr.size()) i = 0; i < amalgRowPtr.size() - 1; i++) {
149 size_t rowLength = amalgRowPtr[i+1] - amalgRowPtr[i];
152 for(
int j = 0; j < maxDofPerNode; j++) {
153 newPRowPtr[i*maxDofPerNode+j] = cnt;
154 if (dofStatus[i*maxDofPerNode+j] ==
's') {
156 for (
size_t k = 0; k < rowLength; k++) {
157 newPCols[cnt ] = amalgCols[k+amalgRowPtr[i]] * maxDofPerNode + j;
158 newPVals[cnt++] = amalgVals[k+amalgRowPtr[i]];
165 newPRowPtr[paddedNrows] = cnt;
166 rowCount = paddedNrows;
174 for (decltype(amalgRowPtr.size()) i = 0; i < amalgRowPtr.size() - 1; i++) {
176 size_t rowLength = amalgRowPtr[i+1] - amalgRowPtr[i];
179 for(
int j = 0; j < maxDofPerNode; j++) {
182 if (dofStatus[i*maxDofPerNode+j] ==
's') {
183 newPRowPtr[rowCount++] = cnt;
185 for (
size_t k = 0; k < rowLength; k++) {
186 newPCols[cnt ] = amalgCols[k+amalgRowPtr[i]] * maxDofPerNode + j;
187 newPVals[cnt++] = amalgVals[k+amalgRowPtr[i]];
191 if (dofStatus[i*maxDofPerNode+j] ==
'd') {
192 newPRowPtr[rowCount++] = cnt;
196 newPRowPtr[rowCount] = cnt;
202 std::vector<size_t> stridingInfo(1,maxDofPerNode);
204 GlobalOrdinal nCoarseDofs = amalgP->getDomainMap()->getNodeNumElements() * maxDofPerNode;
205 GlobalOrdinal indexBase = amalgP->getDomainMap()->getIndexBase();
206 RCP<const Map> coarseDomainMap = StridedMapFactory::Build(amalgP->getDomainMap()->lib(),
211 amalgP->getDomainMap()->getComm(),
215 size_t nColCoarseDofs = Teuchos::as<size_t>(amalgP->getColMap()->getNodeNumElements() * maxDofPerNode);
217 for(
size_t c = 0; c < amalgP->getColMap()->getNodeNumElements(); ++c) {
218 GlobalOrdinal gid = (amalgP->getColMap()->getGlobalElement(c)-indexBase) * maxDofPerNode + indexBase;
220 for(
int i = 0; i < maxDofPerNode; ++i) {
221 unsmooshColMapGIDs[c * maxDofPerNode + i] = gid + i;
224 Teuchos::RCP<Map> coarseColMap = MapFactory::Build(amalgP->getDomainMap()->lib(),
226 unsmooshColMapGIDs(),
228 amalgP->getDomainMap()->getComm());
232 for (decltype(rowCount) i = 0; i < rowCount; i++) {
233 unamalgPCrs->insertLocalValues(i, newPCols.view(newPRowPtr[i],newPRowPtr[i+1]-newPRowPtr[i]),
234 newPVals.view(newPRowPtr[i],newPRowPtr[i+1]-newPRowPtr[i]));
236 unamalgPCrs->fillComplete(coarseDomainMap,unamalgA->getRowMap());
240 Set(coarseLevel,
"P",unamalgP);
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)
void Build(Level &fineLevel, Level &coarseLevel) const
Build an object with this factory.
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
void DeclareInput(Level &fineLevel, Level &coarseLevel) const
Input.
Class that holds all level-specific information.
static Teuchos::ArrayRCP< const bool > DetectDirichletRowsExt(const Xpetra::Matrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > &A, bool &bHasZeroDiagonal, const Magnitude &tol=Teuchos::ScalarTraits< Scalar >::zero())
UnsmooshFactory()
Constructor.
int GetLevelID() const
Return level number.
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.