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
28
29#include <cstddef>
30#include <optional>
31#include <utility>
32#include <type_traits>
33#include <vector>
34
41
42namespace Opm { namespace data {
43
47 {
48 private:
50 using Neighbours = std::vector<int>;
51
53 using Offset = Neighbours::size_type;
54
56 using Start = std::vector<Offset>;
57
59 using RateBuffer = std::vector<float>;
60
62 using Window = InterRegFlow<RateBuffer::iterator>;
63
64 public:
66 using ReadOnlyWindow = InterRegFlow<std::vector<float>::const_iterator>;
67
69 using FlowRates = Window::FlowRates;
70
73 using Component = Window::Component;
74
86 void addConnection(const int r1, const int r2, const FlowRates& rates);
87
96 void compress(const std::size_t numRegions);
97
101 Offset numRegions() const;
102
115 std::optional<std::pair<ReadOnlyWindow, ReadOnlyWindow::ElmT>>
116 getInterRegFlows(const int r1, const int r2) const;
117
118 // MessageBufferType API should be similar to Dune::MessageBufferIF
119 template <class MessageBufferType>
120 void write(MessageBufferType& buffer) const
121 {
122 this->connections_.write(buffer);
123 this->writeVector(this->rates_, buffer);
124 }
125
126 // MessageBufferType API should be similar to Dune::MessageBufferIF
127 template <class MessageBufferType>
128 void read(MessageBufferType& buffer)
129 {
130 this->connections_.read(buffer);
131
132 auto rates = RateBuffer{};
133 this->readVector(buffer, rates);
134 this->appendRates(rates);
135 }
136
138 void clear();
139
140 private:
141 // VertexID = int, TrackCompressedIdx = true.
142 using Graph = utility::CSRGraphFromCoordinates<int, true>;
143
144 Graph connections_{};
145 RateBuffer rates_{};
146
147 template <typename T, class A, class MessageBufferType>
148 void writeVector(const std::vector<T,A>& vec,
149 MessageBufferType& buffer) const
150 {
151 const auto n = vec.size();
152 buffer.write(n);
153
154 for (const auto& x : vec) {
155 buffer.write(x);
156 }
157 }
158
159 template <typename T, class A, class MessageBufferType>
160 void readVector(MessageBufferType& buffer,
161 std::vector<T,A>& vec)
162 {
163 auto n = 0 * vec.size();
164 buffer.read(n);
165
166 vec.resize(n);
167
168 for (auto& x : vec) {
169 buffer.read(x);
170 }
171 }
172
173 template <typename Rates>
174 void appendRates(const Rates& rates)
175 {
176 this->rates_.insert(this->rates_.end(), rates.begin(), rates.end());
177 }
178 };
179
180}} // namespace Opm::data
181
182#endif // OPM_OUTPUT_DATA_INTERREGFLOWMAP_HPP
Facility for converting collection of region ID pairs into a sparse (CSR) adjacency matrix representa...
Form CSR adjacency matrix representation of inter-region flow rate graph provided as a list of connec...
Definition InterRegFlowMap.hpp:47
Offset numRegions() const
Retrieve number of rows (source entities) in input graph.
Definition InterRegFlowMap.cpp:113
void compress(const std::size_t numRegions)
Form CSR adjacency matrix representation of input graph from connections established in previous call...
Definition InterRegFlowMap.cpp:76
void clear()
Clear all internal buffers, but preserve allocated capacity.
Definition InterRegFlowMap.cpp:170
void addConnection(const int r1, const int r2, const FlowRates &rates)
Add flow rate connection between regions.
Definition InterRegFlowMap.cpp:44
Window::Component Component
Client type through which to identify a component flow of a single inter-region connection.
Definition InterRegFlowMap.hpp:73
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.
Definition InterRegFlowMap.cpp:122
Window::FlowRates FlowRates
Client type through which to define a single inter-region connection.
Definition InterRegFlowMap.hpp:69
InterRegFlow< std::vector< float >::const_iterator > ReadOnlyWindow
Client view of flows between specified region pair.
Definition InterRegFlowMap.hpp:66
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition Exceptions.hpp:30