46#ifndef MUELU_UNCOUPLEDINDEXMANAGER_DEF_HPP_
47#define MUELU_UNCOUPLEDINDEXMANAGER_DEF_HPP_
49#include <Xpetra_MapFactory.hpp>
50#include <Teuchos_OrdinalTraits.hpp>
55 template <
class LocalOrdinal,
class GlobalOrdinal,
class Node>
58 const int NumDimensions,
const int interpolationOrder,
59 const int MyRank,
const int NumRanks,
60 const Array<GO> GFineNodesPerDir,
const Array<LO> LFineNodesPerDir,
61 const Array<LO> CoarseRate,
const bool singleCoarsePoint) :
62 IndexManager(comm, coupled, singleCoarsePoint, NumDimensions, interpolationOrder,
63 Array<GO>(3, -1), LFineNodesPerDir),
64 myRank(MyRank), numRanks(NumRanks)
68 for(
int dim = 0; dim < 3; ++dim) {
70 if(CoarseRate.size() == 1) {
72 }
else if(CoarseRate.size() == this->numDimensions) {
85 template <
class LocalOrdinal,
class GlobalOrdinal,
class Node>
88 GO input[1] = {as<GO>(this->lNumCoarseNodes)}, output[1] = {0};
89 Teuchos::reduceAll(*(this->comm_), Teuchos::REDUCE_SUM, 1, input, output);
90 this->gNumCoarseNodes = output[0];
93 template <
class LocalOrdinal,
class GlobalOrdinal,
class Node>
96 Array<LO>& ghostedNodeCoarseLIDs,
97 Array<int>& ghostedNodeCoarsePIDs,
101 ghostedNodeCoarseLIDs.resize(this->getNumLocalGhostedNodes());
102 ghostedNodeCoarsePIDs.resize(this->getNumLocalGhostedNodes());
104 for(LO idx = 0; idx < this->getNumLocalGhostedNodes(); ++idx) {
105 ghostedNodeCoarseLIDs[idx] = idx;
106 ghostedNodeCoarsePIDs[idx] = myRank;
110 template <
class LocalOrdinal,
class GlobalOrdinal,
class Node>
113 Array<GO>& coarseNodeCoarseGIDs,
114 Array<GO>& coarseNodeFineGIDs)
const {
117 coarseNodeCoarseGIDs.resize(this->getNumLocalCoarseNodes());
118 coarseNodeFineGIDs.resize(this->getNumLocalCoarseNodes());
121 ArrayView<const GO> fineNodeGIDs = fineCoordinatesMap->getLocalElementList();
125 for(LO coarseLID = 0; coarseLID < this->getNumLocalCoarseNodes(); ++coarseLID) {
126 Array<LO> coarseIndices(3), fineIndices(3);
127 this->getCoarseNodeLocalTuple(coarseLID,
131 for(
int dim = 0; dim < 3; ++dim) {
132 if(coarseIndices[dim] == this->lCoarseNodesPerDir[dim] - 1) {
133 if(this->lCoarseNodesPerDir[dim] == 1) {
134 fineIndices[dim] = 0;
136 fineIndices[dim] = this->lFineNodesPerDir[dim] - 1;
139 fineIndices[dim] = coarseIndices[dim]*this->coarseRate[dim];
143 fineLID = fineIndices[2]*this->lNumFineNodes10
144 + fineIndices[1]*this->lFineNodesPerDir[0]
146 coarseNodeFineGIDs[coarseLID] = fineNodeGIDs[fineLID];
151 template <
class LocalOrdinal,
class GlobalOrdinal,
class Node>
154 std::vector<std::vector<GO> > coarseMeshData;
155 return coarseMeshData;
158 template <
class LocalOrdinal,
class GlobalOrdinal,
class Node>
163 template <
class LocalOrdinal,
class GlobalOrdinal,
class Node>
167 k = myLID / this->lNumFineNodes10;
168 tmp = myLID % this->lNumFineNodes10;
169 j = tmp / this->lFineNodesPerDir[0];
170 i = tmp % this->lFineNodesPerDir[0];
173 template <
class LocalOrdinal,
class GlobalOrdinal,
class Node>
177 k = myLID / this->lNumFineNodes10;
178 tmp = myLID % this->lNumFineNodes10;
179 j = tmp / this->lFineNodesPerDir[0];
180 i = tmp % this->lFineNodesPerDir[0];
182 k += this->offsets[2];
183 j += this->offsets[1];
184 i += this->offsets[0];
187 template <
class LocalOrdinal,
class GlobalOrdinal,
class Node>
192 template <
class LocalOrdinal,
class GlobalOrdinal,
class Node>
197 template <
class LocalOrdinal,
class GlobalOrdinal,
class Node>
202 template <
class LocalOrdinal,
class GlobalOrdinal,
class Node>
206 k = myLID / this->lNumCoarseNodes10;
207 tmp = myLID % this->lNumCoarseNodes10;
208 j = tmp / this->lCoarseNodesPerDir[0];
209 i = tmp % this->lCoarseNodesPerDir[0];
212 template <
class LocalOrdinal,
class GlobalOrdinal,
class Node>
217 template <
class LocalOrdinal,
class GlobalOrdinal,
class Node>
222 template <
class LocalOrdinal,
class GlobalOrdinal,
class Node>
225 myLID = k*this->numGhostedNodes10 + j*this->ghostedNodesPerDir[0] + i;
228 template <
class LocalOrdinal,
class GlobalOrdinal,
class Node>
233 template <
class LocalOrdinal,
class GlobalOrdinal,
class Node>
238 template <
class LocalOrdinal,
class GlobalOrdinal,
class Node>
Container class for mesh layout and indices calculation.
GO gNumCoarseNodes
global number of nodes remaining after coarsening.
void computeMeshParameters()
Array< int > coarseRate
coarsening rate in each direction
const int numDimensions
Number of spacial dimensions in the problem.
GO gNumCoarseNodes10
global number of nodes per 0-1 slice remaining after coarsening.
void getCoarseNodeGID(const GO i, const GO j, const GO k, GO &myGID) const
void getFineNodeGID(const GO i, const GO j, const GO k, GO &myGID) const
void getCoarseNodesData(const RCP< const Map > fineCoordinatesMap, Array< GO > &coarseNodeCoarseGIDs, Array< GO > &coarseNodeFineGIDs) const
void getGhostedNodeFineLID(const LO i, const LO j, const LO k, LO &myLID) const
void getCoarseNodeFineLID(const LO i, const LO j, const LO k, LO &myLID) const
void getGhostedNodesData(const RCP< const Map > fineMap, Array< LO > &ghostedNodeCoarseLIDs, Array< int > &ghostedNodeCoarsePIDs, Array< GO > &ghostedNodeCoarseGIDs) const
void getGhostedNodeCoarseLID(const LO i, const LO j, const LO k, LO &myLID) const
void getCoarseNodeLID(const LO i, const LO j, const LO k, LO &myLID) const
void getFineNodeLID(const LO i, const LO j, const LO k, LO &myLID) const
UncoupledIndexManager()=default
void getCoarseNodeGhostedLID(const LO i, const LO j, const LO k, LO &myLID) const
std::vector< std::vector< GO > > getCoarseMeshData() const
void getFineNodeLocalTuple(const LO myLID, LO &i, LO &j, LO &k) const
void getCoarseNodeGlobalTuple(const GO myGID, GO &i, GO &j, GO &k) const
void getCoarseNodeLocalTuple(const LO myLID, LO &i, LO &j, LO &k) const
void getFineNodeGhostedTuple(const LO myLID, LO &i, LO &j, LO &k) const
void computeGlobalCoarseParameters()
void getFineNodeGlobalTuple(const GO myGID, GO &i, GO &j, GO &k) const
Namespace for MueLu classes and methods.