28 #ifndef OPM_FLUID_STATE_COMPOSITION_MODULES_HPP
29 #define OPM_FLUID_STATE_COMPOSITION_MODULES_HPP
31 #include <opm/common/Valgrind.hpp>
34 #include <opm/common/ErrorMacros.hpp>
35 #include <opm/common/Exceptions.hpp>
46 template <
class Scalar,
51 enum { numPhases = FluidSystem::numPhases };
52 enum { numComponents = FluidSystem::numComponents };
57 for (
int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
58 for (
int compIdx = 0; compIdx < numComponents; ++compIdx)
59 moleFraction_[phaseIdx][compIdx] = 0.0;
61 Valgrind::SetDefined(moleFraction_);
62 Valgrind::SetUndefined(averageMolarMass_);
63 Valgrind::SetUndefined(sumMoleFractions_);
69 const Scalar&
moleFraction(
unsigned phaseIdx,
unsigned compIdx)
const
70 {
return moleFraction_[phaseIdx][compIdx]; }
78 Opm::abs(sumMoleFractions_[phaseIdx])
79 *moleFraction_[phaseIdx][compIdx]
80 *FluidSystem::molarMass(compIdx)
81 / Opm::max(1e-40, Opm::abs(averageMolarMass_[phaseIdx]));
93 {
return averageMolarMass_[phaseIdx]; }
104 Scalar
molarity(
unsigned phaseIdx,
unsigned compIdx)
const
105 {
return asImp_().molarDensity(phaseIdx)*
moleFraction(phaseIdx, compIdx); }
114 Valgrind::CheckDefined(value);
115 Valgrind::SetUndefined(sumMoleFractions_[phaseIdx]);
116 Valgrind::SetUndefined(averageMolarMass_[phaseIdx]);
117 Valgrind::SetUndefined(moleFraction_[phaseIdx][compIdx]);
119 moleFraction_[phaseIdx][compIdx] = value;
122 sumMoleFractions_[phaseIdx] = 0.0;
123 averageMolarMass_[phaseIdx] = 0.0;
124 for (
unsigned compJIdx = 0; compJIdx < numComponents; ++compJIdx) {
125 sumMoleFractions_[phaseIdx] += moleFraction_[phaseIdx][compJIdx];
126 averageMolarMass_[phaseIdx] += moleFraction_[phaseIdx][compJIdx]*FluidSystem::molarMass(compJIdx);
134 template <
class Flu
idState>
137 for (
unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
138 averageMolarMass_[phaseIdx] = 0;
139 sumMoleFractions_[phaseIdx] = 0;
140 for (
unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
141 moleFraction_[phaseIdx][compIdx] =
142 Opm::decay<Scalar>(fs.moleFraction(phaseIdx, compIdx));
144 averageMolarMass_[phaseIdx] += moleFraction_[phaseIdx][compIdx]*FluidSystem::molarMass(compIdx);
145 sumMoleFractions_[phaseIdx] += moleFraction_[phaseIdx][compIdx];
160 Valgrind::CheckDefined(moleFraction_);
161 Valgrind::CheckDefined(averageMolarMass_);
162 Valgrind::CheckDefined(sumMoleFractions_);
166 const Implementation& asImp_()
const
167 {
return *
static_cast<const Implementation*
>(
this); }
169 Scalar moleFraction_[numPhases][numComponents];
170 Scalar averageMolarMass_[numPhases];
171 Scalar sumMoleFractions_[numPhases];
178 template <
class Scalar,
180 class Implementation>
183 enum { numPhases = FluidSystem::numPhases };
186 enum { numComponents = FluidSystem::numComponents };
187 static_assert(static_cast<int>(numPhases) == static_cast<int>(numComponents),
188 "The number of phases must be the same as the number of (pseudo-) components if you assume immiscibility");
197 {
return (phaseIdx == compIdx)?1.0:0.0; }
203 {
return (phaseIdx == compIdx)?1.0:0.0; }
214 {
return FluidSystem::molarMass(phaseIdx); }
225 Scalar
molarity(
unsigned phaseIdx,
unsigned compIdx)
const
226 {
return asImp_().molarDensity(phaseIdx)*
moleFraction(phaseIdx, compIdx); }
232 template <
class Flu
idState>
248 const Implementation& asImp_()
const
249 {
return *
static_cast<const Implementation*
>(
this); }
256 template <
class Scalar>
260 enum { numComponents = 0 };
269 { OPM_THROW(std::logic_error,
"Mole fractions are not provided by this fluid state"); }
275 { OPM_THROW(std::logic_error,
"Mass fractions are not provided by this fluid state"); }
286 { OPM_THROW(std::logic_error,
"Mean molar masses are not provided by this fluid state"); }
298 { OPM_THROW(std::logic_error,
"Molarities are not provided by this fluid state"); }
Scalar molarity(unsigned phaseIdx, unsigned compIdx) const
The concentration of a component in a phase [mol/m^3].
Definition: FluidStateCompositionModules.hpp:225
void assign(const FluidState &fs)
Retrieve all parameters from an arbitrary fluid state.
Definition: FluidStateCompositionModules.hpp:135
Module for the modular fluid state which does not store the compositions but throws std::logic_error ...
Definition: FluidStateCompositionModules.hpp:257
Module for the modular fluid state which stores the phase compositions explicitly in terms of mole fr...
Definition: FluidStateCompositionModules.hpp:49
Scalar massFraction(unsigned phaseIdx, unsigned compIdx) const
The mass fraction of a component in a phase [].
Definition: FluidStateCompositionModules.hpp:202
const Scalar & moleFraction(unsigned phaseIdx, unsigned compIdx) const
The mole fraction of a component in a phase [].
Definition: FluidStateCompositionModules.hpp:69
Scalar averageMolarMass(unsigned phaseIdx) const
The mean molar mass of a fluid phase [kg/mol].
Definition: FluidStateCompositionModules.hpp:213
Scalar averageMolarMass(unsigned) const
The mean molar mass of a fluid phase [kg/mol].
Definition: FluidStateCompositionModules.hpp:285
Scalar molarity(unsigned phaseIdx, unsigned compIdx) const
The concentration of a component in a phase [mol/m^3].
Definition: FluidStateCompositionModules.hpp:104
Scalar moleFraction(unsigned, unsigned) const
The mole fraction of a component in a phase [].
Definition: FluidStateCompositionModules.hpp:268
void checkDefined() const
Make sure that all attributes are defined.
Definition: FluidStateCompositionModules.hpp:308
void assign(const FluidState &)
Retrieve all parameters from an arbitrary fluid state.
Definition: FluidStateCompositionModules.hpp:233
Scalar massFraction(unsigned phaseIdx, unsigned compIdx) const
The mass fraction of a component in a phase [].
Definition: FluidStateCompositionModules.hpp:75
Module for the modular fluid state which provides the phase compositions assuming immiscibility...
Definition: FluidStateCompositionModules.hpp:181
void checkDefined() const
Make sure that all attributes are defined.
Definition: FluidStateCompositionModules.hpp:158
Scalar moleFraction(unsigned phaseIdx, unsigned compIdx) const
The mole fraction of a component in a phase [].
Definition: FluidStateCompositionModules.hpp:196
const Scalar & averageMolarMass(unsigned phaseIdx) const
The mean molar mass of a fluid phase [kg/mol].
Definition: FluidStateCompositionModules.hpp:92
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:112
void checkDefined() const
Make sure that all attributes are defined.
Definition: FluidStateCompositionModules.hpp:244
Scalar massFraction(unsigned, unsigned) const
The mass fraction of a component in a phase [].
Definition: FluidStateCompositionModules.hpp:274
Scalar molarity(unsigned, unsigned) const
The concentration of a component in a phase [mol/m^3].
Definition: FluidStateCompositionModules.hpp:297