27 #ifndef EWOMS_GLOBAL_INDICES_HH 28 #define EWOMS_GLOBAL_INDICES_HH 30 #include <dune/grid/common/datahandleif.hh> 31 #include <dune/istl/bcrsmatrix.hh> 32 #include <dune/istl/scalarproducts.hh> 33 #include <dune/istl/operators.hh> 54 template <
class ForeignOverlap>
59 typedef std::map<Index, Index> GlobalToDomesticMap;
60 typedef std::map<Index, Index> DomesticToGlobalMap;
64 : foreignOverlap_(foreignOverlap)
72 MPI_Comm_rank(MPI_COMM_WORLD, &tmp);
74 MPI_Comm_size(MPI_COMM_WORLD, &tmp);
75 mpiSize_ =
static_cast<size_t>(tmp);
82 buildGlobalIndices_();
90 assert(domesticToGlobal_.find(domesticIdx) != domesticToGlobal_.end());
92 return domesticToGlobal_.find(domesticIdx)->second;
100 const auto& tmp = globalToDomestic_.find(globalIdx);
102 if (tmp == globalToDomestic_.end())
113 {
return foreignOverlap_.numLocal(); }
122 {
return numDomestic_; }
129 domesticToGlobal_[domesticIdx] = globalIdx;
130 globalToDomestic_[globalIdx] = domesticIdx;
131 numDomestic_ = domesticToGlobal_.size();
133 assert(domesticToGlobal_.size() == globalToDomestic_.size());
143 sendBuf.peerIdx = peerLocalIdx;
148 static_cast<int>(peerRank),
165 static_cast<int>(peerRank),
170 Index domesticIdx = foreignOverlap_.nativeToLocal(recvBuf.peerIdx);
171 if (domesticIdx >= 0) {
172 Index globalIdx = recvBuf.globalIdx;
182 {
return globalToDomestic_.find(globalIdx) != globalToDomestic_.end(); }
190 std::cout <<
"(domestic index, global index, domestic->global->domestic)" 191 <<
" list for rank " << myRank_ <<
"\n";
193 for (
size_t domIdx = 0; domIdx < domesticToGlobal_.size(); ++domIdx)
196 std::cout <<
"\n" << std::flush;
202 void buildGlobalIndices_()
207 numDomestic_ = foreignOverlap_.numLocal();
218 MPI_Recv(&domesticOffset_,
221 static_cast<int>(myRank_ - 1),
230 for (
unsigned i = 0; i < foreignOverlap_.numLocal(); ++i) {
231 if (!foreignOverlap_.iAmMasterOf(static_cast<Index>(i)))
235 static_cast<Index>(domesticOffset_ + numMaster));
239 if (myRank_ < mpiSize_ - 1) {
242 int tmp = domesticOffset_ + numMaster;
246 static_cast<int>(myRank_ + 1),
251 typename PeerSet::const_iterator peerIt;
252 typename PeerSet::const_iterator peerEndIt = peerSet_().end();
254 peerIt = peerSet_().begin();
255 for (; peerIt != peerEndIt; ++peerIt) {
256 if (*peerIt < myRank_)
257 receiveBorderFrom_(*peerIt);
261 peerIt = peerSet_().begin();
262 for (; peerIt != peerEndIt; ++peerIt) {
263 if (*peerIt > myRank_)
264 sendBorderTo_(*peerIt);
268 peerIt = peerSet_().begin();
269 for (; peerIt != peerEndIt; ++peerIt) {
270 if (*peerIt > myRank_)
271 receiveBorderFrom_(*peerIt);
275 peerIt = peerSet_().begin();
276 for (; peerIt != peerEndIt; ++peerIt) {
277 if (*peerIt < myRank_)
278 sendBorderTo_(*peerIt);
283 void sendBorderTo_(ProcessRank peerRank)
288 BorderList::const_iterator borderIt = borderList_().begin();
289 BorderList::const_iterator borderEndIt = borderList_().end();
290 for (; borderIt != borderEndIt; ++borderIt) {
293 if (borderPeer != peerRank || borderDistance != 0)
296 Index localIdx = foreignOverlap_.nativeToLocal(borderIt->localIdx);
297 Index peerIdx = borderIt->peerIdx;
298 assert(localIdx >= 0);
299 if (foreignOverlap_.iAmMasterOf(localIdx)) {
306 void receiveBorderFrom_(ProcessRank peerRank)
311 BorderList::const_iterator borderIt = borderList_().begin();
312 BorderList::const_iterator borderEndIt = borderList_().end();
313 for (; borderIt != borderEndIt; ++borderIt) {
316 if (borderPeer != peerRank || borderDistance != 0)
319 Index nativeIdx = borderIt->localIdx;
320 Index localIdx = foreignOverlap_.nativeToLocal(nativeIdx);
321 if (localIdx >= 0 && foreignOverlap_.masterRank(localIdx) == borderPeer)
327 const PeerSet& peerSet_()
const 328 {
return foreignOverlap_.peerSet(); }
331 {
return foreignOverlap_.borderList(); }
338 const ForeignOverlap& foreignOverlap_;
340 GlobalToDomesticMap globalToDomestic_;
341 DomesticToGlobalMap domesticToGlobal_;
Index domesticToGlobal(Index domesticIdx) const
Converts a domestic index to a global one.
Definition: globalindices.hh:88
void sendBorderIndex(ProcessRank peerRank, Index domesticIdx, Index peerLocalIdx)
Send a border index to a remote process.
Definition: globalindices.hh:139
Definition: baseauxiliarymodule.hh:37
bool hasGlobalIndex(Index globalIdx) const
Return true iff a given global index already exists.
Definition: globalindices.hh:181
size_t numLocal() const
Returns the number of indices which are in the interior or on the border of the current rank...
Definition: globalindices.hh:112
unsigned BorderDistance
The type representing the distance of an index to the border.
Definition: overlaptypes.hh:54
This class maps domestic row indices to and from "global" indices which is used to construct an algeb...
Definition: globalindices.hh:55
int Index
The type of an index of a degree of freedom.
Definition: overlaptypes.hh:44
void addIndex(Index domesticIdx, Index globalIdx)
Add an index to the domestic<->global mapping.
Definition: globalindices.hh:127
void print() const
Prints the global indices of all domestic indices for debugging purposes.
Definition: globalindices.hh:188
This files provides several data structures for storing tuples of indices of remote and/or local proc...
unsigned ProcessRank
The type of the rank of a process.
Definition: overlaptypes.hh:49
size_t numDomestic() const
Returns the number domestic indices.
Definition: globalindices.hh:121
This structure stores a local index on a peer process and a global index.
Definition: overlaptypes.hh:69
void receiveBorderIndex(ProcessRank peerRank)
Receive an index on the border from a remote process and add it the translation maps.
Definition: globalindices.hh:158
std::list< BorderIndex > BorderList
This class managages a list of indices which are on the border of a process' partition of the grid...
Definition: overlaptypes.hh:120
Index globalToDomestic(Index globalIdx) const
Converts a global index to a domestic one.
Definition: globalindices.hh:98