MultiPhase Class Reference
[Equilibrium Solver Capability]

A class for multiphase mixtures. More...

#include <MultiPhase.h>

List of all members.

Public Types

typedef size_t index_t
 Shorthand for an index variable that can't be negative.
typedef ThermoPhase phase_t
 Shorthand for a ThermoPhase.
typedef DenseMatrix array_t
 shorthand for a 2D matrix
typedef std::vector< phase_t * > phase_list
 Shorthand for a vector of pointers to ThermoPhase's.

Public Member Functions

 MultiPhase ()
 Constructor.
virtual ~MultiPhase ()
 Destructor.
void addPhases (phase_list &phases, const vector_fp &phaseMoles)
 Add a vector of phases to the mixture.
void addPhases (MultiPhase &mix)
 Add all phases present in 'mix' to this mixture.
void addPhase (phase_t *p, doublereal moles)
 Add a phase to the mixture.
int nElements () const
 Number of elements.
std::string elementName (int m) const
 Returns the string name of the global element m.
int elementIndex (std::string name) const
 Returns the index of the element with name name.
int nSpecies () const
 Number of species, summed over all phases.
std::string speciesName (const int kGlob) const
 Name of species with global index kGlob.
doublereal nAtoms (const int kGlob, const int mGlob) const
 Returns the Number of atoms of global element mGlob in global species kGlob.
void getMoleFractions (doublereal *const x) const
 Returns the global Species mole fractions.
void init ()
 Process phases and build atomic composition array.
std::string phaseName (const index_t iph) const
 Returns the name of the n'th phase.
int phaseIndex (const std::string &pName) const
 Returns the index, given the phase name.
doublereal phaseMoles (const index_t n) const
 Return the number of moles in phase n.
void setPhaseMoles (const index_t n, const doublereal moles)
 Set the number of moles of phase with index n.
phase_tphase (index_t n)
 Return a ThermoPhase reference to phase n.
doublereal speciesMoles (index_t kGlob) const
 Returns the moles of global species k.
int speciesIndex (index_t k, index_t p) const
 Return the global index of the species belonging to phase number p with local index k within the phase.
int speciesIndex (std::string speciesName, std::string phaseName)
 Return the global index of the species belonging to phase name phaseName with species name speciesName.
doublereal minTemp () const
 Minimum temperature for which all solution phases have valid thermo data.
doublereal maxTemp () const
 Maximum temperature for which all solution phases have valid thermo data.
doublereal charge () const
 Total charge (Coulombs).
doublereal phaseCharge (index_t p) const
 Charge (Coulombs) of phase with index p.
doublereal elementMoles (index_t m) const
 Total moles of global element m, summed over all phases.
void getChemPotentials (doublereal *mu) const
 Returns a vector of Chemical potentials.
void getValidChemPotentials (doublereal not_mu, doublereal *mu, bool standard=false) const
 Returns a vector of Valid chemical potentials.
doublereal temperature () const
 Temperature [K].
doublereal equilibrate (int XY, doublereal err=1.0e-9, int maxsteps=1000, int maxiter=200, int loglevel=-99)
 Set the mixture to a state of chemical equilibrium.
void setTemperature (const doublereal T)
 Set the temperature [K].
void setState_TP (const doublereal T, const doublereal Pres)
 Set the state of the underlying ThermoPhase objects in one call.
void setState_TPMoles (const doublereal T, const doublereal Pres, const doublereal *Moles)
 Set the state of the underlying ThermoPhase objects in one call.
doublereal pressure () const
 Pressure [Pa].
doublereal volume () const
 Volume [m^3].
void setPressure (doublereal P)
 Set the pressure [Pa].
doublereal enthalpy () const
 Enthalpy [J].
doublereal IntEnergy () const
 Enthalpy [J].
doublereal entropy () const
 Entropy [J/K].
doublereal gibbs () const
 Gibbs function [J].
doublereal cp () const
 Heat capacity at constant pressure [J/K].
index_t nPhases () const
 Number of phases.
bool solutionSpecies (index_t kGlob) const
 Return true is species kGlob is a species in a multicomponent solution phase.
int speciesPhaseIndex (const index_t kGlob) const
 Returns the phase index of the Kth "global" species.
doublereal moleFraction (const index_t kGlob) const
 Returns the mole fraction of global species k.
void setPhaseMoleFractions (const index_t n, const doublereal *const x)
 Set the Mole fractions of the nth phase.
void setMolesByName (compositionMap &xMap)
 Set the number numbers of species in the MultiPhase.
void setMolesByName (const std::string &x)
 Set the Moles via a string containing their names.
void getMoles (doublereal *molNum) const
 Return a vector of global species mole numbers.
void setMoles (const doublereal *n)
 Sets all of the global species mole numbers.
void addSpeciesMoles (const int indexS, const doublereal addedMoles)
 Adds moles of a certain species to the mixture.
void getElemAbundances (doublereal *elemAbundances) const
 Retrieves a vector of element abundances.
bool tempOK (index_t p) const
 Return true if the phase p has valid thermo data for the current temperature.
void updateMoleFractions ()
 Update the locally-stored composition within this object to match the current compositions of the phase objects.
void uploadMoleFractionsFromPhases ()
 Update the locally-stored composition within this object to match the current compositions of the phase objects.

Private Member Functions

void updatePhases () const
 Set the states of the phase objects to the locally-stored state within this MultiPhase object.
void calcElemAbundances () const
 Calculate the element abundance vector.

Private Attributes

vector_fp m_moles
 Vector of the number of moles in each phase.
std::vector< phase_t * > m_phase
 Vector of the ThermoPhase Pointers.
array_t m_atoms
 Global Stoichiometric Coefficient array.
vector_fp m_moleFractions
 Locally storred vector of mole fractions of all species comprising the MultiPhase object.
vector_int m_spphase
 Mapping between the global species number and the phase ID.
vector_int m_spstart
 Vector of ints containing of first species index in the global list of species for each phase.
std::vector< std::string > m_enames
 String names of the global elements.
vector_int m_atomicNumber
 Atomic number of each element.
std::vector< std::string > m_snames
 Vector of species names in the problem.
std::map< std::string, int > m_enamemap
 Returns the global element index, given the element string name.
index_t m_np
 Number of phases in the MultiPhase object.
doublereal m_temp
 Current value of the temperature (kelvin).
doublereal m_press
 Current value of the pressure (Pa).
index_t m_nel
 Number of distinct elements in all of the phases.
index_t m_nsp
 Number of distinct species in all of the phases.
bool m_init
 True if the init() routine has been called, and the MultiPhase frozen.
int m_eloc
 Global ID of the element corresponding to the electronic charge.
std::vector< bool > m_temp_OK
 Vector of bools indicating whether temperatures are ok for phases.
doublereal m_Tmin
 Minimum temperature for which thermo parameterizations are valid.
doublereal m_Tmax
 Minimum temperature for which thermo parameterizations are valid.
vector_fp m_elemAbundances
 Vector of element abundances.

Detailed Description

A class for multiphase mixtures.

The mixture can contain any number of phases of any type.

This object is the basic tool used by Cantera for use in Multiphase equilibrium calculations.

It is a container for a set of phases. Each phase has a given number of kmoles. Therefore, MultiPhase may be considered an "extrinsic" thermodynamic object, in contrast to the ThermoPhase object, which is an "intrinsic" thermodynamic object.

MultiPhase may be considered to be "upstream" of the ThermoPhase objects in the sense that setting a property within MultiPhase, such as temperature, pressure, or species mole number, affects the underlying ThermoPhase object, but not the other way around.

All phases have the same temperature and pressure, and a specified number of moles for each phase. The phases do not need to have the same elements. For example, a mixture might consist of a gaseous phase with elements (H, C, O, N), a solid carbon phase containing only element C, etc. A master element set will be constructed for the mixture that is the intersection of the elements of each phase.

Below, reference is made to global species and global elements. These refer to the collective species and elements encompassing all of the phases tracked by the object.

The global element list kept by this object is an intersection of the element lists of all the phases that comprise the MultiPhase.

The global species list kept by this object is a concatenated list of all of the species in all the phases that comprise the MultiPhase. The ordering of species is contiguous with respect to the phase id.

Definition at line 65 of file MultiPhase.h.


Member Typedef Documentation

typedef DenseMatrix array_t

shorthand for a 2D matrix

Definition at line 76 of file MultiPhase.h.

typedef size_t index_t

Shorthand for an index variable that can't be negative.

Definition at line 70 of file MultiPhase.h.

typedef std::vector<phase_t*> phase_list

Shorthand for a vector of pointers to ThermoPhase's.

Definition at line 79 of file MultiPhase.h.

Shorthand for a ThermoPhase.

Definition at line 73 of file MultiPhase.h.


Constructor & Destructor Documentation

MultiPhase (  ) 

Constructor.

The constructor takes no arguments, since phases are added using method addPhase().

Definition at line 25 of file MultiPhase.cpp.

virtual ~MultiPhase (  )  [inline, virtual]

Destructor.

Does nothing. Class MultiPhase does not take "ownership" (i.e. responsibility for destroying) the phase objects.

Definition at line 94 of file MultiPhase.h.


Member Function Documentation

void addPhase ( phase_t p,
doublereal  moles 
)
void addPhases ( MultiPhase mix  ) 

Add all phases present in 'mix' to this mixture.

Parameters:
mix Add all of the phases in another MultiPhase object to the current object.

Definition at line 39 of file MultiPhase.cpp.

References MultiPhase::addPhase(), MultiPhase::m_moles, MultiPhase::m_np, and MultiPhase::m_phase.

void addPhases ( phase_list phases,
const vector_fp &  phaseMoles 
)

Add a vector of phases to the mixture.

See the single addPhases command. This just does a bunch of phases at one time

Parameters:
phases Vector of pointers to phases
phaseMoles Vector of mole numbers in each phase (kmol)

Definition at line 47 of file MultiPhase.cpp.

References MultiPhase::addPhase(), and MultiPhase::init().

void addSpeciesMoles ( const int  indexS,
const doublereal  addedMoles 
)

Adds moles of a certain species to the mixture.

Parameters:
indexS Index of the species in the MultiPhase object
addedMoles Value of the moles that are added to the species.

Definition at line 474 of file MultiPhase.cpp.

References DATA_PTR, MultiPhase::getMoles(), MultiPhase::m_nsp, and MultiPhase::setMoles().

void calcElemAbundances (  )  const [private]
doublereal charge (  )  const

Total charge (Coulombs).

Total charge, summed over all phases.

Definition at line 224 of file MultiPhase.cpp.

References MultiPhase::m_np, and MultiPhase::phaseCharge().

doublereal cp (  )  const

Heat capacity at constant pressure [J/K].

The specific heat at constant pressure and composition (J/K).

Note that this does not account for changes in composition of the mixture with temperature.

Definition at line 367 of file MultiPhase.cpp.

References MultiPhase::m_moles, MultiPhase::m_np, MultiPhase::m_phase, and MultiPhase::updatePhases().

Referenced by MultiPhase::equilibrate().

int elementIndex ( std::string  name  )  const

Returns the index of the element with name name.

Parameters:
name String name of the global element

Definition at line 866 of file MultiPhase.cpp.

References MultiPhase::m_enames, and MultiPhase::m_nel.

doublereal elementMoles ( index_t  m  )  const

Total moles of global element m, summed over all phases.

Total moles of element m, summed over all phases.

Parameters:
m Index of the global element

Definition at line 208 of file MultiPhase.cpp.

References MultiPhase::m_atoms, MultiPhase::m_moleFractions, MultiPhase::m_moles, MultiPhase::m_np, MultiPhase::m_phase, and MultiPhase::speciesIndex().

std::string elementName ( int  m  )  const

Returns the string name of the global element m.

Parameters:
m index of the global element

Definition at line 861 of file MultiPhase.cpp.

References MultiPhase::m_enames.

Referenced by Cantera::BasisOptimize(), and Cantera::ElemRearrange().

doublereal enthalpy (  )  const
doublereal entropy (  )  const
doublereal equilibrate ( int  XY,
doublereal  err = 1.0e-9,
int  maxsteps = 1000,
int  maxiter = 200,
int  loglevel = -99 
)

Set the mixture to a state of chemical equilibrium.

Parameters:
XY Integer flag specifying properties to hold fixed.
err Error tolerance for $\Delta \mu/RT $ for all reactions. Also used as the relative error tolerance for the outer loop.
maxsteps Maximum number of steps to take in solving the fixed TP problem.
maxiter Maximum number of "outer" iterations for problems holding fixed something other than (T,P).
loglevel Level of diagnostic output, written to a file in HTML format.

Definition at line 540 of file MultiPhase.cpp.

References Cantera::addLogEntry(), Cantera::beginLogGroup(), MultiPhase::cp(), Cantera::endLogGroup(), MultiPhase::enthalpy(), MultiPhase::entropy(), Cantera::error(), Cantera::fp2str(), MultiPhase::init(), Cantera::int2str(), MultiPhase::m_init, MultiPhase::m_temp, MultiPhase::m_Tmax, MultiPhase::m_Tmin, MultiPhase::pressure(), MultiPhase::setPressure(), MultiPhase::setTemperature(), MultiPhase::temperature(), Cantera::Undef, and MultiPhase::volume().

void getChemPotentials ( doublereal *  mu  )  const

Returns a vector of Chemical potentials.

Get the chemical potentials of all species in all phases.

Write into array mu the chemical potentials of all species [J/kmol]. The chemical potentials are related to the activities by

$ \mu_k = \mu_k^0(T, P) + RT \ln a_k. $.

Parameters:
mu Chemical potential vector. Length = num global species. Units = J/kmol.

Definition at line 261 of file MultiPhase.cpp.

References MultiPhase::m_np, MultiPhase::m_phase, and MultiPhase::updatePhases().

void getElemAbundances ( doublereal *  elemAbundances  )  const

Retrieves a vector of element abundances.

Parameters:
elemAbundances Vector of element abundances Length = number of elements in the MultiPhase object. Index is the global element index units is in kmol.

Definition at line 498 of file MultiPhase.cpp.

References MultiPhase::calcElemAbundances(), MultiPhase::m_elemAbundances, and MultiPhase::m_nel.

void getMoleFractions ( doublereal *const   x  )  const

Returns the global Species mole fractions.

Write the array of species mole fractions into array x. The mole fractions are normalized to sum to one in each phase.

Parameters:
x vector of mole fractions. Length = number of global species.

Definition at line 884 of file MultiPhase.cpp.

References MultiPhase::m_moleFractions.

void getMoles ( doublereal *  molNum  )  const

Return a vector of global species mole numbers.

Returns a vector of the number of moles of each species in the multiphase object.

Parameters:
molNum Vector of doubles of length nSpecies containing the global mole numbers (kmol).

Definition at line 425 of file MultiPhase.cpp.

References MultiPhase::m_moleFractions, MultiPhase::m_moles, MultiPhase::m_np, MultiPhase::m_phase, and Constituents::nSpecies().

Referenced by MultiPhase::addSpeciesMoles(), and Cantera::BasisOptimize().

void getValidChemPotentials ( doublereal  not_mu,
doublereal *  mu,
bool  standard = false 
) const

Returns a vector of Valid chemical potentials.

Write into array mu the chemical potentials of all species with thermo data valid for the current temperature [J/kmol]. For other species, set the chemical potential to the value not_mu. If standard is set to true, then the values returned are standard chemical potentials.

This method is designed for use in computing chemical equilibrium by Gibbs minimization. For solution phases (more than one species), this does the same thing as getChemPotentials. But for stoichiometric phases, this writes into array mu the user-specified value not_mu instead of the chemical potential if the temperature is outside the range for which the thermo data for the one species in the phase are valid. The need for this arises since many condensed phases have thermo data fit only for the temperature range for which they are stable. For example, in the NASA database, the fits for H2O(s) are only done up to 0 C, the fits for H2O(L) are only done from 0 C to 100 C, etc. Using the polynomial fits outside the range for which the fits were done can result in spurious chemical potentials, and can lead to condensed phases appearing when in fact they should be absent.

By setting not_mu to a large positive value, it is possible to force routines which seek to minimize the Gibbs free energy of the mixture to zero out any phases outside the temperature range for which their thermo data are valid.

Parameters:
not_mu Value of the chemical potential to set species in phases, for which the thermo data is not valid
mu Vector of chemical potentials length = Global species, units = J kmol-1
standard If this method is called with standard set to true, then the composition-independent standard chemical potentials are returned instead of the composition-dependent chemical potentials.

Definition at line 297 of file MultiPhase.cpp.

References MultiPhase::m_np, MultiPhase::m_phase, MultiPhase::nSpecies(), MultiPhase::tempOK(), and MultiPhase::updatePhases().

doublereal gibbs (  )  const

Gibbs function [J].

The Gibbs free energy of the mixture (J).

Definition at line 325 of file MultiPhase.cpp.

References MultiPhase::m_moles, MultiPhase::m_np, MultiPhase::m_phase, and MultiPhase::updatePhases().

void init (  ) 
doublereal IntEnergy (  )  const

Enthalpy [J].

The internal energy of the mixture (J).

Definition at line 345 of file MultiPhase.cpp.

References MultiPhase::m_moles, MultiPhase::m_np, MultiPhase::m_phase, and MultiPhase::updatePhases().

Referenced by vcs_MultiPhaseEquil::equilibrate(), and vcs_MultiPhaseEquil::equilibrate_HP().

doublereal maxTemp (  )  const [inline]

Maximum temperature for which all solution phases have valid thermo data.

Stoichiometric phases are not considered, since they may have thermo data only valid for conditions for which they are stable.

Definition at line 259 of file MultiPhase.h.

References MultiPhase::m_Tmax.

Referenced by vcs_MultiPhaseEquil::equilibrate(), vcs_MultiPhaseEquil::equilibrate_HP(), vcs_MultiPhaseEquil::equilibrate_SP(), vcs_MultiPhaseEquil::equilibrate_TV(), and MultiPhase::updatePhases().

doublereal minTemp (  )  const [inline]

Minimum temperature for which all solution phases have valid thermo data.

Stoichiometric phases are not considered, since they may have thermo data only valid for conditions for which they are stable.

Definition at line 253 of file MultiPhase.h.

References MultiPhase::m_Tmin.

Referenced by vcs_MultiPhaseEquil::equilibrate(), vcs_MultiPhaseEquil::equilibrate_HP(), vcs_MultiPhaseEquil::equilibrate_SP(), vcs_MultiPhaseEquil::equilibrate_TV(), and MultiPhase::updatePhases().

doublereal moleFraction ( const index_t  kGlob  )  const

Returns the mole fraction of global species k.

Parameters:
kGlob Index of the global species.

Definition at line 917 of file MultiPhase.cpp.

References MultiPhase::m_moleFractions.

doublereal nAtoms ( const int  kGlob,
const int  mGlob 
) const

Returns the Number of atoms of global element mGlob in global species kGlob.

Parameters:
kGlob global species index
mGlob global element index
Returns:
returns the number of atoms.

Definition at line 880 of file MultiPhase.cpp.

References MultiPhase::m_atoms.

Referenced by Cantera::BasisOptimize(), and Cantera::ElemRearrange().

int nElements (  )  const [inline]
index_t nPhases (  )  const [inline]
int nSpecies (  )  const [inline]
MultiPhase::phase_t & phase ( index_t  n  ) 

Return a ThermoPhase reference to phase n.

The state of phase n is also updated to match the state stored locally in the mixture object.

Parameters:
n Phase Index
Returns:
Reference to the ThermoPhase object for the phase

Definition at line 192 of file MultiPhase.cpp.

References DATA_PTR, MultiPhase::init(), MultiPhase::m_init, MultiPhase::m_moleFractions, MultiPhase::m_phase, MultiPhase::m_press, MultiPhase::m_spstart, and MultiPhase::m_temp.

Referenced by vcs_MultiPhaseEquil::equilibrate_TP(), Cantera::operator<<(), and vcs_MultiPhaseEquil::reportCSV().

doublereal phaseCharge ( index_t  p  )  const

Charge (Coulombs) of phase with index p.

Net charge of one phase (Coulombs).

Parameters:
p Phase Index

The net charge is computed as

\[ Q_p = N_p \sum_k F z_k X_k \]

where the sum runs only over species in phase p.

Parameters:
p index of the phase for which the charge is desired.

Definition at line 249 of file MultiPhase.cpp.

References MultiPhase::m_moleFractions, MultiPhase::m_moles, MultiPhase::m_phase, and MultiPhase::speciesIndex().

Referenced by MultiPhase::charge().

int phaseIndex ( const std::string &  pName  )  const

Returns the index, given the phase name.

Parameters:
pName Name of the phase
Returns:
returns the index. A value of -1 means the phase isn't in the object.

Definition at line 893 of file MultiPhase.cpp.

References Phase::id(), MultiPhase::m_np, and MultiPhase::m_phase.

Referenced by MultiPhase::speciesIndex().

doublereal phaseMoles ( const index_t  n  )  const

Return the number of moles in phase n.

Parameters:
n Index of the phase.

Definition at line 905 of file MultiPhase.cpp.

References MultiPhase::m_moles.

Referenced by vcs_MultiPhaseEquil::equilibrate_HP(), vcs_MultiPhaseEquil::equilibrate_SP(), and Cantera::operator<<().

std::string phaseName ( const index_t  iph  )  const

Returns the name of the n'th phase.

Parameters:
iph phase Index

Definition at line 888 of file MultiPhase.cpp.

References Phase::id(), and MultiPhase::m_phase.

doublereal pressure (  )  const [inline]
void setMoles ( const doublereal *  n  ) 

Sets all of the global species mole numbers.

Set the species moles to the values in array n.

Sets the number of moles of each species in the multiphase object.

Parameters:
n Vector of doubles of length nSpecies containing the global mole numbers (kmol).

The state of each phase object is also updated to have the specified composition and the mixture temperature and pressure.

Definition at line 445 of file MultiPhase.cpp.

References DATA_PTR, State::getMoleFractions(), MultiPhase::init(), MultiPhase::m_init, MultiPhase::m_moleFractions, MultiPhase::m_moles, MultiPhase::m_np, MultiPhase::m_phase, MultiPhase::m_press, MultiPhase::m_temp, Constituents::nSpecies(), and ThermoPhase::setState_TPX().

Referenced by MultiPhase::addSpeciesMoles(), MultiPhase::setMolesByName(), and MultiPhase::setState_TPMoles().

void setMolesByName ( const std::string &  x  ) 

Set the Moles via a string containing their names.

The string x is in the form of a composition map Species which are not listed by name in the composition map are set to zero.

Parameters:
x string x in the form of a composition map where values are the moles of the species.

Definition at line 406 of file MultiPhase.cpp.

References MultiPhase::nSpecies(), Cantera::parseCompString(), MultiPhase::setMolesByName(), and MultiPhase::speciesName().

void setMolesByName ( compositionMap xMap  ) 

Set the number numbers of species in the MultiPhase.

Parameters:
xMap CompositionMap of the species with nonzero mole numbers units = kmol.

Definition at line 393 of file MultiPhase.cpp.

References DATA_PTR, MultiPhase::nSpecies(), MultiPhase::setMoles(), and MultiPhase::speciesName().

Referenced by MultiPhase::setMolesByName().

void setPhaseMoleFractions ( const index_t  n,
const doublereal *const   x 
)

Set the Mole fractions of the nth phase.

Set the mole fractions of phase n to the values in array x.

This function sets the mole fractions of the nth phase. Note, the mole number of the phase stays constant

Parameters:
n ID of the phase
x Vector of input mole fractions.

Definition at line 380 of file MultiPhase.cpp.

References MultiPhase::m_moleFractions, MultiPhase::m_phase, MultiPhase::m_press, MultiPhase::m_spstart, MultiPhase::m_temp, Constituents::nSpecies(), and ThermoPhase::setState_TPX().

void setPhaseMoles ( const index_t  n,
const doublereal  moles 
)

Set the number of moles of phase with index n.

Parameters:
n Index of the phase
moles Number of moles in the phase (kmol)

Definition at line 909 of file MultiPhase.cpp.

References MultiPhase::m_moles.

Referenced by vcs_MultiPhaseEquil::equilibrate_TP().

void setPressure ( doublereal  P  )  [inline]

Set the pressure [Pa].

Parameters:
P Set the pressure in the MultiPhase object (Pa)

Definition at line 396 of file MultiPhase.h.

References MultiPhase::m_press, and MultiPhase::updatePhases().

Referenced by MultiPhase::equilibrate(), and vcs_MultiPhaseEquil::equilibrate_TV().

void setState_TP ( const doublereal  T,
const doublereal  Pres 
)

Set the state of the underlying ThermoPhase objects in one call.

Parameters:
T Temperature of the system (kelvin)
Pres pressure of the system (pascal) (kmol)

Definition at line 484 of file MultiPhase.cpp.

References MultiPhase::init(), MultiPhase::m_init, MultiPhase::m_press, MultiPhase::m_temp, and MultiPhase::updatePhases().

void setState_TPMoles ( const doublereal  T,
const doublereal  Pres,
const doublereal *  Moles 
)

Set the state of the underlying ThermoPhase objects in one call.

Parameters:
T Temperature of the system (kelvin)
Pres pressure of the system (pascal)
Moles Vector of mole numbers of all the species in all the phases (kmol)

Definition at line 491 of file MultiPhase.cpp.

References MultiPhase::m_press, MultiPhase::m_temp, and MultiPhase::setMoles().

void setTemperature ( const doublereal  T  ) 
bool solutionSpecies ( index_t  kGlob  )  const

Return true is species kGlob is a species in a multicomponent solution phase.

True if species k belongs to a solution phase.

Parameters:
kGlob index of the global species

Definition at line 317 of file MultiPhase.cpp.

References MultiPhase::m_phase, MultiPhase::m_spphase, and MultiPhase::nSpecies().

int speciesIndex ( std::string  speciesName,
std::string  phaseName 
)

Return the global index of the species belonging to phase name phaseName with species name speciesName.

Returns the index of the global species

Parameters:
speciesName Species Name
phaseName Phase Name
Returns:
returns the global index

If the species or phase name is not recognized, this routine throws a CanteraError.

Definition at line 233 of file MultiPhase.cpp.

References MultiPhase::m_phase, MultiPhase::m_spstart, and MultiPhase::phaseIndex().

int speciesIndex ( index_t  k,
index_t  p 
) const [inline]

Return the global index of the species belonging to phase number p with local index k within the phase.

Returns the index of the global species

Parameters:
k local index of the species within the phase
p index of the phase

Definition at line 230 of file MultiPhase.h.

References MultiPhase::m_spstart.

Referenced by MultiPhase::elementMoles(), MultiPhase::phaseCharge(), and vcs_MultiPhaseEquil::reportCSV().

doublereal speciesMoles ( index_t  kGlob  )  const

Returns the moles of global species k.

Moles of species k.

Returns the moles of global species k. units = kmol

Parameters:
kGlob Global species index k

Definition at line 201 of file MultiPhase.cpp.

References MultiPhase::m_moleFractions, MultiPhase::m_moles, and MultiPhase::m_spphase.

std::string speciesName ( const int  kGlob  )  const

Name of species with global index kGlob.

Parameters:
kGlob global species index

Definition at line 876 of file MultiPhase.cpp.

References MultiPhase::m_snames.

Referenced by Cantera::BasisOptimize(), and MultiPhase::setMolesByName().

int speciesPhaseIndex ( const index_t  kGlob  )  const

Returns the phase index of the Kth "global" species.

Parameters:
kGlob Global species index.
Returns:
Returns the index of the owning phase.

Definition at line 913 of file MultiPhase.cpp.

References MultiPhase::m_spphase.

doublereal temperature (  )  const [inline]
bool tempOK ( index_t  p  )  const

Return true if the phase p has valid thermo data for the current temperature.

Parameters:
p Index of the phase.

Definition at line 922 of file MultiPhase.cpp.

References MultiPhase::m_temp_OK.

Referenced by MultiPhase::getValidChemPotentials().

void updateMoleFractions (  ) 

Update the locally-stored composition within this object to match the current compositions of the phase objects.

Update the locally-stored species mole fractions.

Deprecated:
'update' is confusing within this context. Switching to the terminology 'uploadFrom' and 'downloadTo'. uploadFrom means to query the underlying ThermoPhase objects and fill in the resulting information within this object. downloadTo means to take information from this object and put it into the underlying ThermoPhase objects. switch to uploadMoleFractionsFromPhases();

Definition at line 927 of file MultiPhase.cpp.

References MultiPhase::uploadMoleFractionsFromPhases().

Referenced by MultiPhase::init().

void updatePhases (  )  const [private]

Set the states of the phase objects to the locally-stored state within this MultiPhase object.

synchronize the phase objects with the mixture state.

Note that if individual phases have T and P different than that stored locally, the phase T and P will be modified.

This is an download operation in the sense that we are taking upstream object information (MultiPhase object) and applying it to downstrean objects (ThermoPhase object information)

Therefore, the term, "update", is appropriate for a downstream operation.

This method sets each phase to the mixture temperature and pressure, and sets the phase mole fractions based on the mixture mole numbers.

Definition at line 953 of file MultiPhase.cpp.

References DATA_PTR, MultiPhase::m_moleFractions, MultiPhase::m_np, MultiPhase::m_phase, MultiPhase::m_press, MultiPhase::m_temp, MultiPhase::m_temp_OK, MultiPhase::maxTemp(), and MultiPhase::minTemp().

Referenced by MultiPhase::cp(), MultiPhase::enthalpy(), MultiPhase::entropy(), MultiPhase::getChemPotentials(), MultiPhase::getValidChemPotentials(), MultiPhase::gibbs(), MultiPhase::IntEnergy(), MultiPhase::setPressure(), MultiPhase::setState_TP(), and MultiPhase::setTemperature().

void uploadMoleFractionsFromPhases (  ) 

Update the locally-stored composition within this object to match the current compositions of the phase objects.

Update the locally-stored species mole fractions.

Query the underlying ThermoPhase objects for their moel fractions and fill in the mole fraction vector of this current object. Adjust element compositions within this object to match.

This is an upload operation in the sense that we are taking downstream information (ThermoPhase object info) and applying it to an upstream object (MultiPhase object).

Definition at line 931 of file MultiPhase.cpp.

References MultiPhase::calcElemAbundances(), DATA_PTR, State::getMoleFractions(), MultiPhase::m_moleFractions, MultiPhase::m_np, MultiPhase::m_phase, and Constituents::nSpecies().

Referenced by vcs_MultiPhaseEquil::equilibrate_TP(), and MultiPhase::updateMoleFractions().

doublereal volume (  )  const

Volume [m^3].

The total mixture volume [m^3].

Returns the cummulative sum of the volumes of all the phases in the MultiPhase.

Definition at line 531 of file MultiPhase.cpp.

References MultiPhase::m_moles, MultiPhase::m_np, and MultiPhase::m_phase.

Referenced by MultiPhase::equilibrate(), and vcs_MultiPhaseEquil::equilibrate_TV().


Member Data Documentation

vector_int m_atomicNumber [private]

Atomic number of each element.

This is the atomic number of each global element.

Definition at line 625 of file MultiPhase.h.

Referenced by MultiPhase::addPhase(), and MultiPhase::init().

array_t m_atoms [private]

Global Stoichiometric Coefficient array.

This is a two dimensional array m_atoms(m, k). The first index is the global element index. The second index, k, is the global species index. The value is the number of atoms of type m in species k.

Definition at line 592 of file MultiPhase.h.

Referenced by MultiPhase::calcElemAbundances(), MultiPhase::elementMoles(), MultiPhase::init(), and MultiPhase::nAtoms().

vector_fp m_elemAbundances [mutable, private]

Vector of element abundances.

m_elemAbundances[mGlobal] = kmol of element mGlobal summed over all species in all phases.

Definition at line 695 of file MultiPhase.h.

Referenced by MultiPhase::calcElemAbundances(), MultiPhase::getElemAbundances(), and MultiPhase::init().

int m_eloc [private]

Global ID of the element corresponding to the electronic charge.

If there is none, then this is equal to -1

Definition at line 667 of file MultiPhase.h.

Referenced by MultiPhase::addPhase(), and MultiPhase::init().

std::map<std::string, int> m_enamemap [private]

Returns the global element index, given the element string name.

-> used in the construction. However, wonder if it needs to be global.

Definition at line 638 of file MultiPhase.h.

Referenced by MultiPhase::addPhase().

std::vector<std::string> m_enames [private]

String names of the global elements.

This has a length equal to the number of global elements.

Definition at line 619 of file MultiPhase.h.

Referenced by MultiPhase::addPhase(), MultiPhase::elementIndex(), MultiPhase::elementName(), and MultiPhase::init().

bool m_init [private]
vector_fp m_moleFractions [private]
vector_fp m_moles [private]
index_t m_nel [private]
index_t m_np [private]
index_t m_nsp [private]

Number of distinct species in all of the phases.

Definition at line 658 of file MultiPhase.h.

Referenced by MultiPhase::addPhase(), MultiPhase::addSpeciesMoles(), MultiPhase::init(), and MultiPhase::nSpecies().

std::vector<phase_t*> m_phase [private]
doublereal m_press [private]
std::vector<std::string> m_snames [private]

Vector of species names in the problem.

Vector is over all species defined in the object, the global species index.

Definition at line 632 of file MultiPhase.h.

Referenced by MultiPhase::init(), and MultiPhase::speciesName().

vector_int m_spphase [private]

Mapping between the global species number and the phase ID.

m_spphase[kGlobal] = iPhase Length = number of global species

Definition at line 605 of file MultiPhase.h.

Referenced by MultiPhase::init(), MultiPhase::solutionSpecies(), MultiPhase::speciesMoles(), and MultiPhase::speciesPhaseIndex().

vector_int m_spstart [private]

Vector of ints containing of first species index in the global list of species for each phase.

kfirst = m_spstart[ip], kfirst is the index of the first species in the ip'th phase.

Definition at line 613 of file MultiPhase.h.

Referenced by MultiPhase::init(), MultiPhase::phase(), MultiPhase::setPhaseMoleFractions(), and MultiPhase::speciesIndex().

doublereal m_temp [private]
std::vector<bool> m_temp_OK [mutable, private]

Vector of bools indicating whether temperatures are ok for phases.

If the current temperature is outside the range of valid temperatures for the phase thermodynamics, the phase flag is set to false.

Definition at line 674 of file MultiPhase.h.

Referenced by MultiPhase::addPhase(), MultiPhase::tempOK(), and MultiPhase::updatePhases().

doublereal m_Tmax [private]

Minimum temperature for which thermo parameterizations are valid.

Stoichiometric phases are ignored in this determination. units Kelvin

Definition at line 688 of file MultiPhase.h.

Referenced by MultiPhase::addPhase(), MultiPhase::equilibrate(), and MultiPhase::maxTemp().

doublereal m_Tmin [private]

Minimum temperature for which thermo parameterizations are valid.

Stoichiometric phases are ignored in this determination. units Kelvin

Definition at line 681 of file MultiPhase.h.

Referenced by MultiPhase::addPhase(), MultiPhase::equilibrate(), and MultiPhase::minTemp().


The documentation for this class was generated from the following files:
Generated by  doxygen 1.6.3