Loading [MathJax]/extensions/tex2jax.js
Cantera  3.2.0a1
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules Pages
LatticePhase.h
Go to the documentation of this file.
1/**
2 * @file LatticePhase.h Header for a simple thermodynamics model of a bulk
3 * phase derived from ThermoPhase, assuming a lattice of solid atoms (see
4 * @ref thermoprops and class @link Cantera::LatticePhase
5 * LatticePhase@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_LATTICE_H
12#define CT_LATTICE_H
13
14#include "ThermoPhase.h"
15
16namespace Cantera
17{
18
19//! A simple thermodynamic model for a bulk phase, assuming a lattice of solid
20//! atoms
21/*!
22 * The bulk consists of a matrix of equivalent sites whose molar density does
23 * not vary with temperature or pressure. The thermodynamics obeys the ideal
24 * solution laws. The phase and the pure species phases which comprise the
25 * standard states of the species are assumed to have zero volume expansivity
26 * and zero isothermal compressibility.
27 *
28 * The density of matrix sites is given by the variable @f$ C_o @f$, which has
29 * SI units of kmol m-3.
30 *
31 * ## Specification of Species Standard State Properties
32 *
33 * It is assumed that the reference state thermodynamics may be obtained by a
34 * pointer to a populated species thermodynamic property manager class (see
35 * ThermoPhase::m_spthermo). However, how to relate pressure changes to the
36 * reference state thermodynamics is within this class.
37 *
38 * Pressure is defined as an independent variable in this phase. However, it has
39 * no effect on any quantities, as the molar concentration is a constant.
40 *
41 * The standard state enthalpy function is given by the following relation,
42 * which has a weak dependence on the system pressure, @f$ P @f$.
43 *
44 * @f[
45 * h^o_k(T,P) =
46 * h^{ref}_k(T) + \left( \frac{P - P_{ref}}{C_o} \right)
47 * @f]
48 *
49 * For an incompressible substance, the molar internal energy is independent of
50 * pressure. Since the thermodynamic properties are specified by giving the
51 * standard-state enthalpy, the term @f$ \frac{P_{ref}}{C_o} @f$ is subtracted
52 * from the specified reference molar enthalpy to compute the standard state
53 * molar internal energy:
54 *
55 * @f[
56 * u^o_k(T,P) = h^{ref}_k(T) - \frac{P_{ref}}{C_o}
57 * @f]
58 *
59 * The standard state heat capacity, internal energy, and entropy are
60 * independent of pressure. The standard state Gibbs free energy is obtained
61 * from the enthalpy and entropy functions.
62 *
63 * The standard state molar volume is independent of temperature, pressure, and
64 * species identity:
65 *
66 * @f[
67 * V^o_k(T,P) = \frac{1.0}{C_o}
68 * @f]
69 *
70 * ## Specification of Solution Thermodynamic Properties
71 *
72 * The activity of species @f$ k @f$ defined in the phase, @f$ a_k @f$, is given
73 * by the ideal solution law:
74 *
75 * @f[
76 * a_k = X_k ,
77 * @f]
78 *
79 * where @f$ X_k @f$ is the mole fraction of species *k*. The chemical potential
80 * for species *k* is equal to
81 *
82 * @f[
83 * \mu_k(T,P) = \mu^o_k(T, P) + R T \ln X_k
84 * @f]
85 *
86 * The partial molar entropy for species *k* is given by the following relation,
87 *
88 * @f[
89 * \tilde{s}_k(T,P) = s^o_k(T,P) - R \ln X_k = s^{ref}_k(T) - R \ln X_k
90 * @f]
91 *
92 * The partial molar enthalpy for species *k* is
93 *
94 * @f[
95 * \tilde{h}_k(T,P) = h^o_k(T,P) = h^{ref}_k(T) + \left( \frac{P - P_{ref}}{C_o} \right)
96 * @f]
97 *
98 * The partial molar Internal Energy for species *k* is
99 *
100 * @f[
101 * \tilde{u}_k(T,P) = u^o_k(T,P) = u^{ref}_k(T)
102 * @f]
103 *
104 * The partial molar Heat Capacity for species *k* is
105 *
106 * @f[
107 * \tilde{Cp}_k(T,P) = Cp^o_k(T,P) = Cp^{ref}_k(T)
108 * @f]
109 *
110 * The partial molar volume is independent of temperature, pressure, and species
111 * identity:
112 *
113 * @f[
114 * \tilde{V}_k(T,P) = V^o_k(T,P) = \frac{1.0}{C_o}
115 * @f]
116 *
117 * It is assumed that the reference state thermodynamics may be obtained by a
118 * pointer to a populated species thermodynamic property manager class (see
119 * ThermoPhase::m_spthermo). How to relate pressure changes to the reference
120 * state thermodynamics is resolved at this level.
121 *
122 * Pressure is defined as an independent variable in this phase. However, it
123 * only has a weak dependence on the enthalpy, and doesn't effect the molar
124 * concentration.
125 *
126 * ## Application within Kinetics Managers
127 *
128 * @f$ C^a_k @f$ are defined such that @f$ C^a_k = a_k = X_k @f$. @f$ C^s_k @f$,
129 * the standard concentration, is defined to be equal to one. @f$ a_k @f$ are
130 * activities used in the thermodynamic functions. These activity (or
131 * generalized) concentrations are used by kinetics manager classes to compute
132 * the forward and reverse rates of elementary reactions. The activity
133 * concentration,@f$ C^a_k @f$, is given by the following expression.
134 *
135 * @f[
136 * C^a_k = C^s_k X_k = X_k
137 * @f]
138 *
139 * The standard concentration for species *k* is identically one
140 *
141 * @f[
142 * C^s_k = C^s = 1.0
143 * @f]
144 *
145 * For example, a bulk-phase binary gas reaction between species j and k,
146 * producing a new species l would have the following equation for its rate of
147 * progress variable, @f$ R^1 @f$, which has units of kmol m-3 s-1.
148 *
149 * @f[
150 * R^1 = k^1 C_j^a C_k^a = k^1 X_j X_k
151 * @f]
152 *
153 * The reverse rate constant can then be obtained from the law of microscopic
154 * reversibility and the equilibrium expression for the system.
155 *
156 * @f[
157 * \frac{X_j X_k}{ X_l} = K_a^{o,1} = \exp(\frac{\mu^o_l - \mu^o_j - \mu^o_k}{R T} )
158 * @f]
159 *
160 * @f$ K_a^{o,1} @f$ is the dimensionless form of the equilibrium constant,
161 * associated with the pressure dependent standard states @f$ \mu^o_l(T,P) @f$
162 * and their associated activities,
163 * @f$ a_l @f$, repeated here:
164 *
165 * @f[
166 * \mu_l(T,P) = \mu^o_l(T, P) + R T \ln a_l
167 * @f]
168 *
169 * The concentration equilibrium constant, @f$ K_c @f$, may be obtained by
170 * changing over to activity concentrations. When this is done:
171 *
172 * @f[
173 * \frac{C^a_j C^a_k}{ C^a_l} = C^o K_a^{o,1} = K_c^1 =
174 * \exp(\frac{\mu^{o}_l - \mu^{o}_j - \mu^{o}_k}{R T} )
175 * @f]
176 *
177 * %Kinetics managers will calculate the concentration equilibrium constant, @f$
178 * K_c @f$, using the second and third part of the above expression as a
179 * definition for the concentration equilibrium constant.
180 *
181 * @ingroup thermoprops
182 */
184{
185public:
186 //! Full constructor for a lattice phase
187 /*!
188 * @param inputFile String name of the input file. If blank,
189 * an empty phase will be created.
190 * @param id string id of the phase name
191 */
192 explicit LatticePhase(const string& inputFile="", const string& id="");
193
194 string type() const override {
195 return "lattice";
196 }
197
198 bool isCompressible() const override {
199 return false;
200 }
201
202 map<string, size_t> nativeState() const override {
203 return { {"T", 0}, {"P", 1}, {"X", 2} };
204 }
205
206 //! @name Molar Thermodynamic Properties of the Solution
207 //! @{
208
209 //! Return the Molar Enthalpy. Units: J/kmol.
210 /*!
211 * For an ideal solution,
212 *
213 * @f[
214 * \hat h(T,P) = \sum_k X_k \hat h^0_k(T,P),
215 * @f]
216 *
217 * The standard-state pure-species Enthalpies @f$ \hat h^0_k(T,P) @f$ are
218 * computed first by the species reference state thermodynamic property
219 * manager and then a small pressure dependent term is added in.
220 *
221 * @see MultiSpeciesThermo
222 */
223 double enthalpy_mole() const override;
224
225 //! Molar entropy of the solution. Units: J/kmol/K
226 /*!
227 * For an ideal, constant partial molar volume solution mixture with
228 * pure species phases which exhibit zero volume expansivity:
229 * @f[
230 * \hat s(T, P, X_k) = \sum_k X_k \hat s^0_k(T) - \hat R \sum_k X_k \ln(X_k)
231 * @f]
232 * The reference-state pure-species entropies @f$ \hat s^0_k(T,p_{ref}) @f$
233 * are computed by the species thermodynamic property manager. The pure
234 * species entropies are independent of pressure since the volume
235 * expansivities are equal to zero.
236 *
237 * Units: J/kmol/K.
238 *
239 * @see MultiSpeciesThermo
240 */
241 double entropy_mole() const override;
242
243 //! Molar Gibbs function of the solution. Units: J/kmol
244 /*!
245 * Uses thermodynamic property relations.
246 */
247 double gibbs_mole() const override;
248
249 //! Molar heat capacity at constant pressure of the solution.
250 //! Units: J/kmol/K.
251 /*!
252 * For an ideal, constant partial molar volume solution mixture with
253 * pure species phases which exhibit zero volume expansivity:
254 * @f[
255 * \hat c_p(T,P) = \sum_k X_k \hat c^0_{p,k}(T) .
256 * @f]
257 * The heat capacity is independent of pressure. The reference-state pure-
258 * species heat capacities @f$ \hat c^0_{p,k}(T) @f$ are computed by the
259 * species thermodynamic property manager.
260 *
261 * @see MultiSpeciesThermo
262 */
263 double cp_mole() const override;
264
265 //! Molar heat capacity at constant volume of the solution.
266 //! Units: J/kmol/K.
267 /*!
268 * For an ideal, constant partial molar volume solution mixture with
269 * pure species phases which exhibit zero volume expansivity:
270 * @f[
271 * \hat c_v(T,P) = \hat c_p(T,P)
272 * @f]
273 *
274 * The two heat capacities are equal.
275 */
276 double cv_mole() const override;
277
278 //! @}
279 //! @name Mechanical Equation of State Properties
280 //!
281 //! In this equation of state implementation, the density is a function only
282 //! of the mole fractions. Therefore, it can't be an independent variable.
283 //! Instead, the pressure is used as the independent variable. Functions
284 //! which try to set the thermodynamic state by calling setDensity() may
285 //! cause an exception to be thrown.
286 //! @{
287
288 //! Pressure. Units: Pa.
289 /*!
290 * For this incompressible system, we return the internally stored
291 * independent value of the pressure.
292 */
293 double pressure() const override {
294 return m_Pcurrent;
295 }
296
297 //! Set the internally stored pressure (Pa) at constant temperature and
298 //! composition
299 /*!
300 * This method sets the pressure within the object. The mass density is not
301 * a function of pressure.
302 *
303 * @param p Input Pressure (Pa)
304 */
305 void setPressure(double p) override;
306
307 //! Calculate the density of the mixture using the partial molar volumes and
308 //! mole fractions as input
309 /*!
310 * The formula for this is
311 *
312 * @f[
313 * \rho = \frac{\sum_k{X_k W_k}}{\sum_k{X_k V_k}}
314 * @f]
315 *
316 * where @f$ X_k @f$ are the mole fractions, @f$ W_k @f$ are the molecular
317 * weights, and @f$ V_k @f$ are the pure species molar volumes.
318 *
319 * Note, the basis behind this formula is that in an ideal solution the
320 * partial molar volumes are equal to the pure species molar volumes. We
321 * have additionally specified in this class that the pure species molar
322 * volumes are independent of temperature and pressure.
323 */
324 double calcDensity();
325
326 //! @}
327 //! @name Activities, Standard States, and Activity Concentrations
328 //!
329 //! The activity @f$ a_k @f$ of a species in solution is related to the
330 //! chemical potential by @f[ \mu_k = \mu_k^0(T) + \hat R T \ln a_k. @f] The
331 //! quantity @f$ \mu_k^0(T,P) @f$ is the chemical potential at unit activity,
332 //! which depends only on temperature and the pressure. Activity is assumed
333 //! to be molality-based here.
334 //! @{
335
336 Units standardConcentrationUnits() const override;
337 void getActivityConcentrations(double* c) const override;
338
339 //! Return the standard concentration for the kth species
340 /*!
341 * The standard concentration @f$ C^0_k @f$ used to normalize
342 * the activity (that is, generalized) concentration for use
343 *
344 * For the time being, we will use the concentration of pure solvent for the
345 * the standard concentration of all species. This has the effect of making
346 * mass-action reaction rates based on the molality of species proportional
347 * to the molality of the species.
348 *
349 * @param k Optional parameter indicating the species. The default is to
350 * assume this refers to species 0.
351 * @returns the standard Concentration in units of m^3/kmol.
352 */
353 double standardConcentration(size_t k=0) const override;
354 double logStandardConc(size_t k=0) const override;
355
356 //! Get the array of non-dimensional activity coefficients at
357 //! the current solution temperature, pressure, and solution concentration.
358 /*!
359 * For this phase, the activity coefficients are all equal to one.
360 *
361 * @param ac Output vector of activity coefficients. Length: m_kk.
362 */
363 void getActivityCoefficients(double* ac) const override;
364
365 //! @}
366 //! @name Partial Molar Properties of the Solution
367 //! @{
368
369 //! Get the species chemical potentials. Units: J/kmol.
370 /*!
371 * This function returns a vector of chemical potentials of the species in
372 * solid solution at the current temperature, pressure and mole fraction of
373 * the solid solution.
374 *
375 * @param mu Output vector of species chemical
376 * potentials. Length: m_kk. Units: J/kmol
377 */
378 void getChemPotentials(double* mu) const override;
379
380 /**
381 * Returns an array of partial molar enthalpies for the species in the
382 * mixture. Units (J/kmol). For this phase, the partial molar enthalpies are
383 * equal to the pure species enthalpies
384 * @f[
385 * \bar h_k(T,P) = \hat h^{ref}_k(T) + (P - P_{ref}) \hat V^0_k
386 * @f]
387 * The reference-state pure-species enthalpies, @f$ \hat h^{ref}_k(T) @f$,
388 * at the reference pressure,@f$ P_{ref} @f$, are computed by the species
389 * thermodynamic property manager. They are polynomial functions of
390 * temperature.
391 * @see MultiSpeciesThermo
392 *
393 * @param hbar Output vector containing partial molar enthalpies.
394 * Length: m_kk.
395 */
396 void getPartialMolarEnthalpies(double* hbar) const override;
397
398 /**
399 * Returns an array of partial molar entropies of the species in the
400 * solution. Units: J/kmol/K. For this phase, the partial molar entropies
401 * are equal to the pure species entropies plus the ideal solution
402 * contribution.
403 * @f[
404 * \bar s_k(T,P) = \hat s^0_k(T) - R \ln(X_k)
405 * @f]
406 * The reference-state pure-species entropies,@f$ \hat s^{ref}_k(T) @f$, at
407 * the reference pressure, @f$ P_{ref} @f$, are computed by the species
408 * thermodynamic property manager. They are polynomial functions of
409 * temperature.
410 * @see MultiSpeciesThermo
411 *
412 * @param sbar Output vector containing partial molar entropies.
413 * Length: m_kk.
414 */
415 void getPartialMolarEntropies(double* sbar) const override;
416
417 /**
418 * Returns an array of partial molar Heat Capacities at constant pressure of
419 * the species in the solution. Units: J/kmol/K. For this phase, the partial
420 * molar heat capacities are equal to the standard state heat capacities.
421 *
422 * @param cpbar Output vector of partial heat capacities. Length: m_kk.
423 */
424 void getPartialMolarCp(double* cpbar) const override;
425
426 void getPartialMolarVolumes(double* vbar) const override;
427 void getStandardChemPotentials(double* mu) const override;
428 void getPureGibbs(double* gpure) const override;
429
430 //! @}
431 //! @name Properties of the Standard State of the Species in the Solution
432 //! @{
433
434 //! Get the nondimensional Enthalpy functions for the species standard
435 //! states at their standard states at the current *T* and *P* of the
436 //! solution.
437 /*!
438 * A small pressure dependent term is added onto the reference state enthalpy
439 * to get the pressure dependence of this term.
440 *
441 * @f[
442 * h^o_k(T,P) = h^{ref}_k(T) + \left( \frac{P - P_{ref}}{C_o} \right)
443 * @f]
444 *
445 * The reference state thermodynamics is obtained by a pointer to a
446 * populated species thermodynamic property manager class (see
447 * ThermoPhase::m_spthermo). How to relate pressure changes to the reference
448 * state thermodynamics is resolved at this level.
449 *
450 * @param hrt Output vector of nondimensional standard state enthalpies.
451 * Length: m_kk.
452 */
453 void getEnthalpy_RT(double* hrt) const override;
454
455 //! Get the array of nondimensional Entropy functions for the species
456 //! standard states at the current *T* and *P* of the solution.
457 /*!
458 * The entropy of the standard state is defined as independent of
459 * pressure here.
460 *
461 * @f[
462 * s^o_k(T,P) = s^{ref}_k(T)
463 * @f]
464 *
465 * The reference state thermodynamics is obtained by a pointer to a
466 * populated species thermodynamic property manager class (see
467 * ThermoPhase::m_spthermo). How to relate pressure changes to the reference
468 * state thermodynamics is resolved at this level.
469 *
470 * @param sr Output vector of nondimensional standard state entropies.
471 * Length: m_kk.
472 */
473 void getEntropy_R(double* sr) const override;
474
475 //! Get the nondimensional Gibbs functions for the species standard states
476 //! at the current *T* and *P* of the solution.
477 /*!
478 * The standard Gibbs free energies are obtained from the enthalpy and
479 * entropy formulation.
480 *
481 * @f[
482 * g^o_k(T,P) = h^{o}_k(T,P) - T s^{o}_k(T,P)
483 * @f]
484 *
485 * @param grt Output vector of nondimensional standard state Gibbs free
486 * energies. Length: m_kk.
487 */
488 void getGibbs_RT(double* grt) const override;
489
490 //! Get the nondimensional Heat Capacities at constant pressure for the
491 //! species standard states at the current *T* and *P* of the solution
492 /*!
493 * The heat capacity of the standard state is independent of pressure
494 *
495 * @f[
496 * Cp^o_k(T,P) = Cp^{ref}_k(T)
497 * @f]
498 *
499 * The reference state thermodynamics is obtained by a pointer to a
500 * populated species thermodynamic property manager class (see
501 * ThermoPhase::m_spthermo). How to relate pressure changes to the reference
502 * state thermodynamics is resolved at this level.
503 *
504 * @param cpr Output vector of nondimensional standard state heat
505 * capacities. Length: m_kk.
506 */
507 void getCp_R(double* cpr) const override;
508
509 void getStandardVolumes(double* vol) const override;
510
511 //! @}
512 //! @name Thermodynamic Values for the Species Reference States
513 //! @{
514
515 const vector<double>& enthalpy_RT_ref() const;
516
517 //! Returns a reference to the dimensionless reference state Gibbs free
518 //! energy vector.
519 /*!
520 * This function is part of the layer that checks/recalculates the reference
521 * state thermo functions.
522 */
523 const vector<double>& gibbs_RT_ref() const;
524
525 void getGibbs_RT_ref(double* grt) const override;
526 void getGibbs_ref(double* g) const override;
527
528 //! Returns a reference to the dimensionless reference state Entropy vector.
529 /*!
530 * This function is part of the layer that checks/recalculates the reference
531 * state thermo functions.
532 */
533 const vector<double>& entropy_R_ref() const;
534
535 //! Returns a reference to the dimensionless reference state Heat Capacity
536 //! vector.
537 /*!
538 * This function is part of the layer that checks/recalculates the reference
539 * state thermo functions.
540 */
541 const vector<double>& cp_R_ref() const;
542
543 //! @}
544 //! @name Utilities for Initialization of the Object
545 //! @{
546
547 bool addSpecies(shared_ptr<Species> spec) override;
548
549 //! Set the density of lattice sites [kmol/m^3]
550 void setSiteDensity(double sitedens);
551
552 void initThermo() override;
553 void getParameters(AnyMap& phaseNode) const override;
554 void getSpeciesParameters(const string& name, AnyMap& speciesNode) const override;
555
556 //! @}
557
558protected:
559 void compositionChanged() override;
560
561 //! Reference state pressure
562 double m_Pref = OneAtm;
563
564 //! The current pressure
565 /*!
566 * Since the density isn't a function of pressure, but only of the mole
567 * fractions, we need to independently specify the pressure. The density
568 * variable which is inherited as part of the State class, m_dens, is always
569 * kept current whenever T, P, or X[] change.
570 */
572
573 //! Reference state enthalpies / RT
574 mutable vector<double> m_h0_RT;
575
576 //! Temporary storage for the reference state heat capacities
577 mutable vector<double> m_cp0_R;
578
579 //! Temporary storage for the reference state Gibbs energies
580 mutable vector<double> m_g0_RT;
581
582 //! Temporary storage for the reference state entropies at the current
583 //! temperature
584 mutable vector<double> m_s0_R;
585
586 //! Vector of molar volumes for each species in the solution
587 /**
588 * Species molar volumes @f$ m^3 kmol^-1 @f$
589 */
590 vector<double> m_speciesMolarVolume;
591
592 //! Site Density of the lattice solid
593 /*!
594 * Currently, this is imposed as a function of T, P or composition
595 *
596 * units are kmol m-3
597 */
598 double m_site_density = 0.0;
599
600private:
601 //! Update the species reference state thermodynamic functions
602 /*!
603 * The polynomials for the standard state functions are only reevaluated if
604 * the temperature has changed.
605 */
606 void _updateThermo() const;
607};
608}
609
610#endif
Header file for class ThermoPhase, the base class for phases with thermodynamic properties,...
A map of string keys to values whose type can vary at runtime.
Definition AnyMap.h:432
A simple thermodynamic model for a bulk phase, assuming a lattice of solid atoms.
const vector< double > & entropy_R_ref() const
Returns a reference to the dimensionless reference state Entropy vector.
double enthalpy_mole() const override
Return the Molar Enthalpy. Units: J/kmol.
double logStandardConc(size_t k=0) const override
Natural logarithm of the standard concentration of the kth species.
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.
double m_site_density
Site Density of the lattice solid.
double pressure() const override
Pressure. Units: Pa.
vector< double > m_g0_RT
Temporary storage for the reference state Gibbs energies.
bool isCompressible() const override
Return whether phase represents a compressible substance.
void getSpeciesParameters(const string &name, AnyMap &speciesNode) const override
Get phase-specific parameters of a Species object such that an identical one could be reconstructed a...
map< string, size_t > nativeState() const override
Return a map of properties defining the native state of a substance.
void getEntropy_R(double *sr) const override
Get the array of nondimensional Entropy functions for the species standard states at the current T an...
vector< double > m_h0_RT
Reference state enthalpies / RT.
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...
double m_Pref
Reference state pressure.
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...
void getParameters(AnyMap &phaseNode) const override
Store the parameters of a ThermoPhase object such that an identical one could be reconstructed using ...
string type() const override
String indicating the thermodynamic model implemented.
void initThermo() override
Initialize the ThermoPhase object after all species have been set up.
void getActivityConcentrations(double *c) const override
This method returns an array of generalized concentrations.
void setPressure(double p) override
Set the internally stored pressure (Pa) at constant temperature and composition.
const vector< double > & gibbs_RT_ref() const
Returns a reference to the dimensionless reference state Gibbs free energy vector.
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 of the solution.
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 getEnthalpy_RT(double *hrt) const override
Get the nondimensional Enthalpy functions for the species standard states at their standard states at...
vector< double > m_s0_R
Temporary storage for the reference state entropies at the current temperature.
void getGibbs_RT(double *grt) const override
Get the nondimensional Gibbs functions for the species standard states at the current T and P of the ...
double entropy_mole() const override
Molar entropy of the solution. Units: J/kmol/K.
void setSiteDensity(double sitedens)
Set the density of lattice sites [kmol/m^3].
vector< double > m_speciesMolarVolume
Vector of molar volumes for each species in the solution.
void getStandardVolumes(double *vol) const override
Get the molar volumes of the species standard states at the current T and P of the solution.
double cp_mole() const override
Molar heat capacity at constant pressure of the solution.
Units standardConcentrationUnits() const override
Returns the units of the "standard concentration" for this phase.
void getPartialMolarCp(double *cpbar) const override
Returns an array of partial molar Heat Capacities at constant pressure of the species in the solution...
void compositionChanged() override
Apply changes to the state which are needed after the composition changes.
double gibbs_mole() const override
Molar Gibbs function of the solution. Units: J/kmol.
double standardConcentration(size_t k=0) const override
Return the standard concentration for the kth species.
vector< double > m_cp0_R
Temporary storage for the reference state heat capacities.
double calcDensity()
Calculate the density of the mixture using the partial molar volumes and mole fractions as input.
bool addSpecies(shared_ptr< Species > spec) override
Add a Species to this Phase.
void _updateThermo() const
Update the species reference state thermodynamic functions.
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,...
void getPartialMolarEntropies(double *sbar) const override
Returns an array of partial molar entropies of the species in the solution.
const vector< double > & cp_R_ref() const
Returns a reference to the dimensionless reference state Heat Capacity vector.
double m_Pcurrent
The current pressure.
string name() const
Return the name of the phase.
Definition Phase.cpp:20
Base class for a phase with thermodynamic properties.
A representation of the units associated with a dimensional quantity.
Definition Units.h:35
const double OneAtm
One atmosphere [Pa].
Definition ct_defs.h:96
Namespace for the Cantera kernel.
Definition AnyMap.cpp:595