53 #ifndef MUELU_AMALGAMATIONINFO_DEF_HPP_ 54 #define MUELU_AMALGAMATIONINFO_DEF_HPP_ 56 #include <Xpetra_MapFactory.hpp> 57 #include <Xpetra_Vector.hpp> 61 #include "MueLu_Aggregates.hpp" 65 template <
class LocalOrdinal,
class GlobalOrdinal,
class Node>
67 Teuchos::ArrayRCP<LocalOrdinal>& aggStart, Teuchos::ArrayRCP<GlobalOrdinal>& aggToRowMap)
const {
68 int myPid = aggregates.
GetMap()->getComm()->getRank();
69 Teuchos::ArrayView<const GO> nodeGlobalElts = aggregates.
GetMap()->getNodeElementList();
70 Teuchos::ArrayRCP<LO> procWinner = aggregates.
GetProcWinner()->getDataNonConst(0);
71 Teuchos::ArrayRCP<LO> vertex2AggId = aggregates.
GetVertex2AggId()->getDataNonConst(0);
72 LO size = procWinner.size();
75 std::vector<LO> sizes(numAggregates);
76 if (stridedblocksize_ == 1) {
77 for (LO lnode = 0; lnode < size; ++lnode) {
78 LO myAgg = vertex2AggId[lnode];
79 if (procWinner[lnode] == myPid)
83 for (LO lnode = 0; lnode < size; ++lnode) {
84 LO myAgg = vertex2AggId[lnode];
85 if (procWinner[lnode] == myPid) {
86 GO gnodeid = nodeGlobalElts[lnode];
87 for (LocalOrdinal k = 0; k < stridedblocksize_; k++) {
88 GlobalOrdinal gDofIndex = ComputeGlobalDOF(gnodeid,k);
89 if (columnMap_->isNodeGlobalElement(gDofIndex))
95 aggStart = ArrayRCP<LO>(numAggregates+1,0);
97 for (GO i=0; i<numAggregates; ++i) {
98 aggStart[i+1] = aggStart[i] + sizes[i];
100 aggToRowMap = ArrayRCP<GO>(aggStart[numAggregates],0);
103 Array<LO> numDofs(numAggregates, 0);
105 if (stridedblocksize_ == 1) {
106 for (LO lnode = 0; lnode < size; ++lnode) {
107 LO myAgg = vertex2AggId[lnode];
108 if (procWinner[lnode] == myPid) {
109 aggToRowMap[ aggStart[myAgg] + numDofs[myAgg] ] = ComputeGlobalDOF(nodeGlobalElts[lnode]);
114 for (LO lnode = 0; lnode < size; ++lnode) {
115 LO myAgg = vertex2AggId[lnode];
117 if (procWinner[lnode] == myPid) {
118 GO gnodeid = nodeGlobalElts[lnode];
119 for (LocalOrdinal k = 0; k < stridedblocksize_; k++) {
120 GlobalOrdinal gDofIndex = ComputeGlobalDOF(gnodeid,k);
121 if (columnMap_->isNodeGlobalElement(gDofIndex)) {
122 aggToRowMap[ aggStart[myAgg] + numDofs[myAgg] ] = gDofIndex;
133 template <
class LocalOrdinal,
class GlobalOrdinal,
class Node>
135 Teuchos::ArrayRCP<LO>& aggStart, Teuchos::ArrayRCP<LO>& aggToRowMap)
const {
137 int myPid = aggregates.
GetMap()->getComm()->getRank();
138 Teuchos::ArrayView<const GO> nodeGlobalElts = aggregates.
GetMap()->getNodeElementList();
140 Teuchos::ArrayRCP<LO> procWinner = aggregates.
GetProcWinner() ->getDataNonConst(0);
141 Teuchos::ArrayRCP<LO> vertex2AggId = aggregates.
GetVertex2AggId()->getDataNonConst(0);
146 LO size = procWinner.size();
148 std::vector<LO> sizes(numAggregates);
149 if (stridedblocksize_ == 1) {
150 for (LO lnode = 0; lnode < size; lnode++)
151 if (procWinner[lnode] == myPid)
152 sizes[vertex2AggId[lnode]]++;
154 for (LO lnode = 0; lnode < size; lnode++)
155 if (procWinner[lnode] == myPid) {
156 GO nodeGID = nodeGlobalElts[lnode];
158 for (LO k = 0; k < stridedblocksize_; k++) {
159 GO GID = ComputeGlobalDOF(nodeGID, k);
160 if (columnMap_->isNodeGlobalElement(GID))
161 sizes[vertex2AggId[lnode]]++;
166 aggStart = ArrayRCP<LO>(numAggregates+1);
168 for (GO i = 0; i < numAggregates; i++)
169 aggStart[i+1] = aggStart[i] + sizes[i];
171 aggToRowMap = ArrayRCP<LO>(aggStart[numAggregates], 0);
174 Array<LO> numDofs(numAggregates, 0);
175 if (stridedblocksize_ == 1) {
176 for (LO lnode = 0; lnode < size; ++lnode)
177 if (procWinner[lnode] == myPid) {
178 LO myAgg = vertex2AggId[lnode];
179 aggToRowMap[aggStart[myAgg] + numDofs[myAgg]] = lnode;
183 for (LO lnode = 0; lnode < size; ++lnode)
184 if (procWinner[lnode] == myPid) {
185 LO myAgg = vertex2AggId[lnode];
186 GO nodeGID = nodeGlobalElts[lnode];
188 for (LO k = 0; k < stridedblocksize_; k++) {
189 GO GID = ComputeGlobalDOF(nodeGID, k);
190 if (columnMap_->isNodeGlobalElement(GID)) {
191 aggToRowMap[aggStart[myAgg] + numDofs[myAgg]] = lnode*stridedblocksize_ + k;
203 template <
class LocalOrdinal,
class GlobalOrdinal,
class Node>
205 Teuchos::RCP<const Map> nodeMap = aggregates.
GetMap();
207 Teuchos::RCP<std::vector<GO> > myDofGids = Teuchos::rcp(
new std::vector<GO>);
208 Teuchos::ArrayView<const GO> gEltList = nodeMap->getNodeElementList();
209 LO nodeElements = Teuchos::as<LO>(nodeMap->getNodeNumElements());
210 if (stridedblocksize_ == 1) {
211 for (LO n = 0; n<nodeElements; n++) {
212 GlobalOrdinal gDofIndex = ComputeGlobalDOF(gEltList[n]);
213 myDofGids->push_back(gDofIndex);
216 for (LO n = 0; n<nodeElements; n++) {
217 for (LocalOrdinal k = 0; k < stridedblocksize_; k++) {
218 GlobalOrdinal gDofIndex = ComputeGlobalDOF(gEltList[n],k);
219 if (columnMap_->isNodeGlobalElement(gDofIndex))
220 myDofGids->push_back(gDofIndex);
225 Teuchos::ArrayRCP<GO> arr_myDofGids = Teuchos::arcp( myDofGids );
226 Teuchos::RCP<Map> importDofMap = MapFactory::Build(aggregates.
GetMap()->lib(), Teuchos::OrdinalTraits<Xpetra::global_size_t>::invalid(), arr_myDofGids(), aggregates.
GetMap()->getIndexBase(), aggregates.
GetMap()->getComm());
232 template <
class LocalOrdinal,
class GlobalOrdinal,
class Node>
236 GlobalOrdinal gDofIndex = offset_ + (gNodeID-indexBase_)*fullblocksize_ + nStridedOffset_ + k + indexBase_;
Container class for aggregation information.
Namespace for MueLu classes and methods.
Teuchos::RCP< Xpetra::Map< LocalOrdinal, GlobalOrdinal, Node > > ComputeUnamalgamatedImportDofMap(const Aggregates &aggregates) const
ComputeUnamalgamatedImportDofMap build overlapping dof row map from aggregates needed for overlapping...
void UnamalgamateAggregatesLO(const Aggregates &aggregates, Teuchos::ArrayRCP< LocalOrdinal > &aggStart, Teuchos::ArrayRCP< LO > &aggToRowMap) const
const RCP< const Map > GetMap() const
returns (overlapping) map of aggregate/node distribution
void UnamalgamateAggregates(const Aggregates &aggregates, Teuchos::ArrayRCP< LocalOrdinal > &aggStart, Teuchos::ArrayRCP< GlobalOrdinal > &aggToRowMap) const
UnamalgamateAggregates.
LO GetNumAggregates() const
returns the number of aggregates of the current processor. Note: could/should be renamed to GetNumLoc...
const RCP< LOMultiVector > & GetVertex2AggId() const
Returns constant vector that maps local node IDs to local aggregates IDs.
const RCP< LOVector > & GetProcWinner() const
Returns constant vector that maps local node IDs to owning processor IDs.
GO ComputeGlobalDOF(GO const &gNodeID, LO const &k=0) const
ComputeGlobalDOF return global dof id associated with global node id gNodeID and dof index k...