From 3d8b7f857d8382b24d2881b5a46bfad36225d31d Mon Sep 17 00:00:00 2001 From: Harry Moffat Date: Sun, 4 Jan 2009 21:28:02 +0000 Subject: [PATCH] Added a few routines. --- Cantera/src/thermo/Constituents.cpp | 21 +- Cantera/src/thermo/Constituents.h | 613 ++++++++++++------------- Cantera/src/thermo/StoichSubstance.cpp | 5 + Cantera/src/thermo/StoichSubstance.h | 10 + Cantera/src/thermo/SurfPhase.cpp | 463 +++++++++---------- Cantera/src/thermo/SurfPhase.h | 17 +- 6 files changed, 566 insertions(+), 563 deletions(-) diff --git a/Cantera/src/thermo/Constituents.cpp b/Cantera/src/thermo/Constituents.cpp index 894536de8..2d42d9564 100755 --- a/Cantera/src/thermo/Constituents.cpp +++ b/Cantera/src/thermo/Constituents.cpp @@ -4,7 +4,7 @@ * manages a set of elements and species (see \ref phases). */ -/* $Author$ +/* * $Date$ * $Revision$ */ @@ -88,7 +88,7 @@ namespace Cantera { int Constituents::nElements() const { return m_Elements->nElements(); } - /** + /* * Return the Atomic weight of element m. * units = Kg / Kmol */ @@ -96,12 +96,12 @@ namespace Cantera { return m_Elements->atomicWeight(m); } - + doublereal Constituents::entropyElement298(int m) const { return m_Elements->entropyElement298(m); } - /** + /* * returns a reference to the vector of atomic weights pertinent * to this constituents object * units = kg / Kmol @@ -111,7 +111,7 @@ namespace Cantera { } - /** + /* * Return the atomic number of element m. */ int Constituents::atomicNumber(int m) const { @@ -119,7 +119,7 @@ namespace Cantera { } - /** + /* * Add an element to the set. * @param symbol symbol string * @param weight atomic weight in kg/mol. @@ -230,8 +230,7 @@ namespace Cantera { return m_Elements->elementNames(); } - /********************************************************************** - * + /* * molecularWeight() * * Returns the molecular weight of a species given the species index @@ -246,8 +245,7 @@ namespace Cantera { return m_weight[k]; } - /********************************************************************** - * + /* * molecularWeights() * * Returns a const reference to the vector of molecular weights @@ -259,8 +257,7 @@ namespace Cantera { return m_weight; } - /********************************************************************** - * + /* * charge(): * * Electrical charge of one species k molecule, divided by diff --git a/Cantera/src/thermo/Constituents.h b/Cantera/src/thermo/Constituents.h index 8cba4a52f..d5d5fe58d 100755 --- a/Cantera/src/thermo/Constituents.h +++ b/Cantera/src/thermo/Constituents.h @@ -4,10 +4,9 @@ * manages a set of elements and species (see \ref phases). */ -/* $Author$ +/* * $Date$ * $Revision$ - * */ // Copyright 2001 California Institute of Technology @@ -27,378 +26,376 @@ namespace Cantera { - class Elements; + class Elements; - /************** DEFINITIONS OF ERRORS *****************************/ + /************** DEFINITIONS OF ERRORS *****************************/ - //! Specific fatal error indicating that the index of a species is out of range. + //! Specific fatal error indicating that the index of a species is out of range. + /*! + * + * @ingroup errorhandling + */ + class SpeciesRangeError : public CanteraError { + public: + //! Constructor /*! + * @param func Function where the error occurred. + * @param k current species index value + * @param kmax Maximum permissible species index value. The + * minimum permissible species index value is assumed to be 0 * - * @ingroup errorhandling */ - class SpeciesRangeError : public CanteraError { - public: - //! Constructor - /*! - * @param func Function where the error occurred. - * @param k current species index value - * @param kmax Maximum permissible species index value. The - * minimum permissible species index value is assumed to be 0 - * - */ - SpeciesRangeError(std::string func, int k, int kmax) : - CanteraError(func, "Species index " + int2str(k) + - " outside valid range of 0 to " + int2str(kmax-1)) {} - }; + SpeciesRangeError(std::string func, int k, int kmax) : + CanteraError(func, "Species index " + int2str(k) + + " outside valid range of 0 to " + int2str(kmax-1)) {} + }; - /******************************************************************/ + /******************************************************************/ - //! Class %Constituents manages a set of elements and species. + //! Class %Constituents manages a set of elements and species. + /*! + * Class %Constituents is designed to provide information + * about the elements and species in a phase - names, index + * numbers (location in arrays), atomic or molecular weights, + * etc. No computations are performed by the methods of this + * class. The set of elements must include all those that compose + * the species, but may include additional elements. The species + * all must belong to the same phase. + * + * @ingroup phases + */ + class Constituents { + + public: + + //! Constructor. /*! - * Class %Constituents is designed to provide information - * about the elements and species in a phase - names, index - * numbers (location in arrays), atomic or molecular weights, - * etc. No computations are performed by the methods of this - * class. The set of elements must include all those that compose - * the species, but may include additional elements. The species - * all must belong to the same phase. + * Constructor sets all base variable types to zero. Also, it + * sets the pointer to the Elements object for this object. * - * @ingroup phases + * @param ptr_Elements + * The default is that a new Elements object is created, so this + * Constituents object is independent of any other object. But if + * ptr_Elements is supplied, it will be used. This way, a class + * implementing a multi-phase mixture is responsible for + * maintaining the global elements list for the mixture, and no + * static global element list is required. */ - class Constituents { - - public: - - //! Constructor. - /*! - * Constructor sets all base variable types to zero. Also, it - * sets the pointer to the Elements object for this object. - * - * @param ptr_Elements - * The default is that a new Elements object is created, so this - * Constituents object is independent of any other object. But if - * ptr_Elements is supplied, it will be used. This way, a class - * implementing a multi-phase mixture is responsible for - * maintaining the global elements list for the mixture, and no - * static global element list is required. - */ - Constituents(Elements* ptr_Elements = 0); + Constituents(Elements* ptr_Elements = 0); - /// Destructor. - ~Constituents(); + /// Destructor. + ~Constituents(); - /// This copy constructor just calls the assignment operator - /// for this class. - /*! - * @param right reference to the object to be copied. - */ - Constituents(const Constituents& right); + /// This copy constructor just calls the assignment operator + /// for this class. + /*! + * @param right reference to the object to be copied. + */ + Constituents(const Constituents& right); - /// Assignment operator - /*! - * @param right Reference to the object to be copied. - */ - Constituents& operator=(const Constituents& right); + /// Assignment operator + /*! + * @param right Reference to the object to be copied. + */ + Constituents& operator=(const Constituents& right); - /// @name Element Information - // @{ + /// @name Element Information + // @{ - /// Name of the element with index m. - /// This is a passthrough routine to the Element object. - /// \param m Element index. - /// \exception If m < 0 or m >= nElements(), the - /// exception, ElementRangeError, is thrown. - std::string elementName(int m) const; + /// Name of the element with index m. + /// This is a passthrough routine to the Element object. + /// \param m Element index. + /// \exception If m < 0 or m >= nElements(), the + /// exception, ElementRangeError, is thrown. + std::string elementName(int m) const; - /// Index of element named 'name'. - /// The index is an integer - /// assigned to each element in the order it was added, - /// beginning with 0 for the first element. - /// @param name name of the element - /// - /// If 'name' is not - /// the name of an element in the set, then the value -1 is - /// returned. - int elementIndex(std::string name) const; + /// Index of element named 'name'. + /// The index is an integer + /// assigned to each element in the order it was added, + /// beginning with 0 for the first element. + /// @param name name of the element + /// + /// If 'name' is not + /// the name of an element in the set, then the value -1 is + /// returned. + int elementIndex(std::string name) const; - /// Atomic weight of element m. - /*! - * @param m Element index - */ - doublereal atomicWeight(int m) const; + /// Atomic weight of element m. + /*! + * @param m Element index + */ + doublereal atomicWeight(int m) const; - /// Entropy of the element in its standard state at 298 K and 1 bar - /*! - * @param m Element index - */ - doublereal entropyElement298(int m) const; + /// Entropy of the element in its standard state at 298 K and 1 bar + /*! + * @param m Element index + */ + doublereal entropyElement298(int m) const; - /// Atomic number of element m. - /*! - * @param m Element index - */ - int atomicNumber(int m) const; + /// Atomic number of element m. + /*! + * @param m Element index + */ + int atomicNumber(int m) const; - /// Return a read-only reference to the vector of element names. - const std::vector& elementNames() const; + /// Return a read-only reference to the vector of element names. + const std::vector& elementNames() const; - /// Return a read-only reference to the vector of atomic weights. - const vector_fp& atomicWeights() const; + /// Return a read-only reference to the vector of atomic weights. + const vector_fp& atomicWeights() const; - /// Number of elements. - int nElements() const; + /// Number of elements. + int nElements() const; - // @} + // @} - /// @name Adding Elements and Species - /// These methods are used to add new elements or species. - /// These are not usually called by user programs. - /// - /// Since species are checked to insure that they are only - /// composed of declared elements, it is necessary to first - /// add all elements before adding any species. + /// @name Adding Elements and Species + /// These methods are used to add new elements or species. + /// These are not usually called by user programs. + /// + /// Since species are checked to insure that they are only + /// composed of declared elements, it is necessary to first + /// add all elements before adding any species. - //@{ + //@{ - //! Add an element. - /*! - * @param symbol Atomic symbol std::string. - * @param weight Atomic mass in amu. - */ - void addElement(const std::string& symbol, doublereal weight); + //! Add an element. + /*! + * @param symbol Atomic symbol std::string. + * @param weight Atomic mass in amu. + */ + void addElement(const std::string& symbol, doublereal weight); - //! Add an element from an XML specification. - /*! - * @param e Reference to the XML_Node where the element is described. - */ - void addElement(const XML_Node& e); + //! Add an element from an XML specification. + /*! + * @param e Reference to the XML_Node where the element is described. + */ + void addElement(const XML_Node& e); - //! Add an element, checking for uniqueness - /*! - * The uniqueness is checked by comparing the string symbol. If - * not unique, nothing is done. - * - * @param symbol String symbol of the element - * @param weight Atomic weight of the element (kg kmol-1). - * @param atomicNumber Atomic number of the element (unitless) - * @param entropy298 Entropy of the element at 298 K and 1 bar - * in its most stable form. The default is - * the value ENTROPY298_UNKNOWN, which is - * interpreted as an unknown, and if used - * will cause Cantera to throw an error. - */ - void addUniqueElement(const std::string& symbol, doublereal weight, - int atomicNumber = 0, - doublereal entropy298 = ENTROPY298_UNKNOWN); + //! Add an element, checking for uniqueness + /*! + * The uniqueness is checked by comparing the string symbol. If + * not unique, nothing is done. + * + * @param symbol String symbol of the element + * @param weight Atomic weight of the element (kg kmol-1). + * @param atomicNumber Atomic number of the element (unitless) + * @param entropy298 Entropy of the element at 298 K and 1 bar + * in its most stable form. The default is + * the value ENTROPY298_UNKNOWN, which is + * interpreted as an unknown, and if used + * will cause Cantera to throw an error. + */ + void addUniqueElement(const std::string& symbol, doublereal weight, + int atomicNumber = 0, + doublereal entropy298 = ENTROPY298_UNKNOWN); - //! Adde an element, checking for uniqueness - /*! - * The uniqueness is checked by comparing the string symbol. If - * not unique, nothing is done. - * - * @param e Reference to the XML_Node where the element is described. - */ - void addUniqueElement(const XML_Node& e); + //! Adde an element, checking for uniqueness + /*! + * The uniqueness is checked by comparing the string symbol. If + * not unique, nothing is done. + * + * @param e Reference to the XML_Node where the element is described. + */ + void addUniqueElement(const XML_Node& e); - //! Add all elements referenced in an XML_Node tree - /*! - * @param phase Reference to the top XML_Node of a phase - */ - void addElementsFromXML(const XML_Node& phase); + //! Add all elements referenced in an XML_Node tree + /*! + * @param phase Reference to the top XML_Node of a phase + */ + void addElementsFromXML(const XML_Node& phase); - /// Prohibit addition of more elements, and prepare to add species. - void freezeElements(); + /// Prohibit addition of more elements, and prepare to add species. + void freezeElements(); - /// True if freezeElements has been called. - bool elementsFrozen(); + /// True if freezeElements has been called. + bool elementsFrozen(); - //@} + //@} - /// Returns the number of species in the phase - int nSpecies() const { return m_kk; } + /// Returns the number of species in the phase + int nSpecies() const { return m_kk; } - //! Molecular weight of species \c k. - /*! - * @param k index of species \c k - * @return - * Returns the molecular weight of species \c k. - */ - doublereal molecularWeight(int k) const; + //! Molecular weight of species \c k. + /*! + * @param k index of species \c k + * @return + * Returns the molecular weight of species \c k. + */ + doublereal molecularWeight(int k) const; - //! Return the Molar mass of species \c k - /*! - * Preferred name for molecular weight. - * - * @param k index for species - * @return - * Return the molar mass of species k kg/kmol. - */ - doublereal molarMass(int k) const { - return molecularWeight(k); - } + //! Return the Molar mass of species \c k + /*! + * Preferred name for molecular weight. + * + * @param k index for species + * @return + * Return the molar mass of species k kg/kmol. + */ + doublereal molarMass(int k) const { + return molecularWeight(k); + } - /** - * Return a const reference to the vector of molecular weights - * of the species - */ - const vector_fp& molecularWeights() const; + /** + * Return a const reference to the vector of molecular weights + * of the species + */ + const vector_fp& molecularWeights() const; - /*! - * Electrical charge of one species k molecule, divided by - * the magnitude of the electron charge ( \f$ e = 1.602 - * \times 10^{-19}\f$ Coulombs). Dimensionless. - * - * @param k species index - */ - doublereal charge(int k) const; + /*! + * Electrical charge of one species k molecule, divided by + * the magnitude of the electron charge ( \f$ e = 1.602 + * \times 10^{-19}\f$ Coulombs). Dimensionless. + * + * @param k species index + */ + doublereal charge(int k) const; - /** - * @name Adding Species - * These methods are used to add new species. - * They are not usually called by user programs. - */ - //@{ - void addSpecies(const std::string& name, const doublereal* comp, - doublereal charge = 0.0, doublereal size = 1.0); + /** + * @name Adding Species + * These methods are used to add new species. + * They are not usually called by user programs. + */ + //@{ + void addSpecies(const std::string& name, const doublereal* comp, + doublereal charge = 0.0, doublereal size = 1.0); - //! Add a species to the phase, checking for uniqueness of the name - /*! - * This routine checks for uniqueness of the string name. It only - * adds the species if it is unique. - * - * @param name String name of the species - * @param comp Double vector containing the elemental composition of the - * species. - * @param charge Charge of the species. Defaults to zero. - * @param size Size of the species (meters). Defaults to 1 meter. - */ - void addUniqueSpecies(const std::string& name, const doublereal* comp, - doublereal charge = 0.0, - doublereal size = 1.0); + //! Add a species to the phase, checking for uniqueness of the name + /*! + * This routine checks for uniqueness of the string name. It only + * adds the species if it is unique. + * + * @param name String name of the species + * @param comp Double vector containing the elemental composition of the + * species. + * @param charge Charge of the species. Defaults to zero. + * @param size Size of the species (meters). Defaults to 1 meter. + */ + void addUniqueSpecies(const std::string& name, const doublereal* comp, + doublereal charge = 0.0, + doublereal size = 1.0); - //! Index of species named 'name' - /*! - * The first species added - * will have index 0, and the last one index nSpecies() - 1. - * - * @param name String name of the species - * @return - * Returns the index of the species. - */ - int speciesIndex(std::string name) const; + //! Index of species named 'name' + /*! + * The first species added + * will have index 0, and the last one index nSpecies() - 1. + * + * @param name String name of the species + * @return + * Returns the index of the species. + */ + int speciesIndex(std::string name) const; - //! Name of the species with index k - /*! - * @param k index of the species - */ - std::string speciesName(int k) const; + //! Name of the species with index k + /*! + * @param k index of the species + */ + std::string speciesName(int k) const; - /// Return a const referernce to the vector of species names - const std::vector& speciesNames() const; + /// Return a const referernce to the vector of species names + const std::vector& speciesNames() const; - //! This routine returns the size of species k - /*! - * @param k index of the species - * @return - * Returns the size of the species. Units are meters. - */ - doublereal size(int k) const { return m_speciesSize[k]; } + //! This routine returns the size of species k + /*! + * @param k index of the species + * @return + * Returns the size of the species. Units are meters. + */ + doublereal size(int k) const { return m_speciesSize[k]; } - /** - * Prohibit addition of more species, and prepare for - * calculations with this set of elements and species. - */ - void freezeSpecies(); + /** + * Prohibit addition of more species, and prepare for + * calculations with this set of elements and species. + */ + void freezeSpecies(); - /// True if freezeSpecies has been called. - bool speciesFrozen() { return m_speciesFrozen; } + /// True if freezeSpecies has been called. + bool speciesFrozen() { return m_speciesFrozen; } - /// Remove all elements and species - void clear(); + /// Remove all elements and species + void clear(); - //@} + //@} - /// True if both elements and species have been frozen - bool ready() const; + /// True if both elements and species have been frozen + bool ready() const; - //! Number of atoms of element \c m in species \c k. - /*! - * @param k species index - * @param m element index - */ - doublereal nAtoms(int k, int m) const; + //! Number of atoms of element \c m in species \c k. + /*! + * @param k species index + * @param m element index + */ + doublereal nAtoms(int k, int m) const; - - //! Get a vector containing the atomic composition of species k - /*! - * @param k species index - * @param atomArray vector containing the atomic number in the species. - * Length: m_mm - */ - void getAtoms(int k, double *atomArray) const; - - - protected: + //! Get a vector containing the atomic composition of species k + /*! + * @param k species index + * @param atomArray vector containing the atomic number in the species. + * Length: m_mm + */ + void getAtoms(int k, double *atomArray) const; - //! Number of species in the phase. - int m_kk; - //! Vector of molecular weights of the species - /*! - * This vector has length m_kk. - * The units of the vector are kg kmol-1. - */ - vector_fp m_weight; + protected: + + //! Number of species in the phase. + int m_kk; + //! Vector of molecular weights of the species + /*! + * This vector has length m_kk. + * The units of the vector are kg kmol-1. + */ + vector_fp m_weight; - //! Boolean indicating whether the number of species has been frozen. - /*! - * During the construction of the phase, this is false. After - * construction of the the phase, this is true. - */ - bool m_speciesFrozen; + //! Boolean indicating whether the number of species has been frozen. + /*! + * During the construction of the phase, this is false. After + * construction of the the phase, this is true. + */ + bool m_speciesFrozen; - /*! - * Pointer to the element object corresponding to this - * phase. Normally, this will be the default Element object - * common to all phases. - */ - Elements * m_Elements; + /*! + * Pointer to the element object corresponding to this + * phase. Normally, this will be the default Element object + * common to all phases. + */ + Elements * m_Elements; - //! Vector of the species names - std::vector m_speciesNames; + //! Vector of the species names + std::vector m_speciesNames; - //! Atomic composition of the species. - /*! - * the number of atoms of i in species k is equal to - * m_speciesComp[k * m_mm + i] - * The length of this vector is equal to m_kk * m_mm - */ - vector_fp m_speciesComp; + //! Atomic composition of the species. + /*! + * the number of atoms of i in species k is equal to + * m_speciesComp[k * m_mm + i] + * The length of this vector is equal to m_kk * m_mm + */ + vector_fp m_speciesComp; - /** - * m_speciesCharge: Vector of species charges - * length = m_kk - */ - vector_fp m_speciesCharge; + /** + * m_speciesCharge: Vector of species charges + * length = m_kk + */ + vector_fp m_speciesCharge; - /** - * m_speciesSize(): Vector of species sizes. - * length m_kk - * This is used in some equations of state - * which employ the constant partial molar - * volume approximation. It's so fundamental - * we've put it at the Constituents class level - */ - vector_fp m_speciesSize; + /** + * m_speciesSize(): Vector of species sizes. + * length m_kk + * This is used in some equations of state + * which employ the constant partial molar + * volume approximation. It's so fundamental + * we've put it at the Constituents class level + */ + vector_fp m_speciesSize; - private: + private: - }; + }; } // namespace diff --git a/Cantera/src/thermo/StoichSubstance.cpp b/Cantera/src/thermo/StoichSubstance.cpp index 339f935d5..7770c149a 100644 --- a/Cantera/src/thermo/StoichSubstance.cpp +++ b/Cantera/src/thermo/StoichSubstance.cpp @@ -272,6 +272,11 @@ namespace Cantera { er[0] = m_s0_R[0]; } + void StoichSubstance::getCp_R_ref(doublereal* cprt) const { + _updateThermo(); + cprt[0] = m_cp0_R[0]; + } + /* * */ diff --git a/Cantera/src/thermo/StoichSubstance.h b/Cantera/src/thermo/StoichSubstance.h index 022cf5e67..991ef4962 100644 --- a/Cantera/src/thermo/StoichSubstance.h +++ b/Cantera/src/thermo/StoichSubstance.h @@ -372,6 +372,16 @@ namespace Cantera { */ virtual void getEntropy_R_ref(doublereal *er) const; + //! 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; virtual void initThermo(); diff --git a/Cantera/src/thermo/SurfPhase.cpp b/Cantera/src/thermo/SurfPhase.cpp index 6585d2f42..ba72b7f2c 100644 --- a/Cantera/src/thermo/SurfPhase.cpp +++ b/Cantera/src/thermo/SurfPhase.cpp @@ -6,8 +6,12 @@ * \link Cantera::SurfPhase SurfPhase\endlink). */ -// Copyright 2002 California Institute of Technology +/* + * $Revision$ + * $Date$ + */ +// Copyright 2002 California Institute of Technology // turn off warnings under Windows #ifdef WIN32 @@ -17,14 +21,10 @@ #include "SurfPhase.h" #include "EdgePhase.h" -#include "utilities.h" -//#include "importCTML.h" #include "ThermoFactory.h" -#include using namespace std; - /////////////////////////////////////////////////////////// // // class SurfPhase methods @@ -60,7 +60,7 @@ namespace Cantera { XML_Node* xphase = get_XML_NameID("phase", std::string("#")+id, root); if (!xphase) { throw CanteraError("SurfPhase::SurfPhase", - "Couldn't find phase name in file:" + id); + "Couldn't find phase name in file:" + id); } // Check the model name to ensure we have compatibility const XML_Node& th = xphase->child("thermo"); @@ -153,23 +153,23 @@ namespace Cantera { return (ThermoPhase *) igp; } - doublereal SurfPhase:: - enthalpy_mole() const { - if (m_n0 <= 0.0) return 0.0; - _updateThermo(); - return mean_X(DATA_PTR(m_h0)); - } + doublereal SurfPhase:: + enthalpy_mole() const { + if (m_n0 <= 0.0) return 0.0; + _updateThermo(); + return mean_X(DATA_PTR(m_h0)); + } - SurfPhase:: - ~SurfPhase() { } + SurfPhase::~SurfPhase() { + } - /* - * For a surface phase, the pressure is not a relevant - * thermodynamic variable, and so the Enthalpy is equal to the - * internal energy. - */ - doublereal SurfPhase:: - intEnergy_mole() const { return enthalpy_mole(); } + /* + * For a surface phase, the pressure is not a relevant + * thermodynamic variable, and so the Enthalpy is equal to the + * internal energy. + */ + doublereal SurfPhase:: + intEnergy_mole() const { return enthalpy_mole(); } /* * Get the array of partial molar enthalpies of the species @@ -208,49 +208,41 @@ namespace Cantera { cpbar[k] *= GasConstant; } } - void SurfPhase::getPartialMolarVolumes(doublereal* vbar) const { getStandardVolumes(vbar); } - void SurfPhase:: - getStandardChemPotentials(doublereal* mu0) const { - _updateThermo(); - copy(m_mu0.begin(), m_mu0.end(), mu0); - } + void SurfPhase::getStandardChemPotentials(doublereal* mu0) const { + _updateThermo(); + copy(m_mu0.begin(), m_mu0.end(), mu0); + } - void SurfPhase:: - getChemPotentials(doublereal* mu) const { - _updateThermo(); - copy(m_mu0.begin(), m_mu0.end(), mu); - int k; - getActivityConcentrations(DATA_PTR(m_work)); - for (k = 0; k < m_kk; k++) { - mu[k] += GasConstant * temperature() * - (log(m_work[k]) - logStandardConc(k)); - } - } + void SurfPhase::getChemPotentials(doublereal* mu) const { + _updateThermo(); + copy(m_mu0.begin(), m_mu0.end(), mu); + int k; + getActivityConcentrations(DATA_PTR(m_work)); + for (k = 0; k < m_kk; k++) { + mu[k] += GasConstant * temperature() * + (log(m_work[k]) - logStandardConc(k)); + } + } - void SurfPhase:: - getActivityConcentrations(doublereal* c) const { - getConcentrations(c); - } + void SurfPhase::getActivityConcentrations(doublereal* c) const { + getConcentrations(c); + } - doublereal SurfPhase:: - standardConcentration(int k) const { - return m_n0/size(k); - } - - doublereal SurfPhase:: - logStandardConc(int k) const { - return m_logn0 - m_logsize[k]; - } + doublereal SurfPhase::standardConcentration(int k) const { + return m_n0/size(k); + } + doublereal SurfPhase::logStandardConc(int k) const { + return m_logn0 - m_logsize[k]; + } /// The only parameter that can be set is the site density. - void SurfPhase:: - setParameters(int n, doublereal* const c) { + void SurfPhase::setParameters(int n, doublereal* const c) { if (n != 1) { throw CanteraError("SurfPhase::setParameters", "Bad value for number of parameter"); @@ -263,59 +255,55 @@ namespace Cantera { m_logn0 = log(m_n0); } - void SurfPhase:: - getGibbs_RT(doublereal* grt) const { + void SurfPhase::getGibbs_RT(doublereal* grt) const { _updateThermo(); double rrt = 1.0/(GasConstant*temperature()); scale(m_mu0.begin(), m_mu0.end(), grt, rrt); } - void SurfPhase:: - getEnthalpy_RT(doublereal* hrt) const { - _updateThermo(); - double rrt = 1.0/(GasConstant*temperature()); - scale(m_h0.begin(), m_h0.end(), hrt, rrt); - } - - void SurfPhase:: - getEntropy_R(doublereal* sr) const { - _updateThermo(); - double rr = 1.0/GasConstant; - scale(m_s0.begin(), m_s0.end(), sr, rr); - } - void SurfPhase:: - getCp_R(doublereal* cpr) const { + getEnthalpy_RT(doublereal* hrt) const { + _updateThermo(); + double rrt = 1.0/(GasConstant*temperature()); + scale(m_h0.begin(), m_h0.end(), hrt, rrt); + } + + void SurfPhase::getEntropy_R(doublereal* sr) const { + _updateThermo(); + double rr = 1.0/GasConstant; + scale(m_s0.begin(), m_s0.end(), sr, rr); + } + + void SurfPhase::getCp_R(doublereal* cpr) const { _updateThermo(); double rr = 1.0/GasConstant; scale(m_cp0.begin(), m_cp0.end(), cpr, rr); } - void SurfPhase:: - getStandardVolumes(doublereal* vol) const { - _updateThermo(); - for (int k = 0; k < m_kk; k++) { - vol[k] = 1.0/standardConcentration(k); - } + void SurfPhase::getStandardVolumes(doublereal* vol) const { + _updateThermo(); + for (int k = 0; k < m_kk; k++) { + vol[k] = 1.0/standardConcentration(k); + } } - void SurfPhase:: - getGibbs_RT_ref(doublereal* grt) const { - getGibbs_RT(grt); - } + void SurfPhase::getGibbs_RT_ref(doublereal* grt) const { + getGibbs_RT(grt); + } - void SurfPhase:: - getEnthalpy_RT_ref(doublereal* hrt) const { - getEnthalpy_RT(hrt); - } + void SurfPhase::getEnthalpy_RT_ref(doublereal* hrt) const { + getEnthalpy_RT(hrt); + } - void SurfPhase:: - getEntropy_R_ref(doublereal* sr) const { - getEntropy_R(sr); - } + void SurfPhase::getEntropy_R_ref(doublereal* sr) const { + getEntropy_R(sr); + } - void SurfPhase:: - initThermo() { + void SurfPhase::getCp_R_ref(doublereal* cprt) const { + getCp_R(cprt); + } + + void SurfPhase::initThermo() { if (m_kk <= 0) { throw CanteraError("SurfPhase::initThermo", "Number of species is less than or equal to zero"); @@ -334,151 +322,148 @@ namespace Cantera { m_logsize[k] = log(size(k)); } - void SurfPhase:: - setPotentialEnergy(int k, doublereal pe) { - m_pe[k] = pe; - _updateThermo(true); + void SurfPhase::setPotentialEnergy(int k, doublereal pe) { + m_pe[k] = pe; + _updateThermo(true); + } + + void SurfPhase::setSiteDensity(doublereal n0) { + doublereal x = n0; + setParameters(1, &x); + } + + //void SurfPhase:: + //setElectricPotential(doublereal V) { + // for (int k = 0; k < m_kk; k++) { + // m_pe[k] = charge(k)*Faraday*V; + // } + // _updateThermo(true); + //} + + + /** + * Set the coverage fractions to a specified + * state. This routine converts to concentrations + * in kmol/m2, using m_n0, the surface site density, + * and size(k), which is defined to be the number of + * surface sites occupied by the kth molecule. + * It then calls State::setConcentrations to set the + * internal concentration in the object. + */ + void SurfPhase:: + setCoverages(const doublereal* theta) { + double sum = 0.0; + int k; + for (k = 0; k < m_kk; k++) { + sum += theta[k]; } - - void SurfPhase:: - setSiteDensity(doublereal n0) { - doublereal x = n0; - setParameters(1, &x); + if (sum <= 0.0) { + for (k = 0; k < m_kk; k++) { + cout << "theta(" << k << ") = " << theta[k] << endl; + } + throw CanteraError("SurfPhase::setCoverages", + "Sum of Coverage fractions is zero or negative"); } - - - //void SurfPhase:: - //setElectricPotential(doublereal V) { - // for (int k = 0; k < m_kk; k++) { - // m_pe[k] = charge(k)*Faraday*V; - // } - // _updateThermo(true); - //} - - - /** - * Set the coverage fractions to a specified - * state. This routine converts to concentrations - * in kmol/m2, using m_n0, the surface site density, - * and size(k), which is defined to be the number of - * surface sites occupied by the kth molecule. - * It then calls State::setConcentrations to set the - * internal concentration in the object. + for (k = 0; k < m_kk; k++) { + m_work[k] = m_n0*theta[k]/(sum*size(k)); + } + /* + * Call the State:: class function + * setConcentrations. */ - void SurfPhase:: - setCoverages(const doublereal* theta) { - double sum = 0.0; - int k; - for (k = 0; k < m_kk; k++) { - sum += theta[k]; - } - if (sum <= 0.0) { - for (k = 0; k < m_kk; k++) { - cout << "theta(" << k << ") = " << theta[k] << endl; - } - throw CanteraError("SurfPhase::setCoverages", - "Sum of Coverage fractions is zero or negative"); - } - for (k = 0; k < m_kk; k++) { - m_work[k] = m_n0*theta[k]/(sum*size(k)); - } - /* - * Call the State:: class function - * setConcentrations. - */ - setConcentrations(DATA_PTR(m_work)); + setConcentrations(DATA_PTR(m_work)); + } + + void SurfPhase:: + setCoveragesNoNorm(const doublereal* theta) { + for (int k = 0; k < m_kk; k++) { + m_work[k] = m_n0*theta[k]/(size(k)); + } + /* + * Call the State:: class function + * setConcentrations. + */ + setConcentrations(DATA_PTR(m_work)); + } + + void SurfPhase:: + getCoverages(doublereal* theta) const { + getConcentrations(theta); + for (int k = 0; k < m_kk; k++) { + theta[k] *= size(k)/m_n0; + } + } + + void SurfPhase:: + setCoveragesByName(std::string cov) { + int kk = nSpecies(); + int k; + compositionMap cc; + for (k = 0; k < kk; k++) { + cc[speciesName(k)] = -1.0; + } + parseCompString(cov, cc); + doublereal c; + vector_fp cv(kk, 0.0); + bool ifound = false; + for (k = 0; k < kk; k++) { + c = cc[speciesName(k)]; + if (c > 0.0) { + ifound = true; + cv[k] = c; + } + } + if (!ifound) { + throw CanteraError("SurfPhase::setCoveragesByName", + "Input coverages are all zero or negative"); + } + setCoverages(DATA_PTR(cv)); + } + + + void SurfPhase:: + _updateThermo(bool force) const { + doublereal tnow = temperature(); + if (m_tlast != tnow || force) { + m_spthermo->update(tnow, DATA_PTR(m_cp0), DATA_PTR(m_h0), + DATA_PTR(m_s0)); + m_tlast = tnow; + doublereal rt = GasConstant * tnow; + int k; + for (k = 0; k < m_kk; k++) { + m_h0[k] *= rt; + m_s0[k] *= GasConstant; + m_cp0[k] *= GasConstant; + m_mu0[k] = m_h0[k] - tnow*m_s0[k]; + } + m_tlast = tnow; + } + } + + void SurfPhase:: + setParametersFromXML(const XML_Node& eosdata) { + eosdata._require("model","Surface"); + doublereal n = getFloat(eosdata, "site_density", "toSI"); + if (n <= 0.0) + throw CanteraError("SurfPhase::setParametersFromXML", + "missing or negative site density"); + m_n0 = n; + m_logn0 = log(m_n0); + } + + + void SurfPhase::setStateFromXML(const XML_Node& state) { + + if (state.hasChild("temperature")) { + double t = getFloat(state, "temperature", "temperature"); + setTemperature(t); } - void SurfPhase:: - setCoveragesNoNorm(const doublereal* theta) { - for (int k = 0; k < m_kk; k++) { - m_work[k] = m_n0*theta[k]/(size(k)); - } - /* - * Call the State:: class function - * setConcentrations. - */ - setConcentrations(DATA_PTR(m_work)); - } - - void SurfPhase:: - getCoverages(doublereal* theta) const { - getConcentrations(theta); - for (int k = 0; k < m_kk; k++) { - theta[k] *= size(k)/m_n0; - } - } - - void SurfPhase:: - setCoveragesByName(std::string cov) { - int kk = nSpecies(); - int k; - compositionMap cc; - for (k = 0; k < kk; k++) { - cc[speciesName(k)] = -1.0; - } - parseCompString(cov, cc); - doublereal c; - vector_fp cv(kk, 0.0); - bool ifound = false; - for (k = 0; k < kk; k++) { - c = cc[speciesName(k)]; - if (c > 0.0) { - ifound = true; - cv[k] = c; - } - } - if (!ifound) { - throw CanteraError("SurfPhase::setCoveragesByName", - "Input coverages are all zero or negative"); - } - setCoverages(DATA_PTR(cv)); - } - - - void SurfPhase:: - _updateThermo(bool force) const { - doublereal tnow = temperature(); - if (m_tlast != tnow || force) { - m_spthermo->update(tnow, DATA_PTR(m_cp0), DATA_PTR(m_h0), - DATA_PTR(m_s0)); - m_tlast = tnow; - doublereal rt = GasConstant * tnow; - int k; - for (k = 0; k < m_kk; k++) { - m_h0[k] *= rt; - m_s0[k] *= GasConstant; - m_cp0[k] *= GasConstant; - m_mu0[k] = m_h0[k] - tnow*m_s0[k]; - } - m_tlast = tnow; - } - } - - void SurfPhase:: - setParametersFromXML(const XML_Node& eosdata) { - eosdata._require("model","Surface"); - doublereal n = getFloat(eosdata, "site_density", "toSI"); - if (n <= 0.0) - throw CanteraError("SurfPhase::setParametersFromXML", - "missing or negative site density"); - m_n0 = n; - m_logn0 = log(m_n0); - } - - - void SurfPhase::setStateFromXML(const XML_Node& state) { - - if (state.hasChild("temperature")) { - double t = getFloat(state, "temperature", "temperature"); - setTemperature(t); - } - - if (state.hasChild("coverages")) { - string comp = getChildValue(state,"coverages"); - setCoveragesByName(comp); - } + if (state.hasChild("coverages")) { + string comp = getChildValue(state,"coverages"); + setCoveragesByName(comp); } + } // Default constructor EdgePhase::EdgePhase(doublereal n0) : SurfPhase(n0) { @@ -521,16 +506,16 @@ namespace Cantera { return (ThermoPhase *) igp; } - void EdgePhase:: - setParametersFromXML(const XML_Node& eosdata) { - eosdata._require("model","Edge"); - doublereal n = getFloat(eosdata, "site_density", "toSI"); - if (n <= 0.0) - throw CanteraError("EdgePhase::setParametersFromXML", - "missing or negative site density"); - m_n0 = n; - m_logn0 = log(m_n0); - } + void EdgePhase:: + setParametersFromXML(const XML_Node& eosdata) { + eosdata._require("model","Edge"); + doublereal n = getFloat(eosdata, "site_density", "toSI"); + if (n <= 0.0) + throw CanteraError("EdgePhase::setParametersFromXML", + "missing or negative site density"); + m_n0 = n; + m_logn0 = log(m_n0); + } } diff --git a/Cantera/src/thermo/SurfPhase.h b/Cantera/src/thermo/SurfPhase.h index 29a0e9f30..2872c9db4 100644 --- a/Cantera/src/thermo/SurfPhase.h +++ b/Cantera/src/thermo/SurfPhase.h @@ -6,7 +6,7 @@ * (see \ref thermoprops and class \link Cantera::SurfPhase SurfPhase\endlink). */ -/* $Author$ +/* * $Date$ * $Revision$ * @@ -14,16 +14,13 @@ * */ - #ifndef CT_SURFPHASE_H #define CT_SURFPHASE_H #include "mix_defs.h" #include "ThermoPhase.h" - namespace Cantera { - //! A simple thermoydnamics model for a surface phase, //! assuming an ideal solution model. @@ -563,6 +560,18 @@ namespace Cantera { */ virtual void getEntropy_R_ref(doublereal *er) const; + //! 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; + + //------- new methods defined in this class ---------- //! Set the surface site fractions to a specified state.