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