My Project
Loading...
Searching...
No Matches
VFPProdTable.hpp
1/*
2 Copyright 2015 SINTEF ICT, Applied Mathematics.
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_PARSER_ECLIPSE_ECLIPSESTATE_TABLES_VFPPRODTABLE_HPP_
21#define OPM_PARSER_ECLIPSE_ECLIPSESTATE_TABLES_VFPPRODTABLE_HPP_
22
23#include <opm/common/OpmLog/KeywordLocation.hpp>
24
25#include <opm/input/eclipse/Units/Dimension.hpp>
26
27#include <array>
28#include <utility>
29#include <vector>
30
31namespace Opm {
32
33 class DeckItem;
34 class DeckKeyword;
35 class UnitSystem;
36
37} // namespace Opm
38
39namespace Opm {
40
45public:
46
47 enum class FLO_TYPE {
48 FLO_OIL=1,
49 FLO_LIQ,
50 FLO_GAS
51 };
52
53 enum class WFR_TYPE {
54 WFR_WOR=11,
55 WFR_WCT,
56 WFR_WGR
57 };
58
59 enum class GFR_TYPE {
60 GFR_GOR=21,
61 GFR_GLR,
62 GFR_OGR
63 };
64
65 enum class ALQ_TYPE {
66 ALQ_GRAT=31,
67 ALQ_IGLR,
68 ALQ_TGLR,
69 ALQ_PUMP,
70 ALQ_COMP,
71 ALQ_BEAN,
72 ALQ_UNDEF
73 };
74
76 VFPProdTable( const DeckKeyword& table, bool gaslift_opt_active, const UnitSystem& deck_unit_system);
77 VFPProdTable(int table_num,
78 double datum_depth,
79 FLO_TYPE flo_type,
80 WFR_TYPE wfr_type,
81 GFR_TYPE gfr_type,
82 ALQ_TYPE alq_type,
83 const std::vector<double>& flo_data,
84 const std::vector<double>& thp_data,
85 const std::vector<double>& wfr_data,
86 const std::vector<double>& gfr_data,
87 const std::vector<double>& alq_data,
88 const std::vector<double>& data);
89
90 static Dimension ALQDimension(const ALQ_TYPE& alq_type, const UnitSystem& unit_system);
91
92 static std::pair<int, ALQ_TYPE>
93 getALQType(const DeckKeyword& inputTable, bool gaslift_opt_active);
94
95 static VFPProdTable serializationTestObject();
96
97 inline int getTableNum() const {
98 return m_table_num;
99 }
100
101 inline const KeywordLocation& location() const {
102 return this->m_location;
103 }
104
105 // The name() method is added to simplify serialization.
106 inline int name() const {
107 return m_table_num;
108 }
109
110 inline double getDatumDepth() const {
111 return m_datum_depth;
112 }
113
114 inline FLO_TYPE getFloType() const {
115 return m_flo_type;
116 }
117
118 inline WFR_TYPE getWFRType() const {
119 return m_wfr_type;
120 }
121
122 inline GFR_TYPE getGFRType() const {
123 return m_gfr_type;
124 }
125
126 inline ALQ_TYPE getALQType() const {
127 return m_alq_type;
128 }
129
130 inline const std::vector<double>& getFloAxis() const {
131 return m_flo_data;
132 }
133
134 inline const std::vector<double>& getTHPAxis() const {
135 return m_thp_data;
136 }
137
138 inline const std::vector<double>& getWFRAxis() const {
139 return m_wfr_data;
140 }
141
142 inline const std::vector<double>& getGFRAxis() const {
143 return m_gfr_data;
144 }
145
146 inline const std::vector<double>& getALQAxis() const {
147 return m_alq_data;
148 }
149
164 inline const std::vector<double>& getTable() const {
165 return m_data;
166 }
167
168 bool operator==(const VFPProdTable& data) const;
169
170 std::array<size_t,5> shape() const;
171
172 double operator()(size_t thp_idx, size_t wfr_idx, size_t gfr_idx, size_t alq_idx, size_t flo_idx) const;
173
174 template<class Serializer>
175 void serializeOp(Serializer& serializer)
176 {
177 serializer(m_table_num);
178 serializer(m_datum_depth);
179 serializer(m_flo_type);
180 serializer(m_wfr_type);
181 serializer(m_gfr_type);
182 serializer(m_alq_type);
183 serializer(m_flo_data);
184 serializer(m_thp_data);
185 serializer(m_wfr_data);
186 serializer(m_gfr_data);
187 serializer(m_alq_data);
188 serializer(m_data);
189 serializer(m_location);
190 }
191
192private:
193 int m_table_num;
194 double m_datum_depth;
195 FLO_TYPE m_flo_type;
196 WFR_TYPE m_wfr_type;
197 GFR_TYPE m_gfr_type;
198 ALQ_TYPE m_alq_type;
199
200 std::vector<double> m_flo_data;
201 std::vector<double> m_thp_data;
202 std::vector<double> m_wfr_data;
203 std::vector<double> m_gfr_data;
204 std::vector<double> m_alq_data;
205
206 std::vector<double> m_data;
207 KeywordLocation m_location;
208
209 void check();
210
211 double& operator()(size_t thp_idx, size_t wfr_idx, size_t gfr_idx, size_t alq_idx, size_t flo_idx);
212
213 static void scaleValues(std::vector<double>& values,
214 const double& scaling_factor);
215
216 static void convertFloToSI(const FLO_TYPE& type,
217 std::vector<double>& values,
218 const UnitSystem& unit_system);
219 static void convertTHPToSI(std::vector<double>& values,
220 const UnitSystem& unit_system);
221 static void convertWFRToSI(const WFR_TYPE& type,
222 std::vector<double>& values,
223 const UnitSystem& unit_system);
224 static void convertGFRToSI(const GFR_TYPE& type,
225 std::vector<double>& values,
226 const UnitSystem& unit_system);
227 static void convertAlqToSI(const ALQ_TYPE& type,
228 std::vector<double>& values,
229 const UnitSystem& unit_system);
230};
231
232} // namespace Opm
233
234#endif // OPM_PARSER_ECLIPSE_ECLIPSESTATE_TABLES_VFPPRODTABLE_HPP_
Definition DeckKeyword.hpp:36
Definition Dimension.hpp:27
Definition KeywordLocation.hpp:27
Class for (de-)serializing.
Definition Serializer.hpp:94
Definition UnitSystem.hpp:34
Class for reading data from a VFPPROD (vertical flow performance production) table.
Definition VFPProdTable.hpp:44
const std::vector< double > & getTable() const
Returns the data of the table itself.
Definition VFPProdTable.hpp:164
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition Exceptions.hpp:30