28 #ifndef OPM_FLUID_STATE_COMPOSITION_MODULES_HPP
29 #define OPM_FLUID_STATE_COMPOSITION_MODULES_HPP
44 template <
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;
59 Valgrind::SetDefined(moleFraction_);
60 Valgrind::SetUndefined(averageMolarMass_);
61 Valgrind::SetUndefined(sumMoleFractions_);
67 const Scalar&
moleFraction(
unsigned phaseIdx,
unsigned compIdx)
const
68 {
return moleFraction_[phaseIdx][compIdx]; }
76 abs(sumMoleFractions_[phaseIdx])
77 *moleFraction_[phaseIdx][compIdx]
78 *FluidSystem::molarMass(compIdx)
79 / max(1e-40, abs(averageMolarMass_[phaseIdx]));
91 {
return averageMolarMass_[phaseIdx]; }
102 Scalar
molarity(
unsigned phaseIdx,
unsigned compIdx)
const
103 {
return asImp_().molarDensity(phaseIdx)*
moleFraction(phaseIdx, compIdx); }
112 Valgrind::CheckDefined(value);
113 Valgrind::SetUndefined(sumMoleFractions_[phaseIdx]);
114 Valgrind::SetUndefined(averageMolarMass_[phaseIdx]);
115 Valgrind::SetUndefined(moleFraction_[phaseIdx][compIdx]);
117 moleFraction_[phaseIdx][compIdx] = value;
120 sumMoleFractions_[phaseIdx] = 0.0;
121 averageMolarMass_[phaseIdx] = 0.0;
122 for (
unsigned compJIdx = 0; compJIdx < numComponents; ++compJIdx) {
123 sumMoleFractions_[phaseIdx] += moleFraction_[phaseIdx][compJIdx];
124 averageMolarMass_[phaseIdx] += moleFraction_[phaseIdx][compJIdx]*FluidSystem::molarMass(compJIdx);
128 void setCompressFactor(
unsigned phaseIdx,
const Scalar& value) {
129 Valgrind::CheckDefined(value);
130 Z_[phaseIdx] = value;
133 Scalar compressFactor(
unsigned phaseIdx)
const {
141 template <
class Flu
idState>
144 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
145 averageMolarMass_[phaseIdx] = 0;
146 sumMoleFractions_[phaseIdx] = 0;
147 for (
unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
148 moleFraction_[phaseIdx][compIdx] =
149 decay<Scalar>(fs.moleFraction(phaseIdx, compIdx));
151 averageMolarMass_[phaseIdx] += moleFraction_[phaseIdx][compIdx]*FluidSystem::molarMass(compIdx);
152 sumMoleFractions_[phaseIdx] += moleFraction_[phaseIdx][compIdx];
167 Valgrind::CheckDefined(moleFraction_);
168 Valgrind::CheckDefined(averageMolarMass_);
169 Valgrind::CheckDefined(sumMoleFractions_);
170 Valgrind::CheckDefined(K_);
171 Valgrind::CheckDefined(L_);
174 const Scalar& K(
unsigned compIdx)
const
190 const Scalar&
L()
const
209 const auto& acf = FluidSystem::acentricFactor(compIdx);
210 const auto& T_crit = FluidSystem::criticalTemperature(compIdx);
211 const auto& T = asImp_().temperature(0);
212 const auto& p_crit = FluidSystem::criticalPressure(compIdx);
213 const auto& p = asImp_().pressure(0);
215 const auto tmp = exp(5.37 * (1+acf) * (1-T_crit/T)) * (p_crit/p);
220 const Implementation& asImp_()
const
222 return *
static_cast<const Implementation*
>(
this);
225 std::array<std::array<Scalar,numComponents>,numPhases> moleFraction_;
226 std::array<Scalar,numPhases> averageMolarMass_;
227 std::array<Scalar,numPhases> sumMoleFractions_;
228 std::array<Scalar,numPhases> Z_;
229 std::array<Scalar,numComponents> K_;
237 template <
class Scalar,
239 class Implementation>
242 enum { numPhases = FluidSystem::numPhases };
245 enum { numComponents = FluidSystem::numComponents };
246 static_assert(
static_cast<int>(numPhases) ==
static_cast<int>(numComponents),
247 "The number of phases must be the same as the number of (pseudo-) components if you assume immiscibility");
255 {
return (phaseIdx == compIdx)?1.0:0.0; }
261 {
return (phaseIdx == compIdx)?1.0:0.0; }
272 {
return FluidSystem::molarMass(phaseIdx); }
283 Scalar
molarity(
unsigned phaseIdx,
unsigned compIdx)
const
284 {
return asImp_().molarDensity(phaseIdx)*
moleFraction(phaseIdx, compIdx); }
290 template <
class Flu
idState>
306 const Implementation& asImp_()
const
307 {
return *
static_cast<const Implementation*
>(
this); }
314 template <
class Scalar>
318 enum { numComponents = 0 };
326 {
throw std::logic_error(
"Mole fractions are not provided by this fluid state"); }
332 {
throw std::logic_error(
"Mass fractions are not provided by this fluid state"); }
343 {
throw std::logic_error(
"Mean molar masses are not provided by this fluid state"); }
355 {
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:198
Scalar molarity(unsigned phaseIdx, unsigned compIdx) const
The concentration of a component in a phase [mol/m^3].
Definition: FluidStateCompositionModules.hpp:102
const Scalar & moleFraction(unsigned phaseIdx, unsigned compIdx) const
The mole fraction of a component in a phase [].
Definition: FluidStateCompositionModules.hpp:67
Scalar wilsonK_(unsigned compIdx) const
Wilson formula to calculate K.
Definition: FluidStateCompositionModules.hpp:207
void assign(const FluidState &fs)
Retrieve all parameters from an arbitrary fluid state.
Definition: FluidStateCompositionModules.hpp:142
const Scalar & L() const
The L value of a composition [-].
Definition: FluidStateCompositionModules.hpp:190
const Scalar & averageMolarMass(unsigned phaseIdx) const
The mean molar mass of a fluid phase [kg/mol].
Definition: FluidStateCompositionModules.hpp:90
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:110
void setKvalue(unsigned compIdx, const Scalar &value)
Set the K value of a component [-].
Definition: FluidStateCompositionModules.hpp:182
Scalar massFraction(unsigned phaseIdx, unsigned compIdx) const
The mass fraction of a component in a phase [].
Definition: FluidStateCompositionModules.hpp:73
void checkDefined() const
Make sure that all attributes are defined.
Definition: FluidStateCompositionModules.hpp:165
Module for the modular fluid state which provides the phase compositions assuming immiscibility.
Definition: FluidStateCompositionModules.hpp:241
Scalar massFraction(unsigned phaseIdx, unsigned compIdx) const
The mass fraction of a component in a phase [].
Definition: FluidStateCompositionModules.hpp:260
Scalar molarity(unsigned phaseIdx, unsigned compIdx) const
The concentration of a component in a phase [mol/m^3].
Definition: FluidStateCompositionModules.hpp:283
void checkDefined() const
Make sure that all attributes are defined.
Definition: FluidStateCompositionModules.hpp:302
void assign(const FluidState &)
Retrieve all parameters from an arbitrary fluid state.
Definition: FluidStateCompositionModules.hpp:291
Scalar moleFraction(unsigned phaseIdx, unsigned compIdx) const
The mole fraction of a component in a phase [].
Definition: FluidStateCompositionModules.hpp:254
Scalar averageMolarMass(unsigned phaseIdx) const
The mean molar mass of a fluid phase [kg/mol].
Definition: FluidStateCompositionModules.hpp:271
Module for the modular fluid state which does not store the compositions but throws std::logic_error ...
Definition: FluidStateCompositionModules.hpp:316
Scalar moleFraction(unsigned, unsigned) const
The mole fraction of a component in a phase [].
Definition: FluidStateCompositionModules.hpp:325
Scalar averageMolarMass(unsigned) const
The mean molar mass of a fluid phase [kg/mol].
Definition: FluidStateCompositionModules.hpp:342
Scalar molarity(unsigned, unsigned) const
The concentration of a component in a phase [mol/m^3].
Definition: FluidStateCompositionModules.hpp:354
void checkDefined() const
Make sure that all attributes are defined.
Definition: FluidStateCompositionModules.hpp:365
Scalar massFraction(unsigned, unsigned) const
The mass fraction of a component in a phase [].
Definition: FluidStateCompositionModules.hpp:331