98 virtual ~Phase() =
default;
127 void setName(
const string& nm);
212 double nAtoms(
size_t k,
size_t m)
const;
219 void getAtoms(
size_t k,
double* atomArray)
const;
321 void saveState(vector<double>& state)
const;
327 virtual void saveState(
size_t lenstate,
double* state)
const;
336 virtual void restoreState(
size_t lenstate,
const double* state);
373 void setState_TRX(
double t,
double dens,
const double* x);
391 void setState_TRY(
double t,
double dens,
const double* y);
682 "Not implemented for thermo model '{}'",
type());
702 virtual void setDensity(
const double density_);
723 "Not implemented for thermo model '{}'",
type());
733 "temperature must be positive. T = {}", temp);
741 "Not implemented for thermo model '{}'",
type());
754 double mean_X(
const double*
const Q)
const;
757 double mean_X(
const vector<double>& Q)
const;
790 size_t addElement(
const string& symbol,
double weight=-12345.0,
804 virtual bool addSpecies(shared_ptr<Species> spec);
813 virtual void modifySpecies(
size_t k, shared_ptr<Species> spec);
829 virtual vector<string>
findIsomers(
const string& comp)
const;
834 shared_ptr<Species>
species(
const string&
name)
const;
839 shared_ptr<Species>
species(
size_t k)
const;
853 struct UndefElement {
enum behavior {
863 virtual bool ready()
const;
910 "Setter '{}' is not available. Density is not an "
911 "independent \nvariable for "
912 "'{}' ('{}')", setter,
name(),
type());
960 map<string, shared_ptr<Species>> m_species;
998 mutable vector<double>
m_y;
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.
#define ENTROPY298_UNKNOWN
Number indicating we don't know the entropy of the element in its most stable state at 298....
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,...
virtual vector< string > partialStates() const
Return a vector of settable partial property sets within a phase.
void getCharges(double *charges) const
Copy the vector of species charges into array charges.
virtual void getConcentrations(double *const c) const
Get the species concentrations (kmol/m^3).
map< string, size_t > m_speciesLower
Map of lower-case species names to indices.
double massFraction(size_t k) const
Return the mass fraction of a single species.
virtual double molarDensity() const
Molar density (kmol/m^3).
void assignDensity(const double density_)
Set the internally stored constant density (kg/m^3) of the phase.
Phase()=default
Default constructor.
virtual bool addSpecies(shared_ptr< Species > spec)
Add a Species to this Phase.
virtual void setMoleFractions(const double *const x)
Set the mole fractions to the specified values.
int changeElementType(int m, int elem_type)
Change the element type of the mth constraint Reassigns an element type.
const vector< double > & atomicWeights() const
Return a read-only reference to the vector of atomic weights.
virtual vector< string > fullStates() const
Return a vector containing full states defining a phase.
void assertCompressible(const string &setter) const
Ensure that phase is compressible.
void checkSpeciesIndex(size_t k) const
Check that the specified species index is in range.
void restoreState(const vector< double > &state)
Restore a state saved on a previous call to saveState.
vector< double > m_speciesComp
Atomic composition of the species.
ValueCache m_cache
Cached for saved calculations within each ThermoPhase.
void getMolecularWeights(vector< double > &weights) const
Copy the vector of molecular weights into vector weights.
const double * moleFractdivMMW() const
Returns a const pointer to the start of the moleFraction/MW array.
double m_temp
Temperature (K). This is an independent variable.
vector< string > m_speciesNames
Vector of the species names.
void setState_RY(double rho, double *y)
Set the density (kg/m^3) and mass fractions.
size_t nSpecies() const
Returns the number of species in the phase.
virtual void setElectronTemperature(double etemp)
Set the internally stored electron temperature of the phase (K).
bool m_caseSensitiveSpecies
Flag determining whether case sensitive species names are enforced.
vector< string > m_elementNames
element names
void ignoreUndefinedElements()
Set behavior when adding a species containing undefined elements to just skip the species.
UndefElement::behavior m_undefinedElementBehavior
Flag determining behavior when adding species with an undefined element.
virtual void setMassFractions_NoNorm(const double *const y)
Set the mass fractions to the specified values without normalizing.
virtual map< string, size_t > nativeState() const
Return a map of properties defining the native state of a substance.
double chargeDensity() const
Charge density [C/m^3].
void addUndefinedElements()
Set behavior when adding a species containing undefined elements to add those elements to the phase.
virtual void setConcentrationsNoNorm(const double *const conc)
Set the concentrations without ignoring negative concentrations.
virtual string type() const
String indicating the thermodynamic model implemented.
void setNDim(size_t ndim)
Set the number of spatial dimensions (1, 2, or 3).
vector< int > m_atomicNumbers
element atomic numbers
size_t m_kk
Number of species in the phase.
int atomicNumber(size_t m) const
Atomic number of element m.
virtual void modifySpecies(size_t k, shared_ptr< Species > spec)
Modify the thermodynamic data associated with a species.
void setState_TRY(double t, double dens, const double *y)
Set the internally stored temperature (K), density, and mass fractions.
double m_mmw
mean molecular weight of the mixture (kg kmol-1)
double elementalMoleFraction(const size_t m) const
Elemental mole fraction of element m.
size_t m_ndim
Dimensionality of the phase.
void setState_TR(double t, double rho)
Set the internally stored temperature (K) and density (kg/m^3)
size_t nDim() const
Returns the number of spatial dimensions (1, 2, or 3)
bool caseSensitiveSpecies() const
Returns true if case sensitive species names are enforced.
void setCaseSensitiveSpecies(bool cflag=true)
Set flag that determines whether case sensitive species are enforced in look-up operations,...
void setState_TD(double t, double rho)
Set the internally stored temperature (K) and density (kg/m^3)
vector< double > m_rmolwts
inverse of species molecular weights (kmol kg-1)
double temperature() const
Temperature (K).
string speciesSPName(int k) const
Returns the expanded species name of a species, including the phase name This is guaranteed to be uni...
virtual void setPressure(double p)
Set the internally stored pressure (Pa) at constant temperature and composition.
virtual bool isCompressible() const
Return whether phase represents a compressible substance.
double meanMolecularWeight() const
The mean molecular weight. Units: (kg/kmol)
virtual double electronTemperature() const
Electron Temperature (K)
void moleFractionsToMassFractions(const double *X, double *Y) const
Converts a mixture composition from mass fractions to mole fractions.
virtual bool hasPhaseTransition() const
Return whether phase represents a substance with phase transitions.
virtual void setConcentrations(const double *const conc)
Set the concentrations to the specified values within the phase.
void setState_RX(double rho, double *x)
Set the density (kg/m^3) and mole fractions.
virtual void setMolarDensity(const double molarDensity)
Set the internally stored molar density (kmol/m^3) of the phase.
void saveState(vector< double > &state) const
Save the current internal state of the phase.
void setState_TX(double t, double *x)
Set the internally stored temperature (K) and mole fractions.
Composition getMoleFractionsByName(double threshold=0.0) const
Get the mole fractions by name.
size_t elementIndex(const string &name) const
Return the index of element named 'name'.
virtual double concentration(const size_t k) const
Concentration of species k.
double atomicWeight(size_t m) const
Atomic weight of element m.
void checkElementArraySize(size_t mm) const
Check that an array size is at least nElements().
void setMassFractionsByName(const Composition &yMap)
Set the species mass fractions by name.
int elementType(size_t m) const
Return the element constraint type Possible types include:
string speciesName(size_t k) const
Name of the species with index k.
map< string, size_t > m_speciesIndices
Map of species names to indices.
virtual void setDensity(const double density_)
Set the internally stored density (kg/m^3) of the phase.
Composition getMassFractionsByName(double threshold=0.0) const
Get the mass fractions by name.
void setState_TRX(double t, double dens, const double *x)
Set the internally stored temperature (K), density, and mole fractions.
virtual size_t stateSize() const
Return size of vector defining internal state of the phase.
string nativeMode() const
Return string acronym representing the native state of a Phase.
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...
size_t findSpeciesLower(const string &nameStr) const
Find lowercase species name in m_speciesIndices when case sensitive species names are not enforced an...
vector< double > m_molwts
species molecular weights (kg kmol-1)
virtual vector< string > findIsomers(const Composition &compMap) const
Return a vector with isomers names matching a given composition map.
const vector< double > & inverseMolecularWeights() const
Return a const reference to the internal vector of molecular weights.
virtual bool isPure() const
Return whether phase represents a pure (single species) substance.
vector< double > m_y
Mass fractions of the species.
void getMoleFractions(double *const x) const
Get the species mole fraction vector.
void setMoleFractionsByName(const Composition &xMap)
Set the species mole fractions by name.
const double * massFractions() const
Return a const pointer to the mass fraction array.
vector< int > m_elem_type
Vector of element types.
void setState_TNX(double t, double n, const double *x)
Set the internally stored temperature (K), molar density (kmol/m^3), and mole fractions.
double sum_xlogx() const
Evaluate .
string m_name
Name of the phase.
const vector< double > & molecularWeights() const
Return a const reference to the internal vector of molecular weights.
size_t speciesIndex(const string &name) const
Returns the index of a species named 'name' within the Phase object.
double moleFraction(size_t k) const
Return the mole fraction of a single species.
double m_dens
Density (kg m-3).
const vector< string > & elementNames() const
Return a read-only reference to the vector of element names.
virtual double density() const
Density (kg/m^3).
virtual void compositionChanged()
Apply changes to the state which are needed after the composition changes.
vector< double > m_atomicWeights
element atomic weights (kg kmol-1)
void checkSpeciesArraySize(size_t kk) const
Check that an array size is at least nSpecies().
double nAtoms(size_t k, size_t m) const
Number of atoms of element m in species k.
int stateMFNumber() const
Return the State Mole Fraction Number.
void addSpeciesAlias(const string &name, const string &alias)
Add a species alias (that is, a user-defined alternative species name).
virtual void setMolesNoTruncate(const double *const N)
Set the state of the object with moles in [kmol].
virtual void setTemperature(double temp)
Set the internally stored temperature of the phase (K).
size_t nElements() const
Number of elements.
void setMolecularWeight(const int k, const double mw)
Set the molecular weight of a single species to a given value.
double mean_X(const double *const Q) const
Evaluate the mole-fraction-weighted mean of an array Q.
vector< double > m_entropy298
Entropy at 298.15 K and 1 bar of stable state pure elements (J kmol-1)
vector< double > m_ym
m_ym[k] = mole fraction of species k divided by the mean molecular weight of mixture.
virtual void setMassFractions(const double *const y)
Set the mass fractions to the specified values and normalize them.
const vector< string > & speciesNames() const
Return a const reference to the vector of species names.
virtual bool ready() const
Returns a bool indicating whether the object is ready for use.
double molecularWeight(size_t k) const
Molecular weight of species k.
double elementalMassFraction(const size_t m) const
Elemental mass fraction of element m.
shared_ptr< Species > species(const string &name) const
Return the Species object for the named species.
virtual double molarVolume() const
Molar volume (m^3/kmol).
virtual void invalidateCache()
Invalidate any cached values which are normally updated only when a change in state is detected.
void checkElementIndex(size_t m) const
Check that the specified element index is in range.
void getMassFractions(double *const y) const
Get the species mass fractions.
void setState_TY(double t, double *y)
Set the internally stored temperature (K) and mass fractions.
void getAtoms(size_t k, double *atomArray) const
Get a vector containing the atomic composition of species k.
int m_stateNum
State Change variable.
void throwUndefinedElements()
Set the behavior when adding a species containing undefined elements to throw an exception.
void setName(const string &nm)
Sets the string name for the phase.
size_t m_mm
Number of elements.
virtual double pressure() const
Return the thermodynamic pressure (Pa).
string elementName(size_t m) const
Name of the element with index m.
double charge(size_t k) const
Dimensionless electrical charge of a single molecule of species k The charge is normalized by the the...
void massFractionsToMoleFractions(const double *Y, double *X) const
Converts a mixture composition from mole fractions to mass fractions.
double entropyElement298(size_t m) const
Entropy of the element in its standard state at 298 K and 1 bar.
virtual void setMoleFractions_NoNorm(const double *const x)
Set the mole fractions to the specified values without normalizing.
vector< double > m_speciesCharge
Vector of species charges. length m_kk.
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.
string name() const
Return the name of the phase.
Storage for cached values.
Definitions for the classes that are thrown when Cantera experiences an error condition (also contain...
Namespace for the Cantera kernel.
map< string, double > Composition
Map from string names to doubles.