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_TWO_PHASE_IMMISCIBLE_FLUID_SYSTEM_HPP
00028 #define OPM_TWO_PHASE_IMMISCIBLE_FLUID_SYSTEM_HPP
00029
00030 #include <limits>
00031 #include <cassert>
00032
00033 #include <opm/material/fluidsystems/LiquidPhase.hpp>
00034 #include <opm/material/fluidsystems/GasPhase.hpp>
00035 #include <opm/material/fluidstates/ImmiscibleFluidState.hpp>
00036
00037 #include "BaseFluidSystem.hpp"
00038 #include "NullParameterCache.hpp"
00039
00040 namespace Opm {
00041 namespace FluidSystems {
00042
00056 template <class Scalar, class WettingPhase, class NonwettingPhase>
00057 class TwoPhaseImmiscible
00058 : public BaseFluidSystem<Scalar, TwoPhaseImmiscible<Scalar, WettingPhase, NonwettingPhase> >
00059 {
00060
00061 TwoPhaseImmiscible()
00062 {}
00063
00064 typedef TwoPhaseImmiscible<Scalar, WettingPhase, NonwettingPhase> ThisType;
00065 typedef BaseFluidSystem<Scalar, ThisType> Base;
00066
00067 public:
00068 template <class Evaluation>
00069 struct ParameterCache : public Opm::NullParameterCache<Evaluation>
00070 {};
00071
00072
00073
00074
00075
00077 static const int numPhases = 2;
00078
00080 static const int wettingPhaseIdx = 0;
00082 static const int nonWettingPhaseIdx = 1;
00083
00085 static const char* phaseName(unsigned phaseIdx)
00086 {
00087 assert(0 <= phaseIdx && phaseIdx < numPhases);
00088
00089 static const char* name[] = {
00090 "wetting",
00091 "nonwetting"
00092 };
00093 return name[phaseIdx];
00094 }
00095
00097 static bool isLiquid(unsigned phaseIdx)
00098 {
00099
00100 return
00101 (phaseIdx == wettingPhaseIdx)
00102 ? WettingPhase::isLiquid()
00103 : NonwettingPhase::isLiquid();
00104 }
00105
00107 static bool isCompressible(unsigned phaseIdx)
00108 {
00109
00110
00111 return
00112 (phaseIdx == wettingPhaseIdx)
00113 ? WettingPhase::isCompressible()
00114 : NonwettingPhase::isCompressible();
00115 }
00116
00118 static bool isIdealGas(unsigned phaseIdx)
00119 {
00120
00121
00122
00123 return
00124 (phaseIdx == wettingPhaseIdx)
00125 ? WettingPhase::isIdealGas()
00126 : NonwettingPhase::isIdealGas();
00127 }
00128
00130 static bool isIdealMixture(unsigned )
00131 {
00132
00133
00134
00135 return true;
00136 }
00137
00138
00139
00140
00141
00143 static const int numComponents = 2;
00144
00146 static const int wettingCompIdx = 0;
00148 static const int nonWettingCompIdx = 1;
00149
00151 static const char* componentName(unsigned compIdx)
00152 {
00153 assert(0 <= compIdx && compIdx < numComponents);
00154
00155 if (compIdx == wettingCompIdx)
00156 return WettingPhase::name();
00157 return NonwettingPhase::name();
00158 }
00159
00161 static Scalar molarMass(unsigned compIdx)
00162 {
00163
00164
00165
00166 return
00167 (compIdx == wettingCompIdx)
00168 ? WettingPhase::molarMass()
00169 : NonwettingPhase::molarMass();
00170 }
00171
00175 static Scalar criticalTemperature(unsigned compIdx)
00176 {
00177
00178
00179 return
00180 (compIdx == wettingCompIdx)
00181 ? WettingPhase::criticalTemperature()
00182 : NonwettingPhase::criticalTemperature();
00183 }
00184
00188 static Scalar criticalPressure(unsigned compIdx)
00189 {
00190
00191
00192 return
00193 (compIdx == wettingCompIdx)
00194 ? WettingPhase::criticalPressure()
00195 : NonwettingPhase::criticalPressure();
00196 }
00197
00201 static Scalar acentricFactor(unsigned compIdx)
00202 {
00203
00204
00205 return
00206 (compIdx == wettingCompIdx)
00207 ? WettingPhase::acentricFactor()
00208 : NonwettingPhase::acentricFactor();
00209 }
00210
00211
00212
00213
00214
00216 static void init()
00217 {
00218
00219
00220 assert(WettingPhase::isLiquid() || NonwettingPhase::isLiquid());
00221 }
00222
00224 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
00225 static LhsEval density(const FluidState& fluidState,
00226 const ParameterCache<ParamCacheEval>& ,
00227 unsigned phaseIdx)
00228 {
00229 assert(0 <= phaseIdx && phaseIdx < numPhases);
00230
00231 const auto& temperature = Opm::decay<LhsEval>(fluidState.temperature(phaseIdx));
00232 const auto& pressure = Opm::decay<LhsEval>(fluidState.pressure(phaseIdx));
00233 if (phaseIdx == wettingPhaseIdx)
00234 return WettingPhase::density(temperature, pressure);
00235 return NonwettingPhase::density(temperature, pressure);
00236 }
00237
00239 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
00240 static LhsEval viscosity(const FluidState& fluidState,
00241 const ParameterCache<ParamCacheEval>& ,
00242 unsigned phaseIdx)
00243 {
00244 assert(0 <= phaseIdx && phaseIdx < numPhases);
00245
00246 const auto& temperature = Opm::decay<LhsEval>(fluidState.temperature(phaseIdx));
00247 const auto& pressure = Opm::decay<LhsEval>(fluidState.pressure(phaseIdx));
00248 if (phaseIdx == wettingPhaseIdx)
00249 return WettingPhase::viscosity(temperature, pressure);
00250 return NonwettingPhase::viscosity(temperature, pressure);
00251 }
00252
00254 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
00255 static LhsEval fugacityCoefficient(const FluidState& ,
00256 const ParameterCache<ParamCacheEval>& ,
00257 unsigned phaseIdx,
00258 unsigned compIdx)
00259 {
00260 assert(0 <= phaseIdx && phaseIdx < numPhases);
00261 assert(0 <= compIdx && compIdx < numComponents);
00262
00263 if (phaseIdx == compIdx)
00264
00265
00266
00267
00268 return 1.0;
00269 return std::numeric_limits<Scalar>::infinity();
00270 }
00271
00273 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
00274 static LhsEval enthalpy(const FluidState& fluidState,
00275 const ParameterCache<ParamCacheEval>& ,
00276 unsigned phaseIdx)
00277 {
00278 assert(0 <= phaseIdx && phaseIdx < numPhases);
00279
00280 const auto& temperature = Opm::decay<LhsEval>(fluidState.temperature(phaseIdx));
00281 const auto& pressure = Opm::decay<LhsEval>(fluidState.pressure(phaseIdx));
00282 if (phaseIdx == wettingPhaseIdx)
00283 return WettingPhase::enthalpy(temperature, pressure);
00284 return NonwettingPhase::enthalpy(temperature, pressure);
00285 }
00286
00288 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
00289 static LhsEval thermalConductivity(const FluidState& fluidState,
00290 const ParameterCache<ParamCacheEval>& ,
00291 unsigned phaseIdx)
00292 {
00293 assert(0 <= phaseIdx && phaseIdx < numPhases);
00294
00295 const auto& temperature = Opm::decay<LhsEval>(fluidState.temperature(phaseIdx));
00296 const auto& pressure = Opm::decay<LhsEval>(fluidState.pressure(phaseIdx));
00297 if (phaseIdx == wettingPhaseIdx)
00298 return WettingPhase::thermalConductivity(temperature, pressure);
00299 return NonwettingPhase::thermalConductivity(temperature, pressure);
00300 }
00301
00303 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
00304 static LhsEval heatCapacity(const FluidState& fluidState,
00305 const ParameterCache<ParamCacheEval>& ,
00306 unsigned phaseIdx)
00307 {
00308 assert(0 <= phaseIdx && phaseIdx < numPhases);
00309
00310 const auto& temperature = Opm::decay<LhsEval>(fluidState.temperature(phaseIdx));
00311 const auto& pressure = Opm::decay<LhsEval>(fluidState.pressure(phaseIdx));
00312 if (phaseIdx == wettingPhaseIdx)
00313 return WettingPhase::heatCapacity(temperature, pressure);
00314 return NonwettingPhase::heatCapacity(temperature, pressure);
00315 }
00316 };
00317
00318 }
00319
00320 }
00321
00322 #endif