26#ifndef OPM_GENERIC_OIL_GAS_WATER_FLUIDSYSTEM_HPP
27#define OPM_GENERIC_OIL_GAS_WATER_FLUIDSYSTEM_HPP
29#include <opm/common/OpmLog/OpmLog.hpp>
32#include <opm/input/eclipse/EclipseState/EclipseState.hpp>
33#include <opm/input/eclipse/Schedule/Schedule.hpp>
36#include <opm/material/eos/CubicEOS.hpp>
47#include <fmt/format.h>
59 template<
class Scalar,
int NumComp,
bool enableWater>
61 :
public BaseFluidSystem<Scalar, GenericOilGasWaterFluidSystem<Scalar, NumComp, enableWater> > {
64 static constexpr bool waterEnabled = enableWater;
65 static constexpr int numPhases = enableWater ? 3 : 2;
66 static constexpr int numComponents = NumComp;
67 static constexpr int numMisciblePhases = 2;
70 static constexpr int numMiscibleComponents = NumComp;
72 static constexpr int oilPhaseIdx = 0;
73 static constexpr int gasPhaseIdx = 1;
74 static constexpr int waterPhaseIdx = 2;
76 static constexpr int waterCompIdx = -1;
77 static constexpr int oilCompIdx = 0;
78 static constexpr int gasCompIdx = 1;
79 static constexpr int compositionSwitchIdx = -1;
81 template <
class ValueType>
98 molar_mass(molar_mass_),
99 critic_temp(critic_temp_),
100 critic_pres(critic_pres_),
101 critic_vol(critic_vol_),
102 acentric_factor(acentric_factor_)
106 static bool phaseIsActive(
unsigned phaseIdx)
108 if constexpr (enableWater) {
109 assert(phaseIdx < numPhases);
113 if (phaseIdx == waterPhaseIdx)
115 assert(phaseIdx < numPhases);
120 template<
typename Param>
121 static void addComponent(
const Param& param)
124 if (component_param_.size() < numComponents) {
125 component_param_.push_back(param);
128 const std::string msg = fmt::format(
"The fluid system has reached its maximum capacity of {} components,"
129 "the component '{}' will not be added.", NumComp, param.name);
139 static void initFromState(
const EclipseState& eclState,
const Schedule& schedule)
142 const auto& comp_config = eclState.compositionalConfig();
144 using FluidSystem = GenericOilGasWaterFluidSystem<Scalar, NumComp, enableWater>;
145 const std::size_t num_comps = comp_config.numComps();
147 assert(num_comps == NumComp);
148 const auto& names = comp_config.compName();
150 const auto& molar_weight = comp_config.molecularWeights(0);
151 const auto& acentric_factor = comp_config.acentricFactors(0);
152 const auto& critic_pressure = comp_config.criticalPressure(0);
153 const auto& critic_temp = comp_config.criticalTemperature(0);
154 const auto& critic_volume = comp_config.criticalVolume(0);
156 using CompParm =
typename FluidSystem::ComponentParam;
157 for (std::size_t c = 0; c < num_comps; ++c) {
159 FluidSystem::addComponent(CompParm{names[c], molar_weight[c], critic_temp[c], critic_pressure[c],
160 critic_volume[c] * 1.e3, acentric_factor[c]});
162 FluidSystem::printComponentParams();
163 interaction_coefficients_ = comp_config.binaryInteractionCoefficient(0);
166 waterPvt_->initFromState(eclState, schedule);
173 waterPvt_ = std::make_shared<WaterPvt>();
174 component_param_.reserve(numComponents);
181 { waterPvt_ = std::move(pvtObj); }
190 assert(isConsistent());
191 assert(compIdx < numComponents);
193 return component_param_[compIdx].acentric_factor;
202 assert(isConsistent());
203 assert(compIdx < numComponents);
205 return component_param_[compIdx].critic_temp;
214 assert(isConsistent());
215 assert(compIdx < numComponents);
217 return component_param_[compIdx].critic_pres;
226 assert(isConsistent());
227 assert(compIdx < numComponents);
229 return component_param_[compIdx].critic_vol;
235 assert(isConsistent());
236 assert(compIdx < numComponents);
238 return component_param_[compIdx].molar_mass;
247 assert(isConsistent());
248 assert(comp1Idx < numComponents);
249 assert(comp2Idx < numComponents);
250 if (interaction_coefficients_.empty() || comp2Idx == comp1Idx) {
254 const auto [column, row] = std::minmax(comp1Idx, comp2Idx);
255 const unsigned index = (row * (row - 1) / 2 + column);
256 return interaction_coefficients_[index];
262 static const std::string_view name[] = {
"o",
266 assert(phaseIdx < numPhases);
267 return name[phaseIdx];
273 assert(isConsistent());
274 assert(compIdx < numComponents);
276 return component_param_[compIdx].name;
282 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar,
class ParamCacheEval = LhsEval>
283 static LhsEval
density(
const FluidState& fluidState,
287 assert(isConsistent());
288 assert(phaseIdx < numPhases);
290 if (phaseIdx == oilPhaseIdx || phaseIdx == gasPhaseIdx) {
291 return decay<LhsEval>(fluidState.averageMolarMass(phaseIdx) / paramCache.
molarVolume(phaseIdx));
294 const LhsEval& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
295 const LhsEval& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
296 const Scalar& rho_w_ref = waterPvt_->waterReferenceDensity(0);
297 const LhsEval& bw = waterPvt_->inverseFormationVolumeFactor(0, T, p, LhsEval(0.0), LhsEval(0.0));
298 return rho_w_ref * bw;
305 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar,
class ParamCacheEval = LhsEval>
310 assert(isConsistent());
311 assert(phaseIdx < numPhases);
313 if (phaseIdx == oilPhaseIdx || phaseIdx == gasPhaseIdx) {
315 return decay<LhsEval>(ViscosityModel::LBC(fluidState, paramCache, phaseIdx));
318 const LhsEval& p = decay<LhsEval>(fluidState.pressure(phaseIdx));
319 const LhsEval& T = decay<LhsEval>(fluidState.temperature(phaseIdx));
320 return waterPvt_->viscosity(0, T, p, LhsEval(0.0), LhsEval(0.0));
325 template <
class Flu
idState,
class LhsEval =
typename Flu
idState::Scalar,
class ParamCacheEval = LhsEval>
331 if (waterEnabled && phaseIdx ==
static_cast<unsigned int>(waterPhaseIdx))
334 assert(isConsistent());
335 assert(phaseIdx < numPhases);
336 assert(compIdx < numComponents);
338 return decay<LhsEval>(CubicEOS::computeFugacityCoefficient(fluidState, paramCache, phaseIdx, compIdx));
345 assert(phaseIdx < numPhases);
353 assert(phaseIdx < numPhases);
361 assert(phaseIdx < numPhases);
363 return (phaseIdx == oilPhaseIdx);
369 assert(phaseIdx < numPhases);
371 return (phaseIdx == gasPhaseIdx);
375 template <
class LhsEval>
376 static LhsEval convertXwGToxwG(
const LhsEval&,
unsigned)
378 assert(
false &&
"convertXwGToxwG not implemented for GenericOilGasWaterFluidSystem!");
381 template <
class LhsEval>
382 static LhsEval convertXoGToxoG(
const LhsEval&,
unsigned)
384 assert(
false &&
"convertXoGToxoG not implemented for GenericOilGasWaterFluidSystem!");
388 template <
class LhsEval>
389 static LhsEval convertxoGToXoG(
const LhsEval&,
unsigned)
391 assert(
false &&
"convertxoGToXoG not implemented for GenericOilGasWaterFluidSystem!");
395 template <
class LhsEval>
396 static LhsEval convertXgOToxgO(
const LhsEval&,
unsigned)
398 assert(
false &&
"convertXgOToxgO not implemented for GenericOilGasWaterFluidSystem!");
402 template <
class LhsEval>
403 static LhsEval convertRswToXwG(
const LhsEval&,
unsigned)
405 assert(
false &&
"convertRswToXwG not implemented for GenericOilGasWaterFluidSystem!");
409 template <
class LhsEval>
410 static LhsEval convertRvwToXgW(
const LhsEval&,
unsigned)
412 assert(
false &&
"convertRvwToXgW not implemented for GenericOilGasWaterFluidSystem!");
416 template <
class LhsEval>
417 static LhsEval convertXgWToxgW(
const LhsEval&,
unsigned)
419 assert(
false &&
"convertXgWToxgW not implemented for GenericOilGasWaterFluidSystem!");
423 static bool enableDissolvedGas()
428 static bool enableDissolvedGasInWater()
433 static bool enableVaporizedWater()
438 static bool enableVaporizedOil()
444 static bool isConsistent() {
445 return component_param_.size() == NumComp;
448 static std::vector<ComponentParam> component_param_;
449 static std::vector<Scalar> interaction_coefficients_;
450 static std::shared_ptr<WaterPvt> waterPvt_;
453 static std::string printComponentParams() {
454 std::string result =
"Components Information:\n";
455 for (
const auto& param : component_param_) {
456 result += fmt::format(
"Name: {}\n", param.name);
457 result += fmt::format(
"Molar Mass: {} g/mol\n", param.molar_mass);
458 result += fmt::format(
"Critical Temperature: {} K\n", param.critic_temp);
459 result += fmt::format(
"Critical Pressure: {} Pascal\n", param.critic_pres);
460 result += fmt::format(
"Critical Volume: {} m^3/kmol\n", param.critic_vol);
461 result += fmt::format(
"Acentric Factor: {}\n", param.acentric_factor);
462 result +=
"---------------------------------\n";
468 template <
class Scalar,
int NumComp,
bool enableWater>
469 std::vector<typename GenericOilGasWaterFluidSystem<Scalar, NumComp, enableWater>::ComponentParam>
470 GenericOilGasWaterFluidSystem<Scalar, NumComp, enableWater>::component_param_;
472 template <
class Scalar,
int NumComp,
bool enableWater>
474 GenericOilGasWaterFluidSystem<Scalar, NumComp, enableWater>::interaction_coefficients_;
476 template <
class Scalar,
int NumComp,
bool enableWater>
477 std::shared_ptr<WaterPvtMultiplexer<Scalar> >
478 GenericOilGasWaterFluidSystem<Scalar, NumComp, enableWater>::waterPvt_;
The base class for all fluid systems.
Specifies the parameter cache used by the SPE-5 fluid system.
This class represents the Pressure-Volume-Temperature relations of the water phase in the black-oil m...
The base class for all fluid systems.
Definition BaseFluidSystem.hpp:43
Scalar Scalar
The type used for scalar quantities.
Definition BaseFluidSystem.hpp:48
Definition CubicEOS.hpp:34
A two phase system that can contain NumComp components.
Definition GenericOilGasWaterFluidSystem.hpp:61
static bool isCompressible(unsigned phaseIdx)
Returns true if and only if a fluid phase is assumed to be compressible.
Definition GenericOilGasWaterFluidSystem.hpp:343
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 GenericOilGasWaterFluidSystem.hpp:326
static Scalar acentricFactor(unsigned compIdx)
The acentric factor of a component [].
Definition GenericOilGasWaterFluidSystem.hpp:188
static std::string_view phaseName(unsigned phaseIdx)
Return the human readable name of a fluid phase.
Definition GenericOilGasWaterFluidSystem.hpp:260
static bool isIdealGas(unsigned phaseIdx)
Returns true if and only if a fluid phase is assumed to be an ideal gas.
Definition GenericOilGasWaterFluidSystem.hpp:367
static Scalar criticalTemperature(unsigned compIdx)
Critical temperature of a component [K].
Definition GenericOilGasWaterFluidSystem.hpp:200
static LhsEval density(const FluidState &fluidState, const ParameterCache< ParamCacheEval > ¶mCache, unsigned phaseIdx)
Calculate the density [kg/m^3] of a fluid phase.
Definition GenericOilGasWaterFluidSystem.hpp:283
static Scalar interactionCoefficient(unsigned comp1Idx, unsigned comp2Idx)
Returns the interaction coefficient for two components.
Definition GenericOilGasWaterFluidSystem.hpp:245
static Scalar criticalPressure(unsigned compIdx)
Critical pressure of a component [Pa].
Definition GenericOilGasWaterFluidSystem.hpp:212
static bool isIdealMixture(unsigned phaseIdx)
Returns true if and only if a fluid phase is assumed to be an ideal mixture.
Definition GenericOilGasWaterFluidSystem.hpp:351
static Scalar molarMass(unsigned compIdx)
Return the molar mass of a component in [kg/mol].
Definition GenericOilGasWaterFluidSystem.hpp:233
static std::string_view componentName(unsigned compIdx)
Return the human readable name of a component.
Definition GenericOilGasWaterFluidSystem.hpp:271
static bool isLiquid(unsigned phaseIdx)
Return whether a phase is liquid.
Definition GenericOilGasWaterFluidSystem.hpp:359
static Scalar criticalVolume(unsigned compIdx)
Critical volume of a component [m3].
Definition GenericOilGasWaterFluidSystem.hpp:224
static void setWaterPvt(std::shared_ptr< WaterPvt > pvtObj)
Set the pressure-volume-saturation (PVT) relations for the water phase.
Definition GenericOilGasWaterFluidSystem.hpp:180
static LhsEval viscosity(const FluidState &fluidState, const ParameterCache< ParamCacheEval > ¶mCache, unsigned phaseIdx)
Calculate the dynamic viscosity of a fluid phase [Pa*s].
Definition GenericOilGasWaterFluidSystem.hpp:306
Specifies the parameter cache used by the SPE-5 fluid system.
Definition PTFlashParameterCache.hpp:51
Scalar molarVolume(unsigned phaseIdx) const
Returns the molar volume of a phase [m^3/mol].
Definition PTFlashParameterCache.hpp:283
This class represents the Pressure-Volume-Temperature relations of the water phase in the black-oil m...
Definition WaterPvtMultiplexer.hpp:90
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition Exceptions.hpp:30
Definition GenericOilGasWaterFluidSystem.hpp:87