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_N2_HPP
00028 #define OPM_BINARY_COEFF_H2O_N2_HPP
00029
00030 #include "HenryIapws.hpp"
00031 #include "FullerMethod.hpp"
00032
00033 #include <opm/material/components/N2.hpp>
00034 #include <opm/material/components/H2O.hpp>
00035
00036 namespace Opm {
00037 namespace BinaryCoeff {
00038
00043 class H2O_N2
00044 {
00045 public:
00051 template <class Evaluation>
00052 static Evaluation henry(const Evaluation& temperature)
00053 {
00054 const double E = 2388.8777;
00055 const double F = -14.9593;
00056 const double G = 42.0179;
00057 const double H = -29.4396;
00058
00059 return henryIAPWS(E, F, G, H, temperature);
00060 }
00061
00069 template <class Evaluation>
00070 static Evaluation gasDiffCoeff(const Evaluation& temperature, const Evaluation& pressure)
00071 {
00072 typedef Opm::H2O<double> H2O;
00073 typedef Opm::N2<double> N2;
00074
00075
00076 const double SigmaNu[2] = { 13.1 , 18.5 };
00077
00078 const double M[2] = { H2O::molarMass()*1e3, N2::molarMass()*1e3 };
00079
00080 return fullerMethod(M, SigmaNu, temperature, pressure);
00081 }
00082
00101 template <class Evaluation>
00102 static Evaluation liquidDiffCoeff(const Evaluation& temperature, const Evaluation& )
00103 {
00104 const double Texp = 273.15 + 25;
00105 const double Dexp = 2.01e-9;
00106
00107 return Dexp * temperature/Texp;
00108 }
00109 };
00110
00111 }
00112 }
00113
00114 #endif