00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00028 #ifndef OPM_BINARY_COEFF_H2O_CO2_HPP
00029 #define OPM_BINARY_COEFF_H2O_CO2_HPP
00030
00031 #include <opm/material/binarycoefficients/HenryIapws.hpp>
00032 #include <opm/material/binarycoefficients/FullerMethod.hpp>
00033
00034 #include <opm/material/components/H2O.hpp>
00035 #include <opm/material/components/SimpleCO2.hpp>
00036
00037 namespace Opm {
00038 namespace BinaryCoeff {
00039
00044 class H2O_CO2
00045 {
00046 public:
00057 template <class Scalar, class Evaluation = Scalar>
00058 static Evaluation henry(const Evaluation& temperature)
00059 {
00060 const Scalar E = 1672.9376;
00061 const Scalar F = 28.1751;
00062 const Scalar G = -112.4619;
00063 const Scalar H = 85.3807;
00064
00065 return henryIAPWS(E, F, G, H, temperature);
00066 }
00067
00073 template <class Scalar, class Evaluation = Scalar>
00074 static Evaluation gasDiffCoeff(const Evaluation& temperature, const Evaluation& pressure)
00075 {
00076 typedef Opm::H2O<Scalar> H2O;
00077 typedef Opm::SimpleCO2<Scalar> CO2;
00078
00079
00080 const Scalar SigmaNu[2] = { 13.1 , 26.9 };
00081
00082 const Scalar M[2] = { H2O::molarMass()*1e3, CO2::molarMass()*1e3 };
00083
00084 return fullerMethod(M, SigmaNu, temperature, pressure);
00085 }
00086
00090 template <class Scalar, class Evaluation = Scalar>
00091 static Evaluation liquidDiffCoeff(const Evaluation& , const Evaluation& )
00092 { OPM_THROW(std::runtime_error, "Not implemented: Binary liquid diffusion coefficients of CO2 and CH4"); }
00093 };
00094
00095 }
00096 }
00097
00098 #endif