46 #ifndef MUELU_INTERFACEAGGREGATIONFACTORY_DEF_HPP_ 47 #define MUELU_INTERFACEAGGREGATIONFACTORY_DEF_HPP_ 52 #include "MueLu_Aggregates.hpp" 54 #include "Xpetra_MapFactory.hpp" 59 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 A");
65 validParamList->set<RCP<const FactoryBase>>(
"Aggregates", Teuchos::null,
"Generating factory of the Aggregates (for block 0,0)");
66 validParamList->set<RCP<const FactoryBase>>(
"DualNodeID2PrimalNodeID", Teuchos::null,
"Generating factory of the DualNodeID2PrimalNodeID map");
68 validParamList->set<
LocalOrdinal>(
"number of DOFs per dual node", Teuchos::ScalarTraits<LocalOrdinal>::one(),
"Number of DOFs per dual node");
70 return validParamList;
73 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
76 Input(currentLevel,
"A");
77 Input(currentLevel,
"Aggregates");
86 "DualNodeID2PrimalNodeID was not provided by the user on level 0!");
91 Input(currentLevel,
"DualNodeID2PrimalNodeID");
95 template <
class Scalar,
class LocalOrdinal,
class GlobalOrdinal,
class Node>
98 using Map = Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node>;
99 using MapFactory = Xpetra::MapFactory<LocalOrdinal, GlobalOrdinal, Node>;
101 using Dual2Primal_type = std::map<LocalOrdinal, LocalOrdinal>;
103 const char prefix[] =
"MueLu::InterfaceAggregationFactory::Build: ";
104 const ParameterList &pL = GetParameterList();
108 RCP<Matrix> A = Get<RCP<Matrix>>(currentLevel,
"A");
109 RCP<Aggregates> aggs00 = Get<RCP<Aggregates>>(currentLevel,
"Aggregates");
110 ArrayRCP<LocalOrdinal> vertex2AggIdInput = aggs00->GetVertex2AggId()->getDataNonConst(0);
112 ArrayRCP<LocalOrdinal> procWinnerInput = aggs00->GetProcWinner()->getDataNonConst(0);
114 RCP<Dual2Primal_type> lagr2dof;
116 if (currentLevel.GetLevelID() == 0)
117 lagr2dof = currentLevel.Get<RCP<Dual2Primal_type>>(
"DualNodeID2PrimalNodeID",
NoFactory::get());
119 lagr2dof = Get<RCP<Dual2Primal_type>>(currentLevel,
"DualNodeID2PrimalNodeID");
123 RCP<const Map> aggRangeMap = A->getRangeMap();
124 const size_t myRank = aggRangeMap->getComm()->getRank();
126 LocalOrdinal globalNumDualNodes = aggRangeMap->getGlobalNumElements() / nDOFPerDualNode;
127 LocalOrdinal numDualNodes = aggRangeMap->getNodeNumElements() / nDOFPerDualNode;
129 TEUCHOS_TEST_FOR_EXCEPTION(numDualNodes !=
LocalOrdinal(lagr2dof->size()), std::runtime_error, prefix <<
" MueLu requires the range map and the DualNodeID2PrimalNodeID map to be compatible.");
131 RCP<const Map> aggVertexMap;
133 if (nDOFPerDualNode == 1)
134 aggVertexMap = aggRangeMap;
138 auto comm = aggRangeMap->getComm();
139 std::vector<GlobalOrdinal> myDualNodes = {};
141 for (
size_t i = 0; i < aggRangeMap->getNodeNumElements(); i += nDOFPerDualNode)
142 myDualNodes.push_back((aggRangeMap->getGlobalElement(i) - indexBase) / nDOFPerDualNode + indexBase);
144 aggVertexMap = MapFactory::Build(aggRangeMap->lib(), globalNumDualNodes, myDualNodes, indexBase, comm);
147 RCP<Aggregates> aggregates = rcp(
new Aggregates(aggVertexMap));
148 aggregates->setObjectLabel(
"IA");
149 aggregates->AggregatesCrossProcessors(
false);
151 ArrayRCP<LocalOrdinal> vertex2AggId = aggregates->GetVertex2AggId()->getDataNonConst(0);
152 ArrayRCP<LocalOrdinal> procWinner = aggregates->GetProcWinner()->getDataNonConst(0);
154 RCP<Dual2Primal_type> coarseLagr2Dof = rcp(
new Dual2Primal_type());
155 RCP<Dual2Primal_type> coarseDof2Lagr = rcp(
new Dual2Primal_type());
160 for (
LocalOrdinal localDualNodeID = 0; localDualNodeID < numDualNodes; ++localDualNodeID)
163 LocalOrdinal localPrimalNodeID = (*lagr2dof)[localDualNodeID];
166 LocalOrdinal currentAggIdInput = vertex2AggIdInput[localPrimalNodeID];
169 if (coarseDof2Lagr->count(currentAggIdInput) == 0)
172 (*coarseDof2Lagr)[currentAggIdInput] = numAggId;
173 (*coarseLagr2Dof)[numAggId] = currentAggIdInput;
178 vertex2AggId[localDualNodeID] = (*coarseDof2Lagr)[currentAggIdInput];
179 procWinner[localDualNodeID] = myRank;
182 aggregates->SetNumAggregates(numAggId);
183 Set(currentLevel,
"Aggregates", aggregates);
184 Set(currentLevel,
"CoarseDualNodeID2PrimalNodeID", coarseLagr2Dof);
185 GetOStream(
Statistics1) << aggregates->description() << std::endl;
MueLu::DefaultLocalOrdinal LocalOrdinal
Container class for aggregation information.
Timer to be used in factories. Similar to Monitor but with additional timers.
Namespace for MueLu classes and methods.
static const NoFactory * get()
int GetLevelID() const
Return level number.
MueLu::DefaultGlobalOrdinal GlobalOrdinal
Class that holds all level-specific information.
void Build(Level ¤tLevel) const override
Build aggregates.
void DeclareInput(Level ¤tLevel) const override
Input.
RCP< const ParameterList > GetValidParameterList() const override
Return a const parameter list of valid parameters that setParameterList() will accept.
bool IsAvailable(const std::string &ename, const FactoryBase *factory=NoFactory::get()) const
Test whether a need's value has been saved.
Exception throws to report errors in the internal logical of the program.
void DeclareInput(const std::string &ename, const FactoryBase *factory, const FactoryBase *requestedBy=NoFactory::get())
Callback from FactoryBase::CallDeclareInput() and FactoryBase::DeclareInput()