27 #ifndef EWOMS_OVERLAPPING_BLOCK_VECTOR_HH 28 #define EWOMS_OVERLAPPING_BLOCK_VECTOR_HH 33 #include <opm/common/Valgrind.hpp> 35 #include <dune/istl/bvector.hh> 36 #include <dune/common/fvector.hh> 48 template <
class FieldVector,
class Overlap>
51 typedef Dune::BlockVector<FieldVector> ParentType;
52 typedef Dune::BlockVector<FieldVector> BlockVector;
60 : ParentType(overlap.numDomestic()), overlap_(&overlap)
68 , numIndicesSendBuff_(obv.numIndicesSendBuff_)
69 , indicesSendBuff_(obv.indicesSendBuff_)
70 , indicesRecvBuff_(obv.indicesRecvBuff_)
71 , valuesSendBuff_(obv.valuesSendBuff_)
72 , valuesRecvBuff_(obv.valuesRecvBuff_)
73 , overlap_(obv.overlap_)
86 using ParentType::operator=;
94 ParentType::operator=(obv);
95 numIndicesSendBuff_ = obv.numIndicesSendBuff_;
96 indicesSendBuff_ = obv.indicesSendBuff_;
97 indicesRecvBuff_ = obv.indicesRecvBuff_;
98 valuesSendBuff_ = obv.valuesSendBuff_;
99 valuesRecvBuff_ = obv.valuesRecvBuff_;
100 overlap_ = obv.overlap_;
108 template <
class BlockVector>
111 size_t numDomestic = overlap_->numDomestic();
114 for (
unsigned domRowIdx = 0; domRowIdx < numDomestic; ++domRowIdx) {
115 Index nativeRowIdx = overlap_->domesticToNative(static_cast<Index>(domRowIdx));
116 if (nativeRowIdx < 0)
117 (*this)[domRowIdx] = 0.0;
119 (*
this)[domRowIdx] = nativeBlockVector[nativeRowIdx];
131 template <
class NativeBlockVector>
132 void assign(
const NativeBlockVector& nativeBlockVector)
134 Index numDomestic = overlap_->numDomestic();
137 for (
Index domRowIdx = 0; domRowIdx < numDomestic; ++domRowIdx) {
138 Index nativeRowIdx = overlap_->domesticToNative(domRowIdx);
139 if (nativeRowIdx < 0)
140 (*this)[
static_cast<unsigned>(domRowIdx)] = 0.0;
142 (*
this)[
static_cast<unsigned>(domRowIdx)] = nativeBlockVector[static_cast<unsigned>(nativeRowIdx)];
154 template <
class NativeBlockVector>
155 void assignTo(NativeBlockVector& nativeBlockVector)
const 158 size_t numNative = overlap_->numNative();
159 nativeBlockVector.resize(numNative);
160 for (
unsigned nativeRowIdx = 0; nativeRowIdx < numNative; ++nativeRowIdx) {
161 Index domRowIdx = overlap_->nativeToDomestic(static_cast<Index>(nativeRowIdx));
164 nativeBlockVector[nativeRowIdx] = 0.0;
166 nativeBlockVector[nativeRowIdx] = (*this)[
static_cast<unsigned>(domRowIdx)];
176 typename PeerSet::const_iterator peerIt;
177 typename PeerSet::const_iterator peerEndIt = overlap_->peerSet().end();
180 peerIt = overlap_->peerSet().begin();
181 for (; peerIt != peerEndIt; ++peerIt) {
183 sendEntries_(peerRank);
187 peerIt = overlap_->peerSet().begin();
188 for (; peerIt != peerEndIt; ++peerIt) {
190 receiveFromMaster_(peerRank);
203 typename PeerSet::const_iterator peerIt;
204 typename PeerSet::const_iterator peerEndIt = overlap_->peerSet().end();
207 peerIt = overlap_->peerSet().begin();
208 for (; peerIt != peerEndIt; ++peerIt) {
210 sendEntries_(peerRank);
214 peerIt = overlap_->peerSet().begin();
215 for (; peerIt != peerEndIt; ++peerIt) {
217 receiveAdd_(peerRank);
230 typename PeerSet::const_iterator peerIt;
231 typename PeerSet::const_iterator peerEndIt = overlap_->peerSet().end();
234 peerIt = overlap_->peerSet().begin();
235 for (; peerIt != peerEndIt; ++peerIt) {
237 sendEntries_(peerRank);
241 peerIt = overlap_->peerSet().begin();
242 for (; peerIt != peerEndIt; ++peerIt) {
244 receiveAddBorder_(peerRank);
254 for (
unsigned i = 0; i < this->size(); ++i) {
255 std::cout <<
"row " << i << (overlap_->isLocal(i) ?
" " :
"*")
256 <<
": " << (*
this)[i] <<
"\n" << std::flush;
261 void createBuffers_()
265 typename PeerSet::const_iterator peerIt;
266 typename PeerSet::const_iterator peerEndIt = overlap_->peerSet().end();
269 peerIt = overlap_->peerSet().begin();
270 for (; peerIt != peerEndIt; ++peerIt) {
271 ProcessRank peerRank = *peerIt;
273 size_t numEntries = overlap_->foreignOverlapSize(peerRank);
274 numIndicesSendBuff_[peerRank] = std::make_shared<MpiBuffer<unsigned> >(1);
275 indicesSendBuff_[peerRank] = std::make_shared<MpiBuffer<Index> >(numEntries);
276 valuesSendBuff_[peerRank] = std::make_shared<MpiBuffer<FieldVector> >(numEntries);
279 MpiBuffer<Index>& indicesSendBuff = *indicesSendBuff_[peerRank];
280 for (
unsigned i = 0; i < numEntries; ++i) {
281 Index domRowIdx = overlap_->foreignOverlapOffsetToDomesticIdx(peerRank, i);
282 indicesSendBuff[i] = overlap_->domesticToGlobal(domRowIdx);
286 (*numIndicesSendBuff_[peerRank])[0] = static_cast<unsigned>(numEntries);
287 numIndicesSendBuff_[peerRank]->send(peerRank);
290 indicesSendBuff.send(peerRank);
294 peerIt = overlap_->peerSet().begin();
295 for (; peerIt != peerEndIt; ++peerIt) {
299 MpiBuffer<unsigned> numRowsRecvBuff(1);
300 numRowsRecvBuff.receive(peerRank);
301 unsigned numRows = numRowsRecvBuff[0];
304 indicesRecvBuff_[peerRank] = std::shared_ptr<MpiBuffer<Index> >(
305 new MpiBuffer<Index>(numRows));
306 valuesRecvBuff_[peerRank] = std::shared_ptr<MpiBuffer<FieldVector> >(
307 new MpiBuffer<FieldVector>(numRows));
308 MpiBuffer<Index>& indicesRecvBuff = *indicesRecvBuff_[peerRank];
311 indicesRecvBuff.receive(peerRank);
314 for (
unsigned i = 0; i != numRows; ++i) {
315 Index globalRowIdx = indicesRecvBuff[i];
316 Index domRowIdx = overlap_->globalToDomestic(globalRowIdx);
318 indicesRecvBuff[i] = domRowIdx;
323 peerIt = overlap_->peerSet().begin();
324 for (; peerIt != peerEndIt; ++peerIt) {
326 numIndicesSendBuff_[peerRank]->wait();
327 indicesSendBuff_[peerRank]->wait();
331 MpiBuffer<Index>& indicesSendBuff = *indicesSendBuff_[peerRank];
332 for (
unsigned i = 0; i < indicesSendBuff.size(); ++i) {
333 indicesSendBuff[i] = overlap_->globalToDomestic(indicesSendBuff[i]);
339 void sendEntries_(ProcessRank peerRank)
342 const MpiBuffer<Index>& indices = *indicesSendBuff_[peerRank];
343 MpiBuffer<FieldVector>& values = *valuesSendBuff_[peerRank];
344 for (
unsigned i = 0; i < indices.size(); ++i)
345 values[i] = (*
this)[
static_cast<unsigned>(indices[i])];
347 values.send(peerRank);
350 void waitSendFinished_()
352 typename PeerSet::const_iterator peerIt;
353 typename PeerSet::const_iterator peerEndIt = overlap_->peerSet().end();
356 peerIt = overlap_->peerSet().begin();
357 for (; peerIt != peerEndIt; ++peerIt) {
359 valuesSendBuff_[peerRank]->wait();
363 void receiveFromMaster_(ProcessRank peerRank)
365 const MpiBuffer<Index>& indices = *indicesRecvBuff_[peerRank];
366 MpiBuffer<FieldVector>& values = *valuesRecvBuff_[peerRank];
369 values.receive(peerRank);
372 for (
unsigned j = 0; j < indices.size(); ++j) {
373 Index domRowIdx = indices[j];
374 if (overlap_->masterRank(domRowIdx) == peerRank) {
375 (*this)[
static_cast<unsigned>(domRowIdx)] = values[j];
380 void receiveAddBorder_(ProcessRank peerRank)
382 const MpiBuffer<Index>& indices = *indicesRecvBuff_[peerRank];
383 MpiBuffer<FieldVector>& values = *valuesRecvBuff_[peerRank];
386 values.receive(peerRank);
389 for (
unsigned j = 0; j < indices.size(); ++j) {
390 Index domRowIdx = indices[j];
391 if (overlap_->isBorderWith(domRowIdx, peerRank))
392 (*
this)[
static_cast<unsigned>(domRowIdx)] += values[j];
394 (*
this)[
static_cast<unsigned>(domRowIdx)] = values[j];
398 void receiveAdd_(ProcessRank peerRank)
400 const MpiBuffer<Index>& indices = *indicesRecvBuff_[peerRank];
401 MpiBuffer<FieldVector>& values = *valuesRecvBuff_[peerRank];
404 values.receive(peerRank);
407 for (
unsigned j = 0; j < indices.size(); ++j) {
408 Index domRowIdx = indices[j];
409 (*this)[
static_cast<unsigned>(domRowIdx)] += values[j];
413 std::map<ProcessRank, std::shared_ptr<MpiBuffer<unsigned> > > numIndicesSendBuff_;
414 std::map<ProcessRank, std::shared_ptr<MpiBuffer<Index> > > indicesSendBuff_;
415 std::map<ProcessRank, std::shared_ptr<MpiBuffer<Index> > > indicesRecvBuff_;
416 std::map<ProcessRank, std::shared_ptr<MpiBuffer<FieldVector> > > valuesSendBuff_;
417 std::map<ProcessRank, std::shared_ptr<MpiBuffer<FieldVector> > > valuesRecvBuff_;
419 const Overlap *overlap_;
Definition: baseauxiliarymodule.hh:37
Simplifies handling of buffers to be used in conjunction with MPI.
An overlap aware block vector.
Definition: overlappingblockvector.hh:49
void syncAddBorder()
Syncronize all values of the block vector from the master rank, but add up the entries on the border...
Definition: overlappingblockvector.hh:228
OverlappingBlockVector(const OverlappingBlockVector &obv)
Copy constructor.
Definition: overlappingblockvector.hh:66
OverlappingBlockVector & operator=(const OverlappingBlockVector &obv)
Assignment operator.
Definition: overlappingblockvector.hh:92
int Index
The type of an index of a degree of freedom.
Definition: overlaptypes.hh:44
OverlappingBlockVector(const Overlap &overlap)
Given a domestic overlap object, create an overlapping block vector coherent to it.
Definition: overlappingblockvector.hh:59
void sync()
Syncronize all values of the block vector from their master process.
Definition: overlappingblockvector.hh:174
This files provides several data structures for storing tuples of indices of remote and/or local proc...
void assignTo(NativeBlockVector &nativeBlockVector) const
Assign the local values to a non-overlapping block vector.
Definition: overlappingblockvector.hh:155
void assign(const NativeBlockVector &nativeBlockVector)
Assign an overlapping block vector from a non-overlapping one, border entries are assigned using thei...
Definition: overlappingblockvector.hh:132
OverlappingBlockVector()
Default constructor.
Definition: overlappingblockvector.hh:79
unsigned ProcessRank
The type of the rank of a process.
Definition: overlaptypes.hh:49
void assignAddBorder(const BlockVector &nativeBlockVector)
Assign an overlapping block vector from a non-overlapping one, border entries are added...
Definition: overlappingblockvector.hh:109
void syncAdd()
Syncronize all values of the block vector by adding up the values of all peer ranks.
Definition: overlappingblockvector.hh:201