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_HENRY_IAPWS_HPP
00029 #define OPM_HENRY_IAPWS_HPP
00030
00031 #include <opm/material/components/H2O.hpp>
00032
00033 namespace Opm
00034 {
00044 template <class Scalar, class Evaluation>
00045 inline Evaluation henryIAPWS(Scalar E,
00046 Scalar F,
00047 Scalar G,
00048 Scalar H,
00049 const Evaluation& temperature)
00050 {
00051 typedef Opm::H2O<Evaluation> H2O;
00052
00053 Evaluation Tr = temperature/H2O::criticalTemperature();
00054 Evaluation tau = 1 - Tr;
00055
00056 static const Scalar c[6] = {
00057 1.99274064, 1.09965342, -0.510839303,
00058 -1.75493479,-45.5170352, -6.7469445e5
00059 };
00060 static const Scalar d[6] = {
00061 1/3.0, 2/3.0, 5/3.0,
00062 16/3.0, 43/3.0, 110/3.0
00063 };
00064 static const Scalar q = -0.023767;
00065
00066 Evaluation f = 0;
00067 for (int i = 0; i < 6; ++i) {
00068 f += c[i]*Opm::pow(tau, d[i]);
00069 }
00070
00071 const Evaluation& exponent =
00072 q*F +
00073 E/temperature*f +
00074 (F +
00075 G*Opm::pow(tau, 2.0/3) +
00076 H*tau)*
00077 Opm::exp((H2O::tripleTemperature() - temperature)/100);
00078
00079
00080
00081 return Opm::exp(exponent)*H2O::vaporPressure(temperature);
00082 }
00083 }
00084
00085 #endif // OPM_HENRY_IAPWS_HPP