All Classes Files Functions Variables Typedefs Enumerations Enumerator Macros Groups Pages
foreignoverlapfrombcrsmatrix.hh
Go to the documentation of this file.
1 // -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
2 // vi: set et ts=4 sw=4 sts=4:
3 /*
4  This file is part of the Open Porous Media project (OPM).
5 
6  OPM is free software: you can redistribute it and/or modify
7  it under the terms of the GNU General Public License as published by
8  the Free Software Foundation, either version 2 of the License, or
9  (at your option) any later version.
10 
11  OPM is distributed in the hope that it will be useful,
12  but WITHOUT ANY WARRANTY; without even the implied warranty of
13  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14  GNU General Public License for more details.
15 
16  You should have received a copy of the GNU General Public License
17  along with OPM. If not, see <http://www.gnu.org/licenses/>.
18 
19  Consult the COPYING file in the top-level source directory of this
20  module for the precise wording of the license and the list of
21  copyright holders.
22 */
27 #ifndef EWOMS_FOREIGN_OVERLAP_FROM_BCRS_MATRIX_HH
28 #define EWOMS_FOREIGN_OVERLAP_FROM_BCRS_MATRIX_HH
29 
30 #include "overlaptypes.hh"
31 #include "blacklist.hh"
32 
34 
35 #include <opm/common/Unused.hpp>
36 
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>
41 
42 #include <algorithm>
43 #include <iostream>
44 #include <map>
45 #include <vector>
46 
47 #if HAVE_MPI
48 #include <mpi.h>
49 #endif // HAVE_MPI
50 
51 namespace Ewoms {
52 namespace Linear {
53 
62 {
63 public:
64  // overlaps should never be copied!
66 
71  template <class BCRSMatrix>
72  ForeignOverlapFromBCRSMatrix(const BCRSMatrix& A,
73  const BorderList& borderList,
74  const BlackList& blackList,
75  unsigned overlapSize)
76  : borderList_(borderList), blackList_(blackList)
77  {
78  overlapSize_ = overlapSize;
79 
80  myRank_ = 0;
81 #if HAVE_MPI
82  {
83  int tmp;
84  MPI_Comm_rank(MPI_COMM_WORLD, &tmp);
85  myRank_ = static_cast<ProcessRank>(tmp);
86  }
87 #endif
88  numNative_ = A.N();
89 
90  // Computes the local <-> native index maps
91  createLocalIndices_();
92 
93  // calculate the set of local indices on the border (beware:
94  // _not_ the native ones)
95  auto it = borderList.begin();
96  const auto& endIt = borderList.end();
97  for (; it != endIt; ++it) {
98  Index localIdx = nativeToLocal(it->localIdx);
99  if (localIdx < 0)
100  continue;
101 
102  localBorderIndices_.insert(localIdx);
103  }
104 
105  // compute the set of processes which are neighbors of the
106  // local process ...
107  neighborPeerSet_.update(borderList);
108  // ... and the initial set of processes which we will have to
109  // communicate with. We must always communicate with our
110  // neighbors, but depending on the size of the overlap region,
111  // we might have to communicate with additional processes as
112  // well (these will be added later).
113  peerSet_ = neighborPeerSet_;
114 
115  // Create an initial seed list of indices which are in the
116  // overlap.
117  SeedList initialSeedList;
118  initialSeedList.update(borderList);
119 
120  // calculate the minimum distance from the border of the
121  // initial seed list
122  unsigned minBorderDist = overlapSize;
123  auto borderIt = borderList.begin();
124  const auto& borderEndIt = borderList.end();
125  for (; borderIt != borderEndIt; ++borderIt) {
126  minBorderDist = std::min(minBorderDist, borderIt->borderDistance);
127  }
128 
129  // calculate the foreign overlap for the local partition,
130  // i.e. find the distance of each row from the seed set.
131  foreignOverlapByLocalIndex_.resize(numLocal());
132  extendForeignOverlap_(A, initialSeedList, minBorderDist, overlapSize);
133 
134  // computes the process with the lowest rank for all local
135  // indices.
136  computeMasterRanks_();
137 
138  // group foreign overlap by peer process rank
139  groupForeignOverlapByRank_();
140  }
141 
145  unsigned overlapSize() const
146  { return overlapSize_; }
147 
151  bool isBorder(Index localIdx) const
152  { return localBorderIndices_.count(localIdx) > 0; }
153 
158  bool isBorderWith(Index localIdx, ProcessRank peerRank) const
159  {
160  const auto& indexOverlap = foreignOverlapByLocalIndex_[static_cast<unsigned>(localIdx)];
161  const auto& borderDistIt = indexOverlap.find(peerRank);
162  if (borderDistIt == indexOverlap.end())
163  return false;
164 
165  // border distance of the index needs to be 0
166  return borderDistIt->second == 0;
167  }
168 
173  ProcessRank masterRank(Index localIdx) const
174  { return masterRank_[static_cast<unsigned>(localIdx)]; }
175 
184  bool iAmMasterOf(Index localIdx) const
185  { return masterRank_[static_cast<unsigned>(localIdx)] == myRank_; }
186 
191  const BorderList& borderList() const
192  { return borderList_; }
193 
199  const OverlapWithPeer& foreignOverlapWithPeer(ProcessRank peerRank) const
200  {
201  assert(foreignOverlapByRank_.find(peerRank) != foreignOverlapByRank_.end());
202  return foreignOverlapByRank_.find(peerRank)->second;
203  }
204 
209  const std::map<ProcessRank, BorderDistance> &
210  foreignOverlapByLocalIndex(Index localIdx) const
211  {
212  assert(isLocal(localIdx));
213  return foreignOverlapByLocalIndex_[static_cast<unsigned>(localIdx)];
214  }
215 
219  bool peerHasIndex(ProcessRank peerRank, Index localIdx) const
220  {
221  const auto& idxOverlap = foreignOverlapByLocalIndex_[localIdx];
222  return idxOverlap.find(peerRank) != idxOverlap.end();
223  }
224 
229  size_t numFront(ProcessRank peerRank) const
230  {
231  const auto& peerOverlap = foreignOverlapByRank_.find(peerRank)->second;
232 
233  size_t n = 0;
234  auto it = peerOverlap.begin();
235  const auto& endIt = peerOverlap.end();
236  for (; it != endIt; ++it) {
237  if (it->borderDistance == overlapSize_)
238  ++n;
239  }
240  return n;
241  }
242 
247  bool isFrontFor(ProcessRank peerRank, Index localIdx) const
248  {
249  const auto& idxOverlap = foreignOverlapByLocalIndex_[localIdx];
250 
251  auto it = idxOverlap.find(peerRank);
252  if (it == idxOverlap.end())
253  return false; // index is not in overlap
254 
255  return it->second == overlapSize_;
256  }
257 
262  const PeerSet& peerSet() const
263  { return peerSet_; }
264 
269  const PeerSet& neighborPeerSet() const
270  { return neighborPeerSet_; }
271 
275  size_t numNative() const
276  { return numNative_; }
277 
281  size_t numLocal() const
282  { return numLocal_; }
283 
287  bool isLocal(Index domesticIdx) const
288  { return static_cast<unsigned>(domesticIdx) < numLocal(); }
289 
296  Index nativeToLocal(Index nativeIdx) const
297  { return nativeToLocalIndices_[static_cast<unsigned>(nativeIdx)]; }
298 
302  Index localToNative(Index localIdx) const
303  {
304  assert(localIdx < static_cast<Index>(localToNativeIndices_.size()));
305  return localToNativeIndices_[static_cast<unsigned>(localIdx)];
306  }
307 
311  const BlackList& blackList() const
312  { return blackList_; }
313 
318  size_t numPeers(Index localIdx) const
319  { return foreignOverlapByLocalIndex_[static_cast<unsigned>(localIdx)].size(); }
320 
325  bool isInOverlap(Index localIdx) const
326  { return foreignOverlapByLocalIndex_[static_cast<unsigned>(localIdx)].size() > 0; }
327 
331  void print() const
332  {
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 << ": ";
337 
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;
343  }
344  }
345 
346 protected:
347  // extend the foreign overlaps by 'overlapSize' levels. this uses
348  // a greedy algorithm which extends the region by one level and
349  // then calls itself recursively...
350  template <class BCRSMatrix>
351  void extendForeignOverlap_(const BCRSMatrix& A,
352  SeedList& seedList,
353  BorderDistance borderDistance,
354  BorderDistance overlapSize)
355  {
356  // communicate the non-neigbor overlap indices
357  addNonNeighborOverlapIndices_(A, seedList, borderDistance);
358 
359  // add all processes in the seed rows of the current overlap level
360  auto seedIt = seedList.begin();
361  const auto& seedEndIt = seedList.end();
362  for (; seedIt != seedEndIt; ++seedIt) {
363  Index localIdx = nativeToLocal(seedIt->index);
364  ProcessRank peerRank = seedIt->peerRank;
365  unsigned distance = borderDistance;
366  if (localIdx < 0)
367  continue;
368  if (foreignOverlapByLocalIndex_[static_cast<unsigned>(localIdx)].count(peerRank) == 0)
369  foreignOverlapByLocalIndex_[static_cast<unsigned>(localIdx)][peerRank] = distance;
370  }
371 
372  // if we have reached the maximum overlap distance, i.e. we're
373  // finished and break the recursion
374  if (borderDistance >= overlapSize)
375  return;
376 
377  // find the seed list for the next overlap level using the
378  // seed set for the current level
379  SeedList nextSeedList;
380  seedIt = seedList.begin();
381  for (; seedIt != seedEndIt; ++seedIt) {
382  Index nativeRowIdx = seedIt->index;
383  if (nativeToLocal(nativeRowIdx) < 0)
384  continue; // ignore blacklisted indices
385  ProcessRank peerRank = seedIt->peerRank;
386 
387  // find all column indices in the row. The indices of the
388  // columns are the additional indices of the overlap which
389  // we would like to add
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());
395  Index localColIdx = nativeToLocal(nativeColIdx);
396 
397  // ignore if the native index is not a local one
398  if (localColIdx < 0)
399  continue;
400  // if the process is already is in the overlap of the
401  // column index, ignore this column index!
402  else if (foreignOverlapByLocalIndex_[static_cast<unsigned>(localColIdx)].count(peerRank) > 0)
403  continue;
404 
405  // check whether the new index is already in the overlap
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) {
411  hasIndex = true;
412  break;
413  }
414  }
415  if (hasIndex)
416  continue; // we already have this index
417 
418  // add the current processes to the seed list for the
419  // next overlap level
420  IndexRankDist newTuple;
421  newTuple.index = nativeColIdx;
422  newTuple.peerRank = peerRank;
423  newTuple.borderDistance = seedIt->borderDistance + 1;
424  nextSeedList.push_back(newTuple);
425  }
426  }
427 
428  // clear the old seed list to save some memory
429  seedList.clear();
430 
431  // Perform the same excercise for the next overlap distance
432  extendForeignOverlap_(A, nextSeedList, borderDistance + 1, overlapSize);
433  }
434 
435  // Computes the local <-> native index maps
436  void createLocalIndices_()
437  {
438  // create the native <-> local maps
439  Index localIdx = 0;
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));
444  ++nativeIdx;
445  ++localIdx;
446  }
447  else {
448  nativeToLocalIndices_.push_back(-1);
449  ++nativeIdx;
450  }
451  }
452 
453  numLocal_ = localToNativeIndices_.size();
454  }
455 
456  Index localToPeerIdx_(Index localIdx, ProcessRank peerRank) const
457  {
458  auto it = borderList_.begin();
459  const auto& endIt = borderList_.end();
460  for (; it != endIt; ++it) {
461  if (it->localIdx == localIdx && it->peerRank == peerRank)
462  return it->peerIdx;
463  }
464 
465  return -1;
466  }
467 
468  template <class BCRSMatrix>
469  void addNonNeighborOverlapIndices_(const BCRSMatrix& A OPM_UNUSED,
470  SeedList& seedList,
471  BorderDistance borderDist)
472  {
473  // TODO: this probably does not work! (the matrix A is unused, but it is needed
474  // from a logical POV.)
475 #if HAVE_MPI
476  // first, create the buffers which will contain the number of
477  // border indices relevant for a neighbor peer
478  std::map<ProcessRank, std::vector<BorderIndex> > borderIndices;
479 
480  // get all indices in the border which have borderDist as
481  // their distance to the closest border of their local process
482  auto it = seedList.begin();
483  const auto& endIt = seedList.end();
484  for (; it != endIt; ++it) {
485  Index localIdx = nativeToLocal(it->index);
486  if (!isBorder(localIdx))
487  continue;
488  BorderIndex borderHandle;
489  borderHandle.localIdx = localIdx;
490  borderHandle.peerRank = it->peerRank;
491  borderHandle.borderDistance = it->borderDistance;
492 
493  // add the border index to all the neighboring peers
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)
498  // not a border index for the neighbor
499  continue;
500  else if (neighborIt->first == borderHandle.peerRank)
501  // don't communicate the indices which are owned
502  // by the peer to itself
503  continue;
504 
505  Index peerIdx = localToPeerIdx_(localIdx, neighborIt->first);
506  if (peerIdx < 0)
507  // the index is on the border, but is not on the border
508  // with the considered neighboring process. Ignore it!
509  continue;
510  borderHandle.peerIdx = peerIdx;
511  borderIndices[neighborIt->first].push_back(borderHandle);
512  }
513  }
514 
515  // now borderIndices contains the lists of indices which we
516  // would like to send to each neighbor. Let's create the MPI
517  // buffers.
518  std::map<ProcessRank, Ewoms::MpiBuffer<unsigned> > numIndicesSendBufs;
519  std::map<ProcessRank, Ewoms::MpiBuffer<BorderIndex> > indicesSendBufs;
520  auto peerIt = neighborPeerSet().begin();
521  const auto& peerEndIt = neighborPeerSet().end();
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);
527 
528  const auto& peerBorderIndices = borderIndices[peerRank];
529  indicesSendBufs[peerRank].resize(numIndices);
530 
531  auto tmpIt = peerBorderIndices.begin();
532  const auto& tmpEndIt = peerBorderIndices.end();
533  size_t i = 0;
534  for (; tmpIt != tmpEndIt; ++tmpIt, ++i) {
535  indicesSendBufs[peerRank][i] = *tmpIt;
536  }
537  }
538 
539  // now, send all these nice buffers to our neighbors
540  peerIt = neighborPeerSet().begin();
541  for (; peerIt != peerEndIt; ++peerIt) {
542  ProcessRank neighborPeer = *peerIt;
543  numIndicesSendBufs[neighborPeer].send(neighborPeer);
544  indicesSendBufs[neighborPeer].send(neighborPeer);
545  }
546 
547  // receive all data from the neighbors
548  std::map<ProcessRank, MpiBuffer<unsigned> > numIndicesRcvBufs;
549  std::map<ProcessRank, MpiBuffer<BorderIndex> > indicesRcvBufs;
550  peerIt = neighborPeerSet().begin();
551  for (; peerIt != peerEndIt; ++peerIt) {
552  ProcessRank neighborPeer = *peerIt;
553  auto& numIndicesRcvBuf = numIndicesRcvBufs[neighborPeer];
554  auto& indicesRcvBuf = indicesRcvBufs[neighborPeer];
555 
556  numIndicesRcvBuf.resize(1);
557  numIndicesRcvBuf.receive(neighborPeer);
558  unsigned numIndices = numIndicesRcvBufs[neighborPeer][0];
559  indicesRcvBuf.resize(numIndices);
560  indicesRcvBuf.receive(neighborPeer);
561 
562  // filter out all indices which are already in the peer
563  // processes' overlap and add them to the seed list. also
564  // extend the set of peer processes.
565  for (unsigned i = 0; i < numIndices; ++i) {
566  // swap the local and the peer indices, because they were
567  // created with the point view of the sender
568  std::swap(indicesRcvBuf[i].localIdx, indicesRcvBuf[i].peerIdx);
569 
570  ProcessRank peerRank = indicesRcvBuf[i].peerRank;
571  // Index peerIdx = indicesRcvBuf[i].peerIdx;
572  Index localIdx = indicesRcvBuf[i].localIdx;
573 
574  // check if the index is already in the overlap for
575  // the peer
576  const auto& distIt = foreignOverlapByLocalIndex_[static_cast<unsigned>(localIdx)].find(peerRank);
577  if (distIt != foreignOverlapByLocalIndex_[static_cast<unsigned>(localIdx)].end())
578  continue;
579 
580  // make sure the index is not already in the seed list
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) {
586  inSeedList = true;
587  break;
588  }
589  }
590  if (inSeedList)
591  continue;
592 
593  IndexRankDist seedEntry;
594  seedEntry.index = localIdx;
595  seedEntry.peerRank = peerRank;
596  seedEntry.borderDistance = borderDist;
597  seedList.push_back(seedEntry);
598 
599  // update the peer set
600  peerSet_.insert(peerRank);
601  }
602  }
603 
604  // make sure all data was send
605  peerIt = neighborPeerSet().begin();
606  for (; peerIt != peerEndIt; ++peerIt) {
607  ProcessRank neighborPeer = *peerIt;
608  numIndicesSendBufs[neighborPeer].wait();
609  indicesSendBufs[neighborPeer].wait();
610  }
611 #endif // HAVE_MPI
612  }
613 
614  // given a list of border indices and provided that
615  // borderListToSeedList_() was already called, calculate the
616  // master process of each local index.
617  void computeMasterRanks_()
618  {
619  // determine the minimum rank for all indices
620  masterRank_.resize(numLocal_);
621  for (unsigned localIdx = 0; localIdx < numLocal_; ++localIdx) {
622  unsigned masterRank = myRank_;
623  if (isBorder(static_cast<Index>(localIdx))) {
624  // if the local index is a border index, loop over all ranks
625  // for which this index is also a border index. the lowest
626  // rank wins!
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) {
631  // if the border distance is zero, the rank with the
632  // minimum
633  masterRank = std::min<ProcessRank>(masterRank, it->first);
634  }
635  }
636  }
637  masterRank_[static_cast<unsigned>(localIdx)] = masterRank;
638  }
639  }
640 
641  // assuming that the foreign overlap has been created for each
642  // local index, this method groups the foreign overlap by peer
643  // process rank
644  void groupForeignOverlapByRank_()
645  {
646  // loop over all indices which are in the overlap of some
647  // process
648  size_t numLocal = foreignOverlapByLocalIndex_.size();
649  for (unsigned localIdx = 0; localIdx < numLocal; ++localIdx) {
650  // loop over the list of processes for the current index
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);
660  }
661  }
662  }
663 
664  // set of processes with which we have to communicate
665  PeerSet peerSet_;
666 
667  // set of processes which are direct neighbors of us
668  PeerSet neighborPeerSet_;
669 
670  // the list of indices on the border
671  const BorderList& borderList_;
672 
673  // the set of indices which should not be considered
674  const BlackList& blackList_;
675 
676  // local indices are the native indices sans the black listed ones
677  std::vector<Index> nativeToLocalIndices_;
678  std::vector<Index> localToNativeIndices_;
679 
680  // an array which contains the rank of the master process for each
681  // index
682  std::vector<ProcessRank> masterRank_;
683 
684  // set of all local indices which are on the border of some remote
685  // process
686  std::set<Index> localBorderIndices_;
687 
688  // stores the set of process ranks which are in the overlap for a
689  // given row index "owned" by the current rank. The second value
690  // store the distance from the nearest process border.
691  OverlapByIndex foreignOverlapByLocalIndex_;
692 
693  // stores a list of foreign overlap indices for each rank
694  OverlapByRank foreignOverlapByRank_;
695 
696  // size of the overlap region
697  unsigned overlapSize_;
698 
699  // number of local indices
700  size_t numLocal_;
701 
702  // number of native indices
703  size_t numNative_;
704 
705  // the MPI rank of the local process
706  ProcessRank myRank_;
707 };
708 
709 } // namespace Linear
710 } // namespace Ewoms
711 
712 #endif
bool iAmMasterOf(Index localIdx) const
Return true if the current rank is the &quot;master&quot; 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