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_REGION2_HPP
00028 #define OPM_IAPWS_REGION2_HPP
00029
00030 #include <opm/material/common/MathToolbox.hpp>
00031
00032 #include <cmath>
00033
00034 namespace Opm {
00035 namespace IAPWS {
00050 template <class Scalar>
00051 class Region2
00052 {
00053 public:
00061 template <class Evaluation>
00062 static bool isValid(const Evaluation& temperature, const Evaluation& pressure)
00063 {
00064 return
00065 temperature <= 623.15 && pressure <= 100e6;
00066
00067
00068
00069
00070
00071
00072
00073
00074 }
00075
00081 template <class Evaluation>
00082 static Evaluation tau(const Evaluation& temperature)
00083 { return 540.0 / temperature; }
00084
00091 template <class Evaluation>
00092 static Evaluation dtau_dT(const Evaluation& temperature)
00093 { return - 540.0 / (temperature*temperature); }
00094
00100 template <class Evaluation>
00101 static Evaluation pi(const Evaluation& pressure)
00102 { return pressure / 1e6; }
00103
00110 template <class Evaluation>
00111 static Scalar dpi_dp(const Evaluation& )
00112 { return 1.0 / 1e6; }
00113
00120 template <class Evaluation>
00121 static Evaluation dp_dpi(const Evaluation& )
00122 { return 1e6; }
00123
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;
00142
00143
00144 result = Opm::log(pi_);
00145 for (int i = 0; i < 9; ++i)
00146 result += n_g(i)*Opm::pow(tau_, J_g(i));
00147
00148
00149 for (int i = 0; i < 43; ++i)
00150 result +=
00151 n_r(i)*
00152 Opm::pow(pi_, I_r(i))*
00153 Opm::pow(tau_ - 0.5, J_r(i));
00154 return result;
00155 }
00156
00169 template <class Evaluation>
00170 static Evaluation dgamma_dtau(const Evaluation& temperature, const Evaluation& pressure)
00171 {
00172 const Evaluation& tau_ = tau(temperature);
00173 const Evaluation& pi_ = pi(pressure);
00174
00175
00176 Evaluation result = 0.0;
00177 for (int i = 0; i < 9; i++) {
00178 result +=
00179 n_g(i) *
00180 J_g(i) *
00181 Opm::pow(tau_, static_cast<Scalar>(J_g(i) - 1));
00182 }
00183
00184
00185 for (int i = 0; i < 43; i++) {
00186 result +=
00187 n_r(i) *
00188 Opm::pow(pi_, static_cast<Scalar>(I_r(i))) *
00189 J_r(i) *
00190 Opm::pow(tau_ - 0.5, static_cast<Scalar>(J_r(i) - 1));
00191 }
00192
00193 return result;
00194 }
00195
00208 template <class Evaluation>
00209 static Evaluation dgamma_dpi(const Evaluation& temperature, const Evaluation& pressure)
00210 {
00211 const Evaluation& tau_ = tau(temperature);
00212 const Evaluation& pi_ = pi(pressure);
00213
00214
00215 Evaluation result = 1/pi_;
00216
00217
00218 for (int i = 0; i < 43; i++) {
00219 result +=
00220 n_r(i) *
00221 I_r(i) *
00222 Opm::pow(pi_, static_cast<Scalar>(I_r(i) - 1)) *
00223 Opm::pow(tau_ - 0.5, static_cast<Scalar>(J_r(i)));
00224 }
00225
00226 return result;
00227 }
00228
00241 template <class Evaluation>
00242 static Evaluation ddgamma_dtaudpi(const Evaluation& temperature, const Evaluation& pressure)
00243 {
00244 const Evaluation& tau_ = tau(temperature);
00245 const Evaluation& pi_ = pi(pressure);
00246
00247
00248 Evaluation result = 0.0;
00249
00250
00251 for (int i = 0; i < 43; i++) {
00252 result +=
00253 n_r(i) *
00254 I_r(i) *
00255 J_r(i) *
00256 Opm::pow(pi_, static_cast<Scalar>(I_r(i) - 1)) *
00257 Opm::pow(tau_ - 0.5, static_cast<Scalar>(J_r(i) - 1));
00258 }
00259
00260 return result;
00261 }
00262
00275 template <class Evaluation>
00276 static Evaluation ddgamma_ddpi(const Evaluation& temperature, const Evaluation& pressure)
00277 {
00278 const Evaluation& tau_ = tau(temperature);
00279 const Evaluation& pi_ = pi(pressure);
00280
00281
00282 Evaluation result = -1/(pi_*pi_);
00283
00284
00285 for (int i = 0; i < 43; i++) {
00286 result +=
00287 n_r(i) *
00288 I_r(i) *
00289 (I_r(i) - 1) *
00290 Opm::pow(pi_, static_cast<Scalar>(I_r(i) - 2)) *
00291 Opm::pow(tau_ - 0.5, static_cast<Scalar>(J_r(i)));
00292 }
00293
00294 return result;
00295 }
00296
00309 template <class Evaluation>
00310 static Evaluation ddgamma_ddtau(const Evaluation& temperature, const Evaluation& pressure)
00311 {
00312 const Evaluation& tau_ = tau(temperature);
00313 const Evaluation& pi_ = pi(pressure);
00314
00315
00316 Evaluation result = 0.0;
00317 for (int i = 0; i < 9; i++) {
00318 result +=
00319 n_g(i) *
00320 J_g(i) *
00321 (J_g(i) - 1) *
00322 Opm::pow(tau_, static_cast<Scalar>(J_g(i) - 2));
00323 }
00324
00325
00326 for (int i = 0; i < 43; i++) {
00327 result +=
00328 n_r(i) *
00329 Opm::pow(pi_, I_r(i)) *
00330 J_r(i) *
00331 (J_r(i) - 1.) *
00332 Opm::pow(tau_ - 0.5, static_cast<Scalar>(J_r(i) - 2));
00333 }
00334
00335 return result;
00336 }
00337
00338
00339 private:
00340 static Scalar n_g(int i)
00341 {
00342 static const Scalar n[9] = {
00343 -0.96927686500217e1, 0.10086655968018e2, -0.56087911283020e-2,
00344 0.71452738081455e-1, -0.40710498223928, 0.14240819171444e1,
00345 -0.43839511319450e1, -0.28408632460772, 0.21268463753307e-1
00346 };
00347 return n[i];
00348 }
00349
00350 static Scalar n_r(int i)
00351 {
00352 static const Scalar n[43] = {
00353 -0.17731742473213e-2, -0.17834862292358e-1, -0.45996013696365e-1,
00354 -0.57581259083432e-1, -0.50325278727930e-1, -0.33032641670203e-4,
00355 -0.18948987516315e-3, -0.39392777243355e-2, -0.43797295650573e-1,
00356 -0.26674547914087e-4, 0.20481737692309e-7, 0.43870667284435e-6,
00357 -0.32277677238570e-4, -0.15033924542148e-2, -0.40668253562649e-1,
00358 -0.78847309559367e-9, 0.12790717852285e-7, 0.48225372718507e-6,
00359 0.22922076337661e-5, -0.16714766451061e-10, -0.21171472321355e-2,
00360 -0.23895741934104e2, -0.59059564324270e-17, -0.12621808899101e-5,
00361 -0.38946842435739e-1, 0.11256211360459e-10, -0.82311340897998e1,
00362 0.19809712802088e-7, 0.10406965210174e-18, -0.10234747095929e-12,
00363 -0.10018179379511e-8, -0.80882908646985e-10, 0.10693031879409,
00364 -0.33662250574171, 0.89185845355421e-24, 0.30629316876232e-12,
00365 -0.42002467698208e-5, -0.59056029685639e-25, 0.37826947613457e-5,
00366 -0.12768608934681e-14, 0.73087610595061e-28, 0.55414715350778e-16,
00367 -0.94369707241210e-6
00368 };
00369 return n[i];
00370 }
00371
00372 static Scalar I_r(int i)
00373 {
00374 static const short int I[43] = {
00375 1, 1, 1,
00376 1, 1, 2,
00377 2, 2, 2,
00378 2, 3, 3,
00379 3, 3, 3,
00380 4, 4, 4,
00381 5, 6, 6,
00382 6, 7, 7,
00383 7, 8, 8,
00384 9, 10, 10,
00385 10, 16, 16,
00386 18, 20, 20,
00387 20, 21, 22,
00388 23, 24, 24,
00389 24
00390 };
00391 return I[i];
00392 }
00393
00394 static Scalar J_g(int i)
00395 {
00396 static const short int J[9] = {
00397 0, 1, -5,
00398 -4, -3, -2,
00399 -1, 2, 3
00400 };
00401 return J[i];
00402 }
00403
00404 static Scalar J_r(int i)
00405 {
00406 static const short int J[43] = {
00407 0, 1, 2,
00408 3, 6, 1,
00409 2, 4, 7,
00410 36, 0, 1,
00411 3, 6, 35,
00412 1, 2, 3,
00413 7, 3, 16,
00414 35, 0, 11,
00415 25, 8, 36,
00416 13, 4, 10,
00417 14, 29, 50,
00418 57, 20, 35,
00419 48, 21, 53,
00420 39, 26, 40,
00421 58
00422 };
00423 return J[i];
00424 }
00425
00426 };
00427
00428 }
00429 }
00430
00431 #endif