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_BINARY_COEFF_AIR_XYLENE_HPP
00028 #define OPM_BINARY_COEFF_AIR_XYLENE_HPP
00029
00030 #include <opm/material/components/Air.hpp>
00031 #include <opm/material/components/Xylene.hpp>
00032
00033 namespace Opm {
00034 namespace BinaryCoeff {
00035
00039 class Air_Xylene
00040 {
00041 public:
00045 template <class Evaluation>
00046 static Evaluation henry(const Evaluation& )
00047 { OPM_THROW(std::runtime_error, "Not implemented: Henry coefficient of air in xylene"); }
00048
00056 template <class Evaluation>
00057 static Evaluation gasDiffCoeff(Evaluation temperature, Evaluation pressure)
00058 {
00059 typedef Opm::Air<double> Air;
00060 typedef Opm::Xylene<double> Xylene;
00061
00062 temperature = Opm::max(temperature, 1e-9);
00063 temperature = Opm::min(temperature, 500.0);
00064 pressure = Opm::max(pressure, 0.0);
00065 pressure = Opm::min(pressure, 1e8);
00066
00067 const double M_x = 1e3*Xylene::molarMass();
00068 const double M_a = 1e3*Air::molarMass();
00069 const double Tb_x = 412.0;
00070 const double sigma_a = 3.711;
00071 const double T_scal_a = 78.6;
00072 const double V_B_x = 140.4;
00073 const double sigma_x = 1.18*std::pow(V_B_x, 0.333);
00074 const double sigma_ax = 0.5*(sigma_a + sigma_x);
00075 const double T_scal_x = 1.15*Tb_x;
00076 const double T_scal_ax = std::sqrt(T_scal_a*T_scal_x);
00077
00078 const Evaluation& T_star = Opm::max(temperature/T_scal_ax, 1e-5);
00079
00080 const Evaluation& Omega = 1.06036/Opm::pow(T_star, 0.1561) + 0.193/Opm::exp(T_star*0.47635)
00081 + 1.03587/Opm::exp(T_star*1.52996) + 1.76474/Opm::exp(T_star*3.89411);
00082 const double B_ = 0.00217 - 0.0005*std::sqrt(1.0/M_a + 1.0/M_x);
00083 const double Mr = (M_a + M_x)/(M_a*M_x);
00084 return 1e-4
00085 *(B_*Opm::pow(temperature,1.5)*std::sqrt(Mr))
00086 /(1e-5*pressure*std::pow(sigma_ax, 2.0)*Omega);
00087 }
00088
00094 template <class Evaluation>
00095 static Evaluation liquidDiffCoeff(const Evaluation& , const Evaluation& )
00096 { OPM_THROW(std::runtime_error, "Not implemented: Binary liquid diffusion coefficients of air and xylene"); }
00097 };
00098
00099 }
00100 }
00101
00102 #endif