Cantera  3.1.0a1
Loading...
Searching...
No Matches
ThermoPhase.h
Go to the documentation of this file.
1/**
2 * @file ThermoPhase.h
3 * Header file for class ThermoPhase, the base class for phases with
4 * thermodynamic properties, and the text for the Module thermoprops
5 * (see @ref thermoprops and class @link Cantera::ThermoPhase ThermoPhase@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_THERMOPHASE_H
12#define CT_THERMOPHASE_H
13
14#include "Phase.h"
15#include "MultiSpeciesThermo.h"
16#include "cantera/base/Units.h"
17#include "cantera/base/AnyMap.h"
18
19namespace Cantera
20{
21
22/**
23 * @defgroup thermoprops Thermodynamic Properties
24 *
25 * These classes are used to compute the thermodynamic properties of phases of matter.
26 * The main base class for describing thermodynamic properties of phases within %Cantera
27 * is called ThermoPhase. %ThermoPhase is a large class that describes the interface
28 * within %Cantera to thermodynamic functions for a phase.
29 *
30 * ## Categorizing the Different ThermoPhase Objects
31 *
32 * ThermoPhase objects may be cataloged into four general bins.
33 *
34 * The first type are those whose underlying species have a reference state associated
35 * with them. The reference state describes the thermodynamic functions for a species at
36 * a single reference pressure, @f$ p_0 @f$. The thermodynamic functions are specified
37 * via derived objects of the SpeciesThermoInterpType object class, and usually consist
38 * of polynomials in temperature such as the NASA polynomial or the Shomate polynomial.
39 * Calculators for these reference states, which manage the calculation for all of the
40 * species in a phase, are all derived from the virtual base class
41 * SpeciesThermoInterpType. Calculators are needed because the actual calculation of the
42 * reference state thermodynamics has been shown to be relatively expensive. A great
43 * deal of work has gone into devising efficient schemes for calculating the
44 * thermodynamic polynomials of a set of species in a phase, in particular gas species
45 * in ideal gas phases whose reference state thermodynamics is specified by NASA
46 * polynomials.
47 *
48 * The reference state thermodynamics combined with the mixing rules and an assumption
49 * about the pressure dependence yields the thermodynamic functions for the phase.
50 * Expressions involving the specification of the fugacities of species would fall into
51 * this category of %ThermoPhase objects. Note, however, that at this time, we do not
52 * have any nontrivial examples of these types of phases. In general, the independent
53 * variables that completely describe the state of the system for this class are
54 * temperature, the phase density, and @f$ N - 1 @f$ species mole or mass fractions.
55 * Additionally, if the phase involves charged species, the phase electric potential is
56 * an added independent variable. Examples of this first class of %ThermoPhase models,
57 * which includes the IdealGasPhase object, the most commonly used object with %Cantera,
58 * include:
59 *
60 * - IdealGasPhase
61 * - StoichSubstance
62 * - SurfPhase
63 * - EdgePhase
64 * - LatticePhase
65 * - LatticeSolidPhase
66 * - PureFluidPhase
67 * - IdealSolidSolnPhase
68 * - VPStandardStateTP
69 *
70 * The second class of objects are all derivatives of the VPStandardStateTP class listed
71 * above. These classes assume that there exists a standard state for each species in
72 * the phase, where the thermodynamic functions are specified as a function of
73 * temperature and pressure. Standard state objects for each species are all derived
74 * from the PDSS virtual base class. In turn, these standard states may employ reference
75 * state calculation to aid in their calculations. However, there are some PDSS objects
76 * which do not employ reference state calculations. An example of this is real equation
77 * of state for liquid water used within the calculation of brine thermodynamics. In
78 * general, the independent variables that completely describe the state of the system
79 * for this class are temperature, the phase pressure, and @f$ N - 1 @f$ species mole or
80 * mass fractions or molalities. The standard state thermodynamics combined with the
81 * mixing rules yields the thermodynamic functions for the phase. Mixing rules are given
82 * in terms of specifying the molar-base activity coefficients or activities. Lists of
83 * phases which belong to this group are given below
84 *
85 * - IdealSolnGasVPSS
86 * - MolalityVPSSTP
87 *
88 * Note, the ideal gas and ideal solution approximations are lumped together in the
89 * class IdealSolnGasVPSS, because at this level they look alike having the same mixing
90 * rules with respect to the specification of the excess thermodynamic properties.
91 *
92 * The third class of objects are all derivatives of the MolalityVPSSTP object. They
93 * assume that the standard states are temperature and pressure dependent but they also
94 * assume that the standard states are molality-based. In other words, they assume that
95 * the standard state of the solute species are in a pseudo state of 1 molality but at
96 * infinite dilution. A solvent must be specified in these calculations, defined as the
97 * first species in the phase, and its standard state is the pure solvent state. Phases
98 * which belong to this group include:
99 *
100 * - DebyeHuckel
101 * - IdealMolalSoln
102 * - HMWSoln
103 *
104 * The fourth class of %ThermoPhase objects are stoichiometric phases. Stoichiometric
105 * phases are phases which consist of one and only one species. The class
106 * SingleSpeciesTP is the base class for these substances. Within the class, the general
107 * %ThermoPhase interface is dumbed down so that phases consisting of one species may be
108 * succinctly described. These phases may have PDSS classes or SpeciesThermoInterpType
109 * calculators associated with them. In general, the independent variables that
110 * completely describe the state of the system for this class are temperature and either
111 * the phase density or the phase pressure. Classes in this group include:
112 *
113 * - StoichSubstance
114 * - WaterSSTP
115 *
116 * ## Creating ThermoPhase objects
117 *
118 * Instances of subclasses of ThermoPhase should be created using the factory methods
119 * newThermo(const string&, const string&), newThermo(const AnyMap&, const AnyMap&), or
120 * newThermoModel(). This allows new classes to be used with the various %Cantera
121 * language interfaces.
122 *
123 * ## Defining new thermodynamic models
124 *
125 * To implement a new equation of state, derive a class from ThermoPhase or a relevant
126 * existing derived class and overload the virtual methods in ThermoPhase. Methods that
127 * are not needed can be left unimplemented, which will cause an exception to be thrown
128 * if they are called.
129 */
130
131//! @name CONSTANTS - Specification of the Molality convention
132//! @{
133
134//! Standard state uses the molar convention
136//! Standard state uses the molality convention
138
139//! @}
140//! @name CONSTANTS - Specification of the SS convention
141//! @{
142
143//! Standard state uses the molar convention
145//! Standard state uses the molality convention
147//! Standard state thermodynamics is obtained from slave ThermoPhase objects
149//! @}
150
151//! Differentiate between mole fractions and mass fractions for input mixture
152//! composition
153enum class ThermoBasis
154{
155 mass,
156 molar
157};
158
159//! Base class for a phase with thermodynamic properties.
160/*!
161 * Class ThermoPhase is the base class for the family of classes that represent
162 * phases of matter of any type. It defines a common public interface, and
163 * implements a few methods. Most of the methods, however, are declared virtual
164 * and are meant to be overloaded in derived classes. The standard way used
165 * throughout %Cantera to compute properties of phases of matter is through
166 * pointers of type `ThermoPhase*` that point to objects of subclasses of
167 * ThermoPhase.
168 *
169 * Class ThermoPhase extends class Phase by adding methods to compute
170 * thermodynamic properties in addition to the ones that are used to define the
171 * state of a substance (temperature, density/pressure and composition). The
172 * distinction is that the methods declared in ThermoPhase require knowing the
173 * particular equation of state of the phase of interest, while those of class
174 * Phase do not, since they only involve data values stored within the object.
175 * These methods are then implemented by the classes derived from ThermoPhase to
176 * represent a phase with a specific equation of state.
177 *
178 * ## Calculating and accessing thermodynamic properties
179 *
180 * The calculation of thermodynamic functions within %ThermoPhase is broken down roughly
181 * into two or more steps. First, the standard state properties of all of the species
182 * are calculated at the current temperature and at either the current pressure or at a
183 * reference pressure. If the calculation is carried out at a reference pressure instead
184 * of at the current pressure the calculation is called a "reference state properties"
185 * calculation, just to make the distinction (even though it may be considered to be a
186 * fixed-pressure standard-state calculation). The next step is to adjust the reference
187 * state calculation to the current pressure. The thermodynamic functions then are
188 * considered to be at the standard state of each species. Lastly the mixing
189 * contributions are added to arrive at the thermodynamic functions for the solution.
190 *
191 * The %ThermoPhase class provides interfaces to thermodynamic properties calculated for
192 * the reference state of each species, the standard state values for each species, the
193 * thermodynamic functions for solution values, both on a per mole of solution basis
194 * (such as ThermoPhase::enthalpy_mole()), on a per kg of solution basis, and on a
195 * partial molar basis for each species (such as
196 * ThermoPhase::getPartialMolarEnthalpies). At each level, functions for the enthalpy,
197 * entropy, Gibbs free energy, internal energy, and volume are provided. So, 5 levels
198 * (reference state, standard state, partial molar, per mole of solution, and per mass
199 * of solution) and 5 functions multiplied together makes 25 possible functions. That's
200 * why %ThermoPhase is such a large class.
201 *
202 * ## Setting the State of the phase
203 *
204 * Typically, the way the ThermoPhase object works is that there are a set of functions
205 * that set the state of the phase via setting the internal independent variables. Then,
206 * there are another set of functions that query the thermodynamic functions evaluated
207 * at the current %State of the phase. Internally, most of the intermediate work
208 * generally occurs at the point where the internal state of the system is set and not
209 * at the time when individual thermodynamic functions are queried (though the actual
210 * breakdown in work is dependent on the individual derived ThermoPhase object).
211 * Therefore, for efficiency, the user should lump together queries of thermodynamic
212 * functions after setting the state. Moreover, in setting the state, if the density is
213 * the independent variable, the following order should be used:
214 *
215 * - Set the temperature
216 * - Set the mole or mass fractions or set the molalities
217 * - set the pressure.
218 *
219 * For classes which inherit from VPStandardStateTP, the above order may be used, or the
220 * following order may be used. It's not important.
221 *
222 * - Set the temperature
223 * - Set the pressure
224 * - Set the mole or mass fractions or set the molalities
225 *
226 * See the @ref sec-thermophase-set-state "list of methods" that can be used to set
227 * the complete state of ThermoPhase objects.
228 *
229 * ## Treatment of the phase potential and the electrochemical potential of a species
230 *
231 * The electrochemical potential of species k in a phase p, @f$ \zeta_k @f$, is related
232 * to the chemical potential as:
233 *
234 * @f[
235 * \zeta_{k}(T,P) = \mu_{k}(T,P) + z_k \phi_p
236 * @f]
237 *
238 * where @f$ \nu_k @f$ is the charge of species k, and @f$ \phi_p @f$ is the electric
239 * potential of phase p.
240 *
241 * The potential @f$ \phi_p @f$ is tracked and internally stored within the base
242 * ThermoPhase object. It constitutes a specification of the internal state of the
243 * phase; it's the third state variable, the first two being temperature and density
244 * (or, pressure, for incompressible equations of state). It may be set with the
245 * function, setElectricPotential(), and may be queried with the function
246 * electricPotential().
247 *
248 * Note, the overall electrochemical potential of a phase may not be changed by the
249 * potential because many phases enforce charge neutrality:
250 *
251 * @f[
252 * 0 = \sum_k z_k X_k
253 * @f]
254 *
255 * Whether charge neutrality is necessary for a phase is also specified within the
256 * ThermoPhase object, by the function call chargeNeutralityNecessary(). Note, that it
257 * is not necessary for the ideal gas phase, currently. However, it is necessary for
258 * liquid phases such as DebyeHuckel and HMWSoln for the proper specification of the
259 * chemical potentials.
260 *
261 * This equation, when applied to the @f$ \zeta_k @f$ equation described above, results
262 * in a zero net change in the effective Gibbs free energy of the phase. However,
263 * specific charged species in the phase may increase or decrease their electrochemical
264 * potentials, which will have an effect on interfacial reactions involving charged
265 * species, when there is a potential drop between phases. This effect is used within
266 * the InterfaceKinetics and EdgeKinetics classes.
267 *
268 * ## Specification of Activities and Activity Conventions
269 *
270 * The activity @f$ a_k @f$ and activity coefficient @f$ \gamma_k @f$ of a species in
271 * solution is related to the chemical potential by
272 *
273 * @f[
274 * \mu_k = \mu_k^0(T,P) + \hat R T \ln a_k = \mu_k^0(T,P) + \hat R T \ln x_k \gamma_k
275 * @f]
276 *
277 * The quantity @f$ \mu_k^0(T,P) @f$ is the standard chemical potential at unit
278 * activity, which depends on the temperature and pressure, but not on the composition.
279 * The activity is dimensionless. Within liquid electrolytes it's common to use a
280 * molality convention, where solute species employ the molality-based activity
281 * coefficients:
282 *
283 * @f[
284 * \mu_k = \mu_k^\triangle(T,P) + R T \ln a_k^{\triangle} =
285 * \mu_k^\triangle(T,P) + R T \ln \frac{\gamma_k^{\triangle} m_k}{m^\triangle}
286 * @f]
287 *
288 * And the solvent employs the convention
289 * @f[
290 * \mu_o = \mu^o_o(T,P) + RT \ln a_o
291 * @f]
292 *
293 * where @f$ a_o @f$ is often redefined in terms of the osmotic coefficient @f$ \phi
294 * @f$:
295 *
296 * @f[
297 * \phi = \frac{- \ln a_o}{\tilde{M}_o \sum_{i \ne o} m_i}
298 * @f]
299 *
300 * ThermoPhase classes which employ the molality based convention are all derived from
301 * the MolalityVPSSTP class. See the class description for further information on its
302 * capabilities.
303 *
304 * The activity convention used by a ThermoPhase object may be queried via the
305 * activityConvention() function. A zero means molar based, while a one
306 * means molality based.
307 *
308 * The function getActivities() returns a vector of activities. Whether these are
309 * molar-based or molality-based depends on the value of activityConvention().
310 *
311 * The function getActivityCoefficients() always returns molar-based activity
312 * coefficients regardless of the activity convention used. The function
313 * MolalityVPSSTP::getMolalityActivityCoefficients() returns molality
314 * based activity coefficients for those ThermoPhase objects derived
315 * from the MolalityVPSSTP class. The function MolalityVPSSTP::osmoticCoefficient()
316 * returns the osmotic coefficient.
317
318 * ## Activity Concentrations: Relationship of ThermoPhase to Kinetics Expressions
319 *
320 * %Cantera can handle both thermodynamics and kinetics mechanisms. Reversible kinetics
321 * mechanisms within %Cantera must be compatible with thermodynamics in the sense that
322 * at equilibrium, or at infinite times, the concentrations of species must conform to
323 * thermodynamics. This means that for every valid reversible kinetics reaction in a
324 * mechanism, it must be reducible to an expression involving the ratio of the product
325 * activity to the reactant activities being equal to the exponential of the
326 * dimensionless standard state gibbs free energies of reaction. Irreversible kinetics
327 * reactions do not have this requirement; however, their usage can yield unexpected and
328 * inconsistent results in many situations.
329 *
330 * The actual units used in a kinetics expression depend on the context or the relative
331 * field of study. For example, in gas phase kinetics, species in kinetics expressions
332 * are expressed in terms of concentrations, for example, gmol cm-3. In solid phase
333 * studies, however, kinetics is usually expressed in terms of unitless activities,
334 * which most often equate to solid phase mole fractions. In order to accommodate
335 * variability here, %Cantera has come up with the idea of activity concentrations,
336 * @f$ C^a_k @f$. Activity concentrations are the expressions used directly in kinetics
337 * expressions. These activity (or generalized) concentrations are used by kinetics
338 * manager classes to compute the forward and reverse rates of elementary reactions.
339 * Note that they may or may not have units of concentration --- they might be partial
340 * pressures, mole fractions, or surface coverages, The activity concentrations for
341 * species *k*, @f$ C^a_k @f$, are related to the activity for species *k*, @f$ a_k @f$,
342 * via the expression:
343 *
344 * @f[
345 * a_k = C^a_k / C^0_k
346 * @f]
347 *
348 * @f$ C^0_k @f$ are called standard concentrations. They serve as multiplicative
349 * factors between the activities and the generalized concentrations. Standard
350 * concentrations may be different for each species. They may depend on both the
351 * temperature and the pressure. However, they may not depend on the composition of the
352 * phase. For example, for the IdealGasPhase object the standard concentration is
353 * defined as
354 *
355 * @f[
356 * C^0_k = \frac{P}{RT}
357 * @f]
358 *
359 * while in many solid phase kinetics problems,
360 *
361 * @f[
362 * C^0_k = 1.0
363 * @f]
364 *
365 * is employed making the units for activity concentrations in solids unitless.
366 *
367 * ThermoPhase member functions dealing with this concept include
368 * getActivityConcentrations(), which provides a vector of the current activity
369 * concentrations. The function standardConcentration() returns the standard
370 * concentration of the kth species. The function logStandardConc(), returns the natural
371 * log of the kth standard concentration. The function standardConcentrationUnits()
372 * returns the units of the standard concentration.
373 *
374 * ### Equilibrium constants
375 *
376 * - @f$ K_a @f$ is the equilibrium constant defined in terms of the standard state
377 * Gibbs free energy values. It is by definition dimensionless.
378 *
379 * - @f$ K_p @f$ is the equilibrium constant defined in terms of the reference state
380 * Gibbs free energy values. It is by definition dimensionless. The pressure
381 * dependence is handled entirely on the RHS of the equilibrium expression.
382 *
383 * - @f$ K_c @f$ is the equilibrium constant defined in terms of the activity
384 * concentrations. The dimensions depend on the number of products and reactants.
385 *
386 * The kinetics manager requires the calculation of @f$ K_c @f$ for the calculation of
387 * the reverse rate constant.
388 *
389 * @ingroup thermoprops
390 */
391class ThermoPhase : public Phase
392{
393public:
394 //! Constructor. Note that ThermoPhase is meant to be used as a base class,
395 //! so this constructor should not be called explicitly.
396 ThermoPhase() = default;
397
398 //! @name Information Methods
399 //! @{
400
401 string type() const override {
402 return "none";
403 }
404
405 //! Boolean indicating whether phase is ideal
406 virtual bool isIdeal() const {
407 return false;
408 }
409
410 //! String indicating the mechanical phase of the matter in this Phase.
411 /*!
412 * Options for the string are:
413 * * `unspecified`
414 * * `supercritical`
415 * * `gas`
416 * * `liquid`
417 * * `solid`
418 * * `solid-liquid-mix`
419 * * `solid-gas-mix`
420 * * `liquid-gas-mix`
421 * * `solid-liquid-gas-mix`
422 *
423 * `unspecified` is the default and should be used when the Phase does not
424 * distinguish between mechanical phases or does not have enough information to
425 * determine which mechanical phase(s) are present.
426 *
427 * @todo Needs to be implemented for all phase types. Currently only implemented for
428 * PureFluidPhase.
429 */
430 virtual string phaseOfMatter() const {
431 return "unspecified";
432 }
433
434 /**
435 * Returns the reference pressure in Pa. This function is a wrapper
436 * that calls the species thermo refPressure function.
437 */
438 virtual double refPressure() const {
439 return m_spthermo.refPressure();
440 }
441
442 //! Minimum temperature for which the thermodynamic data for the species
443 //! or phase are valid.
444 /*!
445 * If no argument is supplied, the value returned will be the lowest
446 * temperature at which the data for @e all species are valid. Otherwise,
447 * the value will be only for species @e k. This function is a wrapper that
448 * calls the species thermo minTemp function.
449 *
450 * @param k index of the species. Default is -1, which will return the max
451 * of the min value over all species.
452 */
453 virtual double minTemp(size_t k = npos) const {
454 return m_spthermo.minTemp(k);
455 }
456
457 //! Report the 298 K Heat of Formation of the standard state of one species
458 //! (J kmol-1)
459 /*!
460 * The 298K Heat of Formation is defined as the enthalpy change to create
461 * the standard state of the species from its constituent elements in their
462 * standard states at 298 K and 1 bar.
463 *
464 * @param k species index
465 * @returns the current value of the Heat of Formation at 298K
466 * and 1 bar
467 */
468 double Hf298SS(const size_t k) const {
469 return m_spthermo.reportOneHf298(k);
470 }
471
472 //! Modify the value of the 298 K Heat of Formation of one species in the
473 //! phase (J kmol-1)
474 /*!
475 * The 298K heat of formation is defined as the enthalpy change to create
476 * the standard state of the species from its constituent elements in their
477 * standard states at 298 K and 1 bar.
478 *
479 * @param k Species k
480 * @param Hf298New Specify the new value of the Heat of Formation at
481 * 298K and 1 bar
482 */
483 virtual void modifyOneHf298SS(const size_t k, const double Hf298New) {
484 m_spthermo.modifyOneHf298(k, Hf298New);
486 }
487
488 //! Restore the original heat of formation of one or more species
489 /*!
490 * Resets changes made by modifyOneHf298SS(). If the species index is not
491 * specified, the heats of formation for all species are restored.
492 */
493 virtual void resetHf298(const size_t k=npos);
494
495 //! Maximum temperature for which the thermodynamic data for the species
496 //! are valid.
497 /*!
498 * If no argument is supplied, the value returned will be the highest
499 * temperature at which the data for @e all species are valid. Otherwise,
500 * the value will be only for species @e k. This function is a wrapper that
501 * calls the species thermo maxTemp function.
502 *
503 * @param k index of the species. Default is -1, which will return the min
504 * of the max value over all species.
505 */
506 virtual double maxTemp(size_t k = npos) const {
507 return m_spthermo.maxTemp(k);
508 }
509
510 //! Returns the chargeNeutralityNecessity boolean
511 /*!
512 * Some phases must have zero net charge in order for their thermodynamics
513 * functions to be valid. If this is so, then the value returned from this
514 * function is true. If this is not the case, then this is false. Now, ideal
515 * gases have this parameter set to false, while solution with molality-
516 * based activity coefficients have this parameter set to true.
517 */
520 }
521
522 //! @}
523 //! @name Molar Thermodynamic Properties of the Solution
524 //! @{
525
526 //! Molar enthalpy. Units: J/kmol.
527 virtual double enthalpy_mole() const {
528 throw NotImplementedError("ThermoPhase::enthalpy_mole");
529 }
530
531 //! Molar internal energy. Units: J/kmol.
532 virtual double intEnergy_mole() const {
533 return enthalpy_mole() - pressure()* molarVolume();
534 }
535
536 //! Molar entropy. Units: J/kmol/K.
537 virtual double entropy_mole() const {
538 throw NotImplementedError("ThermoPhase::entropy_mole");
539 }
540
541 //! Molar Gibbs function. Units: J/kmol.
542 virtual double gibbs_mole() const {
544 }
545
546 //! Molar heat capacity at constant pressure. Units: J/kmol/K.
547 virtual double cp_mole() const {
548 throw NotImplementedError("ThermoPhase::cp_mole");
549 }
550
551 //! Molar heat capacity at constant volume. Units: J/kmol/K.
552 virtual double cv_mole() const {
553 throw NotImplementedError("ThermoPhase::cv_mole");
554 }
555
556 //! @}
557 //! @name Mechanical Properties
558 //! @{
559
560 //! Returns the isothermal compressibility. Units: 1/Pa.
561 /*!
562 * The isothermal compressibility is defined as
563 * @f[
564 * \kappa_T = -\frac{1}{v}\left(\frac{\partial v}{\partial P}\right)_T
565 * @f]
566 * or
567 * @f[
568 * \kappa_T = \frac{1}{\rho}\left(\frac{\partial \rho}{\partial P}\right)_T
569 * @f]
570 */
571 virtual double isothermalCompressibility() const {
572 throw NotImplementedError("ThermoPhase::isothermalCompressibility");
573 }
574
575 //! Return the volumetric thermal expansion coefficient. Units: 1/K.
576 /*!
577 * The thermal expansion coefficient is defined as
578 * @f[
579 * \beta = \frac{1}{v}\left(\frac{\partial v}{\partial T}\right)_P
580 * @f]
581 */
582 virtual double thermalExpansionCoeff() const {
583 throw NotImplementedError("ThermoPhase::thermalExpansionCoeff");
584 }
585
586 //! Return the speed of sound. Units: m/s.
587 /*!
588 * The speed of sound is defined as
589 * @f[
590 * c = \sqrt{\left(\frac{\partial P}{\partial\rho}\right)_s}
591 * @f]
592 */
593 virtual double soundSpeed() const {
594 throw NotImplementedError("ThermoPhase::soundSpeed");
595 }
596
597 //! @}
598 //! @name Electric Potential
599 //!
600 //! The phase may be at some non-zero electrical potential. These methods
601 //! set or get the value of the electric potential.
602 //! @{
603
604 //! Set the electric potential of this phase (V).
605 /*!
606 * This is used by classes InterfaceKinetics and EdgeKinetics to
607 * compute the rates of charge-transfer reactions, and in computing
608 * the electrochemical potentials of the species.
609 *
610 * Each phase may have its own electric potential.
611 *
612 * @param v Input value of the electric potential in Volts
613 */
614 void setElectricPotential(double v) {
615 m_phi = v;
617 }
618
619 //! Returns the electric potential of this phase (V).
620 /*!
621 * Units are Volts (which are Joules/coulomb)
622 */
623 double electricPotential() const {
624 return m_phi;
625 }
626
627 //! @}
628 //! @name Activities, Standard States, and Activity Concentrations
629 //!
630 //! The activity @f$ a_k @f$ of a species in solution is related to the
631 //! chemical potential by @f[ \mu_k = \mu_k^0(T,P) + \hat R T \ln a_k. @f]
632 //! The quantity @f$ \mu_k^0(T,P) @f$ is the standard chemical potential at
633 //! unit activity, which depends on temperature and pressure, but not on
634 //! composition. The activity is dimensionless.
635 //! @{
636
637 //! This method returns the convention used in specification of the
638 //! activities, of which there are currently two, molar- and molality-based
639 //! conventions.
640 /*!
641 * Currently, there are two activity conventions:
642 * - Molar-based activities
643 * %Unit activity of species at either a hypothetical pure
644 * solution of the species or at a hypothetical
645 * pure ideal solution at infinite dilution
646 * cAC_CONVENTION_MOLAR 0
647 * - default
648 *
649 * - Molality-based activities
650 * (unit activity of solutes at a hypothetical 1 molal
651 * solution referenced to infinite dilution at all
652 * pressures and temperatures).
653 * cAC_CONVENTION_MOLALITY 1
654 */
655 virtual int activityConvention() const;
656
657 //! This method returns the convention used in specification of the standard
658 //! state, of which there are currently two, temperature based, and variable
659 //! pressure based.
660 /*!
661 * Currently, there are two standard state conventions:
662 * - Temperature-based activities
663 * cSS_CONVENTION_TEMPERATURE 0
664 * - default
665 *
666 * - Variable Pressure and Temperature -based activities
667 * cSS_CONVENTION_VPSS 1
668 *
669 * - Thermodynamics is set via slave ThermoPhase objects with
670 * nothing being carried out at this ThermoPhase object level
671 * cSS_CONVENTION_SLAVE 2
672 */
673 virtual int standardStateConvention() const;
674
675 //! Returns the units of the "standard concentration" for this phase
676 /*!
677 * These are the units of the values returned by the functions
678 * getActivityConcentrations() and standardConcentration(), which can
679 * vary between different ThermoPhase-derived classes, or change within
680 * a single class depending on input options. See the documentation for
681 * standardConcentration() for the derived class for specific details.
682 */
683 virtual Units standardConcentrationUnits() const;
684
685 //! This method returns an array of generalized concentrations
686 /*!
687 * @f$ C^a_k @f$ are defined such that @f$ a_k = C^a_k / C^0_k, @f$ where
688 * @f$ C^0_k @f$ is a standard concentration defined below and @f$ a_k @f$
689 * are activities used in the thermodynamic functions. These activity (or
690 * generalized) concentrations are used by kinetics manager classes to
691 * compute the forward and reverse rates of elementary reactions. Note that
692 * they may or may not have units of concentration --- they might be partial
693 * pressures, mole fractions, or surface coverages, for example.
694 *
695 * @param c Output array of generalized concentrations. The units depend
696 * upon the implementation of the reaction rate expressions within
697 * the phase.
698 */
699 virtual void getActivityConcentrations(double* c) const {
700 throw NotImplementedError("ThermoPhase::getActivityConcentrations");
701 }
702
703 //! Return the standard concentration for the kth species
704 /*!
705 * The standard concentration @f$ C^0_k @f$ used to normalize the activity
706 * (that is, generalized) concentration. In many cases, this quantity will be
707 * the same for all species in a phase - for example, for an ideal gas @f$
708 * C^0_k = P/\hat R T @f$. For this reason, this method returns a single
709 * value, instead of an array. However, for phases in which the standard
710 * concentration is species-specific (such as surface species of different
711 * sizes), this method may be called with an optional parameter indicating
712 * the species.
713 *
714 * @param k Optional parameter indicating the species. The default
715 * is to assume this refers to species 0.
716 * @return
717 * Returns the standard concentration. The units are by definition
718 * dependent on the ThermoPhase and kinetics manager representation.
719 */
720 virtual double standardConcentration(size_t k=0) const {
721 throw NotImplementedError("ThermoPhase::standardConcentration");
722 }
723
724 //! Natural logarithm of the standard concentration of the kth species.
725 /*!
726 * @param k index of the species (defaults to zero)
727 */
728 virtual double logStandardConc(size_t k=0) const;
729
730 //! Get the array of non-dimensional activities at the current solution
731 //! temperature, pressure, and solution concentration.
732 /*!
733 * Note, for molality based formulations, this returns the molality based
734 * activities.
735 *
736 * We resolve this function at this level by calling on the
737 * activityConcentration function. However, derived classes may want to
738 * override this default implementation.
739 *
740 * @param a Output vector of activities. Length: m_kk.
741 */
742 virtual void getActivities(double* a) const;
743
744 //! Get the array of non-dimensional molar-based activity coefficients at
745 //! the current solution temperature, pressure, and solution concentration.
746 /*!
747 * @param ac Output vector of activity coefficients. Length: m_kk.
748 */
749 virtual void getActivityCoefficients(double* ac) const {
750 if (m_kk == 1) {
751 ac[0] = 1.0;
752 } else {
753 throw NotImplementedError("ThermoPhase::getActivityCoefficients");
754 }
755 }
756
757 //! Get the array of non-dimensional molar-based ln activity coefficients at
758 //! the current solution temperature, pressure, and solution concentration.
759 /*!
760 * @param lnac Output vector of ln activity coefficients. Length: m_kk.
761 */
762 virtual void getLnActivityCoefficients(double* lnac) const;
763
764 //! @}
765 //! @name Partial Molar Properties of the Solution
766 //! @{
767
768 //! Get the species chemical potentials. Units: J/kmol.
769 /*!
770 * This function returns a vector of chemical potentials of the species in
771 * solution at the current temperature, pressure and mole fraction of the
772 * solution.
773 *
774 * @param mu Output vector of species chemical
775 * potentials. Length: m_kk. Units: J/kmol
776 */
777 virtual void getChemPotentials(double* mu) const {
778 throw NotImplementedError("ThermoPhase::getChemPotentials");
779 }
780
781 //! Get the species electrochemical potentials.
782 /*!
783 * These are partial molar quantities. This method adds a term @f$ F z_k
784 * \phi_p @f$ to each chemical potential. The electrochemical potential of
785 * species k in a phase p, @f$ \zeta_k @f$, is related to the chemical
786 * potential via the following equation,
787 *
788 * @f[
789 * \zeta_{k}(T,P) = \mu_{k}(T,P) + F z_k \phi_p
790 * @f]
791 *
792 * @param mu Output vector of species electrochemical
793 * potentials. Length: m_kk. Units: J/kmol
794 */
795 void getElectrochemPotentials(double* mu) const;
796
797 //! Returns an array of partial molar enthalpies for the species
798 //! in the mixture. Units (J/kmol)
799 /*!
800 * @param hbar Output vector of species partial molar enthalpies.
801 * Length: m_kk. units are J/kmol.
802 */
803 virtual void getPartialMolarEnthalpies(double* hbar) const {
804 throw NotImplementedError("ThermoPhase::getPartialMolarEnthalpies");
805 }
806
807 //! Returns an array of partial molar entropies of the species in the
808 //! solution. Units: J/kmol/K.
809 /*!
810 * @param sbar Output vector of species partial molar entropies.
811 * Length = m_kk. units are J/kmol/K.
812 */
813 virtual void getPartialMolarEntropies(double* sbar) const {
814 throw NotImplementedError("ThermoPhase::getPartialMolarEntropies");
815 }
816
817 //! Return an array of partial molar internal energies for the
818 //! species in the mixture. Units: J/kmol.
819 /*!
820 * @param ubar Output vector of species partial molar internal energies.
821 * Length = m_kk. units are J/kmol.
822 */
823 virtual void getPartialMolarIntEnergies(double* ubar) const {
824 throw NotImplementedError("ThermoPhase::getPartialMolarIntEnergies");
825 }
826
827 //! Return an array of partial molar heat capacities for the
828 //! species in the mixture. Units: J/kmol/K
829 /*!
830 * @param cpbar Output vector of species partial molar heat
831 * capacities at constant pressure.
832 * Length = m_kk. units are J/kmol/K.
833 */
834 virtual void getPartialMolarCp(double* cpbar) const {
835 throw NotImplementedError("ThermoPhase::getPartialMolarCp");
836 }
837
838 //! Return an array of partial molar volumes for the
839 //! species in the mixture. Units: m^3/kmol.
840 /*!
841 * @param vbar Output vector of species partial molar volumes.
842 * Length = m_kk. units are m^3/kmol.
843 */
844 virtual void getPartialMolarVolumes(double* vbar) const {
845 throw NotImplementedError("ThermoPhase::getPartialMolarVolumes");
846 }
847
848 //! @}
849 //! @name Properties of the Standard State of the Species in the Solution
850 //! @{
851
852 //! Get the array of chemical potentials at unit activity for the species at
853 //! their standard states at the current *T* and *P* of the solution.
854 /*!
855 * These are the standard state chemical potentials @f$ \mu^0_k(T,P)
856 * @f$. The values are evaluated at the current temperature and pressure of
857 * the solution
858 *
859 * @param mu Output vector of chemical potentials.
860 * Length: m_kk.
861 */
862 virtual void getStandardChemPotentials(double* mu) const {
863 throw NotImplementedError("ThermoPhase::getStandardChemPotentials");
864 }
865
866 //! Get the nondimensional Enthalpy functions for the species at their
867 //! standard states at the current *T* and *P* of the solution.
868 /*!
869 * @param hrt Output vector of nondimensional standard state enthalpies.
870 * Length: m_kk.
871 */
872 virtual void getEnthalpy_RT(double* hrt) const {
873 throw NotImplementedError("ThermoPhase::getEnthalpy_RT");
874 }
875
876 //! Get the array of nondimensional Entropy functions for the standard state
877 //! species at the current *T* and *P* of the solution.
878 /*!
879 * @param sr Output vector of nondimensional standard state entropies.
880 * Length: m_kk.
881 */
882 virtual void getEntropy_R(double* sr) const {
883 throw NotImplementedError("ThermoPhase::getEntropy_R");
884 }
885
886 //! Get the nondimensional Gibbs functions for the species in their standard
887 //! states at the current *T* and *P* of the solution.
888 /*!
889 * @param grt Output vector of nondimensional standard state Gibbs free
890 * energies. Length: m_kk.
891 */
892 virtual void getGibbs_RT(double* grt) const {
893 throw NotImplementedError("ThermoPhase::getGibbs_RT");
894 }
895
896 //! Get the Gibbs functions for the standard state of the species at the
897 //! current *T* and *P* of the solution
898 /*!
899 * Units are Joules/kmol
900 * @param gpure Output vector of standard state Gibbs free energies.
901 * Length: m_kk.
902 */
903 virtual void getPureGibbs(double* gpure) const {
904 throw NotImplementedError("ThermoPhase::getPureGibbs");
905 }
906
907 //! Returns the vector of nondimensional Internal Energies of the standard
908 //! state species at the current *T* and *P* of the solution
909 /*!
910 * @param urt output vector of nondimensional standard state internal energies
911 * of the species. Length: m_kk.
912 */
913 virtual void getIntEnergy_RT(double* urt) const {
914 throw NotImplementedError("ThermoPhase::getIntEnergy_RT");
915 }
916
917 //! Get the nondimensional Heat Capacities at constant pressure for the
918 //! species standard states at the current *T* and *P* of the
919 //! solution
920 /*!
921 * @param cpr Output vector of nondimensional standard state heat
922 * capacities. Length: m_kk.
923 */
924 virtual void getCp_R(double* cpr) const {
925 throw NotImplementedError("ThermoPhase::getCp_R");
926 }
927
928 //! Get the molar volumes of the species standard states at the current
929 //! *T* and *P* of the solution.
930 /*!
931 * units = m^3 / kmol
932 *
933 * @param vol Output vector containing the standard state volumes.
934 * Length: m_kk.
935 */
936 virtual void getStandardVolumes(double* vol) const {
937 throw NotImplementedError("ThermoPhase::getStandardVolumes");
938 }
939
940 //! @}
941 //! @name Thermodynamic Values for the Species Reference States
942 //! @{
943
944 //! Returns the vector of nondimensional enthalpies of the reference state
945 //! at the current temperature of the solution and the reference pressure
946 //! for the species.
947 /*!
948 * @param hrt Output vector containing the nondimensional reference
949 * state enthalpies. Length: m_kk.
950 */
951 virtual void getEnthalpy_RT_ref(double* hrt) const {
952 throw NotImplementedError("ThermoPhase::getEnthalpy_RT_ref");
953 }
954
955 //! Returns the vector of nondimensional Gibbs Free Energies of the
956 //! reference state at the current temperature of the solution and the
957 //! reference pressure for the species.
958 /*!
959 * @param grt Output vector containing the nondimensional reference state
960 * Gibbs Free energies. Length: m_kk.
961 */
962 virtual void getGibbs_RT_ref(double* grt) const {
963 throw NotImplementedError("ThermoPhase::getGibbs_RT_ref");
964 }
965
966 //! Returns the vector of the Gibbs function of the reference state at the
967 //! current temperature of the solution and the reference pressure for the
968 //! species.
969 /*!
970 * @param g Output vector containing the reference state
971 * Gibbs Free energies. Length: m_kk. Units: J/kmol.
972 */
973 virtual void getGibbs_ref(double* g) const {
974 throw NotImplementedError("ThermoPhase::getGibbs_ref");
975 }
976
977 //! Returns the vector of nondimensional entropies of the reference state at
978 //! the current temperature of the solution and the reference pressure for
979 //! each species.
980 /*!
981 * @param er Output vector containing the nondimensional reference
982 * state entropies. Length: m_kk.
983 */
984 virtual void getEntropy_R_ref(double* er) const {
985 throw NotImplementedError("ThermoPhase::getEntropy_R_ref");
986 }
987
988 //! Returns the vector of nondimensional internal Energies of the reference
989 //! state at the current temperature of the solution and the reference
990 //! pressure for each species.
991 /*!
992 * @param urt Output vector of nondimensional reference state internal
993 * energies of the species. Length: m_kk
994 */
995 virtual void getIntEnergy_RT_ref(double* urt) const {
996 throw NotImplementedError("ThermoPhase::getIntEnergy_RT_ref");
997 }
998
999 //! Returns the vector of nondimensional constant pressure heat capacities
1000 //! of the reference state at the current temperature of the solution and
1001 //! reference pressure for each species.
1002 /*!
1003 * @param cprt Output vector of nondimensional reference state
1004 * heat capacities at constant pressure for the species.
1005 * Length: m_kk
1006 */
1007 virtual void getCp_R_ref(double* cprt) const {
1008 throw NotImplementedError("ThermoPhase::getCp_R_ref");
1009 }
1010
1011 //! Get the molar volumes of the species reference states at the current
1012 //! *T* and *P_ref* of the solution.
1013 /*!
1014 * units = m^3 / kmol
1015 *
1016 * @param vol Output vector containing the standard state volumes.
1017 * Length: m_kk.
1018 */
1019 virtual void getStandardVolumes_ref(double* vol) const {
1020 throw NotImplementedError("ThermoPhase::getStandardVolumes_ref");
1021 }
1022
1023 // The methods below are not virtual, and should not be overloaded.
1024
1025 //! @}
1026 //! @name Specific Properties
1027 //! @{
1028
1029 //! Specific enthalpy. Units: J/kg.
1030 double enthalpy_mass() const {
1032 }
1033
1034 //! Specific internal energy. Units: J/kg.
1035 double intEnergy_mass() const {
1037 }
1038
1039 //! Specific entropy. Units: J/kg/K.
1040 double entropy_mass() const {
1042 }
1043
1044 //! Specific Gibbs function. Units: J/kg.
1045 double gibbs_mass() const {
1047 }
1048
1049 //! Specific heat at constant pressure. Units: J/kg/K.
1050 double cp_mass() const {
1051 return cp_mole()/meanMolecularWeight();
1052 }
1053
1054 //! Specific heat at constant volume. Units: J/kg/K.
1055 double cv_mass() const {
1056 return cv_mole()/meanMolecularWeight();
1057 }
1058 //! @}
1059
1060 //! Return the Gas Constant multiplied by the current temperature
1061 /*!
1062 * The units are Joules kmol-1
1063 */
1064 double RT() const {
1065 return temperature() * GasConstant;
1066 }
1067
1068 //! @name Setting the State
1069 //! @anchor sec-thermophase-set-state
1070 //!
1071 //! These methods set all or part of the thermodynamic state.
1072 //! @{
1073
1074 //! Set the temperature (K), pressure (Pa), and mole fractions.
1075 /*!
1076 * Note, the mole fractions are set first before the pressure is set.
1077 * Setting the pressure may involve the solution of a nonlinear equation.
1078 *
1079 * @param t Temperature (K)
1080 * @param p Pressure (Pa)
1081 * @param x Vector of mole fractions.
1082 * Length is equal to m_kk.
1083 */
1084 virtual void setState_TPX(double t, double p, const double* x);
1085
1086 //! Set the temperature (K), pressure (Pa), and mole fractions.
1087 /*!
1088 * Note, the mole fractions are set first before the pressure is set.
1089 * Setting the pressure may involve the solution of a nonlinear equation.
1090 *
1091 * @param t Temperature (K)
1092 * @param p Pressure (Pa)
1093 * @param x Composition map of mole fractions. Species not in
1094 * the composition map are assumed to have zero mole fraction
1095 */
1096 virtual void setState_TPX(double t, double p, const Composition& x);
1097
1098 //! Set the temperature (K), pressure (Pa), and mole fractions.
1099 /*!
1100 * Note, the mole fractions are set first before the pressure is set.
1101 * Setting the pressure may involve the solution of a nonlinear equation.
1102 *
1103 * @param t Temperature (K)
1104 * @param p Pressure (Pa)
1105 * @param x String containing a composition map of the mole fractions.
1106 * Species not in the composition map are assumed to have zero
1107 * mole fraction
1108 */
1109 virtual void setState_TPX(double t, double p, const string& x);
1110
1111 //! Set the internally stored temperature (K), pressure (Pa), and mass
1112 //! fractions of the phase.
1113 /*!
1114 * Note, the mass fractions are set first before the pressure is set.
1115 * Setting the pressure may involve the solution of a nonlinear equation.
1116 *
1117 * @param t Temperature (K)
1118 * @param p Pressure (Pa)
1119 * @param y Vector of mass fractions.
1120 * Length is equal to m_kk.
1121 */
1122 virtual void setState_TPY(double t, double p, const double* y);
1123
1124 //! Set the internally stored temperature (K), pressure (Pa), and mass
1125 //! fractions of the phase
1126 /*!
1127 * Note, the mass fractions are set first before the pressure is set.
1128 * Setting the pressure may involve the solution of a nonlinear equation.
1129 *
1130 * @param t Temperature (K)
1131 * @param p Pressure (Pa)
1132 * @param y Composition map of mass fractions. Species not in
1133 * the composition map are assumed to have zero mass fraction
1134 */
1135 virtual void setState_TPY(double t, double p, const Composition& y);
1136
1137 //! Set the internally stored temperature (K), pressure (Pa), and mass
1138 //! fractions of the phase
1139 /*!
1140 * Note, the mass fractions are set first before the pressure is set.
1141 * Setting the pressure may involve the solution of a nonlinear equation.
1142 *
1143 * @param t Temperature (K)
1144 * @param p Pressure (Pa)
1145 * @param y String containing a composition map of the mass fractions.
1146 * Species not in the composition map are assumed to have zero
1147 * mass fraction
1148 */
1149 virtual void setState_TPY(double t, double p, const string& y);
1150
1151 //! Set the temperature (K) and pressure (Pa)
1152 /*!
1153 * Setting the pressure may involve the solution of a nonlinear equation.
1154 *
1155 * @param t Temperature (K)
1156 * @param p Pressure (Pa)
1157 */
1158 virtual void setState_TP(double t, double p);
1159
1160 //! Set the internally stored specific enthalpy (J/kg) and pressure (Pa) of
1161 //! the phase.
1162 /*!
1163 * @param h Specific enthalpy (J/kg)
1164 * @param p Pressure (Pa)
1165 * @param tol Optional parameter setting the tolerance of the calculation.
1166 * Important for some applications where numerical Jacobians
1167 * are being calculated.
1168 */
1169 virtual void setState_HP(double h, double p, double tol=1e-9);
1170
1171 //! Set the specific internal energy (J/kg) and specific volume (m^3/kg).
1172 /*!
1173 * This function fixes the internal state of the phase so that the specific
1174 * internal energy and specific volume have the value of the input
1175 * parameters.
1176 *
1177 * @param u specific internal energy (J/kg)
1178 * @param v specific volume (m^3/kg).
1179 * @param tol Optional parameter setting the tolerance of the calculation.
1180 * Important for some applications where numerical Jacobians
1181 * are being calculated.
1182 */
1183 virtual void setState_UV(double u, double v, double tol=1e-9);
1184
1185 //! Set the specific entropy (J/kg/K) and pressure (Pa).
1186 /*!
1187 * This function fixes the internal state of the phase so that the specific
1188 * entropy and the pressure have the value of the input parameters.
1189 *
1190 * @param s specific entropy (J/kg/K)
1191 * @param p specific pressure (Pa).
1192 * @param tol Optional parameter setting the tolerance of the calculation.
1193 * Important for some applications where numerical Jacobians
1194 * are being calculated.
1195 */
1196 virtual void setState_SP(double s, double p, double tol=1e-9);
1197
1198 //! Set the specific entropy (J/kg/K) and specific volume (m^3/kg).
1199 /*!
1200 * This function fixes the internal state of the phase so that the specific
1201 * entropy and specific volume have the value of the input parameters.
1202 *
1203 * @param s specific entropy (J/kg/K)
1204 * @param v specific volume (m^3/kg).
1205 * @param tol Optional parameter setting the tolerance of the calculation.
1206 * Important for some applications where numerical Jacobians
1207 * are being calculated.
1208 */
1209 virtual void setState_SV(double s, double v, double tol=1e-9);
1210
1211 //! Set the specific entropy (J/kg/K) and temperature (K).
1212 /*!
1213 * This function fixes the internal state of the phase so that the specific
1214 * entropy and temperature have the value of the input parameters.
1215 * This base class function will throw an exception if not overridden.
1216 *
1217 * @param s specific entropy (J/kg/K)
1218 * @param t temperature (K)
1219 * @param tol Optional parameter setting the tolerance of the calculation.
1220 * Important for some applications where numerical Jacobians
1221 * are being calculated.
1222 */
1223 virtual void setState_ST(double s, double t, double tol=1e-9) {
1224 throw NotImplementedError("ThermoPhase::setState_ST");
1225 }
1226
1227 //! Set the temperature (K) and specific volume (m^3/kg).
1228 /*!
1229 * This function fixes the internal state of the phase so that the
1230 * temperature and specific volume have the value of the input parameters.
1231 * This base class function will throw an exception if not overridden.
1232 *
1233 * @param t temperature (K)
1234 * @param v specific volume (m^3/kg)
1235 * @param tol Optional parameter setting the tolerance of the calculation.
1236 * Important for some applications where numerical Jacobians
1237 * are being calculated.
1238 */
1239 virtual void setState_TV(double t, double v, double tol=1e-9) {
1240 throw NotImplementedError("ThermoPhase::setState_TV");
1241 }
1242
1243 //! Set the pressure (Pa) and specific volume (m^3/kg).
1244 /*!
1245 * This function fixes the internal state of the phase so that the
1246 * pressure and specific volume have the value of the input parameters.
1247 * This base class function will throw an exception if not overridden.
1248 *
1249 * @param p pressure (Pa)
1250 * @param v specific volume (m^3/kg)
1251 * @param tol Optional parameter setting the tolerance of the calculation.
1252 * Important for some applications where numerical Jacobians
1253 * are being calculated.
1254 */
1255 virtual void setState_PV(double p, double v, double tol=1e-9) {
1256 throw NotImplementedError("ThermoPhase::setState_PV");
1257 }
1258
1259 //! Set the specific internal energy (J/kg) and pressure (Pa).
1260 /*!
1261 * This function fixes the internal state of the phase so that the specific
1262 * internal energy and pressure have the value of the input parameters.
1263 * This base class function will throw an exception if not overridden.
1264 *
1265 * @param u specific internal energy (J/kg)
1266 * @param p pressure (Pa)
1267 * @param tol Optional parameter setting the tolerance of the calculation.
1268 * Important for some applications where numerical Jacobians
1269 * are being calculated.
1270 */
1271 virtual void setState_UP(double u, double p, double tol=1e-9) {
1272 throw NotImplementedError("ThermoPhase::setState_UP");
1273 }
1274
1275 //! Set the specific volume (m^3/kg) and the specific enthalpy (J/kg)
1276 /*!
1277 * This function fixes the internal state of the phase so that the specific
1278 * volume and the specific enthalpy have the value of the input parameters.
1279 * This base class function will throw an exception if not overridden.
1280 *
1281 * @param v specific volume (m^3/kg)
1282 * @param h specific enthalpy (J/kg)
1283 * @param tol Optional parameter setting the tolerance of the calculation.
1284 * Important for some applications where numerical Jacobians
1285 * are being calculated.
1286 */
1287 virtual void setState_VH(double v, double h, double tol=1e-9) {
1288 throw NotImplementedError("ThermoPhase::setState_VH");
1289 }
1290
1291 //! Set the temperature (K) and the specific enthalpy (J/kg)
1292 /*!
1293 * This function fixes the internal state of the phase so that the
1294 * temperature and specific enthalpy have the value of the input parameters.
1295 * This base class function will throw an exception if not overridden.
1296 *
1297 * @param t temperature (K)
1298 * @param h specific enthalpy (J/kg)
1299 * @param tol Optional parameter setting the tolerance of the calculation.
1300 * Important for some applications where numerical Jacobians
1301 * are being calculated.
1302 */
1303 virtual void setState_TH(double t, double h, double tol=1e-9) {
1304 throw NotImplementedError("ThermoPhase::setState_TH");
1305 }
1306
1307 //! Set the specific entropy (J/kg/K) and the specific enthalpy (J/kg)
1308 /*!
1309 * This function fixes the internal state of the phase so that the
1310 * temperature and pressure have the value of the input parameters.
1311 * This base class function will throw an exception if not overridden.
1312 *
1313 * @param s specific entropy (J/kg/K)
1314 * @param h specific enthalpy (J/kg)
1315 * @param tol Optional parameter setting the tolerance of the calculation.
1316 * Important for some applications where numerical Jacobians
1317 * are being calculated.
1318 */
1319 virtual void setState_SH(double s, double h, double tol=1e-9) {
1320 throw NotImplementedError("ThermoPhase::setState_SH");
1321 }
1322
1323 //! Set the density (kg/m**3) and pressure (Pa) at constant composition
1324 /*!
1325 * This method must be reimplemented in derived classes, where it may
1326 * involve the solution of a nonlinear equation. Within %Cantera, the
1327 * independent variable is the density. Therefore, this function solves for
1328 * the temperature that will yield the desired input pressure and density.
1329 * The composition is held constant during this process.
1330 *
1331 * This base class function will print an error, if not overridden.
1332 *
1333 * @param rho Density (kg/m^3)
1334 * @param p Pressure (Pa)
1335 * @since New in %Cantera 3.0.
1336 */
1337 virtual void setState_DP(double rho, double p) {
1338 throw NotImplementedError("ThermoPhase::setState_DP");
1339 }
1340
1341 //! Set the state using an AnyMap containing any combination of properties
1342 //! supported by the thermodynamic model
1343 /*!
1344 * Accepted keys are:
1345 * * `X` (mole fractions)
1346 * * `Y` (mass fractions)
1347 * * `T` or `temperature`
1348 * * `P` or `pressure` [Pa]
1349 * * `H` or `enthalpy` [J/kg]
1350 * * `U` or `internal-energy` [J/kg]
1351 * * `S` or `entropy` [J/kg/K]
1352 * * `V` or `specific-volume` [m^3/kg]
1353 * * `D` or `density` [kg/m^3]
1354 *
1355 * Composition can be specified as either an AnyMap of species names to
1356 * values or as a composition string. All other values can be given as
1357 * floating point values in Cantera's default units, or as strings with the
1358 * units specified, which will be converted using the Units class.
1359 *
1360 * If no thermodynamic property pair is given, or only one of temperature or
1361 * pressure is given, then 298.15 K and 101325 Pa will be used as necessary
1362 * to fully set the state.
1363 */
1364 virtual void setState(const AnyMap& state);
1365
1366 //! @}
1367 //! @name Set Mixture Composition by Mixture Fraction
1368 //! @{
1369
1370 //! Set the mixture composition according to the
1371 //! mixture fraction = kg fuel / (kg oxidizer + kg fuel)
1372 /*!
1373 * Fuel and oxidizer compositions are given either as
1374 * mole fractions or mass fractions (specified by `basis`)
1375 * and do not need to be normalized. Pressure and temperature are
1376 * kept constant. Elements C, S, H and O are considered for the oxidation.
1377 *
1378 * @param mixFrac mixture fraction (between 0 and 1)
1379 * @param fuelComp composition of the fuel
1380 * @param oxComp composition of the oxidizer
1381 * @param basis either ThermoPhase::molar or ThermoPhase::mass.
1382 * Fuel and oxidizer composition are interpreted
1383 * as mole or mass fractions (default: molar)
1384 */
1385 void setMixtureFraction(double mixFrac, const double* fuelComp,
1386 const double* oxComp, ThermoBasis basis=ThermoBasis::molar);
1387 //! @copydoc ThermoPhase::setMixtureFraction
1388 void setMixtureFraction(double mixFrac, const string& fuelComp,
1389 const string& oxComp, ThermoBasis basis=ThermoBasis::molar);
1390 //! @copydoc ThermoPhase::setMixtureFraction
1391 void setMixtureFraction(double mixFrac, const Composition& fuelComp,
1392 const Composition& oxComp, ThermoBasis basis=ThermoBasis::molar);
1393 //! @}
1394 //! @name Compute Mixture Fraction
1395 //! @{
1396
1397 //! Compute the mixture fraction = kg fuel / (kg oxidizer + kg fuel) for
1398 //! the current mixture given fuel and oxidizer compositions.
1399 /*!
1400 * Fuel and oxidizer compositions are given either as
1401 * mole fractions or mass fractions (specified by `basis`)
1402 * and do not need to be normalized.
1403 * The mixture fraction @f$ Z @f$ can be computed from a single element
1404 * @f[ Z_m = \frac{Z_{\mathrm{mass},m}-Z_{\mathrm{mass},m,\mathrm{ox}}}
1405 * {Z_{\mathrm{mass},\mathrm{fuel}}-Z_{\mathrm{mass},m,\mathrm{ox}}} @f] where
1406 * @f$ Z_{\mathrm{mass},m} @f$ is the elemental mass fraction of element m
1407 * in the mixture, and @f$ Z_{\mathrm{mass},m,\mathrm{ox}} @f$ and
1408 * @f$ Z_{\mathrm{mass},m,\mathrm{fuel}} @f$ are the elemental mass fractions
1409 * of the oxidizer and fuel, or from the Bilger mixture fraction,
1410 * which considers the elements C, S, H and O @cite bilger1979
1411 * @f[ Z_{\mathrm{Bilger}} = \frac{\beta-\beta_{\mathrm{ox}}}
1412 * {\beta_{\mathrm{fuel}}-\beta_{\mathrm{ox}}} @f]
1413 * with @f$ \beta = 2\frac{Z_C}{M_C}+2\frac{Z_S}{M_S}+\frac{1}{2}\frac{Z_H}{M_H}
1414 * -\frac{Z_O}{M_O} @f$
1415 * and @f$ M_m @f$ the atomic weight of element @f$ m @f$.
1416 *
1417 * @param fuelComp composition of the fuel
1418 * @param oxComp composition of the oxidizer
1419 * @param basis either ThermoBasis::molar or ThermoBasis::mass.
1420 * Fuel and oxidizer composition are interpreted
1421 * as mole or mass fractions (default: molar)
1422 * @param element either "Bilger" to compute the mixture fraction
1423 * in terms of the Bilger mixture fraction, or
1424 * an element name, to compute the mixture fraction
1425 * based on a single element (default: "Bilger")
1426 * @returns mixture fraction (kg fuel / kg mixture)
1427 */
1428 double mixtureFraction(const double* fuelComp, const double* oxComp,
1429 ThermoBasis basis=ThermoBasis::molar,
1430 const string& element="Bilger") const;
1431 //! @copydoc ThermoPhase::mixtureFraction
1432 double mixtureFraction(const string& fuelComp, const string& oxComp,
1433 ThermoBasis basis=ThermoBasis::molar,
1434 const string& element="Bilger") const;
1435 //! @copydoc ThermoPhase::mixtureFraction
1436 double mixtureFraction(const Composition& fuelComp, const Composition& oxComp,
1437 ThermoBasis basis=ThermoBasis::molar,
1438 const string& element="Bilger") const;
1439 //! @}
1440 //! @name Set Mixture Composition by Equivalence Ratio
1441 //! @{
1442
1443 //! Set the mixture composition according to the equivalence ratio.
1444 /*!
1445 * Fuel and oxidizer compositions are given either as
1446 * mole fractions or mass fractions (specified by `basis`)
1447 * and do not need to be normalized. Pressure and temperature are
1448 * kept constant. Elements C, S, H and O are considered for the oxidation.
1449 *
1450 * @param phi equivalence ratio
1451 * @param fuelComp composition of the fuel
1452 * @param oxComp composition of the oxidizer
1453 * @param basis either ThermoBasis::mole or ThermoBasis::mass.
1454 * Fuel and oxidizer composition are interpreted
1455 * as mole or mass fractions (default: molar)
1456 */
1457 void setEquivalenceRatio(double phi, const double* fuelComp, const double* oxComp,
1458 ThermoBasis basis=ThermoBasis::molar);
1459 //! @copydoc ThermoPhase::setEquivalenceRatio
1460 void setEquivalenceRatio(double phi, const string& fuelComp,
1461 const string& oxComp, ThermoBasis basis=ThermoBasis::molar);
1462 //! @copydoc ThermoPhase::setEquivalenceRatio
1463 void setEquivalenceRatio(double phi, const Composition& fuelComp,
1464 const Composition& oxComp, ThermoBasis basis=ThermoBasis::molar);
1465 //! @}
1466
1467 //! @name Compute Equivalence Ratio
1468 //! @{
1469
1470 //! Compute the equivalence ratio for the current mixture
1471 //! given the compositions of fuel and oxidizer
1472 /*!
1473 * The equivalence ratio @f$ \phi @f$ is computed from
1474 * @f[ \phi = \frac{Z}{1-Z}\frac{1-Z_{\mathrm{st}}}{Z_{\mathrm{st}}} @f]
1475 * where @f$ Z @f$ is the Bilger mixture fraction @cite bilger1979 of the mixture
1476 * given the specified fuel and oxidizer compositions
1477 * @f$ Z_{\mathrm{st}} @f$ is the mixture fraction at stoichiometric
1478 * conditions. Fuel and oxidizer compositions are given either as
1479 * mole fractions or mass fractions (specified by `basis`)
1480 * and do not need to be normalized.
1481 * Elements C, S, H and O are considered for the oxidation.
1482 * If fuel and oxidizer composition are unknown or not specified,
1483 * use the version that takes no arguments.
1484 *
1485 * @param fuelComp composition of the fuel
1486 * @param oxComp composition of the oxidizer
1487 * @param basis either ThermoPhase::mole or ThermoPhase::mass.
1488 * Fuel and oxidizer composition are interpreted
1489 * as mole or mass fractions (default: molar)
1490 * @returns equivalence ratio
1491 * @see mixtureFraction for the definition of the Bilger mixture fraction
1492 * @see equivalenceRatio() for the computation of @f$ \phi @f$ without arguments
1493 */
1494 double equivalenceRatio(const double* fuelComp, const double* oxComp,
1495 ThermoBasis basis=ThermoBasis::molar) const;
1496 //! @copydoc ThermoPhase::equivalenceRatio
1497 double equivalenceRatio(const string& fuelComp, const string& oxComp,
1498 ThermoBasis basis=ThermoBasis::molar) const;
1499 //! @copydoc ThermoPhase::equivalenceRatio
1500 double equivalenceRatio(const Composition& fuelComp,
1501 const Composition& oxComp, ThermoBasis basis=ThermoBasis::molar) const;
1502 //! @}
1503
1504 //! Compute the equivalence ratio for the current mixture
1505 //! from available oxygen and required oxygen
1506 /*!
1507 * Computes the equivalence ratio @f$ \phi @f$ from
1508 * @f[ \phi =
1509 * \frac{Z_{\mathrm{mole},C} + Z_{\mathrm{mole},S} + \frac{1}{4}Z_{\mathrm{mole},H}}
1510 * {\frac{1}{2}Z_{\mathrm{mole},O}} @f]
1511 * where @f$ Z_{\mathrm{mole},m} @f$ is the elemental mole fraction
1512 * of element @f$ m @f$. In this special case, the equivalence ratio
1513 * is independent of a fuel or oxidizer composition because it only
1514 * considers the locally available oxygen compared to the required oxygen
1515 * for complete oxidation. It is the same as assuming that the oxidizer
1516 * only contains O (and inert elements) and the fuel contains only
1517 * H, C and S (and inert elements). If either of these conditions is
1518 * not met, use the version of this functions which takes the fuel and
1519 * oxidizer compositions as input
1520 *
1521 * @returns equivalence ratio
1522 * @see equivalenceRatio compute the equivalence ratio from specific
1523 * fuel and oxidizer compositions
1524 */
1525 double equivalenceRatio() const;
1526
1527 //! @name Compute Stoichiometric Air to Fuel Ratio
1528 //! @{
1529
1530 //! Compute the stoichiometric air to fuel ratio (kg oxidizer / kg fuel)
1531 //! given fuel and oxidizer compositions.
1532 /*!
1533 * Fuel and oxidizer compositions are given either as
1534 * mole fractions or mass fractions (specified by `basis`)
1535 * and do not need to be normalized.
1536 * Elements C, S, H and O are considered for the oxidation.
1537 * Note that the stoichiometric air to fuel ratio @f$ \mathit{AFR}_{\mathrm{st}} @f$
1538 * does not depend on the current mixture composition. The current air to fuel ratio
1539 * can be computed from @f$ \mathit{AFR} = \mathit{AFR}_{\mathrm{st}}/\phi @f$
1540 * where @f$ \phi @f$ is the equivalence ratio of the current mixture
1541 *
1542 * @param fuelComp composition of the fuel
1543 * @param oxComp composition of the oxidizer
1544 * @param basis either ThermoPhase::mole or ThermoPhase::mass.
1545 * Fuel and oxidizer composition are interpreted
1546 * as mole or mass fractions (default: molar)
1547 * @returns Stoichiometric Air to Fuel Ratio (kg oxidizer / kg fuel)
1548 */
1549 double stoichAirFuelRatio(const double* fuelComp, const double* oxComp,
1550 ThermoBasis basis=ThermoBasis::molar) const;
1551 //! @copydoc ThermoPhase::stoichAirFuelRatio
1552 double stoichAirFuelRatio(const string& fuelComp, const string& oxComp,
1553 ThermoBasis basis=ThermoBasis::molar) const;
1554 //! @copydoc ThermoPhase::stoichAirFuelRatio
1555 double stoichAirFuelRatio(const Composition& fuelComp,
1556 const Composition& oxComp, ThermoBasis basis=ThermoBasis::molar) const;
1557 //! @}
1558
1559private:
1560
1561 //! Carry out work in HP and UV calculations.
1562 /*!
1563 * @param h Specific enthalpy or internal energy (J/kg)
1564 * @param p Pressure (Pa) or specific volume (m^3/kg)
1565 * @param tol Optional parameter setting the tolerance of the calculation.
1566 * Important for some applications where numerical Jacobians
1567 * are being calculated.
1568 * @param doUV True if solving for UV, false for HP.
1569 */
1570 void setState_HPorUV(double h, double p, double tol=1e-9, bool doUV = false);
1571
1572 //! Carry out work in SP and SV calculations.
1573 /*!
1574 * @param s Specific entropy (J/kg)
1575 * @param p Pressure (Pa) or specific volume (m^3/kg)
1576 * @param tol Optional parameter setting the tolerance of the calculation.
1577 * Important for some applications where numerical Jacobians
1578 * are being calculated.
1579 * @param doSV True if solving for SV, false for SP.
1580 */
1581 void setState_SPorSV(double s, double p, double tol=1e-9, bool doSV = false);
1582
1583 //! Helper function used by setState_HPorUV and setState_SPorSV.
1584 //! Sets the temperature and (if set_p is true) the pressure.
1585 void setState_conditional_TP(double t, double p, bool set_p);
1586
1587 //! Helper function for computing the amount of oxygen required for complete
1588 //! oxidation.
1589 /*!
1590 * @param y array of (possibly non-normalized) mass fractions (length m_kk)
1591 * @returns amount of required oxygen in kmol O / kg mixture
1592 */
1593 double o2Required(const double* y) const;
1594
1595 //! Helper function for computing the amount of oxygen
1596 //! available in the current mixture.
1597 /*!
1598 * @param y array of (possibly non-normalized) mass fractions (length m_kk)
1599 * @returns amount of O in kmol O / kg mixture
1600 */
1601 double o2Present(const double* y) const;
1602
1603public:
1604 //! @name Chemical Equilibrium
1605 //!
1606 //! Chemical equilibrium.
1607 //! @{
1608
1609 //! Equilibrate a ThermoPhase object
1610 /*!
1611 * Set this phase to chemical equilibrium by calling one of several
1612 * equilibrium solvers. The XY parameter indicates what two thermodynamic
1613 * quantities are to be held constant during the equilibration process.
1614 *
1615 * @param XY String representation of what two properties are being
1616 * held constant
1617 * @param solver Name of the solver to be used to equilibrate the phase.
1618 * If solver = 'element_potential', the ChemEquil element potential
1619 * solver will be used. If solver = 'vcs', the VCS solver will be used.
1620 * If solver = 'gibbs', the MultiPhaseEquil solver will be used. If
1621 * solver = 'auto', the solvers will be tried in order if the initial
1622 * solver(s) fail.
1623 * @param rtol Relative tolerance
1624 * @param max_steps Maximum number of steps to take to find the solution
1625 * @param max_iter For the 'gibbs' and 'vcs' solvers, this is the maximum
1626 * number of outer temperature or pressure iterations to take when T
1627 * and/or P is not held fixed.
1628 * @param estimate_equil For MultiPhaseEquil solver, an integer indicating
1629 * whether the solver should estimate its own initial condition. If 0,
1630 * the initial mole fraction vector in the ThermoPhase object is used
1631 * as the initial condition. If 1, the initial mole fraction vector is
1632 * used if the element abundances are satisfied. If -1, the initial
1633 * mole fraction vector is thrown out, and an estimate is formulated.
1634 * @param log_level loglevel Controls amount of diagnostic output.
1635 * log_level=0 suppresses diagnostics, and increasingly-verbose
1636 * messages are written as loglevel increases.
1637 *
1638 * @ingroup equilGroup
1639 */
1640 void equilibrate(const string& XY, const string& solver="auto",
1641 double rtol=1e-9, int max_steps=50000, int max_iter=100,
1642 int estimate_equil=0, int log_level=0);
1643
1644 //!This method is used by the ChemEquil equilibrium solver.
1645 /*!
1646 * It sets the state such that the chemical potentials satisfy
1647 * @f[ \frac{\mu_k}{\hat R T} = \sum_m A_{k,m}
1648 * \left(\frac{\lambda_m} {\hat R T}\right) @f] where
1649 * @f$ \lambda_m @f$ is the element potential of element m. The
1650 * temperature is unchanged. Any phase (ideal or not) that
1651 * implements this method can be equilibrated by ChemEquil.
1652 *
1653 * @param mu_RT Input vector of dimensionless chemical potentials
1654 * The length is equal to nSpecies().
1655 */
1656 virtual void setToEquilState(const double* mu_RT) {
1657 throw NotImplementedError("ThermoPhase::setToEquilState");
1658 }
1659
1660 //! Indicates whether this phase type can be used with class MultiPhase for
1661 //! equilibrium calculations. Returns `false` for special phase types which
1662 //! already represent multi-phase mixtures, namely PureFluidPhase.
1663 virtual bool compatibleWithMultiPhase() const {
1664 return true;
1665 }
1666
1667 //! @}
1668 //! @name Critical State Properties
1669 //!
1670 //! These methods are only implemented by subclasses that implement
1671 //! liquid-vapor equations of state.
1672 //! @{
1673
1674 //! Critical temperature (K).
1675 virtual double critTemperature() const {
1676 throw NotImplementedError("ThermoPhase::critTemperature");
1677 }
1678
1679 //! Critical pressure (Pa).
1680 virtual double critPressure() const {
1681 throw NotImplementedError("ThermoPhase::critPressure");
1682 }
1683
1684 //! Critical volume (m3/kmol).
1685 virtual double critVolume() const {
1686 throw NotImplementedError("ThermoPhase::critVolume");
1687 }
1688
1689 //! Critical compressibility (unitless).
1690 virtual double critCompressibility() const {
1691 throw NotImplementedError("ThermoPhase::critCompressibility");
1692 }
1693
1694 //! Critical density (kg/m3).
1695 virtual double critDensity() const {
1696 throw NotImplementedError("ThermoPhase::critDensity");
1697 }
1698
1699 //! @}
1700 //! @name Saturation Properties
1701 //!
1702 //! These methods are only implemented by subclasses that implement full
1703 //! liquid-vapor equations of state.
1704 //! @{
1705
1706 //! Return the saturation temperature given the pressure
1707 /*!
1708 * @param p Pressure (Pa)
1709 */
1710 virtual double satTemperature(double p) const {
1711 throw NotImplementedError("ThermoPhase::satTemperature");
1712 }
1713
1714 //! Return the saturation pressure given the temperature
1715 /*!
1716 * @param t Temperature (Kelvin)
1717 */
1718 virtual double satPressure(double t) {
1719 throw NotImplementedError("ThermoPhase::satPressure");
1720 }
1721
1722 //! Return the fraction of vapor at the current conditions
1723 virtual double vaporFraction() const {
1724 throw NotImplementedError("ThermoPhase::vaporFraction");
1725 }
1726
1727 //! Set the state to a saturated system at a particular temperature
1728 /*!
1729 * @param t Temperature (kelvin)
1730 * @param x Fraction of vapor
1731 */
1732 virtual void setState_Tsat(double t, double x) {
1733 throw NotImplementedError("ThermoPhase::setState_Tsat");
1734 }
1735
1736 //! Set the state to a saturated system at a particular pressure
1737 /*!
1738 * @param p Pressure (Pa)
1739 * @param x Fraction of vapor
1740 */
1741 virtual void setState_Psat(double p, double x) {
1742 throw NotImplementedError("ThermoPhase::setState_Psat");
1743 }
1744
1745 //! Set the temperature, pressure, and vapor fraction (quality).
1746 /*!
1747 * An exception is thrown if the thermodynamic state is not consistent.
1748 *
1749 * For temperatures below the critical temperature, if the vapor fraction is
1750 * not 0 or 1, the pressure and temperature must fall on the saturation
1751 * line.
1752 *
1753 * Above the critical temperature, the vapor fraction must be 1 if the
1754 * pressure is less than the critical pressure. Above the critical pressure,
1755 * the vapor fraction is not defined, and its value is ignored.
1756 *
1757 * @param T Temperature (K)
1758 * @param P Pressure (Pa)
1759 * @param Q vapor fraction
1760 */
1761 void setState_TPQ(double T, double P, double Q);
1762
1763 //! @}
1764 //! @name Initialization Methods - For Internal Use (ThermoPhase)
1765 //!
1766 //! The following methods are used in the process of constructing
1767 //! the phase and setting its parameters from a specification in an
1768 //! input file. They are not normally used in application programs.
1769 //! To see how they are used, see importPhase().
1770 //! @{
1771
1772 bool addSpecies(shared_ptr<Species> spec) override;
1773
1774 void modifySpecies(size_t k, shared_ptr<Species> spec) override;
1775
1776 //! Return a changeable reference to the calculation manager for species
1777 //! reference-state thermodynamic properties
1778 /*!
1779 * @param k Species id. The default is -1, meaning return the default
1780 */
1781 virtual MultiSpeciesThermo& speciesThermo(int k = -1);
1782
1783 virtual const MultiSpeciesThermo& speciesThermo(int k = -1) const;
1784
1785 /**
1786 * Initialize a ThermoPhase object using an input file.
1787 *
1788 * Used to implement constructors for derived classes which take a
1789 * file name and phase name as arguments.
1790 *
1791 * @param inputFile Input file containing the description of the phase. If blank,
1792 * no setup will be performed.
1793 * @param id Optional parameter identifying the name of the phase. If
1794 * blank, the first phase definition encountered will be used.
1795 */
1796 void initThermoFile(const string& inputFile, const string& id);
1797
1798 //! Initialize the ThermoPhase object after all species have been set up
1799 /*!
1800 * This method is provided to allow subclasses to perform any initialization
1801 * required after all species have been added. For example, it might be used
1802 * to resize internal work arrays that must have an entry for each species.
1803 * The base class implementation does nothing, and subclasses that do not
1804 * require initialization do not need to overload this method. Derived
1805 * classes which do override this function should call their parent class's
1806 * implementation of this function as their last action.
1807 *
1808 * When importing from an AnyMap phase description (or from a YAML file),
1809 * setupPhase() adds all the species, stores the input data in #m_input, and then
1810 * calls this method to set model parameters from the data stored in #m_input.
1811 */
1812 virtual void initThermo();
1813
1814 //! Set equation of state parameters from an AnyMap phase description.
1815 //! Phases that need additional parameters from the root node should
1816 //! override this method.
1817 virtual void setParameters(const AnyMap& phaseNode,
1818 const AnyMap& rootNode=AnyMap());
1819
1820 //! Returns the parameters of a ThermoPhase object such that an identical
1821 //! one could be reconstructed using the newThermo(AnyMap&) function.
1822 //! @param withInput If true, include additional input data fields associated
1823 //! with the phase description, such as user-defined fields from a YAML input
1824 //! file, as returned by the input() method.
1825 AnyMap parameters(bool withInput=true) const;
1826
1827 //! Get phase-specific parameters of a Species object such that an
1828 //! identical one could be reconstructed and added to this phase.
1829 /*!
1830 * @param name Name of the species
1831 * @param speciesNode Mapping to be populated with parameters
1832 */
1833 virtual void getSpeciesParameters(const string& name, AnyMap& speciesNode) const {}
1834
1835 //! Access input data associated with the phase description
1836 const AnyMap& input() const;
1837 AnyMap& input();
1838
1839 void invalidateCache() override;
1840
1841 //! @}
1842 //! @name Derivatives of Thermodynamic Variables needed for Applications
1843 //!
1844 //! Derivatives of the activity coefficients are needed to evaluate terms arising
1845 //! in multicomponent transport models for non-ideal systems. While %Cantera does
1846 //! not currently implement such models, these derivatives are provided by a few
1847 //! phase models.
1848 //! @{
1849
1850 //! Get the change in activity coefficients wrt changes in state (temp, mole
1851 //! fraction, etc) along a line in parameter space or along a line in
1852 //! physical space
1853 /*!
1854 * @param dTds Input of temperature change along the path
1855 * @param dXds Input vector of changes in mole fraction along the
1856 * path. length = m_kk Along the path length it must
1857 * be the case that the mole fractions sum to one.
1858 * @param dlnActCoeffds Output vector of the directional derivatives of the
1859 * log Activity Coefficients along the path. length =
1860 * m_kk units are 1/units(s). if s is a physical
1861 * coordinate then the units are 1/m.
1862 */
1863 virtual void getdlnActCoeffds(const double dTds, const double* const dXds,
1864 double* dlnActCoeffds) const {
1865 throw NotImplementedError("ThermoPhase::getdlnActCoeffds");
1866 }
1867
1868 //! Get the array of ln mole fraction derivatives of the log activity
1869 //! coefficients - diagonal component only
1870 /*!
1871 * For ideal mixtures (unity activity coefficients), this can return zero.
1872 * Implementations should take the derivative of the logarithm of the
1873 * activity coefficient with respect to the logarithm of the mole fraction
1874 * variable that represents the standard state. This quantity is to be used
1875 * in conjunction with derivatives of that mole fraction variable when the
1876 * derivative of the chemical potential is taken.
1877 *
1878 * units = dimensionless
1879 *
1880 * @param dlnActCoeffdlnX_diag Output vector of derivatives of the log
1881 * Activity Coefficients wrt the mole fractions. length = m_kk
1882 */
1883 virtual void getdlnActCoeffdlnX_diag(double* dlnActCoeffdlnX_diag) const {
1884 throw NotImplementedError("ThermoPhase::getdlnActCoeffdlnX_diag");
1885 }
1886
1887 //! Get the array of log species mole number derivatives of the log activity
1888 //! coefficients
1889 /*!
1890 * For ideal mixtures (unity activity coefficients), this can return zero.
1891 * Implementations should take the derivative of the logarithm of the
1892 * activity coefficient with respect to the logarithm of the concentration-
1893 * like variable (for example, moles) that represents the standard state. This
1894 * quantity is to be used in conjunction with derivatives of that species
1895 * mole number variable when the derivative of the chemical potential is
1896 * taken.
1897 *
1898 * units = dimensionless
1899 *
1900 * @param dlnActCoeffdlnN_diag Output vector of derivatives of the
1901 * log Activity Coefficients. length = m_kk
1902 */
1903 virtual void getdlnActCoeffdlnN_diag(double* dlnActCoeffdlnN_diag) const {
1904 throw NotImplementedError("ThermoPhase::getdlnActCoeffdlnN_diag");
1905 }
1906
1907 //! Get the array of derivatives of the log activity coefficients with
1908 //! respect to the log of the species mole numbers
1909 /*!
1910 * Implementations should take the derivative of the logarithm of the
1911 * activity coefficient with respect to a species log mole number (with all
1912 * other species mole numbers held constant). The default treatment in the
1913 * ThermoPhase object is to set this vector to zero.
1914 *
1915 * units = 1 / kmol
1916 *
1917 * dlnActCoeffdlnN[ ld * k + m] will contain the derivative of log
1918 * act_coeff for the *m*-th species with respect to the number of moles of
1919 * the *k*-th species.
1920 *
1921 * @f[
1922 * \frac{d \ln(\gamma_m) }{d \ln( n_k ) }\Bigg|_{n_i}
1923 * @f]
1924 *
1925 * When implemented, this method is used within the VCS equilibrium solver to
1926 * calculate the Jacobian elements, which accelerates convergence of the algorithm.
1927 *
1928 * @param ld Number of rows in the matrix
1929 * @param dlnActCoeffdlnN Output vector of derivatives of the
1930 * log Activity Coefficients. length = m_kk * m_kk
1931 */
1932 virtual void getdlnActCoeffdlnN(const size_t ld, double* const dlnActCoeffdlnN);
1933
1934 virtual void getdlnActCoeffdlnN_numderiv(const size_t ld,
1935 double* const dlnActCoeffdlnN);
1936
1937 //! @}
1938 //! @name Printing
1939 //! @{
1940
1941 //! returns a summary of the state of the phase as a string
1942 /*!
1943 * @param show_thermo If true, extra information is printed out
1944 * about the thermodynamic state of the system.
1945 * @param threshold Show information about species with mole fractions
1946 * greater than *threshold*.
1947 */
1948 virtual string report(bool show_thermo=true, double threshold=-1e-14) const;
1949
1950 //! @}
1951
1952protected:
1953 //! Store the parameters of a ThermoPhase object such that an identical
1954 //! one could be reconstructed using the newThermo(AnyMap&) function. This
1955 //! does not include user-defined fields available in input().
1956 virtual void getParameters(AnyMap& phaseNode) const;
1957
1958 //! Pointer to the calculation manager for species reference-state
1959 //! thermodynamic properties
1960 /*!
1961 * This class is called when the reference-state thermodynamic properties
1962 * of all the species in the phase needs to be evaluated.
1963 */
1965
1966 //! Data supplied via setParameters. When first set, this may include
1967 //! parameters used by different phase models when initThermo() is called.
1969
1970 //! Stored value of the electric potential for this phase. Units are Volts.
1971 double m_phi = 0.0;
1972
1973 //! Boolean indicating whether a charge neutrality condition is a necessity
1974 /*!
1975 * Note, the charge neutrality condition is not a necessity for ideal gas
1976 * phases. There may be a net charge in those phases, because the NASA
1977 * polynomials for ionized species in Ideal gases take this condition into
1978 * account. However, liquid phases usually require charge neutrality in
1979 * order for their derived thermodynamics to be valid.
1980 */
1982
1983 //! Contains the standard state convention
1985
1986 //! last value of the temperature processed by reference state
1987 mutable double m_tlast = 0.0;
1988};
1989
1990}
1991
1992#endif
Header for a general species thermodynamic property manager for a phase (see MultiSpeciesThermo).
Header file for class Phase.
Header for unit conversion utilities, which are used to translate user input from input files (See In...
A map of string keys to values whose type can vary at runtime.
Definition AnyMap.h:427
A species thermodynamic property manager for a phase.
virtual double refPressure() const
The reference-state pressure (Pa) for all species.
virtual double minTemp(size_t k=npos) const
Minimum temperature.
virtual double maxTemp(size_t k=npos) const
Maximum temperature.
virtual void modifyOneHf298(const size_t k, const double Hf298New)
Modify the value of the 298 K Heat of Formation of the standard state of one species in the phase (J ...
virtual double reportOneHf298(const size_t k) const
Report the 298 K Heat of Formation of the standard state of one species (J kmol-1)
An error indicating that an unimplemented function has been called.
Class Phase is the base class for phases of matter, managing the species and elements in a phase,...
Definition Phase.h:95
size_t m_kk
Number of species in the phase.
Definition Phase.h:842
double temperature() const
Temperature (K).
Definition Phase.h:562
double meanMolecularWeight() const
The mean molecular weight. Units: (kg/kmol)
Definition Phase.h:655
virtual double molarVolume() const
Molar volume (m^3/kmol).
Definition Phase.cpp:581
virtual double pressure() const
Return the thermodynamic pressure (Pa).
Definition Phase.h:580
string name() const
Return the name of the phase.
Definition Phase.cpp:20
Base class for a phase with thermodynamic properties.
int m_ssConvention
Contains the standard state convention.
virtual void getPartialMolarEnthalpies(double *hbar) const
Returns an array of partial molar enthalpies for the species in the mixture.
virtual double critTemperature() const
Critical temperature (K).
virtual void setState_HP(double h, double p, double tol=1e-9)
Set the internally stored specific enthalpy (J/kg) and pressure (Pa) of the phase.
double electricPotential() const
Returns the electric potential of this phase (V).
virtual void getEntropy_R(double *sr) const
Get the array of nondimensional Entropy functions for the standard state species at the current T and...
virtual void setState_UV(double u, double v, double tol=1e-9)
Set the specific internal energy (J/kg) and specific volume (m^3/kg).
bool chargeNeutralityNecessary() const
Returns the chargeNeutralityNecessity boolean.
virtual double cp_mole() const
Molar heat capacity at constant pressure. Units: J/kmol/K.
double equivalenceRatio() const
Compute the equivalence ratio for the current mixture from available oxygen and required oxygen.
virtual void setParameters(const AnyMap &phaseNode, const AnyMap &rootNode=AnyMap())
Set equation of state parameters from an AnyMap phase description.
virtual double thermalExpansionCoeff() const
Return the volumetric thermal expansion coefficient. Units: 1/K.
virtual void getEnthalpy_RT_ref(double *hrt) const
Returns the vector of nondimensional enthalpies of the reference state at the current temperature of ...
virtual void getParameters(AnyMap &phaseNode) const
Store the parameters of a ThermoPhase object such that an identical one could be reconstructed using ...
virtual double enthalpy_mole() const
Molar enthalpy. Units: J/kmol.
virtual void setState_TP(double t, double p)
Set the temperature (K) and pressure (Pa)
virtual double standardConcentration(size_t k=0) const
Return the standard concentration for the kth species.
virtual void getCp_R_ref(double *cprt) const
Returns the vector of nondimensional constant pressure heat capacities of the reference state at the ...
virtual void setState_TV(double t, double v, double tol=1e-9)
Set the temperature (K) and specific volume (m^3/kg).
virtual double logStandardConc(size_t k=0) const
Natural logarithm of the standard concentration of the kth species.
double o2Present(const double *y) const
Helper function for computing the amount of oxygen available in the current mixture.
virtual void setState_PV(double p, double v, double tol=1e-9)
Set the pressure (Pa) and specific volume (m^3/kg).
virtual void setState(const AnyMap &state)
Set the state using an AnyMap containing any combination of properties supported by the thermodynamic...
virtual double minTemp(size_t k=npos) const
Minimum temperature for which the thermodynamic data for the species or phase are valid.
virtual void getdlnActCoeffdlnN_diag(double *dlnActCoeffdlnN_diag) const
Get the array of log species mole number derivatives of the log activity coefficients.
virtual void setState_TPX(double t, double p, const double *x)
Set the temperature (K), pressure (Pa), and mole fractions.
void setState_SPorSV(double s, double p, double tol=1e-9, bool doSV=false)
Carry out work in SP and SV calculations.
double RT() const
Return the Gas Constant multiplied by the current temperature.
virtual void getPartialMolarCp(double *cpbar) const
Return an array of partial molar heat capacities for the species in the mixture.
virtual double critPressure() const
Critical pressure (Pa).
virtual void getGibbs_RT_ref(double *grt) const
Returns the vector of nondimensional Gibbs Free Energies of the reference state at the current temper...
virtual double soundSpeed() const
Return the speed of sound. Units: m/s.
virtual void setState_TPY(double t, double p, const double *y)
Set the internally stored temperature (K), pressure (Pa), and mass fractions of the phase.
double m_tlast
last value of the temperature processed by reference state
virtual void setState_ST(double s, double t, double tol=1e-9)
Set the specific entropy (J/kg/K) and temperature (K).
void setState_HPorUV(double h, double p, double tol=1e-9, bool doUV=false)
Carry out work in HP and UV calculations.
double gibbs_mass() const
Specific Gibbs function. Units: J/kg.
virtual void getActivityConcentrations(double *c) const
This method returns an array of generalized concentrations.
double stoichAirFuelRatio(const double *fuelComp, const double *oxComp, ThermoBasis basis=ThermoBasis::molar) const
Compute the stoichiometric air to fuel ratio (kg oxidizer / kg fuel) given fuel and oxidizer composit...
string type() const override
String indicating the thermodynamic model implemented.
AnyMap parameters(bool withInput=true) const
Returns the parameters of a ThermoPhase object such that an identical one could be reconstructed usin...
bool m_chargeNeutralityNecessary
Boolean indicating whether a charge neutrality condition is a necessity.
virtual void getPureGibbs(double *gpure) const
Get the Gibbs functions for the standard state of the species at the current T and P of the solution.
virtual string report(bool show_thermo=true, double threshold=-1e-14) const
returns a summary of the state of the phase as a string
virtual void getPartialMolarIntEnergies(double *ubar) const
Return an array of partial molar internal energies for the species in the mixture.
virtual void getIntEnergy_RT(double *urt) const
Returns the vector of nondimensional Internal Energies of the standard state species at the current T...
virtual void getCp_R(double *cpr) const
Get the nondimensional Heat Capacities at constant pressure for the species standard states at the cu...
virtual double maxTemp(size_t k=npos) const
Maximum temperature for which the thermodynamic data for the species are valid.
double m_phi
Stored value of the electric potential for this phase. Units are Volts.
virtual double isothermalCompressibility() const
Returns the isothermal compressibility. Units: 1/Pa.
double mixtureFraction(const double *fuelComp, const double *oxComp, ThermoBasis basis=ThermoBasis::molar, const string &element="Bilger") const
Compute the mixture fraction = kg fuel / (kg oxidizer + kg fuel) for the current mixture given fuel a...
double o2Required(const double *y) const
Helper function for computing the amount of oxygen required for complete oxidation.
virtual double satTemperature(double p) const
Return the saturation temperature given the pressure.
virtual void getdlnActCoeffds(const double dTds, const double *const dXds, double *dlnActCoeffds) const
Get the change in activity coefficients wrt changes in state (temp, mole fraction,...
void getElectrochemPotentials(double *mu) const
Get the species electrochemical potentials.
virtual void getdlnActCoeffdlnN(const size_t ld, double *const dlnActCoeffdlnN)
Get the array of derivatives of the log activity coefficients with respect to the log of the species ...
virtual void getGibbs_RT(double *grt) const
Get the nondimensional Gibbs functions for the species in their standard states at the current T and ...
virtual double critVolume() const
Critical volume (m3/kmol).
virtual void getActivityCoefficients(double *ac) const
Get the array of non-dimensional molar-based activity coefficients at the current solution temperatur...
virtual string phaseOfMatter() const
String indicating the mechanical phase of the matter in this Phase.
virtual void getStandardVolumes(double *vol) const
Get the molar volumes of the species standard states at the current T and P of the solution.
virtual void setState_Tsat(double t, double x)
Set the state to a saturated system at a particular temperature.
virtual double entropy_mole() const
Molar entropy. Units: J/kmol/K.
void setElectricPotential(double v)
Set the electric potential of this phase (V).
double cv_mass() const
Specific heat at constant volume. Units: J/kg/K.
virtual int activityConvention() const
This method returns the convention used in specification of the activities, of which there are curren...
virtual void initThermo()
Initialize the ThermoPhase object after all species have been set up.
double entropy_mass() const
Specific entropy. Units: J/kg/K.
virtual double critDensity() const
Critical density (kg/m3).
virtual void getGibbs_ref(double *g) const
Returns the vector of the Gibbs function of the reference state at the current temperature of the sol...
virtual MultiSpeciesThermo & speciesThermo(int k=-1)
Return a changeable reference to the calculation manager for species reference-state thermodynamic pr...
virtual void setState_UP(double u, double p, double tol=1e-9)
Set the specific internal energy (J/kg) and pressure (Pa).
void initThermoFile(const string &inputFile, const string &id)
Initialize a ThermoPhase object using an input file.
virtual void setState_SP(double s, double p, double tol=1e-9)
Set the specific entropy (J/kg/K) and pressure (Pa).
virtual void modifyOneHf298SS(const size_t k, const double Hf298New)
Modify the value of the 298 K Heat of Formation of one species in the phase (J kmol-1)
virtual int standardStateConvention() const
This method returns the convention used in specification of the standard state, of which there are cu...
void modifySpecies(size_t k, shared_ptr< Species > spec) override
Modify the thermodynamic data associated with a species.
virtual void setState_SH(double s, double h, double tol=1e-9)
Set the specific entropy (J/kg/K) and the specific enthalpy (J/kg)
virtual void getdlnActCoeffdlnX_diag(double *dlnActCoeffdlnX_diag) const
Get the array of ln mole fraction derivatives of the log activity coefficients - diagonal component o...
void invalidateCache() override
Invalidate any cached values which are normally updated only when a change in state is detected.
virtual void getActivities(double *a) const
Get the array of non-dimensional activities at the current solution temperature, pressure,...
void setMixtureFraction(double mixFrac, const double *fuelComp, const double *oxComp, ThermoBasis basis=ThermoBasis::molar)
Set the mixture composition according to the mixture fraction = kg fuel / (kg oxidizer + kg fuel)
virtual void getStandardVolumes_ref(double *vol) const
Get the molar volumes of the species reference states at the current T and P_ref of the solution.
virtual double vaporFraction() const
Return the fraction of vapor at the current conditions.
virtual void resetHf298(const size_t k=npos)
Restore the original heat of formation of one or more species.
virtual void getStandardChemPotentials(double *mu) const
Get the array of chemical potentials at unit activity for the species at their standard states at the...
virtual void getEnthalpy_RT(double *hrt) const
Get the nondimensional Enthalpy functions for the species at their standard states at the current T a...
virtual void getEntropy_R_ref(double *er) const
Returns the vector of nondimensional entropies of the reference state at the current temperature of t...
virtual void getChemPotentials(double *mu) const
Get the species chemical potentials. Units: J/kmol.
double cp_mass() const
Specific heat at constant pressure. Units: J/kg/K.
virtual void setState_TH(double t, double h, double tol=1e-9)
Set the temperature (K) and the specific enthalpy (J/kg)
virtual void getLnActivityCoefficients(double *lnac) const
Get the array of non-dimensional molar-based ln activity coefficients at the current solution tempera...
double intEnergy_mass() const
Specific internal energy. Units: J/kg.
virtual void getSpeciesParameters(const string &name, AnyMap &speciesNode) const
Get phase-specific parameters of a Species object such that an identical one could be reconstructed a...
virtual Units standardConcentrationUnits() const
Returns the units of the "standard concentration" for this phase.
virtual void getIntEnergy_RT_ref(double *urt) const
Returns the vector of nondimensional internal Energies of the reference state at the current temperat...
double Hf298SS(const size_t k) const
Report the 298 K Heat of Formation of the standard state of one species (J kmol-1)
ThermoPhase()=default
Constructor.
virtual bool isIdeal() const
Boolean indicating whether phase is ideal.
virtual double cv_mole() const
Molar heat capacity at constant volume. Units: J/kmol/K.
MultiSpeciesThermo m_spthermo
Pointer to the calculation manager for species reference-state thermodynamic properties.
virtual double satPressure(double t)
Return the saturation pressure given the temperature.
virtual double refPressure() const
Returns the reference pressure in Pa.
virtual double critCompressibility() const
Critical compressibility (unitless).
bool addSpecies(shared_ptr< Species > spec) override
Add a Species to this Phase.
AnyMap m_input
Data supplied via setParameters.
virtual double intEnergy_mole() const
Molar internal energy. Units: J/kmol.
virtual void setState_DP(double rho, double p)
Set the density (kg/m**3) and pressure (Pa) at constant composition.
void setEquivalenceRatio(double phi, const double *fuelComp, const double *oxComp, ThermoBasis basis=ThermoBasis::molar)
Set the mixture composition according to the equivalence ratio.
void setState_TPQ(double T, double P, double Q)
Set the temperature, pressure, and vapor fraction (quality).
virtual void setState_VH(double v, double h, double tol=1e-9)
Set the specific volume (m^3/kg) and the specific enthalpy (J/kg)
virtual void getPartialMolarEntropies(double *sbar) const
Returns an array of partial molar entropies of the species in the solution.
virtual double gibbs_mole() const
Molar Gibbs function. Units: J/kmol.
virtual void setState_SV(double s, double v, double tol=1e-9)
Set the specific entropy (J/kg/K) and specific volume (m^3/kg).
const AnyMap & input() const
Access input data associated with the phase description.
virtual void setState_Psat(double p, double x)
Set the state to a saturated system at a particular pressure.
void setState_conditional_TP(double t, double p, bool set_p)
Helper function used by setState_HPorUV and setState_SPorSV.
virtual void getPartialMolarVolumes(double *vbar) const
Return an array of partial molar volumes for the species in the mixture.
double enthalpy_mass() const
Specific enthalpy. Units: J/kg.
A representation of the units associated with a dimensional quantity.
Definition Units.h:35
void equilibrate(const string &XY, const string &solver="auto", double rtol=1e-9, int max_steps=50000, int max_iter=100, int estimate_equil=0, int log_level=0)
Equilibrate a ThermoPhase object.
virtual bool compatibleWithMultiPhase() const
Indicates whether this phase type can be used with class MultiPhase for equilibrium calculations.
virtual void setToEquilState(const double *mu_RT)
This method is used by the ChemEquil equilibrium solver.
const double GasConstant
Universal Gas Constant [J/kmol/K].
Definition ct_defs.h:120
Namespace for the Cantera kernel.
Definition AnyMap.cpp:564
const size_t npos
index returned by functions to indicate "no position"
Definition ct_defs.h:180
const int cSS_CONVENTION_VPSS
Standard state uses the molality convention.
const int cAC_CONVENTION_MOLAR
Standard state uses the molar convention.
const int cSS_CONVENTION_TEMPERATURE
Standard state uses the molar convention.
ThermoBasis
Differentiate between mole fractions and mass fractions for input mixture composition.
const int cSS_CONVENTION_SLAVE
Standard state thermodynamics is obtained from slave ThermoPhase objects.
map< string, double > Composition
Map from string names to doubles.
Definition ct_defs.h:177
const int cAC_CONVENTION_MOLALITY
Standard state uses the molality convention.