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_MESITYLENE_HPP
00028 #define OPM_BINARY_COEFF_H2O_MESITYLENE_HPP
00029
00030 #include <opm/material/components/H2O.hpp>
00031 #include <opm/material/components/Mesitylene.hpp>
00032
00033 namespace Opm
00034 {
00035 namespace BinaryCoeff
00036 {
00037
00041 class H2O_Mesitylene
00042 {
00043 public:
00051 template <class Evaluation>
00052 static Evaluation henry(const Evaluation& )
00053 {
00054
00055 double sanderH = 1.7e-1;
00056
00057 double opmH = sanderH / 101.325;
00058 opmH *= 18.02e-6;
00059 return 1.0/opmH;
00060 }
00061
00066 template <class Evaluation>
00067 static Evaluation gasDiffCoeff(Evaluation temperature, Evaluation pressure)
00068 {
00069 typedef Opm::H2O<double> H2O;
00070 typedef Opm::Mesitylene<double> Mesitylene;
00071
00072 temperature = Opm::max(temperature, 1e-9);
00073 temperature = Opm::min(temperature, 500.0);
00074 pressure = Opm::max(pressure, 0.0);
00075 pressure = Opm::min(pressure, 1e8);
00076
00077 const double M_m = 1e3*Mesitylene::molarMass();
00078 const double M_w = 1e3*H2O::molarMass();
00079 const double Tb_m = 437.9;
00080 const double Tb_w = 373.15;
00081 const double V_B_w = 18.0;
00082 const double sigma_w = 1.18*std::pow(V_B_w, 0.333);
00083 const double T_scal_w = 1.15*Tb_w;
00084 const double V_B_m = 162.6;
00085 const double sigma_m = 1.18*std::pow(V_B_m, 0.333);
00086 const double sigma_wm = 0.5*(sigma_w + sigma_m);
00087 const double T_scal_m = 1.15*Tb_m;
00088 const double T_scal_wm = std::sqrt(T_scal_w*T_scal_m);
00089
00090 const Evaluation T_star = Opm::max(temperature/T_scal_wm, 1e-5);
00091
00092 const Evaluation& Omega = 1.06036/Opm::pow(T_star,0.1561) + 0.193/Opm::exp(T_star*0.47635)
00093 + 1.03587/Opm::exp(T_star*1.52996) + 1.76474/Opm::exp(T_star*3.89411);
00094 const double B_ = 0.00217 - 0.0005*std::sqrt(1.0/M_w + 1.0/M_m);
00095 const double Mr = (M_w + M_m)/(M_w*M_m);
00096 const Evaluation D_wm = (B_*Opm::pow(temperature,1.6)*std::sqrt(Mr))
00097 /(1e-5*pressure*std::pow(sigma_wm, 2.0)*Omega);
00098
00099 return D_wm*1e-4;
00100 }
00101
00105 template <class Evaluation>
00106 static Evaluation liquidDiffCoeff(const Evaluation& , const Evaluation& )
00107 {
00108
00109 return 1e-9;
00110 }
00111 };
00112
00113 }
00114 }
00115
00116 #endif