diff --git a/Cantera/src/NasaPoly1.h b/Cantera/src/NasaPoly1.h index 66efd0744..1e71b17e9 100755 --- a/Cantera/src/NasaPoly1.h +++ b/Cantera/src/NasaPoly1.h @@ -245,6 +245,11 @@ namespace Cantera { #endif } + //! Modify parameters for the standard state + /*! + * @param coeffs Vector of coefficients used to set the + * parameters for the standard state. + */ virtual void modifyParameters(doublereal* coeffs) { m_coeff[0] = coeffs[5]; m_coeff[1] = coeffs[6]; diff --git a/Cantera/src/StoichSubstance.h b/Cantera/src/StoichSubstance.h index 063dfddd4..7b35751af 100644 --- a/Cantera/src/StoichSubstance.h +++ b/Cantera/src/StoichSubstance.h @@ -136,26 +136,30 @@ namespace Cantera { * @{ */ - /** - * Pressure. Units: Pa. - * For an incompressible substance, the density is independent - * of pressure. This method simply returns the stored - * pressure value. - */ - virtual doublereal pressure() const { - return m_press; - } - - /** - * Set the pressure at constant temperature. Units: Pa. - * For an incompressible substance, the density is - * independent of pressure. Therefore, this method only - * stores the specified pressure value. It does not - * modify the density. - */ - virtual void setPressure(doublereal p) { - m_press = p; - } + + //! Report the Pressure. Units: Pa. + /*! + * For an incompressible substance, the density is independent + * of pressure. This method simply returns the storred + * pressure value. + */ + virtual doublereal pressure() const { + return m_press; + } + + + //! Set the pressure at constant temperature. Units: Pa. + /*! + * For an incompressible substance, the density is + * independent of pressure. Therefore, this method only + * stores the specified pressure value. It does not + * modify the density. + * + * @param p Pressure (units - Pa) + */ + virtual void setPressure(doublereal p) { + m_press = p; + } //@} diff --git a/Cantera/src/SurfPhase.h b/Cantera/src/SurfPhase.h index af3d510c0..45d182b38 100644 --- a/Cantera/src/SurfPhase.h +++ b/Cantera/src/SurfPhase.h @@ -35,33 +35,44 @@ namespace Cantera { * * The density of surface sites is given by the variable \f$ n_0 \f$, which has MKS units * of kmol m-2. - + * + * + * Specification of Species Standard State Properties + * + * It is assumed that the reference state thermodynamics may be + * obtained by a pointer to a populated species thermodynamic property + * manager class (see ThermoPhase::m_spthermo). How to relate pressure + * changes to the reference state thermodynamics is resolved at this level. + * + * Pressure is defined as an independent variable in this phase. However, it has + * no effect on any quantities, as the molar concentration is a constant. + * + * Therefore, The standard state internal energy for species k is + * equal to the enthalpy for species k. + * + * \f[ + * u^o_k = h^o_k + * \f] + * + * Also, the standard state chemical potentials, entropy, and heat capacities + * are independent of pressure. The standard state gibbs free energy is obtained + * from the enthalpy and entropy functions. + * + * Specification of Solution Thermodynamic Properties * * The activity of species defined in the phase is given by * \f[ * a_k = \theta_k * \f] * - * The activity concentration,\f$ C^a_k \f$, used by the kinetics manager, is equal to - * the actual concentration, \f$ C^s_k \f$, and is given by the following - * expression. - * \f[ - * C^a_k = C^s_k = \frac{\theta_k n_0}{s_k} - * \f] - * - * The standard concentration for species k is: - * \f[ - * C^0_k = \frac{n_0}{s_k} - * \f] - * - * Pressure is defined as an independent variable in this phase. However, it has - * no effect on any quantities, as the molar concentration is a constant. - * * The chemical potential for species k is equal to * \f[ * \mu_k(T,P) = \mu^o_k(T) + R T \log(\theta_k) * \f] * + * Pressure is defined as an independent variable in this phase. However, it has + * no effect on any quantities, as the molar concentration is a constant. + * * The internal energy for species k is equal to the enthalpy for species k * \f[ * u_k = h_k @@ -74,7 +85,23 @@ namespace Cantera { * s_k(T,P) = s^o_k(T) - R \log(\theta_k) * \f] * - * The constructor for this phase is located in the default ThermoFactory + * Application within %Kinetics Managers + * + * The activity concentration,\f$ C^a_k \f$, used by the kinetics manager, is equal to + * the actual concentration, \f$ C^s_k \f$, and is given by the following + * expression. + * \f[ + * C^a_k = C^s_k = \frac{\theta_k n_0}{s_k} + * \f] + * + * The standard concentration for species k is: + * \f[ + * C^0_k = \frac{n_0}{s_k} + * \f] + * + * Instanteation of the Class + * + * The constructor for this phase is located in the default ThermoFactory * for Cantera. A new SurfPhase may be created by the following code snippet: * * @code @@ -90,10 +117,12 @@ namespace Cantera { * SurfPhase *diamond100TP = new SurfPhase(*xs); * @endcode * + * XML Example + * * An example of an XML Element named phase setting up a SurfPhase object named diamond_100 * is given below. * - * @code + * @verbatim * * H C * c6HH c6H* c6*H c6** c6HM c6HM* c6*M c6B @@ -112,7 +141,7 @@ namespace Cantera { * * * - * @endcode + * @endverbatim * * The model attribute, "Surface", on the thermo element identifies the phase as being * a SurfPhase object. @@ -273,11 +302,11 @@ namespace Cantera { * site density in any convenient form. Internally it is changed * into MKS form. * - * @code + * @verbatim * * 3e-09 * - * @endcode + * @endverbatim */ virtual void setParametersFromXML(const XML_Node& thermoData); @@ -311,12 +340,12 @@ namespace Cantera { * * An example of the XML code block is given below. * - * @code + * @verbatim * * 1200.0 * c6H*:0.1, c6HH:0.9 * - * @endcode + * @endverbatim */ virtual void setStateFromXML(const XML_Node& state); diff --git a/Cantera/src/ThermoPhase.h b/Cantera/src/ThermoPhase.h index a395a2d4d..7e423a17e 100755 --- a/Cantera/src/ThermoPhase.h +++ b/Cantera/src/ThermoPhase.h @@ -306,22 +306,22 @@ namespace Cantera { return err("pressure"); } - //! Set the internally storred pressure (Pa) at constant - //! temperature and composition - /*! - * This method must be reimplemented in derived classes, where it - * may involve the solution of a nonlinear equation. Within %Cantera, - * the independent variable is the density. Therefore, this function - * solves for the density that will yield the desired input pressure. - * The temperature and composition iare held constant during this process. - * - * This base class function will print an error, if not overwritten. - * - * @param p input Pressure (Pa) - */ - virtual void setPressure(doublereal p) { - err("setPressure"); - } + //! Set the internally storred pressure (Pa) at constant + //! temperature and composition + /*! + * This method must be reimplemented in derived classes, where it + * may involve the solution of a nonlinear equation. Within %Cantera, + * the independent variable is the density. Therefore, this function + * solves for the density that will yield the desired input pressure. + * The temperature and composition iare held constant during this process. + * + * This base class function will print an error, if not overwritten. + * + * @param p input Pressure (Pa) + */ + virtual void setPressure(doublereal p) { + err("setPressure"); + } //! Returns the isothermal compressibility. Units: 1/Pa. /*! @@ -467,36 +467,36 @@ namespace Cantera { err("logStandardConc"); return -1.0; } - - /** - * Returns the units of the standard and generalized - * concentrations. Note they have the same units, as their - * ratio is defined to be equal to the activity of the kth - * species in the solution, which is unitless. - * - * This routine is used in print out applications where the - * units are needed. Usually, MKS units are assumed throughout - * the program and in the XML input files. - * - * The base %ThermoPhase class assigns thedefault quantities - * of (kmol/m3) for all species. - * Inherited classes are responsible for overriding the default - * values if necessary. - * - * @param uA Output vector containing the units - * uA[0] = kmol units - default = 1 - * uA[1] = m units - default = -nDim(), the number of spatial - * dimensions in the Phase class. - * uA[2] = kg units - default = 0; - * uA[3] = Pa(pressure) units - default = 0; - * uA[4] = Temperature units - default = 0; - * uA[5] = time units - default = 0 - * @param k species index. Defaults to 0. - * @param sizeUA output int containing the size of the vector. - * Currently, this is equal to 6. - */ - virtual void getUnitsStandardConc(double *uA, int k = 0, - int sizeUA = 6); + + //! Returns the units of the standard and generalized concentrations. + /*! + * Note they have the same units, as their + * ratio is defined to be equal to the activity of the kth + * species in the solution, which is unitless. + * + * This routine is used in print out applications where the + * units are needed. Usually, MKS units are assumed throughout + * the program and in the XML input files. + * + * The base %ThermoPhase class assigns thedefault quantities + * of (kmol/m3) for all species. + * Inherited classes are responsible for overriding the default + * values if necessary. + * + * @param uA Output vector containing the units + * uA[0] = kmol units - default = 1 + * uA[1] = m units - default = -nDim(), the number of spatial + * dimensions in the Phase class. + * uA[2] = kg units - default = 0; + * uA[3] = Pa(pressure) units - default = 0; + * uA[4] = Temperature units - default = 0; + * uA[5] = time units - default = 0 + * @param k species index. Defaults to 0. + * @param sizeUA output int containing the size of the vector. + * Currently, this is equal to 6. + */ + virtual void getUnitsStandardConc(double *uA, int k = 0, + int sizeUA = 6); /** * Get the array of non-dimensional activities at @@ -621,85 +621,85 @@ namespace Cantera { err("getPartialMolarVolumes"); } - //@} - /// @name Properties of the Standard State of the Species in the Solution - //@{ - - //! Get the array of chemical potentials at unit activity for the species - //! at their standard states at the current T and P of the solution. - /*! - * These are the standard state chemical potentials \f$ \mu^0_k(T,P) - * \f$. The values are evaluated at the current - * temperature and pressure of the solution - * - * @param mu Output vector of chemical potentials. - * Length: m_kk. - */ - virtual void getStandardChemPotentials(doublereal* mu) const { - err("getStandardChemPotentials"); + //@} + /// @name Properties of the Standard State of the Species in the Solution + //@{ + + //! Get the array of chemical potentials at unit activity for the species + //! at their standard states at the current T and P of the solution. + /*! + * These are the standard state chemical potentials \f$ \mu^0_k(T,P) + * \f$. The values are evaluated at the current + * temperature and pressure of the solution + * + * @param mu Output vector of chemical potentials. + * Length: m_kk. + */ + virtual void getStandardChemPotentials(doublereal* mu) const { + err("getStandardChemPotentials"); } - //! Get the nondimensional Enthalpy functions for the species - //! at their standard states at the current T and P of the solution. - /*! - * @param hrt Output vector of nondimensional standard state enthalpies. - * Length: m_kk. - */ - virtual void getEnthalpy_RT(doublereal* hrt) const { - err("getEnthalpy_RT"); - } + //! Get the nondimensional Enthalpy functions for the species + //! at their standard states at the current T and P of the solution. + /*! + * @param hrt Output vector of nondimensional standard state enthalpies. + * Length: m_kk. + */ + virtual void getEnthalpy_RT(doublereal* hrt) const { + err("getEnthalpy_RT"); + } - //! Get the array of nondimensional Entropy functions for the - //! standard state species at the current T and P of the solution. - /*! - * @param sr Output vector of nondimensional standard state entropies. - * Length: m_kk. - */ - virtual void getEntropy_R(doublereal* sr) const { - err("getEntropy_R"); - } + //! Get the array of nondimensional Entropy functions for the + //! standard state species at the current T and P of the solution. + /*! + * @param sr Output vector of nondimensional standard state entropies. + * Length: m_kk. + */ + virtual void getEntropy_R(doublereal* sr) const { + err("getEntropy_R"); + } - //! Get the nondimensional Gibbs functions for the species - //! in their standard states at the current T and P of the solution. - /*! - * @param grt Output vector of nondimensional standard state gibbs free energies - * Length: m_kk. - */ - virtual void getGibbs_RT(doublereal* grt) const { - err("getGibbs_RT"); - } + //! Get the nondimensional Gibbs functions for the species + //! in their standard states at the current T and P of the solution. + /*! + * @param grt Output vector of nondimensional standard state gibbs free energies + * Length: m_kk. + */ + virtual void getGibbs_RT(doublereal* grt) const { + err("getGibbs_RT"); + } - //! Get the Gibbs functions for the standard - //! state of the species at the current T and P of the solution - /*! - * Units are Joules/kmol - * @param gpure Output vector of standard state gibbs free energies - * Length: m_kk. - */ - virtual void getPureGibbs(doublereal* gpure) const { - err("getPureGibbs"); - } + //! Get the Gibbs functions for the standard + //! state of the species at the current T and P of the solution + /*! + * Units are Joules/kmol + * @param gpure Output vector of standard state gibbs free energies + * Length: m_kk. + */ + virtual void getPureGibbs(doublereal* gpure) const { + err("getPureGibbs"); + } - //! Returns the vector of nondimensional Internal Energies of the standard - //! state species at the current T and P of the solution - /*! - * @param urt output vector of nondimensional standard state internal energies - * of the species. Length: m_kk. - */ - virtual void getIntEnergy_RT(doublereal *urt) const { - err("getIntEnergy_RT"); - } + //! Returns the vector of nondimensional Internal Energies of the standard + //! state species at the current T and P of the solution + /*! + * @param urt output vector of nondimensional standard state internal energies + * of the species. Length: m_kk. + */ + virtual void getIntEnergy_RT(doublereal *urt) const { + err("getIntEnergy_RT"); + } - //! Get the nondimensional Heat Capacities at constant - //! pressure for the species standard states - //! at the current T and P of the solution - /*! - * @param cpr Output vector of nondimensional standard state heat capacities - * Length: m_kk. - */ - virtual void getCp_R(doublereal* cpr) const { - err("getCp_R"); - } + //! Get the nondimensional Heat Capacities at constant + //! pressure for the species standard states + //! at the current T and P of the solution + /*! + * @param cpr Output vector of nondimensional standard state heat capacities + * Length: m_kk. + */ + virtual void getCp_R(doublereal* cpr) const { + err("getCp_R"); + } //! Get the molar volumes of the species standard states at the current //! T and P of the solution. @@ -757,41 +757,41 @@ namespace Cantera { err("getGibbs_ref"); } - //! Returns the vector of nondimensional - //! entropies of the reference state at the current temperature - //! of the solution and the reference pressure for each species. - /*! - * @param er Output vector containing the nondimensional reference state - * entropies. Length: m_kk. - */ - virtual void getEntropy_R_ref(doublereal *er) const { - err("getEntropy_R_ref"); - } - - //! Returns the vector of nondimensional - //! internal Energies of the reference state at the current temperature - //! of the solution and the reference pressure for each species. - /*! - * @param urt Output vector of nondimensional reference state - * internal energies of the species. - * Length: m_kk - */ - virtual void getIntEnergy_RT_ref(doublereal *urt) const { - err("getIntEnergy_RT_ref"); - } - - //! Returns the vector of nondimensional - //! constant pressure heat capacities of the reference state - //! at the current temperature of the solution - //! and reference pressure for each species. - /*! - * @param cprt Output vector of nondimensional reference state - * heat capacities at constant pressure for the species. - * Length: m_kk - */ - virtual void getCp_R_ref(doublereal *cprt) const { - err("getCp_R_ref()"); - } + //! Returns the vector of nondimensional + //! entropies of the reference state at the current temperature + //! of the solution and the reference pressure for each species. + /*! + * @param er Output vector containing the nondimensional reference state + * entropies. Length: m_kk. + */ + virtual void getEntropy_R_ref(doublereal *er) const { + err("getEntropy_R_ref"); + } + + //! Returns the vector of nondimensional + //! internal Energies of the reference state at the current temperature + //! of the solution and the reference pressure for each species. + /*! + * @param urt Output vector of nondimensional reference state + * internal energies of the species. + * Length: m_kk + */ + virtual void getIntEnergy_RT_ref(doublereal *urt) const { + err("getIntEnergy_RT_ref"); + } + + //! Returns the vector of nondimensional + //! constant pressure heat capacities of the reference state + //! at the current temperature of the solution + //! and reference pressure for each species. + /*! + * @param cprt Output vector of nondimensional reference state + * heat capacities at constant pressure for the species. + * Length: m_kk + */ + virtual void getCp_R_ref(doublereal *cprt) const { + err("getCp_R_ref()"); + } /////////////////////////////////////////////////////// @@ -1283,38 +1283,42 @@ namespace Cantera { void setIndex(int m) { m_index = m; } - /** - * @internal - * Set equation of state parameters. The number and meaning of - * these depends on the subclass. - * @param n number of parameters - * @param c array of \a n coefficients - */ - virtual void setParameters(int n, doublereal* c) {} + //! Set the equation of state parameters + /*! + * @internal + * The number and meaning of these depends on the subclass. + * + * @param n number of parameters + * @param c array of \a n coefficients + */ + virtual void setParameters(int n, doublereal* c) {} - /** - * @internal - * Get equation of state parameters. The number and meaning of - * these depends on the subclass. - * @param n number of parameters - * @param c array of \a n coefficients - */ - virtual void getParameters(int &n, doublereal * const c) {} - /** - * Set equation of state parameter values from XML entries. - * - * This method is called by function importPhase() in - * file importCTML.cpp when processing a phase definition in - * an input file. It should be overloaded in subclasses to set - * any parameters that are specific to that particular phase - * model. Note, this method is called before the phase is - * initialzed with elements and/or species. - * - * @param eosdata An XML_Node object corresponding to - * the "thermo" entry for this phase in the input file. - */ - virtual void setParametersFromXML(const XML_Node& eosdata) {} + //! Get the equation of state parameters in a vector + /*! + * @internal + * The number and meaning of these depends on the subclass. + * + * @param n number of parameters + * @param c array of \a n coefficients + */ + virtual void getParameters(int &n, doublereal * const c) {} + + + //! Set equation of state parameter values from XML entries. + /*! + * + * This method is called by function importPhase() in + * file importCTML.cpp when processing a phase definition in + * an input file. It should be overloaded in subclasses to set + * any parameters that are specific to that particular phase + * model. Note, this method is called before the phase is + * initialzed with elements and/or species. + * + * @param eosdata An XML_Node object corresponding to + * the "thermo" entry for this phase in the input file. + */ + virtual void setParametersFromXML(const XML_Node& eosdata) {} /** * Set the initial state of the phase to the conditions diff --git a/Cantera/src/importCTML.cpp b/Cantera/src/importCTML.cpp index 222b72069..592526a10 100755 --- a/Cantera/src/importCTML.cpp +++ b/Cantera/src/importCTML.cpp @@ -104,7 +104,7 @@ namespace Cantera { * part of the file_ID string. Searches are based on the * ID attribute of the XML element only. * - * @param file_ID This is a concatenation of two strings seperated + * param file_ID This is a concatenation of two strings seperated * by the "#" character. The string before the * pound character is the file name of an xml * file to carry out the search. The string after @@ -113,7 +113,7 @@ namespace Cantera { * The string is interpreted as a file string if * no # character is in the string. * - * @param root If the file string is empty, searches for the + * param root If the file string is empty, searches for the * xml element with matching ID attribute are * carried out from this XML node. */ diff --git a/Cantera/src/importCTML.h b/Cantera/src/importCTML.h index 9194bccb6..3b6b6af68 100755 --- a/Cantera/src/importCTML.h +++ b/Cantera/src/importCTML.h @@ -59,11 +59,11 @@ namespace Cantera { * * will search in the file gri30.xml for an XML element of the following form, where * the XML element name, phase, is an optional hit: - * @code - * + . . . + + * @endverbatim * * It will return a pointer to an xml tree for the XML phase element. * @@ -108,11 +108,11 @@ namespace Cantera { * @endcode * * will search in the file gri30.xml for an XML element of the following form: - * @code + * @verbatim * +#include "importCTML.h" namespace Cantera { - /* - * ---- Constructors ------- - */ + /* + * ---- Constructors ------- + */ - /** - * Default Constructor for the StoichSubstanceSSTP class - */ - StoichSubstanceSSTP::StoichSubstanceSSTP(): - SingleSpeciesTP() - { + /* + * Default Constructor for the StoichSubstanceSSTP class + */ + StoichSubstanceSSTP::StoichSubstanceSSTP(): + SingleSpeciesTP() + { + } + + StoichSubstanceSSTP::StoichSubstanceSSTP(XML_Node& xmlphase, std::string id) { + if (id != "") { + std::string idxml = xmlphase["id"]; + if (id != idxml) { + throw CanteraError("StoichSubstanceSSTP::StoichSubstanceSSTP", + "id's don't match"); + } } - - /** - * Destructor for the routine (virtual) - * - */ - StoichSubstanceSSTP::~StoichSubstanceSSTP() - { + const XML_Node& th = xmlphase.child("thermo"); + std::string model = th["model"]; + if (model != "StoichSubstanceSSTP") { + throw CanteraError("StoichSubstanceSSTP::StoichSubstanceSSTP", + "thermo model attribute must be StoichSubstance"); } + importPhase(xmlphase, this); + } - /* - * ---- Utilities ----- - */ + /* + * Destructor for the routine (virtual) + * + */ + StoichSubstanceSSTP::~StoichSubstanceSSTP() + { + } - /** - * Equation of state flag. Returns the value cStoichSubstance, - * defined in mix_defs.h. - */ - int StoichSubstanceSSTP::eosType() const { - return cStoichSubstance; - } + /* + * ---- Utilities ----- + */ - /* - * ---- Molar Thermodynamic properties of the solution ---- - */ + /* + * Equation of state flag. Returns the value cStoichSubstance, + * defined in mix_defs.h. + */ + int StoichSubstanceSSTP::eosType() const { + return cStoichSubstance; + } - /** - * ----- Mechanical Equation of State ------ - */ + /* + * ---- Molar Thermodynamic properties of the solution ---- + */ - /** - * Pressure. Units: Pa. - * For an incompressible substance, the density is independent - * of pressure. This method simply returns the stored - * pressure value. - */ - doublereal StoichSubstanceSSTP::pressure() const { - return m_press; - } + /** + * ----- Mechanical Equation of State ------ + */ + + /* + * Pressure. Units: Pa. + * For an incompressible substance, the density is independent + * of pressure. This method simply returns the stored + * pressure value. + */ + doublereal StoichSubstanceSSTP::pressure() const { + return m_press; + } - /** - * Set the pressure at constant temperature. Units: Pa. - * For an incompressible substance, the density is - * independent of pressure. Therefore, this method only - * stores the specified pressure value. It does not - * modify the density. - */ - void StoichSubstanceSSTP::setPressure(doublereal p) { - m_press = p; - } + /* + * Set the pressure at constant temperature. Units: Pa. + * For an incompressible substance, the density is + * independent of pressure. Therefore, this method only + * stores the specified pressure value. It does not + * modify the density. + */ + void StoichSubstanceSSTP::setPressure(doublereal p) { + m_press = p; + } - /** - * The isothermal compressibility. Units: 1/Pa. - * The isothermal compressibility is defined as - * \f[ - * \kappa_T = -\frac{1}{v}\left(\frac{\partial v}{\partial P}\right)_T - * \f] - * - * It's equal to zero for this model, since the molar volume - * doesn't change with pressure or temperature. - */ - doublereal StoichSubstanceSSTP::isothermalCompressibility() const { - return 0.0; - } + /* + * The isothermal compressibility. Units: 1/Pa. + * The isothermal compressibility is defined as + * \f[ + * \kappa_T = -\frac{1}{v}\left(\frac{\partial v}{\partial P}\right)_T + * \f] + * + * It's equal to zero for this model, since the molar volume + * doesn't change with pressure or temperature. + */ + doublereal StoichSubstanceSSTP::isothermalCompressibility() const { + return 0.0; + } - /** - * The thermal expansion coefficient. Units: 1/K. - * The thermal expansion coefficient is defined as - * - * \f[ - * \beta = \frac{1}{v}\left(\frac{\partial v}{\partial T}\right)_P - * \f] - * - * It's equal to zero for this model, since the molar volume - * doesn't change with pressure or temperature. - */ - doublereal StoichSubstanceSSTP::thermalExpansionCoeff() const { - return 0.0; - } + /* + * The thermal expansion coefficient. Units: 1/K. + * The thermal expansion coefficient is defined as + * + * \f[ + * \beta = \frac{1}{v}\left(\frac{\partial v}{\partial T}\right)_P + * \f] + * + * It's equal to zero for this model, since the molar volume + * doesn't change with pressure or temperature. + */ + doublereal StoichSubstanceSSTP::thermalExpansionCoeff() const { + return 0.0; + } - /* - * ---- Chemical Potentials and Activities ---- - */ + /* + * ---- Chemical Potentials and Activities ---- + */ - /** - * This method returns the array of generalized - * concentrations. For a stoichiomeetric substance, there is - * only one species, and the generalized concentration is 1.0. - */ - void StoichSubstanceSSTP:: - getActivityConcentrations(doublereal* c) const { - c[0] = 1.0; + /* + * This method returns the array of generalized + * concentrations. For a stoichiomeetric substance, there is + * only one species, and the generalized concentration is 1.0. + */ + void StoichSubstanceSSTP:: + getActivityConcentrations(doublereal* c) const { + c[0] = 1.0; + } + + /* + * The standard concentration. This is defined as the concentration + * by which the generalized concentration is normalized to produce + * the activity. + */ + doublereal StoichSubstanceSSTP::standardConcentration(int k) const { + return 1.0; + } + + /* + * Returns the natural logarithm of the standard + * concentration of the kth species + */ + doublereal StoichSubstanceSSTP::logStandardConc(int k) const { + return 0.0; + } + + /* + * Returns the units of the standard and generalized + * concentrations Note they have the same units, as their + * ratio is defined to be equal to the activity of the kth + * species in the solution, which is unitless. + * + * This routine is used in print out applications where the + * units are needed. Usually, MKS units are assumed throughout + * the program and in the XML input files. + * + * uA[0] = kmol units - default = 1 + * uA[1] = m units - default = -nDim(), the number of spatial + * dimensions in the Phase class. + * uA[2] = kg units - default = 0; + * uA[3] = Pa(pressure) units - default = 0; + * uA[4] = Temperature units - default = 0; + * uA[5] = time units - default = 0 + */ + void StoichSubstanceSSTP:: + getUnitsStandardConc(double *uA, int k, int sizeUA) { + for (int i = 0; i < 6; i++) { + uA[i] = 0; } + } - /** - * The standard concentration. This is defined as the concentration - * by which the generalized concentration is normalized to produce - * the activity. - */ - doublereal StoichSubstanceSSTP::standardConcentration(int k) const { - return 1.0; - } - - /** - * Returns the natural logarithm of the standard - * concentration of the kth species - */ - doublereal StoichSubstanceSSTP::logStandardConc(int k) const { - return 0.0; - } - - /** - * Returns the units of the standard and generalized - * concentrations Note they have the same units, as their - * ratio is defined to be equal to the activity of the kth - * species in the solution, which is unitless. - * - * This routine is used in print out applications where the - * units are needed. Usually, MKS units are assumed throughout - * the program and in the XML input files. - * - * uA[0] = kmol units - default = 1 - * uA[1] = m units - default = -nDim(), the number of spatial - * dimensions in the Phase class. - * uA[2] = kg units - default = 0; - * uA[3] = Pa(pressure) units - default = 0; - * uA[4] = Temperature units - default = 0; - * uA[5] = time units - default = 0 - */ - void StoichSubstanceSSTP:: - getUnitsStandardConc(double *uA, int k, int sizeUA) { - for (int i = 0; i < 6; i++) { - uA[i] = 0; - } - } - - /* - * ---- Partial Molar Properties of the Solution ---- - */ + /* + * ---- Partial Molar Properties of the Solution ---- + */ - /* - * ---- Properties of the Standard State of the Species in the Solution - * ---- - */ + /* + * ---- Properties of the Standard State of the Species in the Solution + * ---- + */ - /** - * Get the array of chemical potentials at unit activity - * \f$ \mu^0_k \f$. - * - * For a stoichiometric substance, there is no activity term in - * the chemical potential expression, and therefore the - * standard chemical potential and the chemical potential - * are both equal to the molar Gibbs function. - */ - void StoichSubstanceSSTP:: - getStandardChemPotentials(doublereal* mu0) const { - getGibbs_RT(mu0); - mu0[0] *= GasConstant * temperature(); - } + /* + * Get the array of chemical potentials at unit activity + * \f$ \mu^0_k \f$. + * + * For a stoichiometric substance, there is no activity term in + * the chemical potential expression, and therefore the + * standard chemical potential and the chemical potential + * are both equal to the molar Gibbs function. + */ + void StoichSubstanceSSTP:: + getStandardChemPotentials(doublereal* mu0) const { + getGibbs_RT(mu0); + mu0[0] *= GasConstant * temperature(); + } - /** - * Get the nondimensional Enthalpy functions for the species - * at their standard states at the current - * T and P of the solution. - * Molar enthalpy. Units: J/kmol. For an incompressible, - * stoichiometric substance, the internal energy is - * independent of pressure, and therefore the molar enthalpy - * is \f[ \hat h(T, P) = \hat u(T) + P \hat v \f], where the - * molar specific volume is constant. - */ - void StoichSubstanceSSTP::getEnthalpy_RT(doublereal* hrt) const { - getEnthalpy_RT_ref(hrt); - double RT = GasConstant * temperature(); - double presCorrect = (m_press - m_p0) / molarDensity(); - hrt[0] += presCorrect / RT; - } + /* + * Get the nondimensional Enthalpy functions for the species + * at their standard states at the current + * T and P of the solution. + * Molar enthalpy. Units: J/kmol. For an incompressible, + * stoichiometric substance, the internal energy is + * independent of pressure, and therefore the molar enthalpy + * is \f[ \hat h(T, P) = \hat u(T) + P \hat v \f], where the + * molar specific volume is constant. + */ + void StoichSubstanceSSTP::getEnthalpy_RT(doublereal* hrt) const { + getEnthalpy_RT_ref(hrt); + double RT = GasConstant * temperature(); + double presCorrect = (m_press - m_p0) / molarDensity(); + hrt[0] += presCorrect / RT; + } - /** - * Get the array of nondimensional Entropy functions for the - * standard state species - * at the current T and P of the solution. - */ - void StoichSubstanceSSTP::getEntropy_R(doublereal* sr) const { - getEntropy_R_ref(sr); - } + /* + * Get the array of nondimensional Entropy functions for the + * standard state species + * at the current T and P of the solution. + */ + void StoichSubstanceSSTP::getEntropy_R(doublereal* sr) const { + getEntropy_R_ref(sr); + } - /** - * Get the nondimensional Gibbs functions for the species - * at their standard states of solution at the current T and P - * of the solution - */ - void StoichSubstanceSSTP::getGibbs_RT(doublereal* grt) const { - getEnthalpy_RT(grt); - grt[0] -= m_s0_R[0]; - } + /* + * Get the nondimensional Gibbs functions for the species + * at their standard states of solution at the current T and P + * of the solution + */ + void StoichSubstanceSSTP::getGibbs_RT(doublereal* grt) const { + getEnthalpy_RT(grt); + grt[0] -= m_s0_R[0]; + } - /** - * Get the nondimensional Gibbs functions for the standard - * state of the species at the current T and P. - */ - void StoichSubstanceSSTP::getCp_R(doublereal* cpr) const { - _updateThermo(); - cpr[0] = m_cp0_R[0]; - } + /* + * Get the nondimensional Gibbs functions for the standard + * state of the species at the current T and P. + */ + void StoichSubstanceSSTP::getCp_R(doublereal* cpr) const { + _updateThermo(); + cpr[0] = m_cp0_R[0]; + } - /** - * Molar internal energy (J/kmol). - * For an incompressible, - * stoichiometric substance, the molar internal energy is - * independent of pressure. Since the thermodynamic properties - * are specified by giving the standard-state enthalpy, the - * term \f$ P_0 \hat v\f$ is subtracted from the specified molar - * enthalpy to compute the molar internal energy. - */ - void StoichSubstanceSSTP::getIntEnergy_RT(doublereal* urt) const { - _updateThermo(); - double RT = GasConstant * temperature(); - double PV = m_press / molarDensity(); - urt[0] = m_h0_RT[0] - PV / RT; - } + /* + * Molar internal energy (J/kmol). + * For an incompressible, + * stoichiometric substance, the molar internal energy is + * independent of pressure. Since the thermodynamic properties + * are specified by giving the standard-state enthalpy, the + * term \f$ P_0 \hat v\f$ is subtracted from the specified molar + * enthalpy to compute the molar internal energy. + */ + void StoichSubstanceSSTP::getIntEnergy_RT(doublereal* urt) const { + _updateThermo(); + double RT = GasConstant * temperature(); + double PV = m_p0 / molarDensity(); + urt[0] = m_h0_RT[0] - PV / RT; + } - /* - * ---- Thermodynamic Values for the Species Reference States ---- - */ - /** - * Molar internal energy or the reference state at the current - * temperature, T (J/kmol). - * For an incompressible, - * stoichiometric substance, the molar internal energy is - * independent of pressure. Since the thermodynamic properties - * are specified by giving the standard-state enthalpy, the - * term \f$ P_0 \hat v\f$ is subtracted from the specified molar - * enthalpy to compute the molar internal energy. - * - * Note, this is equal to the standard state internal energy - * evaluated at the reference pressure. - */ - void StoichSubstanceSSTP::getIntEnergy_RT_ref(doublereal* urt) const { - _updateThermo(); - double RT = GasConstant * temperature(); - double PV = m_p0 / molarDensity(); - urt[0] = m_h0_RT[0] - PV / RT; - } + /* + * ---- Thermodynamic Values for the Species Reference States ---- + */ + /* + * Molar internal energy or the reference state at the current + * temperature, T (J/kmol). + * For an incompressible, + * stoichiometric substance, the molar internal energy is + * independent of pressure. Since the thermodynamic properties + * are specified by giving the standard-state enthalpy, the + * term \f$ P_0 \hat v\f$ is subtracted from the specified molar + * enthalpy to compute the molar internal energy. + * + * Note, this is equal to the standard state internal energy + * evaluated at the reference pressure. + */ + void StoichSubstanceSSTP::getIntEnergy_RT_ref(doublereal* urt) const { + _updateThermo(); + double RT = GasConstant * temperature(); + double PV = m_p0 / molarDensity(); + urt[0] = m_h0_RT[0] - PV / RT; + } + + /* + * ---- Saturation Properties + */ + + + /* + * ---- Initialization and Internal functions + */ + + /** + * @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. When importing a CTML phase + * description, this method is called just prior to returning + * from function importPhase. + * + * @see importCTML.cpp + */ + void StoichSubstanceSSTP::initThermo() { /* - * ---- Critical State Properties + * Make sure there is one and only one species in this phase. */ - /// Critical temperature (K). - doublereal StoichSubstanceSSTP::critTemperature() const { - return -1.0; + m_kk = nSpecies(); + if (m_kk != 1) { + throw CanteraError("initThermo", + "stoichiometric substances may only contain one species."); } - - /// Critical pressure (Pa). - doublereal StoichSubstanceSSTP::critPressure() const { - return -1.0; - } - - /// Critical density (kg/m3). - doublereal StoichSubstanceSSTP::critDensity() const { - return -1.0; - } - + doublereal tmin = m_spthermo->minTemp(); + doublereal tmax = m_spthermo->maxTemp(); + if (tmin > 0.0) m_tmin = tmin; + if (tmax > 0.0) m_tmax = tmax; /* - * ---- Saturation Properties + * Store the reference pressure in the variables for the class. */ - - doublereal StoichSubstanceSSTP::satTemperature(doublereal p) const { - return (-1.0); - } - doublereal StoichSubstanceSSTP::satPressure(doublereal t) const { - return 0.0; - } - doublereal StoichSubstanceSSTP::vaporFraction() const { - return 0.0; - } - void StoichSubstanceSSTP::setState_Tsat(doublereal t, doublereal x) { - setTemperature(t); - } - void StoichSubstanceSSTP::setState_Psat(doublereal p, doublereal x) { - setPressure(p); - } + m_p0 = refPressure(); /* - * ---- Initialization and Internal functions + * Resize temporary arrays. */ - - /** - * @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. When importing a CTML phase - * description, this method is called just prior to returning - * from function importPhase. - * - * @see importCTML.cpp + int leng = 1; + m_h0_RT.resize(leng); + m_cp0_R.resize(leng); + m_s0_R.resize(leng); + /* + * Call the base class thermo initializer */ - void StoichSubstanceSSTP::initThermo() { - /* - * Make sure there is one and only one species in this phase. - */ - m_kk = nSpecies(); - if (m_kk != 1) { - throw CanteraError("initThermo", - "stoichiometric substances may only contain one species."); - } - doublereal tmin = m_spthermo->minTemp(); - doublereal tmax = m_spthermo->maxTemp(); - if (tmin > 0.0) m_tmin = tmin; - if (tmax > 0.0) m_tmax = tmax; - /* - * Store the reference pressure in the variables for the class. - */ - m_p0 = refPressure(); + SingleSpeciesTP::initThermo(); + } - /* - * Resize temporary arrays. - */ - int leng = 1; - m_h0_RT.resize(leng); - m_cp0_R.resize(leng); - m_s0_R.resize(leng); - /* - * Call the base class thermo initializer - */ - SingleSpeciesTP::initThermo(); - } + /** + * setParameters: + * + * Generic routine that is used to set the parameters used + * by this model. + * C[0] = density of phase [ kg/m3 ] + */ + void StoichSubstanceSSTP::setParameters(int n, double * c) { + double rho = c[0]; + setDensity(rho); + } - /** - * setParameters: - * - * Generic routine that is used to set the parameters used - * by this model. - * C[0] = density of phase [ kg/m3 ] - */ - void StoichSubstanceSSTP::setParameters(int n, double * c) { - double rho = c[0]; - setDensity(rho); - } + /** + * getParameters: + * + * Generic routine that is used to get the parameters used + * by this model. + * n = 1 + * C[0] = density of phase [ kg/m3 ] + */ + void StoichSubstanceSSTP::getParameters(int &n, double * const c) { + double rho = density(); + n = 1; + c[0] = rho; + } - /** - * getParameters: - * - * Generic routine that is used to get the parameters used - * by this model. - * n = 1 - * C[0] = density of phase [ kg/m3 ] - */ - void StoichSubstanceSSTP::getParameters(int &n, double * const c) { - double rho = density(); - n = 1; - c[0] = rho; - } - - /** - * Reads an xml data block for the parameters needed by this - * routine. eosdata is a reference to the xml thermo block, and looks - * like this: - * - * - * - * 3.52 - * - * - */ - void StoichSubstanceSSTP::setParametersFromXML(const XML_Node& eosdata) { - eosdata._require("model","StoichSubstanceSSTP"); - doublereal rho = getFloat(eosdata, "density", "-"); - setDensity(rho); - } + /* + * Reads an xml data block for the parameters needed by this + * routine. eosdata is a reference to the xml thermo block, and looks + * like this: + * + * + * + * 3.52 + * + * + */ + void StoichSubstanceSSTP::setParametersFromXML(const XML_Node& eosdata) { + eosdata._require("model","StoichSubstanceSSTP"); + doublereal rho = getFloat(eosdata, "density", "-"); + setDensity(rho); + } } diff --git a/Cantera/src/thermo/StoichSubstanceSSTP.h b/Cantera/src/thermo/StoichSubstanceSSTP.h index 60850ed18..905d7e9c7 100644 --- a/Cantera/src/thermo/StoichSubstanceSSTP.h +++ b/Cantera/src/thermo/StoichSubstanceSSTP.h @@ -26,278 +26,476 @@ namespace Cantera { + /** + * @ingroup thermoprops + * + * Class %StoichSubstanceSSTP represents a stoichiometric (fixed composition) + * incompressible substance. + * This class internally changes the independent degree of freedom from + * density to pressure. This is necessary because the phase is incompressible. + * It uses a constant volume approximation. + * + * + * Specification of Species Standard %State Properties + * + * This class inherits from SingleSpeciesTP. + * It is assumed that the reference state thermodynamics may be + * obtained by a pointer to a populated species thermodynamic property + * manager class (see ThermoPhase::m_spthermo). How to relate pressure + * changes to the reference state thermodynamics is resolved at this level. + * + * For an incompressible, + * stoichiometric substance, the molar internal energy is + * independent of pressure. Since the thermodynamic properties + * are specified by giving the standard-state enthalpy, the + * term \f$ P_0 \hat v\f$ is subtracted from the specified molar + * enthalpy to compute the molar internal energy. The entropy is + * assumed to be independent of the pressure. + * + * The enthalpy function is given by the following relation. + * + * \f[ + * h^o_k(T,P) = h^{ref}_k(T) + \tilde v \left( P - P_{ref} \right) + * \f] + * + * For an incompressible, + * stoichiometric substance, the molar internal energy is + * independent of pressure. Since the thermodynamic properties + * are specified by giving the standard-state enthalpy, the + * term \f$ P_{ref} \tilde v\f$ is subtracted from the specified reference molar + * enthalpy to compute the molar internal energy. + * + * \f[ + * u^o_k(T,P) = h^{ref}_k(T) - P_{ref} \tilde v + * \f] + * + * The standard state heat capacity and entropy are independent + * of pressure. The standard state gibbs free energy is obtained + * from the enthalpy and entropy functions. + * + * + * Specification of Solution Thermodynamic Properties + * + * All solution properties are obtained from the standard state + * species functions, since there is only one species in the phase. + * + * Application within %Kinetics Managers + * + * The standard concentration is equal to 1.0. This means that the + * kinetics operator works on an (activities basis). Since this + * is a stoichiometric substance, this means that the concentration + * of this phase drops out of kinetics expressions. + * + * An example of a reaction using this is a sticking coefficient + * reaction of a substance in an ideal gas phase on a surface with a bulk phase + * species in this phase. In this case, the rate of progress for this + * reaction, \f$ R_s \f$, may be expressed via the following equation: + * \f[ + * R_s = k_s C_{gas} + * \f] + * where the units for \f$ R_s \f$ are kmol m-2 s-1. \f$ C_{gas} \f$ has units + * of kmol m-3. Therefore, the kinetic rate constant, \f$ k_s \f$, has + * units of m s-1. Nowhere does the concentration of the bulk phase + * appear in the rate constant expression, since it's a stoichiometric + * phase and the activity is always equal to 1.0. + * + * Instanteation of the Class + * + * The constructor for this phase is NOT located in the default ThermoFactory + * for %Cantera. However, a new %StoichSubstanceSSTP may be created by + * the following code snippets: + * + * @code + * sprintf(file_ID,"%s#NaCl(S)", iFile); + * XML_Node *xm = get_XML_NameID("phase", file_ID, 0); + * StoichSubstanceSSTP *solid = new StoichSubstanceSSTP(*xm); + * @endcode + * + * or by the following call to importPhase(): + * + * @code + * sprintf(file_ID,"%s#NaCl(S)", iFile); + * XML_Node *xm = get_XML_NameID("phase", file_ID, 0); + * StoichSubstanceSSTP solid; + * importPhase(*xm, &solid); + * @endcode + * + * XML Example + * + * The phase model name for this is called StoichSubstance. It must be supplied + * as the model attribute of the thermo XML element entry. + * Within the phase XML block, + * the density of the phase must be specified. An example of an XML file + * this phase is given below. + * + * @verbatim + + + + Na Cl + + NaCl(S) + + 2.165 + + + + + + + + + + Na:1 Cl:1 + + + + 50.72389, 6.672267, -2.517167, + 10.15934, -0.200675, -427.2115, + 130.3973 + + + + 2.165 + + @endverbatim + * + * The model attribute, "StoichSubstanceSSTP", on the thermo element identifies the phase as being + * a StoichSubstanceSSTP object. + * + */ + class StoichSubstanceSSTP : public SingleSpeciesTP { + + public: /** - * @ingroup thermoprops + * Default Constructor for the StoichSubstanceSSTP class + */ + StoichSubstanceSSTP(); + + //! Constructor. + /*! + * @param phaseRef XML node pointing to a StoichSubstanceSSTP description + * @param id Id of the phase. + */ + StoichSubstanceSSTP(XML_Node& phaseRef, std::string id = ""); + + + /** + * Destructor for the routine (virtual) + * + */ + virtual ~StoichSubstanceSSTP(); + + /** + * + * @name Utilities + * @{ + */ + + /** + * Equation of state flag. * - * Class StoichSubstance represents a stoichiometric (fixed composition) - * incompressible substance. + * Returns the value cStoichSubstance, defined in mix_defs.h. + */ + virtual int eosType() const; + + /** + * @} + * @name Molar Thermodynamic Properties of the Solution + * @{ + */ + + /** + * @} + * @name Mechanical Equation of State + * @{ + */ + + + //! Report the Pressure. Units: Pa. + /*! + * For an incompressible substance, the density is independent + * of pressure. This method simply returns the storred + * pressure value. + */ + virtual doublereal pressure() const; + + //! Set the pressure at constant temperature. Units: Pa. + /*! + * For an incompressible substance, the density is + * independent of pressure. Therefore, this method only + * stores the specified pressure value. It does not + * modify the density. + * + * @param p Pressure (units - Pa) + */ + virtual void setPressure(doublereal p); + + //! Returns the isothermal compressibility. Units: 1/Pa. + /*! + * The isothermal compressibility is defined as + * \f[ + * \kappa_T = -\frac{1}{v}\left(\frac{\partial v}{\partial P}\right)_T + * \f] + */ + virtual doublereal isothermalCompressibility() const; + + //! Return the volumetric thermal expansion coefficient. Units: 1/K. + /*! + * The thermal expansion coefficient is defined as + * \f[ + * \beta = \frac{1}{v}\left(\frac{\partial v}{\partial T}\right)_P + * \f] + */ + virtual doublereal thermalExpansionCoeff() const ; + + /** + * @} + * @name Activities, Standard States, and Activity Concentrations + * + * This section is largely handled by parent classes, since there + * is only one species. Therefore, the activity is equal to one. + * @{ + */ + + //! This method returns an array of generalized concentrations + /*! + * \f$ C^a_k\f$ are defined such that \f$ a_k = C^a_k / + * C^0_k, \f$ where \f$ C^0_k \f$ is a standard concentration + * defined below and \f$ a_k \f$ are activities used in the + * thermodynamic functions. These activity (or generalized) + * concentrations are used + * by kinetics manager classes to compute the forward and + * reverse rates of elementary reactions. + * + * For a stoichiomeetric substance, there is + * only one species, and the generalized concentration is 1.0. + * + * @param c Output array of generalized concentrations. The + * units depend upon the implementation of the + * reaction rate expressions within the phase. + */ + virtual void getActivityConcentrations(doublereal* c) const; + + //! Return the standard concentration for the kth species + /*! + * The standard concentration \f$ C^0_k \f$ used to normalize + * the activity (i.e., generalized) concentration. + * This phase assumes that the kinetics operator works on an + * dimensionless basis. Thus, the standard concentration is + * equal to 1.0. + * + * @param k Optional parameter indicating the species. The default + * is to assume this refers to species 0. + * @return + * Returns The standard Concentration as 1.0 + */ + virtual doublereal standardConcentration(int k=0) const; + + //! Natural logarithm of the standard concentration of the kth species. + /*! + * @param k index of the species (defaults to zero) + */ + virtual doublereal logStandardConc(int k=0) const; + + //! Get the array of chemical potentials at unit activity for the species + //! at their standard states at the current T and P of the solution. + /*! + * For a stoichiometric substance, there is no activity term in + * the chemical potential expression, and therefore the + * standard chemical potential and the chemical potential + * are both equal to the molar Gibbs function. + * + * These are the standard state chemical potentials \f$ \mu^0_k(T,P) + * \f$. The values are evaluated at the current + * temperature and pressure of the solution + * + * @param mu0 Output vector of chemical potentials. + * Length: m_kk. + */ + virtual void getStandardChemPotentials(doublereal* mu0) const; + + //! Returns the units of the standard and generalized concentrations. + /*! + * Note they have the same units, as their + * ratio is defined to be equal to the activity of the kth + * species in the solution, which is unitless. + * + * This routine is used in print out applications where the + * units are needed. Usually, MKS units are assumed throughout + * the program and in the XML input files. + * + * The base %ThermoPhase class assigns thedefault quantities + * of (kmol/m3) for all species. + * Inherited classes are responsible for overriding the default + * values if necessary. + * + * @param uA Output vector containing the units + * uA[0] = kmol units - default = 1 + * uA[1] = m units - default = -nDim(), the number of spatial + * dimensions in the Phase class. + * uA[2] = kg units - default = 0; + * uA[3] = Pa(pressure) units - default = 0; + * uA[4] = Temperature units - default = 0; + * uA[5] = time units - default = 0 + * @param k species index. Defaults to 0. + * @param sizeUA output int containing the size of the vector. + * Currently, this is equal to 6. + */ + virtual void getUnitsStandardConc(double *uA, int k = 0, + int sizeUA = 6); + + //@} + /// @name Partial Molar Properties of the Solution + /// + /// These properties are handled by the parent class, + /// SingleSpeciesTP + //@{ + + + //@} + /// @name Properties of the Standard State of the Species in the Solution + //@{ + + //! Get the nondimensional Enthalpy functions for the species + //! at their standard states at the current T and P of the solution. + /*! + * @param hrt Output vector of nondimensional standard state enthalpies. + * Length: m_kk. + */ + virtual void getEnthalpy_RT(doublereal* hrt) const; + + //! Get the array of nondimensional Entropy functions for the + //! standard state species at the current T and P of the solution. + /*! + * @param sr Output vector of nondimensional standard state entropies. + * Length: m_kk. + */ + virtual void getEntropy_R(doublereal* sr) const; + + //! Get the nondimensional Gibbs functions for the species + //! in their standard states at the current T and P of the solution. + /*! + * @param grt Output vector of nondimensional standard state gibbs free energies + * Length: m_kk. + */ + virtual void getGibbs_RT(doublereal* grt) const; + + //! Get the nondimensional Heat Capacities at constant + //! pressure for the species standard states + //! at the current T and P of the solution + /*! + * @param cpr Output vector of nondimensional standard state heat capacities + * Length: m_kk. + */ + virtual void getCp_R(doublereal* cpr) const; + + //! Returns the vector of nondimensional Internal Energies of the standard + //! state species at the current T and P of the solution + /*! + * For an incompressible, + * stoichiometric substance, the molar internal energy is + * independent of pressure. Since the thermodynamic properties + * are specified by giving the standard-state enthalpy, the + * term \f$ P_{ref} \hat v\f$ is subtracted from the specified reference molar + * enthalpy to compute the standard state molar internal energy. + * + * @param urt output vector of nondimensional standard state internal energies + * of the species. Length: m_kk. + */ + virtual void getIntEnergy_RT(doublereal* urt) const; + + //@} + /// @name Thermodynamic Values for the Species Reference States + //@{ + + //! Returns the vector of nondimensional + //! internal Energies of the reference state at the current temperature + //! of the solution and the reference pressure for each species. + /*! + * @param urt Output vector of nondimensional reference state + * internal energies of the species. + * Length: m_kk + */ + virtual void getIntEnergy_RT_ref(doublereal *urt) const; + + /* + * ---- Critical State Properties + */ + + + /* + * ---- Saturation Properties + */ + + /* + * @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. When importing a CTML phase + * description, this method is called just prior to returning + * from function importPhase. + * + * @see importCTML.cpp + */ + virtual void initThermo(); + + //! Set the equation of state parameters + /*! + * @internal + * The number and meaning of these depends on the subclass. + * + * @param n number of parameters + * @param c array of \a n coefficients + * c[0] = density of phase [ kg/m3 ] + */ + virtual void setParameters(int n, double *c); + + //! Get the equation of state parameters in a vector + /*! + * @internal + * + * @param n number of parameters + * @param c array of \a n coefficients + * + * For this phase: + * - n = 1 + * - c[0] = density of phase [ kg/m3 ] + */ + virtual void getParameters(int &n, double * const c); + + //! Set equation of state parameter values from XML entries. + /*! + * This method is called by function importPhase() in + * file importCTML.cpp when processing a phase definition in + * an input file. It should be overloaded in subclasses to set + * any parameters that are specific to that particular phase + * model. Note, this method is called before the phase is + * initialzed with elements and/or species. + * + * For this phase, the density of the phase is specified in this block. + * + * @param eosdata An XML_Node object corresponding to + * the "thermo" entry for this phase in the input file. + * + * eosdata points to the thermo block, and looks like this: + * + * @verbatim + + + 3.52 + + @endverbatim * */ - class StoichSubstanceSSTP : public SingleSpeciesTP { + virtual void setParametersFromXML(const XML_Node& eosdata); - public: - /** - * Default Constructor for the StoichSubstanceSSTP class - */ - StoichSubstanceSSTP(); + protected: - /** - * Destructor for the routine (virtual) - * - */ - virtual ~StoichSubstanceSSTP(); - - /** - * - * @name Utilities - * @{ - */ - - /** - * Equation of state flag. - * - * Returns the value cStoichSubstance, defined in mix_defs.h. - */ - virtual int eosType() const; - - /** - * @} - * @name Molar Thermodynamic Properties of the Solution - * @{ - */ - - /** - * @} - * @name Mechanical Equation of State - * @{ - */ - - /** - * Pressure. Units: Pa. - * For an incompressible substance, the density is independent - * of pressure. This method simply returns the stored - * pressure value. - */ - virtual doublereal pressure() const; - - /** - * Set the pressure at constant temperature. Units: Pa. - * For an incompressible substance, the density is - * independent of pressure. Therefore, this method only - * stores the specified pressure value. It does not - * modify the density. - */ - virtual void setPressure(doublereal p); - - /** - * The isothermal compressibility. Units: 1/Pa. - * The isothermal compressibility is defined as - * \f[ - * \kappa_T = -\frac{1}{v}\left(\frac{\partial v}{\partial P}\right)_T - * \f] - */ - virtual doublereal isothermalCompressibility() const; - - /** - * The thermal expansion coefficient. Units: 1/K. - * The thermal expansion coefficient is defined as - * - * \f[ - * \beta = \frac{1}{v}\left(\frac{\partial v}{\partial T}\right)_P - * \f] - */ - virtual doublereal thermalExpansionCoeff() const ; - - - /** - * @} - * @name Activities, Standard States, and Activity Concentrations - * - * This section is largely handled by parent classes, since there - * is only one species. Therefore, the activity is equal to one. - * @{ - */ - - /** - * This method returns the array of generalized - * concentrations. For a stoichiomeetric substance, there is - * only one species, and the generalized concentration is 1.0. - */ - virtual void getActivityConcentrations(doublereal* c) const; - - /** - * The standard concentration. This is defined as the concentration - * by which the generalized concentration is normalized to produce - * the activity. - */ - virtual doublereal standardConcentration(int k=0) const; - - /** - * Returns the natural logarithm of the standard - * concentration of the kth species - */ - virtual doublereal logStandardConc(int k=0) const; - - /** - * Get the array of chemical potentials at unit activity - * \f$ \mu^0_k \f$. - * - * For a stoichiometric substance, there is no activity term in - * the chemical potential expression, and therefore the - * standard chemical potential and the chemical potential - * are both equal to the molar Gibbs function. - */ - virtual void getStandardChemPotentials(doublereal* mu0) const; - - /** - * Returns the units of the standard and generalized - * concentrations Note they have the same units, as their - * ratio is defined to be equal to the activity of the kth - * species in the solution, which is unitless. - * - * This routine is used in print out applications where the - * units are needed. Usually, MKS units are assumed throughout - * the program and in the XML input files. - * - * uA[0] = kmol units - default = 0 - * uA[1] = m units - default = 0 - * uA[2] = kg units - default = 0; - * uA[3] = Pa(pressure) units - default = 0; - * uA[4] = Temperature units - default = 0; - * uA[5] = time units - default = 0 - */ - virtual void getUnitsStandardConc(double *uA, int k = 0, - int sizeUA = 6); - - //@} - /// @name Partial Molar Properties of the Solution - /// - /// These properties are handled by the parent class, - /// SingleSpeciesTP - //@{ - - - //@} - /// @name Properties of the Standard State of the Species in the Solution - //@{ - - /** - * Get the nondimensional Enthalpy functions for the species - * at their standard states at the current - * T and P of the solution. - */ - virtual void getEnthalpy_RT(doublereal* hrt) const; - - /** - * Get the array of nondimensional Entropy functions for the - * standard state species - * at the current T and P of the solution. - */ - virtual void getEntropy_R(doublereal* sr) const; - - /** - * Get the nondimensional Gibbs functions for the species - * at their standard states of solution at the current T and P - * of the solution - */ - virtual void getGibbs_RT(doublereal* grt) const; - - /** - * Get the nondimensional Gibbs functions for the standard - * state of the species at the current T and P. - */ - virtual void getCp_R(doublereal* cpr) const; - - - /** - * Molar internal energy. J/kmol. For an incompressible, - * stoichiometric substance, the molar internal energy is - * independent of pressure. Since the thermodynamic properties - * are specified by giving the standard-state enthalpy, the - * term \f$ P_0 \hat v\f$ is subtracted from the specified molar - * enthalpy to compute the molar internal energy. - */ - virtual void getIntEnergy_RT(doublereal* urt) const; - - //@} - /// @name Thermodynamic Values for the Species Reference States - //@{ - - /** - * Returns the vector of nondimensional - * internal Energies of the reference state at the current temperature - * of the solution and the reference pressure for each species. - */ - virtual void getIntEnergy_RT_ref(doublereal *urt) const; - - /* - * ---- Critical State Properties - */ - /// Critical temperature (K). - virtual doublereal critTemperature() const; - /// Critical pressure (Pa). - virtual doublereal critPressure() const; - /// Critical density (kg/m3). - virtual doublereal critDensity() const; - - /* - * ---- Saturation Properties - */ - virtual doublereal satTemperature(doublereal p) const; - virtual doublereal satPressure(doublereal t) const; - virtual doublereal vaporFraction() const; - virtual void setState_Tsat(doublereal t, doublereal x); - virtual void setState_Psat(doublereal p, doublereal x); - - /* - * @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. When importing a CTML phase - * description, this method is called just prior to returning - * from function importPhase. - * - * @see importCTML.cpp - */ - virtual void initThermo(); - - /* - * setParameters: - * - * Generic routine that is used to set the parameters used - * by this model. - * C[0] = density of phase [ kg/m3 ] - */ - virtual void setParameters(int n, double *c); - /* - * getParameters: - * - * Generic routine that is used to get the parameters used - * by this model. - * n = 1 - * C[0] = density of phase [ kg/m3 ] - */ - virtual void getParameters(int &n, double * const c); - - /* - * Reads an xml data block for the parameters needed by this - * routine. eosdata points to the thermo block, and looks - * like this: - * - * - * - * 3.52 - * - * - */ - virtual void setParametersFromXML(const XML_Node& eosdata); - - protected: - - }; + }; } diff --git a/test_problems/cathermo/Makefile.in b/test_problems/cathermo/Makefile.in index 0c5896a9c..c1e8f9bce 100644 --- a/test_problems/cathermo/Makefile.in +++ b/test_problems/cathermo/Makefile.in @@ -10,6 +10,7 @@ test_electrolytes=@COMPILE_ELECTROLYTES@ all: ifeq ($(test_issp),1) cd issp; @MAKE@ all + cd stoichSubSSTP; @MAKE@ all endif ifeq ($(test_electrolytes),1) cd ims; @MAKE@ all @@ -31,6 +32,7 @@ endif test: ifeq ($(test_issp),1) cd issp; @MAKE@ -s test + cd stoichSubSSTP; @MAKE@ -s test endif ifeq ($(test_electrolytes),1) cd ims; @MAKE@ -s test @@ -52,6 +54,7 @@ endif clean: $(RM) *.*~ cd issp; @MAKE@ clean + cd stoichSubSSTP; @MAKE@ clean cd ims; @MAKE@ clean cd testIAPWS; @MAKE@ clean cd testIAPWSPres; @MAKE@ clean @@ -69,7 +72,8 @@ clean: depends: ifeq ($(test_issp),1) - cd issp;@MAKE@ depends + cd issp; @MAKE@ depends + cd stoichSubSSTP; @MAKE@ clean endif ifeq ($(test_electrolytes),1) cd ims; @MAKE@ depends diff --git a/test_problems/cathermo/stoichSubSSTP/.cvsignore b/test_problems/cathermo/stoichSubSSTP/.cvsignore new file mode 100644 index 000000000..b52b3c16a --- /dev/null +++ b/test_problems/cathermo/stoichSubSSTP/.cvsignore @@ -0,0 +1,14 @@ +Makefile +.cvsignore.swp +.depends +Gex_standalone +HMW_graph_GvT +HMW_graph_GvT.d +diff_test.out +output.txt +outputa.txt +sortAlgorithms.d +csvCode.txt +ct2ctml.log +stoichSubSSTP +stoichSubSSTP.d diff --git a/test_problems/cathermo/stoichSubSSTP/Makefile.in b/test_problems/cathermo/stoichSubSSTP/Makefile.in new file mode 100644 index 000000000..3d9d463da --- /dev/null +++ b/test_problems/cathermo/stoichSubSSTP/Makefile.in @@ -0,0 +1,112 @@ +#!/bin/sh + +############################################################################ +# +# Makefile to compile and link a C++ application to +# Cantera. +# +############################################################################# + +# addition to suffixes +.SUFFIXES : .d + +# the name of the executable program to be created +PROG_NAME = stoichSubSSTP + +# the object files to be linked together. List those generated from Fortran +# and from C/C++ separately +OBJS = stoichSubSSTP.o sortAlgorithms.o + +# Location of the current build. Will assume that tests are run +# in the source directory tree location +src_dir_tree = 1 + +# additional flags to be passed to the linker. If your program +# requires other external libraries, put them here +LINK_OPTIONS = @EXTRA_LINK@ + +############################################################################# + +# Check to see whether we are in the msvc++ environment +os_is_win = @OS_IS_WIN@ + +# Fortran libraries +FORT_LIBS = @FLIBS@ + +# the C++ compiler +CXX = @CXX@ + +# C++ compile flags +ifeq ($(src_dir_tree), 1) +CXX_FLAGS = -DSRCDIRTREE @CXXFLAGS@ +else +CXX_FLAGS = @CXXFLAGS@ +endif + +# Ending C++ linking libraries +LCXX_END_LIBS = @LCXX_END_LIBS@ + +# the directory where the Cantera libraries are located +CANTERA_LIBDIR=@buildlib@ + +# required Cantera libraries +CANTERA_LIBS = @LOCAL_LIBS@ -lctcxx + +# the directory where Cantera include files may be found. +ifeq ($(src_dir_tree), 1) +CANTERA_INCDIR=../../../Cantera/src +INCLUDES=-I$(CANTERA_INCDIR) -I$(CANTERA_INCDIR)/thermo +else +CANTERA_INCDIR=@ctroot@/build/include/cantera +INCLUDES=-I$(CANTERA_INCDIR) -I$(CANTERA_INCDIR)/kernel +endif + +# flags passed to the C++ compiler/linker for the linking step +LCXX_FLAGS = -L$(CANTERA_LIBDIR) @LOCAL_LIB_DIRS@ @CXXFLAGS@ + +# How to compile C++ source files to object files +.@CXX_EXT@.@OBJ_EXT@: + $(CXX) -c $< $(INCLUDES) $(CXX_FLAGS) + +# How to compile the dependency file +.cpp.d: + @CXX_DEPENDS@ $(INCLUDES) $(CXX_FLAGS) $*.cpp > $*.d + +# List of dependency files to be created +DEPENDS=$(OBJS:.o=.d) + +# Program Name +PROGRAM = $(PROG_NAME)$(EXE_EXT) + +all: $(PROGRAM) .depends + +$(PROGRAM): $(OBJS) $(CANTERA_LIBDIR)/libcantera.a \ + $(CANTERA_LIBDIR)/libcaThermo.a + $(CXX) -o $(PROGRAM) $(OBJS) $(LCXX_FLAGS) $(LINK_OPTIONS) \ + $(CANTERA_LIBS) @LIBS@ $(FORT_LIBS) \ + $(LCXX_END_LIBS) + +# depends target -> forces recalculation of dependencies +depends: + @MAKE@ .depends + +.depends: $(DEPENDS) + cat $(DEPENDS) > .depends + +# Do the test -> For the windows vc++ environment, we have to skip checking on +# whether the program is uptodate, because we don't utilize make +# in that environment to build programs. +test: +ifeq ($(os_is_win), 1) +else + @ @MAKE@ -s $(PROGRAM) +endif + @ ./runtest + +clean: + $(RM) $(OBJS) $(PROGRAM) $(DEPENDS) .depends *.o + ../../../bin/rm_cvsignore + (if test -d SunWS_cache ; then \ + $(RM) -rf SunWS_cache ; \ + fi ) + diff --git a/test_problems/cathermo/stoichSubSSTP/NaCl_Solid.xml b/test_problems/cathermo/stoichSubSSTP/NaCl_Solid.xml new file mode 100644 index 000000000..d7d1fdc96 --- /dev/null +++ b/test_problems/cathermo/stoichSubSSTP/NaCl_Solid.xml @@ -0,0 +1,39 @@ + + + + + + + + + O H C Fe Ca N Na Cl + + NaCl(S) + + 2.165 + + + + + + + + + + + Na:1 Cl:1 + + + + 50.72389, 6.672267, -2.517167, + 10.15934, -0.200675, -427.2115, + 130.3973 + + + + 2.165 + + + + + diff --git a/test_problems/cathermo/stoichSubSSTP/TemperatureTable.h b/test_problems/cathermo/stoichSubSSTP/TemperatureTable.h new file mode 100644 index 000000000..70cb2bf56 --- /dev/null +++ b/test_problems/cathermo/stoichSubSSTP/TemperatureTable.h @@ -0,0 +1,129 @@ +/* + * $Id$ + */ +/* + * Copywrite 2004 Sandia Corporation. Under the terms of Contract + * DE-AC04-94AL85000, there is a non-exclusive license for use of this + * work by or on behalf of the U.S. Government. Export of this program + * may require a license from the United States Government. + */ + +#ifndef TEMPERATURE_TABLE_H +#define TEMPERATURE_TABLE_H +#include "sortAlgorithms.h" +//#include "mdp_allo.h" +#include +using std::vector; + +/***********************************************************************/ +/***********************************************************************/ +/***********************************************************************/ +/** + * This Class constructs a vector of temperature from which to make + * a table. + */ +class TemperatureTable { + +public: + int NPoints; + bool Include298; + double Tlow; //!< Min temperature for thermo data fit + double Thigh; //!< Max temperature for thermo table + double DeltaT; + vector T; + int numAddedTs; + vector AddedTempVector; +public: + /* + * Default constructor for TemperatureTable() + */ + TemperatureTable(const int nPts = 14, + const bool inc298 = true, + const double tlow = 300., + const double deltaT = 100., + const int numAdded = 0, + const double *addedTempVector = 0) : + NPoints(nPts), + Include298(inc298), + Tlow(tlow), + DeltaT(deltaT), + T(0), + numAddedTs(numAdded) { + /****************************/ + int i; + // AddedTempVector = mdp_alloc_dbl_1(numAdded, 0.0); + AddedTempVector.resize(numAdded, 0.0); + for (int i = 0; i < numAdded; i++) { + AddedTempVector[i] = addedTempVector[i]; + } + //mdp_copy_dbl_1(AddedTempVector, addedTempVector, numAdded); + // T = mdp_alloc_dbl_1(NPoints, 0.0); + T.resize(NPoints, 0.0); + double TCurrent = Tlow; + for (i = 0; i < NPoints; i++) { + T[i] = TCurrent; + TCurrent += DeltaT; + } + if (Include298) { + T.push_back(298.15); + //mdp_realloc_dbl_1(&T, NPoints+1, NPoints, 298.15); + NPoints++; + } + if (numAdded > 0) { + //mdp_realloc_dbl_1(&T, NPoints+numAdded, NPoints, 0.0); + T.resize( NPoints+numAdded, 0.0); + for (i = 0; i < numAdded; i++) { + T[i+NPoints] = addedTempVector[i]; + } + NPoints += numAdded; + } + + sort_dbl_1(DATA_PTR(T), NPoints); + + + } + /***********************************************************************/ + /***********************************************************************/ + /***********************************************************************/ + /* + * Destructor + */ + ~TemperatureTable() { + //mdp_safe_free((void **) &AddedTempVector); + // mdp_safe_free((void **) &T); + } + + /***********************************************************************/ + /***********************************************************************/ + /***********************************************************************/ + /* + * Overloaded operator[] + * + * return the array value in the vector + */ + double operator[](const int i) { + return T[i]; + } + /***********************************************************************/ + /***********************************************************************/ + /***********************************************************************/ + /* + * size() + */ + int size() { + return NPoints; + } +/***********************************************************************/ +/***********************************************************************/ +/***********************************************************************/ + /* + * Block assignment and copy constructors: not needed. + */ +private: + TemperatureTable(const TemperatureTable &); + TemperatureTable& operator=(const TemperatureTable&); +}; +/***********************************************************************/ +/***********************************************************************/ +/***********************************************************************/ +#endif diff --git a/test_problems/cathermo/stoichSubSSTP/output_blessed.txt b/test_problems/cathermo/stoichSubSSTP/output_blessed.txt new file mode 100644 index 000000000..b4a0a1a3a --- /dev/null +++ b/test_problems/cathermo/stoichSubSSTP/output_blessed.txt @@ -0,0 +1,13 @@ + Data from http://webbook.nist.gov + + T, Pres, molarGibbs0, Enthalpy, Entropy, Cp , -(G-H298)/T, H-H298 + Kelvin, bars, kJ/gmol, kJ/gmol, J/gmolK, J/gmolK , J/gmolK, J/gmol + 298.15, 1.01325, -432.62, -411.121, 72.1093, 50.5012, 72.1093, 0 + 300, 1.01325, -432.754, -411.027, 72.4218, 50.5436, 72.1103, 0.0934666 + 400, 1.01325, -440.767, -405.875, 87.2308, 52.386, 74.117, 5.24556 + 500, 1.01325, -450.102, -400.56, 99.0843, 53.8979, 77.9635, 10.5604 + 600, 1.01325, -460.522, -395.094, 109.047, 55.4581, 82.3351, 16.0269 + 700, 1.01325, -471.869, -389.461, 117.726, 57.2362, 86.7837, 21.6593 + 800, 1.01325, -484.037, -383.636, 125.502, 59.3387, 91.1453, 27.485 + 900, 1.01325, -496.948, -377.58, 132.631, 61.8484, 95.3638, 33.5407 + 1000, 1.01325, -510.548, -371.25, 139.298, 64.8377, 99.4271, 39.8707 diff --git a/test_problems/cathermo/stoichSubSSTP/runtest b/test_problems/cathermo/stoichSubSSTP/runtest new file mode 100755 index 000000000..49933c3a1 --- /dev/null +++ b/test_problems/cathermo/stoichSubSSTP/runtest @@ -0,0 +1,42 @@ +#!/bin/sh +# +# +temp_success="1" +/bin/rm -f output.txt outputa.txt + +########################################################################## +prog=stoichSubSSTP +if test ! -x $prog ; then + echo $prog ' does not exist' + exit -1 +fi +########################################################################## +/bin/rm -f test.out test.diff output.txt + +################################################################# +# +CANTERA_DATA=${CANTERA_DATA:=../../../data/inputs}; export CANTERA_DATA +CANTERA_BIN=${CANTERA_BIN:=../../../bin} + +################################################################# + +$prog > output.txt +retnStat=$? +if [ $retnStat != "0" ] +then + temp_success="0" + echo "$prog returned with bad status, $retnStat, check output" +fi + +$CANTERA_BIN/exp3to2.sh output.txt > outputa.txt +diff -w outputa.txt output_blessed.txt > diff_test.out +retnStat=$? +if [ $retnStat = "0" ] +then + echo "successful diff comparison on $prog test" +else + echo "unsuccessful diff comparison on $prog test" + echo "FAILED" > csvCode.txt + temp_success="0" +fi + diff --git a/test_problems/cathermo/stoichSubSSTP/sortAlgorithms.cpp b/test_problems/cathermo/stoichSubSSTP/sortAlgorithms.cpp new file mode 100644 index 000000000..d97e51b40 --- /dev/null +++ b/test_problems/cathermo/stoichSubSSTP/sortAlgorithms.cpp @@ -0,0 +1,54 @@ +/* + * @file sortAlgorithms.h + * + * $Author$ + * $Revision$ + * $Date$ + */ +/* + * Copywrite 2004 Sandia Corporation. Under the terms of Contract + * DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government + * retains certain rights in this software. + * See file License.txt for licensing information. + */ + +#include "sortAlgorithms.h" + +/**************************************************************/ + +void sort_dbl_1(double * const x, const int n) { + double rra; + int ll = n/2; + int iret = n - 1; + while (1 > 0) { + if (ll > 0) { + ll--; + rra = x[ll]; + } else { + rra = x[iret]; + x[iret] = x[0]; + iret--; + if (iret == 0) { + x[0] = rra; + return; + } + } + int i = ll; + int j = ll + ll + 1; + while (j <= iret) { + if (j < iret) { + if (x[j] < x[j+1]) + j++; + } + if (rra < x[j]) { + x[i] = x[j]; + i = j; + j = j + j + 1; + } else { + j = iret + 1; + } + } + x[i] = rra; + } +} +/*****************************************************/ diff --git a/test_problems/cathermo/stoichSubSSTP/sortAlgorithms.h b/test_problems/cathermo/stoichSubSSTP/sortAlgorithms.h new file mode 100644 index 000000000..72a7fc2a2 --- /dev/null +++ b/test_problems/cathermo/stoichSubSSTP/sortAlgorithms.h @@ -0,0 +1,21 @@ +/* + * @file sortAlgorithms.h + * + * $Author$ + * $Revision$ + * $Date$ + */ +/* + * Copywrite 2004 Sandia Corporation. Under the terms of Contract + * DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government + * retains certain rights in this software. + * See file License.txt for licensing information. + */ + +#ifndef SORTALGORITHMS_H +#define SORTALGORITHMS_H + + +void sort_dbl_1(double * const x, const int n); + +#endif diff --git a/test_problems/cathermo/stoichSubSSTP/stoichSubSSTP.cpp b/test_problems/cathermo/stoichSubSSTP/stoichSubSSTP.cpp new file mode 100644 index 000000000..3a89e4782 --- /dev/null +++ b/test_problems/cathermo/stoichSubSSTP/stoichSubSSTP.cpp @@ -0,0 +1,205 @@ +/** + * + * @file HMW_graph_1.cpp + */ + +/* + * $Author$ + * $Date$ + * $Revision$ + */ +#include + +#ifdef SRCDIRTREE +#include "ct_defs.h" +#include "logger.h" +#include "ThermoPhase.h" +#include "StoichSubstanceSSTP.h" +#include "importCTML.h" +#else +#include "ThermoPhase.h" + +#include "cantera/Cantera.h" +#include "cantera/kernel/logger.h" +#include "cantera/thermo.h" +#include "cantera/kernel/thermo/HMWSoln.h" +#endif + +#include "TemperatureTable.h" + +using namespace std; +using namespace Cantera; + +class fileLog: public Logger { +public: + fileLog(string fName) { + m_fName = fName; + m_fs.open(fName.c_str()); + } + + virtual void write(const string& msg) { + m_fs << msg; + m_fs.flush(); + } + + virtual ~fileLog() { + m_fs.close(); + } + + string m_fName; + ofstream m_fs; + +}; + +void printUsage() { + cout << "usage: stoichSubSSTP " << endl; + cout <<" -> Everything is hardwired" << endl; +} + + + +int main(int argc, char **argv) +{ + + int retn = 0; + int i; + + try { + //Cantera::ThermoPhase *tp = 0; + char iFile[80], file_ID[80]; + strcpy(iFile, "NaCl_Solid.xml"); + if (argc > 1) { + strcpy(iFile, argv[1]); + } + + //fileLog *fl = new fileLog("HMW_graph_1.log"); + //setLogger(fl); + sprintf(file_ID,"%s#NaCl(S)", iFile); + XML_Node *xm = get_XML_NameID("phase", file_ID, 0); + StoichSubstanceSSTP *solid = new StoichSubstanceSSTP(*xm); + + + /* + * Load in and initialize the + */ + //string nacl_s = "NaCl_Solid.xml"; + //string id = "NaCl(S)"; + //Cantera::ThermoPhase *solid = Cantera::newPhase(nacl_s, id); + + + int nsp = solid->nSpecies(); + if (nsp != 1) { + throw CanteraError("","Should just be one species"); + } + double acMol[100]; + double act[100]; + double mf[100]; + double moll[100]; + for (i = 0; i < 100; i++) { + acMol[i] = 1.0; + act[i] = 1.0; + mf[i] = 0.0; + moll[i] = 0.0; + } + string sName; + + TemperatureTable TTable(8, true, 300, 100., 0, 0); + + /* + * Set the Pressure + */ + double pres = OneAtm; + double T = 298.15; + solid->setState_TP(T, pres); + + /* + * ThermoUnknowns + */ + double mu0_RT[20], mu[20], cp_r[20];; + double enth_RT[20]; + double entrop_RT[20], intE_RT[20]; + double mu_NaCl, enth_NaCl, entrop_NaCl; + double mu0_NaCl, molarGibbs, intE_NaCl, cp_NaCl; + /* + * Create a Table of NaCl Properties as a Function + * of the Temperature + */ + + double RT = GasConstant * T; + solid->getEnthalpy_RT(enth_RT); + double enth_NaCl_298 = enth_RT[0] * RT * 1.0E-6; + + printf(" Data from http://webbook.nist.gov\n"); + printf("\n"); + + + printf(" T, Pres, molarGibbs0, Enthalpy, Entropy, Cp ," + " -(G-H298)/T, H-H298 "); + printf("\n"); + + printf(" Kelvin, bars, kJ/gmol, kJ/gmol, J/gmolK, J/gmolK ," + " J/gmolK, J/gmol"); + printf("\n"); + + for (i = 0; i < TTable.NPoints; i++) { + T = TTable.T[i]; + + // GasConstant is in J/kmol + RT = GasConstant * T; + + pres = OneAtm; + + + solid->setState_TP(T, pres); + /* + * Get the Standard State DeltaH + */ + solid->getGibbs_RT(mu0_RT); + mu0_NaCl = mu0_RT[0] * RT * 1.0E-6; + + solid->getEnthalpy_RT(enth_RT); + enth_NaCl = enth_RT[0] * RT * 1.0E-6; + + + solid->getChemPotentials(mu); + mu_NaCl = mu[0] * 1.0E-6; + + solid->getEntropy_R(entrop_RT); + entrop_NaCl = entrop_RT[0] * GasConstant * 1.0E-3; + + molarGibbs = solid->gibbs_mole() * 1.0E-6; + + solid->getIntEnergy_RT(intE_RT); + intE_NaCl = intE_RT[0] * RT * 1.0E-6; + + solid->getCp_R(cp_r); + cp_NaCl = cp_r[0] * GasConstant * 1.0E-3; + + /* + * Need the gas constant in kJ/gmolK + */ + // double rgas = 8.314472 * 1.0E-3; + + double pbar = pres * 1.0E-5; + + printf("%10g, %10g, %12g, %12g, %12g, %12g, %12g, %12g", + T, pbar, mu_NaCl, enth_NaCl, entrop_NaCl, cp_NaCl, -1.0E3*(mu_NaCl-enth_NaCl_298)/T, enth_NaCl-enth_NaCl_298); + printf("\n"); + } + + + + delete solid; + solid = 0; + Cantera::appdelete(); + + return retn; + + } catch (CanteraError) { + + showErrors(); + Cantera::appdelete(); + return -1; + } + return 0; +} diff --git a/tools/doc/Cantera.cfg.in b/tools/doc/Cantera.cfg.in index 71d2150b6..1a088bb6e 100755 --- a/tools/doc/Cantera.cfg.in +++ b/tools/doc/Cantera.cfg.in @@ -114,7 +114,8 @@ FILE_PATTERNS = Kinetics.h Kinetics.cpp \ SingleSpeciesTP.h SingleSpeciesTP.cpp \ MolalityVPSSTP.h MolalityVPSSTP.cpp \ IdealMolalSoln.h IdealMolalSoln.cpp \ - IdealSolidSolnPhase.h IdealSolidSolnPhase.cpp + IdealSolidSolnPhase.h IdealSolidSolnPhase.cpp \ + StoichSubstanceSSTP.h StoichSubstanceSSTP.cpp RECURSIVE = NO EXCLUDE = CVS examples converters zeroD EXCLUDE_SYMLINKS = NO