19#ifndef OPM_RUNSPEC_HPP
20#define OPM_RUNSPEC_HPP
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>
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;
42 static Phases serializationTestObject();
44 bool active( Phase )
const noexcept;
45 size_t size()
const noexcept;
47 bool operator==(
const Phases& data)
const;
49 template<
class Serializer>
56 std::bitset< NUM_PHASES_IN_ENUM > bits;
65 static Welldims serializationTestObject();
67 int maxConnPerWell()
const
72 int maxWellsPerGroup()
const
77 int maxGroupsInField()
const
82 int maxWellsInField()
const
87 int maxWellListsPrWell()
const
89 return this->nWlistPrWellMax;
92 int maxDynamicWellLists()
const
94 return this->nDynWlistMax;
97 const std::optional<KeywordLocation>& location()
const
99 return this->m_location;
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();
111 bool operator==(
const Welldims& data)
const {
112 return this->location() == data.location() &&
113 rst_cmp(*
this, data);
116 template<
class Serializer>
123 serializer(nWlistPrWellMax);
124 serializer(nDynWlistMax);
125 serializer(m_location);
133 int nWlistPrWellMax { 1 };
134 int nDynWlistMax { 1 };
135 std::optional<KeywordLocation> m_location;
145 int maxSegmentedWells()
const
147 return this->nSegWellMax;
150 int maxSegmentsPerWell()
const
152 return this->nSegmentMax;
155 int maxLateralBranchesPerWell()
const
157 return this->nLatBranchMax;
160 const std::optional<KeywordLocation>& location()
const
162 return this->location_;
167 template<
class Serializer>
170 serializer(nSegWellMax);
171 serializer(nSegmentMax);
172 serializer(nLatBranchMax);
173 serializer(location_);
180 std::optional<KeywordLocation> location_;
190 int maxNONodes()
const
192 return this->nMaxNoNodes;
195 int maxNoBranches()
const
197 return this->nMaxNoBranches;
200 int maxNoBranchesConToNode()
const
202 return this->nMaxNoBranchesConToNode;
209 template<
class Serializer>
212 serializer(nMaxNoNodes);
213 serializer(nMaxNoBranches);
214 serializer(nMaxNoBranchesConToNode);
220 int nMaxNoBranchesConToNode;
230 int maxAnalyticAquifers()
const
232 return this->maxNumAnalyticAquifers;
235 int maxAnalyticAquiferConnections()
const
237 return this->maxNumAnalyticAquiferConn;
240 template <
class Serializer>
243 serializer(this->maxNumAnalyticAquifers);
244 serializer(this->maxNumAnalyticAquiferConn);
248 int maxNumAnalyticAquifers;
249 int maxNumAnalyticAquiferConn;
309 template<
class Serializer>
312 serializer(activeHyst);
313 serializer(pcHystMod);
314 serializer(krHystMod);
315 serializer(modParamTrappedValue);
316 serializer(curvatureCapPrsValue);
317 serializer(activeWagHyst);
322 bool activeHyst {
false };
325 int pcHystMod { -1 };
326 int krHystMod { -1 };
328 double modParamTrappedValue { 0.1 };
330 double curvatureCapPrsValue { 0.1 };
333 bool activeWagHyst {
false };
338 enum class ThreePhaseOilKrModel {
344 enum class KeywordFamily {
355 const ThreePhaseOilKrModel model,
356 const KeywordFamily family);
360 double minimumRelpermMobilityThreshold()
const
362 return this->tolcrit;
365 ThreePhaseOilKrModel krModel()
const
367 return this->krmodel;
370 KeywordFamily family()
const
372 return this->satfunc_family;
377 template<
class Serializer>
382 serializer(satfunc_family);
387 ThreePhaseOilKrModel krmodel = ThreePhaseOilKrModel::Default;
388 KeywordFamily satfunc_family = KeywordFamily::Undefined;
395 explicit Nupcol(
int min_value);
396 void update(
int value);
399 static Nupcol serializationTestObject();
400 bool operator==(
const Nupcol& data)
const;
402 template<
class Serializer>
404 serializer(this->nupcol_value);
405 serializer(this->min_nupcol);
420 int water_tracers()
const;
422 template<
class 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);
433 static Tracers serializationTestObject();
434 bool operator==(
const Tracers& data)
const;
441 bool diffusion_control;
454 static Runspec serializationTestObject();
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;
462 const Welldims& wellDimensions()
const noexcept;
464 const NetworkDims& networkDimensions()
const noexcept;
466 int eclPhaseMask( )
const noexcept;
468 const Actdims& actdims()
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;
477 bool operator==(
const Runspec& data)
const;
478 static bool rst_cmp(
const Runspec& full_state,
const Runspec& rst_state);
480 template<
class Serializer>
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);
494 serializer(m_actdims);
495 serializer(m_sfuncctrl);
496 serializer(m_nupcol);
497 serializer(m_co2storage);
498 serializer(m_h2storage);
504 std::time_t m_start_time;
Definition Actdims.hpp:30
Definition Runspec.hpp:223
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