24#ifndef OPM_COLLECT_DATA_ON_IO_RANK_IMPL_HPP
25#define OPM_COLLECT_DATA_ON_IO_RANK_IMPL_HPP
27#include <opm/simulators/flow/CollectDataOnIORank.hpp>
29#include <dune/common/version.hh>
30#include <dune/grid/common/gridenums.hh>
31#include <dune/grid/common/mcmgmapper.hh>
32#include <dune/grid/common/partitionset.hh>
34#include <opm/grid/common/CartesianIndexMapper.hpp>
43 std::vector<std::string> toVector(
const std::set<std::string>& string_set)
45 return { string_set.begin(), string_set.end() };
65 { isInterior_ =
false; }
67 void setId(
int globalId)
68 { globalId_ = globalId; }
69 void setIndex(
int localIndex)
70 { localIndex_ = localIndex; }
72 int localIndex ()
const
73 {
return localIndex_; }
76 bool isInterior()
const
77 {
return isInterior_; }
80using IndexMapType = std::vector<int>;
81using IndexMapStorageType = std::vector<IndexMapType>;
82using P2PCommunicatorType = Dune::Point2PointCommunicator<Dune::SimpleMessageBuffer>;
83using MessageBufferType =
typename P2PCommunicatorType::MessageBufferType;
87 const std::vector<int>& distributedGlobalIndex_;
88 IndexMapType& localIndexMap_;
89 IndexMapStorageType& indexMaps_;
90 std::map<int, int> globalPosition_;
92 std::vector<int>& ranks_;
99 std::vector<int>&
ranks,
112 for (std::size_t index = 0; index < size; ++index)
113 globalPosition_.insert(std::make_pair(
globalIndex[index], index));
116 if (!indexMaps_.empty()) {
118 ranks_.resize(size, -1);
120 IndexMapType&
indexMap = indexMaps_.back();
121 std::size_t
localSize = localIndexMap_.size();
123 for (std::size_t i = 0; i <
localSize; ++i)
125 int id = distributedGlobalIndex_[localIndexMap_[i]];
133 if (!indexMaps_.size())
138 auto rankIt = recv_.begin();
143 if (
rankIt != recv_.end())
152 ranks_[
entry] = std::max(ranks_[
entry], rank);
154 if (
rankIt != recv_.end())
158 for (
const auto& rank: ranks_)
164 void pack(
int link, MessageBufferType&
buffer)
168 throw std::logic_error(
"link in method pack is not 0 as execpted");
171 int size = localIndexMap_.size();
174 for (
int index = 0; index < size; ++index) {
175 int globalIdx = distributedGlobalIndex_[localIndexMap_[index]];
180 void unpack(
int link, MessageBufferType&
buffer)
184 assert(!globalPosition_.empty());
190 for (
int index = 0; index < numCells; ++index) {
198template<
class EquilMapper,
class Mapper>
205 using DataType =
int;
206 bool fixedSize(
int ,
int )
212 std::size_t size(
const T&)
216 template<
class B,
class T>
217 void gather(B&
buffer,
const T&
t)
219 buffer.write(sendMapper_.index(
t));
221 template<
class B,
class T>
222 void scatter(B&
buffer,
const T&
t, std::size_t)
224 buffer.read(elementIndices_[recvMapper_.index(
t)]);
227 bool contains(
int dim,
int codim)
229 return dim==3 &&
codim==0;
233 const Mapper& recvMapper_;
234 std::vector<int>& elementIndices_;
238template<
class Mapper>
245 using DataType =
int;
246 bool fixedSize(
int ,
int )
252 std::size_t size(
const T&)
256 template<
class B,
class T>
257 void gather(B&
buffer,
const T&
t)
259 buffer.write(elementIndices_[mapper_.index(
t)]);
261 template<
class B,
class T>
262 void scatter(B&
buffer,
const T&
t, std::size_t)
264 buffer.read(elementIndices_[mapper_.index(
t)]);
267 bool contains(
int dim,
int codim)
269 return dim==3 &&
codim==0;
272 const Mapper& mapper_;
273 std::vector<int>& elementIndices_;
278 const data::Solution& localCellData_;
279 data::Solution& globalCellData_;
281 const IndexMapType& localIndexMap_;
282 const IndexMapStorageType& indexMaps_;
286 data::Solution& globalCellData,
292 , globalCellData_(globalCellData)
298 for (
const auto&
pair : localCellData_) {
299 const std::string& key =
pair.first;
316 void pack(
int link, MessageBufferType&
buffer)
320 throw std::logic_error(
"link in method pack is not 0 as expected");
323 for (
const auto&
pair : localCellData_) {
324 const auto& data =
pair.second.data<
double>();
327 write(
buffer, localIndexMap_, data);
331 void doUnpack(
const IndexMapType&
indexMap, MessageBufferType&
buffer)
335 for (
auto&
pair : localCellData_) {
336 const std::string& key =
pair.first;
337 auto& data = globalCellData_.data<
double>(key);
345 void unpack(
int link, MessageBufferType&
buffer)
349 template <
class Vector>
350 void write(MessageBufferType&
buffer,
353 unsigned int offset = 0,
354 unsigned int stride = 1)
const
359 for (
unsigned int i=0; i<size; ++i)
367 template <
class Vector>
368 void read(MessageBufferType&
buffer,
371 unsigned int offset = 0,
372 unsigned int stride = 1)
const
374 unsigned int size = 0;
377 for (
unsigned int i=0; i<size; ++i) {
387 const data::Wells& localWellData_;
388 data::Wells& globalWellData_;
392 data::Wells& globalWellData,
395 , globalWellData_(globalWellData)
408 void pack(
int link, MessageBufferType&
buffer)
412 throw std::logic_error(
"link in method pack is not 0 as expected");
414 localWellData_.write(
buffer);
418 void unpack(
int , MessageBufferType&
buffer)
419 { globalWellData_.read(
buffer); }
425 const data::GroupAndNetworkValues& localGroupAndNetworkData_;
426 data::GroupAndNetworkValues& globalGroupAndNetworkData_;
430 data::GroupAndNetworkValues& globalGroupAndNetworkData,
433 , globalGroupAndNetworkData_(globalGroupAndNetworkData)
435 if (! isIORank) {
return; }
446 void pack(
int link, MessageBufferType&
buffer)
450 throw std::logic_error {
451 "link in method pack is not 0 as expected"
456 this->localGroupAndNetworkData_.write(
buffer);
460 void unpack(
int , MessageBufferType&
buffer)
461 { this->globalGroupAndNetworkData_.read(
buffer); }
466 const std::map<std::pair<std::string, int>,
double>& localBlockData_;
467 std::map<std::pair<std::string, int>,
double>& globalBlockValues_;
487 void pack(
int link, MessageBufferType&
buffer)
491 throw std::logic_error(
"link in method pack is not 0 as expected");
494 unsigned int size = localBlockData_.size();
496 for (
const auto&
map : localBlockData_) {
504 void unpack(
int , MessageBufferType&
buffer)
507 unsigned int size = 0;
509 for (std::size_t i = 0; i < size; ++i) {
516 globalBlockValues_[std::make_pair(name,
idx)] = data;
523 const data::WellBlockAveragePressures& localWBPData_;
524 data::WellBlockAveragePressures& globalWBPValues_;
546 void pack(
int link, MessageBufferType&
buffer)
550 throw std::logic_error {
551 "link (" + std::to_string(
link) +
552 ") in method pack() is not 0 as expected"
557 this->localWBPData_.write(
buffer);
562 MessageBufferType&
buffer)
564 this->globalWBPValues_.read(
buffer);
587 void pack(
int link, MessageBufferType&
buffer) {
589 throw std::logic_error(
"link in method pack is not 0 as expected");
590 this->local_.pack(
buffer);
593 void unpack(
int, MessageBufferType&
buffer) {
594 this->global_.unpack(
buffer);
604 const data::Aquifers& localAquiferData_;
605 data::Aquifers& globalAquiferData_;
609 data::Aquifers& globalAquiferData,
612 , globalAquiferData_(globalAquiferData)
625 void pack(
int link, MessageBufferType&
buffer)
629 throw std::logic_error(
"link in method pack is not 0 as expected");
631 int size = localAquiferData_.size();
633 for (
const auto& [key, data] : localAquiferData_) {
640 void unpack(
int , MessageBufferType&
buffer)
644 for (
int i = 0; i < size; ++i) {
647 data::AquiferData data;
650 auto&
aq = this->globalAquiferData_[key];
652 this->unpackCommon(data,
aq);
654 if (
auto const*
aquFet = data.typeData.get<data::AquiferType::Fetkovich>();
659 else if (
auto const*
aquCT = data.typeData.get<data::AquiferType::CarterTracy>();
662 this->unpackCarterTracy(*
aquCT,
aq);
664 else if (
auto const*
aquNum = data.typeData.get<data::AquiferType::Numerical>();
667 this->unpackNumericAquifer(*
aquNum,
aq);
673 void unpackCommon(
const data::AquiferData& data, data::AquiferData&
aq)
675 aq.aquiferID = std::max(
aq.aquiferID, data.aquiferID);
676 aq.pressure = std::max(
aq.pressure, data.pressure);
677 aq.initPressure = std::max(
aq.initPressure, data.initPressure);
678 aq.datumDepth = std::max(
aq.datumDepth, data.datumDepth);
679 aq.fluxRate += data.fluxRate;
680 aq.volume += data.volume;
683 void unpackAquFet(
const data::FetkovichData&
aquFet, data::AquiferData&
aq)
685 if (!
aq.typeData.is<data::AquiferType::Fetkovich>()) {
686 auto*
tData =
aq.typeData.create<data::AquiferType::Fetkovich>();
691 auto&
dst = *
aq.typeData.getMutable<data::AquiferType::Fetkovich>();
693 dst.initVolume = std::max(
dst.initVolume ,
src.initVolume);
694 dst.prodIndex = std::max(
dst.prodIndex ,
src.prodIndex);
695 dst.timeConstant = std::max(
dst.timeConstant,
src.timeConstant);
699 void unpackCarterTracy(
const data::CarterTracyData&
aquCT, data::AquiferData&
aq)
701 if (!
aq.typeData.is<data::AquiferType::CarterTracy>()) {
702 auto*
tData =
aq.typeData.create<data::AquiferType::CarterTracy>();
707 auto&
dst = *
aq.typeData.getMutable<data::AquiferType::CarterTracy>();
709 dst.timeConstant = std::max(
dst.timeConstant ,
src.timeConstant);
710 dst.influxConstant = std::max(
dst.influxConstant,
src.influxConstant);
711 dst.waterDensity = std::max(
dst.waterDensity ,
src.waterDensity);
712 dst.waterViscosity = std::max(
dst.waterViscosity,
src.waterViscosity);
714 dst.dimensionless_time = std::max(
dst.dimensionless_time ,
src.dimensionless_time);
715 dst.dimensionless_pressure = std::max(
dst.dimensionless_pressure,
src.dimensionless_pressure);
719 void unpackNumericAquifer(
const data::NumericAquiferData&
aquNum, data::AquiferData&
aq)
721 if (!
aq.typeData.is<data::AquiferType::Numerical>()) {
722 auto*
tData =
aq.typeData.create<data::AquiferType::Numerical>();
727 auto&
dst = *
aq.typeData.getMutable<data::AquiferType::Numerical>();
729 std::transform(
src.initPressure.begin(),
730 src.initPressure.end(),
731 dst.initPressure.begin(),
732 dst.initPressure.begin(),
733 [](
const double p0_1,
const double p0_2)
735 return std::max(p0_1, p0_2);
751 , globalInterRegFlows_(globalInterRegFlows)
753 if (! isIORank) {
return; }
764 void pack(
int link, MessageBufferType&
buffer)
768 throw std::logic_error {
769 "link in method pack is not 0 as expected"
778 void unpack(
int , MessageBufferType&
buffer)
784 const std::array<FlowsData<double>, 3>& localFlows_;
785 std::array<FlowsData<double>, 3>& globalFlows_;
794 if (! isIORank) {
return; }
804 void pack(
int link, MessageBufferType&
buffer)
807 throw std::logic_error(
"link in method pack is not 0 as expected");
808 for (
int i = 0; i < 3; ++i) {
809 const auto& name = localFlows_[i].name;
811 const std::size_t size = localFlows_[i].indices.size();
813 for (std::size_t j = 0; j < size; ++j) {
814 const auto&
nncIdx = localFlows_[i].indices[j];
816 const auto&
flows = localFlows_[i].values[j];
822 void unpack(
int , MessageBufferType&
buffer)
824 for (
int i = 0; i < 3; ++i) {
827 globalFlows_[i].name = name;
828 std::size_t size = 0;
830 for (
unsigned int j = 0; j < size; ++j) {
838 globalFlows_[i].values[
nncIdx] = data;
844template <
class Gr
id,
class EquilGr
id,
class Gr
idView>
851 : toIORankComm_(grid.comm())
855 if (!needsReordering && !isParallel())
858 const CollectiveCommunication& comm = grid.comm();
861 std::set<int> send,
recv;
863 typename std::is_same<Grid, EquilGrid>::type
isSameGrid;
865 typedef Dune::MultipleCodimMultipleGeomTypeMapper<GridView> ElementMapper;
872 sortedCartesianIdx_.push_back(
cartMapper.cartesianIndex(
idx));
875 std::sort(sortedCartesianIdx_.begin(), sortedCartesianIdx_.end());
883 const std::size_t
globalSize = equilGrid->leafGridView().size(0);
893 grid.scatterData(handle);
900 globalCartesianIndex_[elemIdx] = cartElemIdx;
903 for (
int i = 0; i < comm.size(); ++i) {
918 grid.scatterData(handle);
925 grid.communicate(handle, Dune::InteriorBorder_All_Interface,
926 Dune::ForwardCommunication);
929 localIndexMap_.clear();
930 const std::size_t
gridSize = grid.size(0);
943 assert(
elem.partitionType() == Dune::InteriorEntity);
945 localIndexMap_.push_back(elemIdx);
949 toIORankComm_.insertRequest(send,
recv);
953 indexMaps_.resize(comm.size());
967template <
class Gr
id,
class EquilGr
id,
class Gr
idView>
968void CollectDataOnIORank<Grid,EquilGrid,GridView>::
970 const std::map<std::pair<std::string, int>,
double>&
localBlockData,
981 globalCellData_ = {};
982 globalBlockData_.clear();
984 globalWellData_.clear();
985 globalWBPData_.values.clear();
986 globalGroupAndNetworkData_.clear();
987 globalAquiferData_.clear();
988 globalWellTestState_.clear();
989 this->globalInterRegFlows_.clear();
994 if(!needsReordering && !isParallel())
1000 this->globalCellData_,
1001 this->localIndexMap_,
1007 if (! isParallel()) {
1013 for (
int i = 0; i < 3; ++i) {
1015 globalFloresn_[i].resize(
sizeFlr);
1017 globalFlowsn_[i].resize(
sizeFlo);
1022 this->globalWellData_,
1028 this->globalGroupAndNetworkData_,
1034 this->globalBlockData_,
1046 this->globalWBPData_,
1052 this->globalAquiferData_,
1058 this->globalWellTestState_,
1064 this->globalInterRegFlows_,
1070 this->globalFlowsn_,
1076 this->globalFloresn_,
1094 toIORankComm_.barrier();
1100template <
class Gr
id,
class EquilGr
id,
class Gr
idView>
1101int CollectDataOnIORank<Grid,EquilGrid,GridView>::
1102localIdxToGlobalIdx(
unsigned localIdx)
const
1104 if (!isParallel()) {
1108 if (this->localIdxToGlobalIdx_.empty()) {
1109 throw std::logic_error(
"index map is not created on this rank");
1112 if (localIdx >= this->localIdxToGlobalIdx_.size()) {
1113 throw std::logic_error(
"local index is outside map range");
1116 return this->localIdxToGlobalIdx_[localIdx];
1119template <
class Gr
id,
class EquilGr
id,
class Gr
idView>
1120bool CollectDataOnIORank<Grid,EquilGrid,GridView>::
1121isCartIdxOnThisRank(
int cartIdx)
const
1123 if (! this->isParallel()) {
1127 assert (! needsReordering);
1129 auto candidate = std::lower_bound(this->sortedCartesianIdx_.begin(),
1130 this->sortedCartesianIdx_.end(),
1133 return (
candidate != sortedCartesianIdx_.end())
Definition CollectDataOnIORank.hpp:49
Definition CollectDataOnIORank.hpp:56
Definition CollectDataOnIORank_impl.hpp:85
Communication handle to scatter the global index.
Definition CollectDataOnIORank_impl.hpp:240
Communication handle to scatter the global index.
Definition CollectDataOnIORank_impl.hpp:200
Definition CollectDataOnIORank_impl.hpp:53
Inter-region flow accumulation maps for all region definition arrays.
Definition InterRegFlows.hpp:179
void read(MessageBufferType &buffer)
Reconstitute internal object representation from MPI message buffer.
Definition InterRegFlows.hpp:323
void write(MessageBufferType &buffer) const
Serialise internal representation to MPI message buffer.
Definition InterRegFlows.hpp:296
Definition CollectDataOnIORank_impl.hpp:603
Definition CollectDataOnIORank_impl.hpp:465
Definition CollectDataOnIORank_impl.hpp:277
Definition CollectDataOnIORank_impl.hpp:424
Definition CollectDataOnIORank_impl.hpp:522
Definition CollectDataOnIORank_impl.hpp:386
Definition CollectDataOnIORank_impl.hpp:569
Definition CollectDataOnIORank_impl.hpp:783
Definition CollectDataOnIORank_impl.hpp:742
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
Simple container for FLOWS data.
Definition FlowsData.hpp:32