27 #ifndef EWOMS_FOREIGN_OVERLAP_FROM_BCRS_MATRIX_HH
28 #define EWOMS_FOREIGN_OVERLAP_FROM_BCRS_MATRIX_HH
35 #include <opm/common/Unused.hpp>
37 #include <dune/grid/common/datahandleif.hh>
38 #include <dune/istl/bcrsmatrix.hh>
39 #include <dune/istl/scalarproducts.hh>
40 #include <dune/istl/operators.hh>
71 template <
class BCRSMatrix>
76 : borderList_(borderList), blackList_(blackList)
84 MPI_Comm_rank(MPI_COMM_WORLD, &tmp);
85 myRank_ =
static_cast<ProcessRank
>(tmp);
91 createLocalIndices_();
95 auto it = borderList.begin();
96 const auto& endIt = borderList.end();
97 for (; it != endIt; ++it) {
102 localBorderIndices_.insert(localIdx);
107 neighborPeerSet_.update(borderList);
113 peerSet_ = neighborPeerSet_;
118 initialSeedList.update(borderList);
123 auto borderIt = borderList.begin();
124 const auto& borderEndIt = borderList.end();
125 for (; borderIt != borderEndIt; ++borderIt) {
126 minBorderDist = std::min(minBorderDist, borderIt->borderDistance);
131 foreignOverlapByLocalIndex_.resize(
numLocal());
132 extendForeignOverlap_(A, initialSeedList, minBorderDist, overlapSize);
136 computeMasterRanks_();
139 groupForeignOverlapByRank_();
146 {
return overlapSize_; }
152 {
return localBorderIndices_.count(localIdx) > 0; }
160 const auto& indexOverlap = foreignOverlapByLocalIndex_[
static_cast<unsigned>(localIdx)];
161 const auto& borderDistIt = indexOverlap.find(peerRank);
162 if (borderDistIt == indexOverlap.end())
166 return borderDistIt->second == 0;
174 {
return masterRank_[
static_cast<unsigned>(localIdx)]; }
185 {
return masterRank_[
static_cast<unsigned>(localIdx)] == myRank_; }
192 {
return borderList_; }
201 assert(foreignOverlapByRank_.find(peerRank) != foreignOverlapByRank_.end());
202 return foreignOverlapByRank_.find(peerRank)->second;
209 const std::map<ProcessRank, BorderDistance> &
213 return foreignOverlapByLocalIndex_[
static_cast<unsigned>(localIdx)];
221 const auto& idxOverlap = foreignOverlapByLocalIndex_[localIdx];
222 return idxOverlap.find(peerRank) != idxOverlap.end();
231 const auto& peerOverlap = foreignOverlapByRank_.find(peerRank)->second;
234 auto it = peerOverlap.begin();
235 const auto& endIt = peerOverlap.end();
236 for (; it != endIt; ++it) {
237 if (it->borderDistance == overlapSize_)
249 const auto& idxOverlap = foreignOverlapByLocalIndex_[localIdx];
251 auto it = idxOverlap.find(peerRank);
252 if (it == idxOverlap.end())
255 return it->second == overlapSize_;
270 {
return neighborPeerSet_; }
276 {
return numNative_; }
282 {
return numLocal_; }
288 {
return static_cast<unsigned>(domesticIdx) <
numLocal(); }
297 {
return nativeToLocalIndices_[
static_cast<unsigned>(nativeIdx)]; }
304 assert(localIdx < static_cast<Index>(localToNativeIndices_.size()));
305 return localToNativeIndices_[
static_cast<unsigned>(localIdx)];
312 {
return blackList_; }
319 {
return foreignOverlapByLocalIndex_[
static_cast<unsigned>(localIdx)].size(); }
326 {
return foreignOverlapByLocalIndex_[
static_cast<unsigned>(localIdx)].size() > 0; }
333 auto it = foreignOverlapByRank_.begin();
334 const auto& endIt = foreignOverlapByRank_.end();
335 for (; it != endIt; ++it) {
336 std::cout <<
"Overlap rows(distance) for rank " << it->first <<
": ";
338 auto rowIt = it->second.begin();
339 const auto& rowEndIt = it->second.end();
340 for (; rowIt != rowEndIt; ++rowIt)
341 std::cout << rowIt->index <<
"(" << rowIt->borderDistance <<
") ";
342 std::cout <<
"\n" << std::flush;
350 template <
class BCRSMatrix>
351 void extendForeignOverlap_(
const BCRSMatrix& A,
353 BorderDistance borderDistance,
357 addNonNeighborOverlapIndices_(A, seedList, borderDistance);
360 auto seedIt = seedList.begin();
361 const auto& seedEndIt = seedList.end();
362 for (; seedIt != seedEndIt; ++seedIt) {
364 ProcessRank peerRank = seedIt->peerRank;
365 unsigned distance = borderDistance;
368 if (foreignOverlapByLocalIndex_[static_cast<unsigned>(localIdx)].count(peerRank) == 0)
369 foreignOverlapByLocalIndex_[
static_cast<unsigned>(localIdx)][peerRank] = distance;
374 if (borderDistance >= overlapSize)
379 SeedList nextSeedList;
380 seedIt = seedList.begin();
381 for (; seedIt != seedEndIt; ++seedIt) {
382 Index nativeRowIdx = seedIt->index;
385 ProcessRank peerRank = seedIt->peerRank;
390 typedef typename BCRSMatrix::ConstColIterator ColIterator;
391 ColIterator colIt = A[
static_cast<unsigned>(nativeRowIdx)].begin();
392 ColIterator colEndIt = A[
static_cast<unsigned>(nativeRowIdx)].end();
393 for (; colIt != colEndIt; ++colIt) {
394 Index nativeColIdx =
static_cast<Index
>(colIt.index());
402 else if (foreignOverlapByLocalIndex_[static_cast<unsigned>(localColIdx)].count(peerRank) > 0)
406 bool hasIndex =
false;
407 typename SeedList::iterator sIt = nextSeedList.begin();
408 typename SeedList::iterator sEndIt = nextSeedList.end();
409 for (; sIt != sEndIt; ++sIt) {
410 if (sIt->index == nativeColIdx && sIt->peerRank == peerRank) {
420 IndexRankDist newTuple;
421 newTuple.index = nativeColIdx;
422 newTuple.peerRank = peerRank;
423 newTuple.borderDistance = seedIt->borderDistance + 1;
424 nextSeedList.push_back(newTuple);
432 extendForeignOverlap_(A, nextSeedList, borderDistance + 1, overlapSize);
436 void createLocalIndices_()
440 for (
unsigned nativeIdx = 0; nativeIdx < numNative_;) {
441 if (!blackList_.hasIndex(static_cast<Index>(nativeIdx))) {
442 localToNativeIndices_.push_back(static_cast<Index>(nativeIdx));
443 nativeToLocalIndices_.push_back(static_cast<Index>(localIdx));
448 nativeToLocalIndices_.push_back(-1);
453 numLocal_ = localToNativeIndices_.size();
456 Index localToPeerIdx_(Index localIdx, ProcessRank peerRank)
const
458 auto it = borderList_.begin();
459 const auto& endIt = borderList_.end();
460 for (; it != endIt; ++it) {
461 if (it->localIdx == localIdx && it->peerRank == peerRank)
468 template <
class BCRSMatrix>
469 void addNonNeighborOverlapIndices_(
const BCRSMatrix& A OPM_UNUSED,
471 BorderDistance borderDist)
478 std::map<ProcessRank, std::vector<BorderIndex> > borderIndices;
482 auto it = seedList.begin();
483 const auto& endIt = seedList.end();
484 for (; it != endIt; ++it) {
488 BorderIndex borderHandle;
489 borderHandle.localIdx = localIdx;
490 borderHandle.peerRank = it->peerRank;
491 borderHandle.borderDistance = it->borderDistance;
494 auto neighborIt = foreignOverlapByLocalIndex_[
static_cast<unsigned>(localIdx)].begin();
495 const auto& neighborEndIt = foreignOverlapByLocalIndex_[
static_cast<unsigned>(localIdx)].end();
496 for (; neighborIt != neighborEndIt; ++neighborIt) {
497 if (neighborIt->second != 0)
500 else if (neighborIt->first == borderHandle.peerRank)
505 Index peerIdx = localToPeerIdx_(localIdx, neighborIt->first);
510 borderHandle.peerIdx = peerIdx;
511 borderIndices[neighborIt->first].push_back(borderHandle);
518 std::map<ProcessRank, Ewoms::MpiBuffer<unsigned> > numIndicesSendBufs;
519 std::map<ProcessRank, Ewoms::MpiBuffer<BorderIndex> > indicesSendBufs;
522 for (; peerIt != peerEndIt; ++peerIt) {
523 ProcessRank peerRank = *peerIt;
524 size_t numIndices = borderIndices[peerRank].size();
525 numIndicesSendBufs[peerRank].resize(1);
526 numIndicesSendBufs[peerRank][0] =
static_cast<unsigned>(numIndices);
528 const auto& peerBorderIndices = borderIndices[peerRank];
529 indicesSendBufs[peerRank].resize(numIndices);
531 auto tmpIt = peerBorderIndices.begin();
532 const auto& tmpEndIt = peerBorderIndices.end();
534 for (; tmpIt != tmpEndIt; ++tmpIt, ++i) {
535 indicesSendBufs[peerRank][i] = *tmpIt;
541 for (; peerIt != peerEndIt; ++peerIt) {
542 ProcessRank neighborPeer = *peerIt;
543 numIndicesSendBufs[neighborPeer].send(neighborPeer);
544 indicesSendBufs[neighborPeer].send(neighborPeer);
548 std::map<ProcessRank, MpiBuffer<unsigned> > numIndicesRcvBufs;
549 std::map<ProcessRank, MpiBuffer<BorderIndex> > indicesRcvBufs;
551 for (; peerIt != peerEndIt; ++peerIt) {
552 ProcessRank neighborPeer = *peerIt;
553 auto& numIndicesRcvBuf = numIndicesRcvBufs[neighborPeer];
554 auto& indicesRcvBuf = indicesRcvBufs[neighborPeer];
556 numIndicesRcvBuf.resize(1);
557 numIndicesRcvBuf.receive(neighborPeer);
558 unsigned numIndices = numIndicesRcvBufs[neighborPeer][0];
559 indicesRcvBuf.resize(numIndices);
560 indicesRcvBuf.receive(neighborPeer);
565 for (
unsigned i = 0; i < numIndices; ++i) {
568 std::swap(indicesRcvBuf[i].localIdx, indicesRcvBuf[i].peerIdx);
570 ProcessRank peerRank = indicesRcvBuf[i].peerRank;
572 Index localIdx = indicesRcvBuf[i].localIdx;
576 const auto& distIt = foreignOverlapByLocalIndex_[
static_cast<unsigned>(localIdx)].find(peerRank);
577 if (distIt != foreignOverlapByLocalIndex_[static_cast<unsigned>(localIdx)].end())
581 bool inSeedList =
false;
582 auto seedIt = seedList.begin();
583 const auto& seedEndIt = seedList.end();
584 for (; seedIt != seedEndIt; ++seedIt) {
585 if (seedIt->index == localIdx && seedIt->peerRank == peerRank) {
593 IndexRankDist seedEntry;
594 seedEntry.index = localIdx;
595 seedEntry.peerRank = peerRank;
596 seedEntry.borderDistance = borderDist;
597 seedList.push_back(seedEntry);
600 peerSet_.insert(peerRank);
606 for (; peerIt != peerEndIt; ++peerIt) {
607 ProcessRank neighborPeer = *peerIt;
608 numIndicesSendBufs[neighborPeer].wait();
609 indicesSendBufs[neighborPeer].wait();
617 void computeMasterRanks_()
620 masterRank_.resize(numLocal_);
621 for (
unsigned localIdx = 0; localIdx < numLocal_; ++localIdx) {
623 if (
isBorder(static_cast<Index>(localIdx))) {
627 auto it = foreignOverlapByLocalIndex_[
static_cast<unsigned>(localIdx)].begin();
628 const auto& endIt = foreignOverlapByLocalIndex_[
static_cast<unsigned>(localIdx)].end();
629 for (; it != endIt; ++it) {
630 if (it->second == 0) {
633 masterRank = std::min<ProcessRank>(
masterRank, it->first);
637 masterRank_[
static_cast<unsigned>(localIdx)] = masterRank;
644 void groupForeignOverlapByRank_()
648 size_t numLocal = foreignOverlapByLocalIndex_.size();
649 for (
unsigned localIdx = 0; localIdx <
numLocal; ++localIdx) {
651 auto it = foreignOverlapByLocalIndex_[localIdx].begin();
652 const auto& endIt = foreignOverlapByLocalIndex_[localIdx].end();
653 size_t nRanks = foreignOverlapByLocalIndex_[localIdx].size();
654 for (; it != endIt; ++it) {
655 IndexDistanceNpeers tmp;
656 tmp.index =
static_cast<Index
>(localIdx);
657 tmp.borderDistance = it->second;
658 tmp.numPeers =
static_cast<unsigned>(nRanks);
659 foreignOverlapByRank_[it->first].push_back(tmp);
668 PeerSet neighborPeerSet_;
671 const BorderList& borderList_;
674 const BlackList& blackList_;
677 std::vector<Index> nativeToLocalIndices_;
678 std::vector<Index> localToNativeIndices_;
682 std::vector<ProcessRank> masterRank_;
686 std::set<Index> localBorderIndices_;
691 OverlapByIndex foreignOverlapByLocalIndex_;
694 OverlapByRank foreignOverlapByRank_;
697 unsigned overlapSize_;
bool iAmMasterOf(Index localIdx) const
Return true if the current rank is the "master" of an index.
Definition: foreignoverlapfrombcrsmatrix.hh:184
The list of indices which are on the process boundary.
Definition: overlaptypes.hh:125
void print() const
Print the foreign overlap for debugging purposes.
Definition: foreignoverlapfrombcrsmatrix.hh:331
bool isBorderWith(Index localIdx, ProcessRank peerRank) const
Returns true iff a local index is a border index shared with a given peer process.
Definition: foreignoverlapfrombcrsmatrix.hh:158
bool isInOverlap(Index localIdx) const
Returns true if a given local index is in the foreign overlap of any rank.
Definition: foreignoverlapfrombcrsmatrix.hh:325
size_t numFront(ProcessRank peerRank) const
Returns the number of front indices of a peer process in the local partition.
Definition: foreignoverlapfrombcrsmatrix.hh:229
const PeerSet & peerSet() const
Return the set of process ranks which share an overlap with the current process.
Definition: foreignoverlapfrombcrsmatrix.hh:262
Expresses which degrees of freedom are blacklisted for the parallel linear solvers and which domestic...
bool isFrontFor(ProcessRank peerRank, Index localIdx) const
Returns whether a given local index is on the front of a given peer rank.
Definition: foreignoverlapfrombcrsmatrix.hh:247
Simplifies handling of buffers to be used in conjunction with MPI.
ForeignOverlapFromBCRSMatrix(const BCRSMatrix &A, const BorderList &borderList, const BlackList &blackList, unsigned overlapSize)
Constructs the foreign overlap given a BCRS matrix and an initial list of border indices.
Definition: foreignoverlapfrombcrsmatrix.hh:72
const std::map< ProcessRank, BorderDistance > & foreignOverlapByLocalIndex(Index localIdx) const
Return the map of (peer rank, border distance) for a given local index.
Definition: foreignoverlapfrombcrsmatrix.hh:210
This class creates and manages the foreign overlap given an initial list of border indices and a BCRS...
Definition: foreignoverlapfrombcrsmatrix.hh:61
const BorderList & borderList() const
Returns the list of indices which intersect the process border.
Definition: foreignoverlapfrombcrsmatrix.hh:191
unsigned overlapSize() const
Returns the size of the overlap region.
Definition: foreignoverlapfrombcrsmatrix.hh:145
Expresses which degrees of freedom are blacklisted for the parallel linear solvers and which domestic...
Definition: blacklist.hh:47
bool peerHasIndex(ProcessRank peerRank, Index localIdx) const
Returns true iff a local index is seen by a peer rank.
Definition: foreignoverlapfrombcrsmatrix.hh:219
size_t numLocal() const
Returns the number of local indices.
Definition: foreignoverlapfrombcrsmatrix.hh:281
This files provides several data structures for storing tuples of indices of remote and/or local proc...
size_t numPeers(Index localIdx) const
Return the number of peer ranks for which a given local index is visible.
Definition: foreignoverlapfrombcrsmatrix.hh:318
const BlackList & blackList() const
Returns the object which represents the black-listed native indices.
Definition: foreignoverlapfrombcrsmatrix.hh:311
size_t numNative() const
Returns the number of native indices.
Definition: foreignoverlapfrombcrsmatrix.hh:275
Index nativeToLocal(Index nativeIdx) const
Convert a native index to a local one.
Definition: foreignoverlapfrombcrsmatrix.hh:296
ProcessRank masterRank(Index localIdx) const
Return the rank of the master process of an index.
Definition: foreignoverlapfrombcrsmatrix.hh:173
Index localToNative(Index localIdx) const
Convert a local index to a native one.
Definition: foreignoverlapfrombcrsmatrix.hh:302
const OverlapWithPeer & foreignOverlapWithPeer(ProcessRank peerRank) const
Return the list of (local indices, border distance, number of processes) triples which are in the ove...
Definition: foreignoverlapfrombcrsmatrix.hh:199
bool isLocal(Index domesticIdx) const
Returns true iff a domestic index is local.
Definition: foreignoverlapfrombcrsmatrix.hh:287
A set of process ranks.
Definition: overlaptypes.hh:148
const PeerSet & neighborPeerSet() const
Return the set of process ranks which share a border index with the current process.
Definition: foreignoverlapfrombcrsmatrix.hh:269
bool isBorder(Index localIdx) const
Returns true iff a local index is a border index.
Definition: foreignoverlapfrombcrsmatrix.hh:151