Cantera  3.1.0a1
IdealGasPhase.h
Go to the documentation of this file.
1 /**
2  * @file IdealGasPhase.h
3  * ThermoPhase object for the ideal gas equation of
4  * state - workhorse for %Cantera (see @ref thermoprops
5  * and class @link Cantera::IdealGasPhase IdealGasPhase@endlink).
6  */
7 
8 // This file is part of Cantera. See License.txt in the top-level directory or
9 // at https://cantera.org/license.txt for license and copyright information.
10 
11 #ifndef CT_IDEALGASPHASE_H
12 #define CT_IDEALGASPHASE_H
13 
14 #include "ThermoPhase.h"
15 
16 namespace Cantera
17 {
18 
19 //! Class IdealGasPhase represents low-density gases that obey the ideal gas
20 //! equation of state.
21 /*!
22  *
23  * IdealGasPhase derives from class ThermoPhase, and overloads the virtual
24  * methods defined there with ones that use expressions appropriate for ideal
25  * gas mixtures.
26  *
27  * The independent unknowns are density, mass fraction, and temperature. the
28  * #setPressure() function will calculate the density consistent with the
29  * current mass fraction vector and temperature and the desired pressure, and
30  * then set the density.
31  *
32  * ## Specification of Species Standard State Properties
33  *
34  * It is assumed that the reference state thermodynamics may be obtained by a
35  * pointer to a populated species thermodynamic property manager class in the
36  * base class, ThermoPhase::m_spthermo (see the base class @link
37  * Cantera::MultiSpeciesThermo MultiSpeciesThermo @endlink for a description of
38  * the specification of reference state species thermodynamics functions). The
39  * reference state, where the pressure is fixed at a single pressure, is a key
40  * species property calculation for the Ideal Gas Equation of state.
41  *
42  * This class is optimized for speed of execution. All calls to thermodynamic
43  * functions first call internal routines (aka #enthalpy_RT_ref()) which return
44  * references the reference state thermodynamics functions. Within these
45  * internal reference state functions, the function #updateThermo() is called,
46  * that first checks to see whether the temperature has changed. If it has, it
47  * updates the internal reference state thermo functions by calling the
48  * MultiSpeciesThermo object.
49  *
50  * Functions for the calculation of standard state properties for species at
51  * arbitrary pressure are provided in IdealGasPhase. However, they are all
52  * derived from their reference state counterparts.
53  *
54  * The standard state enthalpy is independent of pressure:
55  *
56  * @f[
57  * h^o_k(T,P) = h^{ref}_k(T)
58  * @f]
59  *
60  * The standard state constant-pressure heat capacity is independent of pressure:
61  *
62  * @f[
63  * Cp^o_k(T,P) = Cp^{ref}_k(T)
64  * @f]
65  *
66  * The standard state entropy depends in the following fashion on pressure:
67  *
68  * @f[
69  * S^o_k(T,P) = S^{ref}_k(T) - R \ln(\frac{P}{P_{ref}})
70  * @f]
71  * The standard state Gibbs free energy is obtained from the enthalpy and entropy
72  * functions:
73  *
74  * @f[
75  * \mu^o_k(T,P) = h^o_k(T,P) - S^o_k(T,P) T
76  * @f]
77  *
78  * @f[
79  * \mu^o_k(T,P) = \mu^{ref}_k(T) + R T \ln( \frac{P}{P_{ref}})
80  * @f]
81  *
82  * where
83  * @f[
84  * \mu^{ref}_k(T) = h^{ref}_k(T) - T S^{ref}_k(T)
85  * @f]
86  *
87  * The standard state internal energy is obtained from the enthalpy function also
88  *
89  * @f[
90  * u^o_k(T,P) = h^o_k(T) - R T
91  * @f]
92  *
93  * The molar volume of a species is given by the ideal gas law
94  *
95  * @f[
96  * V^o_k(T,P) = \frac{R T}{P}
97  * @f]
98  *
99  * where R is the molar gas constant. For a complete list of physical constants
100  * used within %Cantera, see @ref physConstants .
101  *
102  * ## Specification of Solution Thermodynamic Properties
103  *
104  * The activity of a species defined in the phase is given by the ideal gas law:
105  * @f[
106  * a_k = X_k
107  * @f]
108  * where @f$ X_k @f$ is the mole fraction of species *k*. The chemical potential
109  * for species *k* is equal to
110  *
111  * @f[
112  * \mu_k(T,P) = \mu^o_k(T, P) + R T \ln X_k
113  * @f]
114  *
115  * In terms of the reference state, the above can be rewritten
116  *
117  * @f[
118  * \mu_k(T,P) = \mu^{ref}_k(T, P) + R T \ln \frac{P X_k}{P_{ref}}
119  * @f]
120  *
121  * The partial molar entropy for species *k* is given by the following relation,
122  *
123  * @f[
124  * \tilde{s}_k(T,P) = s^o_k(T,P) - R \ln X_k = s^{ref}_k(T) - R \ln \frac{P X_k}{P_{ref}}
125  * @f]
126  *
127  * The partial molar enthalpy for species *k* is
128  *
129  * @f[
130  * \tilde{h}_k(T,P) = h^o_k(T,P) = h^{ref}_k(T)
131  * @f]
132  *
133  * The partial molar Internal Energy for species *k* is
134  *
135  * @f[
136  * \tilde{u}_k(T,P) = u^o_k(T,P) = u^{ref}_k(T)
137  * @f]
138  *
139  * The partial molar Heat Capacity for species *k* is
140  *
141  * @f[
142  * \tilde{Cp}_k(T,P) = Cp^o_k(T,P) = Cp^{ref}_k(T)
143  * @f]
144  *
145  * ## Application within Kinetics Managers
146  *
147  * @f$ C^a_k @f$ are defined such that @f$ a_k = C^a_k / C^s_k, @f$ where @f$
148  * C^s_k @f$ is a standard concentration defined below and @f$ a_k @f$ are
149  * activities used in the thermodynamic functions. These activity (or
150  * generalized) concentrations are used by kinetics manager classes to compute
151  * the forward and reverse rates of elementary reactions. The activity
152  * concentration,@f$ C^a_k @f$,is given by the following expression.
153  *
154  * @f[
155  * C^a_k = C^s_k X_k = \frac{P}{R T} X_k
156  * @f]
157  *
158  * The standard concentration for species *k* is independent of *k* and equal to
159  *
160  * @f[
161  * C^s_k = C^s = \frac{P}{R T}
162  * @f]
163  *
164  * For example, a bulk-phase binary gas reaction between species j and k,
165  * producing a new gas species l would have the following equation for its rate
166  * of progress variable, @f$ R^1 @f$, which has units of kmol m-3 s-1.
167  *
168  * @f[
169  * R^1 = k^1 C_j^a C_k^a = k^1 (C^s a_j) (C^s a_k)
170  * @f]
171  * where
172  * @f[
173  * C_j^a = C^s a_j \quad \mbox{and} \quad C_k^a = C^s a_k
174  * @f]
175  *
176  * @f$ C_j^a @f$ is the activity concentration of species j, and
177  * @f$ C_k^a @f$ is the activity concentration of species k. @f$ C^s @f$ is the
178  * standard concentration. @f$ a_j @f$ is the activity of species j which is
179  * equal to the mole fraction of j.
180  *
181  * The reverse rate constant can then be obtained from the law of microscopic
182  * reversibility and the equilibrium expression for the system.
183  *
184  * @f[
185  * \frac{a_j a_k}{ a_l} = K_a^{o,1} = \exp(\frac{\mu^o_l - \mu^o_j - \mu^o_k}{R T} )
186  * @f]
187  *
188  * @f$ K_a^{o,1} @f$ is the dimensionless form of the equilibrium constant,
189  * associated with the pressure dependent standard states @f$ \mu^o_l(T,P) @f$
190  * and their associated activities, @f$ a_l @f$, repeated here:
191  *
192  * @f[
193  * \mu_l(T,P) = \mu^o_l(T, P) + R T \ln a_l
194  * @f]
195  *
196  * We can switch over to expressing the equilibrium constant in terms of the
197  * reference state chemical potentials
198  *
199  * @f[
200  * K_a^{o,1} = \exp(\frac{\mu^{ref}_l - \mu^{ref}_j - \mu^{ref}_k}{R T} ) * \frac{P_{ref}}{P}
201  * @f]
202  *
203  * The concentration equilibrium constant, @f$ K_c @f$, may be obtained by
204  * changing over to activity concentrations. When this is done:
205  *
206  * @f[
207  * \frac{C^a_j C^a_k}{ C^a_l} = C^o K_a^{o,1} = K_c^1 =
208  * \exp(\frac{\mu^{ref}_l - \mu^{ref}_j - \mu^{ref}_k}{R T} ) * \frac{P_{ref}}{RT}
209  * @f]
210  *
211  * %Kinetics managers will calculate the concentration equilibrium constant,
212  * @f$ K_c @f$, using the second and third part of the above expression as a
213  * definition for the concentration equilibrium constant.
214  *
215  * For completeness, the pressure equilibrium constant may be obtained as well
216  *
217  * @f[
218  * \frac{P_j P_k}{ P_l P_{ref}} = K_p^1 =
219  * \exp\left(\frac{\mu^{ref}_l - \mu^{ref}_j - \mu^{ref}_k}{R T} \right)
220  * @f]
221  *
222  * @f$ K_p @f$ is the simplest form of the equilibrium constant for ideal gases.
223  * However, it isn't necessarily the simplest form of the equilibrium constant
224  * for other types of phases; @f$ K_c @f$ is used instead because it is
225  * completely general.
226  *
227  * The reverse rate of progress may be written down as
228  * @f[
229  * R^{-1} = k^{-1} C_l^a = k^{-1} (C^o a_l)
230  * @f]
231  *
232  * where we can use the concept of microscopic reversibility to write the
233  * reverse rate constant in terms of the forward rate constant and the
234  * concentration equilibrium constant, @f$ K_c @f$.
235  *
236  * @f[
237  * k^{-1} = k^1 K^1_c
238  * @f]
239  *
240  * @f$ k^{-1} @f$ has units of s-1.
241  *
242  * ## YAML Example
243  *
244  * An example ideal gas phase definition is given in the
245  * <a href="../../sphinx/html/yaml/phases.html#ideal-gas">YAML API Reference</a>.
246  *
247  * @ingroup thermoprops
248  */
250 {
251 public:
252  //! Construct and initialize an IdealGasPhase ThermoPhase object
253  //! directly from an input file
254  /*!
255  * @param inputFile Name of the input file containing the phase definition
256  * to set up the object. If blank, an empty phase will be
257  * created.
258  * @param id ID of the phase in the input file. Defaults to the
259  * empty string.
260  */
261  explicit IdealGasPhase(const string& inputFile="", const string& id="");
262 
263  string type() const override {
264  return "ideal-gas";
265  }
266 
267  bool isIdeal() const override {
268  return true;
269  }
270 
271  //! String indicating the mechanical phase of the matter in this Phase.
272  /*!
273  * For the `IdealGasPhase`, this string is always `gas`.
274  */
275  string phaseOfMatter() const override {
276  return "gas";
277  }
278 
279  //! @name Molar Thermodynamic Properties of the Solution
280  //! @{
281 
282  //! Return the Molar enthalpy. Units: J/kmol.
283  /*!
284  * For an ideal gas mixture,
285  * @f[
286  * \hat h(T) = \sum_k X_k \hat h^0_k(T),
287  * @f]
288  * and is a function only of temperature. The standard-state pure-species
289  * enthalpies @f$ \hat h^0_k(T) @f$ are computed by the species
290  * thermodynamic property manager.
291  *
292  * \see MultiSpeciesThermo
293  */
294  double enthalpy_mole() const override {
295  return RT() * mean_X(enthalpy_RT_ref());
296  }
297 
298  /**
299  * Molar entropy. Units: J/kmol/K.
300  * For an ideal gas mixture,
301  * @f[
302  * \hat s(T, P) = \sum_k X_k \hat s^0_k(T) - \hat R \ln \frac{P}{P^0}.
303  * @f]
304  * The reference-state pure-species entropies @f$ \hat s^0_k(T) @f$ are
305  * computed by the species thermodynamic property manager.
306  * @see MultiSpeciesThermo
307  */
308  double entropy_mole() const override;
309 
310  /**
311  * Molar heat capacity at constant pressure. Units: J/kmol/K.
312  * For an ideal gas mixture,
313  * @f[
314  * \hat c_p(t) = \sum_k \hat c^0_{p,k}(T).
315  * @f]
316  * The reference-state pure-species heat capacities @f$ \hat c^0_{p,k}(T) @f$
317  * are computed by the species thermodynamic property manager.
318  * @see MultiSpeciesThermo
319  */
320  double cp_mole() const override;
321 
322  /**
323  * Molar heat capacity at constant volume. Units: J/kmol/K.
324  * For an ideal gas mixture,
325  * @f[ \hat c_v = \hat c_p - \hat R. @f]
326  */
327  double cv_mole() const override;
328 
329  //! @}
330  //! @name Mechanical Equation of State
331  //! @{
332 
333  /**
334  * Pressure. Units: Pa.
335  * For an ideal gas mixture,
336  * @f[ P = n \hat R T. @f]
337  */
338  double pressure() const override {
339  return GasConstant * molarDensity() * temperature();
340  }
341 
342  //! Set the pressure at constant temperature and composition.
343  /*!
344  * Units: Pa.
345  * This method is implemented by setting the mass density to
346  * @f[
347  * \rho = \frac{P \overline W}{\hat R T }.
348  * @f]
349  *
350  * @param p Pressure (Pa)
351  */
352  void setPressure(double p) override {
353  setDensity(p * meanMolecularWeight() / RT());
354  }
355 
356  //! Set the density and pressure at constant composition.
357  /*!
358  * Units: kg/m^3, Pa.
359  * This method is implemented by setting the density to the input value and
360  * setting the temperature to
361  * @f[
362  * T = \frac{P \overline W}{\hat R \rho}.
363  * @f]
364  *
365  * @param rho Density (kg/m^3)
366  * @param p Pressure (Pa)
367  * @since New in %Cantera 3.0.
368  */
369  void setState_DP(double rho, double p) override {
370  if (p <= 0) {
371  throw CanteraError("IdealGasPhase::setState_DP",
372  "pressure must be positive");
373  }
374  setDensity(rho);
376  }
377 
378  //! Returns the isothermal compressibility. Units: 1/Pa.
379  /**
380  * The isothermal compressibility is defined as
381  * @f[
382  * \kappa_T = -\frac{1}{v}\left(\frac{\partial v}{\partial P}\right)_T
383  * @f]
384  * For ideal gases it's equal to the inverse of the pressure
385  */
386  double isothermalCompressibility() const override {
387  return 1.0 / pressure();
388  }
389 
390  //! Return the volumetric thermal expansion coefficient. Units: 1/K.
391  /*!
392  * The thermal expansion coefficient is defined as
393  * @f[
394  * \beta = \frac{1}{v}\left(\frac{\partial v}{\partial T}\right)_P
395  * @f]
396  * For ideal gases, it's equal to the inverse of the temperature.
397  */
398  double thermalExpansionCoeff() const override {
399  return 1.0 / temperature();
400  }
401 
402  double soundSpeed() const override;
403 
404  //! @}
405  //! @name Chemical Potentials and Activities
406  //!
407  //! The activity @f$ a_k @f$ of a species in solution is
408  //! related to the chemical potential by
409  //! @f[
410  //! \mu_k(T,P,X_k) = \mu_k^0(T,P) + \hat R T \ln a_k.
411  //! @f]
412  //! The quantity @f$ \mu_k^0(T,P) @f$ is the standard state chemical potential
413  //! at unit activity. It may depend on the pressure and the temperature.
414  //! However, it may not depend on the mole fractions of the species in the
415  //! solution.
416  //!
417  //! The activities are related to the generalized concentrations, @f$ \tilde
418  //! C_k @f$, and standard concentrations, @f$ C^0_k @f$, by the following
419  //! formula:
420  //!
421  //! @f[
422  //! a_k = \frac{\tilde C_k}{C^0_k}
423  //! @f]
424  //! The generalized concentrations are used in the kinetics classes to
425  //! describe the rates of progress of reactions involving the species. Their
426  //! formulation depends upon the specification of the rate constants for
427  //! reaction, especially the units used in specifying the rate constants. The
428  //! bridge between the thermodynamic equilibrium expressions that use a_k and
429  //! the kinetics expressions which use the generalized concentrations is
430  //! provided by the multiplicative factor of the standard concentrations.
431  //! @{
432 
433  //! This method returns the array of generalized concentrations.
434  /*!
435  * For an ideal gas mixture, these are simply the actual concentrations.
436  *
437  * @param c Output array of generalized concentrations. The units depend
438  * upon the implementation of the reaction rate expressions within
439  * the phase.
440  */
441  void getActivityConcentrations(double* c) const override {
443  }
444 
445  //! Returns the standard concentration @f$ C^0_k @f$, which is used to
446  //! normalize the generalized concentration.
447  /*!
448  * This is defined as the concentration by which the generalized
449  * concentration is normalized to produce the activity. In many cases, this
450  * quantity will be the same for all species in a phase. Since the activity
451  * for an ideal gas mixture is simply the mole fraction, for an ideal gas
452  * @f$ C^0_k = P/\hat R T @f$.
453  *
454  * @param k Optional parameter indicating the species. The default
455  * is to assume this refers to species 0.
456  * @return
457  * Returns the standard Concentration in units of m3 kmol-1.
458  */
459  double standardConcentration(size_t k=0) const override;
460 
461  //! Get the array of non-dimensional activity coefficients at the current
462  //! solution temperature, pressure, and solution concentration.
463  /*!
464  * For ideal gases, the activity coefficients are all equal to one.
465  *
466  * @param ac Output vector of activity coefficients. Length: m_kk.
467  */
468  void getActivityCoefficients(double* ac) const override;
469 
470  //! @}
471  //! @name Partial Molar Properties of the Solution
472  //! @{
473 
474  void getChemPotentials(double* mu) const override;
475  void getPartialMolarEnthalpies(double* hbar) const override;
476  void getPartialMolarEntropies(double* sbar) const override;
477  void getPartialMolarIntEnergies(double* ubar) const override;
478  void getPartialMolarCp(double* cpbar) const override;
479  void getPartialMolarVolumes(double* vbar) const override;
480 
481  //! @}
482  //! @name Properties of the Standard State of the Species in the Solution
483  //! @{
484 
485  void getStandardChemPotentials(double* mu) const override;
486  void getEnthalpy_RT(double* hrt) const override;
487  void getEntropy_R(double* sr) const override;
488  void getGibbs_RT(double* grt) const override;
489  void getPureGibbs(double* gpure) const override;
490  void getIntEnergy_RT(double* urt) const override;
491  void getCp_R(double* cpr) const override;
492  void getStandardVolumes(double* vol) const override;
493 
494  //! @}
495  //! @name Thermodynamic Values for the Species Reference States
496  //! @{
497 
498  void getEnthalpy_RT_ref(double* hrt) const override;
499  void getGibbs_RT_ref(double* grt) const override;
500  void getGibbs_ref(double* g) const override;
501  void getEntropy_R_ref(double* er) const override;
502  void getIntEnergy_RT_ref(double* urt) const override;
503  void getCp_R_ref(double* cprt) const override;
504  void getStandardVolumes_ref(double* vol) const override;
505 
506  //! @}
507  //! @name NonVirtual Internal methods to Return References to Reference State Thermo
508  //! @{
509 
510  //! Returns a reference to the dimensionless reference state enthalpy vector.
511  /*!
512  * This function is part of the layer that checks/recalculates the reference
513  * state thermo functions.
514  */
515  const vector<double>& enthalpy_RT_ref() const {
516  updateThermo();
517  return m_h0_RT;
518  }
519 
520  //! Returns a reference to the dimensionless reference state Gibbs free energy vector.
521  /*!
522  * This function is part of the layer that checks/recalculates the reference
523  * state thermo functions.
524  */
525  const vector<double>& gibbs_RT_ref() const {
526  updateThermo();
527  return m_g0_RT;
528  }
529 
530  //! Returns a reference to the dimensionless reference state Entropy vector.
531  /*!
532  * This function is part of the layer that checks/recalculates the reference
533  * state thermo functions.
534  */
535  const vector<double>& entropy_R_ref() const {
536  updateThermo();
537  return m_s0_R;
538  }
539 
540  //! Returns a reference to the dimensionless reference state Heat Capacity vector.
541  /*!
542  * This function is part of the layer that checks/recalculates the reference
543  * state thermo functions.
544  */
545  const vector<double>& cp_R_ref() const {
546  updateThermo();
547  return m_cp0_R;
548  }
549 
550  //! @}
551 
552  bool addSpecies(shared_ptr<Species> spec) override;
553  void setToEquilState(const double* mu_RT) override;
554 
555 protected:
556  //! Reference state pressure
557  /*!
558  * Value of the reference state pressure in Pascals.
559  * All species must have the same reference state pressure.
560  */
561  double m_p0 = -1.0;
562 
563  //! Temporary storage for dimensionless reference state enthalpies
564  mutable vector<double> m_h0_RT;
565 
566  //! Temporary storage for dimensionless reference state heat capacities
567  mutable vector<double> m_cp0_R;
568 
569  //! Temporary storage for dimensionless reference state Gibbs energies
570  mutable vector<double> m_g0_RT;
571 
572  //! Temporary storage for dimensionless reference state entropies
573  mutable vector<double> m_s0_R;
574 
575  mutable vector<double> m_expg0_RT;
576 
577  //! Temporary array containing internally calculated partial pressures
578  mutable vector<double> m_pp;
579 
580  //! Update the species reference state thermodynamic functions
581  /*!
582  * This method is called each time a thermodynamic property is requested,
583  * to check whether the internal species properties within the object
584  * need to be updated. Currently, this updates the species thermo
585  * polynomial values for the current value of the temperature. A check is
586  * made to see if the temperature has changed since the last evaluation.
587  * This object does not contain any persistent data that depends on the
588  * concentration, that needs to be updated. The state object modifies its
589  * concentration dependent information at the time the setMoleFractions()
590  * (or equivalent) call is made.
591  */
592  virtual void updateThermo() const;
593 };
594 
595 }
596 
597 #endif
Header file for class ThermoPhase, the base class for phases with thermodynamic properties,...
Base class for exceptions thrown by Cantera classes.
Definition: ctexceptions.h:66
Class IdealGasPhase represents low-density gases that obey the ideal gas equation of state.
double enthalpy_mole() const override
Return the Molar enthalpy. Units: J/kmol.
double thermalExpansionCoeff() const override
Return the volumetric thermal expansion coefficient. Units: 1/K.
bool isIdeal() const override
Boolean indicating whether phase is ideal.
void getPartialMolarEnthalpies(double *hbar) const override
Returns an array of partial molar enthalpies for the species in the mixture.
void getChemPotentials(double *mu) const override
Get the species chemical potentials. Units: J/kmol.
void setState_DP(double rho, double p) override
Set the density and pressure at constant composition.
double m_p0
Reference state pressure.
double soundSpeed() const override
Return the speed of sound. Units: m/s.
double pressure() const override
Pressure.
vector< double > m_g0_RT
Temporary storage for dimensionless reference state Gibbs energies.
const vector< double > & cp_R_ref() const
Returns a reference to the dimensionless reference state Heat Capacity vector.
const vector< double > & gibbs_RT_ref() const
Returns a reference to the dimensionless reference state Gibbs free energy vector.
void getEntropy_R(double *sr) const override
Get the array of nondimensional Entropy functions for the standard state species at the current T and...
vector< double > m_h0_RT
Temporary storage for dimensionless reference state enthalpies.
void getGibbs_ref(double *g) const override
Returns the vector of the Gibbs function of the reference state at the current temperature of the sol...
vector< double > m_pp
Temporary array containing internally calculated partial pressures.
void getStandardChemPotentials(double *mu) const override
Get the array of chemical potentials at unit activity for the species at their standard states at the...
void getCp_R(double *cpr) const override
Get the nondimensional Heat Capacities at constant pressure for the species standard states at the cu...
string type() const override
String indicating the thermodynamic model implemented.
void getActivityConcentrations(double *c) const override
This method returns the array of generalized concentrations.
void setPressure(double p) override
Set the pressure at constant temperature and composition.
void getStandardVolumes_ref(double *vol) const override
Get the molar volumes of the species reference states at the current T and P_ref of the solution.
void getPartialMolarVolumes(double *vbar) const override
Return an array of partial molar volumes for the species in the mixture.
double cv_mole() const override
Molar heat capacity at constant volume.
virtual void updateThermo() const
Update the species reference state thermodynamic functions.
void getPureGibbs(double *gpure) const override
Get the Gibbs functions for the standard state of the species at the current T and P of the solution.
void getIntEnergy_RT_ref(double *urt) const override
Returns the vector of nondimensional internal Energies of the reference state at the current temperat...
void getEnthalpy_RT(double *hrt) const override
Get the nondimensional Enthalpy functions for the species at their standard states at the current T a...
void getEntropy_R_ref(double *er) const override
Returns the vector of nondimensional entropies of the reference state at the current temperature of t...
vector< double > m_s0_R
Temporary storage for dimensionless reference state entropies.
double isothermalCompressibility() const override
Returns the isothermal compressibility. Units: 1/Pa.
void getGibbs_RT(double *grt) const override
Get the nondimensional Gibbs functions for the species in their standard states at the current T and ...
double entropy_mole() const override
Molar entropy.
const vector< double > & enthalpy_RT_ref() const
Returns a reference to the dimensionless reference state enthalpy vector.
void getCp_R_ref(double *cprt) const override
Returns the vector of nondimensional constant pressure heat capacities of the reference state at the ...
void getStandardVolumes(double *vol) const override
Get the molar volumes of the species standard states at the current T and P of the solution.
void getPartialMolarIntEnergies(double *ubar) const override
Return an array of partial molar internal energies for the species in the mixture.
double cp_mole() const override
Molar heat capacity at constant pressure.
void getIntEnergy_RT(double *urt) const override
Returns the vector of nondimensional Internal Energies of the standard state species at the current T...
void getPartialMolarCp(double *cpbar) const override
Return an array of partial molar heat capacities for the species in the mixture.
double standardConcentration(size_t k=0) const override
Returns the standard concentration , which is used to normalize the generalized concentration.
vector< double > m_cp0_R
Temporary storage for dimensionless reference state heat capacities.
const vector< double > & entropy_R_ref() const
Returns a reference to the dimensionless reference state Entropy vector.
bool addSpecies(shared_ptr< Species > spec) override
Add a Species to this Phase.
void setToEquilState(const double *mu_RT) override
This method is used by the ChemEquil equilibrium solver.
void getGibbs_RT_ref(double *grt) const override
Returns the vector of nondimensional Gibbs Free Energies of the reference state at the current temper...
void getActivityCoefficients(double *ac) const override
Get the array of non-dimensional activity coefficients at the current solution temperature,...
IdealGasPhase(const string &inputFile="", const string &id="")
Construct and initialize an IdealGasPhase ThermoPhase object directly from an input file.
void getPartialMolarEntropies(double *sbar) const override
Returns an array of partial molar entropies of the species in the solution.
string phaseOfMatter() const override
String indicating the mechanical phase of the matter in this Phase.
void getEnthalpy_RT_ref(double *hrt) const override
Returns the vector of nondimensional enthalpies of the reference state at the current temperature of ...
virtual void getConcentrations(double *const c) const
Get the species concentrations (kmol/m^3).
Definition: Phase.cpp:482
virtual double molarDensity() const
Molar density (kmol/m^3).
Definition: Phase.cpp:576
double temperature() const
Temperature (K).
Definition: Phase.h:562
double meanMolecularWeight() const
The mean molecular weight. Units: (kg/kmol)
Definition: Phase.h:655
virtual void setDensity(const double density_)
Set the internally stored density (kg/m^3) of the phase.
Definition: Phase.cpp:586
virtual void setTemperature(double temp)
Set the internally stored temperature of the phase (K).
Definition: Phase.h:623
double mean_X(const double *const Q) const
Evaluate the mole-fraction-weighted mean of an array Q.
Definition: Phase.cpp:616
Base class for a phase with thermodynamic properties.
Definition: ThermoPhase.h:390
double RT() const
Return the Gas Constant multiplied by the current temperature.
Definition: ThermoPhase.h:1062
const double GasConstant
Universal Gas Constant [J/kmol/K].
Definition: ct_defs.h:120
Namespace for the Cantera kernel.
Definition: AnyMap.cpp:564