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_MESITYLENE_HPP
00028 #define OPM_BINARY_COEFF_AIR_MESITYLENE_HPP
00029
00030 #include <opm/material/components/Air.hpp>
00031 #include <opm/material/components/Mesitylene.hpp>
00032
00033 namespace Opm
00034 {
00035 namespace BinaryCoeff
00036 {
00037
00041 class Air_Mesitylene
00042 {
00043 public:
00047 template <class Evaluation>
00048 static Evaluation henry(const Evaluation& )
00049 { OPM_THROW(std::runtime_error, "Not implemented: Henry coefficient of air in mesitylene"); }
00050
00058 template <class Evaluation>
00059 static Evaluation gasDiffCoeff(Evaluation temperature, Evaluation pressure)
00060 {
00061 typedef Opm::Air<double> Air;
00062 typedef Opm::Mesitylene<double> Mesitylene;
00063
00064 temperature = Opm::max(temperature, 1e-9);
00065 temperature = Opm::min(temperature, 500.0);
00066 pressure = Opm::max(pressure, 0.0);
00067 pressure = Opm::min(pressure, 1e8);
00068
00069 const double M_m = 1e3*Mesitylene::molarMass();
00070 const double M_a = 1e3*Air::molarMass();
00071 const double Tb_m = 437.9;
00072 const double sigma_a = 3.711;
00073 const double T_scal_a = 78.6;
00074 const double V_B_m = 162.6;
00075 const double sigma_m = 1.18*std::pow(V_B_m, 0.333);
00076 const double sigma_am = 0.5*(sigma_a + sigma_m);
00077 const double T_scal_m = 1.15*Tb_m;
00078 const double T_scal_am = std::sqrt(T_scal_a*T_scal_m);
00079
00080 Evaluation T_star = temperature/T_scal_am;
00081 T_star = Opm::max(T_star, 1e-5);
00082
00083 const Evaluation Omega = 1.06036/Opm::pow(T_star, 0.1561) + 0.193/Opm::exp(T_star*0.47635)
00084 + 1.03587/Opm::exp(T_star*1.52996) + 1.76474/Opm::exp(T_star*3.89411);
00085 const double B_ = 0.00217 - 0.0005*std::sqrt(1.0/M_a + 1.0/M_m);
00086 const double Mr = (M_a + M_m)/(M_a*M_m);
00087 const Evaluation D_am = (B_*Opm::pow(temperature, 1.5) * std::sqrt(Mr))
00088 /(1e-5*pressure*std::pow(sigma_am, 2.0) * Omega);
00089
00090 return 1e-4*D_am;
00091 }
00092
00096 template <class Evaluation>
00097 static Evaluation liquidDiffCoeff(const Evaluation& , const Evaluation& )
00098 { OPM_THROW(std::runtime_error, "Not implemented: Binary liquid diffusion coefficients of air and mesitylene"); }
00099 };
00100
00101 }
00102 }
00103
00104 #endif