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_FUGACITY_MODULES_HPP
00029 #define OPM_FLUID_STATE_FUGACITY_MODULES_HPP
00030
00031 #include <opm/common/Valgrind.hpp>
00032
00033 #include <opm/common/ErrorMacros.hpp>
00034 #include <opm/common/Exceptions.hpp>
00035
00036 #include <algorithm>
00037 #include <limits>
00038
00039 namespace Opm {
00040
00045 template <class Scalar,
00046 unsigned numPhases,
00047 unsigned numComponents,
00048 class Implementation>
00049 class FluidStateExplicitFugacityModule
00050 {
00051 public:
00052 FluidStateExplicitFugacityModule()
00053 {
00054 Valgrind::SetUndefined(fugacityCoefficient_);
00055 }
00056
00060 const Scalar& fugacityCoefficient(unsigned phaseIdx, unsigned compIdx) const
00061 { return fugacityCoefficient_[phaseIdx][compIdx]; }
00062
00066 Scalar fugacity(unsigned phaseIdx, unsigned compIdx) const
00067 { return asImp_().pressure(phaseIdx)*fugacityCoefficient_[phaseIdx][compIdx]*asImp_().moleFraction(phaseIdx, compIdx); }
00068
00072 void setFugacityCoefficient(unsigned phaseIdx, unsigned compIdx, const Scalar& value)
00073 { fugacityCoefficient_[phaseIdx][compIdx] = value; }
00074
00079 template <class FluidState>
00080 void assign(const FluidState& fs)
00081 {
00082 for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
00083 for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx) {
00084 fugacityCoefficient_[phaseIdx][compIdx] = fs.fugacityCoefficient(phaseIdx, compIdx);
00085 }
00086 }
00087 }
00088
00097 void checkDefined() const
00098 {
00099 Valgrind::CheckDefined(fugacityCoefficient_);
00100 }
00101
00102 protected:
00103 const Implementation& asImp_() const
00104 { return *static_cast<const Implementation*>(this); }
00105
00106 Scalar fugacityCoefficient_[numPhases][numComponents];
00107 };
00108
00113 template <class Scalar,
00114 unsigned numPhases,
00115 unsigned numComponents,
00116 class Implementation>
00117 class FluidStateImmiscibleFugacityModule
00118 {
00119 static_assert(numPhases == numComponents,
00120 "The number of phases must be the same as the number of (pseudo-) components if you assume immiscibility");
00121
00122 public:
00123 FluidStateImmiscibleFugacityModule()
00124 { Valgrind::SetUndefined(fugacityCoefficient_); }
00125
00129 Scalar fugacityCoefficient(unsigned phaseIdx, unsigned compIdx) const
00130 { return (phaseIdx == compIdx)?fugacityCoefficient_[phaseIdx]:std::numeric_limits<Scalar>::infinity(); }
00131
00135 Scalar fugacity(unsigned phaseIdx, unsigned compIdx) const
00136 { return asImp_().pressure(phaseIdx)*fugacityCoefficient(phaseIdx, compIdx)*asImp_().moleFraction(phaseIdx, compIdx); }
00137
00141 void setFugacityCoefficient(unsigned phaseIdx, const Scalar& value)
00142 { fugacityCoefficient_[phaseIdx] = value; }
00143
00148 template <class FluidState>
00149 void assign(const FluidState& fs)
00150 {
00151 for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx) {
00152 fugacityCoefficient_[phaseIdx] = fs.fugacityCoefficient(phaseIdx, phaseIdx);
00153 }
00154 }
00155
00164 void checkDefined() const
00165 {
00166 Valgrind::CheckDefined(fugacityCoefficient_);
00167 }
00168
00169 protected:
00170 const Implementation& asImp_() const
00171 { return *static_cast<const Implementation*>(this); }
00172
00173 Scalar fugacityCoefficient_[numPhases];
00174 };
00175
00180 template <class Scalar>
00181 class FluidStateNullFugacityModule
00182 {
00183 public:
00184 FluidStateNullFugacityModule()
00185 { }
00186
00190 const Scalar& fugacityCoefficient(unsigned , unsigned ) const
00191 { OPM_THROW(std::logic_error, "Fugacity coefficients are not provided by this fluid state"); }
00192
00196 const Scalar& fugacity(unsigned , unsigned ) const
00197 { OPM_THROW(std::logic_error, "Fugacities coefficients are not provided by this fluid state"); }
00198
00207 void checkDefined() const
00208 { }
00209 };
00210
00211
00212 }
00213
00214 #endif