00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00028 #ifndef OPM_FLUID_STATE_COMPOSITION_MODULES_HPP
00029 #define OPM_FLUID_STATE_COMPOSITION_MODULES_HPP
00030
00031 #include <opm/common/Valgrind.hpp>
00032 #include <opm/material/common/MathToolbox.hpp>
00033
00034 #include <opm/common/ErrorMacros.hpp>
00035 #include <opm/common/Exceptions.hpp>
00036
00037 #include <algorithm>
00038 #include <cmath>
00039
00040 namespace Opm {
00041
00046 template <class Scalar,
00047 class FluidSystem,
00048 class Implementation>
00049 class FluidStateExplicitCompositionModule
00050 {
00051 enum { numPhases = FluidSystem::numPhases };
00052 enum { numComponents = FluidSystem::numComponents };
00053
00054 public:
00055 FluidStateExplicitCompositionModule()
00056 {
00057 for (int phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
00058 for (int compIdx = 0; compIdx < numComponents; ++compIdx)
00059 moleFraction_[phaseIdx][compIdx] = 0.0;
00060
00061 Valgrind::SetDefined(moleFraction_);
00062 Valgrind::SetUndefined(averageMolarMass_);
00063 Valgrind::SetUndefined(sumMoleFractions_);
00064 }
00065
00069 const Scalar& moleFraction(unsigned phaseIdx, unsigned compIdx) const
00070 { return moleFraction_[phaseIdx][compIdx]; }
00071
00075 Scalar massFraction(unsigned phaseIdx, unsigned compIdx) const
00076 {
00077 return
00078 Opm::abs(sumMoleFractions_[phaseIdx])
00079 *moleFraction_[phaseIdx][compIdx]
00080 *FluidSystem::molarMass(compIdx)
00081 / Opm::max(1e-40, Opm::abs(averageMolarMass_[phaseIdx]));
00082 }
00083
00092 const Scalar& averageMolarMass(unsigned phaseIdx) const
00093 { return averageMolarMass_[phaseIdx]; }
00094
00104 Scalar molarity(unsigned phaseIdx, unsigned compIdx) const
00105 { return asImp_().molarDensity(phaseIdx)*moleFraction(phaseIdx, compIdx); }
00106
00112 void setMoleFraction(unsigned phaseIdx, unsigned compIdx, const Scalar& value)
00113 {
00114 Valgrind::CheckDefined(value);
00115 Valgrind::SetUndefined(sumMoleFractions_[phaseIdx]);
00116 Valgrind::SetUndefined(averageMolarMass_[phaseIdx]);
00117 Valgrind::SetUndefined(moleFraction_[phaseIdx][compIdx]);
00118
00119 moleFraction_[phaseIdx][compIdx] = value;
00120
00121
00122 sumMoleFractions_[phaseIdx] = 0.0;
00123 averageMolarMass_[phaseIdx] = 0.0;
00124 for (unsigned compJIdx = 0; compJIdx < numComponents; ++compJIdx) {
00125 sumMoleFractions_[phaseIdx] += moleFraction_[phaseIdx][compJIdx];
00126 averageMolarMass_[phaseIdx] += moleFraction_[phaseIdx][compJIdx]*FluidSystem::molarMass(compJIdx);
00127 }
00128 }
00129
00134 template <class FluidState>
00135 void assign(const FluidState& fs)
00136 {
00137 for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
00138 averageMolarMass_[phaseIdx] = 0;
00139 sumMoleFractions_[phaseIdx] = 0;
00140 for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
00141 moleFraction_[phaseIdx][compIdx] =
00142 Opm::decay<Scalar>(fs.moleFraction(phaseIdx, compIdx));
00143
00144 averageMolarMass_[phaseIdx] += moleFraction_[phaseIdx][compIdx]*FluidSystem::molarMass(compIdx);
00145 sumMoleFractions_[phaseIdx] += moleFraction_[phaseIdx][compIdx];
00146 }
00147 }
00148 }
00149
00158 void checkDefined() const
00159 {
00160 Valgrind::CheckDefined(moleFraction_);
00161 Valgrind::CheckDefined(averageMolarMass_);
00162 Valgrind::CheckDefined(sumMoleFractions_);
00163 }
00164
00165 protected:
00166 const Implementation& asImp_() const
00167 { return *static_cast<const Implementation*>(this); }
00168
00169 Scalar moleFraction_[numPhases][numComponents];
00170 Scalar averageMolarMass_[numPhases];
00171 Scalar sumMoleFractions_[numPhases];
00172 };
00173
00178 template <class Scalar,
00179 class FluidSystem,
00180 class Implementation>
00181 class FluidStateImmiscibleCompositionModule
00182 {
00183 enum { numPhases = FluidSystem::numPhases };
00184
00185 public:
00186 enum { numComponents = FluidSystem::numComponents };
00187 static_assert(static_cast<int>(numPhases) == static_cast<int>(numComponents),
00188 "The number of phases must be the same as the number of (pseudo-) components if you assume immiscibility");
00189
00190 FluidStateImmiscibleCompositionModule()
00191 { }
00192
00196 Scalar moleFraction(unsigned phaseIdx, unsigned compIdx) const
00197 { return (phaseIdx == compIdx)?1.0:0.0; }
00198
00202 Scalar massFraction(unsigned phaseIdx, unsigned compIdx) const
00203 { return (phaseIdx == compIdx)?1.0:0.0; }
00204
00213 Scalar averageMolarMass(unsigned phaseIdx) const
00214 { return FluidSystem::molarMass(phaseIdx); }
00215
00225 Scalar molarity(unsigned phaseIdx, unsigned compIdx) const
00226 { return asImp_().molarDensity(phaseIdx)*moleFraction(phaseIdx, compIdx); }
00227
00232 template <class FluidState>
00233 void assign(const FluidState& )
00234 { }
00235
00244 void checkDefined() const
00245 { }
00246
00247 protected:
00248 const Implementation& asImp_() const
00249 { return *static_cast<const Implementation*>(this); }
00250 };
00251
00256 template <class Scalar>
00257 class FluidStateNullCompositionModule
00258 {
00259 public:
00260 enum { numComponents = 0 };
00261
00262 FluidStateNullCompositionModule()
00263 { }
00264
00268 Scalar moleFraction(unsigned , unsigned ) const
00269 { OPM_THROW(std::logic_error, "Mole fractions are not provided by this fluid state"); }
00270
00274 Scalar massFraction(unsigned , unsigned ) const
00275 { OPM_THROW(std::logic_error, "Mass fractions are not provided by this fluid state"); }
00276
00285 Scalar averageMolarMass(unsigned ) const
00286 { OPM_THROW(std::logic_error, "Mean molar masses are not provided by this fluid state"); }
00287
00297 Scalar molarity(unsigned , unsigned ) const
00298 { OPM_THROW(std::logic_error, "Molarities are not provided by this fluid state"); }
00299
00308 void checkDefined() const
00309 { }
00310 };
00311
00312
00313 }
00314
00315 #endif