My Project
Loading...
Searching...
No Matches
InterRegFlowMap.hpp
Go to the documentation of this file.
1/*
2 Copyright 2016 SINTEF ICT, Applied Mathematics.
3 Copyright 2016 Statoil ASA.
4 Copyright 2022 Equinor ASA
5
6 This file is part of the Open Porous Media Project (OPM).
7
8 OPM is free software: you can redistribute it and/or modify
9 it under the terms of the GNU General Public License as published by
10 the Free Software Foundation, either version 3 of the License, or
11 (at your option) any later version.
12
13 OPM is distributed in the hope that it will be useful,
14 but WITHOUT ANY WARRANTY; without even the implied warranty of
15 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16 GNU General Public License for more details.
17
18 You should have received a copy of the GNU General Public License
19 along with OPM. If not, see <http://www.gnu.org/licenses/>.
20*/
21
22#ifndef OPM_OUTPUT_DATA_INTERREGFLOWMAP_HPP
23#define OPM_OUTPUT_DATA_INTERREGFLOWMAP_HPP
24
25#include <opm/output/data/InterRegFlow.hpp>
26
27#include <cstddef>
28#include <optional>
29#include <utility>
30#include <type_traits>
31#include <vector>
32
39
40namespace Opm { namespace data {
41
45 {
46 private:
48 using Neighbours = std::vector<int>;
49
51 using Offset = Neighbours::size_type;
52
54 using Start = std::vector<Offset>;
55
57 using RateBuffer = std::vector<float>;
58
60 using Window = InterRegFlow<RateBuffer::iterator>;
61
62 public:
64 using ReadOnlyWindow = InterRegFlow<std::vector<float>::const_iterator>;
65
67 using FlowRates = Window::FlowRates;
68
71 using Component = Window::Component;
72
84 void addConnection(const int r1, const int r2, const FlowRates& rates);
85
94 void compress(const std::size_t numRegions);
95
99 Offset numRegions() const;
100
113 std::optional<std::pair<ReadOnlyWindow, ReadOnlyWindow::ElmT>>
114 getInterRegFlows(const int r1, const int r2) const;
115
116 // MessageBufferType API should be similar to Dune::MessageBufferIF
117 template <class MessageBufferType>
118 void write(MessageBufferType& buffer) const
119 {
120 this->csr_.write(buffer);
121 }
122
123 // MessageBufferType API should be similar to Dune::MessageBufferIF
124 template <class MessageBufferType>
125 void read(MessageBufferType& buffer)
126 {
127 auto other = CSR{};
128 other.read(buffer);
129
130 this->uncompressed_
131 .add(other.maxRowIdx(),
132 other.maxColIdx(),
133 other.coordinateFormatRowIndices(),
134 other.columnIndices(),
135 other.values());
136 }
137
139 void clear();
140
141 private:
144 class Connections
145 {
146 public:
155 void add(const int r1, const int r2, const FlowRates& rates);
156
173 void add(const int maxRowIdx,
174 const int maxColIdx,
175 const Neighbours& rows,
176 const Neighbours& cols,
177 const RateBuffer& rates);
178
180 void clear();
181
185 bool empty() const;
186
189 bool isValid() const;
190
192 int maxRow() const;
193
195 int maxCol() const;
196
198 Neighbours::size_type numContributions() const;
199
201 const Neighbours& rowIndices() const;
202
204 const Neighbours& columnIndices() const;
205
207 const RateBuffer& values() const;
208
209 private:
211 Neighbours i_{};
212
214 Neighbours j_{};
215
218 RateBuffer v_{};
219
221 int max_i_{ -1 };
222
224 int max_j_{ -1 };
225 };
226
233 class CSR
234 {
235 public:
253 void merge(const Connections& conns,
254 const Offset numRegions);
255
264 std::optional<ReadOnlyWindow> getWindow(const int i, const int j) const;
265
267 Offset numRows() const;
268
270 int maxRowIdx() const;
271
273 int maxColIdx() const;
274
276 const Start& startPointers() const;
277
280 const Neighbours& columnIndices() const;
281
285 const RateBuffer& values() const;
286
289 Neighbours coordinateFormatRowIndices() const;
290
291 // MessageBufferType API should be similar to Dune::MessageBufferIF
292 template <class MessageBufferType>
293 void write(MessageBufferType& buffer) const
294 {
295 this->writeVector(this->ia_, buffer);
296 this->writeVector(this->ja_, buffer);
297 this->writeVector(this->sa_, buffer);
298 this->writeVector(this->compressedIdx_, buffer);
299
300 buffer.write(this->numRows_);
301 buffer.write(this->numCols_);
302 }
303
304 // MessageBufferType API should be similar to Dune::MessageBufferIF
305 template <class MessageBufferType>
306 void read(MessageBufferType& buffer)
307 {
308 this->readVector(buffer, this->ia_);
309 this->readVector(buffer, this->ja_);
310 this->readVector(buffer, this->sa_);
311 this->readVector(buffer, this->compressedIdx_);
312
313 buffer.read(this->numRows_);
314 buffer.read(this->numCols_);
315 }
316
318 void clear();
319
320 private:
322 Start ia_{};
323
326 Neighbours ja_{};
327
331 RateBuffer sa_{};
332
334 Start compressedIdx_{};
335
337 int numRows_{ 0 };
338
341 int numCols_{ 0 };
342
343 // ---------------------------------------------------------
344 // Implementation of read()/write()
345 // ---------------------------------------------------------
346
347 template <typename T, class A, class MessageBufferType>
348 void writeVector(const std::vector<T,A>& vec,
349 MessageBufferType& buffer) const
350 {
351 const auto n = vec.size();
352 buffer.write(n);
353
354 for (const auto& x : vec) {
355 buffer.write(x);
356 }
357 }
358
359 template <class MessageBufferType, typename T, class A>
360 void readVector(MessageBufferType& buffer,
361 std::vector<T,A>& vec)
362 {
363 auto n = 0 * vec.size();
364 buffer.read(n);
365
366 vec.resize(n);
367
368 for (auto& x : vec) {
369 buffer.read(x);
370 }
371 }
372
373 // ---------------------------------------------------------
374 // Implementation of merge()
375 // ---------------------------------------------------------
376
396 void assemble(const Neighbours& rows,
397 const Neighbours& cols,
398 const int maxRowIdx,
399 const int maxColIdx);
400
401
415 void compress(const Offset numRegions,
416 const RateBuffer& rates);
417
424 void sortColumnIndicesPerRow();
425
435 void condenseDuplicates();
436
451 void accumulateFlowRates(const RateBuffer& v);
452
453 // ---------------------------------------------------------
454 // Implementation of assemble()
455 // ---------------------------------------------------------
456
469 void preparePushbackRowGrouping(const int numRows,
470 const Neighbours& rowIdx);
471
483 void groupAndTrackColumnIndicesByRow(const Neighbours& rowIdx,
484 const Neighbours& colIdx);
485
486 // ---------------------------------------------------------
487 // General utilities
488 // ---------------------------------------------------------
489
494 void transpose();
495
512 void condenseAndTrackUniqueColumnsForSingleRow(Neighbours::const_iterator begin,
513 Neighbours::const_iterator end);
514
520 void remapCompressedIndex(Start&& compressedIdx);
521 };
522
525 Connections uncompressed_;
526
528 CSR csr_;
529 };
530
531}} // namespace Opm::data
532
533#endif // OPM_OUTPUT_DATA_INTERREGFLOWMAP_HPP
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...
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.
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
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:30