53 #ifndef MUELU_AMALGAMATIONINFO_KOKKOS_DEF_HPP_ 54 #define MUELU_AMALGAMATIONINFO_KOKKOS_DEF_HPP_ 56 #include <Xpetra_MapFactory.hpp> 57 #include <Xpetra_Vector.hpp> 61 #ifdef HAVE_MUELU_KOKKOS_REFACTOR 63 #include "MueLu_Aggregates_kokkos.hpp" 67 template <
class LocalOrdinal,
class GlobalOrdinal,
class Node>
68 void AmalgamationInfo_kokkos<LocalOrdinal, GlobalOrdinal, Node>::
69 UnamalgamateAggregates(
const Aggregates_kokkos& aggregates,
70 Teuchos::ArrayRCP<LocalOrdinal>& aggStart,
71 Teuchos::ArrayRCP<GlobalOrdinal>& aggToRowMap)
const {
73 int myPid = aggregates.GetMap()->getComm()->getRank();
74 Teuchos::ArrayView<const GO> nodeGlobalElts = aggregates.GetMap()->getNodeElementList();
75 Teuchos::ArrayRCP<LO> procWinner = aggregates.GetProcWinner()->getDataNonConst(0);
76 Teuchos::ArrayRCP<LO> vertex2AggId = aggregates.GetVertex2AggId()->getDataNonConst(0);
77 LO size = procWinner.size();
78 GO numAggregates = aggregates.GetNumAggregates();
80 std::vector<LO> sizes(numAggregates);
81 if (stridedblocksize_ == 1) {
82 for (LO lnode = 0; lnode < size; ++lnode) {
83 LO myAgg = vertex2AggId[lnode];
84 if (procWinner[lnode] == myPid)
88 for (LO lnode = 0; lnode < size; ++lnode) {
89 LO myAgg = vertex2AggId[lnode];
90 if (procWinner[lnode] == myPid) {
91 GO gnodeid = nodeGlobalElts[lnode];
94 if (columnMap_->isNodeGlobalElement(gDofIndex))
100 aggStart = ArrayRCP<LO>(numAggregates+1,0);
102 for (GO i=0; i<numAggregates; ++i) {
103 aggStart[i+1] = aggStart[i] + sizes[i];
105 aggToRowMap = ArrayRCP<GO>(aggStart[numAggregates],0);
108 Array<LO> numDofs(numAggregates, 0);
110 if (stridedblocksize_ == 1) {
111 for (LO lnode = 0; lnode < size; ++lnode) {
112 LO myAgg = vertex2AggId[lnode];
113 if (procWinner[lnode] == myPid) {
114 aggToRowMap[ aggStart[myAgg] + numDofs[myAgg] ] = ComputeGlobalDOF(nodeGlobalElts[lnode]);
119 for (LO lnode = 0; lnode < size; ++lnode) {
120 LO myAgg = vertex2AggId[lnode];
122 if (procWinner[lnode] == myPid) {
123 GO gnodeid = nodeGlobalElts[lnode];
126 if (columnMap_->isNodeGlobalElement(gDofIndex)) {
127 aggToRowMap[ aggStart[myAgg] + numDofs[myAgg] ] = gDofIndex;
138 template <
class LocalOrdinal,
class GlobalOrdinal,
class Node>
139 void AmalgamationInfo_kokkos<LocalOrdinal, GlobalOrdinal, Node>::
140 UnamalgamateAggregatesLO(
const Aggregates_kokkos& aggregates,
141 Teuchos::ArrayRCP<LO>& aggStart,
142 Teuchos::ArrayRCP<LO>& aggToRowMap)
const {
144 int myPid = aggregates.GetMap()->getComm()->getRank();
145 Teuchos::ArrayView<const GO> nodeGlobalElts = aggregates.GetMap()->getNodeElementList();
147 Teuchos::ArrayRCP<LO> procWinner = aggregates.GetProcWinner() ->getDataNonConst(0);
148 Teuchos::ArrayRCP<LO> vertex2AggId = aggregates.GetVertex2AggId()->getDataNonConst(0);
149 const GO numAggregates = aggregates.GetNumAggregates();
153 LO size = procWinner.size();
155 std::vector<LO> sizes(numAggregates);
156 if (stridedblocksize_ == 1) {
157 for (LO lnode = 0; lnode < size; lnode++)
158 if (procWinner[lnode] == myPid)
159 sizes[vertex2AggId[lnode]]++;
161 for (LO lnode = 0; lnode < size; lnode++)
162 if (procWinner[lnode] == myPid) {
163 GO nodeGID = nodeGlobalElts[lnode];
165 for (LO k = 0; k < stridedblocksize_; k++) {
166 GO GID = ComputeGlobalDOF(nodeGID, k);
167 if (columnMap_->isNodeGlobalElement(GID))
168 sizes[vertex2AggId[lnode]]++;
173 aggStart = ArrayRCP<LO>(numAggregates+1);
175 for (GO i = 0; i < numAggregates; i++)
176 aggStart[i+1] = aggStart[i] + sizes[i];
178 aggToRowMap = ArrayRCP<LO>(aggStart[numAggregates], 0);
181 Array<LO> numDofs(numAggregates, 0);
182 if (stridedblocksize_ == 1) {
183 for (LO lnode = 0; lnode < size; ++lnode)
184 if (procWinner[lnode] == myPid) {
185 LO myAgg = vertex2AggId[lnode];
186 aggToRowMap[aggStart[myAgg] + numDofs[myAgg]] = lnode;
190 for (LO lnode = 0; lnode < size; ++lnode)
191 if (procWinner[lnode] == myPid) {
192 LO myAgg = vertex2AggId[lnode];
193 GO nodeGID = nodeGlobalElts[lnode];
195 for (LO k = 0; k < stridedblocksize_; k++) {
196 GO GID = ComputeGlobalDOF(nodeGID, k);
197 if (columnMap_->isNodeGlobalElement(GID)) {
198 aggToRowMap[aggStart[myAgg] + numDofs[myAgg]] = lnode*stridedblocksize_ + k;
210 template <
class LocalOrdinal,
class GlobalOrdinal,
class Node>
211 RCP<Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node> > AmalgamationInfo_kokkos<LocalOrdinal, GlobalOrdinal, Node>::
212 ComputeUnamalgamatedImportDofMap(
const Aggregates_kokkos& aggregates)
const {
214 Teuchos::RCP<const Map> nodeMap = aggregates.GetMap();
216 Teuchos::RCP<std::vector<GO> > myDofGids = Teuchos::rcp(
new std::vector<GO>);
217 Teuchos::ArrayView<const GO> gEltList = nodeMap->getNodeElementList();
218 LO nodeElements = Teuchos::as<LO>(nodeMap->getNodeNumElements());
219 if (stridedblocksize_ == 1) {
220 for (LO n = 0; n<nodeElements; n++) {
222 myDofGids->push_back(gDofIndex);
225 for (LO n = 0; n<nodeElements; n++) {
228 if (columnMap_->isNodeGlobalElement(gDofIndex))
229 myDofGids->push_back(gDofIndex);
234 Teuchos::ArrayRCP<GO> arr_myDofGids = Teuchos::arcp( myDofGids );
235 Teuchos::RCP<Map> importDofMap = MapFactory::Build(aggregates.GetMap()->lib(), Teuchos::OrdinalTraits<Xpetra::global_size_t>::invalid(), arr_myDofGids(), aggregates.GetMap()->getIndexBase(), aggregates.GetMap()->getComm());
241 template <
class LocalOrdinal,
class GlobalOrdinal,
class Node>
242 GlobalOrdinal AmalgamationInfo_kokkos<LocalOrdinal, GlobalOrdinal, Node>::
246 GlobalOrdinal gDofIndex = offset_ + (gNodeID-indexBase_)*fullblocksize_ + nStridedOffset_ + k + indexBase_;
253 #endif // HAVE_MUELU_KOKKOS_REFACTOR MueLu::DefaultLocalOrdinal LocalOrdinal
Namespace for MueLu classes and methods.
MueLu::DefaultGlobalOrdinal GlobalOrdinal