27#ifndef OPM_TWO_PHASE_IMMISCIBLE_FLUID_SYSTEM_HPP
28#define OPM_TWO_PHASE_IMMISCIBLE_FLUID_SYSTEM_HPP
56template <
class Scalar,
class WettingPhase,
class NonwettingPhase>
58 :
public BaseFluidSystem<Scalar, TwoPhaseImmiscibleFluidSystem<Scalar, WettingPhase, NonwettingPhase> >
68 template <
class Evaluation>
85 static std::string_view
phaseName(
unsigned phaseIdx)
89 static const std::string_view name[] = {
93 return name[phaseIdx];
102 ? WettingPhase::isLiquid()
103 : NonwettingPhase::isLiquid();
113 ? WettingPhase::isCompressible()
114 : NonwettingPhase::isCompressible();
125 ? WettingPhase::isIdealGas()
126 : NonwettingPhase::isIdealGas();
156 return WettingPhase::name();
157 return NonwettingPhase::name();
168 ? WettingPhase::molarMass()
169 : NonwettingPhase::molarMass();
181 ? WettingPhase::criticalTemperature()
182 : NonwettingPhase::criticalTemperature();
194 ? WettingPhase::criticalPressure()
195 : NonwettingPhase::criticalPressure();
207 ? WettingPhase::acentricFactor()
208 : NonwettingPhase::acentricFactor();
220 assert(WettingPhase::isLiquid() || NonwettingPhase::isLiquid());
224 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar,
class ParamCacheEval = LhsEval>
225 static LhsEval
density(
const FluidState& fluidState,
231 const auto& temperature = decay<LhsEval>(fluidState.temperature(phaseIdx));
232 const auto& pressure = decay<LhsEval>(fluidState.pressure(phaseIdx));
234 return WettingPhase::density(temperature, pressure);
235 return NonwettingPhase::density(temperature, pressure);
239 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar,
class ParamCacheEval = LhsEval>
246 const auto& temperature = decay<LhsEval>(fluidState.temperature(phaseIdx));
247 const auto& pressure = decay<LhsEval>(fluidState.pressure(phaseIdx));
249 return WettingPhase::viscosity(temperature, pressure);
250 return NonwettingPhase::viscosity(temperature, pressure);
254 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar,
class ParamCacheEval = LhsEval>
263 if (phaseIdx == compIdx)
269 return std::numeric_limits<Scalar>::infinity();
273 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar,
class ParamCacheEval = LhsEval>
274 static LhsEval
enthalpy(
const FluidState& fluidState,
280 const auto& temperature = decay<LhsEval>(fluidState.temperature(phaseIdx));
281 const auto& pressure = decay<LhsEval>(fluidState.pressure(phaseIdx));
283 return WettingPhase::enthalpy(temperature, pressure);
284 return NonwettingPhase::enthalpy(temperature, pressure);
288 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar,
class ParamCacheEval = LhsEval>
295 const auto& temperature = decay<LhsEval>(fluidState.temperature(phaseIdx));
296 const auto& pressure = decay<LhsEval>(fluidState.pressure(phaseIdx));
298 return WettingPhase::thermalConductivity(temperature, pressure);
299 return NonwettingPhase::thermalConductivity(temperature, pressure);
303 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar,
class ParamCacheEval = LhsEval>
310 const auto& temperature = decay<LhsEval>(fluidState.temperature(phaseIdx));
311 const auto& pressure = decay<LhsEval>(fluidState.pressure(phaseIdx));
313 return WettingPhase::heatCapacity(temperature, pressure);
314 return NonwettingPhase::heatCapacity(temperature, pressure);
The base class for all fluid systems.
Represents the gas phase of a single (pseudo-) component.
Represents all relevant thermodynamic quantities of a multi-phase, multi-component fluid system assum...
Represents the liquid phase of a single (pseudo-) component.
A parameter cache which does nothing.
The base class for all fluid systems.
Definition BaseFluidSystem.hpp:47
Scalar Scalar
The type used for scalar quantities.
Definition BaseFluidSystem.hpp:52
A parameter cache which does nothing.
Definition NullParameterCache.hpp:40
A fluid system for two-phase models assuming immiscibility and thermodynamic equilibrium.
Definition TwoPhaseImmiscibleFluidSystem.hpp:59
static bool isIdealMixture(unsigned)
Returns true if and only if a fluid phase is assumed to be an ideal mixture.
Definition TwoPhaseImmiscibleFluidSystem.hpp:130
static Scalar acentricFactor(unsigned compIdx)
The acentric factor of a component [].
Definition TwoPhaseImmiscibleFluidSystem.hpp:201
static LhsEval thermalConductivity(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &, unsigned phaseIdx)
Thermal conductivity of a fluid phase [W/(m K)].
Definition TwoPhaseImmiscibleFluidSystem.hpp:289
static void init()
Initialize the fluid system's static parameters.
Definition TwoPhaseImmiscibleFluidSystem.hpp:216
static bool isLiquid(unsigned phaseIdx)
Return whether a phase is liquid.
Definition TwoPhaseImmiscibleFluidSystem.hpp:97
static LhsEval fugacityCoefficient(const FluidState &, const ParameterCache< ParamCacheEval > &, unsigned phaseIdx, unsigned compIdx)
Calculate the fugacity coefficient [Pa] of an individual component in a fluid phase.
Definition TwoPhaseImmiscibleFluidSystem.hpp:255
static Scalar molarMass(unsigned compIdx)
Return the molar mass of a component in [kg/mol].
Definition TwoPhaseImmiscibleFluidSystem.hpp:161
static std::string_view componentName(unsigned compIdx)
Return the human readable name of a component.
Definition TwoPhaseImmiscibleFluidSystem.hpp:151
static bool isIdealGas(unsigned phaseIdx)
Returns true if and only if a fluid phase is assumed to be an ideal gas.
Definition TwoPhaseImmiscibleFluidSystem.hpp:118
static const int numPhases
Number of fluid phases in the fluid system.
Definition TwoPhaseImmiscibleFluidSystem.hpp:77
static const int wettingCompIdx
Index of the wetting phase's component.
Definition TwoPhaseImmiscibleFluidSystem.hpp:146
static const int nonWettingCompIdx
Index of the non-wetting phase's component.
Definition TwoPhaseImmiscibleFluidSystem.hpp:148
static LhsEval enthalpy(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &, unsigned phaseIdx)
Given a phase's composition, temperature, pressure and density, calculate its specific enthalpy [J/kg...
Definition TwoPhaseImmiscibleFluidSystem.hpp:274
static LhsEval viscosity(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &, unsigned phaseIdx)
Calculate the dynamic viscosity of a fluid phase [Pa*s].
Definition TwoPhaseImmiscibleFluidSystem.hpp:240
static std::string_view phaseName(unsigned phaseIdx)
Return the human readable name of a fluid phase.
Definition TwoPhaseImmiscibleFluidSystem.hpp:85
static Scalar criticalTemperature(unsigned compIdx)
Critical temperature of a component [K].
Definition TwoPhaseImmiscibleFluidSystem.hpp:175
static LhsEval heatCapacity(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &, unsigned phaseIdx)
Specific isobaric heat capacity of a fluid phase [J/kg].
Definition TwoPhaseImmiscibleFluidSystem.hpp:304
static const int nonWettingPhaseIdx
Index of the non-wetting phase.
Definition TwoPhaseImmiscibleFluidSystem.hpp:82
static const int numComponents
Number of chemical species in the fluid system.
Definition TwoPhaseImmiscibleFluidSystem.hpp:143
static LhsEval density(const FluidState &fluidState, const ParameterCache< ParamCacheEval > &, unsigned phaseIdx)
Calculate the density [kg/m^3] of a fluid phase.
Definition TwoPhaseImmiscibleFluidSystem.hpp:225
static const int wettingPhaseIdx
Index of the wetting phase.
Definition TwoPhaseImmiscibleFluidSystem.hpp:80
static Scalar criticalPressure(unsigned compIdx)
Critical pressure of a component [Pa].
Definition TwoPhaseImmiscibleFluidSystem.hpp:188
static bool isCompressible(unsigned phaseIdx)
Returns true if and only if a fluid phase is assumed to be compressible.
Definition TwoPhaseImmiscibleFluidSystem.hpp:107
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition Exceptions.hpp:30
Definition TwoPhaseImmiscibleFluidSystem.hpp:70