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_H2O_XYLENE_HPP
00028 #define OPM_BINARY_COEFF_H2O_XYLENE_HPP
00029
00030 #include <opm/material/components/H2O.hpp>
00031 #include <opm/material/components/Xylene.hpp>
00032
00033 namespace Opm {
00034 namespace BinaryCoeff {
00035
00039 class H2O_Xylene
00040 {
00041 public:
00050 template <class Evaluation>
00051 static Evaluation henry(const Evaluation& )
00052 {
00053
00054 double sanderH = 1.5e-1;
00055
00056 double opmH = sanderH / 101.325;
00057 opmH *= 18.02e-6;
00058 return 1.0/opmH;
00059 }
00060
00065 template <class Evaluation>
00066 static Evaluation gasDiffCoeff(Evaluation temperature, Evaluation pressure)
00067 {
00068 typedef Opm::H2O<double> H2O;
00069 typedef Opm::Xylene<double> Xylene;
00070
00071 temperature = Opm::max(temperature, 1e-9);
00072 temperature = Opm::min(temperature, 500.0);
00073 pressure = Opm::max(pressure, 0.0);
00074 pressure = Opm::min(pressure, 1e8);
00075
00076 const double M_x = 1e3*Xylene::molarMass();
00077 const double M_w = 1e3*H2O::molarMass();
00078 const double Tb_x = 412.9;
00079 const double Tb_w = 373.15;
00080 const double V_B_w = 18.0;
00081 const double sigma_w = 1.18*std::pow(V_B_w, 0.333);
00082 const double T_scal_w = 1.15*Tb_w;
00083 const double V_B_x = 140.4;
00084 const double sigma_x = 1.18*std::pow(V_B_x, 0.333);
00085 const double sigma_wx = 0.5*(sigma_w + sigma_x);
00086 const double T_scal_x = 1.15*Tb_x;
00087 const double T_scal_wx = std::sqrt(T_scal_w*T_scal_x);
00088
00089 const Evaluation& T_star = Opm::max(temperature/T_scal_wx, 1e-5);
00090
00091 const Evaluation& Omega = 1.06036/Opm::pow(T_star, 0.1561) + 0.193/Opm::exp(T_star*0.47635)
00092 + 1.03587/Opm::exp(T_star*1.52996) + 1.76474/Opm::exp(T_star*3.89411);
00093 const double B_ = 0.00217 - 0.0005*std::sqrt(1.0/M_w + 1.0/M_x);
00094 const double Mr = (M_w + M_x)/(M_w*M_x);
00095 return 1e-4
00096 *(B_*Opm::pow(temperature,1.6)*std::sqrt(Mr))
00097 /(1e-5*pressure*std::pow(sigma_wx, 2.0)*Omega);
00098 }
00099
00105 template <class Evaluation>
00106 static Evaluation liquidDiffCoeff(const Evaluation& , const Evaluation& )
00107 {
00108 return 1.e-9;
00109 }
00110 };
00111
00112 }
00113 }
00114
00115 #endif