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