My Project
Runspec.hpp
1/*
2 Copyright 2016 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 it under the terms
7 of the GNU General Public License as published by the Free Software
8 Foundation, either version 3 of the License, or (at your option) any later
9 version.
10
11 OPM is distributed in the hope that it will be useful, but WITHOUT ANY
12 WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
13 A PARTICULAR PURPOSE. See the GNU General Public License for more details.
14
15 You should have received a copy of the GNU General Public License along with
16 OPM. If not, see <http://www.gnu.org/licenses/>.
17*/
18
19#ifndef OPM_RUNSPEC_HPP
20#define OPM_RUNSPEC_HPP
21
22#include <iosfwd>
23#include <string>
24#include <optional>
25
26#include <opm/common/OpmLog/KeywordLocation.hpp>
27#include <opm/input/eclipse/EclipseState/Tables/Tabdims.hpp>
28#include <opm/input/eclipse/EclipseState/Tables/Regdims.hpp>
29#include <opm/input/eclipse/EclipseState/EndpointScaling.hpp>
30#include <opm/input/eclipse/Schedule/UDQ/UDQParams.hpp>
31#include <opm/input/eclipse/Schedule/Action/Actdims.hpp>
32
33namespace Opm {
34class Deck;
35
36
37enum class Phase {
38 OIL = 0,
39 GAS = 1,
40 WATER = 2,
41 SOLVENT = 3,
42 POLYMER = 4,
43 ENERGY = 5,
44 POLYMW = 6,
45 FOAM = 7,
46 BRINE = 8,
47 ZFRACTION = 9
48
49 // If you add more entries to this enum, remember to update NUM_PHASES_IN_ENUM below.
50};
51
52constexpr int NUM_PHASES_IN_ENUM = static_cast<int>(Phase::ZFRACTION) + 1; // Used to get correct size of the bitset in class Phases.
53
54Phase get_phase( const std::string& );
55std::ostream& operator<<( std::ostream&, const Phase& );
56
57class Phases {
58 public:
59 Phases() noexcept = default;
60 Phases( bool oil, bool gas, bool wat, bool solvent = false, bool polymer = false, bool energy = false,
61 bool polymw = false, bool foam = false, bool brine = false, bool zfraction = false ) noexcept;
62
63 static Phases serializationTestObject();
64
65 bool active( Phase ) const noexcept;
66 size_t size() const noexcept;
67
68 bool operator==(const Phases& data) const;
69
70 template<class Serializer>
71 void serializeOp(Serializer& serializer)
72 {
73 serializer(bits);
74 }
75
76 private:
77 std::bitset< NUM_PHASES_IN_ENUM > bits;
78};
79
80
81class Welldims {
82public:
83 Welldims() = default;
84 explicit Welldims(const Deck& deck);
85
86 static Welldims serializationTestObject();
87
88 int maxConnPerWell() const
89 {
90 return this->nCWMax;
91 }
92
93 int maxWellsPerGroup() const
94 {
95 return this->nWGMax;
96 }
97
98 int maxGroupsInField() const
99 {
100 return this->nGMax;
101 }
102
103 int maxWellsInField() const
104 {
105 return this->nWMax;
106 }
107
108 int maxWellListsPrWell() const
109 {
110 return this->nWlistPrWellMax;
111 }
112
113 int maxDynamicWellLists() const
114 {
115 return this->nDynWlistMax;
116 }
117
118 const std::optional<KeywordLocation>& location() const {
119 return this->m_location;
120 }
121
122 static bool rst_cmp(const Welldims& full_dims, const Welldims& rst_dims) {
123 return full_dims.maxConnPerWell() == rst_dims.maxConnPerWell() &&
124 full_dims.maxWellsPerGroup() == rst_dims.maxWellsPerGroup() &&
125 full_dims.maxGroupsInField() == rst_dims.maxGroupsInField() &&
126 full_dims.maxWellsInField() == rst_dims.maxWellsInField() &&
127 full_dims.maxWellListsPrWell() == rst_dims.maxWellListsPrWell() &&
128 full_dims.maxDynamicWellLists() == rst_dims.maxDynamicWellLists();
129 }
130
131 bool operator==(const Welldims& data) const {
132 return this->location() == data.location() &&
133 rst_cmp(*this, data);
134 }
135
136
137 template<class Serializer>
138 void serializeOp(Serializer& serializer)
139 {
140 serializer(nWMax);
141 serializer(nCWMax);
142 serializer(nWGMax);
143 serializer(nGMax);
144 serializer(nWlistPrWellMax);
145 serializer(nDynWlistMax);
146 serializer(m_location);
147 }
148
149private:
150 int nWMax { 0 };
151 int nCWMax { 0 };
152 int nWGMax { 0 };
153 int nGMax { 0 };
154 int nWlistPrWellMax { 1 };
155 int nDynWlistMax { 1 };
156 std::optional<KeywordLocation> m_location;
157};
158
160public:
162 explicit WellSegmentDims(const Deck& deck);
163
164 static WellSegmentDims serializationTestObject();
165
166
167 int maxSegmentedWells() const
168 {
169 return this->nSegWellMax;
170 }
171
172 int maxSegmentsPerWell() const
173 {
174 return this->nSegmentMax;
175 }
176
177 int maxLateralBranchesPerWell() const
178 {
179 return this->nLatBranchMax;
180 }
181
182 bool operator==(const WellSegmentDims& data) const;
183
184 template<class Serializer>
185 void serializeOp(Serializer& serializer)
186 {
187 serializer(nSegWellMax);
188 serializer(nSegmentMax);
189 serializer(nLatBranchMax);
190 }
191
192private:
193 int nSegWellMax;
194 int nSegmentMax;
195 int nLatBranchMax;
196};
197
199public:
200 NetworkDims();
201 explicit NetworkDims(const Deck& deck);
202
203 static NetworkDims serializationTestObject();
204
205 int maxNONodes() const
206 {
207 return this->nMaxNoNodes;
208 }
209
210 int maxNoBranches() const
211 {
212 return this->nMaxNoBranches;
213 }
214
215 int maxNoBranchesConToNode() const
216 {
217 return this->nMaxNoBranchesConToNode;
218 }
219
220 bool active() const;
221
222 bool operator==(const NetworkDims& data) const;
223
224 template<class Serializer>
225 void serializeOp(Serializer& serializer)
226 {
227 serializer(nMaxNoNodes);
228 serializer(nMaxNoBranches);
229 serializer(nMaxNoBranchesConToNode);
230 }
231
232private:
233 int nMaxNoNodes;
234 int nMaxNoBranches;
235 int nMaxNoBranchesConToNode;
236};
237
239public:
241 explicit AquiferDimensions(const Deck& deck);
242
243 static AquiferDimensions serializationTestObject();
244
245 int maxAnalyticAquifers() const
246 {
247 return this->maxNumAnalyticAquifers;
248 }
249
250 int maxAnalyticAquiferConnections() const
251 {
252 return this->maxNumAnalyticAquiferConn;
253 }
254
255 template <class Serializer>
256 void serializeOp(Serializer& serializer)
257 {
258 serializer(this->maxNumAnalyticAquifers);
259 serializer(this->maxNumAnalyticAquiferConn);
260 }
261
262private:
263 int maxNumAnalyticAquifers;
264 int maxNumAnalyticAquiferConn;
265};
266
267bool operator==(const AquiferDimensions& lhs, const AquiferDimensions& rhs);
268
270{
271public:
272 EclHysterConfig() = default;
273 explicit EclHysterConfig(const Deck& deck);
274
275 static EclHysterConfig serializationTestObject();
276
280 //void setActive(bool yesno);
281
285 bool active() const;
286
293 int pcHysteresisModel() const;
294
301 int krHysteresisModel() const;
302
308 double modParamTrapped() const;
309
315 double curvatureCapPrs() const;
316
317 bool operator==(const EclHysterConfig& data) const;
318
319 template<class Serializer>
320 void serializeOp(Serializer& serializer)
321 {
322 serializer(activeHyst);
323 serializer(pcHystMod);
324 serializer(krHystMod);
325 serializer(modParamTrappedValue);
326 serializer(curvatureCapPrsValue);
327 }
328
329private:
330 // enable hysteresis at all
331 bool activeHyst { false };
332
333 // the capillary pressure and the relperm hysteresis models to be used
334 int pcHystMod { -1 };
335 int krHystMod { -1 };
336 // regularisation parameter used for Killough model
337 double modParamTrappedValue { 0.1 };
338 // curvature parameter for capillary pressure
339 double curvatureCapPrsValue { 0.1 };
340};
341
343public:
344 enum class ThreePhaseOilKrModel {
345 Default,
346 Stone1,
347 Stone2
348 };
349
350 enum class KeywordFamily {
351 Family_I, // SGOF, SWOF, SLGOF
352 Family_II, // SGFN, SOF{2,3}, SWFN
353 Undefined,
354 };
355
357 explicit SatFuncControls(const Deck& deck);
358 explicit SatFuncControls(const double tolcritArg,
359 const ThreePhaseOilKrModel model,
360 const KeywordFamily family);
361
362 static SatFuncControls serializationTestObject();
363
364 double minimumRelpermMobilityThreshold() const
365 {
366 return this->tolcrit;
367 }
368
369 ThreePhaseOilKrModel krModel() const
370 {
371 return this->krmodel;
372 }
373
374 KeywordFamily family() const
375 {
376 return this->satfunc_family;
377 }
378
379 bool operator==(const SatFuncControls& rhs) const;
380
381 template<class Serializer>
382 void serializeOp(Serializer& serializer)
383 {
384 serializer(tolcrit);
385 serializer(krmodel);
386 serializer(satfunc_family);
387 }
388
389private:
390 double tolcrit;
391 ThreePhaseOilKrModel krmodel = ThreePhaseOilKrModel::Default;
392 KeywordFamily satfunc_family = KeywordFamily::Undefined;
393};
394
395
396class Nupcol {
397public:
398 Nupcol();
399 explicit Nupcol(int min_value);
400 void update(int value);
401 int value() const;
402
403 static Nupcol serializationTestObject();
404 bool operator==(const Nupcol& data) const;
405
406 template<class Serializer>
407 void serializeOp(Serializer& serializer) {
408 serializer(this->nupcol_value);
409 serializer(this->min_nupcol);
410 }
411
412private:
413 int min_nupcol;
414 int nupcol_value;
415};
416
417
418class Tracers {
419public:
420
421 Tracers() = default;
422
423 explicit Tracers(const Deck& );
424 int water_tracers() const;
425
426 template<class Serializer>
427 void serializeOp(Serializer& serializer) {
428 serializer(this->m_oil_tracers);
429 serializer(this->m_water_tracers);
430 serializer(this->m_gas_tracers);
431 serializer(this->m_env_tracers);
432 serializer(this->diffusion_control);
433 serializer(this->max_iter);
434 serializer(this->min_iter);
435 }
436
437 static Tracers serializationTestObject();
438 bool operator==(const Tracers& data) const;
439
440private:
441 int m_oil_tracers;
442 int m_water_tracers;
443 int m_gas_tracers;
444 int m_env_tracers;
445 bool diffusion_control;
446 int max_iter;
447 int min_iter;
448 // The TRACERS keyword has some additional options which seem quite arcane,
449 // for now not included here.
450};
451
452
453class Runspec {
454public:
455 Runspec() = default;
456 explicit Runspec( const Deck& );
457
458 static Runspec serializationTestObject();
459
460 std::time_t start_time() const noexcept;
461 const UDQParams& udqParams() const noexcept;
462 const Phases& phases() const noexcept;
463 const Tabdims& tabdims() const noexcept;
464 const Regdims& regdims() const noexcept;
465 const EndpointScaling& endpointScaling() const noexcept;
466 const Welldims& wellDimensions() const noexcept;
467 const WellSegmentDims& wellSegmentDimensions() const noexcept;
468 const NetworkDims& networkDimensions() const noexcept;
469 const AquiferDimensions& aquiferDimensions() const noexcept;
470 int eclPhaseMask( ) const noexcept;
471 const EclHysterConfig& hysterPar() const noexcept;
472 const Actdims& actdims() const noexcept;
473 const SatFuncControls& saturationFunctionControls() const noexcept;
474 const Nupcol& nupcol() const noexcept;
475 const Tracers& tracers() const;
476 bool co2Storage() const noexcept;
477 bool micp() const noexcept;
478
479 bool operator==(const Runspec& data) const;
480 static bool rst_cmp(const Runspec& full_state, const Runspec& rst_state);
481
482 template<class Serializer>
483 void serializeOp(Serializer& serializer)
484 {
485 serializer(this->m_start_time);
486 serializer(active_phases);
487 serializer(m_tabdims);
488 serializer(m_regdims);
489 serializer(endscale);
490 serializer(welldims);
491 serializer(wsegdims);
492 serializer(netwrkdims);
493 serializer(aquiferdims);
494 serializer(udq_params);
495 serializer(hystpar);
496 serializer(m_actdims);
497 serializer(m_sfuncctrl);
498 serializer(m_nupcol);
499 serializer(m_co2storage);
500 serializer(m_micp);
501 }
502
503private:
504 std::time_t m_start_time;
505 Phases active_phases;
506 Tabdims m_tabdims;
507 Regdims m_regdims;
508 EndpointScaling endscale;
509 Welldims welldims;
510 WellSegmentDims wsegdims;
511 NetworkDims netwrkdims;
512 AquiferDimensions aquiferdims;
513 UDQParams udq_params;
514 EclHysterConfig hystpar;
515 Actdims m_actdims;
516 SatFuncControls m_sfuncctrl;
517 Nupcol m_nupcol;
518 Tracers m_tracers;
519 bool m_co2storage;
520 bool m_micp;
521};
522
523
524}
525
526#endif // OPM_RUNSPEC_HPP
Definition: Actdims.hpp:30
Definition: Runspec.hpp:238
Definition: Deck.hpp:63
Definition: Runspec.hpp:270
int pcHysteresisModel() const
Return the type of the hysteresis model which is used for capillary pressure.
double modParamTrapped() const
Regularisation parameter used for Killough model.
double curvatureCapPrs() const
Curvature parameter used for capillary pressure hysteresis.
bool active() const
Specify whether hysteresis is enabled or not.
int krHysteresisModel() const
Return the type of the hysteresis model which is used for relative permeability.
Definition: EndpointScaling.hpp:28
Definition: Runspec.hpp:198
Definition: Runspec.hpp:396
Definition: Runspec.hpp:57
Definition: Regdims.hpp:36
Definition: Runspec.hpp:453
Definition: Runspec.hpp:342
Class for (de-)serializing.
Definition: Serializer.hpp:75
Definition: Tabdims.hpp:36
Definition: Runspec.hpp:418
Definition: UDQParams.hpp:31
Definition: Runspec.hpp:159
Definition: Runspec.hpp:81
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition: Exceptions.hpp:29