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_H2O_N2_LIQUIDPHASE_FLUID_SYSTEM_HPP
00028 #define OPM_H2O_N2_LIQUIDPHASE_FLUID_SYSTEM_HPP
00029
00030 #include "BaseFluidSystem.hpp"
00031 #include "NullParameterCache.hpp"
00032
00033 #include <opm/material/IdealGas.hpp>
00034 #include <opm/material/components/N2.hpp>
00035 #include <opm/material/components/H2O.hpp>
00036 #include <opm/material/components/SimpleH2O.hpp>
00037 #include <opm/material/components/TabulatedComponent.hpp>
00038 #include <opm/material/binarycoefficients/H2O_N2.hpp>
00039 #include <opm/common/Valgrind.hpp>
00040
00041 #include <opm/common/Exceptions.hpp>
00042 #include <opm/common/ErrorMacros.hpp>
00043
00044 #include <iostream>
00045 #include <cassert>
00046
00047 namespace Opm {
00048 namespace FluidSystems {
00049
00056 template <class Scalar, bool useComplexRelations = true>
00057 class H2ON2LiquidPhase
00058 : public BaseFluidSystem<Scalar, H2ON2LiquidPhase<Scalar, useComplexRelations> >
00059 {
00060 typedef H2ON2LiquidPhase<Scalar, useComplexRelations> ThisType;
00061 typedef BaseFluidSystem<Scalar, ThisType> Base;
00062
00063
00064 typedef Opm::H2O<Scalar> IapwsH2O;
00065 typedef Opm::TabulatedComponent<Scalar, IapwsH2O > TabulatedH2O;
00066 typedef Opm::N2<Scalar> SimpleN2;
00067
00068 public:
00070 template <class Evaluation>
00071 struct ParameterCache : public Opm::NullParameterCache<Evaluation>
00072 {};
00073
00074
00075
00076
00077
00079 static const int numPhases = 1;
00080
00082 static const int liquidPhaseIdx = 0;
00083
00085 static const char* phaseName(unsigned phaseIdx OPM_OPTIM_UNUSED)
00086 {
00087 assert(phaseIdx == liquidPhaseIdx);
00088
00089 return "liquid";
00090 }
00091
00093 static bool isLiquid(unsigned )
00094 {
00095
00096 return true;
00097 }
00098
00100 static bool isCompressible(unsigned )
00101 {
00102
00103
00104 return H2O::liquidIsCompressible();
00105 }
00106
00108 static bool isIdealGas(unsigned )
00109 {
00110
00111 return false;
00112 }
00113
00115 static bool isIdealMixture(unsigned )
00116 {
00117
00118
00119
00120
00121 return true;
00122 }
00123
00124
00125
00126
00127
00129 static const int numComponents = 2;
00130
00132 static const int H2OIdx = 0;
00134 static const int N2Idx = 1;
00135
00137 typedef TabulatedH2O H2O;
00138
00139
00140
00142 typedef SimpleN2 N2;
00143
00145 static const char* componentName(unsigned compIdx)
00146 {
00147 static const char* name[] = {
00148 H2O::name(),
00149 N2::name()
00150 };
00151
00152 assert(0 <= compIdx && compIdx < numComponents);
00153 return name[compIdx];
00154 }
00155
00157 static Scalar molarMass(unsigned compIdx)
00158 {
00159
00160 return (compIdx == H2OIdx)
00161 ? H2O::molarMass()
00162 : (compIdx == N2Idx)
00163 ? N2::molarMass()
00164 : 1e30;
00165 }
00166
00172 static Scalar criticalTemperature(unsigned compIdx)
00173 {
00174
00175 return (compIdx == H2OIdx)
00176 ? H2O::criticalTemperature()
00177 : (compIdx == N2Idx)
00178 ? N2::criticalTemperature()
00179 : 1e30;
00180 }
00181
00187 static Scalar criticalPressure(unsigned compIdx)
00188 {
00189
00190 return (compIdx == H2OIdx)
00191 ? H2O::criticalPressure()
00192 : (compIdx == N2Idx)
00193 ? N2::criticalPressure()
00194 : 1e30;
00195 }
00196
00202 static Scalar acentricFactor(unsigned compIdx)
00203 {
00204
00205 return (compIdx == H2OIdx)
00206 ? H2O::acentricFactor()
00207 : (compIdx == N2Idx)
00208 ? N2::acentricFactor()
00209 : 1e30;
00210 }
00211
00212
00213
00214
00215
00222 static void init()
00223 {
00224 init(273.15,
00225 623.15,
00226 50,
00227 0.0,
00228 20e6,
00229 50);
00230 }
00231
00243 static void init(Scalar tempMin, Scalar tempMax, unsigned nTemp,
00244 Scalar pressMin, Scalar pressMax, unsigned nPress)
00245 {
00246 if (H2O::isTabulated) {
00247 TabulatedH2O::init(tempMin, tempMax, nTemp,
00248 pressMin, pressMax, nPress);
00249 }
00250 }
00251
00253 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
00254 static LhsEval density(const FluidState& fluidState,
00255 const ParameterCache<ParamCacheEval>& ,
00256 unsigned phaseIdx)
00257 {
00258 assert(0 <= phaseIdx && phaseIdx < numPhases);
00259
00260 const auto& T = Opm::decay<LhsEval>(fluidState.temperature(phaseIdx));
00261 const auto& p = Opm::decay<LhsEval>(fluidState.pressure(phaseIdx));
00262
00263 LhsEval sumMoleFrac = 0;
00264 for (unsigned compIdx = 0; compIdx < numComponents; ++compIdx)
00265 sumMoleFrac += Opm::decay<LhsEval>(fluidState.moleFraction(phaseIdx, compIdx));
00266
00267 assert(phaseIdx == liquidPhaseIdx);
00268
00269 if (!useComplexRelations)
00270
00271 return H2O::liquidDensity(T, p);
00272 else
00273 {
00274
00275 const LhsEval& rholH2O = H2O::liquidDensity(T, p);
00276 const LhsEval& clH2O = rholH2O/H2O::molarMass();
00277
00278 const auto& xlH2O = Opm::decay<LhsEval>(fluidState.moleFraction(liquidPhaseIdx, H2OIdx));
00279 const auto& xlN2 = Opm::decay<LhsEval>(fluidState.moleFraction(liquidPhaseIdx, N2Idx));
00280
00281
00282
00283 return clH2O*(H2O::molarMass()*xlH2O + N2::molarMass()*xlN2)/sumMoleFrac;
00284 }
00285 }
00286
00288 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
00289 static LhsEval viscosity(const FluidState& fluidState,
00290 const ParameterCache<ParamCacheEval>& ,
00291 unsigned phaseIdx)
00292 {
00293 assert(phaseIdx == liquidPhaseIdx);
00294
00295 const auto& T = Opm::decay<LhsEval>(fluidState.temperature(phaseIdx));
00296 const auto& p = Opm::decay<LhsEval>(fluidState.pressure(phaseIdx));
00297
00298
00299 return H2O::liquidViscosity(T, p);
00300 }
00301
00303 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
00304 static LhsEval fugacityCoefficient(const FluidState& fluidState,
00305 const ParameterCache<ParamCacheEval>& ,
00306 unsigned phaseIdx,
00307 unsigned compIdx)
00308 {
00309 assert(phaseIdx == liquidPhaseIdx);
00310 assert(0 <= compIdx && compIdx < numComponents);
00311
00312 const auto& T = Opm::decay<LhsEval>(fluidState.temperature(phaseIdx));
00313 const auto& p = Opm::decay<LhsEval>(fluidState.pressure(phaseIdx));
00314
00315 if (compIdx == H2OIdx)
00316 return H2O::vaporPressure(T)/p;
00317 return Opm::BinaryCoeff::H2O_N2::henry(T)/p;
00318 }
00319
00321 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
00322 static LhsEval diffusionCoefficient(const FluidState& fluidState,
00323 const ParameterCache<ParamCacheEval>& ,
00324 unsigned phaseIdx,
00325 unsigned )
00326
00327 {
00328 assert(phaseIdx == liquidPhaseIdx);
00329
00330 const auto& T = Opm::decay<LhsEval>(fluidState.temperature(phaseIdx));
00331 const auto& p = Opm::decay<LhsEval>(fluidState.pressure(phaseIdx));
00332
00333 return BinaryCoeff::H2O_N2::liquidDiffCoeff(T, p);
00334 }
00335
00337 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
00338 static LhsEval enthalpy(const FluidState& fluidState,
00339 const ParameterCache<ParamCacheEval>& ,
00340 unsigned phaseIdx)
00341 {
00342 assert (phaseIdx == liquidPhaseIdx);
00343
00344 const auto& T = Opm::decay<LhsEval>(fluidState.temperature(phaseIdx));
00345 const auto& p = Opm::decay<LhsEval>(fluidState.pressure(phaseIdx));
00346 Valgrind::CheckDefined(T);
00347 Valgrind::CheckDefined(p);
00348
00349
00350 return H2O::liquidEnthalpy(T, p);
00351 }
00352
00354 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
00355 static LhsEval thermalConductivity(const FluidState& fluidState,
00356 const ParameterCache<ParamCacheEval>& ,
00357 const unsigned phaseIdx)
00358 {
00359 assert(phaseIdx == liquidPhaseIdx);
00360
00361 if(useComplexRelations){
00362 const auto& T = Opm::decay<LhsEval>(fluidState.temperature(phaseIdx));
00363 const auto& p = Opm::decay<LhsEval>(fluidState.pressure(phaseIdx));
00364 return H2O::liquidThermalConductivity(T, p);
00365 }
00366 else
00367 return 0.578078;
00368 }
00369
00371 template <class FluidState, class LhsEval = typename FluidState::Scalar, class ParamCacheEval = LhsEval>
00372 static LhsEval heatCapacity(const FluidState& fluidState,
00373 const ParameterCache<ParamCacheEval>& ,
00374 unsigned phaseIdx)
00375 {
00376 assert (phaseIdx == liquidPhaseIdx);
00377
00378 const auto& T = Opm::decay<LhsEval>(fluidState.temperature(phaseIdx));
00379 const auto& p = Opm::decay<LhsEval>(fluidState.pressure(phaseIdx));
00380
00381 return H2O::liquidHeatCapacity(T, p);
00382 }
00383 };
00384
00385 }
00386
00387 }
00388
00389 #endif