Loading [MathJax]/extensions/tex2jax.js
Cantera  3.2.0a1
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules Pages
IdealSolidSolnPhase.h
Go to the documentation of this file.
1/**
2 * @file IdealSolidSolnPhase.h Header file for an ideal solid solution model
3 * with incompressible thermodynamics (see @ref thermoprops and
4 * @link Cantera::IdealSolidSolnPhase IdealSolidSolnPhase@endlink).
5 *
6 * This class inherits from the %Cantera class ThermoPhase and implements an
7 * ideal solid solution model with incompressible thermodynamics.
8 */
9
10// This file is part of Cantera. See License.txt in the top-level directory or
11// at https://cantera.org/license.txt for license and copyright information.
12
13#ifndef CT_IDEALSOLIDSOLNPHASE_H
14#define CT_IDEALSOLIDSOLNPHASE_H
15
16#include "ThermoPhase.h"
17
18namespace Cantera
19{
20
21/**
22 * Class IdealSolidSolnPhase represents a condensed phase ideal solution
23 * compound. The phase and the pure species phases which comprise the standard
24 * states of the species are assumed to have zero volume expansivity and zero
25 * isothermal compressibility. Each species does, however, have constant but
26 * distinct partial molar volumes equal to their pure species molar volumes. The
27 * class derives from class ThermoPhase, and overloads the virtual methods
28 * defined there with ones that use expressions appropriate for ideal solution
29 * mixtures.
30 *
31 * The generalized concentrations can have three different forms depending on
32 * the value of the member attribute #m_formGC, which is supplied in the
33 * constructor and in the input file. The value and form of the generalized
34 * concentration will affect reaction rate constants involving species in this
35 * phase.
36 *
37 * @ingroup thermoprops
38 */
40{
41public:
42 //! Construct and initialize an IdealSolidSolnPhase ThermoPhase object
43 //! directly from an input file
44 /*!
45 * This constructor will also fully initialize the object.
46 * The generalized concentrations can have three different forms
47 * depending on the value of the member attribute #m_formGC, which
48 * is supplied in the constructor or read from the input file.
49 *
50 * @param infile File name for the input file containing information
51 * for this phase. If blank, an empty phase will be
52 * created.
53 * @param id The name of this phase. This is used to look up
54 * the phase in the input file.
55 */
56 explicit IdealSolidSolnPhase(const string& infile="", const string& id="");
57
58 string type() const override {
59 return "ideal-condensed";
60 }
61
62 bool isIdeal() const override {
63 return true;
64 }
65
66 bool isCompressible() const override {
67 return false;
68 }
69
70 //! @name Molar Thermodynamic Properties of the Solution
71 //! @{
72
73 /**
74 * Molar entropy of the solution. Units: J/kmol/K. For an ideal, constant
75 * partial molar volume solution mixture with pure species phases which
76 * exhibit zero volume expansivity:
77 * @f[
78 * \hat s(T, P, X_k) = \sum_k X_k \hat s^0_k(T) - \hat R \sum_k X_k \ln(X_k)
79 * @f]
80 * The reference-state pure-species entropies
81 * @f$ \hat s^0_k(T,p_{ref}) @f$ are computed by the species thermodynamic
82 * property manager. The pure species entropies are independent of
83 * pressure since the volume expansivities are equal to zero.
84 * @see MultiSpeciesThermo
85 */
86 double entropy_mole() const override;
87
88 /**
89 * Molar Gibbs free energy of the solution. Units: J/kmol. For an ideal,
90 * constant partial molar volume solution mixture with pure species phases
91 * which exhibit zero volume expansivity:
92 * @f[
93 * \hat g(T, P) = \sum_k X_k \hat g^0_k(T,P) + \hat R T \sum_k X_k \ln(X_k)
94 * @f]
95 * The reference-state pure-species Gibbs free energies
96 * @f$ \hat g^0_k(T) @f$ are computed by the species thermodynamic
97 * property manager, while the standard state Gibbs free energies
98 * @f$ \hat g^0_k(T,P) @f$ are computed by the member function, gibbs_RT().
99 * @see MultiSpeciesThermo
100 */
101 double gibbs_mole() const override;
102
103 /**
104 * Molar heat capacity at constant pressure of the solution.
105 * Units: J/kmol/K.
106 * For an ideal, constant partial molar volume solution mixture with
107 * pure species phases which exhibit zero volume expansivity:
108 * @f[
109 * \hat c_p(T,P) = \sum_k X_k \hat c^0_{p,k}(T) .
110 * @f]
111 * The heat capacity is independent of pressure. The reference-state pure-
112 * species heat capacities @f$ \hat c^0_{p,k}(T) @f$ are computed by the
113 * species thermodynamic property manager.
114 * @see MultiSpeciesThermo
115 */
116 double cp_mole() const override;
117
118 /**
119 * Molar heat capacity at constant volume of the solution. Units: J/kmol/K.
120 * For an ideal, constant partial molar volume solution mixture with pure
121 * species phases which exhibit zero volume expansivity:
122 * @f[ \hat c_v(T,P) = \hat c_p(T,P) @f]
123 * The two heat capacities are equal.
124 */
125 double cv_mole() const override {
126 return cp_mole();
127 }
128
129 //! @}
130 //! @name Mechanical Equation of State Properties
131 //!
132 //! In this equation of state implementation, the density is a function only
133 //! of the mole fractions. Therefore, it can't be an independent variable.
134 //! Instead, the pressure is used as the independent variable. Functions
135 //! which try to set the thermodynamic state by calling setDensity() will
136 //! cause an exception to be thrown.
137 //! @{
138
139 /**
140 * Pressure. Units: Pa. For this incompressible system, we return the
141 * internally stored independent value of the pressure.
142 */
143 double pressure() const override {
144 return m_Pcurrent;
145 }
146
147 /**
148 * Set the pressure at constant temperature. Units: Pa. This method sets a
149 * constant within the object. The mass density is not a function of
150 * pressure.
151 *
152 * @param p Input Pressure (Pa)
153 */
154 void setPressure(double p) override;
155
156 /**
157 * Calculate the density of the mixture using the partial molar volumes and
158 * mole fractions as input
159 *
160 * The formula for this is
161 *
162 * @f[
163 * \rho = \frac{\sum_k{X_k W_k}}{\sum_k{X_k V_k}}
164 * @f]
165 *
166 * where @f$ X_k @f$ are the mole fractions, @f$ W_k @f$ are the molecular
167 * weights, and @f$ V_k @f$ are the pure species molar volumes.
168 *
169 * Note, the basis behind this formula is that in an ideal solution the
170 * partial molar volumes are equal to the pure species molar volumes. We
171 * have additionally specified in this class that the pure species molar
172 * volumes are independent of temperature and pressure.
173 */
174 virtual void calcDensity();
175
176 //! @}
177 //! @name Chemical Potentials and Activities
178 //!
179 //! The activity @f$ a_k @f$ of a species in solution is related to the
180 //! chemical potential by
181 //! @f[
182 //! \mu_k(T,P,X_k) = \mu_k^0(T,P) + \hat R T \ln a_k.
183 //! @f]
184 //! The quantity @f$ \mu_k^0(T,P) @f$ is the standard state chemical potential
185 //! at unit activity. It may depend on the pressure and the temperature.
186 //! However, it may not depend on the mole fractions of the species in the
187 //! solid solution.
188 //!
189 //! The activities are related to the generalized concentrations, @f$ \tilde
190 //! C_k @f$, and standard concentrations, @f$ C^0_k @f$, by the following
191 //! formula:
192 //!
193 //! @f[
194 //! a_k = \frac{\tilde C_k}{C^0_k}
195 //! @f]
196 //! The generalized concentrations are used in the kinetics classes to
197 //! describe the rates of progress of reactions involving the species. Their
198 //! formulation depends upon the specification of the rate constants for
199 //! reaction, especially the units used in specifying the rate constants. The
200 //! bridge between the thermodynamic equilibrium expressions that use a_k and
201 //! the kinetics expressions which use the generalized concentrations is
202 //! provided by the multiplicative factor of the standard concentrations.
203 //! @{
204
205 Units standardConcentrationUnits() const override;
206
207 /**
208 * This method returns the array of generalized concentrations. The
209 * generalized concentrations are used in the evaluation of the rates of
210 * progress for reactions involving species in this phase. The generalized
211 * concentration divided by the standard concentration is also equal to the
212 * activity of species.
213 *
214 * For this implementation the activity is defined to be the mole fraction
215 * of the species. The generalized concentration is defined to be equal to
216 * the mole fraction divided by the partial molar volume. The generalized
217 * concentrations for species in this phase therefore have units of
218 * kmol/m^3. Rate constants must reflect this fact.
219 *
220 * On a general note, the following must be true. For an ideal solution, the
221 * generalized concentration must consist of the mole fraction multiplied by
222 * a constant. The constant may be fairly arbitrarily chosen, with
223 * differences adsorbed into the reaction rate expression. 1/V_N, 1/V_k, or
224 * 1 are equally good, as long as the standard concentration is adjusted
225 * accordingly. However, it must be a constant (and not the concentration,
226 * btw, which is a function of the mole fractions) in order for the ideal
227 * solution properties to hold at the same time having the standard
228 * concentration to be independent of the mole fractions.
229 *
230 * In this implementation the form of the generalized concentrations
231 * depend upon the member attribute, #m_formGC.
232 *
233 * HKM Note: We have absorbed the pressure dependence of the pure species
234 * state into the thermodynamics functions. Therefore the standard
235 * state on which the activities are based depend on both temperature
236 * and pressure. If we hadn't, it would have appeared in this
237 * function in a very awkward exp[] format.
238 *
239 * @param c Pointer to array of doubles of length m_kk, which on exit
240 * will contain the generalized concentrations.
241 */
242 void getActivityConcentrations(double* c) const override;
243
244 /**
245 * The standard concentration @f$ C^0_k @f$ used to normalize the
246 * generalized concentration. In many cases, this quantity will be the
247 * same for all species in a phase. However, for this case, we will return
248 * a distinct concentration for each species. This is the inverse of the
249 * species molar volume. Units for the standard concentration are kmol/m^3.
250 *
251 * @param k Species number: this is a require parameter, a change from the
252 * ThermoPhase base class, where it was an optional parameter.
253 */
254 double standardConcentration(size_t k) const override;
255
256 //! Get the array of species activity coefficients
257 /*!
258 * @param ac output vector of activity coefficients. Length: m_kk
259 */
260 void getActivityCoefficients(double* ac) const override;
261
262 /**
263 * Get the species chemical potentials. Units: J/kmol.
264 *
265 * This function returns a vector of chemical potentials of the
266 * species in solution.
267 * @f[
268 * \mu_k = \mu^{ref}_k(T) + V_k * (p - p_o) + R T \ln(X_k)
269 * @f]
270 * or another way to phrase this is
271 * @f[
272 * \mu_k = \mu^o_k(T,p) + R T \ln(X_k)
273 * @f]
274 * where @f$ \mu^o_k(T,p) = \mu^{ref}_k(T) + V_k * (p - p_o) @f$
275 *
276 * @param mu Output vector of chemical potentials.
277 */
278 void getChemPotentials(double* mu) const override;
279
280 //! @}
281 //! @name Partial Molar Properties of the Solution
282 //! @{
283
284 //! Returns an array of partial molar enthalpies for the species in the
285 //! mixture.
286 /*!
287 * Units (J/kmol). For this phase, the partial molar enthalpies are equal to
288 * the pure species enthalpies
289 * @f[
290 * \bar h_k(T,P) = \hat h^{ref}_k(T) + (P - P_{ref}) \hat V^0_k
291 * @f]
292 * The reference-state pure-species enthalpies, @f$ \hat h^{ref}_k(T) @f$,
293 * at the reference pressure,@f$ P_{ref} @f$, are computed by the species
294 * thermodynamic property manager. They are polynomial functions of
295 * temperature.
296 * @see MultiSpeciesThermo
297 *
298 * @param hbar Output vector containing partial molar enthalpies.
299 * Length: m_kk.
300 */
301 void getPartialMolarEnthalpies(double* hbar) const override;
302
303 /**
304 * Returns an array of partial molar entropies of the species in the
305 * solution. Units: J/kmol/K. For this phase, the partial molar entropies
306 * are equal to the pure species entropies plus the ideal solution
307 * contribution.
308 * @f[
309 * \bar s_k(T,P) = \hat s^0_k(T) - R \ln(X_k)
310 * @f]
311 * The reference-state pure-species entropies,@f$ \hat s^{ref}_k(T) @f$, at
312 * the reference pressure, @f$ P_{ref} @f$, are computed by the species
313 * thermodynamic property manager. They are polynomial functions of
314 * temperature.
315 * @see MultiSpeciesThermo
316 *
317 * @param sbar Output vector containing partial molar entropies.
318 * Length: m_kk.
319 */
320 void getPartialMolarEntropies(double* sbar) const override;
321
322 /**
323 * Returns an array of partial molar Heat Capacities at constant pressure of
324 * the species in the solution. Units: J/kmol/K. For this phase, the partial
325 * molar heat capacities are equal to the standard state heat capacities.
326 *
327 * @param cpbar Output vector of partial heat capacities. Length: m_kk.
328 */
329 void getPartialMolarCp(double* cpbar) const override;
330
331 /**
332 * returns an array of partial molar volumes of the species
333 * in the solution. Units: m^3 kmol-1.
334 *
335 * For this solution, the partial molar volumes are equal to the
336 * constant species molar volumes.
337 *
338 * @param vbar Output vector of partial molar volumes. Length: m_kk.
339 */
340 void getPartialMolarVolumes(double* vbar) const override;
341
342 //! @}
343 //! @name Properties of the Standard State of the Species in the Solution
344 //! @{
345
346 /**
347 * Get the standard state chemical potentials of the species. This is the
348 * array of chemical potentials at unit activity @f$ \mu^0_k(T,P) @f$. We
349 * define these here as the chemical potentials of the pure species at the
350 * temperature and pressure of the solution. This function is used in the
351 * evaluation of the equilibrium constant Kc. Therefore, Kc will also depend
352 * on T and P. This is the norm for liquid and solid systems.
353 *
354 * units = J / kmol
355 *
356 * @param mu0 Output vector of standard state chemical potentials.
357 * Length: m_kk.
358 */
359 void getStandardChemPotentials(double* mu0) const override {
360 getPureGibbs(mu0);
361 }
362
363 //! Get the array of nondimensional Enthalpy functions for the standard
364 //! state species at the current *T* and *P* of the solution.
365 /*!
366 * We assume an incompressible constant partial molar volume here:
367 * @f[
368 * h^0_k(T,P) = h^{ref}_k(T) + (P - P_{ref}) * V_k
369 * @f]
370 * where @f$ V_k @f$ is the molar volume of pure species *k*.
371 * @f$ h^{ref}_k(T) @f$ is the enthalpy of the pure species *k* at the
372 * reference pressure, @f$ P_{ref} @f$.
373 *
374 * @param hrt Vector of length m_kk, which on return hrt[k] will contain the
375 * nondimensional standard state enthalpy of species k.
376 */
377 void getEnthalpy_RT(double* hrt) const override;
378
379 //! Get the nondimensional Entropies for the species standard states at the
380 //! current T and P of the solution.
381 /*!
382 * Note, this is equal to the reference state entropies due to the zero
383 * volume expansivity: that is, (dS/dP)_T = (dV/dT)_P = 0.0
384 *
385 * @param sr Vector of length m_kk, which on return sr[k] will contain the
386 * nondimensional standard state entropy for species k.
387 */
388 void getEntropy_R(double* sr) const override;
389
390 /**
391 * Get the nondimensional Gibbs function for the species standard states at
392 * the current T and P of the solution.
393 *
394 * @f[
395 * \mu^0_k(T,P) = \mu^{ref}_k(T) + (P - P_{ref}) * V_k
396 * @f]
397 * where @f$ V_k @f$ is the molar volume of pure species *k*.
398 * @f$ \mu^{ref}_k(T) @f$ is the chemical potential of pure species *k*
399 * at the reference pressure, @f$ P_{ref} @f$.
400 *
401 * @param grt Vector of length m_kk, which on return sr[k] will contain the
402 * nondimensional standard state Gibbs function for species k.
403 */
404 void getGibbs_RT(double* grt) const override;
405
406 /**
407 * Get the Gibbs functions for the pure species at the current *T* and *P*
408 * of the solution. We assume an incompressible constant partial molar
409 * volume here:
410 * @f[
411 * \mu^0_k(T,P) = \mu^{ref}_k(T) + (P - P_{ref}) * V_k
412 * @f]
413 * where @f$ V_k @f$ is the molar volume of pure species *k*.
414 * @f$ \mu^{ref}_k(T) @f$ is the chemical potential of pure species *k* at
415 * the reference pressure, @f$ P_{ref} @f$.
416 *
417 * @param gpure Output vector of Gibbs functions for species. Length: m_kk.
418 */
419 void getPureGibbs(double* gpure) const override;
420
421 void getIntEnergy_RT(double* urt) const override;
422
423 /**
424 * Get the nondimensional heat capacity at constant pressure function for
425 * the species standard states at the current T and P of the solution.
426 * @f[
427 * Cp^0_k(T,P) = Cp^{ref}_k(T)
428 * @f]
429 * where @f$ V_k @f$ is the molar volume of pure species *k*.
430 * @f$ Cp^{ref}_k(T) @f$ is the constant pressure heat capacity of species
431 * *k* at the reference pressure, @f$ p_{ref} @f$.
432 *
433 * @param cpr Vector of length m_kk, which on return cpr[k] will contain the
434 * nondimensional constant pressure heat capacity for species k.
435 */
436 void getCp_R(double* cpr) const override;
437
438 void getStandardVolumes(double* vol) const override;
439
440 //! @}
441 //! @name Thermodynamic Values for the Species Reference States
442 //! @{
443
444 void getEnthalpy_RT_ref(double* hrt) const override;
445 void getGibbs_RT_ref(double* grt) const override;
446 void getGibbs_ref(double* g) const override;
447 void getEntropy_R_ref(double* er) const override;
448 void getIntEnergy_RT_ref(double* urt) const override;
449 void getCp_R_ref(double* cprt) const override;
450
451 /**
452 * Returns a reference to the vector of nondimensional enthalpies of the
453 * reference state at the current temperature. Real reason for its existence
454 * is that it also checks to see if a recalculation of the reference
455 * thermodynamics functions needs to be done.
456 */
457 const vector<double>& enthalpy_RT_ref() const;
458
459 /**
460 * Returns a reference to the vector of nondimensional enthalpies of the
461 * reference state at the current temperature. Real reason for its existence
462 * is that it also checks to see if a recalculation of the reference
463 * thermodynamics functions needs to be done.
464 */
465 const vector<double>& gibbs_RT_ref() const {
467 return m_g0_RT;
468 }
469
470 /**
471 * Returns a reference to the vector of nondimensional enthalpies of the
472 * reference state at the current temperature. Real reason for its existence
473 * is that it also checks to see if a recalculation of the reference
474 * thermodynamics functions needs to be done.
475 */
476 const vector<double>& entropy_R_ref() const;
477
478 /**
479 * Returns a reference to the vector of nondimensional enthalpies of the
480 * reference state at the current temperature. Real reason for its existence
481 * is that it also checks to see if a recalculation of the reference
482 * thermodynamics functions needs to be done.
483 */
484 const vector<double>& cp_R_ref() const {
486 return m_cp0_R;
487 }
488
489 //! @}
490 //! @name Utility Functions
491 //! @{
492
493 bool addSpecies(shared_ptr<Species> spec) override;
494 void initThermo() override;
495 void getParameters(AnyMap& phaseNode) const override;
496 void getSpeciesParameters(const string& name, AnyMap& speciesNode) const override;
497 void setToEquilState(const double* mu_RT) override;
498
499 //! Set the form for the standard and generalized concentrations
500 /*!
501 * Must be one of 'unity', 'species-molar-volume', or
502 * 'solvent-molar-volume'. The default is 'unity'.
503 *
504 * | m_formGC | GeneralizedConc | StandardConc |
505 * | -------------------- | --------------- | ------------ |
506 * | unity | X_k | 1.0 |
507 * | species-molar-volume | X_k / V_k | 1.0 / V_k |
508 * | solvent-molar-volume | X_k / V_N | 1.0 / V_N |
509 *
510 * The value and form of the generalized concentration will affect
511 * reaction rate constants involving species in this phase.
512 */
513 void setStandardConcentrationModel(const string& model);
514
515 /**
516 * Report the molar volume of species k
517 *
518 * units - @f$ m^3 kmol^-1 @f$
519 *
520 * @param k species index
521 */
522 double speciesMolarVolume(int k) const;
523
524 /**
525 * Fill in a return vector containing the species molar volumes.
526 *
527 * units - @f$ m^3 kmol^-1 @f$
528 *
529 * @param smv output vector containing species molar volumes.
530 * Length: m_kk.
531 */
532 void getSpeciesMolarVolumes(double* smv) const;
533
534 //! @}
535
536protected:
537 void compositionChanged() override;
538
539 /**
540 * The standard concentrations can have one of three different forms:
541 * 0 = 'unity', 1 = 'species-molar-volume', 2 = 'solvent-molar-volume'. See
542 * setStandardConcentrationModel().
543 */
544 int m_formGC = 0;
545
546 /**
547 * Value of the reference pressure for all species in this phase. The T
548 * dependent polynomials are evaluated at the reference pressure. Note,
549 * because this is a single value, all species are required to have the same
550 * reference pressure.
551 */
552 double m_Pref = OneAtm;
553
554 /**
555 * m_Pcurrent = The current pressure
556 * Since the density isn't a function of pressure, but only of the
557 * mole fractions, we need to independently specify the pressure.
558 * The density variable which is inherited as part of the State class,
559 * m_dens, is always kept current whenever T, P, or X[] change.
560 */
562
563 //! Vector of molar volumes for each species in the solution
564 /**
565 * Species molar volumes (@f$ m^3 kmol^-1 @f$) at the current mixture state.
566 * For the IdealSolidSolnPhase class, these are constant.
567 */
568 mutable vector<double> m_speciesMolarVolume;
569
570 //! Vector containing the species reference enthalpies at T = m_tlast
571 mutable vector<double> m_h0_RT;
572
573 //! Vector containing the species reference constant pressure heat
574 //! capacities at T = m_tlast
575 mutable vector<double> m_cp0_R;
576
577 //! Vector containing the species reference Gibbs functions at T = m_tlast
578 mutable vector<double> m_g0_RT;
579
580 //! Vector containing the species reference entropies at T = m_tlast
581 mutable vector<double> m_s0_R;
582
583 //! Vector containing the species reference exp(-G/RT) functions at
584 //! T = m_tlast
585 mutable vector<double> m_expg0_RT;
586
587 //! Temporary array used in equilibrium calculations
588 mutable vector<double> m_pp;
589
590private:
591 //! @name Utility Functions
592 //! @{
593 /**
594 * This function gets called for every call to functions in this class. It
595 * checks to see whether the temperature has changed and thus the reference
596 * thermodynamics functions for all of the species must be recalculated. If
597 * the temperature has changed, the species thermo manager is called to
598 * recalculate G, Cp, H, and S at the current temperature.
599 */
600 virtual void _updateThermo() const;
601
602 //! @}
603};
604}
605
606#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
Class IdealSolidSolnPhase represents a condensed phase ideal solution compound.
const vector< double > & entropy_R_ref() const
Returns a reference to the vector of nondimensional enthalpies of the reference state at the current ...
void getStandardChemPotentials(double *mu0) const override
Get the standard state chemical potentials of the species.
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.
double pressure() const override
Pressure.
vector< double > m_g0_RT
Vector containing the species reference Gibbs functions at T = m_tlast.
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...
void getEntropy_R(double *sr) const override
Get the nondimensional Entropies for the species standard states at the current T and P of the soluti...
vector< double > m_h0_RT
Vector containing the species reference enthalpies at T = m_tlast.
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 speciesMolarVolume(int k) const
Report the molar volume of species k.
vector< double > m_pp
Temporary array used in equilibrium calculations.
double m_Pref
Value of the reference pressure for all species in this phase.
void getCp_R(double *cpr) const override
Get the nondimensional heat capacity at constant pressure function for the species standard states at...
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 the array of generalized concentrations.
void getSpeciesMolarVolumes(double *smv) const
Fill in a return vector containing the species molar volumes.
void setPressure(double p) override
Set the pressure at constant temperature.
const vector< double > & gibbs_RT_ref() const
Returns a reference to the vector of nondimensional enthalpies of the reference state at the current ...
void getPartialMolarVolumes(double *vbar) const override
returns an array of partial molar volumes of the species in the solution.
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 pure species at the current T and P of the solution.
double standardConcentration(size_t k) const override
The standard concentration used to normalize the generalized concentration.
void setStandardConcentrationModel(const string &model)
Set the form for the standard and generalized concentrations.
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 array of nondimensional Enthalpy functions for the standard state species at the current T an...
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
Vector containing the species reference entropies at T = m_tlast.
void getGibbs_RT(double *grt) const override
Get the nondimensional Gibbs function for the species standard states at the current T and P of the s...
double entropy_mole() const override
Molar entropy of the solution.
int m_formGC
The standard concentrations can have one of three different forms: 0 = 'unity', 1 = 'species-molar-vo...
vector< double > m_speciesMolarVolume
Vector of molar volumes for each species in the solution.
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.
double cp_mole() const override
Molar heat capacity at constant pressure of the solution.
void getIntEnergy_RT(double *urt) const override
Returns the vector of nondimensional Internal Energies of the standard state species at the current T...
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 free energy of the solution.
vector< double > m_cp0_R
Vector containing the species reference constant pressure heat capacities at T = m_tlast.
bool addSpecies(shared_ptr< Species > spec) override
Add a Species to this Phase.
virtual void _updateThermo() const
This function gets called for every call to functions in this class.
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 species activity coefficients.
vector< double > m_expg0_RT
Vector containing the species reference exp(-G/RT) functions at T = m_tlast.
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 vector of nondimensional enthalpies of the reference state at the current ...
double m_Pcurrent
m_Pcurrent = The current pressure Since the density isn't a function of pressure, but only of the mol...
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 calcDensity()
Calculate the density of the mixture using the partial molar volumes and mole fractions as input.
const vector< double > & enthalpy_RT_ref() const
Returns a reference to the vector of nondimensional enthalpies of the reference state at the current ...
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