Cantera  2.3.0
WaterPropsIAPWS.h
Go to the documentation of this file.
1 /**
2  * @file WaterPropsIAPWS.h
3  * Headers for a class for calculating the equation of state of water
4  * from the IAPWS 1995 Formulation based on the steam tables thermodynamic
5  * basis (See class \link Cantera::WaterPropsIAPWS WaterPropsIAPWS\endlink).
6  */
7 
8 // This file is part of Cantera. See License.txt in the top-level directory or
9 // at http://www.cantera.org/license.txt for license and copyright information.
10 
11 #ifndef WATERPROPSIAPWS_H
12 #define WATERPROPSIAPWS_H
13 
14 #include "WaterPropsIAPWSphi.h"
15 
16 namespace Cantera
17 {
18 /**
19  * @name Names for the phase regions
20  *
21  * These constants are defined and used in the interface to describe the
22  * location of where we are in (T,rho) space.
23  *
24  * WATER_UNSTABLELIQUID indicates that we are in the unstable region, inside the
25  * spinodal curve where dpdrho < 0.0 amongst other properties. The difference
26  * between WATER_UNSTABLELIQUID and WATER_UNSTABLEGAS is that
27  * for WATER_UNSTABLELIQUID d2pdrho2 > 0 and dpdrho < 0.0
28  * for WATER_UNSTABLEGAS d2pdrho2 < 0 and dpdrho < 0.0
29  */
30 //@{
31 #define WATER_GAS 0
32 #define WATER_LIQUID 1
33 #define WATER_SUPERCRIT 2
34 #define WATER_UNSTABLELIQUID 3
35 #define WATER_UNSTABLEGAS 4
36 //@}
37 
38 //! Class for calculating the equation of state of water.
39 /*!
40  * The reference is W. Wagner, A. Pruss, "The IAPWS Formulation 1995 for the
41  * Thermodynamic Properties of Ordinary Water Substance for General and
42  * Scientific Use," J. Phys. Chem. Ref. Dat, 31, 387, 2002.
43  *
44  * This class provides a very complicated polynomial for the specific
45  * Helmholtz free energy of water, as a function of temperature and density.
46  *
47  * \f[
48  * \frac{M\hat{f}(\rho,T)}{R T} = \phi(\delta, \tau) =
49  * \phi^o(\delta, \tau) + \phi^r(\delta, \tau)
50  * \f]
51  *
52  * where
53  *
54  * \f[
55  * \delta = \rho / \rho_c \quad \mathrm{and} \quad \tau = T_c / T
56  * \f]
57  *
58  * The following constants are assumed
59  *
60  * \f[
61  * T_c = 647.096\mathrm{\;K}
62  * \f]
63  * \f[
64  * \rho_c = 322 \mathrm{\;kg\,m^{-3}}
65  * \f]
66  * \f[
67  * R/M = 0.46151805 \mathrm{\;kJ\,kg^{-1}\,K^{-1}}
68  * \f]
69  *
70  * The free energy is a unique single-valued function of the temperature and
71  * density over its entire range.
72  *
73  * Note, the base thermodynamic state for this class is the one used in the
74  * steam tables, i.e., the liquid at the triple point for water has the
75  * following properties:
76  *
77  * - u(273.16, rho) = 0.0
78  * - s(273.16, rho) = 0.0
79  * - psat(273.16) = 611.655 Pascal
80  * - rho(273.16, psat) = 999.793 kg m-3
81  *
82  * Therefore, to use this class within %Cantera, offsets to u() and s() must
83  * be used to put the water class onto the same basis as other thermodynamic
84  * quantities. For example, in the WaterSSTP class, these offsets are
85  * calculated in the following way. The thermodynamic base state for water is
86  * set to the NIST basis here by specifying constants EW_Offset and SW_Offset.
87  * These offsets are calculated on the fly so that the following properties
88  * hold:
89  *
90  * - Delta_Hfo_idealGas(298.15, 1bar) = -241.826 kJ/gmol
91  * - So_idealGas(298.15, 1bar) = 188.835 J/gmolK
92  *
93  * The offsets are calculated by actually computing the above quantities and
94  * then calculating the correction factor.
95  *
96  * This class provides an interface to the WaterPropsIAPWSphi class, which
97  * actually calculates the \f$ \phi^o(\delta, \tau) \f$ and the
98  * \f$ \phi^r(\delta, \tau) \f$ polynomials in dimensionless form.
99  *
100  * All thermodynamic results from this class are returned in dimensional form.
101  * This is because the gas constant (and molecular weight) used within this
102  * class is allowed to be potentially different than that used elsewhere in
103  * %Cantera. Therefore, everything has to be in dimensional units. Note,
104  * however, the thermodynamic basis is set to that used in the steam tables.
105  * (u = s = 0 for liquid water at the triple point).
106  *
107  * This class is not a ThermoPhase. However, it does maintain an internal
108  * state of the object that is dependent on temperature and density. The
109  * internal state is characterized by an internally stored \f$ \tau\f$ and a
110  * \f$ \delta \f$ value, and an iState value, which indicates whether the
111  * point is a liquid, a gas, or a supercritical fluid. Along with that the
112  * \f$ \tau\f$ and a \f$ \delta \f$ values are polynomials of \f$ \tau\f$ and
113  * a \f$ \delta \f$ that are kept by the WaterPropsIAPWSphi class. Therefore,
114  * whenever \f$ \tau\f$ or \f$ \delta \f$ is changed, the function setState()
115  * must be called in order for the internal state to be kept up to date.
116  *
117  * The class is pretty straightforward. However, one function deserves
118  * mention. The density() function calculates the density that is consistent
119  * with a particular value of the temperature and pressure. It may therefore
120  * be multivalued or potentially there may be no answer from this function. It
121  * therefore takes a phase guess and a density guess as optional parameters.
122  * If no guesses are supplied to density(), a gas phase guess is assumed. This
123  * may or may not be what is wanted. Therefore, density() should usually at
124  * least be supplied with a phase guess so that it may manufacture an
125  * appropriate density guess. density() manufactures the initial density
126  * guess, nondimensionalizes everything, and then calls
127  * WaterPropsIAPWSphi::dfind(), which does the iterative calculation to find
128  * the density condition that matches the desired input pressure.
129  *
130  * The phase guess defines are located in the .h file. they are
131  *
132  * - WATER_GAS
133  * - WATER_LIQUID
134  * - WATER_SUPERCRIT
135  *
136  * There are only three functions which actually change the value of the
137  * internal state of this object after it's been instantiated
138  *
139  * - setState_TR(temperature, rho)
140  * - density(temperature, pressure, phase, rhoguess)
141  * - psat(temperature, waterState);
142  *
143  * The setState_TR() is the main function that sets the temperature and rho
144  * value. The density() function serves as a setState_TP() function, in that
145  * it sets internal state to a temperature and pressure. However, note that
146  * this is potentially multivalued. Therefore, we need to supply in addition a
147  * phase guess and a rho guess to the input temperature and pressure. The
148  * psat() function sets the internal state to the saturated liquid or
149  * saturated gas state, depending on the waterState parameter.
150  *
151  * Because the underlying object WaterPropsIAPWSphi is privately held, you can
152  * be sure that the underlying state of this object doesn't change except due
153  * to the three function calls listed above.
154  *
155  * @ingroup thermoprops
156  */
158 {
159 public:
160  //! Base constructor
161  WaterPropsIAPWS();
162 
163  WaterPropsIAPWS(const WaterPropsIAPWS& right);
164  WaterPropsIAPWS& operator=(const WaterPropsIAPWS& right);
165 
166  //! Set the internal state of the object wrt temperature and density
167  /*!
168  * @param temperature temperature (kelvin)
169  * @param rho density (kg m-3)
170  */
171  void setState_TR(doublereal temperature, doublereal rho);
172 
173  //! Calculate the Helmholtz free energy in mks units of J kmol-1 K-1,
174  //! using the last temperature and density
175  doublereal helmholtzFE() const;
176 
177  //! Calculate the Gibbs free energy in mks units of J kmol-1 K-1.
178  //! using the last temperature and density
179  doublereal Gibbs() const;
180 
181  //! Calculate the enthalpy in mks units of J kmol-1
182  //! using the last temperature and density
183  doublereal enthalpy() const;
184 
185  //! Calculate the internal energy in mks units of J kmol-1
186  doublereal intEnergy() const;
187 
188  //! Calculate the entropy in mks units of J kmol-1 K-1
189  doublereal entropy() const;
190 
191  //! Calculate the constant volume heat capacity in mks units of J kmol-1 K-1
192  //! at the last temperature and density
193  doublereal cv() const;
194 
195  //! Calculate the constant pressure heat capacity in mks units of J kmol-1 K-1
196  //! at the last temperature and density
197  doublereal cp() const;
198 
199  //! Calculate the molar volume (kmol m-3) at the last temperature and
200  //! density
201  doublereal molarVolume() const;
202 
203  //! Calculates the pressure (Pascals), given the current value of the
204  //! temperature and density.
205  /*!
206  * The density is an independent variable in the underlying equation of state
207  *
208  * @returns the pressure (Pascal)
209  */
210  doublereal pressure() const;
211 
212  //! Calculates the density given the temperature and the pressure,
213  //! and a guess at the density. Sets the internal state.
214  /*!
215  * Note, below T_c, this is a multivalued function.
216  *
217  * The density() function calculates the density that is consistent with
218  * a particular value of the temperature and pressure. It may therefore be
219  * multivalued or potentially there may be no answer from this function.
220  * It therefore takes a phase guess and a density guess as optional
221  * parameters. If no guesses are supplied to density(), a gas phase guess
222  * is assumed. This may or may not be what is wanted. Therefore, density()
223  * should usually at least be supplied with a phase guess so that it may
224  * manufacture an appropriate density guess. density() manufactures the
225  * initial density guess, nondimensionalizes everything, and then calls
226  * WaterPropsIAPWSphi::dfind(), which does the iterative calculation to
227  * find the density condition that matches the desired input pressure.
228  *
229  * @param temperature Kelvin
230  * @param pressure Pressure in Pascals (Newton/m**2)
231  * @param phase guessed phase of water; -1: no guessed phase
232  * @param rhoguess guessed density of the water; -1.0 no guessed density
233  * @returns the density. If an error is encountered in the calculation the
234  * value of -1.0 is returned.
235  */
236  doublereal density(doublereal temperature, doublereal pressure,
237  int phase = -1, doublereal rhoguess = -1.0);
238 
239  //! Calculates the density given the temperature and the pressure,
240  //! and a guess at the density, while not changing the internal state
241  /*!
242  * Note, below T_c, this is a multivalued function.
243  *
244  * The density() function calculates the density that is consistent with a
245  * particular value of the temperature and pressure. It may therefore be
246  * multivalued or potentially there may be no answer from this function.
247  * It therefore takes a phase guess and a density guess as optional
248  * parameters. If no guesses are supplied to density(), a gas phase guess
249  * is assumed. This may or may not be what is wanted. Therefore, density()
250  * should usually at least be supplied with a phase guess so that it may
251  * manufacture an appropriate density guess. density() manufactures the
252  * initial density guess, nondimensionalizes everything, and then calls
253  * WaterPropsIAPWSphi::dfind(), which does the iterative calculation to
254  * find the density condition that matches the desired input pressure.
255  *
256  * @param pressure Pressure in Pascals (Newton/m**2)
257  * @param phase guessed phase of water; -1: no guessed phase
258  * @param rhoguess guessed density of the water; -1.0: no guessed density
259  * @returns the density. If an error is encountered in the calculation the
260  * value of -1.0 is returned.
261  */
262  doublereal density_const(doublereal pressure, int phase = -1, doublereal rhoguess = -1.0) const;
263 
264  //! Returns the density (kg m-3)
265  /*!
266  * The density is an independent variable in the underlying equation of state
267  *
268  * @returns the density (kg m-3)
269  */
270  doublereal density() const;
271 
272  //! Returns the temperature (Kelvin)
273  /*!
274  * @return s the internally stored temperature
275  */
276  doublereal temperature() const;
277 
278  //! Returns the coefficient of thermal expansion.
279  /*!
280  * alpha = d (ln V) / dT at constant P.
281  *
282  * @returns the coefficient of thermal expansion
283  */
284  doublereal coeffThermExp() const;
285 
286  //! Returns the isochoric pressure derivative wrt temperature
287  /*!
288  * beta = M / (rho * Rgas) (d (pressure) / dT) at constant rho
289  *
290  * Note for ideal gases this is equal to one.
291  *
292  * beta = delta (phi0_d() + phiR_d()) - tau delta (phi0_dt() + phiR_dt())
293  */
294  doublereal coeffPresExp() const;
295 
296  //! Returns the coefficient of isothermal compressibility for the state of
297  //! the object
298  /*!
299  * kappa = - d (ln V) / dP at constant T.
300  *
301  * units - 1/Pascal
302  *
303  * @returns the isothermal compressibility
304  */
305  doublereal isothermalCompressibility() const;
306 
307  //! Returns the value of dp / drho at constant T for the state of the object
308  /*!
309  * units - Joules / kg
310  *
311  * @returns dpdrho
312  */
313  doublereal dpdrho() const;
314 
315  //! This function returns an estimated value for the saturation pressure.
316  /*!
317  * It does this via a polynomial fit of the vapor pressure curve.
318  * units = (Pascals)
319  *
320  * @param temperature Input temperature (Kelvin)
321  *
322  * @returns the estimated saturation pressure
323  */
324  doublereal psat_est(doublereal temperature) const;
325 
326  //! This function returns the saturation pressure given the temperature as
327  //! an input parameter, and sets the internal state to the saturated
328  //! conditions.
329  /*!
330  * Note this function will return the saturation pressure, given the
331  * temperature. It will then set the state of the system to the saturation
332  * condition. The input parameter waterState is used to either specify the
333  * liquid state or the gas state at the desired temperature and saturated
334  * pressure.
335  *
336  * If the input temperature, T, is above T_c, this routine will set the
337  * internal state to T and the pressure to P_c. Then, return P_c.
338  *
339  * @param temperature input temperature (kelvin)
340  * @param waterState integer specifying the water state
341  * @returns the saturation pressure. units = Pascal
342  */
343  doublereal psat(doublereal temperature, int waterState = WATER_LIQUID);
344 
345  //! Return the value of the density at the water spinodal point (on the
346  //! liquid side) for the current temperature.
347  /*!
348  * @returns the density with units of kg m-3
349  */
350  doublereal densSpinodalWater() const;
351 
352  //! Return the value of the density at the water spinodal point (on the gas
353  //! side) for the current temperature.
354  /*!
355  * @returns the density with units of kg m-3
356  */
357  doublereal densSpinodalSteam() const;
358 
359  //! Returns the Phase State flag for the current state of the object
360  /*!
361  * @param checkState If true, this function does a complete check to see
362  * where in parameters space we are
363  *
364  * There are three values:
365  * - WATER_GAS below the critical temperature but below the critical density
366  * - WATER_LIQUID below the critical temperature but above the critical density
367  * - WATER_SUPERCRIT above the critical temperature
368  */
369  int phaseState(bool checkState = false) const;
370 
371  //! Returns the critical temperature of water (Kelvin)
372  /*!
373  * This is hard coded to the value 647.096 Kelvin
374  */
375  doublereal Tcrit() const {
376  return 647.096;
377  }
378 
379  //! Returns the critical pressure of water (22.064E6 Pa)
380  /*!
381  * This is hard coded to the value of 22.064E6 pascals
382  */
383  doublereal Pcrit() const {
384  return 22.064E6;
385  }
386 
387  //! Return the critical density of water (kg m-3)
388  /*!
389  * This is equal to 322 kg m-3.
390  */
391  doublereal Rhocrit() const {
392  return 322.;
393  }
394 
395 private:
396  //! Calculate the dimensionless temp and rho and store internally.
397  /*!
398  * @param temperature input temperature (kelvin)
399  * @param rho density in kg m-3
400  */
401  void calcDim(doublereal temperature, doublereal rho);
402 
403  //! Utility routine in the calculation of the saturation pressure
404  /*!
405  * Calculate the Gibbs free energy in mks units of J kmol-1 K-1.
406  *
407  * @param temperature temperature (kelvin)
408  * @param pressure pressure (Pascal)
409  * @param densLiq Output density of liquid
410  * @param densGas output Density of gas
411  * @param delGRT output delGRT
412  */
413  void corr(doublereal temperature, doublereal pressure, doublereal& densLiq,
414  doublereal& densGas, doublereal& delGRT);
415 
416  //! Utility routine in the calculation of the saturation pressure
417  /*!
418  * @param temperature temperature (kelvin)
419  * @param pressure pressure (Pascal)
420  * @param densLiq Output density of liquid
421  * @param densGas output Density of gas
422  * @param pcorr output corrected pressure
423  */
424  void corr1(doublereal temperature, doublereal pressure, doublereal& densLiq,
425  doublereal& densGas, doublereal& pcorr);
426 
427  //! pointer to the underlying object that does the calculations.
429 
430  //! Dimensionless temperature, tau = T_C / T
431  doublereal tau;
432 
433  //! Dimensionless density, delta = rho / rho_c
434  mutable doublereal delta;
435 
436  //! Current state of the system
437  mutable int iState;
438 };
439 
440 }
441 #endif
doublereal psat_est(doublereal temperature) const
This function returns an estimated value for the saturation pressure.
doublereal dpdrho() const
Returns the value of dp / drho at constant T for the state of the object.
doublereal densSpinodalSteam() const
Return the value of the density at the water spinodal point (on the gas side) for the current tempera...
doublereal delta
Dimensionless density, delta = rho / rho_c.
doublereal enthalpy() const
Calculate the enthalpy in mks units of J kmol-1 using the last temperature and density.
doublereal Pcrit() const
Returns the critical pressure of water (22.064E6 Pa)
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...
void setState_TR(doublereal temperature, doublereal rho)
Set the internal state of the object wrt temperature and density.
doublereal intEnergy() const
Calculate the internal energy in mks units of J kmol-1.
Header for Lowest level of the classes which support a real water model (see class WaterPropsIAPWS an...
doublereal entropy() const
Calculate the entropy in mks units of J kmol-1 K-1.
int phaseState(bool checkState=false) const
Returns the Phase State flag for the current state of the object.
Class for calculating the equation of state of water.
WaterPropsIAPWS()
Base constructor.
doublereal molarVolume() const
Calculate the molar volume (kmol m-3) at the last temperature and density.
void corr(doublereal temperature, doublereal pressure, doublereal &densLiq, doublereal &densGas, doublereal &delGRT)
Utility routine in the calculation of the saturation pressure.
doublereal temperature() const
Returns the temperature (Kelvin)
doublereal pressure() const
Calculates the pressure (Pascals), given the current value of the temperature and density...
void calcDim(doublereal temperature, doublereal rho)
Calculate the dimensionless temp and rho and store internally.
doublereal cp() const
Calculate the constant pressure heat capacity in mks units of J kmol-1 K-1 at the last temperature an...
doublereal density() const
Returns the density (kg m-3)
doublereal Rhocrit() const
Return the critical density of water (kg m-3)
int iState
Current state of the system.
doublereal isothermalCompressibility() const
Returns the coefficient of isothermal compressibility for the state of the object.
doublereal helmholtzFE() const
Calculate the Helmholtz free energy in mks units of J kmol-1 K-1, using the last temperature and dens...
doublereal tau
Dimensionless temperature, tau = T_C / T.
doublereal coeffThermExp() const
Returns the coefficient of thermal expansion.
doublereal coeffPresExp() const
Returns the isochoric pressure derivative wrt temperature.
Low level class for the real description of water.
void corr1(doublereal temperature, doublereal pressure, doublereal &densLiq, doublereal &densGas, doublereal &pcorr)
Utility routine in the calculation of the saturation pressure.
WaterPropsIAPWSphi m_phi
pointer to the underlying object that does the calculations.
doublereal Tcrit() const
Returns the critical temperature of water (Kelvin)
doublereal Gibbs() const
Calculate the Gibbs free energy in mks units of J kmol-1 K-1.
doublereal densSpinodalWater() const
Return the value of the density at the water spinodal point (on the liquid side) for the current temp...
doublereal cv() const
Calculate the constant volume heat capacity in mks units of J kmol-1 K-1 at the last temperature and ...
Namespace for the Cantera kernel.
Definition: application.cpp:29
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.