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_AIR_HPP
00028 #define OPM_AIR_HPP
00029
00030 #include <opm/material/components/Component.hpp>
00031 #include <opm/material/common/MathToolbox.hpp>
00032 #include <opm/material/IdealGas.hpp>
00033
00034 #include <opm/common/Exceptions.hpp>
00035 #include <opm/common/ErrorMacros.hpp>
00036
00037 namespace Opm {
00038
00046 template <class Scalar>
00047 class Air : public Component<Scalar, Air<Scalar> >
00048 {
00049 typedef Opm::IdealGas<Scalar> IdealGas;
00050
00051 public:
00055 static bool liquidIsCompressible()
00056 { OPM_THROW(std::runtime_error, "Not implemented: Component::liquidIsCompressible()"); }
00057
00061 static const char* name()
00062 { return "Air"; }
00063
00067 static bool gasIsCompressible()
00068 { return true; }
00069
00073 static bool gasIsIdeal()
00074 { return true; }
00075
00081 static Scalar molarMass()
00082 { return 0.02896; }
00083
00087 static Scalar criticalTemperature()
00088 { return 132.531 ; }
00089
00093 static Scalar criticalPressure()
00094 { return 37.86e5; }
00095
00102 template <class Evaluation>
00103 static Evaluation gasDensity(const Evaluation& temperature, const Evaluation& pressure)
00104 { return IdealGas::density(Evaluation(molarMass()), temperature, pressure); }
00105
00112 template <class Evaluation>
00113 static Evaluation gasPressure(const Evaluation& temperature, Scalar density)
00114 { return IdealGas::pressure(temperature, density/molarMass()); }
00115
00137 template <class Evaluation>
00138 static Evaluation gasViscosity(const Evaluation& temperature, const Evaluation& )
00139 {
00140 Scalar Tc = criticalTemperature();
00141 Scalar Vc = 84.525138;
00142 Scalar omega = 0.078;
00143 Scalar M = molarMass() * 1e3;
00144 Scalar dipole = 0.0;
00145
00146 Scalar mu_r4 = 131.3 * dipole / std::sqrt(Vc * Tc);
00147 mu_r4 *= mu_r4;
00148 mu_r4 *= mu_r4;
00149
00150 Scalar Fc = 1 - 0.2756*omega + 0.059035*mu_r4;
00151 Evaluation Tstar = 1.2593 * temperature/Tc;
00152 Evaluation Omega_v =
00153 1.16145*Opm::pow(Tstar, -0.14874) +
00154 0.52487*Opm::exp(- 0.77320*Tstar) +
00155 2.16178*Opm::exp(- 2.43787*Tstar);
00156 return 40.7851e-7*Fc*Opm::sqrt(M*temperature)/(std::pow(Vc, 2./3)*Omega_v);
00157 }
00158
00159
00160 template <class Evaluation>
00161 static Evaluation simpleGasViscosity(const Evaluation& temperature, const Evaluation& )
00162 {
00163 if(temperature < 273.15 || temperature > 660.) {
00164 OPM_THROW(NumericalProblem,
00165 "Air: Temperature (" << temperature << "K) out of range");
00166 }
00167 return 1.496e-6*Opm::pow(temperature, 1.5)/(temperature + 120);
00168 }
00169
00181 template <class Evaluation>
00182 static Evaluation gasEnthalpy(const Evaluation& temperature, const Evaluation& )
00183 {
00184 return 1005*(temperature - 273.15);
00185 }
00186
00198 template <class Evaluation>
00199 static Evaluation gasInternalEnergy(const Evaluation& temperature,
00200 const Evaluation& pressure)
00201 {
00202 return
00203 gasEnthalpy(temperature, pressure)
00204 - (IdealGas::R*temperature/molarMass());
00205 }
00206
00218 template <class Evaluation>
00219 static Evaluation gasThermalConductivity(const Evaluation& ,
00220 const Evaluation& )
00221 {
00222
00223
00224
00225
00226
00227
00228 return 0.0255535;
00229 }
00230
00247 template <class Evaluation>
00248 static Evaluation gasHeatCapacity(const Evaluation& temperature,
00249 const Evaluation& )
00250 {
00251
00252 Evaluation phi = temperature/100;
00253
00254 Evaluation c_p =
00255 0.661738E+01
00256 -0.105885E+01 * phi
00257 +0.201650E+00 * Opm::pow(phi,2.)
00258 -0.196930E-01 * Opm::pow(phi,3.)
00259 +0.106460E-02 * Opm::pow(phi,4.)
00260 -0.303284E-04 * Opm::pow(phi,5.)
00261 +0.355861E-06 * Opm::pow(phi,6.);
00262 c_p +=
00263 -0.549169E+01 * Opm::pow(phi,-1.)
00264 +0.585171E+01* Opm::pow(phi,-2.)
00265 -0.372865E+01* Opm::pow(phi,-3.)
00266 +0.133981E+01* Opm::pow(phi,-4.)
00267 -0.233758E+00* Opm::pow(phi,-5.)
00268 +0.125718E-01* Opm::pow(phi,-6.);
00269 c_p *= IdealGas::R / (molarMass() * 1000);
00270
00271 return c_p;
00272 }
00273 };
00274
00275 }
00276
00277 #endif