My Project
Loading...
Searching...
No Matches
EclInterRegFlows.hpp
Go to the documentation of this file.
1// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
2// vi: set et ts=4 sw=4 sts=4:
3/*
4 Copyright 2022 Equinor AS
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 ECL_INTERREG_FLOWS_MODULE_HPP
23#define ECL_INTERREG_FLOWS_MODULE_HPP
24
25#include <opm/output/data/InterRegFlowMap.hpp>
26
27#include <algorithm>
28#include <cstddef>
29#include <string>
30#include <utility>
31#include <vector>
32
38
39namespace Opm {
40
41 class EclInterRegFlowMap;
42
47 {
48 public:
50 struct Cell {
52 int activeIndex{-1};
53
56
58 bool isInterior{true};
59 };
60
61 friend class EclInterRegFlowMap;
62
66 explicit EclInterRegFlowMapSingleFIP(const std::vector<int>& region);
67
82 void addConnection(const Cell& source,
83 const Cell& destination,
84 const data::InterRegFlowMap::FlowRates& rates);
85
91 void compress();
92
94 void clear();
95
99 const data::InterRegFlowMap& getInterRegFlows() const;
100
102 std::size_t getLocalMaxRegionID() const;
103
112 bool assignGlobalMaxRegionID(const std::size_t regID);
113
122 template <class MessageBufferType>
124 {
125 // Note: this->region_ is *intentionally* omitted here.
126 buffer.write(this->maxGlobalRegionID_);
127 this->iregFlow_.write(buffer);
128 }
129
144 template <class MessageBufferType>
146 {
147 auto otherMaxRegionID = 0 * this->maxGlobalRegionID_;
149 this->maxGlobalRegionID_ = std::max(this->maxGlobalRegionID_, otherMaxRegionID);
150
151 this->iregFlow_.read(buffer);
152 this->isReadFromStream_ = true;
153 }
154
155 private:
157 std::vector<int> region_{};
158
160 std::size_t maxLocalRegionID_{0};
161
164 std::size_t maxGlobalRegionID_{0};
165
167 data::InterRegFlowMap iregFlow_{};
168
171 bool isReadFromStream_{false};
172
174 EclInterRegFlowMapSingleFIP() = default;
175 };
176
179 {
180 public:
186 std::string name;
187
189 std::reference_wrapper<const std::vector<int>> definition;
190 };
191
194
197
204 static EclInterRegFlowMap
205 createMapFromNames(std::vector<std::string> names);
206
213 explicit EclInterRegFlowMap(const std::size_t numCells,
214 const std::vector<SingleRegion>& regions);
215
216 EclInterRegFlowMap(const EclInterRegFlowMap& rhs) = default;
217 EclInterRegFlowMap(EclInterRegFlowMap&& rhs) noexcept = default;
218
219 EclInterRegFlowMap& operator=(const EclInterRegFlowMap& rhs) = default;
220 EclInterRegFlowMap& operator=(EclInterRegFlowMap&& rhs) noexcept = default;
221
237 void addConnection(const Cell& source,
238 const Cell& destination,
239 const data::InterRegFlowMap::FlowRates& rates);
240
246 void compress();
247
249 void clear();
250
254 const std::vector<std::string>& names() const;
255
259 std::vector<data::InterRegFlowMap> getInterRegFlows() const;
260
262 std::vector<std::size_t> getLocalMaxRegionID() const;
263
272 bool assignGlobalMaxRegionID(const std::vector<std::size_t>& regID);
273
275 bool readIsConsistent() const;
276
277 bool wantInterRegflowSummary() const
278 {
279 return !this->regionMaps_.empty();
280 }
281
290 template <class MessageBufferType>
292 {
293 buffer.write(this->names_.size());
294 for (const auto& name : this->names_) {
295 buffer.write(name);
296 }
297
298 for (const auto& map : this->regionMaps_) {
299 map.write(buffer);
300 }
301 }
302
317 template <class MessageBufferType>
319 {
320 const auto names = this->readNames(buffer);
321
322 if (names == this->names_) {
323 // Input stream holds the maps in expected order (common
324 // case). Read the maps and merge with internal values.
325 for (auto& map : this->regionMaps_) {
326 map.read(buffer);
327 }
328 }
329 else {
330 // Input stream does not hold the maps in expected order (or
331 // different number of maps). Unexpected. Read the values
332 // from the input stream, but do not merge with internal
333 // values.
335 std::vector<int>(this->numCells_, 1)
336 };
337
338 const auto numMaps = names.size();
339 for (auto n = 0*numMaps; n < numMaps; ++n) {
340 map.read(buffer);
341 }
342
343 this->readIsConsistent_ = false;
344 }
345 }
346
347 private:
350 std::vector<EclInterRegFlowMapSingleFIP> regionMaps_{};
351
354 std::vector<std::string> names_;
355
358 std::size_t numCells_{0};
359
362 bool readIsConsistent_{true};
363
367 template <class MessageBufferType>
368 std::vector<std::string> readNames(MessageBufferType& buffer) const
369 {
370 auto numNames = 0 * this->names_.size();
371 buffer.read(numNames);
372
373 auto names = std::vector<std::string>(numNames);
374 for (auto name = 0*numNames; name < numNames; ++name) {
375 buffer.read(names[name]);
376 }
377
378 return names;
379 }
380 };
381} // namespace Opm
382
383#endif // ECL_INTERREG_FLOWS_MODULE_HPP
Definition AquiferInterface.hpp:35
Form CSR adjacency matrix representation of inter-region flow rate graph provided as a list of connec...
Definition EclInterRegFlows.hpp:47
void clear()
Clear all internal buffers, but preserve allocated capacity.
Definition EclInterRegFlows.cpp:83
void addConnection(const Cell &source, const Cell &destination, const data::InterRegFlowMap::FlowRates &rates)
Add flow rate connection between regions.
Definition EclInterRegFlows.cpp:47
const data::InterRegFlowMap & getInterRegFlows() const
Get read-only access to the underlying CSR representation.
Definition EclInterRegFlows.cpp:90
bool assignGlobalMaxRegionID(const std::size_t regID)
Assign maximum FIP region ID across all MPI ranks.
Definition EclInterRegFlows.cpp:102
std::size_t getLocalMaxRegionID() const
Retrieve maximum FIP region ID on local MPI rank.
Definition EclInterRegFlows.cpp:95
void compress()
Form CSR adjacency matrix representation of input graph from connections established in previous call...
Definition EclInterRegFlows.cpp:78
EclInterRegFlowMapSingleFIP(const std::vector< int > &region)
Constructor.
Definition EclInterRegFlows.cpp:32
void write(MessageBufferType &buffer) const
Serialise internal representation to MPI message buffer.
Definition EclInterRegFlows.hpp:123
void read(MessageBufferType &buffer)
Reconstitute internal object representation from MPI message buffer.
Definition EclInterRegFlows.hpp:145
Inter-region flow accumulation maps for all region definition arrays.
Definition EclInterRegFlows.hpp:179
static EclInterRegFlowMap createMapFromNames(std::vector< std::string > names)
Special purpose constructor for global object being collected on the I/O rank.
Definition EclInterRegFlows.cpp:119
bool assignGlobalMaxRegionID(const std::vector< std::size_t > &regID)
Assign maximum FIP region ID across all MPI ranks.
Definition EclInterRegFlows.cpp:205
std::vector< std::size_t > getLocalMaxRegionID() const
Retrieve maximum FIP region ID on local MPI rank.
Definition EclInterRegFlows.cpp:191
EclInterRegFlowMap()=default
Default constructor.
void read(MessageBufferType &buffer)
Reconstitute internal object representation from MPI message buffer.
Definition EclInterRegFlows.hpp:318
bool readIsConsistent() const
Whether or not previous read() operation succeeded.
Definition EclInterRegFlows.cpp:224
const std::vector< std::string > & names() const
Names of all applicable region definition arrays.
Definition EclInterRegFlows.cpp:172
void compress()
Form CSR adjacency matrix representation of input graph from connections established in previous call...
Definition EclInterRegFlows.cpp:155
std::vector< data::InterRegFlowMap > getInterRegFlows() const
Get read-only access to the underlying CSR representation.
Definition EclInterRegFlows.cpp:178
void addConnection(const Cell &source, const Cell &destination, const data::InterRegFlowMap::FlowRates &rates)
Add flow rate connection between regions for all region definitions.
Definition EclInterRegFlows.cpp:146
void write(MessageBufferType &buffer) const
Serialise internal representation to MPI message buffer.
Definition EclInterRegFlows.hpp:291
void clear()
Clear all internal buffers, but preserve allocated capacity.
Definition EclInterRegFlows.cpp:162
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition BlackoilPhases.hpp:27
Minimal characteristics of a cell from a simulation grid.
Definition EclInterRegFlows.hpp:50
int cartesianIndex
Cell's global cell ID.
Definition EclInterRegFlows.hpp:55
int activeIndex
Cell's active index on local rank.
Definition EclInterRegFlows.hpp:52
bool isInterior
Whether or not cell is interior to local rank.
Definition EclInterRegFlows.hpp:58
Minimal representation of a single named region defintion.
Definition EclInterRegFlows.hpp:184
std::string name
Region definition array name.
Definition EclInterRegFlows.hpp:186
std::reference_wrapper< const std::vector< int > > definition
Region definition array.
Definition EclInterRegFlows.hpp:189