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>
54template <
class ForeignOverlap>
59 using GlobalToDomesticMap = std::map<Index, Index>;
60 using DomesticToGlobalMap = std::map<Index, Index>;
75 mpiSize_ =
static_cast<size_t>(tmp);
82 buildGlobalIndices_();
100 const auto& tmp = globalToDomestic_.find(globalIdx);
102 if (tmp == globalToDomestic_.end())
113 {
return foreignOverlap_.numLocal(); }
122 {
return numDomestic_; }
131 numDomestic_ = domesticToGlobal_.size();
133 assert(domesticToGlobal_.size() == globalToDomestic_.size());
150 static_cast<int>(peerRank),
167 static_cast<int>(peerRank),
184 {
return globalToDomestic_.find(globalIdx) != globalToDomestic_.end(); }
192 std::cout <<
"(domestic index, global index, domestic->global->domestic)"
193 <<
" list for rank " << myRank_ <<
"\n";
198 std::cout <<
"\n" << std::flush;
204 void buildGlobalIndices_()
209 numDomestic_ = foreignOverlap_.numLocal();
223 static_cast<int>(myRank_ - 1),
232 for (
unsigned i = 0; i < foreignOverlap_.numLocal(); ++i) {
233 if (!foreignOverlap_.iAmMasterOf(
static_cast<Index
>(i)))
237 static_cast<Index
>(domesticOffset_ +
numMaster));
241 if (myRank_ < mpiSize_ - 1) {
248 static_cast<int>(myRank_ + 1),
253 typename PeerSet::const_iterator
peerIt;
254 typename PeerSet::const_iterator
peerEndIt = peerSet_().end();
256 peerIt = peerSet_().begin();
259 receiveBorderFrom_(*
peerIt);
263 peerIt = peerSet_().begin();
270 peerIt = peerSet_().begin();
273 receiveBorderFrom_(*
peerIt);
277 peerIt = peerSet_().begin();
285 void sendBorderTo_([[
maybe_unused]] ProcessRank peerRank)
290 BorderList::const_iterator
borderIt = borderList_().begin();
291 BorderList::const_iterator
borderEndIt = borderList_().end();
295 if (
borderPeer != peerRank || borderDistance != 0)
298 Index localIdx = foreignOverlap_.nativeToLocal(
borderIt->localIdx);
301 if (foreignOverlap_.iAmMasterOf(localIdx)) {
308 void receiveBorderFrom_([[
maybe_unused]] ProcessRank peerRank)
313 BorderList::const_iterator
borderIt = borderList_().begin();
314 BorderList::const_iterator
borderEndIt = borderList_().end();
318 if (
borderPeer != peerRank || borderDistance != 0)
323 if (localIdx >= 0 && foreignOverlap_.masterRank(localIdx) ==
borderPeer)
329 const PeerSet& peerSet_()
const
330 {
return foreignOverlap_.peerSet(); }
333 {
return foreignOverlap_.borderList(); }
340 const ForeignOverlap& foreignOverlap_;
342 GlobalToDomesticMap globalToDomestic_;
343 DomesticToGlobalMap domesticToGlobal_;
This class maps domestic row indices to and from "global" indices which is used to construct an algeb...
Definition globalindices.hh:56
void print() const
Prints the global indices of all domestic indices for debugging purposes.
Definition globalindices.hh:190
void addIndex(Index domesticIdx, Index globalIdx)
Add an index to the domestic<->global mapping.
Definition globalindices.hh:127
size_t numDomestic() const
Returns the number domestic indices.
Definition globalindices.hh:121
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
void receiveBorderIndex(ProcessRank peerRank)
Receive an index on the border from a remote process and add it the translation maps.
Definition globalindices.hh:160
Index globalToDomestic(Index globalIdx) const
Converts a global index to a domestic one.
Definition globalindices.hh:98
bool hasGlobalIndex(Index globalIdx) const
Return true iff a given global index already exists.
Definition globalindices.hh:183
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
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
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
This structure stores a local index on a peer process and a global index.
Definition overlaptypes.hh:70