Cantera  3.2.0a4
Loading...
Searching...
No Matches
Phase.h
Go to the documentation of this file.
1/**
2 * @file Phase.h
3 * Header file for class Phase.
4 */
5
6// This file is part of Cantera. See License.txt in the top-level directory or
7// at https://cantera.org/license.txt for license and copyright information.
8
9#ifndef CT_PHASE_H
10#define CT_PHASE_H
11
15
16namespace Cantera
17{
18
19class Solution;
20class Species;
21class Kinetics;
22
23//! Class Phase is the base class for phases of matter, managing the species and
24//! elements in a phase, as well as the independent variables of temperature,
25//! mass density (compressible substances) or pressure (incompressible
26//! substances), species mass/mole fraction, and other generalized forces and
27//! intrinsic properties (such as electric potential) that define the
28//! thermodynamic state.
29/*!
30 *
31 * Class Phase provides information about the elements and species in a
32 * phase - names, index numbers (location in arrays), atomic or molecular
33 * weights, etc. The set of elements must include all those that compose the
34 * species, but may include additional elements.
35 *
36 * It also stores an array of species molecular weights, which are used to
37 * convert between mole and mass representations of the composition. For
38 * efficiency in mass/mole conversion, the vector of mass fractions divided
39 * by molecular weight @f$ Y_k/M_k @f$ is also stored.
40 *
41 * Class Phase is not usually used directly. Its primary use is as a base class
42 * for class ThermoPhase. It is not generally necessary to overloaded any of
43 * class Phase's methods, which handles both compressible and incompressible
44 * phases. For incompressible phases, the density is replaced by the pressure
45 * as the independent variable, and can no longer be set directly. In this case,
46 * the density needs to be calculated from a suitable equation of state, and
47 * assigned to the object using the assignDensity() method. This also applies
48 * for nearly-incompressible phases or phases which utilize standard states
49 * based on a T and P, in which case they need to overload these functions too.
50 *
51 * Class Phase contains a number of utility functions that will set the state
52 * of the phase in its entirety, by first setting the composition, and then
53 * temperature and pressure. An example of this is the function
54 * Phase::setState_TPY(double t, double p, const double* y).
55 *
56 * For bulk (3-dimensional) phases, the mass density has units of kg/m^3, and the molar
57 * density and concentrations have units of kmol/m^3, and the units listed in the
58 * methods of the Phase class assume a bulk phase. However, for surface (2-dimensional)
59 * phases have units of kg/m^2 and kmol/m^2, respectively. And for edge (1-dimensional)
60 * phases, these units kg/m and kmol/m.
61 *
62 * Class Phase contains methods for saving and restoring the full internal state
63 * of a given phase. These are saveState() and restoreState(). These functions
64 * operate on a state vector, which by default uses the first two entries for
65 * temperature and density (compressible substances) or temperature and
66 * pressure (incompressible substances). If the substance is not pure in a
67 * thermodynamic sense (that is, it may contain multiple species), the state also
68 * contains nSpecies() entries that specify the composition by corresponding
69 * mass fractions. Default definitions can be overloaded by derived classes.
70 * For any phase, the native definition of its thermodynamic state is defined
71 * the method nativeState(), with the length of the state vector returned by
72 * by stateSize(). In addition, methods isPure() and isCompressible() provide
73 * information on the implementation of a Phase object.
74 *
75 * A species name is referred to via speciesName(), which is unique within a
76 * given phase. Note that within multiphase mixtures (MultiPhase()), both a
77 * phase name/index as well as species name are required to access information
78 * about a species in a particular phase. For surfaces, the species names are
79 * unique among the phases.
80 *
81 * @todo
82 * - Specify that the input mole, mass, and volume fraction vectors must sum
83 * to one on entry to the set state routines. Non-conforming mole/mass
84 * fraction vectors are not thermodynamically consistent. Moreover, unless
85 * we do this, the calculation of Jacobians will be altered whenever the
86 * treatment of non- conforming mole fractions is changed. Add setState
87 * functions corresponding to specifying mole numbers, which is actually
88 * what is being done (well one of the options, there are many) when non-
89 * conforming mole fractions are input. Note, we realize that most numerical
90 * Jacobian and some analytical Jacobians use non-conforming calculations.
91 * These can easily be changed to the set mole number setState functions.
92 *
93 * @ingroup thermoprops
94 */
95class Phase
96{
97public:
98 Phase() = default; //!< Default constructor.
99 virtual ~Phase() = default;
100
101 // Phase objects are not copyable or assignable
102 Phase(const Phase&) = delete;
103 Phase& operator=(const Phase&) = delete;
104
105 /**
106 * @name Name
107 * Class Phase uses the string name to identify a phase. For phases instantiated
108 * from YAML input files, the name is the value of the corresponding key in the
109 * phase map.
110 *
111 * However, the name field may be changed to another value during the
112 * course of a calculation. For example, if duplicates of a phase object
113 * are instantiated and used in multiple places (such as a ReactorNet), they
114 * will have the same constitutive input, that is, the names of the phases will
115 * be the same. Note that this is not a problem for %Cantera internally;
116 * however, a user may want to rename phase objects in order to clarify.
117 */
118 //!@{
119
120 //! Return the name of the phase.
121 /*!
122 * Names are unique within a %Cantera problem.
123 */
124 string name() const;
125
126 //! Sets the string name for the phase.
127 //! @param nm String name of the phase
128 void setName(const string& nm);
129
130 //! String indicating the thermodynamic model implemented. Usually
131 //! corresponds to the name of the derived class, less any suffixes such as
132 //! "Phase", TP", "VPSS", etc.
133 //! @since Starting in %Cantera 3.0, the name returned by this method corresponds
134 //! to the canonical name used in the YAML input format.
135 virtual string type() const {
136 return "Phase";
137 }
138
139 //! @} end group Name
140
141 //! @name Element and Species Information
142 //! @{
143
144 //! Name of the element with index m.
145 //! @param m Element index.
146 string elementName(size_t m) const;
147
148 //! Return the index of element named 'name'. The index is an integer
149 //! assigned to each element in the order it was added. Returns @ref npos
150 //! if the specified element is not found.
151 //! @param name Name of the element
152 size_t elementIndex(const string& name) const;
153
154 //! Return a read-only reference to the vector of element names.
155 const vector<string>& elementNames() const;
156
157 //! Atomic weight of element m.
158 //! @param m Element index
159 double atomicWeight(size_t m) const;
160
161 //! Entropy of the element in its standard state at 298 K and 1 bar.
162 //! If no entropy value was provided when the phase was constructed,
163 //! returns the value `ENTROPY298_UNKNOWN`.
164 //! @param m Element index
165 double entropyElement298(size_t m) const;
166
167 //! Atomic number of element m.
168 //! @param m Element index
169 int atomicNumber(size_t m) const;
170
171 //! Return the element constraint type
172 //! Possible types include:
173 //!
174 //! - `CT_ELEM_TYPE_TURNEDOFF -1`
175 //! - `CT_ELEM_TYPE_ABSPOS 0`
176 //! - `CT_ELEM_TYPE_ELECTRONCHARGE 1`
177 //! - `CT_ELEM_TYPE_CHARGENEUTRALITY 2`
178 //! - `CT_ELEM_TYPE_LATTICERATIO 3`
179 //! - `CT_ELEM_TYPE_KINETICFROZEN 4`
180 //! - `CT_ELEM_TYPE_SURFACECONSTRAINT 5`
181 //! - `CT_ELEM_TYPE_OTHERCONSTRAINT 6`
182 //!
183 //! The default is `CT_ELEM_TYPE_ABSPOS`.
184 //! @param m Element index
185 //! @returns the element type
186 int elementType(size_t m) const;
187
188 //! Change the element type of the mth constraint
189 //! Reassigns an element type.
190 //! @param m Element index
191 //! @param elem_type New elem type to be assigned
192 //! @returns the old element type
193 int changeElementType(int m, int elem_type);
194
195 //! Return a read-only reference to the vector of atomic weights.
196 const vector<double>& atomicWeights() const;
197
198 //! Number of elements.
199 size_t nElements() const;
200
201 //! Check that the specified element index is in range.
202 //! Throws an exception if m is greater than nElements()-1
203 void checkElementIndex(size_t m) const;
204
205 //! Check that an array size is at least nElements().
206 //! Throws an exception if mm is less than nElements(). Used before calls
207 //! which take an array pointer.
208 void checkElementArraySize(size_t mm) const;
209
210 //! Number of atoms of element @c m in species @c k.
211 //! @param k species index
212 //! @param m element index
213 double nAtoms(size_t k, size_t m) const;
214
215 //! Returns the index of a species named 'name' within the Phase object.
216 //! The first species in the phase will have an index 0, and the last one
217 //! will have an index of nSpecies() - 1.
218 //! @param name String name of the species. It may also be in the form
219 //! phaseName:speciesName
220 //! @return The index of the species. If the name is not found,
221 //! the value @ref npos is returned.
222 size_t speciesIndex(const string& name) const;
223
224 //! Name of the species with index k
225 //! @param k index of the species
226 string speciesName(size_t k) const;
227
228 //! Return a const reference to the vector of species names
229 const vector<string>& speciesNames() const;
230
231 //! Returns the number of species in the phase
232 size_t nSpecies() const {
233 return m_kk;
234 }
235
236 //! Check that the specified species index is in range.
237 //! Throws an exception if k is greater than nSpecies()-1
238 void checkSpeciesIndex(size_t k) const;
239
240 //! Check that an array size is at least nSpecies().
241 //! Throws an exception if kk is less than nSpecies(). Used before calls
242 //! which take an array pointer.
243 void checkSpeciesArraySize(size_t kk) const;
244
245 //! @} end group Element and Species Information
246
247 //! Return whether phase represents a pure (single species) substance
248 virtual bool isPure() const {
249 return false;
250 }
251
252 //! Return whether phase represents a substance with phase transitions
253 virtual bool hasPhaseTransition() const {
254 return false;
255 }
256
257 //! Return whether phase represents a compressible substance
258 virtual bool isCompressible() const {
259 return true;
260 }
261
262 //! Return a map of properties defining the native state of a substance.
263 //! By default, entries include "T", "D", "Y" for a compressible substance
264 //! and "T", "P", "Y" for an incompressible substance, with offsets 0, 1 and
265 //! 2, respectively. Mass fractions "Y" are omitted for pure species.
266 //! In all cases, offsets into the state vector are used by saveState()
267 //! and restoreState().
268 virtual map<string, size_t> nativeState() const;
269
270 //! Return string acronym representing the native state of a Phase.
271 //! Examples: "TP", "TDY", "TPY".
272 //! @see nativeState
273 //! @since New in %Cantera 3.0
274 string nativeMode() const;
275
276 //! Return a vector containing full states defining a phase.
277 //! Full states list combinations of properties that allow for the
278 //! specification of a thermodynamic state based on user input.
279 //! Properties and states are represented by single letter acronyms, and
280 //! combinations of letters, respectively (for example, "TDY", "TPX", "SVX").
281 //! Supported property acronyms are:
282 //! "T": temperature
283 //! "P": pressure
284 //! "D": density
285 //! "X": mole fractions
286 //! "Y": mass fractions
287 //! "T": temperature
288 //! "U": specific internal energy
289 //! "V": specific volume
290 //! "H": specific enthalpy
291 //! "S": specific entropy
292 //! "Q": vapor fraction
293 virtual vector<string> fullStates() const;
294
295 //! Return a vector of settable partial property sets within a phase.
296 //! Partial states encompass all valid combinations of properties that allow
297 //! for the specification of a state while ignoring species concentrations
298 //! (such as "TD", "TP", "SV").
299 virtual vector<string> partialStates() const;
300
301 //! Get the size of the partial state vector of the phase.
302 //! The partial state vector excludes composition. Vectors of this size are used by
303 //! savePartialState() and restorePartialState().
304 //! @since New in %Cantera 3.2
305 virtual size_t partialStateSize() const { return 2; }
306
307 //! Save the current thermodynamic state of the phase, excluding composition.
308 //! The default implementation corresponds to the default implementation of
309 //! nativeState().
310 //! @param lenstate Length of the state array. Must be >= partialStateSize().
311 //! @param[out] state Array of state variables, in the order defined by
312 //! nativeState().
313 //! @since New in %Cantera 3.2
314 virtual void savePartialState(size_t lenstate, double* state) const;
315
316 //! Set the internal thermodynamic state of the phase, excluding composition.
317 //! The default implementation corresponds to the default implementation of
318 //! nativeState().
319 //! @param lenstate Length of the state array. Must be >= partialStateSize().
320 //! @param[in] state Array of state variables, in the order defined by
321 //! nativeState().
322 //! @since New in %Cantera 3.2
323 virtual void restorePartialState(size_t lenstate, const double* state);
324
325 //! Return size of vector defining internal state of the phase.
326 //! Used by saveState() and restoreState().
327 virtual size_t stateSize() const;
328
329 //! Save the current internal state of the phase.
330 //! Write to vector 'state' the current internal state.
331 //! @param state output vector. Will be resized to stateSize().
332 void saveState(vector<double>& state) const;
333
334 //! Write to array 'state' the current internal state.
335 //! @param lenstate length of the state array. Must be >= stateSize()
336 //! @param state output vector. Must be of length stateSizes() or
337 //! greater.
338 virtual void saveState(size_t lenstate, double* state) const;
339
340 //! Restore a state saved on a previous call to saveState.
341 //! @param state State vector containing the previously saved state.
342 void restoreState(const vector<double>& state);
343
344 //! Restore the state of the phase from a previously saved state vector.
345 //! @param lenstate Length of the state vector
346 //! @param state Vector of state conditions.
347 virtual void restoreState(size_t lenstate, const double* state);
348
349 //! @name Set Thermodynamic State
350 //!
351 //! Set the internal thermodynamic state by setting the internally stored
352 //! temperature, density and species composition. Note that the composition
353 //! is always set first.
354 //!
355 //! Temperature and density are held constant if not explicitly set.
356 //! @{
357
358 //! Set the species mole fractions by name.
359 //! Species not listed by name in @c xMap are set to zero.
360 //! @param xMap map from species names to mole fraction values.
361 void setMoleFractionsByName(const Composition& xMap);
362
363 //! Set the mole fractions of a group of species by name. Species which
364 //! are not listed by name in the composition map are set to zero.
365 //! @param x string x in the form of a composition map
366 void setMoleFractionsByName(const string& x);
367
368 //! Set the species mass fractions by name.
369 //! Species not listed by name in @c yMap are set to zero.
370 //! @param yMap map from species names to mass fraction values.
371 void setMassFractionsByName(const Composition& yMap);
372
373 //! Set the species mass fractions by name.
374 //! Species not listed by name in @c x are set to zero.
375 //! @param x String containing a composition map
376 void setMassFractionsByName(const string& x);
377
378 //! Set the internally stored temperature (K) and density (kg/m^3)
379 //! @param t Temperature in kelvin
380 //! @param rho Density (kg/m^3)
381 //! @since New in %Cantera 3.0.
382 void setState_TD(double t, double rho);
383
384 //! @} end group set thermo state
385
386 //! Molecular weight of species @c k.
387 //! @param k index of species @c k
388 //! @returns the molecular weight of species @c k.
389 double molecularWeight(size_t k) const;
390
391 //! Copy the vector of molecular weights into array weights.
392 //! @param weights Output array of molecular weights (kg/kmol)
393 void getMolecularWeights(double* weights) const;
394
395 //! Return a const reference to the internal vector of molecular weights.
396 //! units = kg / kmol
397 const vector<double>& molecularWeights() const;
398
399 //! Return a const reference to the internal vector of molecular weights.
400 //! units = kmol / kg
401 const vector<double>& inverseMolecularWeights() const;
402
403 //! Copy the vector of species charges into array charges.
404 //! @param charges Output array of species charges (elem. charge)
405 void getCharges(double* charges) const;
406
407 //! @name Composition
408 //! @{
409
410 //! Get the mole fractions by name.
411 //! @param threshold Exclude species with mole fractions less than or
412 //! equal to this threshold.
413 //! @return Map of species names to mole fractions
414 Composition getMoleFractionsByName(double threshold=0.0) const;
415
416 //! Return the mole fraction of a single species
417 //! @param k species index
418 //! @return Mole fraction of the species
419 double moleFraction(size_t k) const;
420
421 //! Return the mole fraction of a single species
422 //! @param name String name of the species
423 //! @return Mole fraction of the species
424 double moleFraction(const string& name) const;
425
426 //! Get the mass fractions by name.
427 //! @param threshold Exclude species with mass fractions less than or
428 //! equal to this threshold.
429 //! @return Map of species names to mass fractions
430 Composition getMassFractionsByName(double threshold=0.0) const;
431
432 //! Return the mass fraction of a single species
433 //! @param k species index
434 //! @return Mass fraction of the species
435 double massFraction(size_t k) const;
436
437 //! Return the mass fraction of a single species
438 //! @param name String name of the species
439 //! @return Mass Fraction of the species
440 double massFraction(const string& name) const;
441
442 //! Get the species mole fraction vector.
443 //! @param x On return, x contains the mole fractions. Must have a
444 //! length greater than or equal to the number of species.
445 void getMoleFractions(double* const x) const;
446
447 //! Set the mole fractions to the specified values.
448 //! There is no restriction on the sum of the mole fraction vector.
449 //! Internally, the Phase object will normalize this vector before storing
450 //! its contents.
451 //! @param x Array of unnormalized mole fraction values (input). Must
452 //! have a length greater than or equal to the number of species, m_kk.
453 virtual void setMoleFractions(const double* const x);
454
455 //! Set the mole fractions to the specified values without normalizing.
456 //! This is useful when the normalization condition is being handled by
457 //! some other means, for example by a constraint equation as part of a
458 //! larger set of equations.
459 //! @param x Input vector of mole fractions. Length is m_kk.
460 virtual void setMoleFractions_NoNorm(const double* const x);
461
462 //! Get the species mass fractions.
463 //! @param[out] y Array of mass fractions, length nSpecies()
464 void getMassFractions(double* const y) const;
465
466 //! Return a const pointer to the mass fraction array
467 const double* massFractions() const {
468 return &m_y[0];
469 }
470
471 //! Set the mass fractions to the specified values and normalize them.
472 //! @param[in] y Array of unnormalized mass fraction values. Length
473 //! must be greater than or equal to the number of
474 //! species. The Phase object will normalize this vector
475 //! before storing its contents.
476 virtual void setMassFractions(const double* const y);
477
478 //! Set the mass fractions to the specified values without normalizing.
479 //! This is useful when the normalization condition is being handled by
480 //! some other means, for example by a constraint equation as part of a
481 //! larger set of equations.
482 //! @param y Input vector of mass fractions. Length is m_kk.
483 virtual void setMassFractions_NoNorm(const double* const y);
484
485 //! Get the species concentrations (kmol/m^3).
486 /*!
487 * @param[out] c The vector of species concentrations. Units are
488 * kmol/m^3. The length of the vector must be greater than
489 * or equal to the number of species within the phase.
490 */
491 virtual void getConcentrations(double* const c) const;
492
493 //! Concentration of species k.
494 //! If k is outside the valid range, an exception will be thrown.
495 /*!
496 * @param[in] k Index of the species within the phase.
497 *
498 * @returns the concentration of species k (kmol m-3).
499 */
500 virtual double concentration(const size_t k) const;
501
502 //! Set the concentrations to the specified values within the phase.
503 //! We set the concentrations here and therefore we set the overall density
504 //! of the phase. We hold the temperature constant during this operation.
505 //! Therefore, we have possibly changed the pressure of the phase by
506 //! calling this routine.
507 //! @param[in] conc Array of concentrations in dimensional units. For
508 //! bulk phases c[k] is the concentration of the kth
509 //! species in kmol/m3. For surface phases, c[k] is the
510 //! concentration in kmol/m2. The length of the vector
511 //! is the number of species in the phase.
512 virtual void setConcentrations(const double* const conc);
513
514 //! Set the concentrations without ignoring negative concentrations
515 virtual void setConcentrationsNoNorm(const double* const conc);
516 //! @}
517
518 //! Set the state of the object with moles in [kmol]
519 virtual void setMolesNoTruncate(const double* const N);
520
521 //! Elemental mass fraction of element m
522 /*!
523 * The elemental mass fraction @f$ Z_{\mathrm{mass},m} @f$ of element @f$ m @f$
524 * is defined as
525 * @f[
526 * Z_{\mathrm{mass},m} = \sum_k \frac{a_{m,k} M_m}{M_k} Y_k
527 * @f]
528 * with @f$ a_{m,k} @f$ being the number of atoms of element @f$ m @f$ in
529 * species @f$ k @f$, @f$ M_m @f$ the atomic weight of element @f$ m @f$,
530 * @f$ M_k @f$ the molecular weight of species @f$ k @f$, and @f$ Y_k @f$
531 * the mass fraction of species @f$ k @f$.
532 *
533 * @param[in] m Index of the element within the phase. If m is outside
534 * the valid range, an exception will be thrown.
535 *
536 * @return the elemental mass fraction of element m.
537 */
538 double elementalMassFraction(const size_t m) const;
539
540 //! Elemental mole fraction of element m
541 /*!
542 * The elemental mole fraction @f$ Z_{\mathrm{mole},m} @f$ of element @f$ m @f$
543 * is the number of atoms of element *m* divided by the total number of
544 * atoms. It is defined as:
545 *
546 * @f[
547 * Z_{\mathrm{mole},m} = \frac{\sum_k a_{m,k} X_k}
548 * {\sum_k \sum_j a_{j,k} X_k}
549 * @f]
550 * with @f$ a_{m,k} @f$ being the number of atoms of element @f$ m @f$ in
551 * species @f$ k @f$, @f$ \sum_j @f$ being a sum over all elements, and
552 * @f$ X_k @f$ being the mole fraction of species @f$ k @f$.
553 *
554 * @param[in] m Index of the element within the phase. If m is outside the
555 * valid range, an exception will be thrown.
556 * @return the elemental mole fraction of element m.
557 */
558 double elementalMoleFraction(const size_t m) const;
559
560 //! Dimensionless electrical charge of a single molecule of species k
561 //! The charge is normalized by the the magnitude of the electron charge
562 //! @param k species index
563 double charge(size_t k) const {
564 return m_speciesCharge[k];
565 }
566
567 //! Charge density [C/m^3].
568 double chargeDensity() const;
569
570 //! Returns the number of spatial dimensions (1, 2, or 3)
571 size_t nDim() const {
572 return m_ndim;
573 }
574
575 //! Set the number of spatial dimensions (1, 2, or 3). The number of
576 //! spatial dimensions is used for vector involving directions.
577 //! @param ndim Input number of dimensions.
578 void setNDim(size_t ndim) {
579 m_ndim = ndim;
580 }
581
582 //! @name Thermodynamic Properties
583 //! @{
584
585 //! Temperature (K).
586 //! @return The temperature of the phase
587 double temperature() const {
588 return m_temp;
589 }
590
591 //! Electron Temperature (K)
592 //! @return The electron temperature of the phase
593 virtual double electronTemperature() const {
594 return m_temp;
595 }
596
597 //! Return the thermodynamic pressure (Pa).
598 /*!
599 * This method must be overloaded in derived classes. Within %Cantera, the
600 * independent variable is either density or pressure. If the state is
601 * defined by temperature, density, and mass fractions, this method should
602 * use these values to implement the mechanical equation of state @f$ P(T,
603 * \rho, Y_1, \dots, Y_K) @f$. Alternatively, it returns a stored value.
604 */
605 virtual double pressure() const {
606 throw NotImplementedError("Phase::pressure",
607 "Not implemented for thermo model '{}'", type());
608 }
609
610 //! Density (kg/m^3).
611 //! @return The density of the phase
612 virtual double density() const {
613 return m_dens;
614 }
615
616 //! Molar density (kmol/m^3).
617 //! @return The molar density of the phase
618 virtual double molarDensity() const;
619
620 //! Molar volume (m^3/kmol).
621 //! @return The molar volume of the phase
622 virtual double molarVolume() const;
623
624 //! Set the internally stored density (kg/m^3) of the phase.
625 //! Note the density of a phase is an independent variable.
626 //! @param[in] density_ density (kg/m^3).
627 virtual void setDensity(const double density_);
628
629 //! Set the internally stored pressure (Pa) at constant temperature and
630 //! composition
631 /*!
632 * This method must be reimplemented in derived classes, where it may
633 * involve the solution of a nonlinear equation. Within %Cantera, the
634 * independent variable is either density or pressure. Therefore, this
635 * function may either solve for the density that will yield the desired
636 * input pressure or set an independent variable. The temperature
637 * and composition are held constant during this process.
638 *
639 * @param p input Pressure (Pa)
640 */
641 virtual void setPressure(double p) {
642 throw NotImplementedError("Phase::setPressure",
643 "Not implemented for thermo model '{}'", type());
644 }
645
646 //! Set the internally stored temperature of the phase (K).
647 //! @param temp Temperature in Kelvin
648 virtual void setTemperature(double temp) {
649 if (temp > 0) {
650 m_temp = temp;
651 } else {
652 throw CanteraError("Phase::setTemperature",
653 "temperature must be positive. T = {}", temp);
654 }
655 }
656
657 //! Set the internally stored electron temperature of the phase (K).
658 //! @param etemp Electron temperature in Kelvin
659 virtual void setElectronTemperature(double etemp) {
660 throw NotImplementedError("Phase::setElectronTemperature",
661 "Not implemented for thermo model '{}'", type());
662 }
663
664 //! @}
665
666 //! @name Mean Properties
667 //! @{
668
669 //! Evaluate the mole-fraction-weighted mean of an array Q.
670 //! @f[ \sum_k X_k Q_k. @f]
671 //! Q should contain pure-species molar property values.
672 //! @param[in] Q Array of length m_kk that is to be averaged.
673 //! @return mole-fraction-weighted mean of Q
674 double mean_X(const double* const Q) const;
675
676 //! @copydoc Phase::mean_X(const double* const Q) const
677 double mean_X(const vector<double>& Q) const;
678
679 //! The mean molecular weight. Units: (kg/kmol)
680 double meanMolecularWeight() const {
681 return m_mmw;
682 }
683
684 //! Evaluate @f$ \sum_k X_k \ln X_k @f$.
685 //! @return The indicated sum. Dimensionless.
686 double sum_xlogx() const;
687
688 //! @}
689 //! @name Adding Elements and Species
690 //!
691 //! These methods are used to add new elements or species. These are not
692 //! usually called by user programs.
693 //!
694 //! Since species are checked to insure that they are only composed of
695 //! declared elements, it is necessary to first add all elements before
696 //! adding any species.
697 //! @{
698
699 //! Add an element.
700 //! @param symbol Atomic symbol string.
701 //! @param weight Atomic mass in amu.
702 //! @param atomicNumber Atomic number of the element (unitless)
703 //! @param entropy298 Entropy of the element at 298 K and 1 bar in its
704 //! most stable form. The default is the value ENTROPY298_UNKNOWN,
705 //! which is interpreted as an unknown, and if used will cause
706 //! %Cantera to throw an error.
707 //! @param elem_type Specifies the type of the element constraint
708 //! equation. This defaults to CT_ELEM_TYPE_ABSPOS, that is, an element.
709 //! @return index of the element added
710 size_t addElement(const string& symbol, double weight=-12345.0,
711 int atomicNumber=0, double entropy298=ENTROPY298_UNKNOWN,
712 int elem_type=CT_ELEM_TYPE_ABSPOS);
713
714 //! Add a Species to this Phase. Returns `true` if the species was
715 //! successfully added, or `false` if the species was ignored.
716 //!
717 //! Derived classes which need to size arrays according to the number of
718 //! species should overload this method. The derived class implementation
719 //! should call the base class method, and, if this returns `true`
720 //! (indicating that the species has been added), adjust their array sizes
721 //! accordingly.
722 //!
723 //! @see ignoreUndefinedElements addUndefinedElements throwUndefinedElements
724 virtual bool addSpecies(shared_ptr<Species> spec);
725
726 //! Modify the thermodynamic data associated with a species.
727 /*!
728 * The species name, elemental composition, and type of thermo
729 * parameterization must be unchanged. If there are Kinetics objects that
730 * depend on this phase, Kinetics::invalidateCache() should be called on
731 * those objects after calling this function.
732 */
733 virtual void modifySpecies(size_t k, shared_ptr<Species> spec);
734
735 //! Add a species alias (that is, a user-defined alternative species name).
736 //! Aliases are case-sensitive.
737 //! @param name original species name
738 //! @param alias alternate name
739 void addSpeciesAlias(const string& name, const string& alias);
740
741 //! Lock species list to prevent addition of new species.
742 //! Increments a reference counter used to track whether the Phase is being used by
743 //! a Reactor, Domain1D, or MultiPhase object, which require the number of species
744 //! to remain constant. Should be called in C++ by the object owning the reference.
747 }
748
749 //! Decrement species lock counter.
750 //! Should only be called in C++ by the object owning the reference.
751 void removeSpeciesLock();
752
753 //! Return a vector with isomers names matching a given composition map
754 //! @param compMap Composition of the species.
755 //! @return A vector of species names for matching species.
756 virtual vector<string> findIsomers(const Composition& compMap) const;
757
758 //! Return a vector with isomers names matching a given composition string
759 //! @param comp String containing a composition map
760 //! @return A vector of species names for matching species.
761 virtual vector<string> findIsomers(const string& comp) const;
762
763 //! Return the Species object for the named species. Changes to this object
764 //! do not affect the ThermoPhase object until the #modifySpecies function
765 //! is called.
766 shared_ptr<Species> species(const string& name) const;
767
768 //! Return the Species object for species whose index is *k*. Changes to
769 //! this object do not affect the ThermoPhase object until the
770 //! #modifySpecies function is called.
771 shared_ptr<Species> species(size_t k) const;
772
773 //! Set behavior when adding a species containing undefined elements to just
774 //! skip the species.
776
777 //! Set behavior when adding a species containing undefined elements to add
778 //! those elements to the phase. This is the default behavior.
780
781 //! Set the behavior when adding a species containing undefined elements to
782 //! throw an exception.
784
785 struct UndefElement { enum behavior {
786 error, ignore, add
787 }; };
788
789 //! @} end group adding species and elements
790
791 //! Returns a bool indicating whether the object is ready for use
792 /*!
793 * @returns true if the object is ready for calculation, false otherwise.
794 */
795 virtual bool ready() const;
796
797 //! Return the State Mole Fraction Number
798 int stateMFNumber() const {
799 return m_stateNum;
800 }
801
802 //! Invalidate any cached values which are normally updated only when a
803 //! change in state is detected
804 virtual void invalidateCache();
805
806 //! Returns `true` if case sensitive species names are enforced
807 bool caseSensitiveSpecies() const {
809 }
810
811 //! Set flag that determines whether case sensitive species are enforced
812 //! in look-up operations, for example speciesIndex
813 void setCaseSensitiveSpecies(bool cflag = true) {
815 }
816
817 //! Converts a Composition to a vector with entries for each species
818 //! Species that are not specified are set to zero in the vector
819 /*!
820 * @param[in] comp Composition containing the mixture composition
821 * @return vector with length m_kk
822 */
823 vector<double> getCompositionFromMap(const Composition& comp) const;
824
825 //! Converts a mixture composition from mole fractions to mass fractions
826 //! @param[in] Y mixture composition in mass fractions (length m_kk)
827 //! @param[out] X mixture composition in mole fractions (length m_kk)
828 void massFractionsToMoleFractions(const double* Y, double* X) const;
829
830 //! Converts a mixture composition from mass fractions to mole fractions
831 //! @param[in] X mixture composition in mole fractions (length m_kk)
832 //! @param[out] Y mixture composition in mass fractions (length m_kk)
833 void moleFractionsToMassFractions(const double* X, double* Y) const;
834
835protected:
836 //! Ensure that phase is compressible.
837 //! An error is raised if the state is incompressible
838 //! @param setter name of setter (used for exception handling)
839 void assertCompressible(const string& setter) const {
840 if (!isCompressible()) {
841 throw CanteraError("Phase::assertCompressible",
842 "Setter '{}' is not available. Density is not an "
843 "independent \nvariable for "
844 "'{}' ('{}')", setter, name(), type());
845 }
846 }
847
848 //! Set the internally stored constant density (kg/m^3) of the phase.
849 //! Used for incompressible phases where the density is not an independent
850 //! variable, that is, density does not affect pressure in state calculations.
851 //! @param[in] density_ density (kg/m^3).
852 void assignDensity(const double density_);
853
854 //! Cached for saved calculations within each ThermoPhase.
855 /*!
856 * For more information on how to use this, see examples within the source
857 * code and documentation for this within ValueCache class itself.
858 */
860
861 //! Set the molecular weight of a single species to a given value.
862 //!
863 //! Used by phases where the equation of state is defined for a specific
864 //! value of the molecular weight which may not exactly correspond to the
865 //! value computed from the chemical formula.
866 //! @param k id of the species
867 //! @param mw Molecular Weight (kg kmol-1)
868 void setMolecularWeight(const int k, const double mw);
869
870 //! Apply changes to the state which are needed after the composition
871 //! changes. This function is called after any call to setMassFractions(),
872 //! setMoleFractions(), or similar. For phases which need to execute a
873 //! callback after any change to the composition, it should be done by
874 //! overriding this function rather than overriding all of the composition-
875 //! setting functions. Derived class implementations of compositionChanged()
876 //! should call the parent class method as well.
877 virtual void compositionChanged();
878
879 size_t m_kk = 0; //!< Number of species in the phase.
880
881 //! Dimensionality of the phase. Volumetric phases have dimensionality 3
882 //! and surface phases have dimensionality 2.
883 size_t m_ndim = 3;
884
885 //! Atomic composition of the species. The number of atoms of element i
886 //! in species k is equal to m_speciesComp[k * m_mm + i]
887 //! The length of this vector is equal to m_kk * m_mm
888 vector<double> m_speciesComp;
889
890 vector<double> m_speciesCharge; //!< Vector of species charges. length m_kk.
891
892 map<string, shared_ptr<Species>> m_species; //!< Map of Species objects
893
894 size_t m_nSpeciesLocks = 0; //!< Reference counter preventing species addition
895
896 //! Flag determining behavior when adding species with an undefined element
897 UndefElement::behavior m_undefinedElementBehavior = UndefElement::add;
898
899 //! Flag determining whether case sensitive species names are enforced
901
902 //! Vector of size m_kk, used as a temporary holding area.
903 mutable vector<double> m_workS;
904
905private:
906 //! Find lowercase species name in m_speciesIndices when case sensitive
907 //! species names are not enforced and a user specifies a non-canonical
908 //! species name. Raise exception if lowercase name is not unique.
909 size_t findSpeciesLower(const string& nameStr) const;
910
911 //! Name of the phase.
912 //! Initially, this is the name specified in the YAML input file. It may be changed
913 //! to another value during the course of a calculation.
914 string m_name;
915
916 double m_temp = 0.001; //!< Temperature (K). This is an independent variable
917
918 //! Density (kg m-3). This is an independent variable except in the case
919 //! of incompressible phases, where it has to be changed using the
920 //! assignDensity() method. For compressible substances, the pressure is
921 //! determined from this variable rather than other way round.
922 double m_dens = 0.001;
923
924 double m_mmw = 0.0; //!< mean molecular weight of the mixture (kg kmol-1)
925
926 //! m_ym[k] = mole fraction of species k divided by the mean molecular
927 //! weight of mixture.
928 mutable vector<double> m_ym;
929
930 //! Mass fractions of the species
931 /*!
932 * Note, this vector
933 * Length is m_kk
934 */
935 mutable vector<double> m_y;
936
937 vector<double> m_molwts; //!< species molecular weights (kg kmol-1)
938
939 vector<double> m_rmolwts; //!< inverse of species molecular weights (kmol kg-1)
940
941 //! State Change variable. Whenever the mole fraction vector changes,
942 //! this int is incremented.
943 int m_stateNum = -1;
944
945 //! Vector of the species names
946 vector<string> m_speciesNames;
947
948 //! Map of species names to indices
949 map<string, size_t> m_speciesIndices;
950
951 //! Map of lower-case species names to indices
952 map<string, size_t> m_speciesLower;
953
954 size_t m_mm = 0; //!< Number of elements.
955 vector<double> m_atomicWeights; //!< element atomic weights (kg kmol-1)
956 vector<int> m_atomicNumbers; //!< element atomic numbers
957 vector<string> m_elementNames; //!< element names
958 vector<int> m_elem_type; //!< Vector of element types
959
960 //! Entropy at 298.15 K and 1 bar of stable state pure elements (J kmol-1)
961 vector<double> m_entropy298;
962};
963
964}
965
966#endif
Contains the getElementWeight function and the definitions of element constraint types.
#define CT_ELEM_TYPE_ABSPOS
Normal element constraint consisting of positive coefficients for the formula matrix.
Definition Elements.h:35
#define ENTROPY298_UNKNOWN
Number indicating we don't know the entropy of the element in its most stable state at 298....
Definition Elements.h:85
Base class for exceptions thrown by Cantera classes.
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:96
virtual vector< string > partialStates() const
Return a vector of settable partial property sets within a phase.
Definition Phase.cpp:215
void getCharges(double *charges) const
Copy the vector of species charges into array charges.
Definition Phase.cpp:433
virtual void getConcentrations(double *const c) const
Get the species concentrations (kmol/m^3).
Definition Phase.cpp:510
map< string, size_t > m_speciesLower
Map of lower-case species names to indices.
Definition Phase.h:952
double massFraction(size_t k) const
Return the mass fraction of a single species.
Definition Phase.cpp:483
virtual double molarDensity() const
Molar density (kmol/m^3).
Definition Phase.cpp:604
void assignDensity(const double density_)
Set the internally stored constant density (kg/m^3) of the phase.
Definition Phase.cpp:625
Phase()=default
Default constructor.
virtual bool addSpecies(shared_ptr< Species > spec)
Add a Species to this Phase.
Definition Phase.cpp:729
virtual void setMoleFractions(const double *const x)
Set the mole fractions to the specified values.
Definition Phase.cpp:310
int changeElementType(int m, int elem_type)
Change the element type of the mth constraint Reassigns an element type.
Definition Phase.cpp:96
const vector< double > & atomicWeights() const
Return a read-only reference to the vector of atomic weights.
Definition Phase.cpp:81
virtual vector< string > fullStates() const
Return a vector containing full states defining a phase.
Definition Phase.cpp:197
void assertCompressible(const string &setter) const
Ensure that phase is compressible.
Definition Phase.h:839
vector< double > m_workS
Vector of size m_kk, used as a temporary holding area.
Definition Phase.h:903
void checkSpeciesIndex(size_t k) const
Check that the specified species index is in range.
Definition Phase.cpp:153
void restoreState(const vector< double > &state)
Restore a state saved on a previous call to saveState.
Definition Phase.cpp:286
vector< double > m_speciesComp
Atomic composition of the species.
Definition Phase.h:888
ValueCache m_cache
Cached for saved calculations within each ThermoPhase.
Definition Phase.h:859
virtual void restorePartialState(size_t lenstate, const double *state)
Set the internal thermodynamic state of the phase, excluding composition.
Definition Phase.cpp:243
double m_temp
Temperature (K).
Definition Phase.h:916
size_t m_nSpeciesLocks
Reference counter preventing species addition.
Definition Phase.h:894
vector< string > m_speciesNames
Vector of the species names.
Definition Phase.h:946
size_t nSpecies() const
Returns the number of species in the phase.
Definition Phase.h:232
virtual void setElectronTemperature(double etemp)
Set the internally stored electron temperature of the phase (K).
Definition Phase.h:659
bool m_caseSensitiveSpecies
Flag determining whether case sensitive species names are enforced.
Definition Phase.h:900
vector< string > m_elementNames
element names
Definition Phase.h:957
void ignoreUndefinedElements()
Set behavior when adding a species containing undefined elements to just skip the species.
Definition Phase.cpp:919
UndefElement::behavior m_undefinedElementBehavior
Flag determining behavior when adding species with an undefined element.
Definition Phase.h:897
virtual void setMassFractions_NoNorm(const double *const y)
Set the mass fractions to the specified values without normalizing.
Definition Phase.cpp:376
virtual map< string, size_t > nativeState() const
Return a map of properties defining the native state of a substance.
Definition Phase.cpp:167
double chargeDensity() const
Charge density [C/m^3].
Definition Phase.cpp:635
void addUndefinedElements()
Set behavior when adding a species containing undefined elements to add those elements to the phase.
Definition Phase.cpp:923
virtual void setConcentrationsNoNorm(const double *const conc)
Set the concentrations without ignoring negative concentrations.
Definition Phase.cpp:537
virtual string type() const
String indicating the thermodynamic model implemented.
Definition Phase.h:135
void setNDim(size_t ndim)
Set the number of spatial dimensions (1, 2, or 3).
Definition Phase.h:578
vector< int > m_atomicNumbers
element atomic numbers
Definition Phase.h:956
size_t m_kk
Number of species in the phase.
Definition Phase.h:879
int atomicNumber(size_t m) const
Atomic number of element m.
Definition Phase.cpp:86
virtual void modifySpecies(size_t k, shared_ptr< Species > spec)
Modify the thermodynamic data associated with a species.
Definition Phase.cpp:842
double m_mmw
mean molecular weight of the mixture (kg kmol-1)
Definition Phase.h:924
double elementalMoleFraction(const size_t m) const
Elemental mole fraction of element m.
Definition Phase.cpp:586
size_t m_ndim
Dimensionality of the phase.
Definition Phase.h:883
size_t nDim() const
Returns the number of spatial dimensions (1, 2, or 3)
Definition Phase.h:571
bool caseSensitiveSpecies() const
Returns true if case sensitive species names are enforced.
Definition Phase.h:807
void setCaseSensitiveSpecies(bool cflag=true)
Set flag that determines whether case sensitive species are enforced in look-up operations,...
Definition Phase.h:813
void setState_TD(double t, double rho)
Set the internally stored temperature (K) and density (kg/m^3)
Definition Phase.cpp:398
vector< double > m_rmolwts
inverse of species molecular weights (kmol kg-1)
Definition Phase.h:939
double temperature() const
Temperature (K).
Definition Phase.h:587
virtual void setPressure(double p)
Set the internally stored pressure (Pa) at constant temperature and composition.
Definition Phase.h:641
virtual bool isCompressible() const
Return whether phase represents a compressible substance.
Definition Phase.h:258
double meanMolecularWeight() const
The mean molecular weight. Units: (kg/kmol)
Definition Phase.h:680
virtual double electronTemperature() const
Electron Temperature (K)
Definition Phase.h:593
void moleFractionsToMassFractions(const double *X, double *Y) const
Converts a mixture composition from mass fractions to mole fractions.
Definition Phase.cpp:984
virtual bool hasPhaseTransition() const
Return whether phase represents a substance with phase transitions.
Definition Phase.h:253
virtual void setConcentrations(const double *const conc)
Set the concentrations to the specified values within the phase.
Definition Phase.cpp:515
void removeSpeciesLock()
Decrement species lock counter.
Definition Phase.cpp:875
void addSpeciesLock()
Lock species list to prevent addition of new species.
Definition Phase.h:745
void saveState(vector< double > &state) const
Save the current internal state of the phase.
Definition Phase.cpp:266
Composition getMoleFractionsByName(double threshold=0.0) const
Get the mole fractions by name.
Definition Phase.cpp:438
size_t elementIndex(const string &name) const
Return the index of element named 'name'.
Definition Phase.cpp:55
virtual double concentration(const size_t k) const
Concentration of species k.
Definition Phase.cpp:504
double atomicWeight(size_t m) const
Atomic weight of element m.
Definition Phase.cpp:70
void checkElementArraySize(size_t mm) const
Check that an array size is at least nElements().
Definition Phase.cpp:42
void setMassFractionsByName(const Composition &yMap)
Set the species mass fractions by name.
Definition Phase.cpp:387
int elementType(size_t m) const
Return the element constraint type Possible types include:
Definition Phase.cpp:91
string speciesName(size_t k) const
Name of the species with index k.
Definition Phase.cpp:142
map< string, size_t > m_speciesIndices
Map of species names to indices.
Definition Phase.h:949
virtual void setDensity(const double density_)
Set the internally stored density (kg/m^3) of the phase.
Definition Phase.cpp:614
Composition getMassFractionsByName(double threshold=0.0) const
Get the mass fractions by name.
Definition Phase.cpp:450
virtual size_t stateSize() const
Return size of vector defining internal state of the phase.
Definition Phase.cpp:258
string nativeMode() const
Return string acronym representing the native state of a Phase.
Definition Phase.cpp:184
vector< double > getCompositionFromMap(const Composition &comp) const
Converts a Composition to a vector with entries for each species Species that are not specified are s...
Definition Phase.cpp:955
map< string, shared_ptr< Species > > m_species
Map of Species objects.
Definition Phase.h:892
size_t findSpeciesLower(const string &nameStr) const
Find lowercase species name in m_speciesIndices when case sensitive species names are not enforced an...
Definition Phase.cpp:110
vector< double > m_molwts
species molecular weights (kg kmol-1)
Definition Phase.h:937
virtual vector< string > findIsomers(const Composition &compMap) const
Return a vector with isomers names matching a given composition map.
Definition Phase.cpp:884
const vector< double > & inverseMolecularWeights() const
Return a const reference to the internal vector of molecular weights.
Definition Phase.cpp:428
virtual bool isPure() const
Return whether phase represents a pure (single species) substance.
Definition Phase.h:248
vector< double > m_y
Mass fractions of the species.
Definition Phase.h:935
void getMoleFractions(double *const x) const
Get the species mole fraction vector.
Definition Phase.cpp:462
void setMoleFractionsByName(const Composition &xMap)
Set the species mole fractions by name.
Definition Phase.cpp:351
const double * massFractions() const
Return a const pointer to the mass fraction array.
Definition Phase.h:467
vector< int > m_elem_type
Vector of element types.
Definition Phase.h:958
double sum_xlogx() const
Evaluate .
Definition Phase.cpp:654
string m_name
Name of the phase.
Definition Phase.h:914
const vector< double > & molecularWeights() const
Return a const reference to the internal vector of molecular weights.
Definition Phase.cpp:423
size_t speciesIndex(const string &name) const
Returns the index of a species named 'name' within the Phase object.
Definition Phase.cpp:129
double moleFraction(size_t k) const
Return the mole fraction of a single species.
Definition Phase.cpp:467
double m_dens
Density (kg m-3).
Definition Phase.h:922
const vector< string > & elementNames() const
Return a read-only reference to the vector of element names.
Definition Phase.cpp:65
virtual double density() const
Density (kg/m^3).
Definition Phase.h:612
virtual void compositionChanged()
Apply changes to the state which are needed after the composition changes.
Definition Phase.cpp:951
vector< double > m_atomicWeights
element atomic weights (kg kmol-1)
Definition Phase.h:955
void checkSpeciesArraySize(size_t kk) const
Check that an array size is at least nSpecies().
Definition Phase.cpp:160
void getMolecularWeights(double *weights) const
Copy the vector of molecular weights into array weights.
Definition Phase.cpp:417
double nAtoms(size_t k, size_t m) const
Number of atoms of element m in species k.
Definition Phase.cpp:103
int stateMFNumber() const
Return the State Mole Fraction Number.
Definition Phase.h:798
void addSpeciesAlias(const string &name, const string &alias)
Add a species alias (that is, a user-defined alternative species name).
Definition Phase.cpp:859
virtual void setMolesNoTruncate(const double *const N)
Set the state of the object with moles in [kmol].
Definition Phase.cpp:556
virtual void setTemperature(double temp)
Set the internally stored temperature of the phase (K).
Definition Phase.h:648
size_t nElements() const
Number of elements.
Definition Phase.cpp:30
void setMolecularWeight(const int k, const double mw)
Set the molecular weight of a single species to a given value.
Definition Phase.cpp:941
double mean_X(const double *const Q) const
Evaluate the mole-fraction-weighted mean of an array Q.
Definition Phase.cpp:644
vector< double > m_entropy298
Entropy at 298.15 K and 1 bar of stable state pure elements (J kmol-1)
Definition Phase.h:961
virtual void savePartialState(size_t lenstate, double *state) const
Save the current thermodynamic state of the phase, excluding composition.
Definition Phase.cpp:228
vector< double > m_ym
m_ym[k] = mole fraction of species k divided by the mean molecular weight of mixture.
Definition Phase.h:928
virtual void setMassFractions(const double *const y)
Set the mass fractions to the specified values and normalize them.
Definition Phase.cpp:362
const vector< string > & speciesNames() const
Return a const reference to the vector of species names.
Definition Phase.cpp:148
virtual bool ready() const
Returns a bool indicating whether the object is ready for use.
Definition Phase.cpp:931
double molecularWeight(size_t k) const
Molecular weight of species k.
Definition Phase.cpp:411
double elementalMassFraction(const size_t m) const
Elemental mass fraction of element m.
Definition Phase.cpp:575
shared_ptr< Species > species(const string &name) const
Return the Species object for the named species.
Definition Phase.cpp:902
virtual double molarVolume() const
Molar volume (m^3/kmol).
Definition Phase.cpp:609
virtual void invalidateCache()
Invalidate any cached values which are normally updated only when a change in state is detected.
Definition Phase.cpp:936
void checkElementIndex(size_t m) const
Check that the specified element index is in range.
Definition Phase.cpp:35
void getMassFractions(double *const y) const
Get the species mass fractions.
Definition Phase.cpp:499
virtual size_t partialStateSize() const
Get the size of the partial state vector of the phase.
Definition Phase.h:305
int m_stateNum
State Change variable.
Definition Phase.h:943
void throwUndefinedElements()
Set the behavior when adding a species containing undefined elements to throw an exception.
Definition Phase.cpp:927
void setName(const string &nm)
Sets the string name for the phase.
Definition Phase.cpp:25
size_t m_mm
Number of elements.
Definition Phase.h:954
virtual double pressure() const
Return the thermodynamic pressure (Pa).
Definition Phase.h:605
string elementName(size_t m) const
Name of the element with index m.
Definition Phase.cpp:49
double charge(size_t k) const
Dimensionless electrical charge of a single molecule of species k The charge is normalized by the the...
Definition Phase.h:563
void massFractionsToMoleFractions(const double *Y, double *X) const
Converts a mixture composition from mole fractions to mass fractions.
Definition Phase.cpp:969
double entropyElement298(size_t m) const
Entropy of the element in its standard state at 298 K and 1 bar.
Definition Phase.cpp:75
virtual void setMoleFractions_NoNorm(const double *const x)
Set the mole fractions to the specified values without normalizing.
Definition Phase.cpp:342
vector< double > m_speciesCharge
Vector of species charges.
Definition Phase.h:890
size_t addElement(const string &symbol, double weight=-12345.0, int atomicNumber=0, double entropy298=ENTROPY298_UNKNOWN, int elem_type=CT_ELEM_TYPE_ABSPOS)
Add an element.
Definition Phase.cpp:663
string name() const
Return the name of the phase.
Definition Phase.cpp:20
Storage for cached values.
Definition ValueCache.h:153
Definitions for the classes that are thrown when Cantera experiences an error condition (also contain...
Namespace for the Cantera kernel.
Definition AnyMap.cpp:595
map< string, double > Composition
Map from string names to doubles.
Definition ct_defs.h:177