26#ifndef OPM_GENERICOILGASFLUIDSYSTEM_HPP
27#define OPM_GENERICOILGASFLUIDSYSTEM_HPP
29#include <opm/common/OpmLog/OpmLog.hpp>
39#include <fmt/format.h>
52 template<
class Scalar,
int NumComp>
56 static const int numPhases=2;
57 static const int numComponents = NumComp;
58 static const int numMisciblePhases=2;
61 static const int numMiscibleComponents = NumComp;
63 static constexpr int oilPhaseIdx = 0;
64 static constexpr int gasPhaseIdx = 1;
66 template <
class ValueType>
82 molar_mass(molar_mass_),
83 critic_temp(critic_temp_),
84 critic_pres(critic_pres_),
85 critic_vol(critic_vol_),
86 acentric_factor(acentric_factor_)
90 template<
typename Param>
91 static void addComponent(
const Param& param)
94 if (component_param_.size() < numComponents) {
95 component_param_.push_back(param);
98 const std::string msg = fmt::format(
"The fluid system has reached its maximum capacity of {} components,"
99 "the component '{}' will not be added.", NumComp, param.name);
107 component_param_.reserve(numComponents);
117 assert(isConsistent());
118 assert(compIdx < numComponents);
120 return component_param_[compIdx].acentric_factor;
129 assert(isConsistent());
130 assert(compIdx < numComponents);
132 return component_param_[compIdx].critic_temp;
141 assert(isConsistent());
142 assert(compIdx < numComponents);
144 return component_param_[compIdx].critic_pres;
153 assert(isConsistent());
154 assert(compIdx < numComponents);
156 return component_param_[compIdx].critic_vol;
162 assert(isConsistent());
163 assert(compIdx < numComponents);
165 return component_param_[compIdx].molar_mass;
174 assert(isConsistent());
182 static const std::string_view name[] = {
"o",
185 assert(phaseIdx < numPhases);
186 return name[phaseIdx];
192 assert(isConsistent());
193 assert(compIdx < numComponents);
195 return component_param_[compIdx].name;
201 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar,
class ParamCacheEval = LhsEval>
202 static LhsEval
density(
const FluidState& fluidState,
206 assert(isConsistent());
207 assert(phaseIdx < numPhases);
209 if (phaseIdx == oilPhaseIdx || phaseIdx == gasPhaseIdx) {
210 return decay<LhsEval>(fluidState.averageMolarMass(phaseIdx) / paramCache.
molarVolume(phaseIdx));
217 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar,
class ParamCacheEval = LhsEval>
222 assert(isConsistent());
223 assert(phaseIdx < numPhases);
225 return decay<LhsEval>(ViscosityModel::LBC(fluidState, paramCache, phaseIdx));
229 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar,
class ParamCacheEval = LhsEval>
235 assert(isConsistent());
236 assert(phaseIdx < numPhases);
237 assert(compIdx < numComponents);
245 assert(phaseIdx < numPhases);
253 assert(phaseIdx < numPhases);
261 assert(phaseIdx < numPhases);
263 return (phaseIdx == 0);
269 assert(phaseIdx < numPhases);
271 return (phaseIdx == 1);
274 static bool isConsistent() {
275 return component_param_.size() == NumComp;
278 static std::vector<ComponentParam> component_param_;
281 template <
class Scalar,
int NumComp>
282 std::vector<typename GenericOilGasFluidSystem<Scalar, NumComp>::ComponentParam>
283 GenericOilGasFluidSystem<Scalar, NumComp>::component_param_;
The base class for all fluid systems.
Specifies the parameter cache used by the SPE-5 fluid system.
The base class for all fluid systems.
Definition BaseFluidSystem.hpp:47
Scalar Scalar
The type used for scalar quantities.
Definition BaseFluidSystem.hpp:52
A two phase system that can contain NumComp components.
Definition GenericOilGasFluidSystem.hpp:53
static bool isIdealMixture(unsigned phaseIdx)
Returns true if and only if a fluid phase is assumed to be an ideal mixture.
Definition GenericOilGasFluidSystem.hpp:251
static LhsEval fugacityCoefficient(const FluidState &fluidState, const ParameterCache< ParamCacheEval > ¶mCache, unsigned phaseIdx, unsigned compIdx)
Calculate the fugacity coefficient [Pa] of an individual component in a fluid phase.
Definition GenericOilGasFluidSystem.hpp:230
static Scalar interactionCoefficient(unsigned, unsigned)
Returns the interaction coefficient for two components.
Definition GenericOilGasFluidSystem.hpp:172
static Scalar criticalTemperature(unsigned compIdx)
Critical temperature of a component [K].
Definition GenericOilGasFluidSystem.hpp:127
static bool isCompressible(unsigned phaseIdx)
Returns true if and only if a fluid phase is assumed to be compressible.
Definition GenericOilGasFluidSystem.hpp:243
static Scalar criticalVolume(unsigned compIdx)
Critical volume of a component [m3].
Definition GenericOilGasFluidSystem.hpp:151
static Scalar criticalPressure(unsigned compIdx)
Critical pressure of a component [Pa].
Definition GenericOilGasFluidSystem.hpp:139
static bool isIdealGas(unsigned phaseIdx)
Returns true if and only if a fluid phase is assumed to be an ideal gas.
Definition GenericOilGasFluidSystem.hpp:267
static LhsEval density(const FluidState &fluidState, const ParameterCache< ParamCacheEval > ¶mCache, unsigned phaseIdx)
Calculate the density [kg/m^3] of a fluid phase.
Definition GenericOilGasFluidSystem.hpp:202
static bool isLiquid(unsigned phaseIdx)
Return whether a phase is liquid.
Definition GenericOilGasFluidSystem.hpp:259
static std::string_view phaseName(unsigned phaseIdx)
Return the human readable name of a fluid phase.
Definition GenericOilGasFluidSystem.hpp:180
static Scalar acentricFactor(unsigned compIdx)
The acentric factor of a component [].
Definition GenericOilGasFluidSystem.hpp:115
static Scalar molarMass(unsigned compIdx)
Return the molar mass of a component in [kg/mol].
Definition GenericOilGasFluidSystem.hpp:160
static LhsEval viscosity(const FluidState &fluidState, const ParameterCache< ParamCacheEval > ¶mCache, unsigned phaseIdx)
Calculate the dynamic viscosity of a fluid phase [Pa*s].
Definition GenericOilGasFluidSystem.hpp:218
static std::string_view componentName(unsigned compIdx)
Return the human readable name of a component.
Definition GenericOilGasFluidSystem.hpp:190
Specifies the parameter cache used by the SPE-5 fluid system.
Definition PTFlashParameterCache.hpp:48
Scalar molarVolume(unsigned phaseIdx) const
Returns the molar volume of a phase [m^3/mol].
Definition PTFlashParameterCache.hpp:199
Implements the Peng-Robinson equation of state for a mixture.
Definition PengRobinsonMixture.hpp:41
static LhsEval computeFugacityCoefficient(const FluidState &fs, const Params ¶ms, unsigned phaseIdx, unsigned compIdx)
Returns the fugacity coefficient of an individual component in the phase.
Definition PengRobinsonMixture.hpp:74
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition Exceptions.hpp:30
Definition GenericOilGasFluidSystem.hpp:71