Go to the documentation of this file.
27#ifndef OPM_IAPWS_REGION2_HPP
28#define OPM_IAPWS_REGION2_HPP
50template < class Scalar>
61 template < class Evaluation>
62 static bool isValid( const Evaluation& temperature, const Evaluation& pressure)
65 temperature <= 623.15 && pressure <= 100e6;
81 template < class Evaluation>
82 static Evaluation tau( const Evaluation& temperature)
83 { return 540.0 / temperature; }
91 template < class Evaluation>
92 static Evaluation dtau_dT( const Evaluation& temperature)
93 { return - 540.0 / (temperature*temperature); }
100 template < class Evaluation>
101 static Evaluation pi( const Evaluation& pressure)
102 { return pressure / 1e6; }
110 template < class Evaluation>
112 { return 1.0 / 1e6; }
120 template < class Evaluation>
121 static Evaluation dp_dpi( const Evaluation& )
135 template < class Evaluation>
136 static Evaluation gamma( const Evaluation& temperature, const Evaluation& pressure)
138 const Evaluation& tau_ = tau(temperature);
139 const Evaluation& pi_ = pi(pressure);
145 for ( int i = 0; i < 9; ++i)
146 result += n_g(i)* pow(tau_, J_g(i));
149 for ( int i = 0; i < 43; ++i)
153 pow(tau_ - 0.5, J_r(i));
169 template < class Evaluation>
170 static Evaluation dgamma_dtau( const Evaluation& temperature, const Evaluation& pressure)
172 const Evaluation& tau_ = tau(temperature);
173 const Evaluation& pi_ = pi(pressure);
176 Evaluation result = 0.0;
177 for ( int i = 0; i < 9; i++) {
181 pow(tau_, static_cast<Scalar >(J_g(i) - 1));
185 for ( int i = 0; i < 43; i++) {
188 pow(pi_, static_cast<Scalar >(I_r(i))) *
190 pow(tau_ - 0.5, static_cast<Scalar >(J_r(i) - 1));
208 template < class Evaluation>
209 static Evaluation dgamma_dpi( const Evaluation& temperature, const Evaluation& pressure)
211 const Evaluation& tau_ = tau(temperature);
212 const Evaluation& pi_ = pi(pressure);
215 Evaluation result = 1/pi_;
218 for ( int i = 0; i < 43; i++) {
222 pow(pi_, static_cast<Scalar >(I_r(i) - 1)) *
223 pow(tau_ - 0.5, static_cast<Scalar >(J_r(i)));
241 template < class Evaluation>
242 static Evaluation ddgamma_dtaudpi( const Evaluation& temperature, const Evaluation& pressure)
244 const Evaluation& tau_ = tau(temperature);
245 const Evaluation& pi_ = pi(pressure);
248 Evaluation result = 0.0;
251 for ( int i = 0; i < 43; i++) {
256 pow(pi_, static_cast<Scalar >(I_r(i) - 1)) *
257 pow(tau_ - 0.5, static_cast<Scalar >(J_r(i) - 1));
275 template < class Evaluation>
276 static Evaluation ddgamma_ddpi( const Evaluation& temperature, const Evaluation& pressure)
278 const Evaluation& tau_ = tau(temperature);
279 const Evaluation& pi_ = pi(pressure);
282 Evaluation result = -1/(pi_*pi_);
285 for ( int i = 0; i < 43; i++) {
290 pow(pi_, static_cast<Scalar >(I_r(i) - 2)) *
291 pow(tau_ - 0.5, static_cast<Scalar >(J_r(i)));
309 template < class Evaluation>
310 static Evaluation ddgamma_ddtau( const Evaluation& temperature, const Evaluation& pressure)
312 const Evaluation& tau_ = tau(temperature);
313 const Evaluation& pi_ = pi(pressure);
316 Evaluation result = 0.0;
317 for ( int i = 0; i < 9; i++) {
322 pow(tau_, static_cast<Scalar >(J_g(i) - 2));
326 for ( int i = 0; i < 43; i++) {
332 pow(tau_ - 0.5, static_cast<Scalar >(J_r(i) - 2));
340 static Scalar n_g( int i)
342 static const Scalar n[9] = {
343 -0.96927686500217e1, 0.10086655968018e2, -0.56087911283020e-2,
344 0.71452738081455e-1, -0.40710498223928, 0.14240819171444e1,
345 -0.43839511319450e1, -0.28408632460772, 0.21268463753307e-1
350 static Scalar n_r( int i)
352 static const Scalar n[43] = {
353 -0.17731742473213e-2, -0.17834862292358e-1, -0.45996013696365e-1,
354 -0.57581259083432e-1, -0.50325278727930e-1, -0.33032641670203e-4,
355 -0.18948987516315e-3, -0.39392777243355e-2, -0.43797295650573e-1,
356 -0.26674547914087e-4, 0.20481737692309e-7, 0.43870667284435e-6,
357 -0.32277677238570e-4, -0.15033924542148e-2, -0.40668253562649e-1,
358 -0.78847309559367e-9, 0.12790717852285e-7, 0.48225372718507e-6,
359 0.22922076337661e-5, -0.16714766451061e-10, -0.21171472321355e-2,
360 -0.23895741934104e2, -0.59059564324270e-17, -0.12621808899101e-5,
361 -0.38946842435739e-1, 0.11256211360459e-10, -0.82311340897998e1,
362 0.19809712802088e-7, 0.10406965210174e-18, -0.10234747095929e-12,
363 -0.10018179379511e-8, -0.80882908646985e-10, 0.10693031879409,
364 -0.33662250574171, 0.89185845355421e-24, 0.30629316876232e-12,
365 -0.42002467698208e-5, -0.59056029685639e-25, 0.37826947613457e-5,
366 -0.12768608934681e-14, 0.73087610595061e-28, 0.55414715350778e-16,
372 static Scalar I_r( int i)
374 static const short int I[43] = {
394 static Scalar J_g( int i)
396 static const short int J[9] = {
404 static Scalar J_r( int i)
406 static const short int J[43] = {
Implements the equations for region 2 of the IAPWS '97 formulation. Definition: Region2.hpp:52
static Scalar dpi_dp(const Evaluation &) Returns the derivative of the reduced pressure to the pressure for IAPWS region 2 in . Definition: Region2.hpp:111
static Evaluation ddgamma_ddtau(const Evaluation &temperature, const Evaluation &pressure) The second partial derivative of the Gibbs free energy to the normalized temperature for IAPWS region... Definition: Region2.hpp:310
static Evaluation ddgamma_ddpi(const Evaluation &temperature, const Evaluation &pressure) The second partial derivative of the Gibbs free energy to the normalized pressure for IAPWS region 2 ... Definition: Region2.hpp:276
static Evaluation dtau_dT(const Evaluation &temperature) Returns the derivative of the reduced temperature to the temperature for IAPWS region 2. Definition: Region2.hpp:92
static Evaluation pi(const Evaluation &pressure) Returns the reduced pressure (dimensionless) for IAPWS region 2. Definition: Region2.hpp:101
static Evaluation dp_dpi(const Evaluation &) Returns the derivative of the pressure to the reduced pressure for IAPWS region 2 (dimensionless). Definition: Region2.hpp:121
static Evaluation dgamma_dtau(const Evaluation &temperature, const Evaluation &pressure) The partial derivative of the Gibbs free energy to the normalized temperature for IAPWS region 2 (i.... Definition: Region2.hpp:170
static Evaluation tau(const Evaluation &temperature) Returns the reduced temperature (dimensionless) for IAPWS region 2. Definition: Region2.hpp:82
static Evaluation gamma(const Evaluation &temperature, const Evaluation &pressure) The Gibbs free energy for IAPWS region 2 (i.e. sub-critical steam) (dimensionless). Definition: Region2.hpp:136
static Evaluation ddgamma_dtaudpi(const Evaluation &temperature, const Evaluation &pressure) The partial derivative of the Gibbs free energy to the normalized pressure and to the normalized temp... Definition: Region2.hpp:242
static bool isValid(const Evaluation &temperature, const Evaluation &pressure) Returns true if IAPWS region 2 applies for a (temperature, pressure) pair. Definition: Region2.hpp:62
static Evaluation dgamma_dpi(const Evaluation &temperature, const Evaluation &pressure) The partial derivative of the Gibbs free energy to the normalized pressure for IAPWS region 2 (i.... Definition: Region2.hpp:209
Definition: Air_Mesitylene.hpp:34
Evaluation log(const Evaluation &value) Definition: MathToolbox.hpp:407
ReturnEval_< Evaluation1, Evaluation2 >::type pow(const Evaluation1 &base, const Evaluation2 &exp) Definition: MathToolbox.hpp:416
|