22 #ifndef OPM_OUTPUT_DATA_INTERREGFLOWMAP_HPP
23 #define OPM_OUTPUT_DATA_INTERREGFLOWMAP_HPP
25 #include <opm/output/data/InterRegFlow.hpp>
30 #include <type_traits>
40 namespace Opm {
namespace data {
48 using Neighbours = std::vector<int>;
51 using Offset = Neighbours::size_type;
54 using Start = std::vector<Offset>;
57 using RateBuffer = std::vector<float>;
60 using Window = InterRegFlow<RateBuffer::iterator>;
113 std::optional<std::pair<ReadOnlyWindow, ReadOnlyWindow::ElmT>>
117 template <
class MessageBufferType>
118 void write(MessageBufferType& buffer)
const
120 this->csr_.write(buffer);
124 template <
class MessageBufferType>
125 void read(MessageBufferType& buffer)
131 .add(other.maxRowIdx(),
133 other.coordinateFormatRowIndices(),
134 other.columnIndices(),
155 void add(
const int r1,
const int r2,
const FlowRates& rates);
173 void add(
const int maxRowIdx,
175 const Neighbours& rows,
176 const Neighbours& cols,
177 const RateBuffer& rates);
189 bool isValid()
const;
198 Neighbours::size_type numContributions()
const;
201 const Neighbours& rowIndices()
const;
204 const Neighbours& columnIndices()
const;
207 const RateBuffer& values()
const;
253 void merge(
const Connections& conns,
264 std::optional<ReadOnlyWindow> getWindow(
const int i,
const int j)
const;
267 Offset numRows()
const;
270 int maxRowIdx()
const;
273 int maxColIdx()
const;
276 const Start& startPointers()
const;
280 const Neighbours& columnIndices()
const;
285 const RateBuffer& values()
const;
289 Neighbours coordinateFormatRowIndices()
const;
292 template <
class MessageBufferType>
293 void write(MessageBufferType& buffer)
const
295 this->writeVector(this->ia_, buffer);
296 this->writeVector(this->ja_, buffer);
297 this->writeVector(this->sa_, buffer);
298 this->writeVector(this->compressedIdx_, buffer);
300 buffer.write(this->numRows_);
301 buffer.write(this->numCols_);
305 template <
class MessageBufferType>
306 void read(MessageBufferType& buffer)
308 this->readVector(buffer, this->ia_);
309 this->readVector(buffer, this->ja_);
310 this->readVector(buffer, this->sa_);
311 this->readVector(buffer, this->compressedIdx_);
313 buffer.read(this->numRows_);
314 buffer.read(this->numCols_);
334 Start compressedIdx_{};
347 template <
typename T,
class A,
class MessageBufferType>
348 void writeVector(
const std::vector<T,A>& vec,
349 MessageBufferType& buffer)
const
351 const auto n = vec.size();
354 for (
const auto& x : vec) {
359 template <
class MessageBufferType,
typename T,
class A>
360 void readVector(MessageBufferType& buffer,
361 std::vector<T,A>& vec)
363 auto n = 0 * vec.size();
368 for (
auto& x : vec) {
396 void assemble(
const Neighbours& rows,
397 const Neighbours& cols,
399 const int maxColIdx);
416 const RateBuffer& rates);
424 void sortColumnIndicesPerRow();
435 void condenseDuplicates();
451 void accumulateFlowRates(
const RateBuffer& v);
469 void preparePushbackRowGrouping(
const int numRows,
470 const Neighbours& rowIdx);
483 void groupAndTrackColumnIndicesByRow(
const Neighbours& rowIdx,
484 const Neighbours& colIdx);
512 void condenseAndTrackUniqueColumnsForSingleRow(Neighbours::const_iterator begin,
513 Neighbours::const_iterator end);
520 void remapCompressedIndex(Start&& compressedIdx);
525 Connections uncompressed_;
Form CSR adjacency matrix representation of inter-region flow rate graph provided as a list of connec...
Definition: InterRegFlowMap.hpp:45
void compress(const std::size_t numRegions)
Form CSR adjacency matrix representation of input graph from connections established in previous call...
Offset numRegions() const
Retrieve number of rows (source entities) in input graph.
void clear()
Clear all internal buffers, but preserve allocated capacity.
void addConnection(const int r1, const int r2, const FlowRates &rates)
Add flow rate connection between regions.
Window::Component Component
Client type through which to identify a component flow of a single inter-region connection.
Definition: InterRegFlowMap.hpp:71
std::optional< std::pair< ReadOnlyWindow, ReadOnlyWindow::ElmT > > getInterRegFlows(const int r1, const int r2) const
Retrieve accumulated inter-region flow rates for identified pair of regions.
Window::FlowRates FlowRates
Client type through which to define a single inter-region connection.
Definition: InterRegFlowMap.hpp:67
InterRegFlow< std::vector< float >::const_iterator > ReadOnlyWindow
Client view of flows between specified region pair.
Definition: InterRegFlowMap.hpp:64
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition: Exceptions.hpp:29