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>
39 Phases()
noexcept =
default;
40 Phases(
bool oil,
bool gas,
bool wat,
bool solvent =
false,
bool polymer =
false,
bool energy =
false,
41 bool polymw =
false,
bool foam =
false,
bool brine =
false,
bool zfraction =
false )
noexcept;
43 static Phases serializationTestObject();
45 bool active( Phase )
const noexcept;
46 size_t size()
const noexcept;
48 bool operator==(
const Phases& data)
const;
50 template<
class Serializer>
57 std::bitset< NUM_PHASES_IN_ENUM > bits;
66 static Welldims serializationTestObject();
68 int maxConnPerWell()
const
73 int maxWellsPerGroup()
const
78 int maxGroupsInField()
const
83 int maxWellsInField()
const
88 int maxWellListsPrWell()
const
90 return this->nWlistPrWellMax;
93 int maxDynamicWellLists()
const
95 return this->nDynWlistMax;
98 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);
117 template<
class Serializer>
124 serializer(nWlistPrWellMax);
125 serializer(nDynWlistMax);
126 serializer(m_location);
134 int nWlistPrWellMax { 1 };
135 int nDynWlistMax { 1 };
136 std::optional<KeywordLocation> m_location;
147 int maxSegmentedWells()
const
149 return this->nSegWellMax;
152 int maxSegmentsPerWell()
const
154 return this->nSegmentMax;
157 int maxLateralBranchesPerWell()
const
159 return this->nLatBranchMax;
164 template<
class Serializer>
167 serializer(nSegWellMax);
168 serializer(nSegmentMax);
169 serializer(nLatBranchMax);
185 int maxNONodes()
const
187 return this->nMaxNoNodes;
190 int maxNoBranches()
const
192 return this->nMaxNoBranches;
195 int maxNoBranchesConToNode()
const
197 return this->nMaxNoBranchesConToNode;
204 template<
class Serializer>
207 serializer(nMaxNoNodes);
208 serializer(nMaxNoBranches);
209 serializer(nMaxNoBranchesConToNode);
215 int nMaxNoBranchesConToNode;
225 int maxAnalyticAquifers()
const
227 return this->maxNumAnalyticAquifers;
230 int maxAnalyticAquiferConnections()
const
232 return this->maxNumAnalyticAquiferConn;
235 template <
class Serializer>
238 serializer(this->maxNumAnalyticAquifers);
239 serializer(this->maxNumAnalyticAquiferConn);
243 int maxNumAnalyticAquifers;
244 int maxNumAnalyticAquiferConn;
299 template<
class Serializer>
302 serializer(activeHyst);
303 serializer(pcHystMod);
304 serializer(krHystMod);
305 serializer(modParamTrappedValue);
306 serializer(curvatureCapPrsValue);
311 bool activeHyst {
false };
314 int pcHystMod { -1 };
315 int krHystMod { -1 };
317 double modParamTrappedValue { 0.1 };
319 double curvatureCapPrsValue { 0.1 };
324 enum class ThreePhaseOilKrModel {
330 enum class KeywordFamily {
340 const ThreePhaseOilKrModel model,
341 const KeywordFamily family);
345 double minimumRelpermMobilityThreshold()
const
347 return this->tolcrit;
350 ThreePhaseOilKrModel krModel()
const
352 return this->krmodel;
355 KeywordFamily family()
const
357 return this->satfunc_family;
362 template<
class Serializer>
367 serializer(satfunc_family);
372 ThreePhaseOilKrModel krmodel = ThreePhaseOilKrModel::Default;
373 KeywordFamily satfunc_family = KeywordFamily::Undefined;
380 explicit Nupcol(
int min_value);
381 void update(
int value);
384 static Nupcol serializationTestObject();
385 bool operator==(
const Nupcol& data)
const;
387 template<
class Serializer>
389 serializer(this->nupcol_value);
390 serializer(this->min_nupcol);
405 int water_tracers()
const;
407 template<
class Serializer>
409 serializer(this->m_oil_tracers);
410 serializer(this->m_water_tracers);
411 serializer(this->m_gas_tracers);
412 serializer(this->m_env_tracers);
413 serializer(this->diffusion_control);
414 serializer(this->max_iter);
415 serializer(this->min_iter);
418 static Tracers serializationTestObject();
419 bool operator==(
const Tracers& data)
const;
426 bool diffusion_control;
439 static Runspec serializationTestObject();
441 std::time_t start_time()
const noexcept;
442 const UDQParams& udqParams()
const noexcept;
443 const Phases& phases()
const noexcept;
444 const Tabdims& tabdims()
const noexcept;
445 const Regdims& regdims()
const noexcept;
447 const Welldims& wellDimensions()
const noexcept;
449 const NetworkDims& networkDimensions()
const noexcept;
451 int eclPhaseMask( )
const noexcept;
453 const Actdims& actdims()
const noexcept;
455 const Nupcol& nupcol()
const noexcept;
456 const Tracers& tracers()
const;
457 bool co2Storage()
const noexcept;
458 bool micp()
const noexcept;
460 bool operator==(
const Runspec& data)
const;
461 static bool rst_cmp(
const Runspec& full_state,
const Runspec& rst_state);
463 template<
class Serializer>
466 serializer(this->m_start_time);
467 serializer(active_phases);
468 serializer(m_tabdims);
469 serializer(m_regdims);
470 serializer(endscale);
471 serializer(welldims);
472 serializer(wsegdims);
473 serializer(netwrkdims);
474 serializer(aquiferdims);
475 serializer(udq_params);
477 serializer(m_actdims);
478 serializer(m_sfuncctrl);
479 serializer(m_nupcol);
480 serializer(m_co2storage);
485 std::time_t m_start_time;
Definition: Actdims.hpp:30
Definition: Runspec.hpp:218
Definition: Runspec.hpp:250
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:178
Definition: Runspec.hpp:377
Definition: Runspec.hpp:37
Definition: Regdims.hpp:36
Definition: Runspec.hpp:434
Definition: Runspec.hpp:322
Class for (de-)serializing.
Definition: Serializer.hpp:84
Definition: Tabdims.hpp:36
Definition: Runspec.hpp:399
Definition: UDQParams.hpp:31
Definition: Runspec.hpp:139
Definition: Runspec.hpp:61
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition: Exceptions.hpp:30