research/cantera
1
0
Fork 0
Code Issues Pull requests Projects Releases Packages Wiki Activity Actions

Doxygen update

Browse source 
Finished documenting IdealGasPhase
  Started documenting equilibrium solvers
This commit is contained in:
Harry Moffat 2007-03-13 15:55:14 +00:00
parent 2a6f9647e7
commit a0b5f3241d
4 changed files with 356 additions and 211 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

127
Cantera/src/IdealGasPhase.h
View file

@ -1,7 +1,5 @@
/**
*
* @file IdealGasPhase.h
* `
* @file IdealGasPhase.h
* ThermoPhase object for the ideal gas equation of state.
*/
@ -26,15 +24,18 @@
namespace Cantera {
//!Class %IdealGasPhase represents low-density gases that obey the
//! Class %IdealGasPhase represents low-density gases that obey the
//! ideal gas equation of state.
/*!
*
* %IdealGasPhase derives from class ThermoPhase,
* and overloads the virtual methods defined there with ones that
* use expressions appropriate for ideal gas mixtures.
*
* This class is optimized for speed of execution.
*
* The independent unknowns are density, mass fraction, and temperature.
* the #setPressure() function will calculate the density consistent with
* the current mass fraction vector and temperature and the desired pressure,
* and then set the density in the derived State object.
*
* <HR>
* <H2> Specification of Species Standard %State Properties </H2>
@ -47,9 +48,16 @@ namespace Cantera {
* description of the specification of reference state species thermodynamics functions).
* The reference state,
* where the pressure is fixed at a single pressure,
* is key species property calculation for the Ideal Gas Equation
* is a key species property calculation for the Ideal Gas Equation
* of state.
*
* This class is optimized for speed of execution. All calls to thermodynamic functions
* first call internal routines (aka #enthalpy_RT_ref()) which return references
* the reference state thermodynamics functions. Within these internal reference
* state functions, the function #_updateThermo() is called, that first checks to see
* whether the temperature has changed. If it has, it updates the internal reference
* state thermo functions by calling the SpeciesThermo object.
*
* Functions for the calculation of standard state properties for species
* at arbitray pressure are provided in %IdealGasPhase. However, they
* are all derived from their reference state conterparts.
@ -87,10 +95,10 @@ namespace Cantera {
* \mu^{ref}_k(T) = h^{ref}_k(T) - T S^{ref}_k(T)
* \f]
*
* The standard state internal energy is obtained from the enthalpy function too
* The standard state internal energy is obtained from the enthalpy function also
*
* \f[
* u^o_k(T,P) = h^o_k(T) - R T
* u^o_k(T,P) = h^o_k(T) - R T
* \f]
*
* The molar volume of a species is given by the ideal gas law
@ -124,19 +132,25 @@ namespace Cantera {
* \mu_k(T,P) = \mu^{ref}_k(T, P) + R T \log(\frac{P X_k}{P_{ref}})
* \f]
*
* The partial molar entropy for species k is given by the following relation,
* The partial molar entropy for species <I>k</I> is given by the following relation,
*
* \f[
* \tilde{s}_k(T,P) = s^o_k(T,P) - R \log(X_k) = s^{ref}_k(T) - R \log(\frac{P X_k}{P_{ref}})
* \f]
*
* The partial molar enthalpy for species k is
* The partial molar enthalpy for species <I>k</I> is
*
* \f[
* \tilde{h}_k(T,P) = h^o_k(T,P) = h^{ref}_k(T)
* \f]
*
* The partial molar heat capacity for species k is
* The partial molar Internal Energy for species <I>k</I> is
*
* \f[
* \tilde{u}_k(T,P) = u^o_k(T,P) = u^{ref}_k(T)
* \f]
*
* The partial molar Heat Capacity for species <I>k</I> is
*
* \f[
* \tilde{Cp}_k(T,P) = Cp^o_k(T,P) = Cp^{ref}_k(T)
@ -146,7 +160,7 @@ namespace Cantera {
* <HR>
* <H2> %Application within %Kinetics Managers </H2>
* <HR>
*
* \f$ C^a_k\f$ are defined such that \f$ a_k = C^a_k /
* C^s_k, \f$ where \f$ C^s_k \f$ is a standard concentration
* defined below and \f$ a_k \f$ are activities used in the
@ -246,12 +260,49 @@ namespace Cantera {
* <H2> Instantiation of the Class </H2>
* <HR>
*
*
* The constructor for this phase is located in the default ThermoFactory
* for %Cantera. A new %IdealGasPhase may be created by the following code snippet:
*
* @code
* XML_Node * const xs = xc->findNameID("phase", "silane");
* ThermoPhase *silane_tp = newPhase(*xs);
* IdealGasPhase *silaneGas = dynamic_cast <IdealGasPhase *>(silane_tp);
* @endcode
*
* or by the following constructor:
*
* @code
* XML_Node * const xs = xc->findNameID("phase", "silane");
* IdealGasPhase *silaneGas = new IdealGasPhase(*xs);
* @endcode
* <HR>
* <H2> XML Example </H2>
* <HR>
* An example of an XML Element named phase setting up a IdealGasPhase object named silane
* is given below.
*
* @verbatim
<!-- phase silane -->
<phase dim="3" id="silane">
<elementArray datasrc="elements.xml"> Si H He </elementArray>
<speciesArray datasrc="#species_data">
H2 H HE SIH4 SI SIH SIH2 SIH3 H3SISIH SI2H6
H2SISIH2 SI3H8 SI2 SI3
</speciesArray>
<reactionArray datasrc="#reaction_data"/>
<thermo model="IdealGas"/>
<kinetics model="GasKinetics"/>
<transport model="None"/>
</phase>
@endverbatim
*
* The model attribute "IdealGas" of the thermo XML element identifies the phase as
* being of the type handled by the IdealGasPhase object.
*
* @ingroup thermoprops
*
* @ingroup thermoprops
*/
class IdealGasPhase : public ThermoPhase {
@ -698,19 +749,34 @@ namespace Cantera {
virtual void getStandardVolumes_ref(doublereal *vol) const;
//@}
/// @name New Methods Defined Here -------------------------------------------------
/// @name NonVirtual Internal methods to Return References to Reference State Thermo
//@{
//! Returns a reference to the dimensionless reference state enthalpy vector.
/*!
* This function is part of the layer that checks/recalculates the reference
* state thermo functions.
*/
const array_fp& enthalpy_RT_ref() const {
_updateThermo();
return m_h0_RT;
}
//! Returns a reference to the dimensionless reference state Gibbs free energy vector.
/*!
* This function is part of the layer that checks/recalculates the reference
* state thermo functions.
*/
const array_fp& gibbs_RT_ref() const {
_updateThermo();
return m_g0_RT;
}
//! Returns a reference to the exponent of the dimensionless reference state Gibbs Free energy vector.
/*!
* This function is part of the layer that checks/recalculates the reference
* state thermo functions.
*/
const array_fp& expGibbs_RT_ref() const {
_updateThermo();
int k;
@ -718,26 +784,42 @@ namespace Cantera {
return m_expg0_RT;
}
//! Returns a reference to the dimensionless reference state Entropy vector.
/*!
* This function is part of the layer that checks/recalculates the reference
* state thermo functions.
*/
const array_fp& entropy_R_ref() const {
_updateThermo();
return m_s0_R;
}
//! Returns a reference to the dimensionless reference state Heat Capacity vector.
/*!
* This function is part of the layer that checks/recalculates the reference
* state thermo functions.
*/
const array_fp& cp_R_ref() const {
_updateThermo();
return m_cp0_R;
}
// @}
//@}
/**
* @internal Initialize. This method is provided to allow
//! Initialize the ThermoPhase object after all species have been set up
/*!
* @internal Initialize.
*
* This method is provided to allow
* subclasses to perform any initialization required after all
* species have been added. For example, it might be used to
* resize internal work arrays that must have an entry for
* each species. The base class implementation does nothing,
* and subclasses that do not require initialization do not
* need to overload this method.
* need to overload this method. When importing a CTML phase
* description, this method is called from ThermoPhase::initThermoXML(),
* which is called from importPhase(),
* just prior to returning from function importPhase().
*
* @see importCTML.cpp
*/
@ -829,7 +911,14 @@ namespace Cantera {
private:
//! Update the species reference state thermodynamic functions
/*!
* The polynomials for the standard state functions are only
* reevalulated if the temperature has changed.
*
*/
void _updateThermo() const;
};
}

2
Cantera/src/SurfPhase.h
View file

@ -537,7 +537,7 @@ namespace Cantera {
private:
//! Update the species standard state thermodynamic functions
//! Update the species reference state thermodynamic functions
/*!
* The polynomials for the standard state functions are only
* reevalulated if the temperature has changed.

99
Cantera/src/equil.h
View file

@ -6,7 +6,7 @@
***********************************************************************/
// Copyright 2001 California Institute of Technology
/**
/**
* @file equil.h
* This file contains the definition of some high level general equilibration
* routines and the text for the module \ref equilfunctions.
@ -22,27 +22,84 @@
namespace Cantera {
/*!
* @defgroup equilfunctions Equilibrium Solver Capability
*
* Cantera has several different equilibrium routines.
*/
//-----------------------------------------------------------
// convenience functions
//-----------------------------------------------------------
/**
* Return variable is equal to the number of subroutine attempts
* it took to equilibrate the system.
*/
int equilibrate(thermo_t& s, const char* XY,
int solver = -1, doublereal rtol = 1.0e-9, int maxsteps = 1000,
int maxiter = 100, int loglevel = -99);
doublereal equilibrate(MultiPhase& s, const char* XY,
doublereal tol = 1.0e-9, int maxsteps = 1000, int maxiter = 100,
int loglevel = -99);
/*!
* @defgroup equilfunctions Equilibrium Solver Capability
*
* Cantera has several different equilibrium routines.
*/
//-----------------------------------------------------------
// convenience functions
//-----------------------------------------------------------
//! Equilibrate a ThermoPhase object
/*!
* Set a single-phase chemical solution to chemical equilibrium.
* This is a convenience function that uses one or the other of
* the two chemical equilibrium solvers. The XY parameter indicates what two
* thermodynamic quantities, other than element composition, are to be held
* constant during the equilibration process.
*
* @param s ThermoPhase object that will be equilibrated.
* @param XY String representation of what two properties
* are being held constant
* @param solver ID of the solver to be used to equlibrate the phase.
* If solver = 0, the ChemEquil solver will be used,
* and if solver = 1, the
* MultiPhaseEquil solver will be used (slower than ChemEquil,
* but more stable). If solver < 0 (default, then ChemEquil will
* be tried first, and if it fails MultiPhaseEquil will be tried.
* @param rtol Relative tolerance
* @param maxsteps Maximum number of steps to take to find the solution
* @param maxiter For the MultiPhaseEquil solver only, this is
* the maximum number of outer temperature or pressure iterations
* to take when T and/or P is not held fixed.
* @param loglevel loglevel Controls amount of diagnostic output. loglevel
* = 0 suppresses diagnostics, and increasingly-verbose messages
* are written as loglevel increases. The messages are written to
* a file in HTML format for viewing in a web browser.
* @see HTML_logs
*
* @return
* Return variable is equal to the number of subroutine attempts
* it took to equilibrate the system.
*
*
* @ingroup equilfunctions
* @ingroup equil
*/
int equilibrate(thermo_t& s, const char* XY,
int solver = -1, doublereal rtol = 1.0e-9, int maxsteps = 1000,
int maxiter = 100, int loglevel = -99);
//! Equilibrate a MultiPhase object
/*!
* Equilibrate a MultiPhase object. The XY parameter indicates what two
* thermodynamic quantities, other than element composition, are to be held
* constant during the equilibration process.
*
* This is the top-level driver for multiphase equilibrium. It
* doesn't do much more than call the equilibrate method of class
* MultiPhase, except that it adds some messages to the logfile,
* if loglevel is set > 0.
*
* @param s MultiPhase object that will be equilibrated.
* @param XY String representation of what is being held constant
* @param rtol Relative tolerance
* @param maxsteps Maximum number of steps
* @param maxiter Maximum iterations
* @param loglevel loglevel
*
* @return
* Return variable is equal to the number of subroutine attempts
* it took to equilibrate the system.
*
* @ingroup equilfunctions
* @ingroup equil
*/
doublereal equilibrate(MultiPhase& s, const char* XY,
doublereal rtol = 1.0e-9, int maxsteps = 1000, int maxiter = 100,
int loglevel = -99);
}
#endif

339
Cantera/src/equilibrate.cpp
View file

@ -1,6 +1,5 @@
/**
* @file equilibrate.cpp
*
* Driver routines for the chemical equilibrium solvers.
*
*/
@ -12,177 +11,177 @@
namespace Cantera {
/**
* Set a multiphase mixture to a state of chemical equilibrium.
* This is the top-level driver for multiphase equilibrium. It
* doesn't do much more than call the equilibrate method of class
* MultiPhase, except that it adds some messages to the logfile,
* if loglevel is set > 0.
*
* @ingroup equil
*/
doublereal equilibrate(MultiPhase& s, const char* XY,
doublereal tol, int maxsteps, int maxiter,
int loglevel) {
/*
* Set a multiphase mixture to a state of chemical equilibrium.
* This is the top-level driver for multiphase equilibrium. It
* doesn't do much more than call the equilibrate method of class
* MultiPhase, except that it adds some messages to the logfile,
* if loglevel is set > 0.
*
* @ingroup equil
*/
doublereal equilibrate(MultiPhase& s, const char* XY,
doublereal tol, int maxsteps, int maxiter,
int loglevel) {
beginLogGroup("equilibrate",loglevel);
addLogEntry("multiphase equilibrate function");
beginLogGroup("arguments");
addLogEntry("XY",XY);
addLogEntry("tol",tol);
addLogEntry("maxsteps",maxsteps);
addLogEntry("maxiter",maxiter);
addLogEntry("loglevel",loglevel);
endLogGroup("arguments");
beginLogGroup("equilibrate",loglevel);
addLogEntry("multiphase equilibrate function");
beginLogGroup("arguments");
addLogEntry("XY",XY);
addLogEntry("tol",tol);
addLogEntry("maxsteps",maxsteps);
addLogEntry("maxiter",maxiter);
addLogEntry("loglevel",loglevel);
endLogGroup("arguments");
s.init();
int ixy = _equilflag(XY);
if (ixy == TP || ixy == HP || ixy == SP || ixy == TV) {
try {
double err = s.equilibrate(ixy, tol, maxsteps, maxiter);
addLogEntry("Success. Error",err);
endLogGroup("equilibrate");
return err;
}
catch (CanteraError e) {
addLogEntry("Failure.",lastErrorMessage());
endLogGroup("equilibrate");
throw e;
}
}
else {
addLogEntry("multiphase equilibrium can be done only for TP, HP, SP, or TV");
endLogGroup("equilibrate");
throw CanteraError("equilibrate","unsupported option");
return -1.0;
}
s.init();
int ixy = _equilflag(XY);
if (ixy == TP || ixy == HP || ixy == SP || ixy == TV) {
try {
double err = s.equilibrate(ixy, tol, maxsteps, maxiter);
addLogEntry("Success. Error",err);
endLogGroup("equilibrate");
return err;
}
catch (CanteraError e) {
addLogEntry("Failure.",lastErrorMessage());
endLogGroup("equilibrate");
throw e;
}
}
/// Set a single-phase chemical solution to chemical equilibrium.
/// This is a convenience function that uses one or the other of
/// the two chemical equilibrium solvers.
///
/// @param s The object to set to an equilibrium state
///
/// @param XY An integer specifying the two properties to be held
/// constant.
///
/// @param solver The equilibrium solver to use. If solver = 0,
/// the ChemEquil solver will be used, and if solver = 1, the
/// MultiPhaseEquil solver will be used (slower than ChemEquil,
/// but more stable). If solver < 0 (default, then ChemEquil will
/// be tried first, and if it fails MultiPhaseEquil will be tried.
///
/// @param maxsteps The maximum number of steps to take to find
/// the solution.
///
/// @param maxiter For the MultiPhaseEquil solver only, this is
/// the maximum number of outer temperature or pressure iterations
/// to take when T and/or P is not held fixed.
///
/// @param loglevel Controls amount of diagnostic output. loglevel
/// = 0 suppresses diagnostics, and increasingly-verbose messages
/// are written as loglevel increases. The messages are written to
/// a file in HTML format for viewing in a web browser.
/// @see HTML_logs
///
/// @ingroup equil
int equilibrate(thermo_t& s, const char* XY, int solver,
doublereal rtol, int maxsteps, int maxiter, int loglevel) {
MultiPhase* m = 0;
ChemEquil* e = 0;
bool redo = true;
int retn = -1;
int nAttempts = 0;
int retnSub = 0;
beginLogGroup("equilibrate", loglevel);
addLogEntry("Single-phase equilibrate function");
{
beginLogGroup("arguments");
addLogEntry("phase",s.id());
addLogEntry("XY",XY);
addLogEntry("solver",solver);
addLogEntry("rtol",rtol);
addLogEntry("maxsteps",maxsteps);
addLogEntry("maxiter",maxiter);
addLogEntry("loglevel",loglevel);
endLogGroup("arguments");
}
while (redo) {
if (solver > 0) {
m = new MultiPhase;
try {
m->addPhase(&s, 1.0);
m->init();
nAttempts++;
(void) equilibrate(*m, XY, rtol, maxsteps, maxiter, loglevel);
redo = false;
addLogEntry("MultiPhaseEquil solver succeeded.");
delete m;
retn = nAttempts;
}
catch (CanteraError err) {
addLogEntry("MultiPhaseEquil solver failed.");
delete m;
if (nAttempts < 2) {
addLogEntry("Trying single phase ChemEquil solver.");
solver = -1;
} else {
endLogGroup("equilibrate");
throw err;
}
}
}
else { // solver <= 0
/*
* Call the element potential solver
*/
e = new ChemEquil;
try {
e->options.maxIterations = maxsteps;
e->options.relTolerance = rtol;
nAttempts++;
retnSub = e->equilibrate(s,XY);
if (retnSub < 0) {
addLogEntry("ChemEquil solver failed.");
if (nAttempts < 2) {
addLogEntry("Trying MultiPhaseEquil solver.");
solver = 1;
} else {
throw CanteraError("equilibrate",
"Both equilibrium solvers failed");
}
}
retn = nAttempts;
s.setElementPotentials(e->elementPotentials());
redo = false;
delete e;
addLogEntry("ChemEquil solver succeeded.");
}
catch (CanteraError err) {
delete e;
addLogEntry("ChemEquil solver failed.");
// If ChemEquil fails, try the MultiPhase solver
if (solver < 0) {
addLogEntry("Trying MultiPhaseEquil solver.");
solver = 1;
}
else {
redo = false;
endLogGroup("equilibrate");
throw err;
}
}
}
} // while (redo)
/*
* We are here only for a success
*/
endLogGroup("equilibrate");
return retn;
else {
addLogEntry("multiphase equilibrium can be done only for TP, HP, SP, or TV");
endLogGroup("equilibrate");
throw CanteraError("equilibrate","unsupported option");
return -1.0;
}
}
/*
* Set a single-phase chemical solution to chemical equilibrium.
* This is a convenience function that uses one or the other of
* the two chemical equilibrium solvers.
*
* @param s The object to set to an equilibrium state
*
* @param XY An integer specifying the two properties to be held
* constant.
*
* @param solver The equilibrium solver to use. If solver = 0,
* the ChemEquil solver will be used, and if solver = 1, the
* MultiPhaseEquil solver will be used (slower than ChemEquil,
* but more stable). If solver < 0 (default, then ChemEquil will
* be tried first, and if it fails MultiPhaseEquil will be tried.
*
* @param maxsteps The maximum number of steps to take to find
* the solution.
*
* @param maxiter For the MultiPhaseEquil solver only, this is
* the maximum number of outer temperature or pressure iterations
* to take when T and/or P is not held fixed.
*
* @param loglevel Controls amount of diagnostic output. loglevel
* = 0 suppresses diagnostics, and increasingly-verbose messages
* are written as loglevel increases. The messages are written to
* a file in HTML format for viewing in a web browser.
* @see HTML_logs
*
* @ingroup equil
*/
int equilibrate(thermo_t& s, const char* XY, int solver,
doublereal rtol, int maxsteps, int maxiter, int loglevel) {
MultiPhase* m = 0;
ChemEquil* e = 0;
bool redo = true;
int retn = -1;
int nAttempts = 0;
int retnSub = 0;
beginLogGroup("equilibrate", loglevel);
addLogEntry("Single-phase equilibrate function");
{
beginLogGroup("arguments");
addLogEntry("phase",s.id());
addLogEntry("XY",XY);
addLogEntry("solver",solver);
addLogEntry("rtol",rtol);
addLogEntry("maxsteps",maxsteps);
addLogEntry("maxiter",maxiter);
addLogEntry("loglevel",loglevel);
endLogGroup("arguments");
}
while (redo) {
if (solver > 0) {
m = new MultiPhase;
try {
m->addPhase(&s, 1.0);
m->init();
nAttempts++;
(void) equilibrate(*m, XY, rtol, maxsteps, maxiter, loglevel);
redo = false;
addLogEntry("MultiPhaseEquil solver succeeded.");
delete m;
retn = nAttempts;
}
catch (CanteraError err) {
addLogEntry("MultiPhaseEquil solver failed.");
delete m;
if (nAttempts < 2) {
addLogEntry("Trying single phase ChemEquil solver.");
solver = -1;
} else {
endLogGroup("equilibrate");
throw err;
}
}
}
else { // solver <= 0
/*
* Call the element potential solver
*/
e = new ChemEquil;
try {
e->options.maxIterations = maxsteps;
e->options.relTolerance = rtol;
nAttempts++;
retnSub = e->equilibrate(s,XY);
if (retnSub < 0) {
addLogEntry("ChemEquil solver failed.");
if (nAttempts < 2) {
addLogEntry("Trying MultiPhaseEquil solver.");
solver = 1;
} else {
throw CanteraError("equilibrate",
"Both equilibrium solvers failed");
}
}
retn = nAttempts;
s.setElementPotentials(e->elementPotentials());
redo = false;
delete e;
addLogEntry("ChemEquil solver succeeded.");
}
catch (CanteraError err) {
delete e;
addLogEntry("ChemEquil solver failed.");
// If ChemEquil fails, try the MultiPhase solver
if (solver < 0) {
addLogEntry("Trying MultiPhaseEquil solver.");
solver = 1;
}
else {
redo = false;
endLogGroup("equilibrate");
throw err;
}
}
}
} // while (redo)
/*
* We are here only for a success
*/
endLogGroup("equilibrate");
return retn;
}
}