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