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_N2_HPP
00028 #define OPM_N2_HPP
00029
00030 #include <opm/material/IdealGas.hpp>
00031 #include <opm/material/common/MathToolbox.hpp>
00032
00033 #include "Component.hpp"
00034
00035 #include <cmath>
00036
00037 namespace Opm
00038 {
00039
00047 template <class Scalar>
00048 class N2 : public Component<Scalar, N2<Scalar> >
00049 {
00050 typedef Opm::IdealGas<Scalar> IdealGas;
00051
00052 public:
00056 static const char* name()
00057 { return "N2"; }
00058
00062 static Scalar molarMass()
00063 { return 28.0134e-3;}
00064
00068 static Scalar criticalTemperature()
00069 { return 126.192; }
00070
00074 static Scalar criticalPressure()
00075 { return 3.39858e6; }
00076
00080 static Scalar tripleTemperature()
00081 { return 63.151; }
00082
00086 static Scalar triplePressure()
00087 { return 12.523e3; }
00088
00103 template <class Evaluation>
00104 static Evaluation vaporPressure(const Evaluation& temperature)
00105 {
00106 if (temperature > criticalTemperature())
00107 return criticalPressure();
00108 if (temperature < tripleTemperature())
00109 return 0;
00110
00111
00112
00113 const Evaluation& sigma = 1.0 - temperature/criticalTemperature();
00114 const Evaluation& sqrtSigma = Opm::sqrt(sigma);
00115 const Scalar N1 = -6.12445284;
00116 const Scalar N2 = 1.26327220;
00117 const Scalar N3 = -0.765910082;
00118 const Scalar N4 = -1.77570564;
00119 return
00120 criticalPressure() *
00121 Opm::exp(criticalTemperature()/temperature*
00122 (sigma*(N1 +
00123 sqrtSigma*N2 +
00124 sigma*(sqrtSigma*N3 +
00125 sigma*sigma*sigma*N4))));
00126 }
00127
00134 template <class Evaluation>
00135 static Evaluation gasDensity(const Evaluation& temperature, const Evaluation& pressure)
00136 {
00137
00138 return IdealGas::density(Evaluation(molarMass()), temperature, pressure);
00139 }
00140
00144 static bool gasIsCompressible()
00145 { return true; }
00146
00150 static bool gasIsIdeal()
00151 { return true; }
00152
00159 template <class Evaluation>
00160 static Evaluation gasPressure(const Evaluation& temperature, const Evaluation& density)
00161 {
00162
00163 return IdealGas::pressure(temperature, density/molarMass());
00164 }
00165
00175 template <class Evaluation>
00176 static Evaluation gasEnthalpy(const Evaluation& temperature,
00177 const Evaluation& )
00178 {
00179
00180 const Scalar cpVapA = 31.15;
00181 const Scalar cpVapB = -0.01357;
00182 const Scalar cpVapC = 2.680e-5;
00183 const Scalar cpVapD = -1.168e-8;
00184
00185
00186 return
00187 1/molarMass()*
00188
00189 temperature*(cpVapA + temperature*
00190 (cpVapB/2 + temperature*
00191 (cpVapC/3 + temperature*
00192 (cpVapD/4))));
00193
00194
00195
00196
00197
00198
00199
00200
00201 }
00202
00216 template <class Evaluation>
00217 static Evaluation gasInternalEnergy(const Evaluation& temperature,
00218 const Evaluation& pressure)
00219 {
00220 return
00221 gasEnthalpy(temperature, pressure) -
00222 1/molarMass()*
00223 IdealGas::R*temperature;
00224 }
00225
00233 template <class Evaluation>
00234 static Evaluation gasHeatCapacity(const Evaluation& temperature,
00235 const Evaluation& )
00236 {
00237
00238 const Scalar cpVapA = 31.15;
00239 const Scalar cpVapB = -0.01357;
00240 const Scalar cpVapC = 2.680e-5;
00241 const Scalar cpVapD = -1.168e-8;
00242
00243 return
00244 1/molarMass()*
00245
00246 cpVapA + temperature*
00247 (cpVapB + temperature*
00248 (cpVapC + temperature*
00249 (cpVapD)));
00250 }
00251
00265 template <class Evaluation>
00266 static Evaluation gasViscosity(const Evaluation& temperature, const Evaluation& )
00267 {
00268 const Scalar Tc = criticalTemperature();
00269 const Scalar Vc = 90.1;
00270 const Scalar omega = 0.037;
00271 const Scalar M = molarMass() * 1e3;
00272 const Scalar dipole = 0.0;
00273
00274 Scalar mu_r4 = 131.3 * dipole / std::sqrt(Vc * Tc);
00275 mu_r4 *= mu_r4;
00276 mu_r4 *= mu_r4;
00277
00278 Scalar Fc = 1 - 0.2756*omega + 0.059035*mu_r4;
00279 const Evaluation& Tstar = 1.2593 * temperature/Tc;
00280 const Evaluation& Omega_v =
00281 1.16145*Opm::pow(Tstar, -0.14874) +
00282 0.52487*Opm::exp(- 0.77320*Tstar) +
00283 2.16178*Opm::exp(- 2.43787*Tstar);
00284 const Evaluation& mu = 40.785*Fc*Opm::sqrt(M*temperature)/(std::pow(Vc, 2./3)*Omega_v);
00285
00286
00287 return mu/1e6 / 10;
00288 }
00289
00301 template <class Evaluation>
00302 static Evaluation gasThermalConductivity(const Evaluation& ,
00303 const Evaluation& )
00304 { return 0.024572; }
00305 };
00306
00307 }
00308
00309 #endif