All Classes Files Functions Variables Typedefs Enumerations Enumerator Macros Groups Pages
domesticoverlapfrombcrsmatrix.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_DOMESTIC_OVERLAP_FROM_BCRS_MATRIX_HH
28 #define EWOMS_DOMESTIC_OVERLAP_FROM_BCRS_MATRIX_HH
29 
31 #include "blacklist.hh"
32 #include "globalindices.hh"
33 
35 
36 #include <algorithm>
37 #include <limits>
38 #include <set>
39 #include <map>
40 #include <vector>
41 
42 namespace Ewoms {
43 namespace Linear {
44 
53 {
56 
57 public:
58  // overlaps should never be copied!
60 
65  template <class BCRSMatrix>
66  DomesticOverlapFromBCRSMatrix(const BCRSMatrix& A,
67  const BorderList& borderList,
68  const BlackList& blackList,
69  unsigned overlapSize)
70  : foreignOverlap_(A, borderList, blackList, overlapSize)
71  , blackList_(blackList)
72  , globalIndices_(foreignOverlap_)
73  {
74  myRank_ = 0;
75  worldSize_ = 1;
76 
77 #if HAVE_MPI
78  int tmp;
79  MPI_Comm_rank(MPI_COMM_WORLD, &tmp);
80  myRank_ = static_cast<ProcessRank>(tmp);
81  MPI_Comm_size(MPI_COMM_WORLD, &tmp);
82  worldSize_ = static_cast<unsigned>(tmp);
83 #endif // HAVE_MPI
84 
85  buildDomesticOverlap_();
86  updateMasterRanks_();
87  blackList_.updateNativeToDomesticMap(*this);
88 
89  setupDebugMapping_();
90  }
91 
92  void check() const
93  {
94 #ifndef NDEBUG
95  // check consistency of global indices
96  for (unsigned domIdx = 0; domIdx < numDomestic(); ++domIdx) {
97  assert(globalToDomestic(domesticToGlobal(domIdx)) == static_cast<Index>(domIdx));
98  }
99 #endif // NDEBUG
100 
101  // send the foreign overlap for which we are master to the
102  // peers
103  std::map<int, MpiBuffer<unsigned> *> sizeBufferMap;
104 
105  auto peerIt = peerSet_.begin();
106  const auto& peerEndIt = peerSet_.end();
107  for (; peerIt != peerEndIt; ++peerIt) {
108  auto& buffer = *(new MpiBuffer<unsigned>(1));
109  sizeBufferMap[*peerIt] = &buffer;
110  buffer[0] = foreignOverlap_.foreignOverlapWithPeer(*peerIt).size();
111  buffer.send(*peerIt);
112  }
113 
114  peerIt = peerSet_.begin();
115  for (; peerIt != peerEndIt; ++peerIt) {
116  MpiBuffer<unsigned> rcvBuffer(1);
117  rcvBuffer.receive(*peerIt);
118 
119  assert(rcvBuffer[0] == domesticOverlapWithPeer_.find(*peerIt)->second.size());
120  }
121 
122  peerIt = peerSet_.begin();
123  for (; peerIt != peerEndIt; ++peerIt) {
124  sizeBufferMap[*peerIt]->wait();
125  delete sizeBufferMap[*peerIt];
126  }
127  }
128 
132  ProcessRank myRank() const
133  { return myRank_; }
134 
138  unsigned worldSize() const
139  { return worldSize_; }
140 
145  const PeerSet& peerSet() const
146  { return peerSet_; }
147 
151  bool isBorder(Index domesticIdx) const
152  {
153  return isLocal(domesticIdx)
154  && foreignOverlap_.isBorder(mapExternalToInternal_(domesticIdx));
155  }
156 
161  bool isBorderWith(Index domesticIdx, ProcessRank peerRank) const
162  {
163  return isLocal(domesticIdx)
164  && foreignOverlap_.isBorderWith(mapExternalToInternal_(domesticIdx),
165  peerRank);
166  }
167 
172  size_t numFront(ProcessRank peerRank) const
173  { return foreignOverlap_.numFront(peerRank); }
174 
178  bool isFront(Index domesticIdx) const
179  {
180  if (isLocal(domesticIdx))
181  return false;
182  Index internalDomesticIdx = mapExternalToInternal_(domesticIdx);
183 
184  // check wether the border distance of the domestic overlap is
185  // maximal for the index
186  const auto& domOverlap = domesticOverlapByIndex_[internalDomesticIdx];
187  return domOverlap.size() > 0
188  && domOverlap.begin()->second == foreignOverlap_.overlapSize();
189  }
190 
194  const BlackList& blackList() const
195  { return blackList_; }
196 
201  size_t numPeers(Index domesticIdx) const
202  { return domesticOverlapByIndex_[mapExternalToInternal_(domesticIdx)].size(); }
203 
207  unsigned overlapSize() const
208  { return foreignOverlap_.overlapSize(); }
209 
217  size_t numNative() const
218  { return foreignOverlap_.numNative(); }
219 
226  size_t numLocal() const
227  { return foreignOverlap_.numLocal(); }
228 
236  size_t numDomestic() const
237  { return globalIndices_.numDomestic(); }
238 
245  bool isLocal(Index domesticIdx) const
246  { return mapExternalToInternal_(domesticIdx) < static_cast<Index>(numLocal()); }
247 
252  bool iAmMasterOf(Index domesticIdx) const
253  {
254  if (!isLocal(domesticIdx))
255  return false;
256  return foreignOverlap_.iAmMasterOf(mapExternalToInternal_(domesticIdx));
257  }
258 
262  ProcessRank masterRank(Index domesticIdx) const
263  { return masterRank_[static_cast<unsigned>(mapExternalToInternal_(domesticIdx))]; }
264 
268  void print() const
269  { globalIndices_.print(); }
270 
274  Index globalToDomestic(Index globalIdx) const
275  {
276  Index internalIdx = globalIndices_.globalToDomestic(globalIdx);
277  if (internalIdx < 0)
278  return -1;
279  return mapInternalToExternal_(internalIdx);
280  }
281 
285  Index domesticToGlobal(Index domIdx) const
286  { return globalIndices_.domesticToGlobal(mapExternalToInternal_(domIdx)); }
287 
291  Index domesticToNative(Index domIdx) const
292  {
293  Index internalIdx = mapExternalToInternal_(domIdx);
294  if (internalIdx >= static_cast<Index>(numLocal()))
295  return -1;
296  return foreignOverlap_.localToNative(internalIdx);
297  }
298 
302  Index nativeToDomestic(Index nativeIdx) const
303  {
304  Index localIdx = foreignOverlap_.nativeToLocal(nativeIdx);
305  if (localIdx < 0)
306  return localIdx;
307  return mapInternalToExternal_(localIdx);
308  }
309 
314  bool isInOverlap(Index domesticIdx) const
315  {
316  return !this->isLocal(domesticIdx)
317  || this->foreignOverlap_.isInOverlap(mapExternalToInternal_(domesticIdx));
318  }
319 
324  bool isFrontFor(ProcessRank peerRank, Index domesticIdx) const
325  {
326  Index internalIdx = mapExternalToInternal_(domesticIdx);
327  return this->foreignOverlap_.isFrontFor(peerRank, internalIdx);
328  }
329 
333  bool peerHasIndex(int peerRank, Index domesticIdx) const
334  {
335  return foreignOverlap_.peerHasIndex(peerRank,
336  mapExternalToInternal_(domesticIdx));
337  }
338 
343  size_t foreignOverlapSize(ProcessRank peerRank) const
344  { return foreignOverlap_.foreignOverlapWithPeer(peerRank).size(); }
345 
351  Index foreignOverlapOffsetToDomesticIdx(ProcessRank peerRank, unsigned overlapOffset) const
352  {
353  Index internalIdx =
354  foreignOverlap_.foreignOverlapWithPeer(peerRank)[overlapOffset].index;
355  return mapInternalToExternal_(internalIdx);
356  }
357 
362  size_t domesticOverlapSize(ProcessRank peerRank) const
363  { return domesticOverlapWithPeer_.at(peerRank).size(); }
364 
370  Index domesticOverlapOffsetToDomesticIdx(ProcessRank peerRank, Index overlapOffset) const
371  {
372  Index internalIdx = domesticOverlapWithPeer_.at(peerRank)[overlapOffset];
373  return mapInternalToExternal_(internalIdx);
374  }
375 
376 protected:
377  void buildDomesticOverlap_()
378  {
379  // copy the set of peers from the foreign overlap
380  peerSet_ = foreignOverlap_.peerSet();
381 
382  // resize the array which stores the number of peers for
383  // each entry.
384  domesticOverlapByIndex_.resize(numLocal());
385  borderDistance_.resize(numLocal(), 0);
386 
387  PeerSet::const_iterator peerIt;
388  PeerSet::const_iterator peerEndIt = peerSet_.end();
389 
390  // send the overlap indices to all peer processes
391  peerIt = peerSet_.begin();
392  for (; peerIt != peerEndIt; ++peerIt) {
393  ProcessRank peerRank = *peerIt;
394  sendIndicesToPeer_(peerRank);
395  }
396 
397  // receive our overlap from the processes to all peer processes
398  peerIt = peerSet_.begin();
399  for (; peerIt != peerEndIt; ++peerIt) {
400  ProcessRank peerRank = *peerIt;
401  receiveIndicesFromPeer_(peerRank);
402  }
403 
404  // wait until all send operations complete
405  peerIt = peerSet_.begin();
406  for (; peerIt != peerEndIt; ++peerIt) {
407  ProcessRank peerRank = *peerIt;
408  waitSendIndices_(peerRank);
409  }
410  }
411 
412  void updateMasterRanks_()
413  {
414  size_t nLocal = numLocal();
415  size_t nDomestic = numDomestic();
416  masterRank_.resize(nDomestic);
417 
418  // take the master ranks for the local indices from the
419  // foreign overlap
420  for (unsigned i = 0; i < nLocal; ++i) {
421  masterRank_[i] = foreignOverlap_.masterRank(static_cast<Index>(i));
422  }
423 
424  // for non-local indices, initially use INT_MAX as their master
425  // rank
426  for (size_t i = nLocal; i < nDomestic; ++i)
427  masterRank_[i] = std::numeric_limits<ProcessRank>::max();
428 
429  // for the non-local indices, take the peer process for which
430  // a given local index is in the interior
431  auto peerIt = peerSet_.begin();
432  const auto& peerEndIt = peerSet_.end();
433  for (; peerIt != peerEndIt; ++peerIt) {
434  const auto& overlapWithPeer = domesticOverlapWithPeer_.find(*peerIt)->second;
435 
436  auto idxIt = overlapWithPeer.begin();
437  const auto& idxEndIt = overlapWithPeer.end();
438  for (; idxIt != idxEndIt; ++idxIt) {
439  if (*idxIt >= 0 && foreignOverlap_.isLocal(*idxIt))
440  continue; // ignore border indices
441 
442  masterRank_[static_cast<unsigned>(*idxIt)] = std::min(masterRank_[static_cast<unsigned>(*idxIt)], *peerIt);
443  }
444  }
445  }
446 
447  void sendIndicesToPeer_(ProcessRank peerRank)
448  {
449 #if HAVE_MPI
450  const auto& foreignOverlap = foreignOverlap_.foreignOverlapWithPeer(peerRank);
451 
452  // first, send a message containing the number of additional
453  // indices stemming from the overlap (i.e. without the border
454  // indices)
455  size_t numIndices = foreignOverlap.size();
456  numIndicesSendBuffer_[peerRank] = new MpiBuffer<size_t>(1);
457  (*numIndicesSendBuffer_[peerRank])[0] = numIndices;
458  numIndicesSendBuffer_[peerRank]->send(peerRank);
459 
460  // create MPI buffers
461  indicesSendBuffer_[peerRank] = new MpiBuffer<IndexDistanceNpeers>(numIndices);
462 
463  // then send the additional indices themselfs
464  auto overlapIt = foreignOverlap.begin();
465  const auto& overlapEndIt = foreignOverlap.end();
466  for (unsigned i = 0; overlapIt != overlapEndIt; ++overlapIt, ++i) {
467  Index localIdx = overlapIt->index;
468  BorderDistance borderDistance = overlapIt->borderDistance;
469  size_t numPeers = foreignOverlap_.foreignOverlapByLocalIndex(localIdx).size();
470 
471  IndexDistanceNpeers tmp;
472  tmp.index = globalIndices_.domesticToGlobal(localIdx);
473  tmp.borderDistance = borderDistance;
474  tmp.numPeers = static_cast<unsigned>(numPeers);
475 
476  (*indicesSendBuffer_[peerRank])[i] = tmp;
477  }
478 
479  indicesSendBuffer_[peerRank]->send(peerRank);
480 #endif // HAVE_MPI
481  }
482 
483  void waitSendIndices_(ProcessRank peerRank)
484  {
485  numIndicesSendBuffer_[peerRank]->wait();
486  delete numIndicesSendBuffer_[peerRank];
487 
488  indicesSendBuffer_[peerRank]->wait();
489  delete indicesSendBuffer_[peerRank];
490  }
491 
492  void receiveIndicesFromPeer_(ProcessRank peerRank)
493  {
494 #if HAVE_MPI
495  // receive the number of additional indices
496  int numIndices = -1;
497  MpiBuffer<size_t> numIndicesRecvBuff(1);
498  numIndicesRecvBuff.receive(peerRank);
499  numIndices = static_cast<int>(numIndicesRecvBuff[0]);
500 
501  // receive the additional indices themselfs
502  MpiBuffer<IndexDistanceNpeers> recvBuff(static_cast<size_t>(numIndices));
503  recvBuff.receive(peerRank);
504  for (unsigned i = 0; i < static_cast<unsigned>(numIndices); ++i) {
505  Index globalIdx = recvBuff[i].index;
506  BorderDistance borderDistance = recvBuff[i].borderDistance;
507 
508  // if the index is not already known, add it to the
509  // domestic indices
510  if (!globalIndices_.hasGlobalIndex(globalIdx)) {
511  Index newDomesticIdx = static_cast<Index>(globalIndices_.numDomestic());
512  globalIndices_.addIndex(newDomesticIdx, globalIdx);
513 
514  size_t newSize = globalIndices_.numDomestic();
515  borderDistance_.resize(newSize, std::numeric_limits<int>::max());
516  domesticOverlapByIndex_.resize(newSize);
517  }
518 
519  // convert the global index into a domestic one
520  Index domesticIdx = globalIndices_.globalToDomestic(globalIdx);
521 
522  // extend the domestic overlap
523  domesticOverlapByIndex_[static_cast<unsigned>(domesticIdx)][static_cast<unsigned>(peerRank)] = borderDistance;
524  domesticOverlapWithPeer_[static_cast<unsigned>(peerRank)].push_back(domesticIdx);
525 
526  //assert(borderDistance >= 0);
527  assert(globalIdx >= 0);
528  assert(domesticIdx >= 0);
529  assert(!(borderDistance == 0 && !foreignOverlap_.isLocal(domesticIdx)));
530  assert(!(borderDistance > 0 && foreignOverlap_.isLocal(domesticIdx)));
531 
532  borderDistance_[static_cast<unsigned>(domesticIdx)] = std::min(borderDistance, borderDistance_[static_cast<unsigned>(domesticIdx)]);
533  }
534 #endif // HAVE_MPI
535  }
536 
537  // this method is intended to set up the code mapping code for
538  // mapping domestic indices to the same ones used by a sequential
539  // grid. this requires detailed knowledge about how a grid
540  // distributes the degrees of freedom over multiple processes, but
541  // it can simplify debugging considerably because the indices can
542  // be made identical for the parallel and the sequential
543  // computations.
544  //
545  // by default, this method does nothing
546  void setupDebugMapping_()
547  {}
548 
549  // this method is intended to map domestic indices to the ones
550  // used by a sequential grid.
551  //
552  // by default, this method does nothing
553  Index mapInternalToExternal_(Index internalIdx) const
554  { return internalIdx; }
555 
556  // this method is intended to map the indices used by a sequential
557  // to grid domestic indices ones.
558  //
559  // by default, this method does nothing
560  Index mapExternalToInternal_(Index externalIdx) const
561  { return externalIdx; }
562 
563  ProcessRank myRank_;
564  unsigned worldSize_;
565  ForeignOverlap foreignOverlap_;
566 
567  BlackList blackList_;
568 
569  DomesticOverlapByRank domesticOverlapWithPeer_;
570  OverlapByIndex domesticOverlapByIndex_;
571  std::vector<BorderDistance> borderDistance_;
572  std::vector<ProcessRank> masterRank_;
573 
574  std::map<ProcessRank, MpiBuffer<size_t> *> numIndicesSendBuffer_;
575  std::map<ProcessRank, MpiBuffer<IndexDistanceNpeers> *> indicesSendBuffer_;
576  GlobalIndices globalIndices_;
577  PeerSet peerSet_;
578 };
579 
580 } // namespace Linear
581 } // namespace Ewoms
582 
583 #endif
Simplifies handling of buffers to be used in conjunction with MPI.
Definition: mpibuffer.hh:45
bool iAmMasterOf(Index localIdx) const
Return true if the current rank is the &quot;master&quot; of an index.
Definition: foreignoverlapfrombcrsmatrix.hh:184
bool isFront(Index domesticIdx) const
Returns true iff a domestic index is on the front.
Definition: domesticoverlapfrombcrsmatrix.hh:178
Index domesticToGlobal(Index domIdx) const
Returns a global index given a domestic one.
Definition: domesticoverlapfrombcrsmatrix.hh:285
bool peerHasIndex(int peerRank, Index domesticIdx) const
Returns true iff a domestic index is seen by a peer rank.
Definition: domesticoverlapfrombcrsmatrix.hh:333
Index foreignOverlapOffsetToDomesticIdx(ProcessRank peerRank, unsigned overlapOffset) const
Returns the domestic index given an offset in the foreign overlap of a peer process with the local pr...
Definition: domesticoverlapfrombcrsmatrix.hh:351
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
size_t numDomestic() const
Returns the number domestic indices.
Definition: globalindices.hh:121
unsigned overlapSize() const
Returns the size of the overlap region.
Definition: domesticoverlapfrombcrsmatrix.hh:207
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
Index domesticOverlapOffsetToDomesticIdx(ProcessRank peerRank, Index overlapOffset) const
Returns the domestic index given an offset in the domestic overlap of a peer process with the local p...
Definition: domesticoverlapfrombcrsmatrix.hh:370
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 isBorderWith(Index domesticIdx, ProcessRank peerRank) const
Returns true iff a domestic index is on the border with a given peer process.
Definition: domesticoverlapfrombcrsmatrix.hh:161
void print() const
Print the foreign overlap for debugging purposes.
Definition: domesticoverlapfrombcrsmatrix.hh:268
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.
size_t numLocal() const
Returns the number local indices.
Definition: domesticoverlapfrombcrsmatrix.hh:226
DomesticOverlapFromBCRSMatrix(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: domesticoverlapfrombcrsmatrix.hh:66
Index globalToDomestic(Index globalIdx) const
Returns a domestic index given a global one.
Definition: domesticoverlapfrombcrsmatrix.hh:274
const BlackList & blackList() const
Returns the object which represents the black-listed native indices.
Definition: domesticoverlapfrombcrsmatrix.hh:194
size_t numNative() const
Returns the number native indices.
Definition: domesticoverlapfrombcrsmatrix.hh:217
bool isInOverlap(Index domesticIdx) const
Returns true if a given domestic index is either in the foreign or in the domestic overlap...
Definition: domesticoverlapfrombcrsmatrix.hh:314
This class maps domestic row indices to and from &quot;global&quot; indices which is used to construct an algeb...
Index nativeToDomestic(Index nativeIdx) const
Returns a domestic index given a native one.
Definition: domesticoverlapfrombcrsmatrix.hh:302
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
This class maps domestic row indices to and from &quot;global&quot; indices which is used to construct an algeb...
Definition: globalindices.hh:55
bool iAmMasterOf(Index domesticIdx) const
Return true iff the current process is the master of a given domestic index.
Definition: domesticoverlapfrombcrsmatrix.hh:252
size_t foreignOverlapSize(ProcessRank peerRank) const
Returns number of indices which are contained in the foreign overlap with a peer. ...
Definition: domesticoverlapfrombcrsmatrix.hh:343
unsigned overlapSize() const
Returns the size of the overlap region.
Definition: foreignoverlapfrombcrsmatrix.hh:145
void addIndex(Index domesticIdx, Index globalIdx)
Add an index to the domestic&lt;-&gt;global mapping.
Definition: globalindices.hh:127
bool hasGlobalIndex(Index globalIdx) const
Return true iff a given global index already exists.
Definition: globalindices.hh:181
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
bool isBorder(Index domesticIdx) const
Returns true iff a domestic index is a border index.
Definition: domesticoverlapfrombcrsmatrix.hh:151
unsigned worldSize() const
Returns the number of processes in the global MPI communicator.
Definition: domesticoverlapfrombcrsmatrix.hh:138
const PeerSet & peerSet() const
Return the set of process ranks which share an overlap with the current process.
Definition: domesticoverlapfrombcrsmatrix.hh:145
This class creates and manages the foreign overlap given an initial list of border indices and a BCRS...
size_t numDomestic() const
Returns the number domestic indices.
Definition: domesticoverlapfrombcrsmatrix.hh:236
Index domesticToGlobal(Index domesticIdx) const
Converts a domestic index to a global one.
Definition: globalindices.hh:88
ProcessRank myRank() const
Returns the rank of the current process.
Definition: domesticoverlapfrombcrsmatrix.hh:132
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 globalToDomestic(Index globalIdx) const
Converts a global index to a domestic one.
Definition: globalindices.hh:98
size_t numFront(ProcessRank peerRank) const
Returns the number of indices on the front within a given peer rank&#39;s grid partition.
Definition: domesticoverlapfrombcrsmatrix.hh:172
Index localToNative(Index localIdx) const
Convert a local index to a native one.
Definition: foreignoverlapfrombcrsmatrix.hh:302
This class creates and manages the foreign overlap given an initial list of border indices and a BCRS...
Definition: domesticoverlapfrombcrsmatrix.hh:52
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
Index domesticToNative(Index domIdx) const
Returns a native index given a domestic one.
Definition: domesticoverlapfrombcrsmatrix.hh:291
bool isLocal(Index domesticIdx) const
Return true if a domestic index is local for the process.
Definition: domesticoverlapfrombcrsmatrix.hh:245
size_t domesticOverlapSize(ProcessRank peerRank) const
Returns number of indices which are contained in the domestic overlap with a peer.
Definition: domesticoverlapfrombcrsmatrix.hh:362
size_t numPeers(Index domesticIdx) const
Returns the number of processes which &quot;see&quot; a given index.
Definition: domesticoverlapfrombcrsmatrix.hh:201
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
bool isBorder(Index localIdx) const
Returns true iff a local index is a border index.
Definition: foreignoverlapfrombcrsmatrix.hh:151
bool isFrontFor(ProcessRank peerRank, Index domesticIdx) const
Returns true if a given domestic index is a front index for a peer rank.
Definition: domesticoverlapfrombcrsmatrix.hh:324
ProcessRank masterRank(Index domesticIdx) const
Return the rank of a master process for a domestic index.
Definition: domesticoverlapfrombcrsmatrix.hh:262
void print() const
Prints the global indices of all domestic indices for debugging purposes.
Definition: globalindices.hh:188