Cantera  2.0
Public Member Functions | Private Member Functions | Private Attributes | List of all members
WaterPropsIAPWS Class Reference

Class for calculating the equation of state of water. More...

#include <WaterPropsIAPWS.h>

Collaboration diagram for WaterPropsIAPWS:
[legend]

Public Member Functions

 WaterPropsIAPWS ()
 Base constructor.
 
 WaterPropsIAPWS (const WaterPropsIAPWS &right)
 Copy constructor.
 
WaterPropsIAPWSoperator= (const WaterPropsIAPWS &right)
 assignment constructor
 
 ~WaterPropsIAPWS ()
 destructor
 
void setState_TR (doublereal temperature, doublereal rho)
 Set the internal state of the object wrt temperature and density.
 
doublereal helmholtzFE () const
 Calculate the Helmholtz free energy in mks units of J kmol-1 K-1, using the last temperature and density.
 
doublereal Gibbs () const
 Calculate the Gibbs free energy in mks units of J kmol-1 K-1.
 
doublereal enthalpy () const
 Calculate the enthalpy in mks units of J kmol-1 using the last temperature and density.
 
doublereal intEnergy () const
 Calculate the internal energy in mks units of J kmol-1.
 
doublereal entropy () const
 Calculate the entropy in mks units of J kmol-1 K-1.
 
doublereal cv () const
 Calculate the constant volume heat capacity in mks units of J kmol-1 K-1 at the last temperature and density.
 
doublereal cp () const
 Calculate the constant pressure heat capacity in mks units of J kmol-1 K-1 at the last temperature and density.
 
doublereal molarVolume () const
 Calculate the molar volume (kmol m-3) at the last temperature and density.
 
doublereal pressure () const
 Calculates the pressure (Pascals), given the current value of the temperature and density.
 
doublereal density (doublereal temperature, doublereal pressure, int phase=-1, doublereal rhoguess=-1.0)
 Calculates the density given the temperature and the pressure, and a guess at the density.
 
doublereal density_const (doublereal pressure, int phase=-1, doublereal rhoguess=-1.0) const
 Calculates the density given the temperature and the pressure, and a guess at the density, while not changing the internal state.
 
doublereal density () const
 Returns the density (kg m-3)
 
doublereal temperature () const
 Returns the temperature (Kelvin)
 
doublereal coeffThermExp () const
 Returns the coefficient of thermal expansion.
 
doublereal coeffPresExp () const
 Returns the isochoric pressure derivative wrt temperature.
 
doublereal isothermalCompressibility () const
 Returns the coefficient of isothermal compressibility for the state of the object.
 
doublereal dpdrho () const
 Returns the value of dp / drho at constant T for the state of the object.
 
doublereal psat_est (doublereal temperature) const
 This function returns an estimated value for the saturation pressure.
 
doublereal psat (doublereal temperature, int waterState=WATER_LIQUID)
 This function returns the saturation pressure given the temperature as an input parameter, and sets the internal state to the saturated conditions.
 
doublereal densSpinodalWater () const
 Return the value of the density at the water spinodal point (on the liquid side) for the current temperature.
 
doublereal densSpinodalSteam () const
 Return the value of the density at the water spinodal point (on the gas side) for the current temperature.
 
int phaseState (bool checkState=false) const
 Returns the Phase State flag for the current state of the object.
 
doublereal Tcrit () const
 Returns the critical temperature of water (Kelvin)
 
doublereal Pcrit () const
 Returns the critical pressure of water (22.064E6 Pa)
 
doublereal Rhocrit () const
 Return the critical density of water (kg m-3)
 

Private Member Functions

void calcDim (doublereal temperature, doublereal rho)
 Calculate the dimensionless temp and rho and store internally.
 
void corr (doublereal temperature, doublereal pressure, doublereal &densLiq, doublereal &densGas, doublereal &delGRT)
 Utility routine in the calculation of the saturation pressure.
 
void corr1 (doublereal temperature, doublereal pressure, doublereal &densLiq, doublereal &densGas, doublereal &pcorr)
 Utility routine in the calculation of the saturation pressure.
 

Private Attributes

WaterPropsIAPWSphim_phi
 pointer to the underlying object that does the calculations.
 
doublereal tau
 Dimensionless temperature.
 
doublereal delta
 Dimensionless density.
 
int iState
 Current state of the system.
 

Detailed Description

Class for calculating the equation of state of water.

The reference is W. Wagner, A. Prub, "The IAPWS Formulation 1995 for the Thermodynamic Properties of Ordinary Water Substance for General and Scientific Use," J. Phys. Chem. Ref. Dat, 31, 387, 2002.

This class provides a very complicated polynomial for the specific helmholtz free energy of water, as a function of temperature and density.

\[ \frac{M\hat{f}(\rho,T)}{R T} = \phi(\delta, \tau) = \phi^o(\delta, \tau) + \phi^r(\delta, \tau) \]

where

\[ \delta = \rho / \rho_c \mbox{\qquad and \qquad} \tau = T_c / T \]

The following constants are assumed

\[ T_c = 647.096\mbox{\ K} \]

\[ \rho_c = 322 \mbox{\ kg\ m$^{-3}$} \]

\[ R/M = 0.46151805 \mbox{\ kJ\ kg$^{-1}$\ K$^{-1}$} \]

The free energy is a unique single-valued function of the temperature and density over its entire range.

Note, the base thermodynamic state for this class is the one used in the steam tables, i.e., the liquid at the triple point for water has the following properties:

The class is pretty straightforward. However, one function deserves mention. the density() function calculates the density that is consistent with a particular value of the temperature and pressure. It may therefore be multivalued or potentially there may be no answer from this function. It therefore takes a phase guess and a density guess as optional parameters. If no guesses are supplied to density(), a gas phase guess is assumed. This may or may not be what is wanted. Therefore, density() should usually at least be supplied with a phase guess so that it may manufacture an appropriate density guess. density() manufactures the initial density guess, nondimensionalizes everything, and then calls WaterPropsIAPWSphi::dfind(), which does the iterative calculation to find the density condition that matches the desired input pressure.

The phase guess defines are located in the .h file. they are

There are only three functions which actually change the value of the internal state of this object after it's been instantiated

The setState_TR() is the main function that sets the temperature and rho value. The density() function serves as a setState_TP() function, in that it sets internal state to a temperature and pressure. However, note that this is potentially multivalued. Therefore, we need to supply in addition a phase guess and a rho guess to the input temperature and pressure. The psat() function sets the internal state to the saturated liquid or saturated gas state, depending on the waterState parameter.

Because the underlying object WaterPropsIAPWSphi is privately held, you can be sure that the underlying state of this object doesn't change except due to the three function calls listed above.

Definition at line 158 of file WaterPropsIAPWS.h.

Constructor & Destructor Documentation

Base constructor.

Definition at line 43 of file WaterPropsIAPWS.cpp.

References WaterPropsIAPWS::m_phi.

WaterPropsIAPWS ( const WaterPropsIAPWS right)

Copy constructor.

Parameters
rightObject to be copied

Definition at line 56 of file WaterPropsIAPWS.cpp.

References WaterPropsIAPWS::delta, WaterPropsIAPWS::m_phi, WaterPropsIAPWS::tau, and WaterPropsIAPWSphi::tdpolycalc().

destructor

Definition at line 83 of file WaterPropsIAPWS.cpp.

References WaterPropsIAPWS::m_phi.

Member Function Documentation

WaterPropsIAPWS & operator= ( const WaterPropsIAPWS right)

assignment constructor

Parameters
rightObject to be copied

Definition at line 70 of file WaterPropsIAPWS.cpp.

References WaterPropsIAPWS::delta, WaterPropsIAPWS::iState, WaterPropsIAPWS::m_phi, WaterPropsIAPWS::tau, and WaterPropsIAPWSphi::tdpolycalc().

void setState_TR ( doublereal  temperature,
doublereal  rho 
)
doublereal helmholtzFE ( ) const

Calculate the Helmholtz free energy in mks units of J kmol-1 K-1, using the last temperature and density.

Definition at line 117 of file WaterPropsIAPWS.cpp.

References WaterPropsIAPWS::delta, WaterPropsIAPWS::m_phi, WaterPropsIAPWSphi::phi(), Cantera::Rgas, Cantera::T_c, WaterPropsIAPWS::tau, and WaterPropsIAPWS::temperature().

doublereal Gibbs ( ) const

Calculate the Gibbs free energy in mks units of J kmol-1 K-1.

using the last temperature and density

Definition at line 427 of file WaterPropsIAPWS.cpp.

References WaterPropsIAPWSphi::gibbs_RT(), WaterPropsIAPWS::m_phi, Cantera::Rgas, Cantera::T_c, WaterPropsIAPWS::tau, and WaterPropsIAPWS::temperature().

Referenced by WaterSSTP::getGibbs_RT(), WaterSSTP::getGibbs_RT_ref(), WaterSSTP::getStandardChemPotentials(), and PDSS_Water::gibbs_mole().

doublereal enthalpy ( ) const
doublereal intEnergy ( ) const
doublereal entropy ( ) const

Calculate the entropy in mks units of J kmol-1 K-1.

Definition at line 861 of file WaterPropsIAPWS.cpp.

References WaterPropsIAPWSphi::entropy_R(), WaterPropsIAPWS::m_phi, and Cantera::Rgas.

Referenced by PDSS_Water::entropy_mole(), PDSS_Water::entropy_R_ref(), WaterSSTP::getEntropy_R(), and WaterSSTP::getEntropy_R_ref().

doublereal cv ( ) const

Calculate the constant volume heat capacity in mks units of J kmol-1 K-1 at the last temperature and density.

Definition at line 871 of file WaterPropsIAPWS.cpp.

References WaterPropsIAPWSphi::cv_R(), WaterPropsIAPWS::m_phi, and Cantera::Rgas.

Referenced by WaterSSTP::cv_mole(), and PDSS_Water::cv_mole().

doublereal cp ( ) const

Calculate the constant pressure heat capacity in mks units of J kmol-1 K-1 at the last temperature and density.

Definition at line 879 of file WaterPropsIAPWS.cpp.

References WaterPropsIAPWSphi::cp_R(), WaterPropsIAPWS::m_phi, and Cantera::Rgas.

Referenced by PDSS_Water::cp_mole(), PDSS_Water::cp_R_ref(), WaterSSTP::getCp_R(), and WaterSSTP::getCp_R_ref().

doublereal molarVolume ( ) const

Calculate the molar volume (kmol m-3) at the last temperature and density.

Definition at line 887 of file WaterPropsIAPWS.cpp.

References WaterPropsIAPWS::delta, Cantera::M_water, and Cantera::Rho_c.

Referenced by PDSS_Water::molarVolume(), and PDSS_Water::molarVolume_ref().

doublereal pressure ( ) const

Calculates the pressure (Pascals), given the current value of the temperature and density.

The density is an independent variable in the underlying equation of state

Returns
returns the pressure (Pascal)

Definition at line 131 of file WaterPropsIAPWS.cpp.

References WaterPropsIAPWS::delta, WaterPropsIAPWS::m_phi, Cantera::M_water, WaterPropsIAPWSphi::pressureM_rhoRT(), Cantera::Rgas, Cantera::Rho_c, Cantera::T_c, WaterPropsIAPWS::tau, and WaterPropsIAPWS::temperature().

Referenced by WaterSSTP::pressure(), and PDSS_Water::pressure().

doublereal density ( doublereal  temperature,
doublereal  pressure,
int  phase = -1,
doublereal  rhoguess = -1.0 
)

Calculates the density given the temperature and the pressure, and a guess at the density.

Sets the internal state.

Note, below T_c, this is a multivalued function.

The density() function calculates the density that is consistent with a particular value of the temperature and pressure. It may therefore be multivalued or potentially there may be no answer from this function. It therefore takes a phase guess and a density guess as optional parameters. If no guesses are supplied to density(), a gas phase guess is assumed. This may or may not be what is wanted. Therefore, density() should usually at least be supplied with a phase guess so that it may manufacture an appropriate density guess. density() manufactures the initial density guess, nondimensionalizes everything, and then calls WaterPropsIAPWSphi::dfind(), which does the iterative calculation to find the density condition that matches the desired input pressure.

Parameters
temperature,:Kelvin
pressure: Pressure in Pascals (Newton/m**2)
phase: guessed phase of water : -1: no guessed phase
rhoguess: guessed density of the water : -1.0 no guessed density
Returns
Returns the density. If an error is encountered in the calculation the value of -1.0 is returned.

Definition at line 154 of file WaterPropsIAPWS.cpp.

References WaterPropsIAPWS::delta, WaterPropsIAPWSphi::dfind(), Cantera::int2str(), WaterPropsIAPWS::m_phi, Cantera::M_water, Cantera::Rgas, Cantera::Rho_c, WaterPropsIAPWS::setState_TR(), Cantera::T_c, WaterPropsIAPWS::tau, and WaterPropsIAPWS::temperature().

Referenced by WaterProps::coeffThermalExp_IAPWS(), PDSS_Water::constructSet(), PDSS_Water::cp_R_ref(), WaterProps::density_IAPWS(), WaterSSTP::dthermalExpansionCoeffdT(), PDSS_Water::dthermalExpansionCoeffdT(), PDSS_Water::enthalpy_RT_ref(), PDSS_Water::entropy_R_ref(), WaterSSTP::getCp_R_ref(), WaterSSTP::getEnthalpy_RT_ref(), WaterSSTP::getEntropy_R_ref(), WaterSSTP::getGibbs_RT_ref(), WaterSSTP::getStandardVolumes_ref(), PDSS_Water::gibbs_RT_ref(), WaterSSTP::initThermoXML(), WaterProps::isothermalCompressibility_IAPWS(), PDSS_Water::molarVolume_ref(), PDSS_Water::satPressure(), WaterSSTP::setPressure(), PDSS_Water::setPressure(), WaterProps::thermalConductivityWater(), and WaterProps::viscosityWater().

doublereal density_const ( doublereal  pressure,
int  phase = -1,
doublereal  rhoguess = -1.0 
) const

Calculates the density given the temperature and the pressure, and a guess at the density, while not changing the internal state.

Note, below T_c, this is a multivalued function.

The density() function calculates the density that is consistent with a particular value of the temperature and pressure. It may therefore be multivalued or potentially there may be no answer from this function. It therefore takes a phase guess and a density guess as optional parameters. If no guesses are

supplied to density(), a gas phase guess is assumed. This may or may not be what is wanted. Therefore, density() should usually at least be supplied with a phase guess so that it may manufacture an appropriate density guess. density() manufactures the initial density guess, nondimensionalizes everything, and then calls WaterPropsIAPWSphi::dfind(), which does the iterative calculation to find the density condition that matches the desired input pressure.

Parameters
pressure: Pressure in Pascals (Newton/m**2)
phase: guessed phase of water : -1: no guessed phase
rhoguess: guessed density of the water : -1.0 no guessed density
Returns
Returns the density. If an error is encountered in the calculation the value of -1.0 is returned.

Definition at line 240 of file WaterPropsIAPWS.cpp.

References WaterPropsIAPWS::delta, WaterPropsIAPWSphi::dfind(), Cantera::int2str(), WaterPropsIAPWS::m_phi, Cantera::M_water, Cantera::Rgas, Cantera::Rho_c, Cantera::T_c, WaterPropsIAPWS::tau, WaterPropsIAPWSphi::tdpolycalc(), and WaterPropsIAPWS::temperature().

Referenced by WaterPropsIAPWS::densSpinodalSteam(), and WaterPropsIAPWS::densSpinodalWater().

doublereal density ( ) const

Returns the density (kg m-3)

The density is an independent variable in the underlying equation of state

Returns
Returns the density (kg m-3)

Definition at line 308 of file WaterPropsIAPWS.cpp.

References WaterPropsIAPWS::delta, and Cantera::Rho_c.

Referenced by WaterPropsIAPWS::corr(), WaterPropsIAPWS::corr1(), and WaterPropsIAPWS::psat().

doublereal temperature ( ) const
doublereal coeffThermExp ( ) const
doublereal coeffPresExp ( ) const

Returns the isochoric pressure derivative wrt temperature.

beta = M / (rho * Rgas) (d (pressure) / dT) at constant rho

Note for ideal gases this is equal to one.

beta = delta (phi0_d() + phiR_d())

  • tau delta (phi0_dt() + phiR_dt())

Definition at line 405 of file WaterPropsIAPWS.cpp.

References WaterPropsIAPWS::delta, WaterPropsIAPWSphi::dimdpdT(), WaterPropsIAPWS::m_phi, and WaterPropsIAPWS::tau.

Referenced by WaterPropsIAPWS::coeffThermExp(), and WaterProps::thermalConductivityWater().

doublereal isothermalCompressibility ( ) const

Returns the coefficient of isothermal compressibility for the state of the object.

kappa = - d (ln V) / dP at constant T.

units - 1/Pascal

Returns
returns the isothermal compressibility

Definition at line 374 of file WaterPropsIAPWS.cpp.

References WaterPropsIAPWS::delta, WaterPropsIAPWS::dpdrho(), and Cantera::Rho_c.

Referenced by WaterPropsIAPWS::coeffThermExp(), WaterSSTP::isothermalCompressibility(), PDSS_Water::isothermalCompressibility(), WaterProps::isothermalCompressibility_IAPWS(), and WaterPropsIAPWS::phaseState().

doublereal dpdrho ( ) const
doublereal psat_est ( doublereal  temperature) const

This function returns an estimated value for the saturation pressure.

It does this via a polynomial fit of the vapor pressure curve. units = (Pascals)

Parameters
temperatureInput temperature (Kelvin)
Returns
Returns the estimated saturation pressure

Definition at line 333 of file WaterPropsIAPWS.cpp.

Referenced by WaterPropsIAPWS::densSpinodalSteam(), WaterPropsIAPWS::densSpinodalWater(), PDSS_Water::pref_safe(), and WaterPropsIAPWS::psat().

doublereal psat ( doublereal  temperature,
int  waterState = WATER_LIQUID 
)

This function returns the saturation pressure given the temperature as an input parameter, and sets the internal state to the saturated conditions.

Note this function will return the saturation pressure, given the temperature. It will then set the state of the system to the saturation condition. The input parameter waterState is used to either specify the liquid state or the gas state at the desired temperature and saturated pressure.

If the input temperature, T, is above T_c, this routine will set the internal state to T and the pressure to P_c. Then, return P_c.

Parameters
temperatureinput temperature (kelvin)
waterStateinteger specifying the water state
Returns
Returns the saturation pressure units = Pascal

Definition at line 532 of file WaterPropsIAPWS.cpp.

References WaterPropsIAPWS::corr(), WaterPropsIAPWS::corr1(), WaterPropsIAPWS::density(), Cantera::int2str(), Cantera::M_water, Cantera::P_c, WaterPropsIAPWS::psat_est(), Cantera::Rgas, WaterPropsIAPWS::setState_TR(), and Cantera::T_c.

Referenced by WaterProps::satPressure(), PDSS_Water::satPressure(), and WaterSSTP::satPressure().

doublereal densSpinodalWater ( ) const

Return the value of the density at the water spinodal point (on the liquid side) for the current temperature.

Returns
returns the density with units of kg m-3

Definition at line 634 of file WaterPropsIAPWS.cpp.

References WaterPropsIAPWS::delta, WaterPropsIAPWS::density_const(), WaterPropsIAPWS::dpdrho(), WaterPropsIAPWS::m_phi, ckr::max(), ckr::min(), WaterPropsIAPWS::psat_est(), Cantera::Rho_c, Cantera::T_c, WaterPropsIAPWS::tau, WaterPropsIAPWSphi::tdpolycalc(), and WaterPropsIAPWS::temperature().

doublereal densSpinodalSteam ( ) const

Return the value of the density at the water spinodal point (on the gas side) for the current temperature.

Returns
returns the density with units of kg m-3

Definition at line 731 of file WaterPropsIAPWS.cpp.

References WaterPropsIAPWS::delta, WaterPropsIAPWS::density_const(), WaterPropsIAPWS::dpdrho(), WaterPropsIAPWS::m_phi, ckr::max(), ckr::min(), WaterPropsIAPWS::psat_est(), Cantera::Rho_c, Cantera::T_c, WaterPropsIAPWS::tau, WaterPropsIAPWSphi::tdpolycalc(), and WaterPropsIAPWS::temperature().

int phaseState ( bool  checkState = false) const

Returns the Phase State flag for the current state of the object.

Parameters
checkStateIf true, this function does a complete check to see where in parameters space we are

There are three values: WATER_GAS below the critical temperature but below the critical density WATER_LIQUID below the critical temperature but above the critical density WATER_SUPERCRIT above the critical temperature

Definition at line 584 of file WaterPropsIAPWS.cpp.

References WaterPropsIAPWS::delta, WaterPropsIAPWS::isothermalCompressibility(), WaterPropsIAPWS::iState, WaterPropsIAPWS::m_phi, Cantera::M_water, Cantera::OneAtm, Cantera::Rgas, Cantera::Rho_c, Cantera::T_c, WaterPropsIAPWS::tau, and WaterPropsIAPWSphi::tdpolycalc().

Referenced by PDSS_Water::setPressure().

doublereal Tcrit ( ) const
inline

Returns the critical temperature of water (Kelvin)

This is hard coded to the value 647.096 Kelvin

Definition at line 403 of file WaterPropsIAPWS.h.

Referenced by WaterSSTP::critTemperature(), PDSS_Water::critTemperature(), PDSS_Water::setPressure(), and WaterSSTP::vaporFraction().

doublereal Pcrit ( ) const
inline

Returns the critical pressure of water (22.064E6 Pa)

This is hard coded to the value of 22.064E6 pascals

Definition at line 411 of file WaterPropsIAPWS.h.

Referenced by WaterSSTP::critPressure(), PDSS_Water::critPressure(), and PDSS_Water::pref_safe().

doublereal Rhocrit ( ) const
inline
void calcDim ( doublereal  temperature,
doublereal  rho 
)
private

Calculate the dimensionless temp and rho and store internally.

Private routine

Parameters
temperatureinput temperature (kelvin)
rhodensity in kg m-3

Definition at line 97 of file WaterPropsIAPWS.cpp.

References WaterPropsIAPWS::delta, WaterPropsIAPWS::iState, Cantera::Rho_c, Cantera::T_c, WaterPropsIAPWS::tau, and WaterPropsIAPWS::temperature().

Referenced by WaterPropsIAPWS::setState_TR().

void corr ( doublereal  temperature,
doublereal  pressure,
doublereal &  densLiq,
doublereal &  densGas,
doublereal &  delGRT 
)
private

Utility routine in the calculation of the saturation pressure.

Private routine

Parameters
temperaturetemperature (kelvin)
pressurepressure (Pascal)
densLiqOutput density of liquid
densGasoutput Density of gas
delGRToutput delGRT

Definition at line 449 of file WaterPropsIAPWS.cpp.

References WaterPropsIAPWS::density(), Cantera::fp2str(), WaterPropsIAPWSphi::gibbs_RT(), WaterPropsIAPWS::m_phi, and WaterPropsIAPWS::setState_TR().

Referenced by WaterPropsIAPWS::psat().

void corr1 ( doublereal  temperature,
doublereal  pressure,
doublereal &  densLiq,
doublereal &  densGas,
doublereal &  pcorr 
)
private

Utility routine in the calculation of the saturation pressure.

Private routine

Parameters
temperaturetemperature (kelvin)
pressurepressure (Pascal)
densLiqOutput density of liquid
densGasoutput Density of gas
pcorroutput corrected pressure

Definition at line 485 of file WaterPropsIAPWS.cpp.

References WaterPropsIAPWS::density(), Cantera::fp2str(), WaterPropsIAPWS::m_phi, Cantera::M_water, WaterPropsIAPWSphi::phiR(), Cantera::Rgas, and WaterPropsIAPWS::setState_TR().

Referenced by WaterPropsIAPWS::psat().

Member Data Documentation

WaterPropsIAPWSphi* m_phi
private
doublereal tau
private
doublereal delta
mutableprivate
int iState
mutableprivate

Current state of the system.

Definition at line 477 of file WaterPropsIAPWS.h.

Referenced by WaterPropsIAPWS::calcDim(), WaterPropsIAPWS::operator=(), and WaterPropsIAPWS::phaseState().


The documentation for this class was generated from the following files: