27#ifndef EWOMS_FOREIGN_OVERLAP_FROM_BCRS_MATRIX_HH
28#define EWOMS_FOREIGN_OVERLAP_FROM_BCRS_MATRIX_HH
35#include <dune/grid/common/datahandleif.hh>
36#include <dune/istl/bcrsmatrix.hh>
37#include <dune/istl/scalarproducts.hh>
38#include <dune/istl/operators.hh>
69 template <
class BCRSMatrix>
89 createLocalIndices_();
95 for (; it !=
endIt; ++it) {
100 localBorderIndices_.insert(localIdx);
111 peerSet_ = neighborPeerSet_;
129 foreignOverlapByLocalIndex_.resize(
numLocal());
134 computeMasterRanks_();
137 groupForeignOverlapByRank_();
144 {
return overlapSize_; }
150 {
return localBorderIndices_.count(localIdx) > 0; }
158 const auto&
indexOverlap = foreignOverlapByLocalIndex_[
static_cast<unsigned>(localIdx)];
172 {
return masterRank_[
static_cast<unsigned>(localIdx)]; }
183 {
return masterRank_[
static_cast<unsigned>(localIdx)] == myRank_; }
190 {
return borderList_; }
199 assert(foreignOverlapByRank_.find(peerRank) != foreignOverlapByRank_.end());
200 return foreignOverlapByRank_.find(peerRank)->second;
207 const std::map<ProcessRank, BorderDistance> &
211 return foreignOverlapByLocalIndex_[
static_cast<unsigned>(localIdx)];
219 const auto&
idxOverlap = foreignOverlapByLocalIndex_[localIdx];
229 const auto&
peerOverlap = foreignOverlapByRank_.find(peerRank)->second;
234 for (; it !=
endIt; ++it) {
235 if (it->borderDistance == overlapSize_)
247 const auto&
idxOverlap = foreignOverlapByLocalIndex_[localIdx];
253 return it->second == overlapSize_;
268 {
return neighborPeerSet_; }
274 {
return numNative_; }
280 {
return numLocal_; }
295 {
return nativeToLocalIndices_[
static_cast<unsigned>(
nativeIdx)]; }
302 assert(localIdx <
static_cast<Index>(localToNativeIndices_.size()));
303 return localToNativeIndices_[
static_cast<unsigned>(localIdx)];
310 {
return blackList_; }
317 {
return foreignOverlapByLocalIndex_[
static_cast<unsigned>(localIdx)].size(); }
324 {
return foreignOverlapByLocalIndex_[
static_cast<unsigned>(localIdx)].size() > 0; }
331 auto it = foreignOverlapByRank_.begin();
332 const auto&
endIt = foreignOverlapByRank_.end();
333 for (; it !=
endIt; ++it) {
334 std::cout <<
"Overlap rows(distance) for rank " << it->first <<
": ";
336 auto rowIt = it->second.begin();
337 const auto&
rowEndIt = it->second.end();
339 std::cout <<
rowIt->index <<
"(" <<
rowIt->borderDistance <<
") ";
340 std::cout <<
"\n" << std::flush;
348 template <
class BCRSMatrix>
349 void extendForeignOverlap_(
const BCRSMatrix&
A,
355 addNonNeighborOverlapIndices_(
A,
seedList, borderDistance);
366 if (foreignOverlapByLocalIndex_[
static_cast<unsigned>(localIdx)].count(peerRank) == 0)
367 foreignOverlapByLocalIndex_[
static_cast<unsigned>(localIdx)][peerRank] =
distance;
383 ProcessRank peerRank =
seedIt->peerRank;
388 using ColIterator =
typename BCRSMatrix::ConstColIterator;
400 else if (foreignOverlapByLocalIndex_[
static_cast<unsigned>(
localColIdx)].count(peerRank) > 0)
404 bool hasIndex =
false;
434 void createLocalIndices_()
441 nativeToLocalIndices_.push_back(
static_cast<Index>(localIdx));
446 nativeToLocalIndices_.push_back(-1);
451 numLocal_ = localToNativeIndices_.size();
454 Index localToPeerIdx_(Index localIdx, ProcessRank peerRank)
const
456 auto it = borderList_.begin();
457 const auto&
endIt = borderList_.end();
458 for (; it !=
endIt; ++it) {
459 if (it->localIdx == localIdx && it->peerRank == peerRank)
466 template <
class BCRSMatrix>
467 void addNonNeighborOverlapIndices_(
const BCRSMatrix&,
476 std::map<ProcessRank, std::vector<BorderIndex> >
borderIndices;
482 for (; it !=
endIt; ++it) {
492 auto neighborIt = foreignOverlapByLocalIndex_[
static_cast<unsigned>(localIdx)].begin();
493 const auto&
neighborEndIt = foreignOverlapByLocalIndex_[
static_cast<unsigned>(localIdx)].end();
574 const auto&
distIt = foreignOverlapByLocalIndex_[
static_cast<unsigned>(localIdx)].find(peerRank);
575 if (
distIt != foreignOverlapByLocalIndex_[
static_cast<unsigned>(localIdx)].end())
583 if (
seedIt->index == localIdx &&
seedIt->peerRank == peerRank) {
598 peerSet_.insert(peerRank);
615 void computeMasterRanks_()
618 masterRank_.resize(numLocal_);
619 for (
unsigned localIdx = 0; localIdx < numLocal_; ++localIdx) {
625 auto it = foreignOverlapByLocalIndex_[
static_cast<unsigned>(localIdx)].begin();
626 const auto&
endIt = foreignOverlapByLocalIndex_[
static_cast<unsigned>(localIdx)].end();
627 for (; it !=
endIt; ++it) {
628 if (it->second == 0) {
635 masterRank_[
static_cast<unsigned>(localIdx)] =
masterRank;
642 void groupForeignOverlapByRank_()
646 size_t numLocal = foreignOverlapByLocalIndex_.size();
647 for (
unsigned localIdx = 0; localIdx <
numLocal; ++localIdx) {
649 auto it = foreignOverlapByLocalIndex_[localIdx].begin();
650 const auto&
endIt = foreignOverlapByLocalIndex_[localIdx].end();
651 size_t nRanks = foreignOverlapByLocalIndex_[localIdx].size();
652 for (; it !=
endIt; ++it) {
653 IndexDistanceNpeers tmp;
654 tmp.index =
static_cast<Index>(localIdx);
655 tmp.borderDistance = it->second;
656 tmp.numPeers =
static_cast<unsigned>(
nRanks);
657 foreignOverlapByRank_[it->first].push_back(tmp);
666 PeerSet neighborPeerSet_;
672 const BlackList& blackList_;
675 std::vector<Index> nativeToLocalIndices_;
676 std::vector<Index> localToNativeIndices_;
680 std::vector<ProcessRank> masterRank_;
684 std::set<Index> localBorderIndices_;
695 unsigned overlapSize_;
Expresses which degrees of freedom are blacklisted for the parallel linear solvers and which domestic...
Expresses which degrees of freedom are blacklisted for the parallel linear solvers and which domestic...
Definition blacklist.hh:49
This class creates and manages the foreign overlap given an initial list of border indices and a BCRS...
Definition foreignoverlapfrombcrsmatrix.hh:60
const PeerSet & peerSet() const
Return the set of process ranks which share an overlap with the current process.
Definition foreignoverlapfrombcrsmatrix.hh:260
size_t numLocal() const
Returns the number of local indices.
Definition foreignoverlapfrombcrsmatrix.hh:279
bool isInOverlap(Index localIdx) const
Returns true if a given local index is in the foreign overlap of any rank.
Definition foreignoverlapfrombcrsmatrix.hh:323
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:197
void print() const
Print the foreign overlap for debugging purposes.
Definition foreignoverlapfrombcrsmatrix.hh:329
unsigned overlapSize() const
Returns the size of the overlap region.
Definition foreignoverlapfrombcrsmatrix.hh:143
bool isBorder(Index localIdx) const
Returns true iff a local index is a border index.
Definition foreignoverlapfrombcrsmatrix.hh:149
Index nativeToLocal(Index nativeIdx) const
Convert a native index to a local one.
Definition foreignoverlapfrombcrsmatrix.hh:294
const PeerSet & neighborPeerSet() const
Return the set of process ranks which share a border index with the current process.
Definition foreignoverlapfrombcrsmatrix.hh:267
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:70
bool peerHasIndex(ProcessRank peerRank, Index localIdx) const
Returns true iff a local index is seen by a peer rank.
Definition foreignoverlapfrombcrsmatrix.hh:217
size_t numNative() const
Returns the number of native indices.
Definition foreignoverlapfrombcrsmatrix.hh:273
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:208
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:245
const BorderList & borderList() const
Returns the list of indices which intersect the process border.
Definition foreignoverlapfrombcrsmatrix.hh:189
size_t numFront(ProcessRank peerRank) const
Returns the number of front indices of a peer process in the local partition.
Definition foreignoverlapfrombcrsmatrix.hh:227
Index localToNative(Index localIdx) const
Convert a local index to a native one.
Definition foreignoverlapfrombcrsmatrix.hh:300
size_t numPeers(Index localIdx) const
Return the number of peer ranks for which a given local index is visible.
Definition foreignoverlapfrombcrsmatrix.hh:316
ProcessRank masterRank(Index localIdx) const
Return the rank of the master process of an index.
Definition foreignoverlapfrombcrsmatrix.hh:171
const BlackList & blackList() const
Returns the object which represents the black-listed native indices.
Definition foreignoverlapfrombcrsmatrix.hh:309
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:156
bool iAmMasterOf(Index localIdx) const
Return true if the current rank is the "master" of an index.
Definition foreignoverlapfrombcrsmatrix.hh:182
bool isLocal(Index domesticIdx) const
Returns true iff a domestic index is local.
Definition foreignoverlapfrombcrsmatrix.hh:285
A set of process ranks.
Definition overlaptypes.hh:149
The list of indices which are on the process boundary.
Definition overlaptypes.hh:126
Simplifies handling of buffers to be used in conjunction with MPI.
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition blackoilboundaryratevector.hh:37
constexpr auto getPropValue()
get the value data member of a property
Definition propertysystem.hh:242
This files provides several data structures for storing tuples of indices of remote and/or local proc...
unsigned BorderDistance
The type representing the distance of an index to the border.
Definition overlaptypes.hh:54
unsigned ProcessRank
The type of the rank of a process.
Definition overlaptypes.hh:49
std::vector< std::map< ProcessRank, BorderDistance > > OverlapByIndex
Maps each index to a list of processes .
Definition overlaptypes.hh:176
std::vector< IndexDistanceNpeers > OverlapWithPeer
The list of indices which overlap with a peer rank.
Definition overlaptypes.hh:165
int Index
The type of an index of a degree of freedom.
Definition overlaptypes.hh:44
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
std::map< ProcessRank, OverlapWithPeer > OverlapByRank
A type mapping the process rank to the list of indices shared with this peer.
Definition overlaptypes.hh:171