00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00027 #ifndef OPM_SINGLE_PHASE_FLUIDSYSTEM_HPP
00028 #define OPM_SINGLE_PHASE_FLUIDSYSTEM_HPP
00029
00030 #include "BaseFluidSystem.hpp"
00031 #include "NullParameterCache.hpp"
00032
00033 #include <opm/material/fluidsystems/LiquidPhase.hpp>
00034 #include <opm/material/fluidsystems/GasPhase.hpp>
00035 #include <opm/material/components/SimpleH2O.hpp>
00036 #include <opm/material/components/H2O.hpp>
00037 #include <opm/material/components/N2.hpp>
00038 #include <opm/material/components/TabulatedComponent.hpp>
00039
00040 #include <opm/common/Unused.hpp>
00041
00042 #include <limits>
00043 #include <cassert>
00044
00045 namespace Opm {
00046 namespace FluidSystems {
00047
00057 template <class Scalar, class Fluid>
00058 class SinglePhase
00059 : public BaseFluidSystem<Scalar, SinglePhase<Scalar, Fluid> >
00060 {
00061 typedef SinglePhase<Scalar, Fluid> ThisType;
00062 typedef BaseFluidSystem<Scalar, ThisType> Base;
00063
00064 public:
00066 template <class Evaluation>
00067 struct ParameterCache : public Opm::NullParameterCache<Evaluation>
00068 {};
00069
00070
00071
00072
00073
00075 static const int numPhases = 1;
00076
00078 static const char* phaseName(unsigned phaseIdx OPM_OPTIM_UNUSED)
00079 {
00080 assert(0 <= phaseIdx && phaseIdx < numPhases);
00081
00082 return Fluid::name();
00083 }
00084
00086 static bool isLiquid(unsigned )
00087 {
00088
00089
00090 return Fluid::isLiquid();
00091 }
00092
00094 static bool isCompressible(unsigned )
00095 {
00096
00097
00098
00099 return Fluid::isCompressible();
00100 }
00101
00103 static bool isIdealMixture(unsigned )
00104 {
00105
00106
00107
00108 return true;
00109 }
00110
00112 static bool isIdealGas(unsigned )
00113 {
00114
00115
00116
00117 return Fluid::isIdealGas();
00118 }
00119
00120
00121
00122
00123
00125 static const int numComponents = 1;
00126
00128 static const char* componentName(unsigned compIdx OPM_OPTIM_UNUSED)
00129 {
00130 assert(0 <= compIdx && compIdx < numComponents);
00131
00132 return Fluid::name();
00133 }
00134
00136 static Scalar molarMass(unsigned )
00137 {
00138
00139
00140 return Fluid::molarMass();
00141 }
00142
00148 static Scalar criticalTemperature(unsigned )
00149 {
00150
00151
00152 return Fluid::criticalTemperature();
00153 }
00154
00160 static Scalar criticalPressure(unsigned )
00161 {
00162
00163
00164 return Fluid::criticalPressure();
00165 }
00166
00172 static Scalar acentricFactor(unsigned )
00173 {
00174
00175
00176 return Fluid::acentricFactor();
00177 }
00178
00179
00180
00181
00182
00184 static void init()
00185 { }
00186
00188 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
00189 static LhsEval density(const FluidState& fluidState,
00190 const ParameterCache<ParamCacheEval>& ,
00191 unsigned phaseIdx)
00192 {
00193 assert(0 <= phaseIdx && phaseIdx < numPhases);
00194
00195 const auto& T = Opm::decay<LhsEval>(fluidState.temperature(phaseIdx));
00196 const auto& p = Opm::decay<LhsEval>(fluidState.pressure(phaseIdx));
00197 return Fluid::density(T, p);
00198 }
00199
00201 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
00202 static LhsEval viscosity(const FluidState& fluidState,
00203 const ParameterCache<ParamCacheEval>& ,
00204 unsigned phaseIdx)
00205 {
00206 assert(0 <= phaseIdx && phaseIdx < numPhases);
00207
00208 const auto& T = Opm::decay<LhsEval>(fluidState.temperature(phaseIdx));
00209 const auto& p = Opm::decay<LhsEval>(fluidState.pressure(phaseIdx));
00210 return Fluid::viscosity(T, p);
00211 }
00212
00214 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
00215 static LhsEval fugacityCoefficient(const FluidState& ,
00216 const ParameterCache<ParamCacheEval>& ,
00217 unsigned phaseIdx,
00218 unsigned compIdx)
00219 {
00220 assert(0 <= phaseIdx && phaseIdx < numPhases);
00221 assert(0 <= compIdx && compIdx < numComponents);
00222
00223 if (phaseIdx == compIdx)
00224
00225
00226
00227
00228 return 1.0;
00229 return std::numeric_limits<Scalar>::infinity();
00230 }
00231
00233 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
00234 static LhsEval enthalpy(const FluidState& fluidState,
00235 const ParameterCache<ParamCacheEval>& ,
00236 unsigned phaseIdx)
00237 {
00238 assert(0 <= phaseIdx && phaseIdx < numPhases);
00239
00240 const auto& T = Opm::decay<LhsEval>(fluidState.temperature(phaseIdx));
00241 const auto& p = Opm::decay<LhsEval>(fluidState.pressure(phaseIdx));
00242 return Fluid::enthalpy(T, p);
00243 }
00244
00246 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
00247 static LhsEval thermalConductivity(const FluidState& fluidState,
00248 const ParameterCache<ParamCacheEval>& ,
00249 unsigned phaseIdx)
00250 {
00251 assert(0 <= phaseIdx && phaseIdx < numPhases);
00252
00253 const auto& T = Opm::decay<LhsEval>(fluidState.temperature(phaseIdx));
00254 const auto& p = Opm::decay<LhsEval>(fluidState.pressure(phaseIdx));
00255 return Fluid::thermalConductivity(T, p);
00256 }
00257
00259 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
00260 static LhsEval heatCapacity(const FluidState& fluidState,
00261 const ParameterCache<ParamCacheEval>& ,
00262 unsigned phaseIdx)
00263 {
00264 assert(0 <= phaseIdx && phaseIdx < numPhases);
00265
00266 const auto& T = Opm::decay<LhsEval>(fluidState.temperature(phaseIdx));
00267 const auto& p = Opm::decay<LhsEval>(fluidState.pressure(phaseIdx));
00268 return Fluid::heatCapacity(T, p);
00269 }
00270 };
00271
00272 }}
00273
00274 #endif