IAPWS: International Association for the Properties of Water and Steam (chemicals.iapws)¶
This module contains the core of the IAPWS-95 and IAPWS-97 standards. The objective of this module is to contain extremely fast functions to calculate several basic properties of water.
The simplest interfaces are iapws95_rho
for density calculation only and
iapws95_properties
for some basic properties.
For reporting bugs, adding feature requests, or submitting pull requests, please use the GitHub issue tracker.
IAPWS-95 Basic Solvers¶
- chemicals.iapws.iapws95_rho(T, P)[source]¶
Calculate the density of water according to the IAPWS-95 standard given a temperature T and pressure P. The phase is determined in this calculation.
- Parameters
- Returns
- rho
float
Mass density of water, [kg/m^3]
- rho
See also
Notes
There is a sudden transition at the saturation pressure between liquid and vapor density, by design.
This solution is iterative due to the nature of the equation. The solution procedure begins with IAPWS-97’s explicit equations as an initial guess, extrapolating when out of range. If the temperature is under the critical temperature, the saturation density is calculated, and used to ensure the solver begins in the feasible region. Newton’s method converges extremely, normally after 2 or 3 iterations.
Temperatures under 273.15 K are not officially supported by [1], but a solution is still attempted down to 235 K.
References
- 1(1,2)
Wagner, Wolfgang, and Andreas Pruß. “The IAPWS Formulation 1995 for the Thermodynamic Properties of Ordinary Water Substance for General and Scientific Use.” Journal of Physical and Chemical Reference Data 31, no. 2 (2002): 387-535.
Examples
>>> iapws95_rho(T=300.0, P=1e6) 996.96002269499
1 GPa and 5000 K are suggested as upper limits of [1] although there are no hardcoded limits for temperature and pressure.
>>> iapws95_rho(T=5000.0, P=1e9) 326.79451662743
- chemicals.iapws.iapws95_P(T, rho)[source]¶
Calculate the pressure of water according to the IAPWS-95 standard given a temperature T and mass density rho.
- Parameters
- Returns
- P
float
Pressure, [Pa]
- P
Notes
The IAPWS-95 model is explicit with inputs of temperature and density, so this is a direct calculation with no iteration required.
References
- 1
Wagner, Wolfgang, and Andreas Pruß. “The IAPWS Formulation 1995 for the Thermodynamic Properties of Ordinary Water Substance for General and Scientific Use.” Journal of Physical and Chemical Reference Data 31, no. 2 (2002): 387-535.
Examples
>>> iapws95_P(330.0, iapws95_rho(T=330.0, P=8e5)) 8e5 >>> iapws95_P(823.0, 40.393893559703734) 14e6
Not all temperature and density inputs provide a stable solution; for example anything between the vapor and gas saturation curves. In some but not all of these cases a negative pressure is returned:
>>> iapws95_P(T=300, rho=300) -1.526394720e+23
- chemicals.iapws.iapws95_T(P, rho)[source]¶
Calculate the temperature of water according to the IAPWS-95 standard given a density rho and pressure P.
- Parameters
- Returns
- T
float
Temperature, [K]
- T
Notes
This solution is iterative due to the nature of the equation. The solution procedure begins with IAPWS-97’s equations as an initial guess, extrapolating when out of range. Newton’s method converges extremely, normally after 2 or 3 iterations.
Due to water’s unique density curve, there is a temperature region spanning 273.15 K to 280.005 K where there are two solutions. No guarantee is made as to which solution will be returned.
References
- 1
Wagner, Wolfgang, and Andreas Pruß. “The IAPWS Formulation 1995 for the Thermodynamic Properties of Ordinary Water Substance for General and Scientific Use.” Journal of Physical and Chemical Reference Data 31, no. 2 (2002): 387-535.
Examples
>>> iapws95_T(P=1e6, rho=995.0) 306.461547194
IAPWS-97 Basic Solvers¶
- chemicals.iapws.iapws97_rho(T, P, use_95_boundary=False)[source]¶
Calculate the density of water in kg/m^3 according to the IAPWS-97 standard.
- Parameters
- Returns
- rho
float
Mass density of water, [kg/m^3]
- rho
Notes
The range of validity of this formulation is as follows:
For :
For :
A ValueError is raised if the temperature or the pressure is out of bounds.
IAPWS is implemented in four regions in the T-P domain: Region 1 (liquid), region 2 (gas and supercritical gas), region 5 (high temperature gas), and region 3 (near-critical). Significant discontinuities exist between the transitions of each regions. In region 3, there are 26 sub-regions and the correlation has the least accuracy.
For many applications, the discontinuities in IF-97 can be problematic and the slower IAPWS-95 must be used. IAPWS-95 also has a wider range of applicability.
References
- 1
Cooper, JR, and RB Dooley. “Revised Release on the IAPWS Industrial Formulation 1997 for the Thermodynamic Properties of Water and Steam.” The International Association for the Properties of Water and Steam 1 (2007): 48.
Examples
>>> iapws97_rho(648.6, 22.5e6) 353.06081088726 >>> iapws97_rho(330.0, 8e5) 985.10498080770 >>> iapws97_rho(823.0, 14e6) 40.39293607288123 >>> iapws97_rho(2000.0, 3e7) 32.11456228328856
- chemicals.iapws.iapws97_P(T, rho)[source]¶
Calculate the pressure of water according to the IAPWS-97 standard given a temperature T and mass density rho.
- Parameters
- Returns
- P
float
Pressure, [Pa]
- P
Notes
The range of validity of this formulation is as follows:
For :
For :
A ValueError is raised if the temperature or density is out of bounds.
Newton’s method with analytical derivatives is used here to solve these equations. The solver tolerance is as tight as it can be without causing wasted iterations that do not improve the result at all. Pressure changes quickly with density however, and some discrepancy between solvers is to be expected.
For region 3, there are really two formulations present in IAPWS-97. There is a Helmholtz energy equation (Temperature and density dependent), and also 26 separate backwards equations for rho which depend on T and P. The Helmholtz energy equation is much more accurate and does not have discontinuities. The two sets of equations agree closely not not perfectly. By design,
iapws97_rho
implements the 26 T-P equations and this implements the Helmholtz energy equation. This means that in region 3 solutions will not be consistent. For consistency requirements, IAPWS-95 is recommended.This solver does not have any issues with multiple solutions. The solvers have been checked to achieve a relative solution tolerance of 5e-9 on 100 million points.
References
- 1
Cooper, JR, and RB Dooley. “Revised Release on the IAPWS Industrial Formulation 1997 for the Thermodynamic Properties of Water and Steam.” The International Association for the Properties of Water and Steam 1 (2007): 48.
Examples
>>> iapws97_P(330.0, iapws97_rho(T=330.0, P=8e5)) 8e5 >>> iapws97_P(823.0, 40.39293607288123) 14e6 >>> iapws97_P(T=2000.0, rho=32.11456228328856) 3e7
Region 3 point - does not implement the same equations as
iapws97_rho
!>>> iapws97_P(648.6, iapws97_rho(T=648.6, P=22.5e6)) 22499974.093936257
- chemicals.iapws.iapws97_T(P, rho)[source]¶
Calculate the temperature of water according to the IAPWS-97 standard given a pressure P and mass density rho.
- Parameters
- Returns
- T
float
Temperature, [K]
- T
Notes
The range of validity of this formulation is as follows:
For :
For :
A ValueError is raised if the pressure or density is out of bounds.
Newton’s method with analytical derivatives is used here to solve these equations. The solver tolerance is as tight as it can be without causing wasted iterations that do not improve the result at all.
Due to water’s unique density curve, there is a temperature region spanning 273.15 K to 280.005 K where there are two solutions. No guarantee is made as to which solution will be returned.
References
- 1
Cooper, JR, and RB Dooley. “Revised Release on the IAPWS Industrial Formulation 1997 for the Thermodynamic Properties of Water and Steam.” The International Association for the Properties of Water and Steam 1 (2007): 48.
Examples
>>> iapws97_T(8e5, iapws97_rho(T=330.0, P=8e5)) 330.0 >>> iapws97_T(14e6, 40.39293607288123) 823.0 >>> iapws97_T(P=3e7, rho=32.11456228328856) 2000.0
IAPWS-95 Properties¶
- chemicals.iapws.iapws95_properties(T, P)[source]¶
Calculate some basic properties of water according to the IAPWS-95 standard given a temperature T and pressure P.
The properties are density rho, internal energy U, entropy S, enthalpy H, isochoric heat capacity Cv, isobaric heat capacity Cp, speed of sound w, Joule-Thomson coefficient JT, isothermal throttling coefficient delta_T, isentropic temperature-pressure coefficient beta_s, and the derivative of mass density with respect to pressure at constant temperature drho_dP.
This function is intended as a demonstration of how to use the IAPWS-95 equations. For that reason, mass-units are used in all returned variables.
- Parameters
- Returns
- rho
float
Mass density of water, [kg/m^3]
- U
float
Internal energy of water, [J/(kg)]
- S
float
Entropy of water, [J/(kg*K)]
- H
float
Enthalpy of water, [J/(kg)]
- Cv
float
Isochoric heat capacity, [J/(kg*K)]
- Cp
float
Isobaric heat capacity, [J/(kg*K)]
- w
float
Speed of sound, [m/s]
- JT
float
Joule-Thomson coefficient, [K/Pa]
- delta_T
float
Isothermal throttling coefficient, [J/(kg*Pa)]
- beta_s
float
Isentropic temperature-pressure coefficient, [K/Pa]
- drho_dP
float
Derivative of mass density with respect to pressure at constant temperature, [kg/(m^3*Pa)]
- rho
Notes
Hundreds of useful properties can be obtained from the IAPWS-95 model. It is intended for this function to serve as a useful starting point to those. Calculating every property with every set of units is beyond the scope of chemicals. The functions like
iapws95_dAr_ddelta
can be used directly in your own implementation - where you can calculate only those properties which are necessary, for maximum speed.The formulas are as follows:
This derivative isn’t part of the same table of properties, but it is needed by the transport calculation routines:
References
- 1
Wagner, Wolfgang, and Andreas Pruß. “The IAPWS Formulation 1995 for the Thermodynamic Properties of Ordinary Water Substance for General and Scientific Use.” Journal of Physical and Chemical Reference Data 31, no. 2 (2002): 387-535.
Examples
>>> iapws95_properties(T=300.0, P=1e6) (996.96002269, 112478.998245, 392.813902893, 113482.047492, 4127.21730497, 4178.103605593, 1503.035983829, -2.202166728257e-07, 0.000920088074745, 1.985617879134e-08, 4.48108429028e-07)
>>> rho, U, S, H, Cv, Cp, w, JT, delta_T, beta_s, drho_dP = iapws95_properties(T=500.0, P=1e5) >>> w 548.3138393244
IAPWS Saturation Pressure/Temperature¶
- chemicals.iapws.iapws95_Psat(T)[source]¶
Compute the saturation pressure of the IAPWS-95 equation using high- fidelity polynomial fits. These have a relative accuracy of under 1e-12, and are generated by solving the saturation equations under the high-precision environment of mpmath. The range of the fit is 235 K to 647.096 K, the critical point.
- Parameters
- T
float
Temperature at which to calculate the saturation condition, [K]
- T
- Returns
- Psat
float
Saturation vapor pressure, [Pa]
- Psat
See also
Notes
This method should be used in preference to
iapws95_saturation
. Although using mpmath generates slightly different results than using plain floating point numbers, the requirement for the saturation curve is to be smooth, and continuous; mpmath makes this easy and the saturation equations were solved extremely high precision, well under a floating point’s error.The polynomial coefficients have been carefully chosen to be able to be evaluated accurately with horner’s method, although they are derived as a Chebyshev approximation originally.
Examples
>>> iapws95_Psat(400.0) 245769.3455
- chemicals.iapws.iapws95_dPsat_dT(T)[source]¶
Compute the temperature derivative of saturation pressure of the IAPWS-95 equation using high- fidelity polynomial fits. The range of the fit is 235 K to 647.096 K, the critical point.
- Parameters
- T
float
Temperature at which to calculate the saturation condition and its temperature derivative, [K]
- T
- Returns
Notes
Psat must be calculated in the calculation of the derivative, so it is returned as well which may be useful in some applications.
Examples
>>> iapws95_dPsat_dT(400.0) (7483.62075827, 245769.3455657)
- chemicals.iapws.iapws92_Psat(T)[source]¶
Compute the saturation pressure of the IAPWS-92 equation.
- Parameters
- T
float
Temperature at which to calculate the saturation condition and its temperature derivative, [K]
- T
- Returns
- Psat
float
Saturation vapor pressure, [Pa]
- Psat
Notes
The coefficients are [-7.85951783, 1.84408259, -11.7866497, 22.6807411, -15.9618719, 1.80122502]
Examples
>>> iapws92_Psat(400.0) 245765.2635418
- chemicals.iapws.iapws92_dPsat_dT(T)[source]¶
Compute the temperature derivative of saturation pressure of the IAPWS-92 equation.
- Parameters
- T
float
Temperature at which to calculate the saturation condition and its temperature derivative, [K]
- T
- Returns
Notes
The coefficients are [-7.85951783, 1.84408259, -11.7866497, 22.6807411, -15.9618719, 1.80122502]
Examples
>>> iapws92_dPsat_dT(400.0) (7483.47094105, 245765.263541)
- chemicals.iapws.iapws95_Tsat(Psat)[source]¶
Compute the saturation temperature of the IAPWS-95 equation. The range of the fit is 235 K to 647.096 K, the critical point.
- Parameters
- Psat
float
Saturation vapor pressure specified, [Pa]
- Psat
- Returns
- T
float
Temperature at which the saturation pressure occurs, [K]
- T
See also
iapws95_Psat
Tsat_IAPWS
Notes
This method is quite fast and precise because it starts with great initial guesses and the equation is well-bounded. The precision of this calculation should be the same as
iapws95_Psat
.Examples
>>> iapws95_Tsat(iapws95_Psat(400.0)) 400.0
- chemicals.iapws.iapws95_saturation(T, xtol=1e-05, rhol_guess=None, rhog_guess=None)[source]¶
Solve the vapor-liquid saturation equations of IAPWS-95 given a specified temperature. With floating point numbers, the achievable tolerance is somewhat low so xtol is exposed as a setting - it can be adjusted somewhat. Density guesses may be provided, otherwise they will be estimated.
- Parameters
- Returns
Notes
This is not a perfect function.
With mpmath multiple precision, the equation can be solved down to 233.6 K and up to 647.095999995 K - within 10 parts in a billion of the critical point exactly.
Reasons for non-convergence include floating point issues as delta becomes 1, and zero division errors in the matrix inverse.
Examples
>>> iapws95_saturation(400.0, xtol=1e-6) (245769.345, 937.4860, 1.3694075) >>> iapws95_saturation(647.0955, xtol=1e-7) (22063866.35, 325.70, 318.277)
- chemicals.iapws.iapws11_Psub(T)[source]¶
Compute the sublimation pressure of the frozen water using the IAPWS-11 equation from the Revised Release on the Pressure along the Melting and Sublimation Curves of Ordinary Water Substance.
- Parameters
- T
float
Temperature at which to calculate the sublimation condition [K]
- T
- Returns
- Psub
float
Sublimation vapor pressure, [Pa]
- Psub
Notes
The triple temperature is 273.16 K, and triple pressure 611.657 Pa.
The coefficients are as follows:
ais = [-0.212144006E2, 0.273203819E2, -0.610598130E1]
bis = [0.333333333E-2, 0.120666667E1, 0.170333333E1]
The equation is valid from 50 K to the triple temperature.
Examples
>>> iapws11_Psub(230.0) 8.947352740189151
IAPWS Saturation Density¶
- chemicals.iapws.iapws95_rhol_sat(T)[source]¶
Compute the saturation liquid density of the IAPWS-95 equation using high- fidelity polynomial fits. These have a relative accuracy of under 1e-13, except near the critical point where it rises to 1e-10, and are generated by solving the saturation equations under the high-precision environment of mpmath. The range of the fit is 235 K to 647.096 K, the critical point.
- Parameters
- T
float
Temperature at which to calculate the saturation condition, [K]
- T
- Returns
- rhol
float
Saturation liquid density, [kg/m^3]
- rhol
See also
Notes
This method should be used in preference to
iapws92_rhol_sat
.Examples
>>> iapws95_rhol_sat(400.0) 937.48603939
- chemicals.iapws.iapws95_rhog_sat(T)[source]¶
Compute the saturation vapor density of the IAPWS-95 equation using high- fidelity polynomial fits. These have a relative accuracy of under 1e-13, except near the critical point where it rises to 1e-10, and are generated by solving the saturation equations under the high-precision environment of mpmath. The range of the fit is 235 K to 647.096 K, the critical point.
- Parameters
- T
float
Temperature at which to calculate the saturation condition, [K]
- T
- Returns
- rhol
float
Saturation vapor density, [kg/m^3]
- rhol
See also
Notes
This method should be used in preference to
iapws92_rhog_sat
.Examples
>>> iapws95_rhog_sat(400.0) 1.3694075410
- chemicals.iapws.iapws95_drhol_sat_dT(T)[source]¶
Compute the first temperature derivative of saturation liquid density of the IAPWS-95 equation using high-fidelity polynomial fits. The actual saturated liquid density is returned as well.
The range of the fit is 235 K to 647.096 K, the critical point.
- Parameters
- T
float
Temperature at which to calculate the saturation condition and its derivative, [K]
- T
- Returns
Examples
>>> iapws95_drhol_sat_dT(400.0) (-0.835194603380, 937.486039392)
- chemicals.iapws.iapws92_rhol_sat(T)[source]¶
Calculates saturation liquid mass density of water using the IAPWS SR1-86(1992) [1] [2] explicit equation.
- Parameters
- T
float
Temperature of water, [K]
- T
- Returns
- rhol_sat
float
Saturation liquid mass density of water [kg/m^3]
- rhol_sat
See also
Notes
This equation is fit to experimental data to within its accuracy. It does not satisfy the equilibrium conditions for the IAPWS-95 or IAPWS-97 formulations.
The values of the constants are as follows:
b1 = 1.99274064; b2 = 1.09965342; b3 = -0.510839303; b4 = -1.75493479; b5 = -45.5170352; b6 = -6.74694450e5
References
- 1
IAPWS, Secretariat, B Dooley, and EPRI. “Revised Supplementary Release on Saturation Properties of Ordinary Water Substance”, 1992.
- 2
Wagner, Wolfgang, and A. Pruss. “International Equations for the Saturation Properties of Ordinary Water Substance. Revised According to the International Temperature Scale of 1990. Addendum to J. Phys. Chem. Ref. Data 16, 893 (1987).” Journal of Physical and Chemical Reference Data 22, no. 3 (May 1, 1993): 783-87. https://doi.org/10.1063/1.555926.
Examples
>>> iapws92_rhol_sat(300.) 996.5089712803
- chemicals.iapws.iapws92_rhog_sat(T)[source]¶
Calculates saturation vapor mass density of water using the IAPWS SR1-86(1992) [1] [2] explicit equation.
- Parameters
- T
float
Temperature of water, [K]
- T
- Returns
- rhog_sat
float
Saturation vapor mass density of water [kg/m^3]
- rhog_sat
See also
Notes
This equation is fit to experimental data to within its accuracy. It does not satisfy the equilibrium conditions for the IAPWS-95 or IAPWS-97 formulations.
The values of the constants are as follows:
c1 = -2.03150240; c2 = -2.68302940; c3 = -5.38626492; c4 = -17.2991605; c5 = -44.7586581; c6 = -63.9201063
References
- 1
IAPWS, Secretariat, B Dooley, and EPRI. “Revised Supplementary Release on Saturation Properties of Ordinary Water Substance”, 1992.
- 2
Wagner, Wolfgang, and A. Pruss. “International Equations for the Saturation Properties of Ordinary Water Substance. Revised According to the International Temperature Scale of 1990. Addendum to J. Phys. Chem. Ref. Data 16, 893 (1987).” Journal of Physical and Chemical Reference Data 22, no. 3 (May 1, 1993): 783-87. https://doi.org/10.1063/1.555926.
Examples
>>> iapws92_rhog_sat(300.) 0.0255887212886
IAPWS Constants¶
- chemicals.iapws.iapws95_Tc = 647.096¶
Critical temperature of water in K according to IAPWS-95, also used in IAPWS-97
- chemicals.iapws.iapws95_Pc = 22064000.0¶
Critical pressure of water in Pa according to IAPWS-95, also used in IAPWS-97
- chemicals.iapws.iapws95_rhoc = 322.0¶
Critical density of water in kg/m^3 according to IAPWS-95, also used in IAPWS-97
- chemicals.iapws.iapws95_MW = 18.015268¶
Molecular weight of water in g/mol according to IAPWS-95, also used in IAPWS-97
- chemicals.iapws.iapws95_R = 461.51805¶
Specific gas constant in J/(kg*K) according to IAPWS-95
- chemicals.iapws.iapws97_R = 461.526¶
Specific gas constant in J/(kg*K) according to IAPWS-97
- chemicals.iapws.iapws95_Tt = 273.16¶
Triple temperature of water in K according to IAPWS
IAPWS-97 Region 1¶
- chemicals.iapws.iapws97_G_region1(tau, pi)[source]¶
Calculates the dimensionless Gibbs free energy for water according to the IAPWS-97 standard (for region 1).
- Parameters
- Returns
- G
float
Dimensionless Gibbs energy G/(RT), [-]
- G
Examples
>>> iapws97_G_region1(1386/277.15, 101325/16.53E6) -0.00016341033954414
- chemicals.iapws.iapws97_dG_dpi_region1(tau, pi)[source]¶
Calculates the derivative of dimensionless Gibbs free energy with respect to pi for water according to the IAPWS-97 standard (for region 1).
- Parameters
- Returns
- dG_dpi
float
Derivative of dimensionless Gibbs energy G/(RT) with respect to pi, [-]
- dG_dpi
Notes
Used in density solution. This contains a hand-optimized implementation with a single division, no power operations, 65 multiplications, 16 local variables, and a minimum number of additions.
Examples
>>> iapws97_dG_dpi_region1(1386/277.15, 101325/16.53E6) 0.1292327182544
- chemicals.iapws.iapws97_d2G_dpi2_region1(tau, pi)[source]¶
Calculates the second derivative of dimensionless Gibbs free energy with respect to pi for water according to the IAPWS-97 standard (for region 1).
- Parameters
- Returns
- d2G_dpi2
float
Second Derivative of dimensionless Gibbs energy G/(RT) with respect to pi, [-]
- d2G_dpi2
Examples
>>> iapws97_d2G_dpi2_region1(1386/277.15, 101325/16.53E6) -0.0010570100274769
- chemicals.iapws.iapws97_dG_dtau_region1(tau, pi)[source]¶
Calculates the derivative of dimensionless Gibbs free energy with respect to tau for water according to the IAPWS-97 standard (for region 1).
- Parameters
- Returns
- dG_dtau
float
Derivative of dimensionless Gibbs energy G/(RT) with respect to tau, [-]
- dG_dtau
Examples
>>> iapws97_dG_dtau_region1(1386/277.15, 101325/16.53E6) 0.026440334282967
- chemicals.iapws.iapws97_d2G_dtau2_region1(tau, pi)[source]¶
Calculates the second derivative of dimensionless Gibbs free energy with respect to tau for water according to the IAPWS-97 standard (for region 1).
- Parameters
- Returns
- d2G_dtau2
float
Second Derivative of dimensionless Gibbs energy G/(RT) with respect to tau, [-]
- d2G_dtau2
Examples
>>> iapws97_d2G_dtau2_region1(1386/277.15, 101325/16.53E6) -0.3645169808573
- chemicals.iapws.iapws97_d2G_dpidtau_region1(tau, pi)[source]¶
Calculates the second derivative of dimensionless Gibbs free energy with respect to tau and pi for water according to the IAPWS-97 standard (for region 1).
- Parameters
- Returns
- d2G_dpidtau
float
Second Derivative of dimensionless Gibbs energy G/(RT) with respect to tau and pi, [-]
- d2G_dpidtau
Examples
>>> iapws97_d2G_dpidtau_region1(1386/277.15, 101325/16.53E6) 0.025837659858819
IAPWS-97 Region 2¶
- chemicals.iapws.iapws97_G0_region2(tau, pi)[source]¶
Calculates the dimensionless ideal gas Gibbs free energy for water according to the IAPWS-97 standard (for region 2).
- Parameters
- Returns
- G0
float
Dimensionless ideal gas Gibbs energy G0/(RT), [-]
- G0
Examples
>>> iapws97_G0_region2(540/300.0, 101325/1e6) 3.3180953922351
- chemicals.iapws.iapws97_dG0_dtau_region2(tau, pi)[source]¶
Calculates the first derivative of dimensionless ideal gas Gibbs free energy with respect to tau for water according to the IAPWS-97 standard (for region 2).
- Parameters
- Returns
- dG0_dtau
float
First derivative of dimensionless ideal gas Gibbs energy G0/(RT) with respect to tau, [-]
- dG0_dtau
Notes
This function does not depend on pi but it is accepted for consistency.
Examples
>>> iapws97_dG0_dtau_region2(540/300.0, 101325/1e6) 10.2374188173906
- chemicals.iapws.iapws97_d2G0_dtau2_region2(tau, pi)[source]¶
Calculates the second derivative of dimensionless ideal gas Gibbs free energy with respect to tau for water according to the IAPWS-97 standard (for region 2).
- Parameters
- Returns
- d2G0_dtau2
float
Second derivative of dimensionless ideal gas Gibbs energy G0/(RT) with respect to tau, [-]
- d2G0_dtau2
Notes
This function does not depend on pi but it is accepted for consistency.
Examples
>>> iapws97_d2G0_dtau2_region2(540/300.0, 101325/1e6) -1.2472096479372
- chemicals.iapws.iapws97_Gr_region2(tau, pi)[source]¶
Calculates the dimensionless residual Gibbs free energy for water according to the IAPWS-97 standard (for region 2).
- Parameters
- Returns
- Gr
float
Dimensionless residual Gibbs energy Gr/(RT), [-]
- Gr
Examples
>>> iapws97_Gr_region2(540/300.0, 101325/1e6) -0.71851548053980
- chemicals.iapws.iapws97_dGr_dpi_region2(tau, pi)[source]¶
Calculates the first derivative of dimensionless residual Gibbs free energy with respect to pi for water according to the IAPWS-97 standard (for region 2).
- Parameters
- Returns
- dGr_dpi
float
Derivative of dimensionless residual Gibbs energy Gr/(RT) with respect to pi, [-]
- dGr_dpi
Notes
Used in density solution.
Examples
>>> iapws97_dGr_dpi_region2(540/300.0, 101325/1e6) -27.7714056629532
- chemicals.iapws.iapws97_d2Gr_dpi2_region2(tau, pi)[source]¶
Calculates the second derivative of dimensionless residual Gibbs free energy with respect to pi for water according to the IAPWS-97 standard (for region 2).
- Parameters
- Returns
- d2Gr_dpi2
float
Second Derivative of dimensionless residual Gibbs energy Gr/(RT) with respect to pi, [-]
- d2Gr_dpi2
Examples
>>> iapws97_d2Gr_dpi2_region2(540/300.0, 101325/1e6) -983.15187604898
- chemicals.iapws.iapws97_dGr_dtau_region2(tau, pi)[source]¶
Calculates the first derivative of dimensionless residual Gibbs free energy with respect to tau for water according to the IAPWS-97 standard (for region 2).
- Parameters
- Returns
- dGr_dtau
float
Derivative of dimensionless residual Gibbs energy Gr/(RT) with respect to tau, [-]
- dGr_dtau
Examples
>>> iapws97_dGr_dtau_region2(540/300.0, 101325/1e6) -18.1535856049444
- chemicals.iapws.iapws97_d2Gr_dtau2_region2(tau, pi)[source]¶
Calculates the second derivative of dimensionless residual Gibbs free energy with respect to tau for water according to the IAPWS-97 standard (for region 2).
- Parameters
- Returns
- d2Gr_dtau2
float
Second derivative of dimensionless residual Gibbs energy Gr/(RT) with respect to tau, [-]
- d2Gr_dtau2
Examples
>>> iapws97_d2Gr_dtau2_region2(540/300.0, 101325/1e6) -470.9302933324787
- chemicals.iapws.iapws97_d2Gr_dpidtau_region2(tau, pi)[source]¶
Calculates the second derivative of dimensionless residual Gibbs free energy with respect to tau and pi for water according to the IAPWS-97 standard (for region 2).
- Parameters
- Returns
- d2Gr_dpidtau_
float
Second derivative of dimensionless residual Gibbs energy Gr/(RT) with respect to tau and pi, [-]
- d2Gr_dpidtau_
Examples
>>> iapws97_d2Gr_dpidtau_region2(540/300.0, 101325/1e6) -735.391845360247
IAPWS-97 Region 3¶
- chemicals.iapws.iapws97_A_region3(tau, delta)[source]¶
Calculates the dimensionless Helmholtz free energy for water according to the IAPWS-97 standard (for region 3).
- Parameters
- Returns
- A
float
Helmholtz free energy A/(RT), [-]
- A
Examples
>>> iapws97_A_region3(647.096/500.0, 400.0/322.0) -3.0336402168865
- chemicals.iapws.iapws97_dA_ddelta_region3(tau, delta)[source]¶
Calculates the derivative of dimensionless Helmholtz free energy with respect to delta for water according to the IAPWS-97 standard (for region 3).
- Parameters
- Returns
- dA_ddelta
float
Derivative of dimensionless Helmholtz free energy with respect to delta, [-]
- dA_ddelta
Examples
>>> iapws97_dA_ddelta_region3(647.096/500.0, 400.0/322.0) 7.35562435092
- chemicals.iapws.iapws97_d2A_ddelta2_region3(tau, delta)[source]¶
Calculates the second derivative of dimensionless Helmholtz free energy with respect to delta for water according to the IAPWS-97 standard (for region 3).
- Parameters
- Returns
- d2A_ddelta2
float
Second derivative of dimensionless Helmholtz free energy with respect to delta, [-]
- d2A_ddelta2
Examples
>>> iapws97_d2A_ddelta2_region3(647.096/500.0, 400.0/322.0) -2.2858869882497
- chemicals.iapws.iapws97_dA_dtau_region3(tau, delta)[source]¶
Calculates the derivative of dimensionless Helmholtz free energy with respect to tau for water according to the IAPWS-97 standard (for region 3).
- Parameters
- Returns
- dA_dtau
float
Derivative of dimensionless Helmholtz free energy with respect to tau, [-]
- dA_dtau
Examples
>>> iapws97_dA_dtau_region3(647.096/500.0, 400.0/322.0) -24.9687028688
- chemicals.iapws.iapws97_d2A_dtau2_region3(tau, delta)[source]¶
Calculates the second derivative of dimensionless Helmholtz free energy with respect to tau for water according to the IAPWS-97 standard (for region 3).
- Parameters
- Returns
- d2A_dtau2
float
Second derivative of dimensionless Helmholtz free energy with respect to tau, [-]
- d2A_dtau2
Examples
>>> iapws97_d2A_dtau2_region3(647.096/500.0, 400.0/322.0) -373.6565823701
- chemicals.iapws.iapws97_d2A_ddeltadtau_region3(tau, delta)[source]¶
Calculates the second derivative of dimensionless Helmholtz free energy with respect to tau and delta for water according to the IAPWS-97 standard (for region 3).
- Parameters
- Returns
- d2A_ddeltadtau
float
Second derivative of dimensionless Helmholtz free energy with respect to tau and delta, [-]
- d2A_ddeltadtau
Examples
>>> iapws97_d2A_ddeltadtau_region3(647.096/500.0, 400.0/322.0) 145.85190014717
IAPWS-97 Region 3 PT Backwards Equation Boundaries¶
- chemicals.iapws.iapws97_boundary_3uv(P)[source]¶
Calculates the transition temperature for a region 3 PT backwards equation transition.
Examples
>>> iapws97_boundary_3uv(22.3E6) 647.7996121480069
- chemicals.iapws.iapws97_boundary_3ef(P)[source]¶
Calculates the transition temperature for a region 3 PT backwards equation transition.
Examples
>>> iapws97_boundary_3ef(40E6) 713.959399239744
- chemicals.iapws.iapws97_boundary_3cd(P)[source]¶
Calculates the transition temperature for a region 3 PT backwards equation transition.
Examples
>>> iapws97_boundary_3cd(25E6) 649.3659208321279
- chemicals.iapws.iapws97_boundary_3gh(P)[source]¶
Calculates the transition temperature for a region 3 PT backwards equation transition.
Examples
>>> iapws97_boundary_3gh(25E6) 656.69805722612
- chemicals.iapws.iapws97_boundary_3ij(P)[source]¶
Calculates the transition temperature for a region 3 PT backwards equation transition.
Examples
>>> iapws97_boundary_3ij(25E6) 660.7865756716819
- chemicals.iapws.iapws97_boundary_3jk(P)[source]¶
Calculates the transition temperature for a region 3 PT backwards equation transition.
Examples
>>> iapws97_boundary_3jk(25E6) 668.1915358826951
- chemicals.iapws.iapws97_boundary_3mn(P)[source]¶
Calculates the transition temperature for a region 3 PT backwards equation transition.
Examples
>>> iapws97_boundary_3mn(22.8E6) 649.6054132953997
- chemicals.iapws.iapws97_boundary_3qu(P)[source]¶
Calculates the transition temperature for a region 3 PT backwards equation transition.
Examples
>>> iapws97_boundary_3qu(22E6) 645.6355027340121
- chemicals.iapws.iapws97_boundary_3rx(P)[source]¶
Calculates the transition temperature for a region 3 PT backwards equation transition.
Examples
>>> iapws97_boundary_3rx(22E6) 648.26227536701
- chemicals.iapws.iapws97_boundary_3wx(logP_MPa, logP_MPa_inv)[source]¶
Calculates the transition temperature for a region 3 PT backwards equation transition (for one of “wx”, “ab”, or “op”; the others do not use a log fit). The parameters are provided in the specific units for speed savings only.
- Parameters
- Returns
- T_trans
float
Transition temperature [K]
- T_trans
Examples
>>> iapws97_boundary_3wx(log(22.3), 1/log(22.3)) 648.204947950734
- chemicals.iapws.iapws97_boundary_3ab(logP_MPa, logP_MPa_inv)[source]¶
Calculates the transition temperature for a region 3 PT backwards equation transition (for one of “wx”, “ab”, or “op”; the others do not use a log fit). The parameters are provided in the specific units for speed savings only.
- Parameters
- Returns
- T_trans
float
Transition temperature [K]
- T_trans
Examples
>>> iapws97_boundary_3ab(log(40), 1/log(40)) 693.0341408296053
- chemicals.iapws.iapws97_boundary_3op(logP_MPa, logP_MPa_inv)[source]¶
Calculates the transition temperature for a region 3 PT backwards equation transition (for one of “wx”, “ab”, or “op”; the others do not use a log fit). The parameters are provided in the specific units for speed savings only.
- Parameters
- Returns
- T_trans
float
Transition temperature [K]
- T_trans
Examples
>>> iapws97_boundary_3op(log(22.8), 1/log(22.8)) 650.010694314133
IAPWS-97 Region 3 PT Backwards Equations¶
- chemicals.iapws.iapws97_region3_a(T, P)[source]¶
Calculate the mass density water in one of the 26 region 3 backwards regions of the IAPWS-97 standard.
- Parameters
- Returns
- rho
float
Mass density of water in region 3, [kg/m^3]
- rho
Notes
Significant discontinuities exist between each region. These functions are automatically generated and are not to be edited directly.
- chemicals.iapws.iapws97_region3_b(T, P)[source]¶
Calculate the mass density water in one of the 26 region 3 backwards regions of the IAPWS-97 standard.
- Parameters
- Returns
- rho
float
Mass density of water in region 3, [kg/m^3]
- rho
Notes
Significant discontinuities exist between each region. These functions are automatically generated and are not to be edited directly.
- chemicals.iapws.iapws97_region3_c(T, P)[source]¶
Calculate the mass density water in one of the 26 region 3 backwards regions of the IAPWS-97 standard.
- Parameters
- Returns
- rho
float
Mass density of water in region 3, [kg/m^3]
- rho
Notes
Significant discontinuities exist between each region. These functions are automatically generated and are not to be edited directly.
- chemicals.iapws.iapws97_region3_d(T, P)[source]¶
Calculate the mass density water in one of the 26 region 3 backwards regions of the IAPWS-97 standard.
- Parameters
- Returns
- rho
float
Mass density of water in region 3, [kg/m^3]
- rho
Notes
Significant discontinuities exist between each region. These functions are automatically generated and are not to be edited directly.
- chemicals.iapws.iapws97_region3_e(T, P)[source]¶
Calculate the mass density water in one of the 26 region 3 backwards regions of the IAPWS-97 standard.
- Parameters
- Returns
- rho
float
Mass density of water in region 3, [kg/m^3]
- rho
Notes
Significant discontinuities exist between each region. These functions are automatically generated and are not to be edited directly.
- chemicals.iapws.iapws97_region3_f(T, P)[source]¶
Calculate the mass density water in one of the 26 region 3 backwards regions of the IAPWS-97 standard.
- Parameters
- Returns
- rho
float
Mass density of water in region 3, [kg/m^3]
- rho
Notes
Significant discontinuities exist between each region. These functions are automatically generated and are not to be edited directly.
- chemicals.iapws.iapws97_region3_g(T, P)[source]¶
Calculate the mass density water in one of the 26 region 3 backwards regions of the IAPWS-97 standard.
- Parameters
- Returns
- rho
float
Mass density of water in region 3, [kg/m^3]
- rho
Notes
Significant discontinuities exist between each region. These functions are automatically generated and are not to be edited directly.
- chemicals.iapws.iapws97_region3_h(T, P)[source]¶
Calculate the mass density water in one of the 26 region 3 backwards regions of the IAPWS-97 standard.
- Parameters
- Returns
- rho
float
Mass density of water in region 3, [kg/m^3]
- rho
Notes
Significant discontinuities exist between each region. These functions are automatically generated and are not to be edited directly.
- chemicals.iapws.iapws97_region3_i(T, P)[source]¶
Calculate the mass density water in one of the 26 region 3 backwards regions of the IAPWS-97 standard.
- Parameters
- Returns
- rho
float
Mass density of water in region 3, [kg/m^3]
- rho
Notes
Significant discontinuities exist between each region. These functions are automatically generated and are not to be edited directly.
- chemicals.iapws.iapws97_region3_j(T, P)[source]¶
Calculate the mass density water in one of the 26 region 3 backwards regions of the IAPWS-97 standard.
- Parameters
- Returns
- rho
float
Mass density of water in region 3, [kg/m^3]
- rho
Notes
Significant discontinuities exist between each region. These functions are automatically generated and are not to be edited directly.
- chemicals.iapws.iapws97_region3_k(T, P)[source]¶
Calculate the mass density water in one of the 26 region 3 backwards regions of the IAPWS-97 standard.
- Parameters
- Returns
- rho
float
Mass density of water in region 3, [kg/m^3]
- rho
Notes
Significant discontinuities exist between each region. These functions are automatically generated and are not to be edited directly.
- chemicals.iapws.iapws97_region3_l(T, P)[source]¶
Calculate the mass density water in one of the 26 region 3 backwards regions of the IAPWS-97 standard.
- Parameters
- Returns
- rho
float
Mass density of water in region 3, [kg/m^3]
- rho
Notes
Significant discontinuities exist between each region. These functions are automatically generated and are not to be edited directly.
- chemicals.iapws.iapws97_region3_m(T, P)[source]¶
Calculate the mass density water in one of the 26 region 3 backwards regions of the IAPWS-97 standard.
- Parameters
- Returns
- rho
float
Mass density of water in region 3, [kg/m^3]
- rho
Notes
Significant discontinuities exist between each region. These functions are automatically generated and are not to be edited directly.
- chemicals.iapws.iapws97_region3_n(T, P)[source]¶
Calculate the mass density water in one of the 26 region 3 backwards regions of the IAPWS-97 standard.
- Parameters
- Returns
- rho
float
Mass density of water in region 3, [kg/m^3]
- rho
Notes
Significant discontinuities exist between each region. These functions are automatically generated and are not to be edited directly.
- chemicals.iapws.iapws97_region3_o(T, P)[source]¶
Calculate the mass density water in one of the 26 region 3 backwards regions of the IAPWS-97 standard.
- Parameters
- Returns
- rho
float
Mass density of water in region 3, [kg/m^3]
- rho
Notes
Significant discontinuities exist between each region. These functions are automatically generated and are not to be edited directly.
- chemicals.iapws.iapws97_region3_p(T, P)[source]¶
Calculate the mass density water in one of the 26 region 3 backwards regions of the IAPWS-97 standard.
- Parameters
- Returns
- rho
float
Mass density of water in region 3, [kg/m^3]
- rho
Notes
Significant discontinuities exist between each region. These functions are automatically generated and are not to be edited directly.
- chemicals.iapws.iapws97_region3_q(T, P)[source]¶
Calculate the mass density water in one of the 26 region 3 backwards regions of the IAPWS-97 standard.
- Parameters
- Returns
- rho
float
Mass density of water in region 3, [kg/m^3]
- rho
Notes
Significant discontinuities exist between each region. These functions are automatically generated and are not to be edited directly.
- chemicals.iapws.iapws97_region3_r(T, P)[source]¶
Calculate the mass density water in one of the 26 region 3 backwards regions of the IAPWS-97 standard.
- Parameters
- Returns
- rho
float
Mass density of water in region 3, [kg/m^3]
- rho
Notes
Significant discontinuities exist between each region. These functions are automatically generated and are not to be edited directly.
- chemicals.iapws.iapws97_region3_s(T, P)[source]¶
Calculate the mass density water in one of the 26 region 3 backwards regions of the IAPWS-97 standard.
- Parameters
- Returns
- rho
float
Mass density of water in region 3, [kg/m^3]
- rho
Notes
Significant discontinuities exist between each region. These functions are automatically generated and are not to be edited directly.
- chemicals.iapws.iapws97_region3_t(T, P)[source]¶
Calculate the mass density water in one of the 26 region 3 backwards regions of the IAPWS-97 standard.
- Parameters
- Returns
- rho
float
Mass density of water in region 3, [kg/m^3]
- rho
Notes
Significant discontinuities exist between each region. These functions are automatically generated and are not to be edited directly.
- chemicals.iapws.iapws97_region3_u(T, P)[source]¶
Calculate the mass density water in one of the 26 region 3 backwards regions of the IAPWS-97 standard.
- Parameters
- Returns
- rho
float
Mass density of water in region 3, [kg/m^3]
- rho
Notes
Significant discontinuities exist between each region. These functions are automatically generated and are not to be edited directly.
- chemicals.iapws.iapws97_region3_v(T, P)[source]¶
Calculate the mass density water in one of the 26 region 3 backwards regions of the IAPWS-97 standard.
- Parameters
- Returns
- rho
float
Mass density of water in region 3, [kg/m^3]
- rho
Notes
Significant discontinuities exist between each region. These functions are automatically generated and are not to be edited directly.
- chemicals.iapws.iapws97_region3_w(T, P)[source]¶
Calculate the mass density water in one of the 26 region 3 backwards regions of the IAPWS-97 standard.
- Parameters
- Returns
- rho
float
Mass density of water in region 3, [kg/m^3]
- rho
Notes
Significant discontinuities exist between each region. These functions are automatically generated and are not to be edited directly.
- chemicals.iapws.iapws97_region3_x(T, P)[source]¶
Calculate the mass density water in one of the 26 region 3 backwards regions of the IAPWS-97 standard.
- Parameters
- Returns
- rho
float
Mass density of water in region 3, [kg/m^3]
- rho
Notes
Significant discontinuities exist between each region. These functions are automatically generated and are not to be edited directly.
- chemicals.iapws.iapws97_region3_y(T, P)[source]¶
Calculate the mass density water in one of the 26 region 3 backwards regions of the IAPWS-97 standard.
- Parameters
- Returns
- rho
float
Mass density of water in region 3, [kg/m^3]
- rho
Notes
Significant discontinuities exist between each region. These functions are automatically generated and are not to be edited directly.
- chemicals.iapws.iapws97_region3_z(T, P)[source]¶
Calculate the mass density water in one of the 26 region 3 backwards regions of the IAPWS-97 standard.
- Parameters
- Returns
- rho
float
Mass density of water in region 3, [kg/m^3]
- rho
Notes
Significant discontinuities exist between each region. These functions are automatically generated and are not to be edited directly.
IAPWS-97 Region 5¶
- chemicals.iapws.iapws97_G0_region5(tau, pi)[source]¶
Calculates the dimensionless ideal gas Gibbs free energy for water according to the IAPWS-97 standard (for region 5).
- Parameters
- Returns
- G0
float
Dimensionless ideal gas Gibbs energy G/(RT), [-]
- G0
Examples
>>> iapws97_G0_region5(1000.0/1500, 101325/1e6) -14.9741430290056
- chemicals.iapws.iapws97_dG0_dtau_region5(tau, pi)[source]¶
Calculates the first derivative of dimensionless ideal gas Gibbs free energy with respect to tau for water according to the IAPWS-97 standard (for region 5).
- Parameters
- Returns
- dG0_dtau
float
First derivative of dimensionless ideal gas Gibbs energy G/(RT) with respect to tau, [-]
- dG0_dtau
Notes
This function does not depend on pi but it is accepted for consistency.
Examples
>>> iapws97_dG0_dtau_region5(1000.0/1500, 101325/1e6) 11.311766995978
- chemicals.iapws.iapws97_d2G0_dtau2_region5(tau, pi)[source]¶
Calculates the second derivative of dimensionless ideal gas Gibbs free energy with respect to tau for water according to the IAPWS-97 standard (for region 5).
- Parameters
- Returns
- d2G0_dtau2
float
Second derivative of dimensionless ideal gas Gibbs energy G/(RT) with respect to tau, [-]
- d2G0_dtau2
Notes
This function does not depend on pi but it is accepted for consistency.
Examples
>>> iapws97_d2G0_dtau2_region5(1000.0/1500, 101325/1e6) -12.744650271463655
- chemicals.iapws.iapws97_Gr_region5(tau, pi)[source]¶
Calculates the dimensionless residual Gibbs free energy for water according to the IAPWS-97 standard (for region 5).
- Parameters
- Returns
- Gr
float
Dimensionless residual Gibbs energy Gr/(RT), [-]
- Gr
Examples
>>> iapws97_Gr_region5(1000/300.0, 101325/1e6) -0.0194648291645718
- chemicals.iapws.iapws97_dGr_dpi_region5(tau, pi)[source]¶
Calculates the first derivative of dimensionless residual Gibbs free energy with respect to pi for water according to the IAPWS-97 standard (for region 5).
- Parameters
- Returns
- dGr_dpi
float
Derivative of dimensionless residual Gibbs energy Gr/(RT) with respect to pi, [-]
- dGr_dpi
Notes
Used in density solution.
Examples
>>> iapws97_dGr_dpi_region5(1000/300.0, 101325/1e6) -0.213281155629998
- chemicals.iapws.iapws97_d2Gr_dpi2_region5(tau, pi)[source]¶
Calculates the second derivative of dimensionless residual Gibbs free energy with respect to pi for water according to the IAPWS-97 standard (for region 5).
- Parameters
- Returns
- d2Gr_dpi2
float
Second derivative of dimensionless residual Gibbs energy Gr/(RT) with respect to pi, [-]
- d2Gr_dpi2
Examples
>>> iapws97_d2Gr_dpi2_region5(1000/300.0, 101325/1e6) -0.4179905782304291
- chemicals.iapws.iapws97_dGr_dtau_region5(tau, pi)[source]¶
Calculates the first derivative of dimensionless residual Gibbs free energy with respect to tau for water according to the IAPWS-97 standard (for region 5).
- Parameters
- Returns
- dGr_dtau
float
Derivative of dimensionless residual Gibbs energy Gr/(RT) with respect to tau, [-]
- dGr_dtau
Examples
>>> iapws97_dGr_dtau_region5(1000/300.0, 101325/1e6) -0.02200629869194
- chemicals.iapws.iapws97_d2Gr_dtau2_region5(tau, pi)[source]¶
Calculates the second derivative of dimensionless residual Gibbs free energy with respect to tau for water according to the IAPWS-97 standard (for region 5).
- Parameters
- Returns
- d2Gr_dtau2
float
Second derivative of dimensionless residual Gibbs energy Gr/(RT) with respect to tau, [-]
- d2Gr_dtau2
Examples
>>> iapws97_d2Gr_dtau2_region5(1000/300.0, 101325/1e6) -0.0239165867999155
- chemicals.iapws.iapws97_d2Gr_dpidtau_region5(tau, pi)[source]¶
Calculates the second derivative of dimensionless residual Gibbs free energy with respect to tau and pi for water according to the IAPWS-97 standard (for region 5).
- Parameters
- Returns
- d2Gr_dpidtau
float
Second derivative of dimensionless residual Gibbs energy Gr/(RT) with respect to tau and pi, [-]
- d2Gr_dpidtau
Examples
>>> iapws97_d2Gr_dpidtau_region5(1000/300.0, 101325/1e6) -0.27438379131103097
IAPWS-95 Ideal Gas Terms¶
- chemicals.iapws.iapws95_A0(tau, delta)[source]¶
Calculates the ideal gas Helmholtz energy of water according to the IAPWS-95 standard.
- Parameters
- Returns
- A0
float
Ideal gas dimensionless Helmholtz energy A/(RT) [-]
- A0
Notes
This implementation is checked to have a relative error always under 1e-15.
Examples
>>> iapws95_A0(647.096/300.0, 999.0/322) 9.537075529761053
- chemicals.iapws.iapws95_dA0_dtau(tau, delta)[source]¶
Calculates the first derivative of ideal gas Helmholtz energy of water with respect to tau according to the IAPWS-95 standard.
- Parameters
- Returns
- dA0_dtau
float
First derivative of ideal gas dimensionless Helmholtz energy A/(RT) with respect to tau [-]
- dA0_dtau
Notes
This implementation is checked to have a relative error always under 1e-15.
Examples
>>> iapws95_dA0_dtau(647.096/300.0, 999.0/322) 8.079705548882
- chemicals.iapws.iapws95_d2A0_dtau2(tau, delta)[source]¶
Calculates the second derivative of ideal gas Helmholtz energy of water with respect to tau according to the IAPWS-95 standard.
- Parameters
- Returns
- d2A0_dtau2
float
Second derivative of ideal gas dimensionless Helmholtz energy A/(RT) with respect to tau [-]
- d2A0_dtau2
Notes
This implementation is checked to have a relative error always under 1e-15.
Examples
>>> iapws95_d2A0_dtau2(647.096/300.0, 999.0/322) -0.653543047751809
- chemicals.iapws.iapws95_d3A0_dtau3(tau, delta)[source]¶
Calculates the third derivative of ideal gas Helmholtz energy of water with respect to tau according to the IAPWS-95 standard.
- Parameters
- Returns
- d3A0_dtau3
float
Third derivative of ideal gas dimensionless Helmholtz energy A/(RT) with respect to tau [-]
- d3A0_dtau3
Notes
This implementation is checked to have a relative error always under 1e-15. This equation is not explicitly in IAPWS-95, but is needed to compute some second derivatives.
Examples
>>> iapws95_d3A0_dtau3(647.096/300.0, 999.0/322) 0.6222542507278
- chemicals.iapws.iapws95_A0_tau_derivatives(tau, delta)[source]¶
Calculates the ideal gas Helmholtz energy of water and its first three derivatives with respect to tau according to the IAPWS-95 standard. As each of those calls spends most of their time computing exponentials which are the same for each function, function offers a time saving.
- Parameters
- Returns
- A0
float
Ideal gas dimensionless Helmholtz energy A/(RT) [-]
- dA0_dtau
float
First derivative of ideal gas dimensionless Helmholtz energy A/(RT) with respect to tau [-]
- d2A0_dtau2
float
Second derivative of ideal gas dimensionless Helmholtz energy A/(RT) with respect to tau [-]
- d3A0_dtau3
float
Third derivative of ideal gas dimensionless Helmholtz energy A/(RT) with respect to tau [-]
- A0
Notes
The extra cost of calling this function vs
iapws95_A0
alone is ~15% with numba, ~40% with PyPy, and 120% with CPython.Examples
>>> iapws95_A0_tau_derivatives(647.096/300.0, 999.0/322) (9.53707552976, 8.0797055488, -0.65354304775, 0.62225425072)
IAPWS-95 Residual Terms¶
- chemicals.iapws.iapws95_Ar(tau, delta)[source]¶
Calculates the residual Helmholtz energy of water according to the IAPWS-95 standard.
- Parameters
- Returns
- Ar
float
Residual Helmholtz energy A/(RT) [-]
- Ar
Notes
This is an optimized implementatation taking 9 exp calls, 4 sqrts, and 3 powers. It was generated using SymPy’s CSE functionality, with select polynomial optimizations by hand as well. It is over 10x faster than a naive implementation.
This implementation has been tested against a straightforward implementation with the equations given in IAPWS-95.
Over a linear temperature range of 200 K to 5000 K and a logarithmic density range of 1E-10 kg/m^3 to 5000 kg/m^3, 4E6 points were evaluated. The mean relative error was 5.0416E-15, with a maximum relative error of 1.118E-9 and a standard deviation of 5.773e-13.
Over the same range, the model was evaluated to a precision of 50 decimal places with mpmath, and on 90000 points, the mean relative error was 3.14E-15, with a maximum relative error of 3.54e-12 and a standard deviation of 3.017E-14.
This comparison indicates that this implementation is more accurate than the straightforward implementation.
Examples
>>> iapws95_Ar(647.096/300.0, 999.0/322) -9.57577716026768
- chemicals.iapws.iapws95_dAr_ddelta(tau, delta)[source]¶
Calculates the first derivative of residual Helmholtz energy of water with respect to delta according to the IAPWS-95 standard.
- Parameters
- Returns
- dAr_ddelta
float
First derivative of residual Helmholtz energy A/(RT) with respect to delta, [-]
- dAr_ddelta
Notes
This is an optimized implementatation taking 8 exp calls, 4 sqrts, and 2 powers. It was generated using SymPy’s CSE functionality, with select polynomial optimizations by hand as well. It is over 10x faster than a naive implementation.
This implementation has been tested against a straightforward implementation with the equations given in IAPWS-95.
Over a linear temperature range of 200 K to 5000 K and a logarithmic density range of 1E-10 kg/m^3 to 5000 kg/m^3, 4E6 points were evaluated. The mean relative error was 4.033E-15, with a maximum relative error of 3.8765e-10 and a standard deviation of 3.189e-13.
Over the same range, the model was evaluated to a precision of 50 decimal places with mpmath, and on 90000 points, the mean relative error was 6.046E-15, with a maximum relative error of 3.39E-10 and a standard deviation of 7.056E-13.
There was a singularity at tau = delta = 1, but the limit is correctly returned.
Examples
>>> iapws95_dAr_ddelta(647.096/300.0, 999.0/322) -0.3093321202374
- chemicals.iapws.iapws95_d2Ar_ddelta2(tau, delta)[source]¶
Calculates the second derivative of residual Helmholtz energy of water with respect to delta according to the IAPWS-95 standard.
- Parameters
- Returns
- d2Ar_ddelta2
float
Second derivative of residual Helmholtz energy A/(RT) with respect to delta, [-]
- d2Ar_ddelta2
Notes
This is an optimized implementatation taking 4 exp calls, 4 sqrts, and 2 powers. It was generated using SymPy’s CSE functionality, with select polynomial optimizations by hand as well. It is over 10x faster than a naive implementation.
This implementation has been tested against a straightforward implementation with the equations given in IAPWS-95.
Over a linear temperature range of 200 K to 5000 K and a logarithmic density range of 1E-10 kg/m^3 to 5000 kg/m^3, 4E6 points were evaluated. The mean relative error was 9.566e-16, with a maximum relative error of 1.0518E-10 and a standard deviation of 6.20265E-14.
Over the same range, the model was evaluated to a precision of 50 decimal places with mpmath, and on 250000 points, the mean relative error was 1.039E-15, with a maximum relative error of 2.431E-11 and a standard deviation of 5.31708E-14.
Examples
>>> iapws95_d2Ar_ddelta2(647.096/300.0, 999.0/322) 1.7862535141735987
- chemicals.iapws.iapws95_d3Ar_ddelta3(tau, delta)[source]¶
Calculates the third derivative of residual Helmholtz energy of water with respect to delta according to the IAPWS-95 standard.
- Parameters
- Returns
- d3Ar_ddelta3
float
Third derivative of residual Helmholtz energy A/(RT) with respect to delta, [-]
- d3Ar_ddelta3
Notes
No equation is given for this in IAPWS-95, and the derivative was symbolically computed with SymPy.
This is an optimized implementatation. It was generated using SymPy’s CSE functionality.
Over a linear temperature range of 200 K to 5000 K and a logarithmic density range of 1E-4 kg/m^3 to 5000 kg/m^3, 90000 points were evaluated. The mean relative error was 5.41E-13, with a maximum relative error of 6.3957e-11 and a standard deviation of 3.346e-12.
90000 points were also evaluated with mpmath. The mean relative error was 1.41959E-14, with a maximum relative error of 5.8878E-10 and a standard deviation of 1.978E-12.
Evaluating 10000 points in the 1e-10 to 1e-4 range, the mean relative error was 1.2E-16, maximum relative error 1.2e-16, and standard deviation 6.66e-16.
Examples
>>> iapws95_d3Ar_ddelta3(647.096/300.0, 999.0/322) 0.33621190578
- chemicals.iapws.iapws95_dAr_dtau(tau, delta)[source]¶
Calculates the first derivative of residual Helmholtz energy of water with respect to tau according to the IAPWS-95 standard.
- Parameters
- Returns
- dAr_dtau
float
Derivative of residual Helmholtz energy A/(RT) with respect to tau, [-]
- dAr_dtau
Notes
This is an optimized implementatation taking 9 exp calls, 4 sqrts, and 2 powers. It was generated using SymPy’s CSE functionality, with a limited amount of horner polynomial optimizations as well. It is over 10x faster than a naive implementation.
This implementation has been tested against a straightforward implementation with the equations given in IAPWS-95.
Over a linear temperature range of 200 K to 5000 K and a logarithmic density range of 1E-10 kg/m^3 to 5000 kg/m^3, 250000 points were evaluated. The mean relative error was 5.68E-14, with a maximum relative error of 6.73E-9 and a standard deviation of 1.35E-11.
Over the same range, the model was evaluated to a precision of 50 decimal places with mpmath, and on 90000 points, the mean relative error was 4.66E-14, with a maximum relative error of 4.25E-10 and a standard deviation of 1.77E-12.
The maximum error ocurs in the extremely low density regime, .
Examples
>>> iapws95_dAr_dtau(647.096/300.0, 999.0/322) -7.7043336309570
- chemicals.iapws.iapws95_d2Ar_dtau2(tau, delta)[source]¶
Calculates the second derivative of residual Helmholtz energy of water with respect to tau according to the IAPWS-95 standard.
- Parameters
- Returns
- d2Ar_dtau2
float
Second derivative of residual Helmholtz energy A/(RT) with respect to tau, [-]
- d2Ar_dtau2
Notes
This is an optimized implementatation taking 9 exp calls, 4 sqrts, and 2 powers. It was generated using SymPy’s CSE functionality, with a limited amount of horner polynomial optimizations as well. It is over 10x faster than a naive implementation.
This implementation has been tested against a straightforward implementation with the equations given in IAPWS-95.
Over a linear temperature range of 200 K to 5000 K and a logarithmic density range of 1E-10 kg/m^3 to 5000 kg/m^3, 4E6 points were evaluated. The mean relative error was 4.595E-16, with a maximum relative error of 1.835e-10 and a standard deviation of 1.209E-13.
Over the same range, the model was evaluated to a precision of 50 decimal places with mpmath, and on 250000 points, the mean relative error was 2.6026E-16, with a maximum relative error of 2.36E-12 and a standard deviation of 8.055E-15.
This comparison indicates this implementation is more accurate than the straightforward implementation.
Examples
>>> iapws95_d2Ar_dtau2(647.096/300.0, 999.0/322) -1.2616419775539
- chemicals.iapws.iapws95_d2Ar_ddeltadtau(tau, delta)[source]¶
Calculates the second derivative of residual Helmholtz energy of water with respect to tau and also delta according to the IAPWS-95 standard.
- Parameters
- Returns
- d2Ar_ddeltadtau
float
Second derivative of residual Helmholtz energy A/(RT) with respect to tau and delta, [-]
- d2Ar_ddeltadtau
Notes
This is an optimized implementatation taking 11 exp calls, 4 sqrts, and 3 powers. It was generated using SymPy’s CSE functionality, with select polynomial optimizations by hand as well. It is over 10x faster than a naive implementation.
This implementation has been tested against a straightforward implementation with the equations given in IAPWS-95.
Over a linear temperature range of 200 K to 5000 K and a logarithmic density range of 1E-10 kg/m^3 to 5000 kg/m^3, 4E6 points were evaluated. The mean relative error was 2.82E-14, with a maximum relative error of 8.404E-9 and a standard deviation of 5.166e-12.
Over the same range, the model was evaluated to a precision of 50 decimal places with mpmath, and on 90000 points, the mean relative error was 6.974E-14, with a maximum relative error of 4.286E-9 and a standard deviation of 4.286E-12.
Examples
>>> iapws95_d2Ar_ddeltadtau(647.096/300.0, 999.0/322) -0.198403562385
- chemicals.iapws.iapws95_d3Ar_ddeltadtau2(tau, delta)[source]¶
Calculates the third derivative of residual Helmholtz energy of water with respect to delta once and tau twice according to the IAPWS-95 standard.
- Parameters
- Returns
- d3Ar_ddeltadtau2
float
Third derivative of residual Helmholtz energy A/(RT) with respect to delta once and tau twice, [-]
- d3Ar_ddeltadtau2
Notes
This is an optimized implementatation. It was generated using SymPy’s CSE functionality.
No equation is given for this in IAPWS-95, and the derivative was symbolically computed with SymPy.
Like many higher-order derivatives of functions with exponentials, this one balloons to use many, many terms.
Over a linear temperature range of 200 K to 5000 K and a logarithmic density range of 1E-10 kg/m^3 to 5000 kg/m^3, 250000 points were evaluated. The mean relative error was 7.936e-16, with a maximum relative error of 1.965E-11 and a standard deviation of 4.7938E-14.
Over the same range, the model was evaluated to a precision of 50 decimal places with mpmath, and on 90000 points, the mean relative error was 6.08E-16, with a maximum relative error of 3.537E-12 and a standard deviation of 1.85197E-14.
Examples
>>> iapws95_d3Ar_ddeltadtau2(647.096/300.0, 999.0/322) 1.081479970332
- chemicals.iapws.iapws95_d3Ar_ddelta2dtau(tau, delta)[source]¶
Calculates the third derivative of residual Helmholtz energy of water with respect to delta twice and tau one according to the IAPWS-95 standard.
- Parameters
- Returns
- d3Ar_ddeltadtau2
float
Third derivative of residual Helmholtz energy A/(RT) with respect to delta twice and tau once, [-]
- d3Ar_ddeltadtau2
Notes
This is an optimized implementatation. It was generated using SymPy’s CSE functionality.
No equation is given for this in IAPWS-95, and the derivative was symbolically computed with SymPy.
Like many higher-order derivatives of functions with exponentials, this one balloons to use many, many terms.
Over a linear temperature range of 200 K to 5000 K and a logarithmic density range of 1E-10 kg/m^3 to 5000 kg/m^3, 250000 points were evaluated. The mean relative error was 3.629e-15, with a maximum relative error of 8.38E-11 and a standard deviation of 2.1214E-13.
Over the same range, the model was evaluated to a precision of 50 decimal places with mpmath, and on 10000 points, the mean relative error was 2.4e-15, with a maximum relative error of 7.62E-12 and a standard deviation of 7.818E-14.
Examples
>>> iapws95_d3Ar_ddelta2dtau(647.096/300.0, 999.0/322) 0.015646982949077
- chemicals.iapws.iapws95_d4Ar_ddelta2dtau2(tau, delta)[source]¶
Calculates the fourth derivative of residual Helmholtz energy of water with respect to tau twice and delta twice according to the IAPWS-95 standard.
- Parameters
- Returns
- d4Ar_ddelta2dtau2
float
Fourth derivative of residual Helmholtz energy A/(RT) with respect to tau and delta, [-]
- d4Ar_ddelta2dtau2
Examples
>>> iapws95_d4Ar_ddelta2dtau2(647.096/300.0, 999.0/322) -2.656422915480