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_XYLENE_HPP
00028 #define OPM_XYLENE_HPP
00029
00030 #include <opm/material/IdealGas.hpp>
00031 #include <opm/material/components/Component.hpp>
00032 #include <opm/material/Constants.hpp>
00033
00034 #include <opm/material/common/MathToolbox.hpp>
00035
00036 #include <cmath>
00037
00038 namespace Opm {
00045 template <class Scalar>
00046 class Xylene : public Component<Scalar, Xylene<Scalar> >
00047 {
00048 typedef Constants<Scalar> Consts;
00049 typedef Opm::IdealGas<Scalar> IdealGas;
00050
00051 public:
00055 static const char* name()
00056 { return "xylene"; }
00057
00061 static Scalar molarMass()
00062 { return 0.106; }
00063
00067 static Scalar criticalTemperature()
00068 { return 617.1; }
00069
00073 static Scalar criticalPressure()
00074 { return 35.4e5; }
00075
00079 static Scalar tripleTemperature()
00080 { OPM_THROW(std::runtime_error, "Not implemented: tripleTemperature for xylene"); }
00081
00085 static Scalar triplePressure()
00086 { OPM_THROW(std::runtime_error, "Not implemented: triplePressure for xylene"); }
00087
00094 template <class Evaluation>
00095 static Evaluation vaporPressure(const Evaluation& temperature)
00096 {
00097 const Scalar A = 7.00909;;
00098 const Scalar B = 1462.266;;
00099 const Scalar C = 215.110;;
00100
00101 return 100*1.334*Opm::pow(10.0, (A - (B/(temperature - 273.15 + C))));
00102 }
00103
00112 template <class Evaluation>
00113 static Evaluation spHeatCapLiquidPhase(const Evaluation& temperature, const Evaluation& )
00114 {
00115 Evaluation CH3,C6H5,H;
00116
00117
00118
00119
00120 if(temperature < 298.0){
00121 H = 13.4 + 1.2*(temperature - 273.0)/25.0;
00122 CH3 = 40.0 + 1.6*(temperature - 273.0)/25.0;
00123 C6H5 = 113.0 + 4.2*(temperature - 273.0)/25.0;
00124 }
00125 else if(temperature < 323.0){
00126 H = 14.6 + 0.9*(temperature - 298.0)/25.0;
00127 CH3 = 41.6 + 1.9*(temperature - 298.0)/25.0;
00128 C6H5 = 117.2 + 6.2*(temperature - 298.0)/25.0;
00129 }
00130 else if(temperature < 348.0){
00131 H = 15.5 + 1.2*(temperature - 323.0)/25.0;
00132 CH3 = 43.5 + 2.3*(temperature - 323.0)/25.0;
00133 C6H5 = 123.4 + 6.3*(temperature - 323.0)/25.0;
00134 }
00135 else {
00136 H = 16.7 + 2.1*(temperature - 348.0)/25.0;
00137 CH3 = 45.8 + 2.5*(temperature - 348.0)/25.0;
00138 C6H5 = 129.7 + 6.3*(temperature - 348.0)/25.0;
00139 }
00140
00141 return (C6H5 + 2*CH3 - H)/molarMass();
00142 }
00143
00144
00148 template <class Evaluation>
00149 static Evaluation liquidEnthalpy(const Evaluation& temperature, const Evaluation& pressure)
00150 {
00151
00152
00153
00154
00155
00156
00157
00158
00159
00160
00161
00162 return 0.5*temperature*(spHeatCapLiquidPhase(Evaluation(0.2113*temperature),pressure)
00163 + spHeatCapLiquidPhase(Evaluation(0.7887*temperature),pressure));
00164 }
00165
00169 static Scalar boilingTemperature()
00170 { return 412.3; }
00171
00180 template <class Evaluation>
00181 static Evaluation heatVap(Evaluation temperature,
00182 const Evaluation& )
00183 {
00184 temperature = Opm::min(temperature, criticalTemperature());
00185 temperature = Opm::max(temperature, 0.0);
00186
00187 const Scalar T_crit = criticalTemperature();
00188 const Scalar Tr1 = boilingTemperature()/criticalTemperature();
00189 const Scalar p_crit = criticalPressure();
00190
00191
00192 const Scalar DH_v_boil = Consts::R * T_crit * Tr1
00193 * (3.978 * Tr1 - 3.958 + 1.555*std::log(p_crit * 1e-5 ) )
00194 / (1.07 - Tr1);
00195
00196
00197 const Evaluation& Tr2 = temperature/criticalTemperature();
00198 const Scalar n = 0.375;
00199 const Evaluation& DH_vap = DH_v_boil * Opm::pow(((1.0 - Tr2)/(1.0 - Tr1)), n);
00200
00201 return (DH_vap/molarMass());
00202 }
00203
00210 template <class Evaluation>
00211 static Evaluation gasEnthalpy(const Evaluation& temperature, const Evaluation& pressure)
00212 {
00213 return liquidEnthalpy(temperature, pressure) + heatVap(temperature, pressure);
00214 }
00215
00219 template <class Evaluation>
00220 static Evaluation gasDensity(const Evaluation& temperature, const Evaluation& pressure)
00221 {
00222 return IdealGas::density(Evaluation(molarMass()), temperature, pressure);
00223 }
00224
00231 template <class Evaluation>
00232 static Evaluation molarGasDensity(const Evaluation& temperature, const Evaluation& pressure)
00233 {
00234 return gasDensity(temperature, pressure) / molarMass();
00235 }
00236
00246 template <class Evaluation>
00247 static Evaluation molarLiquidDensity(Evaluation temperature, const Evaluation& )
00248 {
00249
00250
00251
00252 temperature = Opm::min(temperature, 500.0);
00253 temperature = Opm::max(temperature, 250.0);
00254
00255 const Scalar A1 = 0.25919;
00256 const Scalar A2 = 0.0014569;
00257 const Evaluation& expo = 1.0 + Opm::pow((1.0 - temperature/criticalTemperature()), (2.0/7.0));
00258 return 1.0/(A2*Opm::pow(A1, expo));
00259 }
00260
00264 template <class Evaluation>
00265 static Evaluation liquidDensity(const Evaluation& temperature, const Evaluation& pressure)
00266 { return molarLiquidDensity(temperature, pressure)*molarMass(); }
00267
00271 static bool gasIsCompressible()
00272 { return true; }
00273
00277 static bool gasIsIdeal()
00278 { return true; }
00279
00283 static bool liquidIsCompressible()
00284 { return false; }
00285
00289 template <class Evaluation>
00290 static Evaluation gasViscosity(Evaluation temperature, const Evaluation& )
00291 {
00292 temperature = Opm::min(temperature, 500.0);
00293 temperature = Opm::max(temperature, 250.0);
00294
00295 const Evaluation& Tr = Opm::max(temperature/criticalTemperature(), 1e-10);
00296 const Scalar Fp0 = 1.0;
00297 const Scalar xi = 0.004623;
00298 const Evaluation& eta_xi = Fp0*(0.807*Opm::pow(Tr, 0.618)
00299 - 0.357*Opm::exp(-0.449*Tr)
00300 + 0.34*Opm::exp(-4.058*Tr)
00301 + 0.018);
00302 return eta_xi/xi / 1e7;
00303 }
00304
00308 template <class Evaluation>
00309 static Evaluation liquidViscosity(Evaluation temperature, const Evaluation& )
00310 {
00311 temperature = Opm::min(temperature, 500.0);
00312 temperature = Opm::max(temperature, 250.0);
00313
00314 const Scalar A = -3.82;
00315 const Scalar B = 1027.0;
00316 const Scalar C = -6.38e-4;
00317 const Scalar D = 4.52e-7;
00318
00319 return 1e-3*Opm::exp(A
00320 + B/temperature
00321 + C*temperature
00322 + D*temperature*temperature);
00323 }
00324 };
00325
00326 }
00327
00328 #endif