25 #ifndef EWOMS_OVERLAPPING_BLOCK_VECTOR_HH
26 #define EWOMS_OVERLAPPING_BLOCK_VECTOR_HH
31 #include <opm/material/common/Valgrind.hpp>
33 #include <dune/istl/bvector.hh>
34 #include <dune/common/fvector.hh>
46 template <
class FieldVector,
class Overlap>
49 typedef Dune::BlockVector<FieldVector> ParentType;
50 typedef Dune::BlockVector<FieldVector> BlockVector;
58 : ParentType(overlap.numDomestic()), overlap_(&overlap)
66 , numIndicesSendBuff_(obv.numIndicesSendBuff_)
67 , indicesSendBuff_(obv.indicesSendBuff_)
68 , indicesRecvBuff_(obv.indicesRecvBuff_)
69 , valuesSendBuff_(obv.valuesSendBuff_)
70 , valuesRecvBuff_(obv.valuesRecvBuff_)
71 , overlap_(obv.overlap_)
84 using ParentType::operator=;
92 ParentType::operator=(obv);
93 numIndicesSendBuff_ = obv.numIndicesSendBuff_;
94 indicesSendBuff_ = obv.indicesSendBuff_;
95 indicesRecvBuff_ = obv.indicesRecvBuff_;
96 valuesSendBuff_ = obv.valuesSendBuff_;
97 valuesRecvBuff_ = obv.valuesRecvBuff_;
98 overlap_ = obv.overlap_;
108 int numDomestic = overlap_->numDomestic();
111 for (
int domRowIdx = 0; domRowIdx < numDomestic; ++domRowIdx) {
112 int nativeRowIdx = overlap_->domesticToNative(domRowIdx);
113 if (nativeRowIdx < 0)
114 (*this)[domRowIdx] = 0.0;
116 (*
this)[domRowIdx] = nativeBlockVector[nativeRowIdx];
128 void assign(
const BlockVector &nativeBlockVector)
130 int numDomestic = overlap_->numDomestic();
133 for (
int domRowIdx = 0; domRowIdx < numDomestic; ++domRowIdx) {
134 int nativeRowIdx = overlap_->domesticToNative(domRowIdx);
135 if (nativeRowIdx < 0)
136 (*this)[domRowIdx] = 0.0;
138 (*
this)[domRowIdx] = nativeBlockVector[nativeRowIdx];
150 void assignTo(BlockVector &nativeBlockVector)
const
153 int numNative = overlap_->numNative();
154 nativeBlockVector.resize(numNative);
155 for (
int nativeRowIdx = 0; nativeRowIdx < numNative; ++nativeRowIdx) {
156 int domRowIdx = overlap_->nativeToDomestic(nativeRowIdx);
159 nativeBlockVector[nativeRowIdx] = 0.0;
161 nativeBlockVector[nativeRowIdx] = (*this)[domRowIdx];
171 typename PeerSet::const_iterator peerIt;
172 typename PeerSet::const_iterator peerEndIt = overlap_->peerSet().end();
175 peerIt = overlap_->peerSet().begin();
176 for (; peerIt != peerEndIt; ++peerIt) {
177 int peerRank = *peerIt;
178 sendEntries_(peerRank);
182 peerIt = overlap_->peerSet().begin();
183 for (; peerIt != peerEndIt; ++peerIt) {
184 int peerRank = *peerIt;
185 receiveFromMaster_(peerRank);
198 typename PeerSet::const_iterator peerIt;
199 typename PeerSet::const_iterator peerEndIt = overlap_->peerSet().end();
202 peerIt = overlap_->peerSet().begin();
203 for (; peerIt != peerEndIt; ++peerIt) {
204 int peerRank = *peerIt;
205 sendEntries_(peerRank);
209 peerIt = overlap_->peerSet().begin();
210 for (; peerIt != peerEndIt; ++peerIt) {
211 int peerRank = *peerIt;
212 receiveAdd_(peerRank);
225 typename PeerSet::const_iterator peerIt;
226 typename PeerSet::const_iterator peerEndIt = overlap_->peerSet().end();
229 peerIt = overlap_->peerSet().begin();
230 for (; peerIt != peerEndIt; ++peerIt) {
231 int peerRank = *peerIt;
232 sendEntries_(peerRank);
236 peerIt = overlap_->peerSet().begin();
237 for (; peerIt != peerEndIt; ++peerIt) {
238 int peerRank = *peerIt;
239 receiveAddBorder_(peerRank);
249 for (
int i = 0; i < this->size(); ++i) {
250 std::cout <<
"row " << i << (overlap_->isLocal(i) ?
" " :
"*")
251 <<
": " << (*
this)[i] <<
"\n" << std::flush;
256 void createBuffers_()
260 typename PeerSet::const_iterator peerIt;
261 typename PeerSet::const_iterator peerEndIt = overlap_->peerSet().end();
264 peerIt = overlap_->peerSet().begin();
265 for (; peerIt != peerEndIt; ++peerIt) {
266 int peerRank = *peerIt;
268 int numEntries = overlap_->foreignOverlapSize(peerRank);
269 numIndicesSendBuff_[peerRank] = std::make_shared<MpiBuffer<int> >(1);
270 indicesSendBuff_[peerRank] = std::make_shared<MpiBuffer<Index> >(numEntries);
271 valuesSendBuff_[peerRank] = std::make_shared<MpiBuffer<FieldVector> >(numEntries);
275 for (
int i = 0; i < numEntries; ++i) {
276 int domRowIdx = overlap_->foreignOverlapOffsetToDomesticIdx(peerRank, i);
277 indicesSendBuff[i] = overlap_->domesticToGlobal(domRowIdx);
281 (*numIndicesSendBuff_[peerRank])[0] = numEntries;
282 numIndicesSendBuff_[peerRank]->send(peerRank);
285 indicesSendBuff.send(peerRank);
289 peerIt = overlap_->peerSet().begin();
290 for (; peerIt != peerEndIt; ++peerIt) {
291 int peerRank = *peerIt;
294 MpiBuffer<int> numRowsRecvBuff(1);
295 numRowsRecvBuff.receive(peerRank);
296 int numRows = numRowsRecvBuff[0];
299 indicesRecvBuff_[peerRank] = std::shared_ptr<MpiBuffer<Index> >(
300 new MpiBuffer<Index>(numRows));
301 valuesRecvBuff_[peerRank] = std::shared_ptr<MpiBuffer<FieldVector> >(
302 new MpiBuffer<FieldVector>(numRows));
303 MpiBuffer<Index> &indicesRecvBuff = *indicesRecvBuff_[peerRank];
306 indicesRecvBuff.receive(peerRank);
309 for (
int i = 0; i != numRows; ++i) {
310 int globalRowIdx = indicesRecvBuff[i];
311 int domRowIdx = overlap_->globalToDomestic(globalRowIdx);
313 indicesRecvBuff[i] = domRowIdx;
318 peerIt = overlap_->peerSet().begin();
319 for (; peerIt != peerEndIt; ++peerIt) {
320 int peerRank = *peerIt;
321 numIndicesSendBuff_[peerRank]->wait();
322 indicesSendBuff_[peerRank]->wait();
326 MpiBuffer<Index> &indicesSendBuff = *indicesSendBuff_[peerRank];
327 for (
unsigned i = 0; i < indicesSendBuff.size(); ++i) {
328 indicesSendBuff[i] = overlap_->globalToDomestic(indicesSendBuff[i]);
334 void sendEntries_(
int peerRank)
337 const MpiBuffer<Index> &indices = *indicesSendBuff_[peerRank];
338 MpiBuffer<FieldVector> &values = *valuesSendBuff_[peerRank];
339 for (
unsigned i = 0; i < indices.size(); ++i)
340 values[i] = (*
this)[indices[i]];
342 values.send(peerRank);
345 void waitSendFinished_()
347 typename PeerSet::const_iterator peerIt;
348 typename PeerSet::const_iterator peerEndIt = overlap_->peerSet().end();
351 peerIt = overlap_->peerSet().begin();
352 for (; peerIt != peerEndIt; ++peerIt) {
353 int peerRank = *peerIt;
354 valuesSendBuff_[peerRank]->wait();
358 void receiveFromMaster_(
int peerRank)
360 const MpiBuffer<Index> &indices = *indicesRecvBuff_[peerRank];
361 MpiBuffer<FieldVector> &values = *valuesRecvBuff_[peerRank];
364 values.receive(peerRank);
367 for (
unsigned j = 0; j < indices.size(); ++j) {
368 Index domRowIdx = indices[j];
369 if (overlap_->masterRank(domRowIdx) == peerRank) {
370 (*this)[domRowIdx] = values[j];
375 void receiveAddBorder_(
int peerRank)
377 const MpiBuffer<Index> &indices = *indicesRecvBuff_[peerRank];
378 MpiBuffer<FieldVector> &values = *valuesRecvBuff_[peerRank];
381 values.receive(peerRank);
384 for (
unsigned j = 0; j < indices.size(); ++j) {
385 int domRowIdx = indices[j];
386 if (overlap_->isBorderWith(domRowIdx, peerRank))
387 (*
this)[domRowIdx] += values[j];
389 (*
this)[domRowIdx] = values[j];
393 void receiveAdd_(
int peerRank)
395 const MpiBuffer<Index> &indices = *indicesRecvBuff_[peerRank];
396 MpiBuffer<FieldVector> &values = *valuesRecvBuff_[peerRank];
399 values.receive(peerRank);
402 for (
int j = 0; j < indices.size(); ++j) {
403 int domRowIdx = indices[j];
404 (*this)[domRowIdx] += values[j];
408 std::map<ProcessRank, std::shared_ptr<MpiBuffer<int> > > numIndicesSendBuff_;
409 std::map<ProcessRank, std::shared_ptr<MpiBuffer<Index> > > indicesSendBuff_;
410 std::map<ProcessRank, std::shared_ptr<MpiBuffer<Index> > > indicesRecvBuff_;
411 std::map<ProcessRank, std::shared_ptr<MpiBuffer<FieldVector> > > valuesSendBuff_;
412 std::map<ProcessRank, std::shared_ptr<MpiBuffer<FieldVector> > > valuesRecvBuff_;
414 const Overlap *overlap_;
void assignTo(BlockVector &nativeBlockVector) const
Assign the local values to a non-overlapping block vector.
Definition: overlappingblockvector.hh:150
This files provides several data structures for storing tuples of indices of remote and/or local proc...
Simplifies handling of buffers to be used in conjunction with MPI.
OverlappingBlockVector(const OverlappingBlockVector &obv)
Copy constructor.
Definition: overlappingblockvector.hh:64
An overlap aware block vector.
Definition: overlappingblockvector.hh:47
void print() const
Definition: overlappingblockvector.hh:247
OverlappingBlockVector()
Default constructor.
Definition: overlappingblockvector.hh:77
void syncAdd()
Syncronize all values of the block vector by adding up the values of all peer ranks.
Definition: overlappingblockvector.hh:196
void assignAddBorder(const BlockVector &nativeBlockVector)
Assign an overlapping block vector from a non-overlapping one, border entries are added...
Definition: overlappingblockvector.hh:106
Definition: baseauxiliarymodule.hh:35
OverlappingBlockVector(const Overlap &overlap)
Given a domestic overlap object, create an overlapping block vector coherent to it.
Definition: overlappingblockvector.hh:57
OverlappingBlockVector & operator=(const OverlappingBlockVector &obv)
Assignment operator.
Definition: overlappingblockvector.hh:90
void sync()
Syncronize all values of the block vector from their master process.
Definition: overlappingblockvector.hh:169
int Index
The type of an index of a degree of freedom.
Definition: overlaptypes.hh:43
void syncAddBorder()
Syncronize all values of the block vector from the master rank, but add up the entries on the border...
Definition: overlappingblockvector.hh:223
void assign(const BlockVector &nativeBlockVector)
Assign an overlapping block vector from a non-overlapping one, border entries are assigned using thei...
Definition: overlappingblockvector.hh:128