28#ifndef OPM_FLUID_STATE_COMPOSITION_MODULES_HPP
29#define OPM_FLUID_STATE_COMPOSITION_MODULES_HPP
44template <
class Scalar,
49 enum { numPhases = FluidSystem::numPhases };
50 enum { numComponents = FluidSystem::numComponents };
55 for (
int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
56 for (
int compIdx = 0; compIdx < numComponents; ++compIdx)
57 moleFraction_[phaseIdx][compIdx] = 0.0;
62 Valgrind::SetDefined(moleFraction_);
63 Valgrind::SetUndefined(averageMolarMass_);
64 Valgrind::SetUndefined(sumMoleFractions_);
70 const Scalar&
moleFraction(
unsigned phaseIdx,
unsigned compIdx)
const
71 {
return moleFraction_[phaseIdx][compIdx]; }
77 {
return totalModelFractions_[compIdx]; }
85 abs(sumMoleFractions_[phaseIdx])
86 *moleFraction_[phaseIdx][compIdx]
87 *FluidSystem::molarMass(compIdx)
88 / max(1e-40, abs(averageMolarMass_[phaseIdx]));
100 {
return averageMolarMass_[phaseIdx]; }
111 Scalar
molarity(
unsigned phaseIdx,
unsigned compIdx)
const
112 {
return asImp_().molarDensity(phaseIdx)*
moleFraction(phaseIdx, compIdx); }
121 Valgrind::CheckDefined(value);
122 Valgrind::SetUndefined(sumMoleFractions_[phaseIdx]);
123 Valgrind::SetUndefined(averageMolarMass_[phaseIdx]);
124 Valgrind::SetUndefined(moleFraction_[phaseIdx][compIdx]);
126 moleFraction_[phaseIdx][compIdx] = value;
129 sumMoleFractions_[phaseIdx] = 0.0;
130 averageMolarMass_[phaseIdx] = 0.0;
131 for (
unsigned compJIdx = 0; compJIdx < numComponents; ++compJIdx) {
132 sumMoleFractions_[phaseIdx] += moleFraction_[phaseIdx][compJIdx];
133 averageMolarMass_[phaseIdx] += moleFraction_[phaseIdx][compJIdx]*FluidSystem::molarMass(compJIdx);
142 Valgrind::CheckDefined(value);
143 totalModelFractions_[compIdx] = value;
146 void setCompressFactor(
unsigned phaseIdx,
const Scalar& value) {
147 Valgrind::CheckDefined(value);
148 Z_[phaseIdx] = value;
151 Scalar compressFactor(
unsigned phaseIdx)
const {
159 template <
class Flu
idState>
162 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
163 averageMolarMass_[phaseIdx] = 0;
164 sumMoleFractions_[phaseIdx] = 0;
165 for (
unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
166 moleFraction_[phaseIdx][compIdx] =
167 decay<Scalar>(fs.moleFraction(phaseIdx, compIdx));
169 averageMolarMass_[phaseIdx] += moleFraction_[phaseIdx][compIdx]*FluidSystem::molarMass(compIdx);
170 sumMoleFractions_[phaseIdx] += moleFraction_[phaseIdx][compIdx];
185 Valgrind::CheckDefined(moleFraction_);
186 Valgrind::CheckDefined(averageMolarMass_);
187 Valgrind::CheckDefined(sumMoleFractions_);
188 Valgrind::CheckDefined(K_);
189 Valgrind::CheckDefined(L_);
192 const Scalar& K(
unsigned compIdx)
const
208 const Scalar&
L()
const
227 const auto& acf = FluidSystem::acentricFactor(compIdx);
228 const auto& T_crit = FluidSystem::criticalTemperature(compIdx);
229 const auto& T = asImp_().temperature(0);
230 const auto& p_crit = FluidSystem::criticalPressure(compIdx);
231 const auto& p = asImp_().pressure(0);
233 const auto tmp = exp(5.37 * (1+acf) * (1-T_crit/T)) * (p_crit/p);
238 const Implementation& asImp_()
const
240 return *
static_cast<const Implementation*
>(
this);
243 std::array<std::array<Scalar,numComponents>,numPhases> moleFraction_{};
244 std::array<Scalar,numPhases> averageMolarMass_{};
245 std::array<Scalar,numPhases> sumMoleFractions_{};
246 std::array<Scalar, numComponents> totalModelFractions_{};
247 std::array<Scalar,numPhases> Z_{};
248 std::array<Scalar,numComponents> K_{};
256template <
class Scalar,
258 class Implementation>
261 enum { numPhases = FluidSystem::numPhases };
264 enum { numComponents = FluidSystem::numComponents };
265 static_assert(
static_cast<int>(numPhases) ==
static_cast<int>(numComponents),
266 "The number of phases must be the same as the number of (pseudo-) components if you assume immiscibility");
274 {
return (phaseIdx == compIdx)?1.0:0.0; }
280 {
return (phaseIdx == compIdx)?1.0:0.0; }
291 {
return FluidSystem::molarMass(phaseIdx); }
302 Scalar
molarity(
unsigned phaseIdx,
unsigned compIdx)
const
303 {
return asImp_().molarDensity(phaseIdx)*
moleFraction(phaseIdx, compIdx); }
309 template <
class Flu
idState>
325 const Implementation& asImp_()
const
326 {
return *
static_cast<const Implementation*
>(
this); }
333template <
class Scalar>
337 enum { numComponents = 0 };
345 {
throw std::logic_error(
"Mole fractions are not provided by this fluid state"); }
351 {
throw std::logic_error(
"Mass fractions are not provided by this fluid state"); }
362 {
throw std::logic_error(
"Mean molar masses are not provided by this fluid state"); }
374 {
throw std::logic_error(
"Molarities are not provided by this fluid state"); }
Some templates to wrap the valgrind client request macros.
Module for the modular fluid state which stores the phase compositions explicitly in terms of mole fr...
Definition FluidStateCompositionModules.hpp:48
void setLvalue(const Scalar &value)
Set the L value [-].
Definition FluidStateCompositionModules.hpp:216
Scalar molarity(unsigned phaseIdx, unsigned compIdx) const
The concentration of a component in a phase [mol/m^3].
Definition FluidStateCompositionModules.hpp:111
const Scalar & L() const
The L value of a composition [-].
Definition FluidStateCompositionModules.hpp:208
Scalar wilsonK_(unsigned compIdx) const
Wilson formula to calculate K.
Definition FluidStateCompositionModules.hpp:225
void assign(const FluidState &fs)
Retrieve all parameters from an arbitrary fluid state.
Definition FluidStateCompositionModules.hpp:160
const Scalar & averageMolarMass(unsigned phaseIdx) const
The mean molar mass of a fluid phase [kg/mol].
Definition FluidStateCompositionModules.hpp:99
const Scalar & moleFraction(unsigned phaseIdx, unsigned compIdx) const
The mole fraction of a component in a phase [].
Definition FluidStateCompositionModules.hpp:70
const Scalar & moleFraction(unsigned compIdx) const
The total mole fraction of a component [].
Definition FluidStateCompositionModules.hpp:76
void setMoleFraction(unsigned compIdx, const Scalar &value)
Set the total mole fraction of a component.
Definition FluidStateCompositionModules.hpp:140
void setMoleFraction(unsigned phaseIdx, unsigned compIdx, const Scalar &value)
Set the mole fraction of a component in a phase [] and update the average molar mass [kg/mol] accordi...
Definition FluidStateCompositionModules.hpp:119
void setKvalue(unsigned compIdx, const Scalar &value)
Set the K value of a component [-].
Definition FluidStateCompositionModules.hpp:200
Scalar massFraction(unsigned phaseIdx, unsigned compIdx) const
The mass fraction of a component in a phase [].
Definition FluidStateCompositionModules.hpp:82
void checkDefined() const
Make sure that all attributes are defined.
Definition FluidStateCompositionModules.hpp:183
Module for the modular fluid state which provides the phase compositions assuming immiscibility.
Definition FluidStateCompositionModules.hpp:260
Scalar massFraction(unsigned phaseIdx, unsigned compIdx) const
The mass fraction of a component in a phase [].
Definition FluidStateCompositionModules.hpp:279
Scalar molarity(unsigned phaseIdx, unsigned compIdx) const
The concentration of a component in a phase [mol/m^3].
Definition FluidStateCompositionModules.hpp:302
void checkDefined() const
Make sure that all attributes are defined.
Definition FluidStateCompositionModules.hpp:321
void assign(const FluidState &)
Retrieve all parameters from an arbitrary fluid state.
Definition FluidStateCompositionModules.hpp:310
Scalar moleFraction(unsigned phaseIdx, unsigned compIdx) const
The mole fraction of a component in a phase [].
Definition FluidStateCompositionModules.hpp:273
Scalar averageMolarMass(unsigned phaseIdx) const
The mean molar mass of a fluid phase [kg/mol].
Definition FluidStateCompositionModules.hpp:290
Module for the modular fluid state which does not store the compositions but throws std::logic_error ...
Definition FluidStateCompositionModules.hpp:335
Scalar moleFraction(unsigned, unsigned) const
The mole fraction of a component in a phase [].
Definition FluidStateCompositionModules.hpp:344
Scalar averageMolarMass(unsigned) const
The mean molar mass of a fluid phase [kg/mol].
Definition FluidStateCompositionModules.hpp:361
Scalar molarity(unsigned, unsigned) const
The concentration of a component in a phase [mol/m^3].
Definition FluidStateCompositionModules.hpp:373
void checkDefined() const
Make sure that all attributes are defined.
Definition FluidStateCompositionModules.hpp:384
Scalar massFraction(unsigned, unsigned) const
The mass fraction of a component in a phase [].
Definition FluidStateCompositionModules.hpp:350
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition Exceptions.hpp:30