My Project
Loading...
Searching...
No Matches
AggregateAquiferData.hpp
1/*
2 Copyright (c) 2021 Equinor ASA
3
4 This file is part of the Open Porous Media project (OPM).
5
6 OPM is free software: you can redistribute it and/or modify
7 it under the terms of the GNU General Public License as published by
8 the Free Software Foundation, either version 3 of the License, or
9 (at your option) any later version.
10
11 OPM is distributed in the hope that it will be useful,
12 but WITHOUT ANY WARRANTY; without even the implied warranty of
13 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 GNU General Public License for more details.
15
16 You should have received a copy of the GNU General Public License
17 along with OPM. If not, see <http://www.gnu.org/licenses/>.
18*/
19
20#ifndef OPM_AGGREGATE_AQUIFER_DATA_HPP
21#define OPM_AGGREGATE_AQUIFER_DATA_HPP
22
23#include <opm/output/eclipse/InteHEAD.hpp>
25
26#include <opm/output/data/Aquifer.hpp>
27
28#include <vector>
29
30namespace Opm {
31 class AquiferConfig;
32 class EclipseGrid;
33 class ScheduleState;
34 class SummaryState;
35 class UnitSystem;
36} // Opm
37
38namespace Opm { namespace RestartIO { namespace Helpers {
39
41 {
42 public:
57 explicit AggregateAquiferData(const InteHEAD::AquiferDims& aqDims,
58 const AquiferConfig& aqConfig,
59 const EclipseGrid& grid);
60
89 const AquiferConfig& aqConfig,
90 const ScheduleState& sched,
91 const data::Aquifers& aquData,
92 const SummaryState& summaryState,
93 const UnitSystem& usys);
94
101 {
102 return this->maxActiveAnalyticAquiferID_;
103 }
104
106 const std::vector<int>& getIntegerAquiferData() const
107 {
108 return this->integerAnalyticAq_.data();
109 }
110
112 const std::vector<float>& getSinglePrecAquiferData() const
113 {
114 return this->singleprecAnalyticAq_.data();
115 }
116
118 const std::vector<double>& getDoublePrecAquiferData() const
119 {
120 return this->doubleprecAnalyticAq_.data();
121 }
122
124 const std::vector<int>& getNumericAquiferIntegerData() const
125 {
126 return this->integerNumericAq_.data();
127 }
128
130 const std::vector<double>& getNumericAquiferDoublePrecData() const
131 {
132 return this->doubleprecNumericAq_.data();
133 }
134
140 const std::vector<int>& getIntegerAquiferConnectionData(const int aquiferID) const
141 {
142 return this->integerAnalyticAquiferConn_[aquiferID - 1].data();
143 }
144
151 const std::vector<float>& getSinglePrecAquiferConnectionData(const int aquiferID) const
152 {
153 return this->singleprecAnalyticAquiferConn_[aquiferID - 1].data();
154 }
155
163 const std::vector<double>& getDoublePrecAquiferConnectionData(const int aquiferID) const
164 {
165 return this->doubleprecAnalyticAquiferConn_[aquiferID - 1].data();
166 }
167
168 private:
169 int maxActiveAnalyticAquiferID_{0};
170
171 std::vector<int> numActiveConn_{};
172 std::vector<double> totalInflux_{};
173
175 WindowedArray<int> integerAnalyticAq_;
176
178 WindowedArray<float> singleprecAnalyticAq_;
179
181 WindowedArray<double> doubleprecAnalyticAq_;
182
184 WindowedArray<int> integerNumericAq_;
185
187 WindowedArray<double> doubleprecNumericAq_;
188
191 std::vector<WindowedArray<int>> integerAnalyticAquiferConn_;
192
195 std::vector<WindowedArray<float>> singleprecAnalyticAquiferConn_;
196
199 std::vector<WindowedArray<double>> doubleprecAnalyticAquiferConn_;
200
201 void allocateDynamicBackingStorage(const InteHEAD::AquiferDims& aqDims);
202
203 void handleCarterTracy(const AquiferConfig& aqConfig,
204 const data::Aquifers& aquData,
205 const SummaryState& summaryState,
206 const UnitSystem& usys);
207
208 void handleConstantFlux(const AquiferConfig& aqConfig,
209 const ScheduleState& sched,
210 const SummaryState& summaryState,
211 const UnitSystem& usys);
212
213 void handleFetkovich(const AquiferConfig& aqConfig,
214 const data::Aquifers& aquData,
215 const SummaryState& summaryState,
216 const UnitSystem& usys);
217
218 void handleNumeric(const AquiferConfig& aqConfig,
219 const data::Aquifers& aquData,
220 const SummaryState& summaryState,
221 const UnitSystem& usys);
222 };
223
224}}} // Opm::RestartIO::Helpers
225
226#endif // OPM_AGGREGATE_WELL_DATA_HPP
Provide facilities to simplify constructing restart vectors such as IWEL or RSEG.
Definition AquiferConfig.hpp:46
About cell information and dimension: The actual grid information is held in a pointer to an ERT ecl_...
Definition EclipseGrid.hpp:55
Definition AggregateAquiferData.hpp:41
const std::vector< double > & getDoublePrecAquiferData() const
Retrieve Floating-Point (Double Precision) Aquifer Data Array.
Definition AggregateAquiferData.hpp:118
const std::vector< float > & getSinglePrecAquiferConnectionData(const int aquiferID) const
Retrieve Floating-Point (Real) Aquifer Connection Data Array (analytic aquifers)
Definition AggregateAquiferData.hpp:151
const std::vector< float > & getSinglePrecAquiferData() const
Retrieve Floating-Point (Real) Aquifer Data Array.
Definition AggregateAquiferData.hpp:112
const std::vector< int > & getIntegerAquiferData() const
Retrieve Integer Aquifer Data Array.
Definition AggregateAquiferData.hpp:106
const std::vector< int > & getNumericAquiferIntegerData() const
Retrieve Integer Aquifer Data Array for Numeric Aquifers.
Definition AggregateAquiferData.hpp:124
const std::vector< int > & getIntegerAquiferConnectionData(const int aquiferID) const
Retrieve Integer Aquifer Connection Data Array (analytic aquifers)
Definition AggregateAquiferData.hpp:140
const std::vector< double > & getNumericAquiferDoublePrecData() const
Retrieve Double Precision Aquifer Data Array for Numeric Aquifers.
Definition AggregateAquiferData.hpp:130
void captureDynamicAquiferData(const InteHEAD::AquiferDims &aqDims, const AquiferConfig &aqConfig, const ScheduleState &sched, const data::Aquifers &aquData, const SummaryState &summaryState, const UnitSystem &usys)
Linearise dynamic information pertinent to analytic aquifers into internal arrays.
Definition AggregateAquiferData.cpp:648
int maximumActiveAnalyticAquiferID() const
Retrieve the maximum active aquifer ID over all analytic aquifers.
Definition AggregateAquiferData.hpp:100
const std::vector< double > & getDoublePrecAquiferConnectionData(const int aquiferID) const
Retrieve Floating-Point (Double Precision) Aquifer Connection Data Array (analytic aquifers)
Definition AggregateAquiferData.hpp:163
const std::vector< T > & data() const
Get read-only access to full, linearised data items for all windows.
Definition WindowedArray.hpp:131
Definition ScheduleState.hpp:91
Definition SummaryState.hpp:68
Definition UnitSystem.hpp:34
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition Exceptions.hpp:30
Definition InteHEAD.hpp:172