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_IAPWS_REGION4_HPP
00028 #define OPM_IAPWS_REGION4_HPP
00029
00030 #include <opm/material/common/MathToolbox.hpp>
00031
00032 #include <cmath>
00033
00034 namespace Opm {
00035 namespace IAPWS {
00036
00050 template <class Scalar>
00051 class Region4
00052 {
00053 public:
00062 template <class Evaluation>
00063 static Evaluation saturationPressure(const Evaluation& temperature)
00064 {
00065 static const Scalar n[10] = {
00066 0.11670521452767e4, -0.72421316703206e6, -0.17073846940092e2,
00067 0.12020824702470e5, -0.32325550322333e7, 0.14915108613530e2,
00068 -0.48232657361591e4, 0.40511340542057e6, -0.23855557567849,
00069 0.65017534844798e3
00070 };
00071
00072 const Evaluation& sigma = temperature + n[8]/(temperature - n[9]);
00073
00074 const Evaluation& A = (sigma + n[0])*sigma + n[1];
00075 const Evaluation& B = (n[2]*sigma + n[3])*sigma + n[4];
00076 const Evaluation& C = (n[5]*sigma + n[6])*sigma + n[7];
00077
00078 Evaluation tmp = 2*C/(Opm::sqrt(B*B - 4*A*C) - B);
00079 tmp *= tmp;
00080 tmp *= tmp;
00081
00082 return 1e6*tmp;
00083 }
00084
00093 template <class Evaluation>
00094 static Evaluation vaporTemperature(const Evaluation& pressure)
00095 {
00096 static const Scalar n[10] = {
00097 0.11670521452767e4, -0.72421316703206e6, -0.17073846940092e2,
00098 0.12020824702470e5, -0.32325550322333e7, 0.14915108613530e2,
00099 -0.48232657361591e4, 0.40511340542057e6, -0.23855557567849,
00100 0.65017534844798e3
00101 };
00102 const Evaluation& beta = pow((pressure/1e6 ), (1./4.));
00103 const Evaluation& beta2 = pow(beta, 2.);
00104 const Evaluation& E = beta2 + n[2] * beta + n[5];
00105 const Evaluation& F = n[0]*beta2 + n[3]*beta + n[6];
00106 const Evaluation& G = n[1]*beta2 + n[4]*beta + n[7];
00107
00108 const Evaluation& D = ( 2.*G)/(-F -Opm::sqrt(pow(F,2.) - 4.*E*G));
00109
00110 const Evaluation& temperature = (n[9] + D - Opm::sqrt(pow(n[9]+D , 2.) - 4.* (n[8] + n[9]*D)) ) * 0.5;
00111
00112 return temperature;
00113 }
00114 };
00115
00116 }
00117 }
00118
00119 #endif