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_REGION1_HPP
00028 #define OPM_IAPWS_REGION1_HPP
00029
00030 #include <opm/material/common/MathToolbox.hpp>
00031
00032 #include <cmath>
00033
00034 namespace Opm {
00035 namespace IAPWS {
00049 template <class Scalar>
00050 class Region1
00051 {
00052 public:
00060 template <class Evaluation>
00061 static bool isValid(const Evaluation& temperature, const Evaluation& pressure)
00062 {
00063 return temperature <= 623.15 && pressure <= 100e6;
00064
00065
00066
00067
00068
00069
00070
00071
00072
00073
00074
00075 }
00076
00082 template <class Evaluation>
00083 static Evaluation tau(const Evaluation& temperature)
00084 { return 1386.0 / temperature; }
00085
00092 template <class Evaluation>
00093 static Evaluation dtau_dT(const Evaluation& temperature)
00094 { return - 1386.0 / (temperature*temperature); }
00095
00101 template <class Evaluation>
00102 static Evaluation pi(const Evaluation& pressure)
00103 { return pressure / 16.53e6; }
00104
00111 template <class Evaluation>
00112 static Scalar dpi_dp(const Evaluation& )
00113 { return 1.0 / 16.53e6; }
00114
00121 template <class Evaluation>
00122 static Scalar dp_dpi(const Evaluation& )
00123 { return 16.53e6; }
00124
00135 template <class Evaluation>
00136 static Evaluation gamma(const Evaluation& temperature, const Evaluation& pressure)
00137 {
00138 const Evaluation tau_ = tau(temperature);
00139 const Evaluation pi_ = pi(pressure);
00140
00141 Evaluation result = 0;
00142 for (int i = 0; i < 34; ++i) {
00143 result += n(i)*Opm::pow(7.1 - pi_, I(i))*Opm::pow(tau_ - 1.222, J(i));
00144 }
00145
00146 return result;
00147 }
00148
00149
00161 template <class Evaluation>
00162 static Evaluation dgamma_dtau(const Evaluation& temperature, const Evaluation& pressure)
00163 {
00164 const Evaluation tau_ = tau(temperature);
00165 const Evaluation pi_ = pi(pressure);
00166
00167 Evaluation result = 0.0;
00168 for (int i = 0; i < 34; i++) {
00169 result +=
00170 n(i) *
00171 Opm::pow(7.1 - pi_, static_cast<Scalar>(I(i))) *
00172 Opm::pow(tau_ - 1.222, static_cast<Scalar>(J(i)-1)) *
00173 J(i);
00174 }
00175
00176 return result;
00177 }
00178
00190 template <class Evaluation>
00191 static Evaluation dgamma_dpi(const Evaluation& temperature, const Evaluation& pressure)
00192 {
00193 const Evaluation tau_ = tau(temperature);
00194 const Evaluation pi_ = pi(pressure);
00195
00196 Evaluation result = 0.0;
00197 for (int i = 0; i < 34; i++) {
00198 result +=
00199 -n(i) *
00200 I(i) *
00201 Opm::pow(7.1 - pi_, static_cast<Scalar>(I(i) - 1)) *
00202 Opm::pow(tau_ - 1.222, static_cast<Scalar>(J(i)));
00203 }
00204
00205 return result;
00206 }
00207
00220 template <class Evaluation>
00221 static Evaluation ddgamma_dtaudpi(const Evaluation& temperature, const Evaluation& pressure)
00222 {
00223 const Evaluation tau_ = tau(temperature);
00224 const Evaluation pi_ = pi(pressure);
00225
00226 Evaluation result = 0.0;
00227 for (int i = 0; i < 34; i++) {
00228 result +=
00229 -n(i) *
00230 I(i) *
00231 J(i) *
00232 Opm::pow(7.1 - pi_, static_cast<Scalar>(I(i) - 1)) *
00233 Opm::pow(tau_ - 1.222, static_cast<Scalar>(J(i) - 1));
00234 }
00235
00236 return result;
00237 }
00238
00251 template <class Evaluation>
00252 static Evaluation ddgamma_ddpi(const Evaluation& temperature, const Evaluation& pressure)
00253 {
00254 const Evaluation tau_ = tau(temperature);
00255 const Evaluation pi_ = pi(pressure);
00256
00257 Evaluation result = 0.0;
00258 for (int i = 0; i < 34; i++) {
00259 result +=
00260 n(i) *
00261 I(i) *
00262 (I(i) - 1) *
00263 Opm::pow(7.1 - pi_, I(i) - 2) *
00264 Opm::pow(tau_ - 1.222, J(i));
00265 }
00266
00267 return result;
00268 }
00269
00281 template <class Evaluation>
00282 static Evaluation ddgamma_ddtau(const Evaluation& temperature, const Evaluation& pressure)
00283 {
00284 const Evaluation tau_ = tau(temperature);
00285 const Evaluation pi_ = pi(pressure);
00286
00287 Evaluation result = 0.0;
00288 for (int i = 0; i < 34; i++) {
00289 result +=
00290 n(i) *
00291 Opm::pow(7.1 - pi_, I(i)) *
00292 J(i) *
00293 (J(i) - 1) *
00294 Opm::pow(tau_ - 1.222, J(i) - 2);
00295 }
00296
00297 return result;
00298 }
00299
00300 private:
00301 static Scalar n(int i)
00302 {
00303 static const Scalar n[34] = {
00304 0.14632971213167, -0.84548187169114, -0.37563603672040e1,
00305 0.33855169168385e1, -0.95791963387872, 0.15772038513228,
00306 -0.16616417199501e-1, 0.81214629983568e-3, 0.28319080123804e-3,
00307 -0.60706301565874e-3, -0.18990068218419e-1, -0.32529748770505e-1,
00308 -0.21841717175414e-1, -0.52838357969930e-4, -0.47184321073267e-3,
00309 -0.30001780793026e-3, 0.47661393906987e-4, -0.44141845330846e-5,
00310 -0.72694996297594e-15,-0.31679644845054e-4, -0.28270797985312e-5,
00311 -0.85205128120103e-9, -0.22425281908000e-5, -0.65171222895601e-6,
00312 -0.14341729937924e-12,-0.40516996860117e-6, -0.12734301741641e-8,
00313 -0.17424871230634e-9, -0.68762131295531e-18, 0.14478307828521e-19,
00314 0.26335781662795e-22,-0.11947622640071e-22, 0.18228094581404e-23,
00315 -0.93537087292458e-25
00316 };
00317 return n[i];
00318 }
00319
00320 static Scalar I(int i)
00321 {
00322 static const short int I[34] = {
00323 0, 0, 0,
00324 0, 0, 0,
00325 0, 0, 1,
00326 1, 1, 1,
00327 1, 1, 2,
00328 2, 2, 2,
00329 2, 3, 3,
00330 3, 4, 4,
00331 4, 5, 8,
00332 8, 21, 23,
00333 29, 30, 31,
00334 32
00335 };
00336 return I[i];
00337 }
00338
00339 static Scalar J(int i)
00340 {
00341 static const short int J[34] = {
00342 -2, -1, 0,
00343 1, 2, 3,
00344 4, 5, -9,
00345 -7, -1, 0,
00346 1, 3, -3,
00347 0, 1, 3,
00348 17, -4, 0,
00349 6, -5, -2,
00350 10, -8, -11,
00351 -6, -29, -31,
00352 -38, -39, -40,
00353 -41
00354 };
00355 return J[i];
00356 }
00357
00358 };
00359
00360 }
00361 }
00362
00363 #endif