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);
132 template <
class Flu
idState>
135 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
136 averageMolarMass_[phaseIdx] = 0;
137 sumMoleFractions_[phaseIdx] = 0;
138 for (
unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
139 moleFraction_[phaseIdx][compIdx] =
140 decay<Scalar>(fs.moleFraction(phaseIdx, compIdx));
142 averageMolarMass_[phaseIdx] += moleFraction_[phaseIdx][compIdx]*FluidSystem::molarMass(compIdx);
143 sumMoleFractions_[phaseIdx] += moleFraction_[phaseIdx][compIdx];
158 Valgrind::CheckDefined(moleFraction_);
159 Valgrind::CheckDefined(averageMolarMass_);
160 Valgrind::CheckDefined(sumMoleFractions_);
164 const Implementation& asImp_()
const
165 {
return *
static_cast<const Implementation*
>(
this); }
167 std::array<std::array<Scalar,numComponents>,numPhases> moleFraction_;
168 std::array<Scalar,numPhases> averageMolarMass_;
169 std::array<Scalar,numPhases> sumMoleFractions_;
176 template <
class Scalar,
178 class Implementation>
181 enum { numPhases = FluidSystem::numPhases };
184 enum { numComponents = FluidSystem::numComponents };
185 static_assert(
static_cast<int>(numPhases) ==
static_cast<int>(numComponents),
186 "The number of phases must be the same as the number of (pseudo-) components if you assume immiscibility");
195 {
return (phaseIdx == compIdx)?1.0:0.0; }
201 {
return (phaseIdx == compIdx)?1.0:0.0; }
212 {
return FluidSystem::molarMass(phaseIdx); }
223 Scalar
molarity(
unsigned phaseIdx,
unsigned compIdx)
const
224 {
return asImp_().molarDensity(phaseIdx)*
moleFraction(phaseIdx, compIdx); }
230 template <
class Flu
idState>
246 const Implementation& asImp_()
const
247 {
return *
static_cast<const Implementation*
>(
this); }
254 template <
class Scalar>
258 enum { numComponents = 0 };
267 {
throw std::logic_error(
"Mole fractions are not provided by this fluid state"); }
273 {
throw std::logic_error(
"Mass fractions are not provided by this fluid state"); }
284 {
throw std::logic_error(
"Mean molar masses are not provided by this fluid state"); }
296 {
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
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
void assign(const FluidState &fs)
Retrieve all parameters from an arbitrary fluid state.
Definition: FluidStateCompositionModules.hpp:133
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
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:156
Module for the modular fluid state which provides the phase compositions assuming immiscibility.
Definition: FluidStateCompositionModules.hpp:180
Scalar massFraction(unsigned phaseIdx, unsigned compIdx) const
The mass fraction of a component in a phase [].
Definition: FluidStateCompositionModules.hpp:200
Scalar molarity(unsigned phaseIdx, unsigned compIdx) const
The concentration of a component in a phase [mol/m^3].
Definition: FluidStateCompositionModules.hpp:223
void checkDefined() const
Make sure that all attributes are defined.
Definition: FluidStateCompositionModules.hpp:242
void assign(const FluidState &)
Retrieve all parameters from an arbitrary fluid state.
Definition: FluidStateCompositionModules.hpp:231
Scalar moleFraction(unsigned phaseIdx, unsigned compIdx) const
The mole fraction of a component in a phase [].
Definition: FluidStateCompositionModules.hpp:194
Scalar averageMolarMass(unsigned phaseIdx) const
The mean molar mass of a fluid phase [kg/mol].
Definition: FluidStateCompositionModules.hpp:211
Module for the modular fluid state which does not store the compositions but throws std::logic_error ...
Definition: FluidStateCompositionModules.hpp:256
Scalar moleFraction(unsigned, unsigned) const
The mole fraction of a component in a phase [].
Definition: FluidStateCompositionModules.hpp:266
Scalar averageMolarMass(unsigned) const
The mean molar mass of a fluid phase [kg/mol].
Definition: FluidStateCompositionModules.hpp:283
Scalar molarity(unsigned, unsigned) const
The concentration of a component in a phase [mol/m^3].
Definition: FluidStateCompositionModules.hpp:295
void checkDefined() const
Make sure that all attributes are defined.
Definition: FluidStateCompositionModules.hpp:306
Scalar massFraction(unsigned, unsigned) const
The mass fraction of a component in a phase [].
Definition: FluidStateCompositionModules.hpp:272