My Project
Loading...
Searching...
No Matches
AggregateUDQData.hpp
1/*
2 Copyright (c) 2018 Statoil 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_UDQ_DATA_HPP
21#define OPM_AGGREGATE_UDQ_DATA_HPP
22
24
25#include <opm/io/eclipse/PaddedOutputString.hpp>
26
27#include <cstddef>
28#include <optional>
29#include <string>
30#include <vector>
31
32namespace Opm {
33 class Group;
34 class Schedule;
35 class UDQActive;
36 class UDQConfig;
37 class UDQDims;
38 class UDQInput;
39 class UDQState;
40} // Opm
41
42namespace Opm::RestartIO::Helpers {
43
45{
46public:
47 explicit AggregateUDQData(const UDQDims& udqDims);
48
49 void captureDeclaredUDQData(const Schedule& sched,
50 const std::size_t simStep,
51 const UDQState& udqState,
52 const std::vector<int>& inteHead);
53
54 const std::vector<int>& getIUDQ() const
55 {
56 return this->iUDQ_.data();
57 }
58
60 const std::optional<WindowedArray<int>>& getIUAD() const
61 {
62 return this->iUAD_;
63 }
64
65 const std::vector<EclIO::PaddedOutputString<8>>& getZUDN() const
66 {
67 return this->zUDN_.data();
68 }
69
70 const std::vector<EclIO::PaddedOutputString<8>>& getZUDL() const
71 {
72 return this->zUDL_.data();
73 }
74
77 const std::optional<WindowedArray<int>>& getIGPH() const
78 {
79 return this->iGPH_;
80 }
81
83 const std::optional<WindowedArray<int>>& getIUAP() const
84 {
85 return this->iUAP_;
86 }
87
89 const std::optional<WindowedArray<double>>& getDUDF() const
90 {
91 return this->dUDF_;
92 }
93
95 const std::optional<WindowedArray<double>>& getDUDG() const
96 {
97 return this->dUDG_;
98 }
99
101 const std::optional<WindowedMatrix<double>>& getDUDS() const
102 {
103 return this->dUDS_;
104 }
105
107 const std::optional<WindowedArray<double>>& getDUDW() const
108 {
109 return this->dUDW_;
110 }
111
112private:
116 WindowedArray<int> iUDQ_;
117
122 std::optional<WindowedArray<int>> iUAD_{};
123
127 WindowedArray<EclIO::PaddedOutputString<8>> zUDN_;
128
132 WindowedArray<EclIO::PaddedOutputString<8>> zUDL_;
133
138 std::optional<WindowedArray<int>> iGPH_{};
139
143 std::optional<WindowedArray<int>> iUAP_{};
144
148 std::optional<WindowedArray<double>> dUDF_{};
149
154 std::optional<WindowedArray<double>> dUDG_{};
155
161 std::optional<WindowedMatrix<double>> dUDS_{};
162
167 std::optional<WindowedArray<double>> dUDW_{};
168
169 void collectUserDefinedQuantities(const std::vector<UDQInput>& udqInput,
170 const std::vector<int>& inteHead);
171
172 void collectUserDefinedArguments(const Schedule& sched,
173 const std::size_t simStep,
174 const std::vector<int>& inteHead);
175
176 void collectFieldUDQValues(const std::vector<UDQInput>& udqInput,
177 const UDQState& udq_state,
178 const int expectNumFieldUDQs);
179
180 void collectGroupUDQValues(const std::vector<UDQInput>& udqInput,
181 const UDQState& udqState,
182 const std::size_t ngmax,
183 const std::vector<const Group*>& groups,
184 const int expectedNumGroupUDQs);
185
186 void collectSegmentUDQValues(const std::vector<UDQInput>& udqInput,
187 const UDQState& udqState,
188 const std::vector<std::string>& msWells);
189
190 void collectWellUDQValues(const std::vector<UDQInput>& udqInput,
191 const UDQState& udqState,
192 const std::size_t nwmax,
193 const std::vector<std::string>& wells,
194 const int expectedNumWellUDQs);
195
202 void collectIUAD(const UDQActive& udqActive,
203 const std::size_t expectNumIUAD);
204
211 void collectIUAP(const std::vector<int>& wgIndex,
212 const std::size_t expectNumIUAP);
213
221 void collectIGPH(const std::vector<int>& phase_vector,
222 const std::size_t expectNumIGPH);
223};
224
225} // Opm::RestartIO::Helpers
226
227#endif // OPM_AGGREGATE_UDQ_DATA_HPP
Provide facilities to simplify constructing restart vectors such as IWEL or RSEG.
Definition AggregateUDQData.hpp:45
const std::optional< WindowedArray< int > > & getIUAP() const
Associate well/group IDs for IUAD. Nullopt if no UDAs in use.
Definition AggregateUDQData.hpp:83
const std::optional< WindowedArray< int > > & getIUAD() const
Retrieve UDA descriptive data. Nullopt if no UDAs in use.
Definition AggregateUDQData.hpp:60
const std::optional< WindowedArray< int > > & getIGPH() const
Retrive group level injection phase UDAs.
Definition AggregateUDQData.hpp:77
const std::optional< WindowedArray< double > > & getDUDW() const
Retrieve values of well level UDQs. Nullopt if no such UDQs exist.
Definition AggregateUDQData.hpp:107
const std::optional< WindowedArray< double > > & getDUDG() const
Retrieve values of group level UDQs. Nullopt if no such UDQs exist.
Definition AggregateUDQData.hpp:95
const std::optional< WindowedMatrix< double > > & getDUDS() const
Retrieve values of segment level UDQs. Nullopt if no such UDQs exist.
Definition AggregateUDQData.hpp:101
const std::optional< WindowedArray< double > > & getDUDF() const
Retrieve values of field level UDQs. Nullopt if no such UDQs exist.
Definition AggregateUDQData.hpp:89
Provide read-only and read/write access to constantly sized portions/windows of a linearised buffer w...
Definition WindowedArray.hpp:50
const std::vector< T > & data() const
Get read-only access to full, linearised data items for all windows.
Definition WindowedArray.hpp:136
Definition Schedule.hpp:101
Collection of UDQ and UDA related dimension queries.
Definition UDQDims.hpp:40
Definition UDQState.hpp:40
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition Exceptions.hpp:30