My Project
Loading...
Searching...
No Matches
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 <opm/common/OpmLog/KeywordLocation.hpp>
23#include <opm/input/eclipse/EclipseState/Phase.hpp>
24#include <opm/input/eclipse/EclipseState/Tables/Tabdims.hpp>
25#include <opm/input/eclipse/EclipseState/Tables/Regdims.hpp>
26#include <opm/input/eclipse/EclipseState/EndpointScaling.hpp>
27#include <opm/input/eclipse/Schedule/UDQ/UDQParams.hpp>
28#include <opm/input/eclipse/Schedule/Action/Actdims.hpp>
29
30#include <optional>
31#include <string>
32
33namespace Opm {
34class Deck;
35
36class Phases {
37public:
38 Phases() noexcept = default;
39 Phases( bool oil, bool gas, bool wat, bool solvent = false, bool polymer = false, bool energy = false,
40 bool polymw = false, bool foam = false, bool brine = false, bool zfraction = false ) noexcept;
41
42 static Phases serializationTestObject();
43
44 bool active( Phase ) const noexcept;
45 size_t size() const noexcept;
46
47 bool operator==(const Phases& data) const;
48
49 template<class Serializer>
50 void serializeOp(Serializer& serializer)
51 {
52 serializer(bits);
53 }
54
55private:
56 std::bitset< NUM_PHASES_IN_ENUM > bits;
57};
58
59
60class Welldims {
61public:
62 Welldims() = default;
63 explicit Welldims(const Deck& deck);
64
65 static Welldims serializationTestObject();
66
67 int maxConnPerWell() const
68 {
69 return this->nCWMax;
70 }
71
72 int maxWellsPerGroup() const
73 {
74 return this->nWGMax;
75 }
76
77 int maxGroupsInField() const
78 {
79 return this->nGMax;
80 }
81
82 int maxWellsInField() const
83 {
84 return this->nWMax;
85 }
86
87 int maxWellListsPrWell() const
88 {
89 return this->nWlistPrWellMax;
90 }
91
92 int maxDynamicWellLists() const
93 {
94 return this->nDynWlistMax;
95 }
96
97 const std::optional<KeywordLocation>& location() const
98 {
99 return this->m_location;
100 }
101
102 static bool rst_cmp(const Welldims& full_dims, const Welldims& rst_dims) {
103 return full_dims.maxConnPerWell() == rst_dims.maxConnPerWell() &&
104 full_dims.maxWellsPerGroup() == rst_dims.maxWellsPerGroup() &&
105 full_dims.maxGroupsInField() == rst_dims.maxGroupsInField() &&
106 full_dims.maxWellsInField() == rst_dims.maxWellsInField() &&
107 full_dims.maxWellListsPrWell() == rst_dims.maxWellListsPrWell() &&
108 full_dims.maxDynamicWellLists() == rst_dims.maxDynamicWellLists();
109 }
110
111 bool operator==(const Welldims& data) const {
112 return this->location() == data.location() &&
113 rst_cmp(*this, data);
114 }
115
116 template<class Serializer>
117 void serializeOp(Serializer& serializer)
118 {
119 serializer(nWMax);
120 serializer(nCWMax);
121 serializer(nWGMax);
122 serializer(nGMax);
123 serializer(nWlistPrWellMax);
124 serializer(nDynWlistMax);
125 serializer(m_location);
126 }
127
128private:
129 int nWMax { 0 };
130 int nCWMax { 0 };
131 int nWGMax { 0 };
132 int nGMax { 0 };
133 int nWlistPrWellMax { 1 };
134 int nDynWlistMax { 1 };
135 std::optional<KeywordLocation> m_location;
136};
137
139public:
141 explicit WellSegmentDims(const Deck& deck);
142
143 static WellSegmentDims serializationTestObject();
144
145 int maxSegmentedWells() const
146 {
147 return this->nSegWellMax;
148 }
149
150 int maxSegmentsPerWell() const
151 {
152 return this->nSegmentMax;
153 }
154
155 int maxLateralBranchesPerWell() const
156 {
157 return this->nLatBranchMax;
158 }
159
160 const std::optional<KeywordLocation>& location() const
161 {
162 return this->location_;
163 }
164
165 bool operator==(const WellSegmentDims& data) const;
166
167 template<class Serializer>
168 void serializeOp(Serializer& serializer)
169 {
170 serializer(nSegWellMax);
171 serializer(nSegmentMax);
172 serializer(nLatBranchMax);
173 serializer(location_);
174 }
175
176private:
177 int nSegWellMax;
178 int nSegmentMax;
179 int nLatBranchMax;
180 std::optional<KeywordLocation> location_;
181};
182
184public:
185 NetworkDims();
186 explicit NetworkDims(const Deck& deck);
187
188 static NetworkDims serializationTestObject();
189
190 int maxNONodes() const
191 {
192 return this->nMaxNoNodes;
193 }
194
195 int maxNoBranches() const
196 {
197 return this->nMaxNoBranches;
198 }
199
200 int maxNoBranchesConToNode() const
201 {
202 return this->nMaxNoBranchesConToNode;
203 }
204
205 bool active() const;
206
207 bool operator==(const NetworkDims& data) const;
208
209 template<class Serializer>
210 void serializeOp(Serializer& serializer)
211 {
212 serializer(nMaxNoNodes);
213 serializer(nMaxNoBranches);
214 serializer(nMaxNoBranchesConToNode);
215 }
216
217private:
218 int nMaxNoNodes;
219 int nMaxNoBranches;
220 int nMaxNoBranchesConToNode;
221};
222
224public:
226 explicit AquiferDimensions(const Deck& deck);
227
228 static AquiferDimensions serializationTestObject();
229
230 int maxAnalyticAquifers() const
231 {
232 return this->maxNumAnalyticAquifers;
233 }
234
235 int maxAnalyticAquiferConnections() const
236 {
237 return this->maxNumAnalyticAquiferConn;
238 }
239
240 template <class Serializer>
241 void serializeOp(Serializer& serializer)
242 {
243 serializer(this->maxNumAnalyticAquifers);
244 serializer(this->maxNumAnalyticAquiferConn);
245 }
246
247private:
248 int maxNumAnalyticAquifers;
249 int maxNumAnalyticAquiferConn;
250};
251
252bool operator==(const AquiferDimensions& lhs, const AquiferDimensions& rhs);
253
255{
256public:
257 EclHysterConfig() = default;
258 explicit EclHysterConfig(const Deck& deck);
259
260 static EclHysterConfig serializationTestObject();
261
265 //void setActive(bool yesno);
266
270 bool active() const;
271
278 int pcHysteresisModel() const;
279
286 int krHysteresisModel() const;
287
293 double modParamTrapped() const;
294
300 double curvatureCapPrs() const;
301
305 bool activeWag() const;
306
307 bool operator==(const EclHysterConfig& data) const;
308
309 template<class Serializer>
310 void serializeOp(Serializer& serializer)
311 {
312 serializer(activeHyst);
313 serializer(pcHystMod);
314 serializer(krHystMod);
315 serializer(modParamTrappedValue);
316 serializer(curvatureCapPrsValue);
317 serializer(activeWagHyst);
318 }
319
320private:
321 // enable hysteresis at all
322 bool activeHyst { false };
323
324 // the capillary pressure and the relperm hysteresis models to be used
325 int pcHystMod { -1 };
326 int krHystMod { -1 };
327 // regularisation parameter used for Killough model
328 double modParamTrappedValue { 0.1 };
329 // curvature parameter for capillary pressure
330 double curvatureCapPrsValue { 0.1 };
331
332 // enable WAG hysteresis
333 bool activeWagHyst { false };
334};
335
337public:
338 enum class ThreePhaseOilKrModel {
339 Default,
340 Stone1,
341 Stone2
342 };
343
344 enum class KeywordFamily {
345 Family_I, // SGOF, SWOF, SLGOF
346 Family_II, // SGFN, SOF{2,3}, SWFN, SGWFN
347 Family_III, // GSF, WSF
348
349 Undefined,
350 };
351
353 explicit SatFuncControls(const Deck& deck);
354 explicit SatFuncControls(const double tolcritArg,
355 const ThreePhaseOilKrModel model,
356 const KeywordFamily family);
357
358 static SatFuncControls serializationTestObject();
359
360 double minimumRelpermMobilityThreshold() const
361 {
362 return this->tolcrit;
363 }
364
365 ThreePhaseOilKrModel krModel() const
366 {
367 return this->krmodel;
368 }
369
370 KeywordFamily family() const
371 {
372 return this->satfunc_family;
373 }
374
375 bool operator==(const SatFuncControls& rhs) const;
376
377 template<class Serializer>
378 void serializeOp(Serializer& serializer)
379 {
380 serializer(tolcrit);
381 serializer(krmodel);
382 serializer(satfunc_family);
383 }
384
385private:
386 double tolcrit;
387 ThreePhaseOilKrModel krmodel = ThreePhaseOilKrModel::Default;
388 KeywordFamily satfunc_family = KeywordFamily::Undefined;
389};
390
391
392class Nupcol {
393public:
394 Nupcol();
395 explicit Nupcol(int min_value);
396 void update(int value);
397 int value() const;
398
399 static Nupcol serializationTestObject();
400 bool operator==(const Nupcol& data) const;
401
402 template<class Serializer>
403 void serializeOp(Serializer& serializer) {
404 serializer(this->nupcol_value);
405 serializer(this->min_nupcol);
406 }
407
408private:
409 int min_nupcol;
410 int nupcol_value;
411};
412
413
414class Tracers {
415public:
416
417 Tracers() = default;
418
419 explicit Tracers(const Deck& );
420 int water_tracers() const;
421
422 template<class Serializer>
423 void serializeOp(Serializer& serializer) {
424 serializer(this->m_oil_tracers);
425 serializer(this->m_water_tracers);
426 serializer(this->m_gas_tracers);
427 serializer(this->m_env_tracers);
428 serializer(this->diffusion_control);
429 serializer(this->max_iter);
430 serializer(this->min_iter);
431 }
432
433 static Tracers serializationTestObject();
434 bool operator==(const Tracers& data) const;
435
436private:
437 int m_oil_tracers;
438 int m_water_tracers;
439 int m_gas_tracers;
440 int m_env_tracers;
441 bool diffusion_control;
442 int max_iter;
443 int min_iter;
444 // The TRACERS keyword has some additional options which seem quite arcane,
445 // for now not included here.
446};
447
448
449class Runspec {
450public:
451 Runspec() = default;
452 explicit Runspec( const Deck& );
453
454 static Runspec serializationTestObject();
455
456 std::time_t start_time() const noexcept;
457 const UDQParams& udqParams() const noexcept;
458 const Phases& phases() const noexcept;
459 const Tabdims& tabdims() const noexcept;
460 const Regdims& regdims() const noexcept;
461 const EndpointScaling& endpointScaling() const noexcept;
462 const Welldims& wellDimensions() const noexcept;
463 const WellSegmentDims& wellSegmentDimensions() const noexcept;
464 const NetworkDims& networkDimensions() const noexcept;
465 const AquiferDimensions& aquiferDimensions() const noexcept;
466 int eclPhaseMask( ) const noexcept;
467 const EclHysterConfig& hysterPar() const noexcept;
468 const Actdims& actdims() const noexcept;
469 const SatFuncControls& saturationFunctionControls() const noexcept;
470 const Nupcol& nupcol() const noexcept;
471 const Tracers& tracers() const;
472 bool co2Storage() const noexcept;
473 bool h2Storage() const noexcept;
474 bool micp() const noexcept;
475 bool mech() const noexcept;
476
477 bool operator==(const Runspec& data) const;
478 static bool rst_cmp(const Runspec& full_state, const Runspec& rst_state);
479
480 template<class Serializer>
481 void serializeOp(Serializer& serializer)
482 {
483 serializer(this->m_start_time);
484 serializer(active_phases);
485 serializer(m_tabdims);
486 serializer(m_regdims);
487 serializer(endscale);
488 serializer(welldims);
489 serializer(wsegdims);
490 serializer(netwrkdims);
491 serializer(aquiferdims);
492 serializer(udq_params);
493 serializer(hystpar);
494 serializer(m_actdims);
495 serializer(m_sfuncctrl);
496 serializer(m_nupcol);
497 serializer(m_co2storage);
498 serializer(m_h2storage);
499 serializer(m_micp);
500 serializer(m_mech);
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_h2storage;
521 bool m_micp;
522 bool m_mech;
523};
524
525
526}
527
528#endif // OPM_RUNSPEC_HPP
Definition Actdims.hpp:30
Definition Runspec.hpp:223
Definition Deck.hpp:49
Definition Runspec.hpp:255
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 activeWag() const
Wag 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:183
Definition Runspec.hpp:392
Definition Runspec.hpp:36
Definition Regdims.hpp:36
Definition Runspec.hpp:449
Definition Runspec.hpp:336
Class for (de-)serializing.
Definition Serializer.hpp:84
Definition Tabdims.hpp:36
Definition Runspec.hpp:414
Definition UDQParams.hpp:31
Definition Runspec.hpp:138
Definition Runspec.hpp:60
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition Exceptions.hpp:30