diff --git a/include/cantera/thermo/ConstDensityThermo.h b/include/cantera/thermo/ConstDensityThermo.h index 2c59a9a0c..dcb3a5222 100644 --- a/include/cantera/thermo/ConstDensityThermo.h +++ b/include/cantera/thermo/ConstDensityThermo.h @@ -22,31 +22,15 @@ namespace Cantera //! Overloads the virtual methods of class ThermoPhase to implement the //! incompressible equation of state. /** - * - * - * Specification of Species Standard State Properties - * - * * Specification of Solution Thermodynamic Properties * - * The density is assumed to be constant, no matter what the concentration of the solution. - * - * - * Application within %Kinetics Managers - * - * - * XML Example - * - * An example of an XML Element named phase setting up a SurfPhase object named diamond_100 - * is given below. + * The density is assumed to be constant, no matter what the concentration of the solution. * * @ingroup thermoprops */ class ConstDensityThermo : public ThermoPhase { - public: - //! Constructor. ConstDensityThermo(); @@ -67,13 +51,13 @@ public: //! Duplication routine for objects which inherit from %ThermoPhase /*! - * This virtual routine can be used to duplicate %ThermoPhase objects - * inherited from %ThermoPhase even if the application only has + * This virtual routine can be used to duplicate objects + * derived from %ThermoPhase even if the application only has * a pointer to %ThermoPhase to work with. */ virtual ThermoPhase* duplMyselfAsThermoPhase() const; - //! overloaded methods of class ThermoPhase + //! Returns a constant corresponding to this class's equation of state virtual int eosType() const; /// Molar enthalpy. Units: J/kmol. @@ -95,26 +79,11 @@ public: virtual doublereal cv_mole() const; //! Return the thermodynamic pressure (Pa). - /*! - * This method must be overloaded in derived classes. Since the - * mass density, temperature, and mass fractions are stored, - * this method should use these values to implement the - * mechanical equation of state \f$ P(T, \rho, Y_1, \dots, - * Y_K) \f$. - */ virtual doublereal pressure() const; //! Set the internally stored 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); @@ -250,9 +219,6 @@ public: std::copy(_cpr.begin(), _cpr.end(), cpr); } - - // new methods defined here - //! Returns a reference to the vector of nondimensional //! enthalpies of the reference state at the current temperature //! of the solution and the reference pressure for the species. @@ -319,7 +285,6 @@ public: */ virtual void setToEquilState(const doublereal* lambda_RT); - //! Set the equation of state parameters /*! * @internal @@ -384,7 +349,6 @@ protected: doublereal m_press; private: - //! Function to update the reference state thermo functions void _updateThermo() const; }; diff --git a/include/cantera/thermo/IdealGasPhase.h b/include/cantera/thermo/IdealGasPhase.h index 629c3669e..961064569 100644 --- a/include/cantera/thermo/IdealGasPhase.h +++ b/include/cantera/thermo/IdealGasPhase.h @@ -97,10 +97,10 @@ namespace Cantera * The molar volume of a species is given by the ideal gas law * * \f[ - * V^o_k(T,P) = \frac{R T}{P} \mbox{\quad where} + * V^o_k(T,P) = \frac{R T}{P} * \f] * - * R is the molar gas constant. For a complete list of physical constants + * where R is the molar gas constant. For a complete list of physical constants * used within %Cantera, see \ref physConstants . * *
@@ -183,7 +183,7 @@ namespace Cantera * \f] * where * \f[ - * C_j^a = C^s a_j \mbox{\quad and \quad} C_k^a = C^s a_k + * C_j^a = C^s a_j \quad \mbox{and} \quad C_k^a = C^s a_k * \f] * * \f$ C_j^a \f$ is the activity concentration of species j, and @@ -228,7 +228,8 @@ namespace Cantera * For completeness, the pressure equilibrium constant may be obtained as well * * \f[ - * \frac{P_j P_k}{ P_l P_{ref}} = K_p^1 = \exp(\frac{\mu^{ref}_l - \mu^{ref}_j - \mu^{ref}_k}{R T} ) + * \frac{P_j P_k}{ P_l P_{ref}} = K_p^1 = + \exp\left(\frac{\mu^{ref}_l - \mu^{ref}_j - \mu^{ref}_k}{R T} \right) * \f] * * \f$ K_p \f$ is the simplest form of the equilibrium constant for ideal gases. However, it isn't @@ -255,7 +256,6 @@ namespace Cantera *

Instantiation of the Class

*
* - * * The constructor for this phase is located in the default ThermoFactory * for %Cantera. A new %IdealGasPhase may be created by the following code * snippet: @@ -281,20 +281,20 @@ namespace Cantera * An example of an XML Element named phase setting up a IdealGasPhase * object named silane is given below. * - * @verbatim - - - Si H He - - H2 H HE SIH4 SI SIH SIH2 SIH3 H3SISIH SI2H6 - H2SISIH2 SI3H8 SI2 SI3 - - - - - - - @endverbatim + * @code + * + * + * Si H He + * + * H2 H HE SIH4 SI SIH SIH2 SIH3 H3SISIH SI2H6 + * H2SISIH2 SI3H8 SI2 SI3 + * + * + * + * + * + * + * @endcode * * The model attribute "IdealGas" of the thermo XML element identifies the phase as * being of the type handled by the IdealGasPhase object. @@ -304,22 +304,13 @@ namespace Cantera */ class IdealGasPhase: public ThermoPhase { - public: - //! Default empty Constructor IdealGasPhase(); //! Construct and initialize an IdealGasPhase ThermoPhase object //! directly from an ASCII input file /*! - * Working constructors - * - * The two constructors below are a direct way that - * the phase can initialize itself. They are shells that call - * the routine initThermo(), with a reference to the - * XML database to get the info for the phase. - * * @param inputFile Name of the input file containing the phase XML data * to set up the object * @param id ID of the phase in the input file. Defaults to the @@ -360,10 +351,10 @@ public: //! Destructor virtual ~IdealGasPhase(); - //! Duplicator from the %ThermoPhase parent class + //! Duplicator from the ThermoPhase parent class /*! - * Given a pointer to a %ThermoPhase object, this function will - * duplicate the %ThermoPhase object and all underlying structures. + * Given a pointer to a ThermoPhase object, this function will + * duplicate the ThermoPhase object and all underlying structures. * This is basically a wrapper around the inherited copy constructor. * * @return returns a pointer to a %ThermoPhase object, containing @@ -379,10 +370,8 @@ public: return cIdealGas; } - /** - * @name Molar Thermodynamic Properties of the Solution ------------------------------ - * @{ - */ + //! @name Molar Thermodynamic Properties of the Solution + //! @{ //! Return the Molar enthalpy. Units: J/kmol. /*! @@ -496,19 +485,25 @@ public: /** * @returns species vibrational specific heat at - * constant volume. This is defined as + * constant volume, * \f[ - * C^{vib}_{v,s} \equiv \frac{\partial e^{vib}_s}{\partial T_V} = \frac{R_s \theta_{vs}^2 \exp\left(\theta_{vs}/T_V\right)}{\left[\left(\exp\left(\theta_{vs}/T_V\right)-1\right)T_V\right]^2} + * C^{vib}_{v,s} = \frac{\partial e^{vib}_{v,s} }{\partial T} * \f] + * where the species vibration energy \f$ e^{vib}_{v,s} \f$ is + * - atom: + * 0 + * - Diatomic: + * \f[ \frac{R_s \theta_{v,s}}{e^{\theta_{v,s}/T}-1} \f] + * - General Molecule: + * \f[ + * \sum_i \frac{R_s \theta_{v,s,i}}{e^{\theta_{v,s,i}/T}-1} + * \f] */ virtual doublereal cv_vib(int k, doublereal T) const; - //@} - - /** - * @name Mechanical Equation of State ------------------------------------------------ - * @{ - */ + //! @} + //! @name Mechanical Equation of State + //! @{ /** * Pressure. Units: Pa. @@ -560,8 +555,7 @@ public: //@} /** - * @name Chemical Potentials and Activities ------------------------------------------ - * + * @name Chemical Potentials and Activities * * The activity \f$a_k\f$ of a species in solution is * related to the chemical potential by @@ -639,7 +633,7 @@ public: virtual void getActivityCoefficients(doublereal* ac) const; //@} - /// @name Partial Molar Properties of the Solution ---------------------------------- + /// @name Partial Molar Properties of the Solution //@{ //! Get the species chemical potentials. Units: J/kmol. @@ -689,7 +683,7 @@ public: virtual void getPartialMolarVolumes(doublereal* vbar) const; //@} - /// @name Properties of the Standard State of the Species in the Solution ---------- + /// @name Properties of the Standard State of the Species in the Solution //@{ //! Get the array of chemical potentials at unit activity for the @@ -765,7 +759,7 @@ public: virtual void getStandardVolumes(doublereal* vol) const; //@} - /// @name Thermodynamic Values for the Species Reference States --------------------- + /// @name Thermodynamic Values for the Species Reference States //@{ //! Returns the vector of nondimensional @@ -906,11 +900,9 @@ public: */ virtual void initThermo(); - //!This method is used by the ChemEquil equilibrium solver. + //! Method used by the ChemEquil equilibrium solver. /*! * @internal - * @name Chemical Equilibrium - * @{ * * Set mixture to an equilibrium state consistent with specified * element potentials and temperature. @@ -926,8 +918,6 @@ public: */ virtual void setToEquilState(const doublereal* lambda_RT); - //@} - protected: //! Reference state pressure /*! @@ -960,15 +950,19 @@ protected: mutable vector_fp m_pp; private: - //! Update the species reference state thermodynamic functions /*! - * The polynomials for the standard state functions are only - * reevaluated if the temperature has changed. - * + * This method is called each time a thermodynamic property is requested, + * to check whether the internal species properties within the object + * need to be updated. Currently, this updates the species thermo + * polynomial values for the current value of the temperature. A check is + * made to see if the temperature has changed since the last evaluation. + * This object does not contain any persistent data that depends on the + * concentration, that needs to be updated. The state object modifies its + * concentration dependent information at the time the setMoleFractions() + * (or equivalent) call is made. */ void _updateThermo() const; - }; } diff --git a/include/cantera/thermo/IdealSolidSolnPhase.h b/include/cantera/thermo/IdealSolidSolnPhase.h index 4f24aff6f..c9d59bc5c 100644 --- a/include/cantera/thermo/IdealSolidSolnPhase.h +++ b/include/cantera/thermo/IdealSolidSolnPhase.h @@ -22,7 +22,6 @@ #include "ThermoFactory.h" #include "SpeciesThermo.h" - namespace Cantera { @@ -45,17 +44,10 @@ const int cIdealSolidSolnPhase2 = 5012; * The class derives from class ThermoPhase, * and overloads the virtual methods defined there with ones that * use expressions appropriate for ideal solution mixtures. - * File name for the XML datafile containing information - * for this phase + * * The generalized concentrations can have three different forms - * depending on the value of the member attribute m_formGC, which + * depending on the value of the member attribute #m_formGC, which * is supplied in the constructor and in the XML file. - * - * - * - * - * - *
m_formGC GeneralizedConc StandardConc
0 X_k 1.0
1 X_k / V_k 1.0 / V_k
2 X_k / V_N 1.0 / V_N
* The value and form of the generalized concentration will affect * reaction rate constants involving species in this phase. * @@ -63,49 +55,30 @@ const int cIdealSolidSolnPhase2 = 5012; */ class IdealSolidSolnPhase : public ThermoPhase { - public: - /** * Constructor for IdealSolidSolnPhase. * The generalized concentrations can have three different forms - * depending on the value of the member attribute m_formGC, which + * depending on the value of the member attribute #m_formGC, which * is supplied in the constructor or read from the xml data file. - * - * - * - * - * - *
m_formGC GeneralizedConc StandardConc
0 X_k 1.0
1 X_k / V_k 1.0 / V_k
2 X_k / V_N 1.0 / V_N
* - * @param formCG This parameter initializes the m_formGC variable. The default - * is a value of 0. + * @param formCG This parameter initializes the #m_formGC variable. */ IdealSolidSolnPhase(int formCG=0); - //! Construct and initialize an IdealSolidSolnPhase ThermoPhase object //! directly from an ASCII input file /*! - * * This constructor will also fully initialize the object. * The generalized concentrations can have three different forms - * depending on the value of the member attribute m_formGC, which + * depending on the value of the member attribute #m_formGC, which * is supplied in the constructor or read from the xml data file. * - * - * - * - * - * - *
m_formGC GeneralizedConc StandardConc
0 X_k 1.0
1 X_k / V_k 1.0 / V_k
2 X_k / V_N 1.0 / V_N
- * * @param infile File name for the XML datafile containing information * for this phase * @param id The name of this phase. This is used to look up * the phase in the XML datafile. - * @param formCG This parameter initializes the m_formGC variable. The default - * is a value of 0. + * @param formCG This parameter initializes the #m_formGC variable. */ IdealSolidSolnPhase(const std::string& infile, const std::string& id="", int formCG=0); @@ -113,42 +86,28 @@ public: //! directly from an XML database /*! * The generalized concentrations can have three different forms - * depending on the value of the member attribute m_formGC, which + * depending on the value of the member attribute #m_formGC, which * is supplied in the constructor and/or read from the data file. * - * - * - * - * - * - *
m_formGC GeneralizedConc StandardConc
0 X_k 1.0
1 X_k / V_k 1.0 / V_k
2 X_k / V_N 1.0 / V_N
- * * @param root XML tree containing a description of the phase. * The tree must be positioned at the XML element * named phase with id, "id", on input to this routine. * @param id The name of this phase. This is used to look up * the phase in the XML datafile. - * @param formCG This parameter initializes the m_formGC variable. The default - * is a value of 0. + * @param formCG This parameter initializes the #m_formGC variable. */ IdealSolidSolnPhase(XML_Node& root, const std::string& id="", int formCG=0); - /*! - * Copy Constructor - */ + //! Copy Constructor IdealSolidSolnPhase(const IdealSolidSolnPhase&); - /*! - * Assignment operator - */ + //! Assignment operator IdealSolidSolnPhase& operator=(const IdealSolidSolnPhase&); /*! * Base Class Duplication Function - * -> given a pointer to ThermoPhase, this function can - * duplicate the object. (note has to be a separate function - * not the copy constructor, because it has to be - * a virtual function) + * + * Given a pointer to ThermoPhase, this function can duplicate the object. */ virtual ThermoPhase* duplMyselfAsThermoPhase() const; @@ -157,14 +116,12 @@ public: /** * Equation of state flag. Returns a value depending upon the value of - * m_formGC, which is defined at instantiation. + * #m_formGC, which is defined at instantiation. */ virtual int eosType() const; - /** - * @name Molar Thermodynamic Properties of the Solution ------------------------ - * @{ - */ + //! @name Molar Thermodynamic Properties of the Solution + //! @{ /** * Molar enthalpy of the solution. Units: J/kmol. @@ -257,7 +214,7 @@ public: } //@} - /** @name Mechanical Equation of State Properties ------------------------------------ + /** @name Mechanical Equation of State Properties * * In this equation of state implementation, the density is a * function only of the mole fractions. Therefore, it can't be @@ -305,9 +262,6 @@ public: * species molar volumes. We have additionally specified * in this class that the pure species molar volumes are * independent of temperature and pressure. - * - * NOTE: This is a non-virtual function, which is not a - * member of the ThermoPhase base class. */ void calcDensity(); @@ -322,11 +276,6 @@ public: * to create a condition where the density is a function of * the pressure. * - * This function will now throw an error condition. - * - * NOTE: This is a virtual function that overwrites the State.h - * class - * * @param rho Input density */ virtual void setDensity(const doublereal rho); @@ -337,9 +286,6 @@ public: * * This function will now throw an error condition. * - * NOTE: This is virtual function that overwrites the State.h - * class - * * @param rho Input Density */ virtual void setMolarDensity(const doublereal rho); @@ -379,11 +325,10 @@ public: */ virtual void setConcentrations(const doublereal* const c); - //@} /** - * @name Chemical Potentials and Activities ----------------------------------------- + * @name Chemical Potentials and Activities * * The activity \f$a_k\f$ of a species in solution is * related to the chemical potential by @@ -443,14 +388,7 @@ public: * standard concentration to be independent of the mole fractions. * * In this implementation the form of the generalized concentrations - * depend upon the member attribute, m_formGC: - * - * - * - * - * - * - *
m_formGC GeneralizedConc StandardConc
0 X_k 1.0
1 X_k / V_k 1.0 / V_k
2 X_k / V_N 1.0 / V_N
+ * depend upon the member attribute, #m_formGC. * * HKM Note: We have absorbed the pressure dependence of the pure species * state into the thermodynamics functions. Therefore the @@ -464,14 +402,12 @@ public: virtual void getActivityConcentrations(doublereal* c) const; /** - * The standard concentration \f$ C^0_k \f$ used to normalize - * the generalized concentration. - * In many cases, this quantity - * will be the same for all species in a phase. - * However, for this case, we will return a distinct concentration - * for each species. This is the inverse of the species molar - * volume. Units for the standard concentration are - * kmol m-3. + * The standard concentration \f$ C^0_k \f$ used to normalize the + * generalized concentration. In many cases, this quantity will be the + * same for all species in a phase. However, for this case, we will return + * a distinct concentration for each species. This is the inverse of the + * species molar volume. Units for the standard concentration are kmol + * m-3. * * @param k Species number: this is a require parameter, * a change from the ThermoPhase base class, where it was @@ -511,14 +447,16 @@ public: * units are needed. Usually, MKS units are assumed throughout * the program and in the XML input files. * - * @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 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. @@ -531,7 +469,6 @@ public: virtual void getUnitsStandardConc(double* uA, int k = 0, int sizeUA = 6) const; - //! Get the array of species activity coefficients /*! * @param ac output vector of activity coefficients. Length: m_kk @@ -572,10 +509,9 @@ public: virtual void getChemPotentials_RT(doublereal* mu) const; //@} - /// @name Partial Molar Properties of the Solution ----------------------------- + /// @name Partial Molar Properties of the Solution //@{ - //! Returns an array of partial molar enthalpies for the species in the mixture. /*! * Units (J/kmol) @@ -637,10 +573,9 @@ public: virtual void getPartialMolarVolumes(doublereal* vbar) const; //@} - /// @name Properties of the Standard State of the Species in the Solution ------------------------------------- + /// @name Properties of the Standard State of the Species in the Solution //@{ - /** * Get the standard state chemical potentials of the species. * This is the array of chemical potentials at unit activity @@ -660,7 +595,6 @@ public: getPureGibbs(mu0); } - //! Get the array of nondimensional Enthalpy functions for the standard state species //! at the current T and P of the solution. /*! @@ -679,7 +613,6 @@ public: */ void getEnthalpy_RT(doublereal* hrt) const; - //! Get the nondimensional Entropies for the species //! standard states at the current T and P of the solution. /*! @@ -727,21 +660,18 @@ public: */ virtual void getPureGibbs(doublereal* gpure) const; - //! Returns the vector of nondimensional //! internal Energies of the standard state at the current //! temperature and pressure of the solution for each species. /*! - * * @param urt Output vector of standard state nondimensional internal energies. * Length: m_kk. */ virtual void getIntEnergy_RT(doublereal* urt) const; /** - * Get the nondimensional heat capacity at constant pressure - * function for the species - * standard states at the current T and P of the solution. + * Get the nondimensional heat capacity at constant pressure function for + * the species standard states at the current T and P of the solution. * \f[ * Cp^0_k(T,P) = Cp^{ref}_k(T) * \f] @@ -757,8 +687,7 @@ public: /** * Get the molar volumes of each species in their standard - * states at the current - * T and P of the solution. + * states at the current T and P of the solution. * units = m^3 / kmol * * @param vol Output vector of standard state volumes. @@ -766,12 +695,10 @@ public: */ virtual void getStandardVolumes(doublereal* vol) const; - //@} - /// @name Thermodynamic Values for the Species Reference States ------ + /// @name Thermodynamic Values for the Species Reference States //@{ - /** * Returns the vector of nondimensional * enthalpies of the reference state at the current temperature @@ -884,8 +811,9 @@ public: virtual doublereal potentialEnergy(int k) const { return m_pe[k]; } + //@} - /// @name Utility Functions ----------------------------------------------- + /// @name Utility Functions //@{ /** @@ -896,45 +824,34 @@ public: virtual void initThermo(); /** - * @internal - * Import and initialize a ThermoPhase object - * using an XML tree. - * Here we read extra information about the XML description - * of a phase. Regular information about elements and species - * and their reference state thermodynamic information - * have already been read at this point. - * For example, we do not need to call this function for - * ideal gas equations of state. - * This function is called from importPhase() - * after the elements and the - * species are initialized with default ideal solution - * level data. + * @internal Import and initialize a ThermoPhase object using an XML + * tree. Here we read extra information about the XML description of a + * phase. Regular information about elements and species and their + * reference state thermodynamic information have already been read at + * this point. For example, we do not need to call this function for + * ideal gas equations of state. This function is called from + * importPhase() after the elements and the species are initialized + * with default ideal solution level data. * - * @param phaseNode This object must be the phase node of a - * complete XML tree - * description of the phase, including all of the - * species data. In other words while "phase" must - * point to an XML phase object, it must have - * sibling nodes "speciesData" that describe - * the species in the phase. - * @param id ID of the phase. If nonnull, a check is done - * to see if phaseNode is pointing to the phase - * with the correct id. + * @param phaseNode This object must be the phase node of a complete XML + * tree description of the phase, including all of the + * species data. In other words while "phase" must point to + * an XML phase object, it must have sibling nodes + * "speciesData" that describe the species in the phase. + * @param id ID of the phase. If nonnull, a check is done to see if + * phaseNode is pointing to the phase with the correct id. */ virtual void initThermoXML(XML_Node& phaseNode, const std::string& id); - /** * Set mixture to an equilibrium state consistent with specified * element potentials and the temperature. * * @param lambda_RT vector of non-dimensional element potentials * \f$ \lambda_m/RT \f$. - * */ virtual void setToEquilState(const doublereal* lambda_RT); - /** * Report the molar volume of species k * @@ -957,12 +874,18 @@ public: //@} protected: - /** - * Format for the generalized concentrations - * 0 = C_k = X_k. (default) - * 1 = C_k = X_k / V_k - * 2 = C_k = X_k / V_N + * Format for the generalized concentrations. + * + * + * + * + * + * + *
m_formGC GeneralizedConc StandardConc
0 (default) X_k 1.0
1 X_k / V_k 1.0 / V_k
2 X_k / V_N 1.0 / V_N
+ * + * The value and form of the generalized concentration will affect + * reaction rate constants involving species in this phase. */ int m_formGC; @@ -995,9 +918,7 @@ protected: */ mutable doublereal m_tlast; - /** - * Vector containing the species reference enthalpies at T = m_tlast - */ + //! Vector containing the species reference enthalpies at T = m_tlast mutable vector_fp m_h0_RT; /** @@ -1006,16 +927,10 @@ protected: */ mutable vector_fp m_cp0_R; - /** - * Vector containing the species reference Gibbs functions - * at T = m_tlast - */ + //! Vector containing the species reference Gibbs functions at T = m_tlast mutable vector_fp m_g0_RT; - /** - * Vector containing the species reference entropies - * at T = m_tlast - */ + //! Vector containing the species reference entropies at T = m_tlast mutable vector_fp m_s0_R; /** @@ -1024,18 +939,14 @@ protected: */ mutable vector_fp m_expg0_RT; - /** - * Vector of potential energies for the species. - */ + //! Vector of potential energies for the species. mutable vector_fp m_pe; - /** - * Temporary array used in equilibrium calculations - */ + //! Temporary array used in equilibrium calculations mutable vector_fp m_pp; private: - /// @name Utility Functions ------------------------------------------ + /// @name Utility Functions //@{ /** * This function gets called for every call to functions in this @@ -1047,9 +958,7 @@ private: */ void _updateThermo() const; - /** - * This internal function adjusts the lengths of arrays - */ + //! This internal function adjusts the lengths of arrays void initLengths(); //@} @@ -1057,8 +966,3 @@ private: } #endif - - - - - diff --git a/include/cantera/thermo/LatticePhase.h b/include/cantera/thermo/LatticePhase.h index 9bff67b83..93cb9cc0e 100644 --- a/include/cantera/thermo/LatticePhase.h +++ b/include/cantera/thermo/LatticePhase.h @@ -24,17 +24,16 @@ namespace Cantera //! A simple thermodynamic model for a bulk phase, //! assuming a lattice of solid atoms /*! - * The bulk consists of a matrix of equivalent sites whose molar density - * does not vary with temperature or pressure. The thermodynamics - * obeys the ideal solution laws. The phase and the pure species phases which + * The bulk consists of a matrix of equivalent sites whose molar density + * does not vary with temperature or pressure. The thermodynamics + * obeys the ideal solution laws. The phase and the pure species phases which * comprise the standard states of the species are assumed to have * zero volume expansivity and zero isothermal compressibility. * * The density of matrix sites is given by the variable \f$ C_o \f$, * which has SI units of kmol m-3. * - * - * Specification of Species Standard %State Properties + * 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 @@ -48,7 +47,7 @@ namespace Cantera * which has a weak dependence on the system pressure, \f$P\f$. * * \f[ - * \raggedright h^o_k(T,P) = + * h^o_k(T,P) = * h^{ref}_k(T) + \left( \frac{P - P_{ref}}{C_o} \right) * \f] * @@ -73,7 +72,6 @@ namespace Cantera * V^o_k(T,P) = \frac{1.0}{C_o} * \f] * - * *
*

Specification of Solution Thermodynamic Properties

*
@@ -179,24 +177,23 @@ namespace Cantera * \mu_l(T,P) = \mu^o_l(T, P) + R T \log(a_l) * \f] * - * The concentration equilibrium constant, \f$ K_c \f$, may be obtained by changing over - * to activity concentrations. When this is done: + * The concentration equilibrium constant, \f$ K_c \f$, may be obtained by changing over + * to activity concentrations. When this is done: * - * \f[ - * \frac{C^a_j C^a_k}{ C^a_l} = C^o K_a^{o,1} = K_c^1 = - * \exp(\frac{\mu^{o}_l - \mu^{o}_j - \mu^{o}_k}{R T} ) - * \f] + * \f[ + * \frac{C^a_j C^a_k}{ C^a_l} = C^o K_a^{o,1} = K_c^1 = + * \exp(\frac{\mu^{o}_l - \mu^{o}_j - \mu^{o}_k}{R T} ) + * \f] * * - * %Kinetics managers will calculate the concentration equilibrium constant, \f$ K_c \f$, - * using the second and third part of the above expression as a definition for the concentration - * equilibrium constant. + * %Kinetics managers will calculate the concentration equilibrium constant, \f$ K_c \f$, + * using the second and third part of the above expression as a definition for the concentration + * equilibrium constant. * *
*

Instantiation of the Class

*
* - * * The constructor for this phase is located in the default ThermoFactory * for %Cantera. A new %LatticePhase object may be created by the following code snippet: * @@ -224,34 +221,31 @@ namespace Cantera * An example of an XML Element named phase setting up a LatticePhase object named "O_lattice_SiO2" * is given below. * - * @verbatim - - - Si H He - - O_O Vac_O - - - - 73.159 - Vac_O - - - - - @endverbatim + * @code + * + * + * Si H He + * + * O_O Vac_O + * + * + * + * 73.159 + * Vac_O + * + * + * + * + * @endcode * * The model attribute "Lattice" of the thermo XML element identifies the phase as * being of the type handled by the LatticePhase object. * * @ingroup thermoprops - * */ class LatticePhase : public ThermoPhase { - public: - //! Base Empty constructor LatticePhase(); @@ -289,8 +283,6 @@ public: * This virtual function is used to create a duplicate of the * current phase. It's used to duplicate the phase when given * a ThermoPhase pointer to the phase. - * - * @return It returns a ThermoPhase pointer. */ ThermoPhase* duplMyselfAsThermoPhase() const; @@ -299,10 +291,8 @@ public: return cLattice; } - /** - * @name Molar Thermodynamic Properties of the Solution ------------------------ - * @{ - */ + //! @name Molar Thermodynamic Properties of the Solution + //! @{ //! Return the Molar Enthalpy. Units: J/kmol. /*! @@ -405,8 +395,7 @@ public: virtual doublereal cv_mole() const; //@} - /// @name Mechanical Equation of State Properties ------------------------------------ - //@{ + /// @name Mechanical Equation of State Properties /** * In this equation of state implementation, the density is a * function only of the mole fractions. Therefore, it can't be @@ -417,7 +406,6 @@ public: */ //@{ - //! Pressure. Units: Pa. /*! * For this incompressible system, we return the internally stored @@ -455,9 +443,6 @@ public: * species molar volumes. We have additionally specified * in this class that the pure species molar volumes are * independent of temperature and pressure. - * - * NOTE: This is a non-virtual function, which is not a - * member of the ThermoPhase base class. */ doublereal calcDensity(); @@ -496,11 +481,9 @@ public: */ virtual void setConcentrations(const doublereal* const c); - //@} /// @name Activities, Standard States, and Activity Concentrations /** - * * The activity \f$a_k\f$ of a species in solution is * related to the chemical potential by \f[ \mu_k = \mu_k^0(T) * + \hat R T \log a_k. \f] The quantity \f$\mu_k^0(T,P)\f$ is @@ -564,7 +547,6 @@ public: //@} /// @name Partial Molar Properties of the Solution - /// //@{ //! Get the species chemical potentials. Units: J/kmol. @@ -578,12 +560,10 @@ public: */ virtual void getChemPotentials(doublereal* mu) const; - //@} - /// @name Partial Molar Properties of the Solution ----------------------------- + /// @name Partial Molar Properties of the Solution //@{ - /** * Returns an array of partial molar enthalpies for the species * in the mixture. @@ -642,7 +622,6 @@ public: */ virtual void getPartialMolarVolumes(doublereal* vbar) const; - //! Get the array of chemical potentials at unit activity for the //! species standard states at the current T and P of the solution. /*! @@ -664,7 +643,6 @@ public: */ virtual void getPureGibbs(doublereal* gpure) const; - //@} /// @name Properties of the Standard State of the Species in the Solution //@{ @@ -844,7 +822,6 @@ public: */ virtual void initThermo(); - //! Import and initialize a ThermoPhase object using an XML tree. /*! * Here we read extra information about the XML description @@ -912,17 +889,16 @@ public: * * eosdata points to the thermo block, and looks like this: * - * @verbatim - - - 73.159 - "O_vacancy" - - @endverbatim - * + * @code + * + * + * 73.159 + * "O_vacancy" + * + * + * @endcode */ virtual void setParametersFromXML(const XML_Node& eosdata); - //@} protected: @@ -953,7 +929,6 @@ protected: //! Temporary storage for the reference state entropies at the current temperature mutable vector_fp m_s0_R; - //! String name for the species which represents a vacancy //! in the lattice /*! @@ -975,10 +950,7 @@ protected: */ doublereal m_site_density; - // doublereal m_molar_lattice_volume; - private: - //! Update the species reference state thermodynamic functions /*! * The polynomials for the standard state functions are only diff --git a/include/cantera/thermo/LatticeSolidPhase.h b/include/cantera/thermo/LatticeSolidPhase.h index 915e9013d..387308e87 100644 --- a/include/cantera/thermo/LatticeSolidPhase.h +++ b/include/cantera/thermo/LatticeSolidPhase.h @@ -32,7 +32,6 @@ namespace Cantera * lattice phase and the molar densities of the sublattice and the molar density within the LatticeSolidPhase * have the same values. * - * * The mole fraction vector is redefined witin the the LatticeSolidPhase object. Each of the mole * fractions sum to one on each of the sublattices. The routine getMoleFraction() and setMoleFraction() * have been redefined to use this convention. @@ -65,7 +64,6 @@ namespace Cantera * solid. \f$ \theta_i \f$ is a fixed weighting factor for the ith lattice representing the lattice * stoichiometric coefficient. For this object the \f$ \theta_i \f$ values are fixed. * - * * Let's take FeS2 as an example, which may be thought of as a combination of two lattices: Fe and S lattice. * The Fe sublattice has a molar density of 1 gmol cm-3. The S sublattice has a molar density of 2 gmol cm-3. * We then define the LatticeSolidPhase object as having a nominal composition of FeS2, and having a @@ -73,7 +71,6 @@ namespace Cantera * associated with the sublattices. The Fe sublattice will have a weight of 1.0 associated with it. The * S sublattice will have a weight of 2.0 associated with it. * - * *
*

Specification of Solution Density Properties

*
@@ -108,9 +105,7 @@ namespace Cantera */ class LatticeSolidPhase : public ThermoPhase { - public: - //! Base empty constructor LatticeSolidPhase(); @@ -141,7 +136,7 @@ public: //! Equation of state type flag. /*! - * Redefine this to return cLatticeSolid, listed in mix_defs.h. + * Returns cLatticeSolid, listed in mix_defs.h. */ virtual int eosType() const { return cLatticeSolid; @@ -175,7 +170,6 @@ public: */ virtual doublereal maxTemp(size_t k = npos) const; - //! Returns the reference pressure in Pa. This function is a wrapper //! that calls the species thermo refPressure function. virtual doublereal refPressure() const ; @@ -205,7 +199,6 @@ public: */ virtual doublereal enthalpy_mole() const; - //! Return the Molar Internal Energy. Units: J/kmol. /*! * The molar enthalpy is determined by the following formula, where \f$ \theta_n \f$ is the @@ -309,8 +302,6 @@ public: * \f] * * where \f$ \rho_n \f$ is the density of the nth sublattice - * - * Note this is a nonvirtual function. */ doublereal calcDensity(); @@ -370,7 +361,6 @@ public: return err("not implemented"); } - //! Set the mass fractions to the specified values, and then //! normalize them so that they sum to 1.0. /*! @@ -386,7 +376,6 @@ public: err("not implemented"); } - //! Set the mass fractions to the specified values without normalizing. /*! * This is useful when the normalization @@ -452,7 +441,6 @@ public: */ virtual void getChemPotentials(doublereal* mu) const; - //! Returns an array of partial molar enthalpies for the species in the mixture. /*! * Units (J/kmol) @@ -520,7 +508,6 @@ public: * \f$. The values are evaluated at the current * temperature and pressure of the solution. * - * * This returns the underlying lattice standard chemical potentials, as the units are kmol-1 of * the sublattice species. * @@ -555,7 +542,7 @@ public: */ virtual doublereal logStandardConc(size_t k=0) const; //@} - /// @name Thermodynamic Values for the Species Reference States -------------------- + /// @name Thermodynamic Values for the Species Reference States //@{ //! Returns the vector of nondimensional enthalpies of the reference state at the current @@ -571,7 +558,6 @@ public: */ virtual void getGibbs_RT_ref(doublereal* grt) const; - //! Returns the vector of the gibbs function of the reference state at the current //! temperatureof the solution and the reference pressure for the species. /*! @@ -586,7 +572,6 @@ public: */ virtual void getGibbs_ref(doublereal* g) const; - //! Initialize the ThermoPhase object after all species have been set up /*! * @internal Initialize. @@ -617,10 +602,8 @@ public: */ virtual void installSlavePhases(Cantera::XML_Node* phaseNode); - //! 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 @@ -633,10 +616,8 @@ public: */ virtual void setParametersFromXML(const XML_Node& eosdata); - //! Set the Lattice mole fractions using a string /*! - * * @param n Integer value of the lattice whose mole fractions are being set * @param x string containing Name:value pairs that will specify the mole fractions * of species on a particular lattice @@ -697,7 +678,6 @@ protected: std::vector lkstart_; private: - //! Update the reference thermodynamic functions void _updateThermo() const; }; diff --git a/include/cantera/thermo/MixtureFugacityTP.h b/include/cantera/thermo/MixtureFugacityTP.h index f838de033..8dc581061 100644 --- a/include/cantera/thermo/MixtureFugacityTP.h +++ b/include/cantera/thermo/MixtureFugacityTP.h @@ -3,7 +3,6 @@ * Header file for a derived class of ThermoPhase that handles * non-ideal mixtures based on the fugacity models (see \ref thermoprops and * class \link Cantera::MixtureFugacityTP MixtureFugacityTP\endlink). - * */ /* * Copyright (2005) Sandia Corporation. Under the terms of @@ -72,19 +71,14 @@ class PDSS; * * This class is usually used for non-ideal gases. * - * * @nosubgrouping */ class MixtureFugacityTP : public ThermoPhase { - public: + //! @name Constructors and Duplicators for %MixtureFugacityTP + //! @{ - /*! - * - * @name Constructors and Duplicators for %MixtureFugacityTP - * - */ //! Constructor. MixtureFugacityTP(); @@ -103,19 +97,15 @@ public: //! Destructor. virtual ~MixtureFugacityTP(); - //! Duplication routine /*! * @return Returns a duplication */ virtual ThermoPhase* duplMyselfAsThermoPhase() const; - //@} - - /** - * @name Utilities (MixtureFugacityTP) - */ - //@{ + //! @} + //! @name Utilities + //! @{ /** * Equation of state type flag. The base class returns * zero. Subclasses should define this to return a unique @@ -131,12 +121,10 @@ public: //! temperature based, and variable pressure based. /*! * Currently, there are two standard state conventions: - * - Temperature-based activities - * cSS_CONVENTION_TEMPERATURE 0 - * - default - * - * - Variable Pressure and Temperature -based activities - * cSS_CONVENTION_VPSS 1 + * - Temperature-based activities, + * `cSS_CONVENTION_TEMPERATURE 0` (default) + * - Variable Pressure and Temperature based activities, + * `cSS_CONVENTION_VPSS 1` */ virtual int standardStateConvention() const; @@ -164,17 +152,15 @@ public: */ virtual int reportSolnBranchActual() const; - - //! Get the array of log concentration-like derivatives of the //! log activity coefficients /*! - * This function is a virtual method. For ideal mixtures - * (unity activity coefficients), this can return zero. - * Implementations should take the derivative of the - * logarithm of the activity coefficient with respect to the - * logarithm of the concentration-like variable (i.e. moles) - * that represents the standard state. + * For ideal mixtures (unity activity coefficients), this can return zero. + * Implementations should take the derivative of the logarithm of the + * activity coefficient with respect to the logarithm of the + * concentration-like variable (i.e. moles) that represents the standard + * state. + * * This quantity is to be used in conjunction with derivatives of * that concentration-like variable when the derivative of the chemical * potential is taken. @@ -188,12 +174,10 @@ public: err("getdlnActCoeffdlnN_diag"); } - //@} - /// @name Partial Molar Properties of the Solution (MixtureFugacityTP) + /// @name Partial Molar Properties of the Solution //@{ - //! Get the array of non-dimensional species chemical potentials //! These are partial molar Gibbs free energies. /*! @@ -210,10 +194,8 @@ public: void getChemPotentials_RT(doublereal* mu) const; //@} - /*! * @name Properties of the Standard State of the Species in the Solution - * (MixtureFugacityTP) * * Within MixtureFugacityTP, these properties are calculated via a common routine, * _updateStandardStateThermo(), @@ -249,7 +231,6 @@ public: */ virtual void getEnthalpy_RT(doublereal* hrt) const; - //! Get the array of nondimensional Enthalpy functions for the standard state species /*! * at the current T and P of the solution. @@ -297,7 +278,6 @@ public: */ virtual void getIntEnergy_RT(doublereal* urt) const; - //! Get the nondimensional Heat Capacities at constant //! pressure for the standard state of the species at the current T and P. /*! @@ -311,7 +291,6 @@ public: */ virtual void getCp_R(doublereal* cpr) const; - //! Get the molar volumes of each species in their standard //! states at the current T and P of the solution. /*! @@ -324,7 +303,7 @@ public: * units = m^3 / kmol */ virtual void getStandardVolumes(doublereal* vol) const; - + // @} //! Set the temperature of the phase /*! @@ -336,7 +315,6 @@ public: */ virtual void setTemperature(const doublereal temp); - //! Set the internally stored pressure (Pa) at constant //! temperature and composition /*! @@ -348,8 +326,6 @@ public: */ virtual void setPressure(doublereal p); - - protected: /** * Calculate the density of the mixture using the partial @@ -370,9 +346,6 @@ protected: * species molar volumes. We have additionally specified * in this class that the pure species molar volumes are * independent of temperature and pressure. - * - * NOTE: This is a non-virtual function, which is not a - * member of the ThermoPhase base class. */ virtual void calcDensity(); @@ -401,12 +374,10 @@ public: * * @param t Temperature (K) * @param p Pressure (Pa) - * @param x Vector of mole fractions. - * Length is equal to m_kk. + * @param x Vector of mole fractions. Length is equal to m_kk. */ virtual void setState_TPX(doublereal t, doublereal p, const doublereal* x); - //! Set the mass fractions to the specified values, and then //! normalize them so that they sum to 1.0. /*! @@ -415,7 +386,6 @@ public: */ virtual void setMassFractions(const doublereal* const y); - //!Set the mass fractions to the specified values without normalizing. /*! * This is useful when the normalization @@ -423,13 +393,10 @@ public: * by a constraint equation as part of a larger set of * equations. * - * @param y Input vector of mass fractions. - * Length is m_kk. + * @param y Input vector of mass fractions. Length is m_kk. */ virtual void setMassFractions_NoNorm(const doublereal* const y); - - //! Set the mole fractions to the specified values, and then //! normalize them so that they sum to 1.0. /*! @@ -438,19 +405,16 @@ public: */ virtual void setMoleFractions(const doublereal* const x); - //! Set the mole fractions to the specified values without normalizing. /*! * This is useful when the normalization * condition is being handled by some other means, for example * by a constraint equation as part of a larger set of equations. * - * @param x Input vector of mole fractions. - * Length is m_kk. + * @param x Input vector of mole fractions. Length is m_kk. */ virtual void setMoleFractions_NoNorm(const doublereal* const x); - //! Set the concentrations to the specified values within the phase. /*! * @param c The input vector to this routine is in dimensional @@ -465,7 +429,6 @@ public: protected: void setMoleFractions_NoState(const doublereal* const x); - public: //! Returns the current pressure of the phase /*! @@ -478,38 +441,24 @@ public: return m_Pcurrent; } - - - - protected: - //! Updates the reference state thermodynamic functions at the current T of the solution. /*! - * - * If m_useTmpStandardStateStorage is true, - * this function must be called for every call to functions in this - * class. It checks to see whether the temperature or pressure has changed and + * This function must be called for every call to functions in this + * class. It checks to see whether the temperature has changed and * thus the ss thermodynamics functions for all of the species * must be recalculated. * - * This function is responsible for updating the following internal members, - * when m_useTmpStandardStateStorage is true. - * - * - m_hss_RT; - * - m_cpss_R; - * - m_gss_RT; - * - m_sss_R; - * - m_Vss - * - * If m_useTmpStandardStateStorage is not true, this function may be - * required to be called by child classes to update internal member data. + * This function is responsible for updating the following internal members: * + * - m_h0_RT; + * - m_cp0_R; + * - m_g0_RT; + * - m_s0_R; */ virtual void _updateReferenceStateThermo() const; public: - //@} /// @name Thermodynamic Values for the Species Reference States (MixtureFugacityTP) /*! * There are also temporary @@ -520,7 +469,6 @@ public: */ //@{ - //! Returns the vector of nondimensional //! enthalpies of the reference state at the current temperature //! of the solution and the reference pressure for the species. @@ -549,7 +497,6 @@ public: //! Gibbs free energies of the reference state at the current temperature //! of the solution and the reference pressure for the species. /*! - * * @param grt Output vector contains the nondimensional Gibbs free energies * of the reference state of the species * length = m_kk, units = dimensionless. @@ -613,16 +560,8 @@ public: */ virtual void getStandardVolumes_ref(doublereal* vol) const; -protected: - - - //@} - - -public: - - //! @name Initialization Methods - For Internal use (VPStandardState) + //! @name Initialization Methods - For Internal use /*! * The following methods are used in the process of constructing * the phase and setting its parameters from a specification in an @@ -632,13 +571,11 @@ public: */ //@{ - //! Set the initial state of the phase to the conditions specified in the state XML element. /*! - * * This method sets the temperature, pressure, and mole fraction vector to a set default value. * - * @param state AN XML_Node object corresponding to + * @param state An XML_Node object corresponding to * the "state" entry for this phase in the input file. */ virtual void setStateFromXML(const XML_Node& state); @@ -690,17 +627,14 @@ public: */ virtual void initThermoXML(XML_Node& phaseNode, const std::string& id); - private: //! @internal Initialize the internal lengths in this object. - /*! - * Note this is not a virtual function. - */ void initLengths(); + //@} protected: - // Special Functions for fugacity classes - + //! @name Special Functions for fugacity classes + //! @{ //! Calculate the value of z /*! @@ -729,7 +663,6 @@ protected: */ virtual doublereal hresid() const; - //! Estimate for the saturation pressure /*! * Note: this is only used as a starting guess for later routines that actually calculate an @@ -740,6 +673,7 @@ protected: * @return returns the estimated saturation pressure at the given temperature */ virtual doublereal psatEst(doublereal TKelvin) const; + public: //! Estimate for the molar volume of the liquid /*! @@ -757,7 +691,6 @@ public: */ virtual doublereal liquidVolEst(doublereal TKelvin, doublereal& pres) const; -public: //! Calculates the density given the temperature and the pressure and a guess at the density. /*! * Note, below T_c, this is a multivalued function. We do not cross the vapor dome in this. @@ -777,26 +710,26 @@ public: * @param rhoguess Guessed density of the fluid. A value of -1.0 indicates that there * is no guessed density * - * * @return We return the density of the fluid at the requested phase. If we have not found any * acceptable density we return a -1. If we have found an acceptable density at a * different phase, we return a -2. */ virtual doublereal densityCalc(doublereal TKelvin, doublereal pressure, int phaseRequested, doublereal rhoguess); + protected: //! Utility routine in the calculation of the saturation pressure /*! - * Private routine - * * @param TKelvin temperature (kelvin) * @param pres pressure (Pascal) - * @param densLiq Output density of liquid - * @param densGas output density of gas - * @param gasGRT output delGRT + * @param[out] densLiq density of liquid + * @param[out] densGas density of gas + * @param[out] liqGRT deltaG/RT of liquid + * @param[out] gasGRT deltaG/RT of gas */ int corr0(doublereal TKelvin, doublereal pres, doublereal& densLiq, doublereal& densGas, doublereal& liqGRT, doublereal& gasGRT); + public: //! Returns the Phase State flag for the current state of the object /*! @@ -804,9 +737,9 @@ public: * in parameters space we are * * There are three values: - * WATER_GAS below the critical temperature but below the critical density - * WATER_LIQUID below the critical temperature but above the critical density - * WATER_SUPERCRIT above the critical temperature + * - WATER_GAS below the critical temperature but below the critical density + * - WATER_LIQUID below the critical temperature but above the critical density + * - WATER_SUPERCRIT above the critical temperature */ int phaseState(bool checkState = false) const ; @@ -817,7 +750,6 @@ public: */ virtual doublereal densSpinodalLiquid() const; - //! Return the value of the density at the gas spinodal point (on the gas side) //! for the current temperature. /*! @@ -825,9 +757,6 @@ public: */ virtual doublereal densSpinodalGas() const; - - - public: //! Calculate the saturation pressure at the current mixture content for the given temperature /*! @@ -839,6 +768,7 @@ public: */ doublereal calculatePsat(doublereal TKelvin, doublereal& molarVolGas, doublereal& molarVolLiquid); + protected: //! Calculate the pressure given the temperature and the molar volume /*! @@ -851,7 +781,6 @@ protected: */ virtual doublereal pressureCalc(doublereal TKelvin, doublereal molarVol) const; - //! Calculate the pressure and the pressure derivative given the temperature and the molar volume /*! * Temperature and mole number are held constant @@ -865,28 +794,19 @@ protected: */ virtual doublereal dpdVCalc(doublereal TKelvin, doublereal molarVol, doublereal& presCalc) const; - - virtual void updateMixingExpressions(); - //@} - class spinodalFunc : public Cantera::ResidEval { public: - spinodalFunc(MixtureFugacityTP* tp); - virtual int evalSS(const doublereal t, const doublereal* const y, doublereal* const r); - MixtureFugacityTP* m_tp; }; - protected: - //! Current value of the pressures /*! * Because the pressure is now a calculation, we store the result of the calculation whenever @@ -896,29 +816,21 @@ protected: */ doublereal m_Pcurrent; - //! Storage for the current values of the mole fractions of the species /*! * This vector is kept up-to-date when some the setState functions are called. - * - * The State object is allowed to com - * - * Therefore, it may be considered to be an independent variable. - * - */ std::vector moleFractions_; //! Current state of the fluid /*! - * There are three possible states of the fluid - * FLUID_GAS - * FLUID_LIQUID - * FLUID_SUPERCRIT + * There are three possible states of the fluid: + * - FLUID_GAS + * - FLUID_LIQUID + * - FLUID_SUPERCRIT */ int iState_; - //! Force the system to be on a particular side of the spinodal curve int forcedState_; @@ -941,14 +853,13 @@ protected: mutable vector_fp m_s0_R; spinodalFunc* fdpdv_; -private: +private: //! MixtureFugacityTP has its own err routine /*! * @param msg Error message string */ doublereal err(const std::string& msg) const; - }; } diff --git a/include/cantera/thermo/PureFluidPhase.h b/include/cantera/thermo/PureFluidPhase.h index 038470f21..7be1bf37b 100644 --- a/include/cantera/thermo/PureFluidPhase.h +++ b/include/cantera/thermo/PureFluidPhase.h @@ -19,7 +19,6 @@ namespace Cantera { - //! This phase object consists of a single component that can be a //! gas, a liquid, a mixed gas-liquid fluid, or a fluid beyond its //! critical point @@ -27,23 +26,10 @@ namespace Cantera * The object inherits from ThermoPhase. However, it's built on top * of the tpx package. * - * - *

Specification of Species Standard State Properties

- * - * - *

Application within %Kinetics Managers

- * - * - *

XML Example

- * - * - *

Instantiation of the Class

- * * @ingroup thermoprops */ class PureFluidPhase : public ThermoPhase { - public: //! Empty Base Constructor @@ -126,8 +112,6 @@ public: mu[0] = gibbs_mole(); } - - //! Get the species electrochemical potentials. /*! * These are partial molar quantities. This method adds a term \f$ F z_k @@ -159,7 +143,6 @@ public: */ virtual void getPartialMolarEnthalpies(doublereal* hbar) const; - //! Returns an array of partial molar entropies of the species in the //! solution. Units: J/kmol/K. /*! @@ -344,9 +327,6 @@ public: //! Returns the vector of nondimensional enthalpies of the reference state at the current temperature //! of the solution and the reference pressure for the species. /*! - * This base function will throw a Cantera exception unless - * it is overwritten in a derived class. - * * @param hrt Output vector containing the nondimensional reference state enthalpies * Length: m_kk. */ @@ -370,8 +350,6 @@ public: */ virtual void getGibbs_ref(doublereal* g) const; - - //! Returns the vector of nondimensional entropies of the reference state at the current temperature //! of the solution and the reference pressure for each species. /*! @@ -435,13 +413,9 @@ public: */ virtual void setState_SP(doublereal s, doublereal p, doublereal tol = 1.e-8); - //@} //! @name Critical State Properties - /*! - * Critical properties for the pure fluid - */ //@{ //! critical temperature @@ -456,14 +430,8 @@ public: //@} //! @name Saturation properties. - /*! - * These methods are only implemented by subclasses that - * implement full liquid-vapor equations of state. They may be - * moved out of ThermoPhase at a later date. - */ //@{ - //! saturation temperature /*! * @param p Pressure (Pa) @@ -515,7 +483,6 @@ public: //! 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 @@ -528,7 +495,6 @@ public: */ virtual void setParametersFromXML(const XML_Node& eosdata); - //! returns a summary of the state of the phase as a string /*! * @param show_thermo If true, extra information is printed out @@ -550,7 +516,6 @@ protected: void setTPXState() const; private: - //! Pointer to the underlying tpx object Substance that does the work mutable tpx::Substance* m_sub; diff --git a/include/cantera/thermo/RedlichKwongMFTP.h b/include/cantera/thermo/RedlichKwongMFTP.h index 09d243e87..14cf84f73 100644 --- a/include/cantera/thermo/RedlichKwongMFTP.h +++ b/include/cantera/thermo/RedlichKwongMFTP.h @@ -21,12 +21,6 @@ namespace Cantera { class XML_Node; -class PDSS; - -/*! - * @name CONSTANTS - Models for the Standard State of IdealSolnPhase's - */ -//@{ /** * @ingroup thermoprops @@ -34,39 +28,28 @@ class PDSS; * This class can handle either an ideal solution or an ideal gas approximation * of a phase. * - * * @nosubgrouping */ class RedlichKwongMFTP : public MixtureFugacityTP { - public: + //! @name Constructors and Duplicators + //! @{ - /*! - * - * @name Constructors and Duplicators for %RedlichKwongMFTP - * - */ //! Base constructor. RedlichKwongMFTP(); - //! Construct and initialize a RedlichKwongMFTP ThermoPhase object - //! directly from an ASCII input file + //! Construct and initialize a RedlichKwongMFTP object directly from an + //! ASCII input file /*! - * Working constructors - * - * The two constructors below are the normal way the phase initializes itself. They are shells that call - * the routine initThermo(), with a reference to the - * XML database to get the info for the phase. - * * @param infile Name of the input file containing the phase XML data * to set up the object * @param id ID of the phase in the input file. Defaults to the empty string. */ RedlichKwongMFTP(const std::string& infile, std::string id=""); - //! Construct and initialize a RedlichKwongMFTP ThermoPhase object - //! directly from an XML database + //! Construct and initialize a RedlichKwongMFTP object directly from an + //! XML database /*! * @param phaseRef XML phase node containing the description of the phase * @param id id attribute containing the name of the phase. (default is the empty string) @@ -76,7 +59,6 @@ public: //! This is a special constructor, used to replicate test problems //! during the initial verification of the object /*! - * * test problems: * 1: Pure CO2 problem * input file = CO2_RedlickKwongMFTP.xml @@ -107,7 +89,6 @@ public: //! Destructor. virtual ~RedlichKwongMFTP(); - //! Duplicator from the ThermoPhase parent class /*! * Given a pointer to a ThermoPhase object, this function will @@ -118,12 +99,6 @@ public: */ virtual ThermoPhase* duplMyselfAsThermoPhase() const; - //@} - - /** - * @name Utilities (RedlichKwongMFTP) - */ - //@{ /** * Equation of state type flag. The base class returns * zero. Subclasses should define this to return a unique @@ -132,7 +107,8 @@ public: */ virtual int eosType() const; - //@} + //! @name Molar Thermodynamic properties + //! @{ /// Molar enthalpy. Units: J/kmol. virtual doublereal enthalpy_mole() const; @@ -152,11 +128,9 @@ public: /// Molar heat capacity at constant volume. Units: J/kmol/K. virtual doublereal cv_mole() const; - /** - * @} - * @name Mechanical Properties - * @{ - */ + //! @} + //! @name Mechanical Properties + //! @{ //! Return the thermodynamic pressure (Pa). /*! @@ -167,7 +141,6 @@ public: * \f[ * P = \frac{RT}{v-b_{mix}} - \frac{a_{mix}}{T^{0.5} v \left( v + b_{mix} \right) } * \f] - * */ virtual doublereal pressure() const; @@ -179,6 +152,7 @@ public: * \f] */ virtual doublereal isothermalCompressibility() const; + // @} protected: /** @@ -200,9 +174,6 @@ protected: * species standard state molar volumes. * The species molar volumes may be functions * of temperature and pressure. - * - * NOTE: This is a non-virtual function, which is not a - * member of the ThermoPhase base class. */ virtual void calcDensity(); @@ -213,8 +184,6 @@ protected: * function sets the temperature, and makes sure that * the value propagates to underlying objects * - * @todo Make Phase::setTemperature a virtual function - * * @param temp Temperature in kelvin */ virtual void setTemperature(const doublereal temp); @@ -234,8 +203,7 @@ protected: * by a constraint equation as part of a larger set of * equations. * - * @param y Input vector of mass fractions. - * Length is m_kk. + * @param y Input vector of mass fractions. Length is m_kk. */ virtual void setMassFractions_NoNorm(const doublereal* const y); @@ -253,12 +221,10 @@ protected: * condition is being handled by some other means, for example * by a constraint equation as part of a larger set ofequations. * - * @param x Input vector of mole fractions. - * Length is m_kk. + * @param x Input vector of mole fractions. Length is m_kk. */ virtual void setMoleFractions_NoNorm(const doublereal* const x); - //! Set the concentrations to the specified values within the phase. /*! * @param c The input vector to this routine is in dimensional @@ -270,9 +236,7 @@ protected: */ virtual void setConcentrations(const doublereal* const c); - public: - //! This method returns an array of generalized concentrations /*! * \f$ C^a_k\f$ are defined such that \f$ a_k = C^a_k / @@ -330,14 +294,16 @@ public: * 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 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. @@ -355,8 +321,7 @@ public: */ virtual void getActivityCoefficients(doublereal* ac) const; - - /// @name Partial Molar Properties of the Solution (RedlichKwongMFTP) + /// @name Partial Molar Properties of the Solution //@{ //! Get the array of non-dimensional species chemical potentials. @@ -422,35 +387,9 @@ public: virtual void getPartialMolarVolumes(doublereal* vbar) const; //@} - - /*! - * @name Properties of the Standard State of the Species in the Solution - * - * Properties of the standard states are delegated to the VPSSMgr object. - * The values are cached within this object, and are not recalculated unless - * the temperature or pressure changes. - */ - //@{ - - //@} - - /// @name Thermodynamic Values for the Species Reference States (RedlichKwongMFTP) - /*! - * Properties of the reference states are delegated to the VPSSMgr object. - * The values are cached within this object, and are not recalculated unless - * the temperature or pressure changes. - */ - //@{ - - //@} - - - - //--------------------------------------------------------- /// @name Critical State Properties. /// These methods are only implemented by some subclasses, and may /// be moved out of ThermoPhase at a later date. - //@{ /// Critical temperature (K). @@ -461,14 +400,11 @@ public: /// Critical density (kg/m3). virtual doublereal critDensity() const; - //@} - - - public: - //! @name Initialization Methods - For Internal use (VPStandardState) + //@} + //! @name Initialization Methods - For Internal use /*! * The following methods are used in the process of constructing * the phase and setting its parameters from a specification in an @@ -478,9 +414,7 @@ public: */ //@{ - - //! Set equation of state parameter values from XML - //! entries. + //! 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 @@ -509,7 +443,6 @@ public: */ virtual void initThermo(); - //!This method is used by the ChemEquil equilibrium solver. /*! * It sets the state such that the chemical potentials satisfy @@ -560,32 +493,25 @@ private: */ void readXMLPureFluid(XML_Node& pureFluidParam); - //! Apply mixing rules for a coefficients void applyStandardMixingRules(); - //! Read the cross species RedlichKwong input parameters /*! * @param pureFluidParam XML_Node for the cross fluid parameters */ void readXMLCrossFluid(XML_Node& pureFluidParam); - - - //============================================================================== private: //! @internal Initialize the internal lengths in this object. /*! - * Note this is not a virtual function and only handles - * this object + * Note this is not a virtual function and only handles this object */ void initLengths(); - - //============================================================================== - // Special functions inherited from MixtureFugacityTP + // @} protected: + // Special functions inherited from MixtureFugacityTP //! Calculate the deviation terms for the total entropy of the mixture from the //! ideal gas mixture @@ -647,7 +573,6 @@ public: */ virtual doublereal densityCalc(doublereal TKelvin, doublereal pressure, int phase, doublereal rhoguess); -public: //! Return the value of the density at the liquid spinodal point (on the liquid side) //! for the current temperature. /*! @@ -655,7 +580,6 @@ public: */ virtual doublereal densSpinodalLiquid() const; - //! Return the value of the density at the gas spinodal point (on the gas side) //! for the current temperature. /*! @@ -663,8 +587,6 @@ public: */ virtual doublereal densSpinodalGas() const; - - //! Calculate the pressure given the temperature and the molar volume /*! * Calculate the pressure given the temperature and the molar volume @@ -676,7 +598,6 @@ public: */ virtual doublereal pressureCalc(doublereal TKelvin, doublereal molarVol) const; - //! Calculate the pressure and the pressure derivative given the temperature and the molar volume /*! * Temperature and mole number are held constant @@ -690,17 +611,14 @@ public: */ virtual doublereal dpdVCalc(doublereal TKelvin, doublereal molarVol, doublereal& presCalc) const; - //! Calculate dpdV and dpdT at the current conditions /*! * These are stored internally. */ void pressureDerivatives() const; - virtual void updateMixingExpressions(); - //! Update the a and b parameters /*! * The a and the b parameters depend on the mole fraction and the temperature. @@ -708,10 +626,8 @@ public: */ void updateAB(); - //! Calculate the a and the b parameters given the temperature /*! - * * This function doesn't change the internal state of the object, so it is a const * function. It does use the stored mole fractions in the object. * @@ -722,24 +638,27 @@ public: */ void calculateAB(doublereal temp, doublereal& aCalc, doublereal& bCalc) const; - - //========================================================================================= - // Special functions not inherited from MixtureFugacityTP + // Special functions not inherited from MixtureFugacityTP doublereal da_dt() const; void calcCriticalConditions(doublereal a, doublereal b, doublereal a0_coeff, doublereal aT_coeff, doublereal& pc, doublereal& tc, doublereal& vc) const; - - + //! Solve the cubic equation of state + /*! + * The R-K equation of state may be solved via the following formula: + * + * V**3 - V**2(RT/P) - V(RTb/P - a/(P T**.5) + b*b) - (a b / (P T**.5)) = 0 + * + * Returns the number of solutions found. If it only finds the liquid + * branch solution, it will return a -1 or a -2 instead of 1 or 2. If it + * returns 0, then there is an error. + */ int NicholsSolve(double TKelvin, double pres, doublereal a, doublereal b, doublereal Vroot[3]) const; - //@} - //============================================================================== protected: - //! boolean indicating whether standard mixing rules are applied /*! * - 1 = Yes, there are standard cross terms in the a coefficient matrices. @@ -754,7 +673,6 @@ protected: */ int m_formTempParam; - //! Value of b in the equation of state /*! * m_b is a function of the temperature and the mole fraction. @@ -767,13 +685,11 @@ protected: */ doublereal m_a_current; - vector_fp a_vec_Curr_; vector_fp b_vec_Curr_; Array2D a_coeff_vec; - vector_fp m_pc_Species; vector_fp m_tc_Species; vector_fp m_vc_Species; @@ -782,8 +698,6 @@ protected: doublereal Vroot_[3]; - - //! Temporary storage - length = m_kk. mutable vector_fp m_pp; @@ -795,8 +709,6 @@ protected: // Partial molar volumes of the species mutable vector_fp m_partialMolarVolumes; - - //! The derivative of the pressure wrt the volume /*! * Calculated at the current conditions @@ -830,8 +742,6 @@ public: //! Omega constant for the critical molar volume static const doublereal omega_vc; - - }; } diff --git a/src/thermo/IdealGasPhase.cpp b/src/thermo/IdealGasPhase.cpp index 6ac1b45b7..626534d58 100644 --- a/src/thermo/IdealGasPhase.cpp +++ b/src/thermo/IdealGasPhase.cpp @@ -14,7 +14,7 @@ using namespace std; namespace Cantera { -// Default empty Constructor + IdealGasPhase::IdealGasPhase() : m_p0(-1.0), m_tlast(0.0), @@ -42,7 +42,6 @@ IdealGasPhase::~IdealGasPhase() { } -// Copy Constructor IdealGasPhase::IdealGasPhase(const IdealGasPhase& right) : m_p0(right.m_p0), m_tlast(right.m_tlast), @@ -55,14 +54,6 @@ IdealGasPhase::IdealGasPhase(const IdealGasPhase& right) : *this = right; } -// Assignment operator -/* - * Assignment operator for the object. Constructed - * object will be a clone of this object, but will - * also own all of its data. - * - * @param right Object to be copied. - */ IdealGasPhase& IdealGasPhase::operator=(const IdealGasPhase& right) { if (&right != this) { @@ -80,14 +71,6 @@ IdealGasPhase& IdealGasPhase::operator=(const IdealGasPhase& right) return *this; } -// Duplicator from the %ThermoPhase parent class -/* - * Given a pointer to a %ThermoPhase object, this function will - * duplicate the %ThermoPhase object and all underlying structures. - * This is basically a wrapper around the copy constructor. - * - * @return returns a pointer to a %ThermoPhase - */ ThermoPhase* IdealGasPhase::duplMyselfAsThermoPhase() const { return new IdealGasPhase(*this); @@ -95,80 +78,31 @@ ThermoPhase* IdealGasPhase::duplMyselfAsThermoPhase() const // Molar Thermodynamic Properties of the Solution ------------------ -/* - * Molar internal energy. J/kmol. For an ideal gas mixture, - * \f[ - * \hat u(T) = \sum_k X_k \hat h^0_k(T) - \hat R T, - * \f] - * and is a function only of temperature. - * The reference-state pure-species enthalpies - * \f$ \hat h^0_k(T) \f$ are computed by the species thermodynamic - * property manager. - * @see SpeciesThermo - */ doublereal IdealGasPhase::intEnergy_mole() const { return GasConstant * temperature() * (mean_X(&enthalpy_RT_ref()[0]) - 1.0); } -/* - * Molar entropy. Units: J/kmol/K. - * For an ideal gas mixture, - * \f[ - * \hat s(T, P) = \sum_k X_k \hat s^0_k(T) - \hat R \log (P/P^0). - * \f] - * The reference-state pure-species entropies - * \f$ \hat s^0_k(T) \f$ are computed by the species thermodynamic - * property manager. - * @see SpeciesThermo - */ doublereal IdealGasPhase::entropy_mole() const { return GasConstant * (mean_X(&entropy_R_ref()[0]) - sum_xlogx() - std::log(pressure() / m_spthermo->refPressure())); } -/* - * Molar Gibbs free Energy for an ideal gas. - * Units = J/kmol. - */ doublereal IdealGasPhase::gibbs_mole() const { return enthalpy_mole() - temperature() * entropy_mole(); } -/* - * Molar heat capacity at constant pressure. Units: J/kmol/K. - * For an ideal gas mixture, - * \f[ - * \hat c_p(t) = \sum_k \hat c^0_{p,k}(T). - * \f] - * The reference-state pure-species heat capacities - * \f$ \hat c^0_{p,k}(T) \f$ are computed by the species thermodynamic - * property manager. - * @see SpeciesThermo - */ doublereal IdealGasPhase::cp_mole() const { return GasConstant * mean_X(&cp_R_ref()[0]); } -/* - * Molar heat capacity at constant volume. Units: J/kmol/K. - * For an ideal gas mixture, - * \f[ \hat c_v = \hat c_p - \hat R. \f] - */ doublereal IdealGasPhase::cv_mole() const { return cp_mole() - GasConstant; } -/** - * @returns species translational/rotational specific heat at - * constant volume. - * - * Either: $5/2 R_s$ or $3/2 R_s$ for molecules/atoms. - * - */ doublereal IdealGasPhase::cv_tr(doublereal atomicity) const { // k is the species number @@ -189,44 +123,16 @@ doublereal IdealGasPhase::cv_tr(doublereal atomicity) const return c[3]; } -/** - * @returns species translational specific heat at constant volume. - */ doublereal IdealGasPhase::cv_trans() const { return 1.5 * GasConstant; } -/** - * @returns species rotational specific heat at constant volume. - * - */ doublereal IdealGasPhase::cv_rot(double atom) const { return std::max(cv_tr(atom) - cv_trans(), 0.); } -/** - * @returns species vibrational specific heat at - * constant volume. - * - * C^{vib}_{v,s} = \frac{\partial e^{vib}_{v,s} }{\partial T} - * - * The species vibration energy ($e^{vib}_{v,s}$) is: - * - * 0: atom - * - * Diatomic: - * \f[ - * \frac{R_s \theta_{v,s}}{e^{\theta_{v,s}/T}-1} - * \f] - * - * General Molecules: - * \f[ - * \sum_i \frac{R_s \theta_{v,s,i}}{e^{\theta_{v,s,i}/T}-1} - * \f] - * - */ doublereal IdealGasPhase::cv_vib(const int k, const doublereal T) const { @@ -253,23 +159,12 @@ doublereal IdealGasPhase::cv_vib(const int k, const doublereal T) const } -// Mechanical Equation of State ---------------------------- -// Chemical Potentials and Activities ---------------------- - -/* - * Returns the standard concentration \f$ C^0_k \f$, which is used to normalize - * the generalized concentration. - */ doublereal IdealGasPhase::standardConcentration(size_t k) const { double p = pressure(); return p / (GasConstant * temperature()); } -/* - * Returns the natural logarithm of the standard - * concentration of the kth species - */ doublereal IdealGasPhase::logStandardConc(size_t k) const { _updateThermo(); @@ -278,9 +173,6 @@ doublereal IdealGasPhase::logStandardConc(size_t k) const return lc; } -/* - * Get the array of non-dimensional activity coefficients - */ void IdealGasPhase::getActivityCoefficients(doublereal* ac) const { for (size_t k = 0; k < m_kk; k++) { @@ -288,10 +180,6 @@ void IdealGasPhase::getActivityCoefficients(doublereal* ac) const } } -/* - * Get the array of chemical potentials at unit activity \f$ - * \mu^0_k(T,P) \f$. - */ void IdealGasPhase::getStandardChemPotentials(doublereal* muStar) const { const vector_fp& gibbsrt = gibbs_RT_ref(); @@ -318,10 +206,6 @@ void IdealGasPhase::getChemPotentials(doublereal* mu) const } } -/* - * Get the array of partial molar enthalpies of the species - * units = J / kmol - */ void IdealGasPhase::getPartialMolarEnthalpies(doublereal* hbar) const { const vector_fp& _h = enthalpy_RT_ref(); @@ -329,10 +213,6 @@ void IdealGasPhase::getPartialMolarEnthalpies(doublereal* hbar) const scale(_h.begin(), _h.end(), hbar, rt); } -/* - * Get the array of partial molar entropies of the species - * units = J / kmol / K - */ void IdealGasPhase::getPartialMolarEntropies(doublereal* sbar) const { const vector_fp& _s = entropy_R_ref(); @@ -345,10 +225,6 @@ void IdealGasPhase::getPartialMolarEntropies(doublereal* sbar) const } } -/* - * Get the array of partial molar internal energies of the species - * units = J / kmol - */ void IdealGasPhase::getPartialMolarIntEnergies(doublereal* ubar) const { const vector_fp& _h = enthalpy_RT_ref(); @@ -358,19 +234,12 @@ void IdealGasPhase::getPartialMolarIntEnergies(doublereal* ubar) const } } -/* - * Get the array of partial molar heat capacities - */ void IdealGasPhase::getPartialMolarCp(doublereal* cpbar) const { const vector_fp& _cp = cp_R_ref(); scale(_cp.begin(), _cp.end(), cpbar, GasConstant); } -/* - * Get the array of partial molar volumes - * units = m^3 / kmol - */ void IdealGasPhase::getPartialMolarVolumes(doublereal* vbar) const { double vol = 1.0 / molarDensity(); @@ -381,22 +250,12 @@ void IdealGasPhase::getPartialMolarVolumes(doublereal* vbar) const // 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 - */ void IdealGasPhase::getEnthalpy_RT(doublereal* hrt) const { const vector_fp& _h = enthalpy_RT_ref(); copy(_h.begin(), _h.end(), hrt); } -/* - * Get the array of nondimensional entropy functions for the - * standard state species - * at the current T and P of the solution. - */ void IdealGasPhase::getEntropy_R(doublereal* sr) const { const vector_fp& _s = entropy_R_ref(); @@ -407,10 +266,6 @@ void IdealGasPhase::getEntropy_R(doublereal* sr) const } } -/* - * Get the nondimensional gibbs function for the species - * standard states at the current T and P of the solution. - */ void IdealGasPhase::getGibbs_RT(doublereal* grt) const { const vector_fp& gibbsrt = gibbs_RT_ref(); @@ -421,11 +276,6 @@ void IdealGasPhase::getGibbs_RT(doublereal* grt) const } } -/* - * get the pure Gibbs free energies of each species assuming - * it is in its standard state. This is the same as - * getStandardChemPotentials(). - */ void IdealGasPhase::getPureGibbs(doublereal* gpure) const { const vector_fp& gibbsrt = gibbs_RT_ref(); @@ -437,11 +287,6 @@ void IdealGasPhase::getPureGibbs(doublereal* gpure) const } } -/* - * Returns the vector of nondimensional - * internal Energies of the standard state at the current temperature - * and pressure of the solution for each species. - */ void IdealGasPhase::getIntEnergy_RT(doublereal* urt) const { const vector_fp& _h = enthalpy_RT_ref(); @@ -450,25 +295,12 @@ void IdealGasPhase::getIntEnergy_RT(doublereal* urt) const } } -/* - * Get the nondimensional heat capacity at constant pressure - * function for the species - * standard states at the current T and P of the solution. - */ void IdealGasPhase::getCp_R(doublereal* cpr) const { const vector_fp& _cpr = cp_R_ref(); copy(_cpr.begin(), _cpr.end(), cpr); } -/* - * Get the molar volumes of the species standard states at the current - * T and P of the solution. - * units = m^3 / kmol - * - * @param vol Output vector containing the standard state volumes. - * Length: m_kk. - */ void IdealGasPhase::getStandardVolumes(doublereal* vol) const { double tmp = 1.0 / molarDensity(); @@ -479,56 +311,30 @@ void IdealGasPhase::getStandardVolumes(doublereal* vol) const // Thermodynamic Values for the Species Reference States --------- -/* - * Returns the vector of nondimensional - * enthalpies of the reference state at the current temperature - * and reference pressure. - */ void IdealGasPhase::getEnthalpy_RT_ref(doublereal* hrt) const { const vector_fp& _h = enthalpy_RT_ref(); copy(_h.begin(), _h.end(), hrt); } -/* - * Returns the vector of nondimensional - * enthalpies of the reference state at the current temperature - * and reference pressure. - */ void IdealGasPhase::getGibbs_RT_ref(doublereal* grt) const { const vector_fp& gibbsrt = gibbs_RT_ref(); copy(gibbsrt.begin(), gibbsrt.end(), grt); } -/* - * Returns the vector of the - * gibbs function of the reference state at the current temperature - * and reference pressure. - * units = J/kmol - */ void IdealGasPhase::getGibbs_ref(doublereal* g) const { const vector_fp& gibbsrt = gibbs_RT_ref(); scale(gibbsrt.begin(), gibbsrt.end(), g, _RT()); } -/* - * Returns the vector of nondimensional - * entropies of the reference state at the current temperature - * and reference pressure. - */ void IdealGasPhase::getEntropy_R_ref(doublereal* er) const { const vector_fp& _s = entropy_R_ref(); copy(_s.begin(), _s.end(), er); } -/* - * 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. - */ void IdealGasPhase::getIntEnergy_RT_ref(doublereal* urt) const { const vector_fp& _h = enthalpy_RT_ref(); @@ -537,11 +343,6 @@ void IdealGasPhase::getIntEnergy_RT_ref(doublereal* urt) const } } -/* - * Returns the vector of nondimensional - * constant pressure heat capacities of the reference state - * at the current temperature and reference pressure. - */ void IdealGasPhase::getCp_R_ref(doublereal* cprt) const { const vector_fp& _cpr = cp_R_ref(); @@ -556,8 +357,6 @@ void IdealGasPhase::getStandardVolumes_ref(doublereal* vol) const } } -// new methods defined here ------------------------------- - void IdealGasPhase::initThermo() { m_p0 = refPressure(); @@ -569,11 +368,6 @@ void IdealGasPhase::initThermo() m_pp.resize(m_kk); } -/* - * Set mixture to an equilibrium state consistent with specified - * chemical potentials and temperature. This method is needed by - * the ChemEquil equilibrium solver. - */ void IdealGasPhase::setToEquilState(const doublereal* mu_RT) { double tmp, tmp2; @@ -605,17 +399,6 @@ void IdealGasPhase::setToEquilState(const doublereal* mu_RT) setState_PX(pres, &m_pp[0]); } -/// This method is called each time a thermodynamic property is -/// requested, to check whether the internal species properties -/// within the object need to be updated. -/// Currently, this updates the species thermo polynomial values -/// for the current value of the temperature. A check is made -/// to see if the temperature has changed since the last -/// evaluation. This object does not contain any persistent -/// data that depends on the concentration, that needs to be -/// updated. The state object modifies its concentration -/// dependent information at the time the setMoleFractions() -/// (or equivalent) call is made. void IdealGasPhase::_updateThermo() const { doublereal tnow = temperature(); @@ -635,4 +418,3 @@ void IdealGasPhase::_updateThermo() const } } } - diff --git a/src/thermo/IdealSolidSolnPhase.cpp b/src/thermo/IdealSolidSolnPhase.cpp index 622365b46..2b99b5047 100644 --- a/src/thermo/IdealSolidSolnPhase.cpp +++ b/src/thermo/IdealSolidSolnPhase.cpp @@ -21,11 +21,6 @@ using namespace std; namespace Cantera { -/* - * Constructor for IdealSolidSolnPhase class: - * The default form for the generalized concentrations is 0 - * i.e., unity. - */ IdealSolidSolnPhase::IdealSolidSolnPhase(int formGC) : ThermoPhase(), m_formGC(formGC), @@ -53,7 +48,7 @@ IdealSolidSolnPhase::IdealSolidSolnPhase(const std::string& inputFile, } initThermoFile(inputFile, id); } -//==================================================================================================================== + IdealSolidSolnPhase::IdealSolidSolnPhase(XML_Node& root, const std::string& id, int formGC) : ThermoPhase(), @@ -68,12 +63,11 @@ IdealSolidSolnPhase::IdealSolidSolnPhase(XML_Node& root, const std::string& id, } importPhase(*findXMLPhase(&root, id), this); } -//==================================================================================================================== + IdealSolidSolnPhase::IdealSolidSolnPhase(const IdealSolidSolnPhase& b) { *this = b; } -//==================================================================================================================== IdealSolidSolnPhase& IdealSolidSolnPhase:: operator=(const IdealSolidSolnPhase& b) @@ -97,22 +91,11 @@ operator=(const IdealSolidSolnPhase& b) return *this; } -/* - * Base Class Duplication Function - * -> given a pointer to ThermoPhase, this function can - * duplicate the object. (note has to be a separate function - * not the copy constructor, because it has to be - * a virtual function) - */ ThermoPhase* IdealSolidSolnPhase::duplMyselfAsThermoPhase() const { return new IdealSolidSolnPhase(*this); } -//==================================================================================================================== -/** - * Equation of state flag. Returns the value cIdealGas, defined - * in mix_defs.h. - */ + int IdealSolidSolnPhase::eosType() const { integer res; @@ -136,19 +119,7 @@ int IdealSolidSolnPhase::eosType() const /******************************************************************** * Molar Thermodynamic Properties of the Solution ********************************************************************/ -/** - * Molar enthalpy of the solution. Units: J/kmol. - * For an ideal, constant partial molar volume solution mixture with - * pure species phases which exhibit zero volume expansivity and - * zero isothermal compressibility: - * \f[ - * \hat h(T,P) = \sum_k X_k \hat h^0_k(T) + (P - P_{ref}) (\sum_k X_k \hat V^0_k) - * \f] - * The reference-state pure-species enthalpies at the reference pressure Pref - * \f$ \hat h^0_k(T) \f$, are computed by the species thermodynamic - * property manager. They are polynomial functions of temperature. - * @see SpeciesThermo - */ + doublereal IdealSolidSolnPhase:: enthalpy_mole() const { @@ -157,21 +128,6 @@ enthalpy_mole() const return (htp + (pressure() - m_Pref)/molarDensity()); } -/** - * Molar internal energy of the solution. J/kmol. - * For an ideal, constant partial molar volume solution mixture with - * pure species phases which exhibit zero volume expansivity and - * zero isothermal compressibility: - * \f[ - * \hat u(T) = \hat h(T,P) - p \hat V = \sum_k X_k \hat h^0_k(T) - * - P_{ref} (\sum_k X_k \hat V^0_k) - * \f] - * and is a function only of temperature. - * The reference-state pure-species enthalpies - * \f$ \hat h^0_k(T) \f$ are computed by the species thermodynamic - * property manager. - * @see SpeciesThermo - */ doublereal IdealSolidSolnPhase::intEnergy_mole() const { const double* eptr = DATA_PTR(enthalpy_RT_ref().begin()); @@ -180,39 +136,12 @@ doublereal IdealSolidSolnPhase::intEnergy_mole() const return (htp - m_Pref / molarDensity()); } -/** - * Molar entropy of the solution. Units: J/kmol/K. - * For an ideal, constant partial molar volume solution mixture with - * pure species phases which exhibit zero volume expansivity: - * \f[ - * \hat s(T, P, X_k) = \sum_k X_k \hat s^0_k(T) - * - \hat R \sum_k X_k log(X_k) - * \f] - * The reference-state pure-species entropies - * \f$ \hat s^0_k(T,p_{ref}) \f$ are computed by the species thermodynamic - * property manager. The pure species entropies are independent of - * temperature since the volume expansivities are equal to zero. - * @see SpeciesThermo - */ doublereal IdealSolidSolnPhase::entropy_mole() const { const double* dptr = DATA_PTR(entropy_R_ref()); return GasConstant * (mean_X(dptr) - sum_xlogx()); } -/** - * Molar gibbs free energy of the solution. Units: J/kmol. - * For an ideal, constant partial molar volume solution mixture with - * pure species phases which exhibit zero volume expansivity: - * \f[ - * \hat g(T, P) = \sum_k X_k \hat g^0_k(T,P) + \hat R T \sum_k X_k log(X_k) - * \f] - * The reference-state pure-species gibbs free energies - * \f$ \hat g^0_k(T) \f$ are computed by the species thermodynamic - * property manager, while the standard state gibbs free energies - * \f$ \hat g^0_k(T,P) \f$ are computed by the member function, gibbs_RT(). - * @see SpeciesThermo - */ doublereal IdealSolidSolnPhase::gibbs_mole() const { const double* dptr = DATA_PTR(gibbs_RT_ref()); @@ -220,20 +149,6 @@ doublereal IdealSolidSolnPhase::gibbs_mole() const return (GasConstant * temperature() * (g + sum_xlogx())); } -/** - * Molar heat capacity at constant pressure of the solution. - * Units: J/kmol/K. - * For an ideal, constant partial molar volume solution mixture with - * pure species phases which exhibit zero volume expansivity: - * \f[ - * \hat c_p(T,P) = \sum_k X_k \hat c^0_{p,k}(T) . - * \f] - * The heat capacity is independent of pressure. - * The reference-state pure-species heat capacities - * \f$ \hat c^0_{p,k}(T) \f$ are computed by the species thermodynamic - * property manager. - * @see SpeciesThermo - */ doublereal IdealSolidSolnPhase::cp_mole() const { const double* dptr = DATA_PTR(cp_R_ref()); @@ -243,27 +158,7 @@ doublereal IdealSolidSolnPhase::cp_mole() const /******************************************************************** * Mechanical Equation of State ********************************************************************/ -/** - * - * Calculate the density of the mixture using the partial - * molar volumes and mole fractions as input - * - * The formula for this is - * - * \f[ - * \rho = \frac{\sum_k{X_k W_k}}{\sum_k{X_k V_k}} - * \f] - * - * where \f$ X_k \f$ are the mole fractions, \f$W_k\f$ are - * the molecular weights, and \f$V_k\f$ are the pure species - * molar volumes. - * - * Note, the basis behind this formula is that in an ideal - * solution the partial molar volumes are equal to the pure - * species molar volumes. We have additionally specified that - * in this class that the pure species molar volumes are - * independent of temperature and pressure. - */ + void IdealSolidSolnPhase::calcDensity() { /* @@ -280,22 +175,6 @@ void IdealSolidSolnPhase::calcDensity() Phase::setDensity(dens); } -/** - * Overwritten setDensity() function is necessary because the - * density is not an independent variable. - * - * This function will now throw an error condition - * - * @internal May have to adjust the strategy here to make - * the eos for these materials slightly compressible, in order - * to create a condition where the density is a function of - * the pressure. - * - * This function will now throw an error condition. - * - * NOTE: This is a virtual function that overwrites the State.h - * class - */ void IdealSolidSolnPhase:: setDensity(const doublereal rho) { @@ -312,90 +191,42 @@ setDensity(const doublereal rho) } } -/* - * setPressure(double) (virtual from ThermoPhase) - * - * Set the pressure at constant temperature. Units: Pa. - * This method sets a constant within the object. - * The mass density is not a function of pressure. - * Note: This function overrides the setPressure() function - * in the ThermoPhase object. - * We calculate the density and store it in the - * State object, because this density is supposed to - * be current after setting the pressure, and is now - * a dependent variable. - */ void IdealSolidSolnPhase::setPressure(doublereal p) { m_Pcurrent = p; calcDensity(); } -/* - * setMolarDensity() (virtual from State) - * Overwritten setMolarDensity() function is necessary because the - * density is not an independent variable. - * - * This function will now throw an error condition. - * - * NOTE: This is a virtual function that overrides the State.h - * class - */ void IdealSolidSolnPhase::setMolarDensity(const doublereal n) { throw CanteraError("IdealSolidSolnPhase::setMolarDensity", "Density is not an independent variable"); } -/* - * setMoleFractions() (virtual from State) - * - * Sets the mole fractions and adjusts the internal density. - */ void IdealSolidSolnPhase::setMoleFractions(const doublereal* const x) { Phase::setMoleFractions(x); calcDensity(); } -/** - * setMoleFractions_NoNorm() (virtual from State) - * - * Sets the mole fractions and adjusts the internal density. - */ void IdealSolidSolnPhase::setMoleFractions_NoNorm(const doublereal* const x) { Phase::setMoleFractions(x); calcDensity(); } -/* - * setMassFractions() (virtual from State) - * - * Sets the mass fractions and adjusts the internal density. - */ void IdealSolidSolnPhase::setMassFractions(const doublereal* const y) { Phase::setMassFractions(y); calcDensity(); } -/* - * setMassFractions_NoNorm() (virtual from State) - * - * Sets the mass fractions and adjusts the internal density. - */ void IdealSolidSolnPhase::setMassFractions_NoNorm(const doublereal* const y) { Phase::setMassFractions_NoNorm(y); calcDensity(); } -/* - * setConcentrations (virtual from State) - * - * Sets the concentrations and adjusts the internal density - */ void IdealSolidSolnPhase::setConcentrations(const doublereal* const c) { Phase::setConcentrations(c); @@ -406,54 +237,6 @@ void IdealSolidSolnPhase::setConcentrations(const doublereal* const c) * Chemical Potentials and Activities ********************************************************************/ -/******************************************************************** - * - * getActivitConcentrations(): - * - * This method returns the array of generalized - * concentrations. The generalized concentrations are used - * in the evaluation of the rates of progress for reactions - * involving species in this phase. The generalized - * concentration divided by the standard concentration is also - * equal to the activity of species. - * - * For this implementation the activity is defined to be the - * mole fraction of the species. The generalized concentration - * is defined to be equal to the mole fraction divided by - * the partial molar volume. The generalized concentrations - * for species in this phase therefore have units of - * kmol m-3. Rate constants must reflect this fact. - * - * On a general note, the following must be true. - * For an ideal solution, the generalized concentration must consist - * of the mole fraction multiplied by a constant. The constant may be - * fairly arbitrarily chosen, with differences adsorbed into the - * reaction rate expression. 1/V_N, 1/V_k, or 1 are equally good, - * as long as the standard concentration is adjusted accordingly. - * However, it must be a constant (and not the concentration, btw, - * which is a function of the mole fractions) in order for the - * ideal solution properties to hold at the same time having the - * standard concentration to be independent of the mole fractions. - * - * In this implementation the form of the generalized concentrations - * depend upon the member attribute, m_formGC: - * - * - * - * - * - * - *
m_formGC GeneralizedConc StandardConc
0 X_k 1.0
1 X_k / V_k 1.0 / V_k
2 X_k / V_N 1.0 / V_N
- * - * HKM Note: We have absorbed the pressure dependence of the pure species - * state into the thermodynamics functions. Therefore the - * standard state on which the activities are based depend - * on both temperature and pressure. If we hadn't, it would have - * appeared in this function in a very awkward exp[] format. - * - * @param c[] Pointer to array of doubles of length m_kk, which on exit - * will contain the generalized concentrations. - */ void IdealSolidSolnPhase:: getActivityConcentrations(doublereal* c) const { @@ -479,23 +262,6 @@ getActivityConcentrations(doublereal* c) const } } -/********************************************************************* - * - * standardConcentration() - * - * The standard concentration \f$ C^0_k \f$ used to normalize - * the generalized concentration. - * In many cases, this quantity - * will be the same for all species in a phase. - * However, for this case, we will return a distinct concentration - * for each species. This is the inverse of the species molar - * volume. Units are m3 kmol-1. - * - * - * @param k Species number: this is a require parameter, - * a change from the ThermoPhase base class, where it was - * an optional parameter. - */ doublereal IdealSolidSolnPhase:: standardConcentration(size_t k) const { @@ -523,16 +289,6 @@ referenceConcentration(int k) const return 0.0; } -/********************************************************************* - * - * logStandardConc() - * - * Returns the log of the standard concentration - * - * @param k Species number: this is a require parameter, - * a change from the ThermoPhase base class, where it was - * an optional parameter. - */ doublereal IdealSolidSolnPhase:: logStandardConc(size_t k) const { @@ -555,32 +311,6 @@ logStandardConc(size_t k) const return res; } -/*********************************************************************** - * - * getUnitsStandardConcentration() - * - * Returns the units of the standard and general concentrations - * Note they have the same units, as their divisor 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 - * - * For EOS types other than cIdealSolidSolnPhase1, the default - * kmol/m3 holds for standard concentration units. For - * cIdealSolidSolnPhase0 type, the standard concentration is - * unitless. - */ void IdealSolidSolnPhase:: getUnitsStandardConc(double* uA, int, int sizeUA) const { @@ -613,10 +343,6 @@ getUnitsStandardConc(double* uA, int, int sizeUA) const } } -/* - * getActivityCoefficients(): - * - */ void IdealSolidSolnPhase:: getActivityCoefficients(doublereal* ac) const { @@ -624,23 +350,7 @@ getActivityCoefficients(doublereal* ac) const ac[k] = 1.0; } } -//================================================================================================ -/* - * - * getChemPotentials(): - * - * This function returns a vector of chemical potentials of the - * species. - * \f[ - * \mu_k = \mu^o_k(T) + V_k * (p - p_o) + R T ln(X_k) - * \f] - * or another way to phrase this is - * \f[ - * \mu_k = \mu^o_k(T,p) + R T ln(X_k) - * \f] - * where \f$ \mu^o_k(T,p) = \mu^o_k(T) + V_k * (p - p_o)\f$ - * - */ + void IdealSolidSolnPhase:: getChemPotentials(doublereal* mu) const { @@ -654,24 +364,7 @@ getChemPotentials(doublereal* mu) const + delta_p * m_speciesMolarVolume[k]; } } -//================================================================================================ -/* - * - * getChemPotentials_RT() - * - * Get the array of non-dimensional chemical potentials \f$ - * \mu_k / \hat R T \f$, where - * - * \f[ - * \mu_k = \mu^o_k(T) + V_k * (p - p_o) + R T ln(X_k) - * \f] - * or another way to phrase this is - * \f[ - * \mu_k = \mu^o_k(T,p) + R T ln(X_k) - * \f] - * where \f$ \mu^o_k(T,p) = \mu^o_k(T) + V_k * (p - p_o)\f$ - * - */ + void IdealSolidSolnPhase:: getChemPotentials_RT(doublereal* mu) const { @@ -690,21 +383,6 @@ getChemPotentials_RT(doublereal* mu) const * Partial Molar Properties ********************************************************************/ -/******************************************************************** - * - * getPartialMolarEnthalpies() - * - * For this phase, the partial molar enthalpies are equal to the - * pure species enthalpies. - * \f[ - * \hat h_k(T,P) = \sum_k X_k \hat h^0_k(T) + (p - p_{ref}) (\sum_k X_k \hat V^0_k) - * \f] - * The reference-state pure-species enthalpies at the reference - * pressure p_ref - * \f$ \hat h^0_k(T) \f$, are computed by the species thermodynamic - * property manager. They are polynomial functions of temperature. - * @see SpeciesThermo - */ void IdealSolidSolnPhase::getPartialMolarEnthalpies(doublereal* hbar) const { const vector_fp& _h = enthalpy_RT_ref(); @@ -712,23 +390,6 @@ void IdealSolidSolnPhase::getPartialMolarEnthalpies(doublereal* hbar) const scale(_h.begin(), _h.end(), hbar, rt); } -/******************************************************************** - * - * getPartialMolarEntropies() - * - * Returns an array of partial molar entropies of the species in the - * solution. Units: J/kmol. - * For this phase, the partial molar entropies are equal to the - * pure species entropies plus the ideal solution contribution. - * \f[ - * \bar s_k(T,P) = \hat s^0_k(T) - R log(X_k) - * \f] - * The reference-state pure-species entropies,\f$ \hat s^0_k(T) \f$, - * at the reference pressure, \f$ P_{ref} \f$, are computed by the - * species thermodynamic - * property manager. They are polynomial functions of temperature. - * @see SpeciesThermo - */ void IdealSolidSolnPhase:: getPartialMolarEntropies(doublereal* sbar) const { @@ -741,14 +402,6 @@ getPartialMolarEntropies(doublereal* sbar) const } } -/******************************************************************** - * - * getPartialMolarCp() - * - * For this phase, the partial molar heat capacities are equal - * to the standard state heat capacities. - - */ void IdealSolidSolnPhase:: getPartialMolarCp(doublereal* cpbar) const { @@ -758,16 +411,6 @@ getPartialMolarCp(doublereal* cpbar) const } } -/****************************************************************** - * - * getPartialMolarVolumes() - * - * returns an array of partial molar volumes of the species - * in the solution. Units: m^3 kmol-1. - * - * For this solution, thepartial molar volumes are equal to the - * constant species molar volumes. - */ void IdealSolidSolnPhase:: getPartialMolarVolumes(doublereal* vbar) const { @@ -775,25 +418,9 @@ getPartialMolarVolumes(doublereal* vbar) const } /***************************************************************** - * Properties of the Standard State of the Species - * in the Solution + * Properties of the Standard State of the Species in the Solution *****************************************************************/ -/****************************************************************** - * - * getPureGibbs() - * - * Get the Gibbs functions for the pure species - * at the current T and P of the solution. - * We assume an incompressible constant partial molar - * volume here: - * \f[ - * \mu^0_k(T,p) = \mu^{ref}_k(T) + (P - P_{ref}) * V_k - * \f] - * where \f$V_k\f$ is the molar volume of pure species k<\I>. - * \f$ u^{ref}_k(T)\f$ is the chemical potential of pure - * species k<\I> at the reference pressure, \f$P_{ref}\f$. - */ void IdealSolidSolnPhase:: getPureGibbs(doublereal* gpure) const { @@ -806,21 +433,6 @@ getPureGibbs(doublereal* gpure) const } } -/** - * Get the nondimensional gibbs function for the species - * standard states at the current T and P of the solution. - * - * \f[ - * \mu^0_k(T,P) = \mu^{ref}_k(T) + (P - P_{ref}) * V_k - * \f] - * where \f$V_k\f$ is the molar volume of pure species k. - * \f$ \mu^{ref}_k(T)\f$ is the chemical potential of pure - * species k at the reference pressure, \f$P_{ref}\f$. - * - * @param grt Vector of length m_kk, which on return sr[k] - * will contain the nondimensional - * standard state gibbs function for species k. - */ void IdealSolidSolnPhase:: getGibbs_RT(doublereal* grt) const { @@ -833,21 +445,6 @@ getGibbs_RT(doublereal* grt) const } } -/******************************************************************** - * - * getEnthalpy_RT() - * - * Get the array of nondimensional Enthalpy functions for the ss - * species at the current T and P of the solution. - * We assume an incompressible constant partial molar - * volume here: - * \f[ - * h^0_k(T,P) = h^{ref}_k(T) + (P - P_{ref}) * V_k - * \f] - * where \f$V_k\f$ is the molar volume of pure species k<\I>. - * \f$ h^{ref}_k(T)\f$ is the enthalpy of the pure - * species k<\I> at the reference pressure, \f$P_{ref}\f$. - */ void IdealSolidSolnPhase:: getEnthalpy_RT(doublereal* hrt) const { @@ -859,36 +456,12 @@ getEnthalpy_RT(doublereal* hrt) const } } -/** - * Get the nondimensional Entropies for the species - * standard states at the current T and P of the solution. - * - * Note, this is equal to the reference state entropies - * due to the zero volume expansivity: - * i.e., (dS/dp)_T = (dV/dT)_P = 0.0 - * - * @param sr Vector of length m_kk, which on return sr[k] - * will contain the nondimensional - * standard state entropy of species k. - */ void IdealSolidSolnPhase::getEntropy_R(doublereal* sr) const { const vector_fp& _s = entropy_R_ref(); copy(_s.begin(), _s.end(), sr); } -/* - * Returns the vector of nondimensional - * internal Energies of the standard state at the current temperature - * of the solution and current pressure for each species. - * \f[ - * u^0_k(T,P) = h^{ref}_k(T) - P_{ref} * V_k - * \f] - * - * The standard state internal energy is independent of - * pressure in this equation of state. - * (inherited from ThermoPhase.h) - */ void IdealSolidSolnPhase::getIntEnergy_RT(doublereal* urt) const { const vector_fp& _h = enthalpy_RT_ref(); @@ -898,50 +471,21 @@ void IdealSolidSolnPhase::getIntEnergy_RT(doublereal* urt) const } } -/* - * Get the nondimensional heat capacity at constant pressure - * function for the species - * standard states at the current T and P of the solution. - * - * \f[ - * Cp^0_k(T,P) = Cp^{ref}_k(T) - * \f] - * where \f$V_k\f$ is the molar volume of pure species k<\I>. - * \f$ Cp^{ref}_k(T)\f$ is the constant pressure heat capacity - * of species k<\I> at the reference pressure, \f$P_{ref}\f$. - * - * @param cpr Vector of length m_kk, which on return cpr[k] - * will contain the nondimensional - * constant pressure heat capacity for species k. - */ void IdealSolidSolnPhase::getCp_R(doublereal* cpr) const { const vector_fp& _cpr = cp_R_ref(); copy(_cpr.begin(), _cpr.end(), cpr); } -/* - * Get the molar volumes of each species in their standard - * states at the current - * T and P of the solution. - * units = m^3 / kmol - */ void IdealSolidSolnPhase::getStandardVolumes(doublereal* vol) const { copy(m_speciesMolarVolume.begin(), m_speciesMolarVolume.end(), vol); } - /********************************************************************* * Thermodynamic Values for the Species Reference States *********************************************************************/ -/* - * Returns the vector of non-dimensional Enthalpy function - * of the reference state at the current temperature - * of the solution and the reference pressure for the species. - * Units = unitless - */ void IdealSolidSolnPhase::getEnthalpy_RT_ref(doublereal* hrt) const { _updateThermo(); @@ -950,12 +494,6 @@ void IdealSolidSolnPhase::getEnthalpy_RT_ref(doublereal* hrt) const } } -/* - * Returns the vector of non-dimensional Gibbs function - * of the reference state at the current temperature - * of the solution and the reference pressure for the species. - * Units = unitless - */ void IdealSolidSolnPhase::getGibbs_RT_ref(doublereal* grt) const { _updateThermo(); @@ -964,12 +502,6 @@ void IdealSolidSolnPhase::getGibbs_RT_ref(doublereal* grt) const } } -/* - * Returns the vector of Gibbs function - * of the reference state at the current temperature - * of the solution and the reference pressure for the species. - * Units = J / kmol - */ void IdealSolidSolnPhase::getGibbs_ref(doublereal* g) const { _updateThermo(); @@ -979,12 +511,6 @@ void IdealSolidSolnPhase::getGibbs_ref(doublereal* g) const } } -/* - * Returns the vector of nondimensional - * internal Energies of the standard state at the current temperature - * of the solution and current pressure for each species. - * (inherited from ThermoPhase.h) - */ void IdealSolidSolnPhase::getIntEnergy_RT_ref(doublereal* urt) const { const vector_fp& _h = enthalpy_RT_ref(); @@ -994,12 +520,6 @@ void IdealSolidSolnPhase::getIntEnergy_RT_ref(doublereal* urt) const } } -/* - * Returns the vector of non-dimensional Entropy function - * of the reference state at the current temperature - * of the solution and the reference pressure for the species. - * Units = unitless - */ void IdealSolidSolnPhase::getEntropy_R_ref(doublereal* er) const { _updateThermo(); @@ -1008,12 +528,6 @@ void IdealSolidSolnPhase::getEntropy_R_ref(doublereal* er) const } } -/* - * Returns the vector of non-dimensional Entropy function - * of the reference state at the current temperature - * of the solution and the reference pressure for the species. - * Units = unitless - */ void IdealSolidSolnPhase::getCp_R_ref(doublereal* cpr) const { _updateThermo(); @@ -1022,26 +536,12 @@ void IdealSolidSolnPhase::getCp_R_ref(doublereal* cpr) const } } -/* - * Returns a reference to the vector of nondimensional - * enthalpies of the reference state at the current temperature. - * Real reason for its existence is that it also checks - * to see if a recalculation of the reference thermodynamics - * functions needs to be done. - */ const vector_fp& IdealSolidSolnPhase::enthalpy_RT_ref() const { _updateThermo(); return m_h0_RT; } -/* - * Returns a reference to the vector of nondimensional - * enthalpies of the reference state at the current temperature. - * Real reason for its existence is that it also checks - * to see if a recalculation of the reference thermodynamics - * functions needs to be done. - */ const vector_fp& IdealSolidSolnPhase::entropy_R_ref() const { _updateThermo(); @@ -1051,41 +551,11 @@ const vector_fp& IdealSolidSolnPhase::entropy_R_ref() const /********************************************************************* * Utility Functions *********************************************************************/ -/* - * initThermo() function initializes the object for use. - * - * Before its invocation, the class isn't ready for calculation. - */ + void IdealSolidSolnPhase::initThermo() { } -/* - * @internal - * Import and initialize a ThermoPhase object - * using an XML tree. - * Here we read extra information about the XML description - * of a phase. Regular information about elements and species - * and their reference state thermodynamic information - * have already been read at this point. - * For example, we do not need to call this function for - * ideal gas equations of state. - * This function is called from importPhase() - * after the elements and the - * species are initialized with default ideal solution - * level data. - * - * @param phaseNode This object must be the phase node of a - * complete XML tree - * description of the phase, including all of the - * species data. In other words while "phase" must - * point to an XML phase object, it must have - * sibling nodes "speciesData" that describe - * the species in the phase. - * @param id ID of the phase. If nonnull, a check is done - * to see if phaseNode is pointing to the phase - * with the correct id. - */ void IdealSolidSolnPhase::initThermoXML(XML_Node& phaseNode, const std::string& id) { string subname = "IdealSolidSolnPhase::initThermoXML"; @@ -1166,9 +636,6 @@ void IdealSolidSolnPhase::initThermoXML(XML_Node& phaseNode, const std::string& ThermoPhase::initThermoXML(phaseNode, id); } -/* - * This internal function adjusts the lengths of arrays - */ void IdealSolidSolnPhase:: initLengths() { @@ -1190,14 +657,6 @@ initLengths() m_speciesMolarVolume.resize(m_kk); } -/* - * Set mixture to an equilibrium state consistent with specified - * element potentials and temperature. - * - * @param lambda_RT vector of non-dimensional element potentials - * \f$ \lambda_m/RT \f$. - * - */ void IdealSolidSolnPhase:: setToEquilState(const doublereal* lambda_RT) { @@ -1217,45 +676,19 @@ setToEquilState(const doublereal* lambda_RT) doublereal* dptr = DATA_PTR(m_pp); setState_PX(pres, dptr); } -//================================================================================================ -/* - * - * speciesMolarVolume() - * - * Report the molar volume of species k - * - * units - \f$ m^3 kmol^-1 \f$ - */ + double IdealSolidSolnPhase:: speciesMolarVolume(int k) const { return m_speciesMolarVolume[k]; } -/* - * - * getSpeciesMolarVolumes(): - * - * Fill in a return vector containing the species molar volumes - * units - \f$ m^3 kmol^-1 \f$ - */ void IdealSolidSolnPhase:: getSpeciesMolarVolumes(doublereal* smv) const { copy(m_speciesMolarVolume.begin(), m_speciesMolarVolume.end(), smv); } -//================================================================================================ -/* - * - * _updateThermo() - * - * This function gets called for every call to functions in this - * class. It checks to see whether the temperature has changed and - * thus the reference thermodynamics functions for all of the species - * must be recalculated. - * If the temperature has changed, the species thermo manager is called - * to recalculate G, Cp, H, and S at the current temperature. - */ + void IdealSolidSolnPhase:: _updateThermo() const { @@ -1277,6 +710,5 @@ _updateThermo() const m_tlast = tnow; } } -//================================================================================================ + } // end namespace Cantera -//================================================================================================== diff --git a/src/thermo/LatticePhase.cpp b/src/thermo/LatticePhase.cpp index c8d4edfc5..66ada2edd 100644 --- a/src/thermo/LatticePhase.cpp +++ b/src/thermo/LatticePhase.cpp @@ -21,7 +21,6 @@ namespace Cantera { -// Base Empty constructor LatticePhase::LatticePhase() : m_Pref(OneAtm), m_Pcurrent(OneAtm), @@ -31,10 +30,6 @@ LatticePhase::LatticePhase() : { } -// Copy Constructor -/* - * @param right Object to be copied - */ LatticePhase::LatticePhase(const LatticePhase& right) : m_Pref(OneAtm), m_Pcurrent(OneAtm), @@ -45,10 +40,6 @@ LatticePhase::LatticePhase(const LatticePhase& right) : *this = operator=(right); } -// Assignment operator -/* - * @param right Object to be copied - */ LatticePhase& LatticePhase::operator=(const LatticePhase& right) { if (&right != this) { @@ -67,41 +58,20 @@ LatticePhase& LatticePhase::operator=(const LatticePhase& right) return *this; } -// Destructor LatticePhase::~LatticePhase() { } - -// Full constructor for a lattice phase -/* - * @param inputFile String name of the input file - * @param id string id of the phase name - */ LatticePhase::LatticePhase(const std::string& inputFile, const std::string& id) { initThermoFile(inputFile, id); } -// Full constructor for a water phase -/* - * @param phaseRef XML node referencing the lattice phase. - * @param id string id of the phase name - */ LatticePhase::LatticePhase(XML_Node& phaseRef, const std::string& id) { importPhase(*findXMLPhase(&phaseRef, id), this); } - -// Duplication function -/* - * This virtual function is used to create a duplicate of the - * current phase. It's used to duplicate the phase when given - * a ThermoPhase pointer to the phase. - * - * @return It returns a ThermoPhase pointer. - */ ThermoPhase* LatticePhase::duplMyselfAsThermoPhase() const { return new LatticePhase(*this); @@ -128,22 +98,22 @@ doublereal LatticePhase::entropy_mole() const return GasConstant * (mean_X(&entropy_R_ref()[0]) - sum_xlogx()); } -//==================================================================================================================== + doublereal LatticePhase::gibbs_mole() const { return enthalpy_mole() - temperature() * entropy_mole(); } -//==================================================================================================================== + doublereal LatticePhase::cp_mole() const { return GasConstant * mean_X(&cp_R_ref()[0]); } -//==================================================================================================================== + doublereal LatticePhase::cv_mole() const { return cp_mole(); } -//==================================================================================================================== + doublereal LatticePhase::calcDensity() { setMolarDensity(m_site_density); @@ -162,65 +132,65 @@ doublereal LatticePhase::calcDensity() // Phase::setDensity(dens); return dens; } -//==================================================================================================================== + void LatticePhase::setPressure(doublereal p) { m_Pcurrent = p; calcDensity(); } -//==================================================================================================================== + void LatticePhase::setMoleFractions(const doublereal* const x) { Phase::setMoleFractions(x); calcDensity(); } -//==================================================================================================================== + void LatticePhase::setMoleFractions_NoNorm(const doublereal* const x) { Phase::setMoleFractions(x); calcDensity(); } -//==================================================================================================================== + void LatticePhase::setMassFractions(const doublereal* const y) { Phase::setMassFractions(y); calcDensity(); } -//==================================================================================================================== + void LatticePhase::setMassFractions_NoNorm(const doublereal* const y) { Phase::setMassFractions_NoNorm(y); calcDensity(); } -//==================================================================================================================== + void LatticePhase::setConcentrations(const doublereal* const c) { Phase::setConcentrations(c); calcDensity(); } -//==================================================================================================================== + void LatticePhase::getActivityConcentrations(doublereal* c) const { getMoleFractions(c); } -//==================================================================================================================== + void LatticePhase::getActivityCoefficients(doublereal* ac) const { for (size_t k = 0; k < m_kk; k++) { ac[k] = 1.0; } } -//==================================================================================================================== + doublereal LatticePhase::standardConcentration(size_t k) const { return 1.0; } -//==================================================================================================================== + doublereal LatticePhase::logStandardConc(size_t k) const { return 0.0; } -//==================================================================================================================== + void LatticePhase::getChemPotentials(doublereal* mu) const { doublereal delta_p = m_Pcurrent - m_Pref; @@ -234,14 +204,14 @@ void LatticePhase::getChemPotentials(doublereal* mu) const } } -//==================================================================================================================== + void LatticePhase::getPartialMolarEnthalpies(doublereal* hbar) const { const vector_fp& _h = enthalpy_RT_ref(); doublereal rt = GasConstant * temperature(); scale(_h.begin(), _h.end(), hbar, rt); } -//==================================================================================================================== + void LatticePhase::getPartialMolarEntropies(doublereal* sbar) const { const vector_fp& _s = entropy_R_ref(); @@ -252,7 +222,7 @@ void LatticePhase::getPartialMolarEntropies(doublereal* sbar) const sbar[k] = r * (_s[k] - log(xx)); } } -//==================================================================================================================== + void LatticePhase::getPartialMolarCp(doublereal* cpbar) const { getCp_R(cpbar); @@ -260,18 +230,18 @@ void LatticePhase::getPartialMolarCp(doublereal* cpbar) const cpbar[k] *= GasConstant; } } -//==================================================================================================================== + void LatticePhase::getPartialMolarVolumes(doublereal* vbar) const { getStandardVolumes(vbar); } -//==================================================================================================================== + void LatticePhase::getStandardChemPotentials(doublereal* mu0) const { const vector_fp& gibbsrt = gibbs_RT_ref(); scale(gibbsrt.begin(), gibbsrt.end(), mu0, _RT()); } -//==================================================================================================================== + void LatticePhase::getPureGibbs(doublereal* gpure) const { const vector_fp& gibbsrt = gibbs_RT_ref(); @@ -281,7 +251,7 @@ void LatticePhase::getPureGibbs(doublereal* gpure) const gpure[k] = RT * gibbsrt[k] + delta_p * m_speciesMolarVolume[k]; } } -//==================================================================================================================== + void LatticePhase::getEnthalpy_RT(doublereal* hrt) const { const vector_fp& _h = enthalpy_RT_ref(); @@ -290,13 +260,13 @@ void LatticePhase::getEnthalpy_RT(doublereal* hrt) const hrt[k] = _h[k] + delta_prt * m_speciesMolarVolume[k]; } } -//==================================================================================================================== + void LatticePhase::getEntropy_R(doublereal* sr) const { const vector_fp& _s = entropy_R_ref(); std::copy(_s.begin(), _s.end(), sr); } -//==================================================================================================================== + void LatticePhase::getGibbs_RT(doublereal* grt) const { const vector_fp& gibbsrt = gibbs_RT_ref(); @@ -306,7 +276,7 @@ void LatticePhase::getGibbs_RT(doublereal* grt) const grt[k] = gibbsrt[k] + delta_prt * m_speciesMolarVolume[k]; } } -//==================================================================================================================== + void LatticePhase::getGibbs_ref(doublereal* g) const { getGibbs_RT_ref(g); @@ -314,42 +284,30 @@ void LatticePhase::getGibbs_ref(doublereal* g) const g[k] *= GasConstant * temperature(); } } -//=================================================================================================================== void LatticePhase::getCp_R(doublereal* cpr) const { const vector_fp& _cpr = cp_R_ref(); std::copy(_cpr.begin(), _cpr.end(), cpr); } -//=================================================================================================================== + void LatticePhase::getStandardVolumes(doublereal* vbar) const { copy(m_speciesMolarVolume.begin(), m_speciesMolarVolume.end(), vbar); } -//======================================================================================================= -// Returns the vector of nondimensional Enthalpies of the reference state at the current temperature -// of the solution and the reference pressure for the phase. -/* - * @return Output vector of nondimensional reference state Enthalpies of the species. - * Length: m_kk - */ + const vector_fp& LatticePhase::enthalpy_RT_ref() const { _updateThermo(); return m_h0_RT; } -//======================================================================================================= -// Returns a reference to the dimensionless reference state Gibbs free energy vector. -/* - * This function is part of the layer that checks/recalculates the reference - * state thermo functions. - */ + const vector_fp& LatticePhase::gibbs_RT_ref() const { _updateThermo(); return m_g0_RT; } -//==================================================================================================================== + void LatticePhase::getGibbs_RT_ref(doublereal* grt) const { _updateThermo(); @@ -357,43 +315,19 @@ void LatticePhase::getGibbs_RT_ref(doublereal* grt) const grt[k] = m_g0_RT[k]; } } -//======================================================================================================= -// Returns a reference to the dimensionless reference state Entropy vector. -/* - * This function is part of the layer that checks/recalculates the reference - * state thermo functions. - */ + const vector_fp& LatticePhase::entropy_R_ref() const { _updateThermo(); return m_s0_R; } -//======================================================================================================= -// Returns a reference to the dimensionless reference state Heat Capacity vector. -/* - * This function is part of the layer that checks/recalculates the reference - * state thermo functions. - */ + const vector_fp& LatticePhase::cp_R_ref() const { _updateThermo(); return m_cp0_R; } -//==================================================================================================================== -// Initialize the ThermoPhase object after all species have been set up -/* - * @internal Initialize. - * - * This method performs any initialization required after all - * species have been added. For example, it is used to - * resize internal work arrays that must have an entry for - * each species. - * This method is called from ThermoPhase::initThermoXML(), - * which is called from importPhase(), - * just prior to returning from the function, importPhase(). - * - * @see importCTML.cpp - */ + void LatticePhase::initThermo() { m_Pref = refPressure(); @@ -406,7 +340,7 @@ void LatticePhase::initThermo() ThermoPhase::initThermo(); } -//==================================================================================================================== + void LatticePhase::initThermoXML(XML_Node& phaseNode, const std::string& id) { std::string idattrib = phaseNode.id(); @@ -459,12 +393,7 @@ void LatticePhase::initThermoXML(XML_Node& phaseNode, const std::string& id) */ ThermoPhase::initThermoXML(phaseNode, id); } -//===================================================================================================== -// Update the species reference state thermodynamic functions -/* - * The polynomials for the standard state functions are only - * reevaluated if the temperature has changed. - */ + void LatticePhase::_updateThermo() const { doublereal tnow = temperature(); @@ -477,26 +406,25 @@ void LatticePhase::_updateThermo() const m_tlast = tnow; } } -//===================================================================================================== + void LatticePhase::setParameters(int n, doublereal* const c) { m_site_density = c[0]; setMolarDensity(m_site_density); } -//===================================================================================================== + void LatticePhase::getParameters(int& n, doublereal* const c) const { double d = molarDensity(); c[0] = d; n = 1; } -//===================================================================================================== + void LatticePhase::setParametersFromXML(const XML_Node& eosdata) { eosdata._require("model", "Lattice"); m_site_density = ctml::getFloat(eosdata, "site_density", "toSI"); m_vacancy = ctml::getChildValue(eosdata, "vacancy_species"); } -//===================================================================================================== + } -//======================================================================================================= diff --git a/src/thermo/LatticeSolidPhase.cpp b/src/thermo/LatticeSolidPhase.cpp index fdc8a683d..d2de017a6 100644 --- a/src/thermo/LatticeSolidPhase.cpp +++ b/src/thermo/LatticeSolidPhase.cpp @@ -18,12 +18,9 @@ #include using namespace std; -//====================================================================================================================== + namespace Cantera { - -//==================================================================================================================== -// Base empty constructor LatticeSolidPhase::LatticeSolidPhase() : m_tlast(0.0), m_press(-1.0), @@ -35,11 +32,7 @@ LatticeSolidPhase::LatticeSolidPhase() : tmpV_(0) { } -//==================================================================================================================== -// Copy Constructor -/* - * @param right Object to be copied - */ + LatticeSolidPhase::LatticeSolidPhase(const LatticeSolidPhase& right) : m_tlast(0.0), m_press(-1.0), @@ -52,11 +45,7 @@ LatticeSolidPhase::LatticeSolidPhase(const LatticeSolidPhase& right) : { *this = operator=(right); } -//==================================================================================================================== -// Assignment operator -/* - * @param right Object to be copied - */ + LatticeSolidPhase& LatticeSolidPhase::operator=(const LatticeSolidPhase& right) { @@ -73,8 +62,7 @@ LatticeSolidPhase::operator=(const LatticeSolidPhase& right) } return *this; } -//==================================================================================================================== -// Destructor + LatticeSolidPhase::~LatticeSolidPhase() { // We own the sublattices. So we have to delete the sublattices @@ -83,33 +71,12 @@ LatticeSolidPhase::~LatticeSolidPhase() m_lattice[n] = 0; } } -//==================================================================================================================== -// Duplication function -/* - * This virtual function is used to create a duplicate of the - * current phase. It's used to duplicate the phase when given - * a ThermoPhase pointer to the phase. - * - * @return It returns a %ThermoPhase pointer. - */ + ThermoPhase* LatticeSolidPhase::duplMyselfAsThermoPhase() const { return new LatticeSolidPhase(*this); } -//==================================================================================================================== -// Minimum temperature for which the thermodynamic data for the species -// or phase are valid. -/* - * If no argument is supplied, the - * value returned will be the lowest temperature at which the - * data for \e all species are valid. Otherwise, the value - * will be only for species \a k. This function is a wrapper - * that calls the species thermo minTemp function. - * - * @param k index of the species. Default is -1, which will return the max of the min value - * over all species. - */ doublereal LatticeSolidPhase::minTemp(size_t k) const { if (k != npos) { @@ -127,19 +94,7 @@ doublereal LatticeSolidPhase::minTemp(size_t k) const } return mm; } -//==================================================================================================================== -// Maximum temperature for which the thermodynamic data for the species -// or phase are valid. -/* - * If no argument is supplied, the - * value returned will be the lowest temperature at which the - * data for \e all species are valid. Otherwise, the value - * will be only for species \a k. This function is a wrapper - * that calls the species thermo minTemp function. - * - * @param k index of the species. Default is -1, which will return the max of the min value - * over all species. - */ + doublereal LatticeSolidPhase::maxTemp(size_t k) const { if (k != npos) { @@ -157,18 +112,13 @@ doublereal LatticeSolidPhase::maxTemp(size_t k) const } return mm; } -//==================================================================================================================== -/* - * Returns the reference pressure in Pa. This function is a wrapper - * that calls the species thermo refPressure function. - */ + doublereal LatticeSolidPhase::refPressure() const { return m_lattice[0]->refPressure(); } -//==================================================================================================================== -doublereal LatticeSolidPhase:: -enthalpy_mole() const + +doublereal LatticeSolidPhase::enthalpy_mole() const { _updateThermo(); doublereal sum = 0.0; @@ -177,7 +127,7 @@ enthalpy_mole() const } return sum; } -//==================================================================================================================== + doublereal LatticeSolidPhase::intEnergy_mole() const { _updateThermo(); @@ -187,7 +137,7 @@ doublereal LatticeSolidPhase::intEnergy_mole() const } return sum; } -//==================================================================================================================== + doublereal LatticeSolidPhase::entropy_mole() const { _updateThermo(); @@ -197,7 +147,7 @@ doublereal LatticeSolidPhase::entropy_mole() const } return sum; } -//==================================================================================================================== + doublereal LatticeSolidPhase::gibbs_mole() const { _updateThermo(); @@ -207,7 +157,7 @@ doublereal LatticeSolidPhase::gibbs_mole() const } return sum; } -//==================================================================================================================== + doublereal LatticeSolidPhase::cp_mole() const { _updateThermo(); @@ -217,7 +167,7 @@ doublereal LatticeSolidPhase::cp_mole() const } return sum; } -//==================================================================================================================== + void LatticeSolidPhase::getActivityConcentrations(doublereal* c) const { _updateThermo(); @@ -227,30 +177,24 @@ void LatticeSolidPhase::getActivityConcentrations(doublereal* c) const strt += m_lattice[n]->nSpecies(); } } -//==================================================================================================================== + void LatticeSolidPhase::getActivityCoefficients(doublereal* ac) const { for (size_t k = 0; k < m_kk; k++) { ac[k] = 1.0; } } -//==================================================================================================================== + doublereal LatticeSolidPhase::standardConcentration(size_t k) const { return 1.0; } -//==================================================================================================================== + doublereal LatticeSolidPhase::logStandardConc(size_t k) const { return 0.0; } -//==================================================================================================================== -// Set the pressure at constant temperature. Units: Pa. -/* - * - * @param p Pressure (units - Pa) - */ void LatticeSolidPhase::setPressure(doublereal p) { m_press = p; @@ -259,19 +203,7 @@ void LatticeSolidPhase::setPressure(doublereal p) } calcDensity(); } -//==================================================================================================================== -// Calculate the density of the solid mixture -/* - * The formula for this is - * - * \f[ - * \rho = \sum_n{ \rho_n \theta_n } - * \f] - * - * where \f$ \rho_n \f$ is the density of the nth sublattice - * - * Note this is a nonvirtual function. - */ + doublereal LatticeSolidPhase::calcDensity() { double sum = 0.0; @@ -281,20 +213,7 @@ doublereal LatticeSolidPhase::calcDensity() Phase::setDensity(sum); return sum; } -//==================================================================================================================== -// Set the mole fractions to the specified values, and then -// normalize them so that they sum to 1.0 for each of the subphases -/* - * On input, the mole fraction vector is assumed to sum to one for each of the sublattices. The sublattices - * are updated with this mole fraction vector. The mole fractions are also stored within this object, after - * they are normalized to one by dividing by the number of sublattices. - * - * @param x Input vector of mole fractions. There is no restriction - * on the sum of the mole fraction vector. Internally, - * this object will pass portions of this vector to the sublattices which assume that the portions - * individually sum to one. - * Length is m_kk. - */ + void LatticeSolidPhase::setMoleFractions(const doublereal* const x) { size_t nsp, strt = 0; @@ -309,14 +228,7 @@ void LatticeSolidPhase::setMoleFractions(const doublereal* const x) Phase::setMoleFractions(DATA_PTR(m_x)); calcDensity(); } -//==================================================================================================================== -// Get the species mole fraction vector. -/* - * On output the mole fraction vector will sum to one for each of the subphases which make up this phase. - * - * @param x On return, x contains the mole fractions. Must have a - * length greater than or equal to the number of species. - */ + void LatticeSolidPhase::getMoleFractions(doublereal* const x) const { size_t nsp, strt = 0; @@ -348,18 +260,7 @@ void LatticeSolidPhase::getMoleFractions(doublereal* const x) const strt += nsp; } } -//==================================================================================================================== -// Get the species chemical potentials. Units: J/kmol. -/* - * This function returns a vector of chemical potentials of the - * species in solution at the current temperature, pressure - * and mole fraction of the solution. - * - * This returns the underlying lattice chemical potentials - * - * @param mu Output vector of species chemical - * potentials. Length: m_kk. Units: J/kmol - */ + void LatticeSolidPhase::getChemPotentials(doublereal* mu) const { _updateThermo(); @@ -370,7 +271,7 @@ void LatticeSolidPhase::getChemPotentials(doublereal* mu) const strt += nlsp; } } -//==================================================================================================================== + void LatticeSolidPhase::getPartialMolarEnthalpies(doublereal* hbar) const { _updateThermo(); @@ -381,7 +282,7 @@ void LatticeSolidPhase::getPartialMolarEnthalpies(doublereal* hbar) const strt += nlsp; } } -//==================================================================================================================== + void LatticeSolidPhase::getPartialMolarEntropies(doublereal* sbar) const { _updateThermo(); @@ -392,7 +293,7 @@ void LatticeSolidPhase::getPartialMolarEntropies(doublereal* sbar) const strt += nlsp; } } -//==================================================================================================================== + void LatticeSolidPhase::getPartialMolarCp(doublereal* cpbar) const { _updateThermo(); @@ -403,7 +304,7 @@ void LatticeSolidPhase::getPartialMolarCp(doublereal* cpbar) const strt += nlsp; } } -//==================================================================================================================== + void LatticeSolidPhase::getPartialMolarVolumes(doublereal* vbar) const { _updateThermo(); @@ -414,20 +315,7 @@ void LatticeSolidPhase::getPartialMolarVolumes(doublereal* vbar) const strt += nlsp; } } -//==================================================================================================================== -// Get the array of standard state 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. - * - * This returns the underlying lattice standard chemical potentials, as the units are kmol-1 of - * the sublattice species. - * - * @param mu0 Output vector of chemical potentials. - * Length: m_kk. Units: J/kmol - */ + void LatticeSolidPhase::getStandardChemPotentials(doublereal* mu0) const { _updateThermo(); @@ -437,7 +325,7 @@ void LatticeSolidPhase::getStandardChemPotentials(doublereal* mu0) const strt += m_lattice[n]->nSpecies(); } } -//==================================================================================================================== + void LatticeSolidPhase::getGibbs_RT_ref(doublereal* grt) const { _updateThermo(); @@ -445,7 +333,7 @@ void LatticeSolidPhase::getGibbs_RT_ref(doublereal* grt) const m_lattice[n]->getGibbs_RT_ref(grt + lkstart_[n]); } } -//==================================================================================================================== + void LatticeSolidPhase::getGibbs_ref(doublereal* g) const { getGibbs_RT_ref(g); @@ -453,13 +341,7 @@ void LatticeSolidPhase::getGibbs_ref(doublereal* g) const g[k] *= GasConstant * temperature(); } } -//==================================================================================================================== -// Add in species from Slave phases -/* - * This hook is used for cSS_CONVENTION_SLAVE phases - * - * @param phaseNode XML_Node for the current phase - */ + void LatticeSolidPhase::installSlavePhases(Cantera::XML_Node* phaseNode) { size_t kk = 0; @@ -537,22 +419,6 @@ void LatticeSolidPhase::installSlavePhases(Cantera::XML_Node* phaseNode) } } -//==================================================================================================================== -// Initialize the ThermoPhase object after all species have been set up -/* - * @internal Initialize. - * - * This method is provided to allow subclasses to perform any initialization required after all - * species have been added. For example, it might be used to - * resize internal work arrays that must have an entry for - * each species. The base class implementation does nothing, - * and subclasses that do not require initialization do not - * need to overload this method. When importing a CTML phase - * description, this method is called from ThermoPhase::initThermoXML(), - * which is called from importPhase(), just prior to returning from function importPhase(). - * - * @see importCTML.cpp - */ void LatticeSolidPhase::initThermo() { initLengths(); @@ -569,11 +435,7 @@ void LatticeSolidPhase::initThermo() setMoleFractions(DATA_PTR(m_x)); ThermoPhase::initThermo(); } -//==================================================================================================================== -// Initialize vectors that depend on the number of species and sublattices -/* - * - */ + void LatticeSolidPhase::initLengths() { theta_.resize(m_nlattice,0); @@ -581,7 +443,7 @@ void LatticeSolidPhase::initLengths() m_x.resize(m_kk, 0.0); tmpV_.resize(m_kk, 0.0); } -//==================================================================================================================== + void LatticeSolidPhase::_updateThermo() const { doublereal tnow = temperature(); @@ -601,7 +463,7 @@ void LatticeSolidPhase::_updateThermo() const m_tlast = tnow; } } -//==================================================================================================================== + void LatticeSolidPhase::setLatticeMoleFractionsByName(int nn, const std::string& x) { m_lattice[nn]->setMoleFractionsByName(x); @@ -617,15 +479,7 @@ void LatticeSolidPhase::setLatticeMoleFractionsByName(int nn, const std::string& } setMoleFractions(DATA_PTR(m_x)); } -//==================================================================================================================== - -//==================================================================================================================== -// Set the parameters from the XML file -/*! - * Currently, this is the spot that we read in all of the sublattice phases. - * The SetParametersFromXML() call is carried out at - */ void LatticeSolidPhase::setParametersFromXML(const XML_Node& eosdata) { eosdata._require("model","LatticeSolid"); @@ -661,18 +515,8 @@ void LatticeSolidPhase::setParametersFromXML(const XML_Node& eosdata) } } -//==================================================================================================================== #ifdef H298MODIFY_CAPABILITY - -//! Modify the value of the 298 K Heat of Formation of one species in the phase (J kmol-1) -/*! - * The 298K heat of formation is defined as the enthalpy change to create the standard state - * of the species from its constituent elements in their standard states at 298 K and 1 bar. - * - * @param k Species k - * @param Hf298New Specify the new value of the Heat of Formation at 298K and 1 bar - */ void LatticeSolidPhase::modifyOneHf298SS(const size_t& k, const doublereal Hf298New) { for (size_t n = 0; n < m_nlattice; n++) { @@ -686,7 +530,6 @@ void LatticeSolidPhase::modifyOneHf298SS(const size_t& k, const doublereal Hf298 _updateThermo(); } #endif -//==================================================================================================================== doublereal LatticeSolidPhase::err(const std::string& msg) const { @@ -695,4 +538,3 @@ doublereal LatticeSolidPhase::err(const std::string& msg) const } } // End namespace Cantera -//====================================================================================================================== diff --git a/src/thermo/MixtureFugacityTP.cpp b/src/thermo/MixtureFugacityTP.cpp index e81051abf..4a7652f8d 100644 --- a/src/thermo/MixtureFugacityTP.cpp +++ b/src/thermo/MixtureFugacityTP.cpp @@ -3,7 +3,6 @@ * Methods file for a derived class of ThermoPhase that handles * non-ideal mixtures based on the fugacity models (see \ref thermoprops and * class \link Cantera::MixtureFugacityTP MixtureFugacityTP\endlink). - * */ /* * Copyright (2005) Sandia Corporation. Under the terms of @@ -20,10 +19,7 @@ using namespace std; namespace Cantera { -//==================================================================================================================== -/* - * Default constructor - */ + MixtureFugacityTP::MixtureFugacityTP() : ThermoPhase(), m_Pcurrent(-1.0), @@ -38,16 +34,7 @@ MixtureFugacityTP::MixtureFugacityTP() : m_s0_R(0) { } -//==================================================================================================================== -/* - * Copy Constructor: - * - * Note this stuff will not work until the underlying phase - * has a working copy constructor. - * - * The copy constructor just calls the assignment operator - * to do the heavy lifting. - */ + MixtureFugacityTP::MixtureFugacityTP(const MixtureFugacityTP& b) : ThermoPhase(), m_Pcurrent(-1.0), @@ -63,13 +50,7 @@ MixtureFugacityTP::MixtureFugacityTP(const MixtureFugacityTP& b) : { MixtureFugacityTP::operator=(b); } -//==================================================================================================================== -/* - * operator=() - * - * Note this stuff will not work until the underlying phase - * has a working assignment operator - */ + MixtureFugacityTP& MixtureFugacityTP::operator=(const MixtureFugacityTP& b) { @@ -92,102 +73,43 @@ MixtureFugacityTP::operator=(const MixtureFugacityTP& b) m_cp0_R = b.m_cp0_R; m_g0_RT = b.m_g0_RT; m_s0_R = b.m_s0_R; - /* - * The VPSSMgr object contains shallow pointers. Whenever you have shallow - * pointers, they have to be fixed up to point to the correct objects referring - * back to this ThermoPhase's properties. - */ - //m_VPSS_ptr->initAllPtrs(this, m_spthermo); - /* - * The PDSS objects contains shallow pointers. Whenever you have shallow - * pointers, they have to be fixed up to point to the correct objects referring - * back to this ThermoPhase's properties. This function also sets m_VPSS_ptr - * so it occurs after m_VPSS_ptr is set. - */ - - /* - * Ok, the VPSSMgr object is ready for business. - * We need to resync the temperature and the pressure of the new standard states - * with what is stored in this object. - */ - // m_VPSS_ptr->setState_TP(m_Tlast_ss, m_Plast_ss); } return *this; } -//==================================================================================================================== -/* - * ~MixtureFugacityTP(): (virtual) - * - */ + MixtureFugacityTP::~MixtureFugacityTP() { - } -/* - * Duplication function. - * This calls the copy constructor for this object. - */ + ThermoPhase* MixtureFugacityTP::duplMyselfAsThermoPhase() const { return new MixtureFugacityTP(*this); } -//==================================================================================================================== -// This method returns the convention used in specification -// of the standard state, of which there are currently two, -// temperature based, and variable pressure based. -/* - * Currently, there are two standard state conventions: - * - Temperature-based activities - * cSS_CONVENTION_TEMPERATURE 0 - * - default - * - * - Variable Pressure and Temperature -based activities - * cSS_CONVENTION_VPSS 1 - */ + int MixtureFugacityTP::standardStateConvention() const { return cSS_CONVENTION_TEMPERATURE; } -//==================================================================================================================== -// Set the solution branch to force the ThermoPhase to exist on one branch or another -/* - * @param solnBranch Branch that the solution is restricted to. - * the value -1 means gas. The value -2 means unrestricted. - * Values of zero or greater refer to species dominated condensed phases. - */ + void MixtureFugacityTP::setForcedSolutionBranch(int solnBranch) { forcedState_ = solnBranch; } -//==================================================================================================================== -// Report the solution branch which the solution is restricted to -/* - * @return Branch that the solution is restricted to. - * the value -1 means gas. The value -2 means unrestricted. - * Values of zero or greater refer to species dominated condensed phases. - */ + int MixtureFugacityTP::forcedSolutionBranch() const { return forcedState_; } -//==================================================================================================================== -// Report the solution branch which the solution is actually on -/* - * @return Branch that the solution is restricted to. - * the value -1 means gas. The value -2 means superfluid.. - * Values of zero or greater refer to species dominated condensed phases. - */ + int MixtureFugacityTP::reportSolnBranchActual() const { return iState_; } -//==================================================================================================================== /* * ------------Molar Thermodynamic Properties ------------------------- */ -//==================================================================================================================== doublereal MixtureFugacityTP::err(const std::string& msg) const { @@ -195,20 +117,11 @@ doublereal MixtureFugacityTP::err(const std::string& msg) const +msg+" called. Equation of state type: "+int2str(eosType())); return 0; } -//==================================================================================================================== + /* * ---- Partial Molar Properties of the Solution ----------------- */ -//==================================================================================================================== -/* - * Get the array of non-dimensional species chemical potentials - * These are partial molar Gibbs free energies. - * \f$ \mu_k / \hat R T \f$. - * Units: unitless - * - * We close the loop on this function, here, calling - * getChemPotentials() and then dividing by RT. - */ + void MixtureFugacityTP::getChemPotentials_RT(doublereal* muRT) const { getChemPotentials(muRT); @@ -217,10 +130,11 @@ void MixtureFugacityTP::getChemPotentials_RT(doublereal* muRT) const muRT[k] *= invRT; } } -//==================================================================================================================== + /* * ----- Thermodynamic Values for the Species Standard States States ---- */ + void MixtureFugacityTP::getStandardChemPotentials(doublereal* g) const { _updateReferenceStateThermo(); @@ -231,33 +145,20 @@ void MixtureFugacityTP::getStandardChemPotentials(doublereal* g) const g[k] = RT * (g[k] + tmp); } } -//==================================================================================================================== + void MixtureFugacityTP::getEnthalpy_RT(doublereal* hrt) const { getEnthalpy_RT_ref(hrt); } -//================================================================================================ + #ifdef H298MODIFY_CAPABILITY -// Modify the value of the 298 K Heat of Formation of one species in the phase (J kmol-1) -/* - * The 298K heat of formation is defined as the enthalpy change to create the standard state - * of the species from its constituent elements in their standard states at 298 K and 1 bar. - * - * @param k Species k - * @param Hf298New Specify the new value of the Heat of Formation at 298K and 1 bar - */ void MixtureFugacityTP::modifyOneHf298SS(const int k, const doublereal Hf298New) { m_spthermo->modifyOneHf298(k, Hf298New); m_Tlast_ref += 0.0001234; } #endif -//==================================================================================================================== -/* - * Get the array of nondimensional entropy functions for the - * standard state species - * at the current T and P of the solution. - */ + void MixtureFugacityTP::getEntropy_R(doublereal* sr) const { _updateReferenceStateThermo(); @@ -267,11 +168,7 @@ void MixtureFugacityTP::getEntropy_R(doublereal* sr) const sr[k] -= tmp; } } -//==================================================================================================================== -/* - * Get the nondimensional gibbs function for the species - * standard states at the current T and P of the solution. - */ + void MixtureFugacityTP::getGibbs_RT(doublereal* grt) const { _updateReferenceStateThermo(); @@ -281,12 +178,7 @@ void MixtureFugacityTP::getGibbs_RT(doublereal* grt) const grt[k] += tmp; } } -//==================================================================================================================== -/* - * get the pure Gibbs free energies of each species assuming - * it is in its standard state. This is the same as - * getStandardChemPotentials(). - */ + void MixtureFugacityTP::getPureGibbs(doublereal* g) const { _updateReferenceStateThermo(); @@ -297,12 +189,7 @@ void MixtureFugacityTP::getPureGibbs(doublereal* g) const g[k] += tmp; } } -//==================================================================================================================== -/* - * Returns the vector of nondimensional - * internal Energies of the standard state at the current temperature - * and pressure of the solution for each species. - */ + void MixtureFugacityTP::getIntEnergy_RT(doublereal* urt) const { _updateReferenceStateThermo(); @@ -314,26 +201,13 @@ void MixtureFugacityTP::getIntEnergy_RT(doublereal* urt) const urt[i] -= tmp * v0; } } -//==================================================================================================================== -/* - * Get the nondimensional heat capacity at constant pressure - * function for the species - * standard states at the current T and P of the solution. - */ + void MixtureFugacityTP::getCp_R(doublereal* cpr) const { _updateReferenceStateThermo(); copy(m_cp0_R.begin(), m_cp0_R.end(), cpr); } -//==================================================================================================================== -/* - * Get the molar volumes of the species standard states at the current - * T and P of the solution. - * units = m^3 / kmol - * - * @param vol Output vector containing the standard state volumes. - * Length: m_kk. - */ + void MixtureFugacityTP::getStandardVolumes(doublereal* vol) const { _updateReferenceStateThermo(); @@ -342,84 +216,49 @@ void MixtureFugacityTP::getStandardVolumes(doublereal* vol) const vol[i]= v0; } } -//==================================================================================================================== + /* * ----- Thermodynamic Values for the Species Reference States ---- */ -/* - * Returns the vector of nondimensional enthalpies of the - * reference state at the current temperature of the solution and - * the reference pressure for the species. - */ + void MixtureFugacityTP::getEnthalpy_RT_ref(doublereal* hrt) const { _updateReferenceStateThermo(); copy(m_h0_RT.begin(), m_h0_RT.end(), hrt); } -//==================================================================================================================== -/* - * Returns the vector of nondimensional - * enthalpies of the reference state at the current temperature - * of the solution and the reference pressure for the species. - */ + void MixtureFugacityTP::getGibbs_RT_ref(doublereal* grt) const { _updateReferenceStateThermo(); copy(m_g0_RT.begin(), m_g0_RT.end(), grt); } -//==================================================================================================================== -/* - * Returns the vector of the - * gibbs function of the reference state at the current temperature - * of the solution and the reference pressure for the species. - * units = J/kmol - * - * This is filled in here so that derived classes don't have to - * take care of it. - */ + void MixtureFugacityTP::getGibbs_ref(doublereal* g) const { const vector_fp& gibbsrt = gibbs_RT_ref(); scale(gibbsrt.begin(), gibbsrt.end(), g, _RT()); } -//==================================================================================================================== + const vector_fp& MixtureFugacityTP::gibbs_RT_ref() const { _updateReferenceStateThermo(); return m_g0_RT; } -//==================================================================================================================== -/* - * Returns the vector of nondimensional - * entropies of the reference state at the current temperature - * of the solution and the reference pressure for the species. - */ + void MixtureFugacityTP::getEntropy_R_ref(doublereal* er) const { _updateReferenceStateThermo(); copy(m_s0_R.begin(), m_s0_R.end(), er); return; } -//==================================================================================================================== -/* - * Returns the vector of nondimensional - * constant pressure heat capacities of the reference state - * at the current temperature of the solution - * and reference pressure for the species. - */ + void MixtureFugacityTP::getCp_R_ref(doublereal* cpr) const { _updateReferenceStateThermo(); copy(m_cp0_R.begin(), m_cp0_R.end(), cpr); } -//==================================================================================================================== -/* - * Get the molar volumes of the species reference states at the current - * T and reference pressure of the solution. - * - * units = m^3 / kmol - */ + void MixtureFugacityTP::getStandardVolumes_ref(doublereal* vol) const { _updateReferenceStateThermo(); @@ -429,17 +268,7 @@ void MixtureFugacityTP::getStandardVolumes_ref(doublereal* vol) const vol[i]= v0; } } -//==================================================================================================================== -// Set the initial state of the phase to the conditions specified in the state XML element. -/* - * - * This method sets the temperature, pressure, and mole fraction vector to a set default value. - * We modify the default behavior here so that TP is evaluated at the same time. - * - * @param state AN XML_Node object corresponding to - * the "state" entry for this phase in the - * input file. - */ + void MixtureFugacityTP::setStateFromXML(const XML_Node& state) { int doTP = 0; @@ -472,21 +301,13 @@ void MixtureFugacityTP::setStateFromXML(const XML_Node& state) setState_TR(t, rho); } } -//==================================================================================================================== -/* - * Perform initializations after all species have been - * added. - */ + void MixtureFugacityTP::initThermo() { initLengths(); ThermoPhase::initThermo(); } -//==================================================================================================================== -/* - * Initialize the internal lengths. - * (this is not a virtual function) - */ + void MixtureFugacityTP::initLengths() { m_kk = nSpecies(); @@ -497,13 +318,13 @@ void MixtureFugacityTP::initLengths() m_g0_RT.resize(m_kk, 0.0); m_s0_R.resize(m_kk, 0.0); } -//==================================================================================================================== + void MixtureFugacityTP::setTemperature(const doublereal temp) { _updateReferenceStateThermo(); setState_TR(temperature(), density()); } -//==================================================================================================================== + void MixtureFugacityTP::setPressure(doublereal p) { setState_TP(temperature(), p); @@ -514,49 +335,48 @@ void MixtureFugacityTP::setPressure(doublereal p) // printf(" MixFug:setPres: mu(%d = %g) = %18.8g\n", i, mf[i], chemPot[i]); // } } -//==================================================================================================================== + void MixtureFugacityTP::setMassFractions(const doublereal* const y) { Phase::setMassFractions(y); getMoleFractions(DATA_PTR(moleFractions_)); } -//==================================================================================================================== + void MixtureFugacityTP::setMassFractions_NoNorm(const doublereal* const y) { Phase::setMassFractions_NoNorm(y); getMoleFractions(DATA_PTR(moleFractions_)); } -//==================================================================================================================== + void MixtureFugacityTP::setMoleFractions(const doublereal* const x) { Phase::setMoleFractions(x); getMoleFractions(DATA_PTR(moleFractions_)); } -//==================================================================================================================== + void MixtureFugacityTP::setMoleFractions_NoNorm(const doublereal* const x) { Phase::setMoleFractions_NoNorm(x); getMoleFractions(DATA_PTR(moleFractions_)); } -//==================================================================================================================== + void MixtureFugacityTP::setConcentrations(const doublereal* const c) { Phase::setConcentrations(c); getMoleFractions(DATA_PTR(moleFractions_)); } -//==================================================================================================================== + void MixtureFugacityTP::setMoleFractions_NoState(const doublereal* const x) { Phase::setMoleFractions(x); getMoleFractions(DATA_PTR(moleFractions_)); updateMixingExpressions(); } -//==================================================================================================================== + void MixtureFugacityTP::calcDensity() { err("MixtureFugacityTP::calcDensity() called, but EOS for phase is not known"); } -//==================================================================================================================== void MixtureFugacityTP::setState_TP(doublereal t, doublereal pres) { @@ -647,15 +467,7 @@ void MixtureFugacityTP::setState_TP(doublereal t, doublereal pres) //setPressure(pres); //calcDensity(); } -//==================================================================================================================== -// Set the internally stored temperature (K) and density (kg/m^3) -/* - * This overrides the default behavior. In addition to just storing the state in the object, we need to do - * an equation of state calculation and figure out what phase state we are in. - * - * @param t Temperature in kelvin - * @param rho Density (kg/m^3) - */ + void MixtureFugacityTP::setState_TR(doublereal T, doublereal rho) { getMoleFractions(DATA_PTR(moleFractions_)); @@ -672,40 +484,12 @@ void MixtureFugacityTP::setState_TR(doublereal T, doublereal rho) // printf("setState_TR: state at T = %g, rho = %g, mv = %g, P = %20.13g, iState = %d\n", T, rho, mv, m_Pcurrent, iState_); } -//==================================================================================================================== -// Set the temperature (K), pressure (Pa), and mole fractions. -/* - * Note, the mole fractions are set first before the pressure is set. - * Setting the pressure may involve the solution of a nonlinear equation. - * - * @param t Temperature (K) - * @param p Pressure (Pa) - * @param x Vector of mole fractions. - * Length is equal to m_kk. - */ void MixtureFugacityTP::setState_TPX(doublereal t, doublereal p, const doublereal* x) { setMoleFractions_NoState(x); setState_TP(t,p); } -//==================================================================================================================== -/* - * Import and initialize a ThermoPhase object - * - * param phaseNode This object must be the phase node of a - * complete XML tree - * description of the phase, including all of the - * species data. In other words while "phase" must - * point to an XML phase object, it must have - * sibling nodes "speciesData" that describe - * the species in the phase. - * param id ID of the phase. If nonnull, a check is done - * to see if phaseNode is pointing to the phase - * with the correct id. - * - * This routine initializes the lengths in the current object and - * then calls the parent routine. - */ + void MixtureFugacityTP::initThermoXML(XML_Node& phaseNode, const std::string& id) { MixtureFugacityTP::initLengths(); @@ -715,7 +499,7 @@ void MixtureFugacityTP::initThermoXML(XML_Node& phaseNode, const std::string& id // m_VPSS_ptr->initThermoXML(phaseNode, id); ThermoPhase::initThermoXML(phaseNode, id); } -//==================================================================================================================== + doublereal MixtureFugacityTP::z() const { doublereal p = pressure(); @@ -726,19 +510,19 @@ doublereal MixtureFugacityTP::z() const doublereal zz = p * molarV / rt; return zz; } -//==================================================================================================================== + doublereal MixtureFugacityTP::sresid() const { throw CanteraError("MixtureFugacityTP::sresid()", "Base Class: not implemented"); return 0.0; } -//==================================================================================================================== + doublereal MixtureFugacityTP::hresid() const { throw CanteraError("MixtureFugacityTP::hresid()", "Base Class: not implemented"); return 0.0; } -//==================================================================================================================== + doublereal MixtureFugacityTP::psatEst(doublereal TKelvin) const { doublereal tcrit = critTemperature(); @@ -750,32 +534,13 @@ doublereal MixtureFugacityTP::psatEst(doublereal TKelvin) const doublereal lpr = -0.8734*tt*tt - 3.4522*tt + 4.2918; return pcrit*exp(lpr); } -//==================================================================================================================== + doublereal MixtureFugacityTP::liquidVolEst(doublereal TKelvin, doublereal& pres) const { throw CanteraError("MixtureFugacityTP::liquidVolEst()", "unimplemented"); return 0.0; } -//==================================================================================================================== -/* - * Calculates the density given the temperature and the pressure, - * and a guess at the density. Note, below T_c, this is a - * multivalued function. This function assumes that the phase is on one side of the vapor dome - * or the other. It does not allow for crosses of the vapor dome. - * - * parameters: - * temperature: Kelvin - * pressure : Pressure in Pascals (Newton/m**2) - * phase : guessed phase of water - * : -1: no guessed phase - * rhoguess : guessed density of the water - * - * -1.0 no guessed density - * - * If a problem is encountered, a negative 1 is returned. - * - * @TODO make this a const function - */ + doublereal MixtureFugacityTP::densityCalc(doublereal TKelvin, doublereal presPa, int phase, doublereal rhoguess) { @@ -956,18 +721,17 @@ doublereal MixtureFugacityTP::densityCalc(doublereal TKelvin, doublereal presPa, } return densBase; } -//==================================================================================================================== + void MixtureFugacityTP::updateMixingExpressions() { - } -//==================================================================================================================== + MixtureFugacityTP::spinodalFunc::spinodalFunc(MixtureFugacityTP* tp) : ResidEval(), m_tp(tp) { } -//==================================================================================================================== + int MixtureFugacityTP::spinodalFunc::evalSS(const doublereal t, const doublereal* const y, doublereal* const r) { @@ -979,7 +743,6 @@ int MixtureFugacityTP::spinodalFunc::evalSS(const doublereal t, const doublereal r[0] = val; return status; } -//==================================================================================================================== int MixtureFugacityTP::corr0(doublereal TKelvin, doublereal pres, doublereal& densLiqGuess, doublereal& densGasGuess, doublereal& liqGRT, doublereal& gasGRT) @@ -1017,17 +780,7 @@ int MixtureFugacityTP::corr0(doublereal TKelvin, doublereal pres, doublereal& de // delGRT = gibbsLiqRT - gibbsGasRT; return retn; } -//==================================================================================================================== -// Returns the Phase State flag for the current state of the object -/* - * @param checkState If true, this function does a complete check to see where - * in parameter space we are - * - * There are three values: - * WATER_GAS below the critical temperature but below the critical density - * WATER_LIQUID below the critical temperature but above the critical density - * WATER_SUPERCRIT above the critical temperature - */ + int MixtureFugacityTP::phaseState(bool checkState) const { int state = iState_; @@ -1070,62 +823,40 @@ int MixtureFugacityTP::phaseState(bool checkState) const } return state; } -//==================================================================================================================== -// Return the value of the density at the liquid spinodal point (on the liquid side) -// for the current temperature. -/* - * @return returns the density with units of kg m-3 - */ + doublereal MixtureFugacityTP::densSpinodalLiquid() const { throw CanteraError("", "unimplmented"); return 0.0; } -//==================================================================================================================== -// Return the value of the density at the gas spinodal point (on the gas side) -// for the current temperature. -/* - * @return returns the density with units of kg m-3 - */ + doublereal MixtureFugacityTP::densSpinodalGas() const { throw CanteraError("", "unimplmented"); return 0.0; } -//==================================================================================================================== -// Calculate the saturation pressure at the current mixture content for the given temperature -/* - * This is a non-const routine that is public. - * - * The algorithm for this routine has undergone quite a bit of work. It probably needs more work. - * However, it seems now to be fairly robust. - * The key requirement is to find an initial pressure where both the liquid and the gas exist. This - * is not as easy as it sounds, and it gets exceedingly hard as the critical temperature is approached - * from below. - * Once we have this initial state, then we seek to equilibrate the gibbs free energies of the - * gas and liquid and use the formula - * - * dp = VdG - * - * to create an update condition for deltaP using - * - * - (Gliq - Ggas) = (Vliq - Vgas) (deltaP) - * - * - * - * @param TKelvin (input) Temperature (Kelvin) - * @param molarVolGas (return) Molar volume of the gas - * @param molarVolLiquid (return) Molar volume of the liquid - * - * @return Returns the saturation pressure at the given temperature - * - * @TODO Suggestions for the future would be to switch it to an algorithm that uses the gas molar volume - * and the liquid molar volumes as the fundamental unknowns. - * - */ doublereal MixtureFugacityTP::calculatePsat(doublereal TKelvin, doublereal& molarVolGas, doublereal& molarVolLiquid) { + /* + * The algorithm for this routine has undergone quite a bit of work. It probably needs more work. + * However, it seems now to be fairly robust. + * The key requirement is to find an initial pressure where both the liquid and the gas exist. This + * is not as easy as it sounds, and it gets exceedingly hard as the critical temperature is approached + * from below. + * Once we have this initial state, then we seek to equilibrate the gibbs free energies of the + * gas and liquid and use the formula + * + * dp = VdG + * + * to create an update condition for deltaP using + * + * - (Gliq - Ggas) = (Vliq - Vgas) (deltaP) + * + * @TODO Suggestions for the future would be to switch it to an algorithm that uses the gas molar volume + * and the liquid molar volumes as the fundamental unknowns. + */ + // we need this because this is a non-const routine that is public setTemperature(TKelvin); double tcrit = critTemperature(); @@ -1316,36 +1047,18 @@ doublereal MixtureFugacityTP::calculatePsat(doublereal TKelvin, doublereal& mola return pres; } -//==================================================================================================================== -// Calculate the pressure given the temperature and the molar volume doublereal MixtureFugacityTP::pressureCalc(doublereal TKelvin, doublereal molarVol) const { throw CanteraError("MixtureFugacityTP::pressureCalc", "unimplemented"); return 0.0; } -//==================================================================================================================== -// Calculate the pressure given the temperature and the molar volume + doublereal MixtureFugacityTP::dpdVCalc(doublereal TKelvin, doublereal molarVol, doublereal& presCalc) const { throw CanteraError("MixtureFugacityTP::dpdVCalc", "unimplemented"); return 0.0; } -//==================================================================================================================== -/* - * void _updateStandardStateThermo() (protected, virtual, const) - * - * If m_useTmpStandardStateStorage is true, - * This function must be called for every call to functions in this - * class that need standard state properties. - * Child classes may require that it be called even if m_useTmpStandardStateStorage - * is not true. - * It checks to see whether the temperature has changed and - * thus the ss thermodynamics functions for all of the species - * must be recalculated. - * - * This - */ void MixtureFugacityTP::_updateReferenceStateThermo() const { double Tnow = temperature(); @@ -1367,9 +1080,5 @@ void MixtureFugacityTP::_updateReferenceStateThermo() const m_logc0 = log(pref/(GasConstant * Tnow)); } } -//==================================================================================================================== - } - - diff --git a/src/thermo/PureFluidPhase.cpp b/src/thermo/PureFluidPhase.cpp index 3b1ca464d..50bd57246 100644 --- a/src/thermo/PureFluidPhase.cpp +++ b/src/thermo/PureFluidPhase.cpp @@ -22,7 +22,6 @@ using std::setw; namespace Cantera { -// Base Constructor PureFluidPhase::PureFluidPhase() : ThermoPhase(), m_sub(0), @@ -32,7 +31,6 @@ PureFluidPhase::PureFluidPhase() : { } -// CopyConstructor PureFluidPhase::PureFluidPhase(const PureFluidPhase& right) : ThermoPhase(), m_sub(0), @@ -43,10 +41,6 @@ PureFluidPhase::PureFluidPhase(const PureFluidPhase& right) : *this = right; } -//! Assignment operator -/*! - * @param right Object to be copied - */ PureFluidPhase& PureFluidPhase::operator=(const PureFluidPhase& right) { if (&right != this) { @@ -62,22 +56,11 @@ PureFluidPhase& PureFluidPhase::operator=(const PureFluidPhase& right) return *this; } -// Duplicator from the %ThermoPhase parent class -/* - * Given a pointer to a %ThermoPhase object, this function will - * duplicate the %ThermoPhase object and all underlying structures. - * This is basically a wrapper around the copy constructor. - * - * @return returns a pointer to a %ThermoPhase - */ ThermoPhase* PureFluidPhase::duplMyselfAsThermoPhase() const { return new PureFluidPhase(*this); } - - - PureFluidPhase::~PureFluidPhase() { delete m_sub; @@ -175,178 +158,108 @@ pressure() const setTPXState(); return m_sub->P(); } -//==================================================================================================================== -void PureFluidPhase:: -setPressure(doublereal p) + +void PureFluidPhase::setPressure(doublereal p) { Set(tpx::PropertyPair::TP, temperature(), p); setDensity(1.0/m_sub->v()); } -//==================================================================================================================== + void PureFluidPhase::Set(tpx::PropertyPair::type n, double x, double y) const { m_sub->Set(n, x, y); } -//==================================================================================================================== + void PureFluidPhase::setTPXState() const { Set(tpx::PropertyPair::TV, temperature(), 1.0/density()); } -//==================================================================================================================== doublereal PureFluidPhase::isothermalCompressibility() const { return m_sub->isothermalCompressibility(); } -//==================================================================================================================== + doublereal PureFluidPhase::thermalExpansionCoeff() const { return m_sub->thermalExpansionCoeff(); } -//==================================================================================================================== + tpx::Substance& PureFluidPhase::TPX_Substance() { return *m_sub; } -//==================================================================================================================== -// Returns an array of partial molar enthalpies for the species -// in the mixture. Units (J/kmol) -/* - * @param hbar Output vector of species partial molar enthalpies. - * Length: m_kk. units are J/kmol. - */ + void PureFluidPhase::getPartialMolarEnthalpies(doublereal* hbar) const { hbar[0] = enthalpy_mole(); } -//==================================================================================================================== -// Returns an array of partial molar entropies of the species in the -// solution. Units: J/kmol/K. -/* - * @param sbar Output vector of species partial molar entropies. - * Length = m_kk. units are J/kmol/K. - */ + void PureFluidPhase::getPartialMolarEntropies(doublereal* sbar) const { sbar[0] = entropy_mole(); } -//==================================================================================================================== -// Return an array of partial molar internal energies for the -// species in the mixture. Units: J/kmol. -/* - * @param ubar Output vector of species partial molar internal energies. - * Length = m_kk. units are J/kmol. - */ + void PureFluidPhase::getPartialMolarIntEnergies(doublereal* ubar) const { ubar[0] = intEnergy_mole(); } -//==================================================================================================================== -// Return an array of partial molar heat capacities for the -// species in the mixture. Units: J/kmol/K -/* - * @param cpbar Output vector of species partial molar heat - * capacities at constant pressure. - * Length = m_kk. units are J/kmol/K. - */ + void PureFluidPhase::getPartialMolarCp(doublereal* cpbar) const { cpbar[0] = cp_mole(); } -//==================================================================================================================== -// Return an array of partial molar volumes for the -// species in the mixture. Units: m^3/kmol. -/* - * @param vbar Output vector of species partial molar volumes. - * Length = m_kk. units are m^3/kmol. - */ + void PureFluidPhase::getPartialMolarVolumes(doublereal* vbar) const { vbar[0] = 1.0 / molarDensity(); } -//==================================================================================================================== + int PureFluidPhase::standardStateConvention() const { return cSS_CONVENTION_TEMPERATURE; } -//==================================================================================================================== + void PureFluidPhase::getActivityConcentrations(doublereal* c) const { c[0] = 1.0; } -//==================================================================================================================== + doublereal PureFluidPhase::standardConcentration(size_t k) const { return 1.0; } -//==================================================================================================================== + void PureFluidPhase::getActivities(doublereal* a) const { a[0] = 1.0; } -//==================================================================================================================== -// 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. - */ + void PureFluidPhase::getStandardChemPotentials(doublereal* mu) const { mu[0] = gibbs_mole(); } -//==================================================================================================================== -// 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. - */ + void PureFluidPhase::getEnthalpy_RT(doublereal* hrt) const { doublereal rt = _RT(); doublereal h = enthalpy_mole(); hrt[0] = h / 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. - */ + void PureFluidPhase::getEntropy_R(doublereal* sr) const { doublereal s = entropy_mole(); sr[0] = s / GasConstant; } -//==================================================================================================================== -// 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. - */ + void PureFluidPhase::getGibbs_RT(doublereal* grt) const { doublereal rt = _RT(); doublereal g = gibbs_mole(); grt[0] = g / rt; } -//==================================================================================================================== -// Returns the vector of nondimensional enthalpies of the reference state at the current temperature -// of the solution and the reference pressure for the species. -/* - * This base function will throw a CanteraException unless - * it is overwritten in a derived class. - * - * @param hrt Output vector containing the nondimensional reference state enthalpies - * Length: m_kk. - */ + void PureFluidPhase::getEnthalpy_RT_ref(doublereal* hrt) const { double psave = pressure(); @@ -358,13 +271,7 @@ void PureFluidPhase::getEnthalpy_RT_ref(doublereal* hrt) const Set(tpx::PropertyPair::TP, t, psave); } -//==================================================================================================================== -// Returns the vector of nondimensional Gibbs Free Energies of the reference state at the current temperature -// of the solution and the reference pressure for the species. -/* - * @param grt Output vector containing the nondimensional reference state - * Gibbs Free energies. Length: m_kk. - */ + void PureFluidPhase::getGibbs_RT_ref(doublereal* grt) const { double psave = pressure(); @@ -376,27 +283,13 @@ void PureFluidPhase::getGibbs_RT_ref(doublereal* grt) const grt[0] += log(pref/plow); Set(tpx::PropertyPair::TP, t, psave); } -//==================================================================================================================== -// Returns the vector of the gibbs function of the reference state at the current temperature -// of the solution and the reference pressure for the species. -/* - * units = J/kmol - * - * @param g Output vector containing the reference state - * Gibbs Free energies. Length: m_kk. Units: J/kmol. - */ + void PureFluidPhase::getGibbs_ref(doublereal* g) const { getGibbs_RT_ref(g); g[0] *= (GasConstant * temperature()); } -//==================================================================================================================== -// 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. - */ + void PureFluidPhase::getEntropy_R_ref(doublereal* er) const { double psave = pressure(); @@ -408,62 +301,56 @@ void PureFluidPhase::getEntropy_R_ref(doublereal* er) const er[0] -= log(pref/plow); Set(tpx::PropertyPair::TP, t, psave); } -//==================================================================================================================== -// critical temperature + doublereal PureFluidPhase::critTemperature() const { return m_sub->Tcrit(); } -//==================================================================================================================== -/// critical pressure + doublereal PureFluidPhase::critPressure() const { return m_sub->Pcrit(); } -//==================================================================================================================== -/// critical density + doublereal PureFluidPhase::critDensity() const { return 1.0/m_sub->Vcrit(); } -//==================================================================================================================== -/// saturation temperature doublereal PureFluidPhase::satTemperature(doublereal p) const { doublereal ts = m_sub->Tsat(p); return ts; } -//==================================================================================================================== + void PureFluidPhase::setState_HP(doublereal h, doublereal p, doublereal tol) { Set(tpx::PropertyPair::HP, h, p); setState_TR(m_sub->Temp(), 1.0/m_sub->v()); } -//==================================================================================================================== + void PureFluidPhase::setState_UV(doublereal u, doublereal v, doublereal tol) { Set(tpx::PropertyPair::UV, u, v); setState_TR(m_sub->Temp(), 1.0/m_sub->v()); } -//==================================================================================================================== + void PureFluidPhase::setState_SV(doublereal s, doublereal v, doublereal tol) { Set(tpx::PropertyPair::SV, s, v); setState_TR(m_sub->Temp(), 1.0/m_sub->v()); } -//==================================================================================================================== + void PureFluidPhase::setState_SP(doublereal s, doublereal p, doublereal tol) { Set(tpx::PropertyPair::SP, s, p); setState_TR(m_sub->Temp(), 1.0/m_sub->v()); } -//==================================================================================================================== -// saturation pressure + doublereal PureFluidPhase::satPressure(doublereal t) const { doublereal vsv = m_sub->v(); @@ -471,13 +358,13 @@ doublereal PureFluidPhase::satPressure(doublereal t) const doublereal ps = m_sub->Ps(); return ps; } -//==================================================================================================================== + doublereal PureFluidPhase::vaporFraction() const { setTPXState(); return m_sub->x(); } -//==================================================================================================================== + void PureFluidPhase::setState_Tsat(doublereal t, doublereal x) { setTemperature(t); @@ -485,7 +372,7 @@ void PureFluidPhase::setState_Tsat(doublereal t, doublereal x) Set(tpx::PropertyPair::TX, t, x); setDensity(1.0/m_sub->v()); } -//==================================================================================================================== + void PureFluidPhase::setState_Psat(doublereal p, doublereal x) { setTPXState(); @@ -494,10 +381,6 @@ void PureFluidPhase::setState_Psat(doublereal p, doublereal x) setDensity(1.0/m_sub->v()); } -//==================================================================================================================== -/** - * Format a summary of the mixture state for output. - */ std::string PureFluidPhase::report(bool show_thermo) const { char p[800]; diff --git a/src/thermo/RedlichKwongMFTP.cpp b/src/thermo/RedlichKwongMFTP.cpp index d1f6c2df9..2f7e436b5 100644 --- a/src/thermo/RedlichKwongMFTP.cpp +++ b/src/thermo/RedlichKwongMFTP.cpp @@ -28,10 +28,6 @@ const doublereal RedlichKwongMFTP::omega_a = 4.27480233540E-01; const doublereal RedlichKwongMFTP::omega_b = 8.66403499650E-02; const doublereal RedlichKwongMFTP::omega_vc = 3.33333333333333E-01; -//==================================================================================================================== -/* - * Default constructor - */ RedlichKwongMFTP::RedlichKwongMFTP() : MixtureFugacityTP(), m_standardMixingRules(0), @@ -56,7 +52,7 @@ RedlichKwongMFTP::RedlichKwongMFTP() : Vroot_[1] = 0.0; Vroot_[2] = 0.0; } -//==================================================================================================================== + RedlichKwongMFTP::RedlichKwongMFTP(const std::string& infile, std::string id) : MixtureFugacityTP(), m_standardMixingRules(0), @@ -91,7 +87,7 @@ RedlichKwongMFTP::RedlichKwongMFTP(const std::string& infile, std::string id) : } importPhase(*xphase, this); } -//==================================================================================================================== + RedlichKwongMFTP::RedlichKwongMFTP(XML_Node& phaseRefRoot, const std::string& id) : MixtureFugacityTP(), m_standardMixingRules(0), @@ -122,7 +118,6 @@ RedlichKwongMFTP::RedlichKwongMFTP(XML_Node& phaseRefRoot, const std::string& id importPhase(*xphase, this); } -//==================================================================================================================== RedlichKwongMFTP::RedlichKwongMFTP(int testProb) : MixtureFugacityTP(), m_standardMixingRules(0), @@ -161,16 +156,7 @@ RedlichKwongMFTP::RedlichKwongMFTP(int testProb) : } importPhase(*xphase, this); } -//==================================================================================================================== -/* - * Copy Constructor: - * - * Note this stuff will not work until the underlying phase - * has a working copy constructor. - * - * The copy constructor just calls the assignment operator - * to do the heavy lifting. - */ + RedlichKwongMFTP::RedlichKwongMFTP(const RedlichKwongMFTP& b) : MixtureFugacityTP(), m_standardMixingRules(0), @@ -194,13 +180,6 @@ RedlichKwongMFTP::RedlichKwongMFTP(const RedlichKwongMFTP& b) : *this = b; } -//==================================================================================================================== -/* - * operator=() - * - * Note this stuff will not work until the underlying phase - * has a working assignment operator - */ RedlichKwongMFTP& RedlichKwongMFTP:: operator=(const RedlichKwongMFTP& b) { @@ -237,35 +216,25 @@ operator=(const RedlichKwongMFTP& b) } return *this; } -//==================================================================================================================== -/* - * ~RedlichKwongMFTP(): (virtual) - * - */ + RedlichKwongMFTP::~RedlichKwongMFTP() { } -//==================================================================================================================== -/* - * Duplication function. - * This calls the copy constructor for this object. - */ + ThermoPhase* RedlichKwongMFTP::duplMyselfAsThermoPhase() const { return new RedlichKwongMFTP(*this); } -//==================================================================================================================== + int RedlichKwongMFTP::eosType() const { return cRedlichKwongMFTP; } -//==================================================================================================================== /* * ------------Molar Thermodynamic Properties ------------------------- */ -//==================================================================================================================== -// Molar enthalpy. Units: J/kmol. + doublereal RedlichKwongMFTP::enthalpy_mole() const { _updateReferenceStateThermo(); @@ -274,16 +243,14 @@ doublereal RedlichKwongMFTP::enthalpy_mole() const doublereal h_nonideal = hresid(); return (h_ideal + h_nonideal); } -//==================================================================================================================== -// Molar internal energy. Units: J/kmol. + doublereal RedlichKwongMFTP::intEnergy_mole() const { doublereal p0 = pressure(); doublereal md = molarDensity(); return (enthalpy_mole() - p0 / md); } -//==================================================================================================================== -// Molar entropy. Units: J/kmol/K. + doublereal RedlichKwongMFTP::entropy_mole() const { _updateReferenceStateThermo(); @@ -292,14 +259,12 @@ doublereal RedlichKwongMFTP::entropy_mole() const doublereal sr_nonideal = sresid(); return (sr_ideal + sr_nonideal); } -//==================================================================================================================== -// Molar Gibbs function. Units: J/kmol. + doublereal RedlichKwongMFTP::gibbs_mole() const { return enthalpy_mole() - temperature() * entropy_mole(); } -//==================================================================================================================== -/// Molar heat capacity at constant pressure. Units: J/kmol/K. + doublereal RedlichKwongMFTP::cp_mole() const { _updateReferenceStateThermo(); @@ -316,29 +281,15 @@ doublereal RedlichKwongMFTP::cp_mole() const double cp = dHdT_V - (mv + TKelvin * dpdT_ / dpdV_) * dpdT_; return cp; } -//==================================================================================================================== -/// Molar heat capacity at constant volume. Units: J/kmol/K. + doublereal RedlichKwongMFTP::cv_mole() const { throw CanteraError("", "unimplemented"); return cp_mole() - GasConstant; } -//==================================================================================================================== -// Return the thermodynamic pressure (Pa). -/* - * Since the mass density, temperature, and mass fractions are stored, - * this method uses these values to implement the - * mechanical equation of state \f$ P(T, \rho, Y_1, \dots, Y_K) \f$. - * - * \f[ - * P = \frac{RT}{v-b_{mix}} - \frac{a_{mix}}{T^{0.5} v \left( v + b_{mix} \right) } - * \f] - * - */ + doublereal RedlichKwongMFTP::pressure() const { - - #ifdef DEBUG_MODE _updateReferenceStateThermo(); @@ -354,10 +305,9 @@ doublereal RedlichKwongMFTP::pressure() const throw CanteraError(" RedlichKwongMFTP::pressure()", "setState broken down, maybe"); } #endif - return m_Pcurrent; } -//==================================================================================================================== + void RedlichKwongMFTP::calcDensity() { /* @@ -373,58 +323,52 @@ void RedlichKwongMFTP::calcDensity() */ double dens = 1.0/invDens; Phase::setDensity(dens); - } -//==================================================================================================================== void RedlichKwongMFTP::setTemperature(const doublereal temp) { Phase::setTemperature(temp); _updateReferenceStateThermo(); updateAB(); } -//==================================================================================================================== + void RedlichKwongMFTP::setMassFractions(const doublereal* const x) { MixtureFugacityTP::setMassFractions(x); updateAB(); } -//==================================================================================================================== + void RedlichKwongMFTP::setMassFractions_NoNorm(const doublereal* const x) { MixtureFugacityTP::setMassFractions_NoNorm(x); updateAB(); } -//==================================================================================================================== + void RedlichKwongMFTP::setMoleFractions(const doublereal* const x) { MixtureFugacityTP::setMoleFractions(x); updateAB(); } -//==================================================================================================================== + void RedlichKwongMFTP::setMoleFractions_NoNorm(const doublereal* const x) { MixtureFugacityTP::setMoleFractions(x); updateAB(); } -//==================================================================================================================== + void RedlichKwongMFTP::setConcentrations(const doublereal* const c) { MixtureFugacityTP::setConcentrations(c); updateAB(); } -//==================================================================================================================== doublereal RedlichKwongMFTP::isothermalCompressibility() const { - - throw CanteraError("RedlichKwongMFTP::isothermalCompressibility() ", "not implemented"); - return 0.0; } -//==================================================================================================================== + void RedlichKwongMFTP::getActivityConcentrations(doublereal* c) const { getPartialMolarVolumes(DATA_PTR(m_partialMolarVolumes)); @@ -432,59 +376,20 @@ void RedlichKwongMFTP::getActivityConcentrations(doublereal* c) const c[k] = moleFraction(k) / m_partialMolarVolumes[k]; } } -//==================================================================================================================== -/* - * Returns the standard concentration \f$ C^0_k \f$, which is used to normalize - * the generalized concentration. - */ + doublereal RedlichKwongMFTP::standardConcentration(size_t k) const { - getStandardVolumes(DATA_PTR(m_tmpV)); - return 1.0 / m_tmpV[k]; - - - } -//==================================================================================================================== -/* - * Returns the natural logarithm of the standard - * concentration of the kth species - */ + doublereal RedlichKwongMFTP::logStandardConc(size_t k) const { double c = standardConcentration(k); double lc = std::log(c); return lc; } -//==================================================================================================================== -/* - * - * getUnitsStandardConcentration() - * - * Returns the units of the standard and general concentrations - * Note they have the same units, as their divisor 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 - * - * For EOS types other than cIdealSolidSolnPhase1, the default - * kmol/m3 holds for standard concentration units. For - * cIdealSolidSolnPhase0 type, the standard concentration is - * unitless. - */ + void RedlichKwongMFTP::getUnitsStandardConc(double* uA, int, int sizeUA) const { //int eos = eosType(); @@ -512,14 +417,6 @@ void RedlichKwongMFTP::getUnitsStandardConc(double* uA, int, int sizeUA) const } -//==================================================================================================================== -//! Get the array of non-dimensional activity coefficients at -//! the current solution temperature, pressure, and solution concentration. -/*! - * For ideal gases, the activity coefficients are all equal to one. - * - * @param ac Output vector of activity coefficients. Length: m_kk. - */ void RedlichKwongMFTP::getActivityCoefficients(doublereal* ac) const { doublereal TKelvin = temperature(); @@ -551,20 +448,11 @@ void RedlichKwongMFTP::getActivityCoefficients(doublereal* ac) const ac[k] = exp(ac[k]/rt); } } -//==================================================================================================================== + /* * ---- Partial Molar Properties of the Solution ----------------- */ -//==================================================================================================================== -/* - * Get the array of non-dimensional species chemical potentials - * These are partial molar Gibbs free energies. - * \f$ \mu_k / \hat R T \f$. - * Units: unitless - * - * We close the loop on this function, here, calling - * getChemPotentials() and then dividing by RT. - */ + void RedlichKwongMFTP::getChemPotentials_RT(doublereal* muRT) const { getChemPotentials(muRT); @@ -573,7 +461,7 @@ void RedlichKwongMFTP::getChemPotentials_RT(doublereal* muRT) const muRT[k] *= invRT; } } -//==================================================================================================================== + void RedlichKwongMFTP::getChemPotentials(doublereal* mu) const { getGibbs_ref(mu); @@ -610,7 +498,7 @@ void RedlichKwongMFTP::getChemPotentials(doublereal* mu) const ); } } -//==================================================================================================================== + void RedlichKwongMFTP::getPartialMolarEnthalpies(doublereal* hbar) const { /* @@ -669,7 +557,7 @@ void RedlichKwongMFTP::getPartialMolarEnthalpies(doublereal* hbar) const } } -//==================================================================================================================== + void RedlichKwongMFTP::getPartialMolarEntropies(doublereal* sbar) const { getEntropy_R_ref(sbar); @@ -726,21 +614,21 @@ void RedlichKwongMFTP::getPartialMolarEntropies(doublereal* sbar) const sbar[k] -= -m_partialMolarVolumes[k] * dpdT_; } } -//==================================================================================================================== + void RedlichKwongMFTP::getPartialMolarIntEnergies(doublereal* ubar) const { getIntEnergy_RT(ubar); doublereal rt = GasConstant * temperature(); scale(ubar, ubar+m_kk, ubar, rt); } -//==================================================================================================================== + void RedlichKwongMFTP::getPartialMolarCp(doublereal* cpbar) const { getCp_R(cpbar); doublereal r = GasConstant; scale(cpbar, cpbar+m_kk, cpbar, r); } -//==================================================================================================================== + void RedlichKwongMFTP::getPartialMolarVolumes(doublereal* vbar) const { // getStandardVolumes(vbar); @@ -786,7 +674,7 @@ void RedlichKwongMFTP::getPartialMolarVolumes(doublereal* vbar) const } } -//==================================================================================================================== + doublereal RedlichKwongMFTP::critTemperature() const { double pc, tc, vc; @@ -802,7 +690,7 @@ doublereal RedlichKwongMFTP::critTemperature() const calcCriticalConditions(m_a_current, m_b_current, a0, aT, pc, tc, vc); return tc; } -//==================================================================================================================== + doublereal RedlichKwongMFTP::critPressure() const { double pc, tc, vc; @@ -819,7 +707,7 @@ doublereal RedlichKwongMFTP::critPressure() const return pc; } -//==================================================================================================================== + doublereal RedlichKwongMFTP::critDensity() const { double pc, tc, vc; @@ -837,26 +725,13 @@ doublereal RedlichKwongMFTP::critDensity() const double mmw = meanMolecularWeight(); return mmw / vc; } -//==================================================================================================================== -/* - * ----- Thermodynamic Values for the Species Reference States ---- - */ - - -//==================================================================================================================== - -/* - * Perform initializations after all species have been - * added. - */ void RedlichKwongMFTP::initThermo() { initLengths(); MixtureFugacityTP::initThermo(); } -//==================================================================================================================== void RedlichKwongMFTP::setToEquilState(const doublereal* mu_RT) { double tmp, tmp2; @@ -891,11 +766,7 @@ void RedlichKwongMFTP::setToEquilState(const doublereal* mu_RT) // set state setState_PX(pres, &m_pp[0]); } -//==================================================================================================================== -/* - * Initialize the internal lengths. - * (this is not a virtual function) - */ + void RedlichKwongMFTP::initLengths() { @@ -916,24 +787,7 @@ void RedlichKwongMFTP::initLengths() m_partialMolarVolumes.resize(m_kk, 0.0); dpdni_.resize(m_kk, 0.0); } -//==================================================================================================================== -/* - * Import and initialize a ThermoPhase object - * - * param phaseNode This object must be the phase node of a - * complete XML tree - * description of the phase, including all of the - * species data. In other words while "phase" must - * point to an XML phase object, it must have - * sibling nodes "speciesData" that describe - * the species in the phase. - * param id ID of the phase. If nonnull, a check is done - * to see if phaseNode is pointing to the phase - * with the correct id. - * - * This routine initializes the lengths in the current object and - * then calls the parent routine. - */ + void RedlichKwongMFTP::initThermoXML(XML_Node& phaseNode, const std::string& id) { RedlichKwongMFTP::initLengths(); @@ -1018,7 +872,6 @@ void RedlichKwongMFTP::initThermoXML(XML_Node& phaseNode, const std::string& id) MixtureFugacityTP::initThermoXML(phaseNode, id); } -//==================================================================================================================== void RedlichKwongMFTP::readXMLPureFluid(XML_Node& pureFluidParam) { @@ -1083,7 +936,7 @@ void RedlichKwongMFTP::readXMLPureFluid(XML_Node& pureFluidParam) } } } -//==================================================================================================================== + void RedlichKwongMFTP::applyStandardMixingRules() { int nParam = 2; @@ -1100,7 +953,6 @@ void RedlichKwongMFTP::applyStandardMixingRules() } } } -//==================================================================================================================== void RedlichKwongMFTP::readXMLCrossFluid(XML_Node& CrossFluidParam) { @@ -1168,22 +1020,13 @@ void RedlichKwongMFTP::readXMLCrossFluid(XML_Node& CrossFluidParam) } } } -//==================================================================================================================== + void RedlichKwongMFTP::setParametersFromXML(const XML_Node& thermoNode) { MixtureFugacityTP::setParametersFromXML(thermoNode); std::string model = thermoNode["model"]; - - } -//==================================================================================================================== -// Calculate the deviation terms for the total entropy of the mixture from the -// ideal gas mixture -/* - * Here we use the current state conditions - * - * @return Returns the change in entropy in units of J kmol-1 K-1. - */ + doublereal RedlichKwongMFTP::sresid() const { // note this agrees with tpx @@ -1200,14 +1043,7 @@ doublereal RedlichKwongMFTP::sresid() const double sp = GasConstant * sresid_mol_R; return sp; } -//==================================================================================================================== -// Calculate the deviation terms for the total enthalpy of the mixture from the -// ideal gas mixture -/* - * Here we use the current state conditions - * - * @return Returns the change in entropy in units of J kmol-1. - */ + doublereal RedlichKwongMFTP::hresid() const { // note this agrees with tpx @@ -1223,21 +1059,7 @@ doublereal RedlichKwongMFTP::hresid() const double hresid_mol = GasConstant * T * (zz - 1.0) + fac * log(1.0 + hh) / (sqT * m_b_current); return hresid_mol; } -//==================================================================================================================== -// Estimate for the molar volume of the liquid -/* - * Note: this is only used as a starting guess for later routines that actually calculate an - * accurate value for the liquid molar volume. - * This routine doesn't change the state of the system. - * - * @param TKelvin temperature in kelvin - * @param pres Pressure in Pa. This is used as an initial guess. If the routine - * needs to change the pressure to find a stable liquid state, the - * new pressure is returned in this variable. - * - * @return Returns the estimate of the liquid volume. If the liquid can't be found, this - * routine returns -1. - */ + doublereal RedlichKwongMFTP::liquidVolEst(doublereal TKelvin, doublereal& presGuess) const { double v = m_b_current * 1.1; @@ -1282,31 +1104,7 @@ doublereal RedlichKwongMFTP::liquidVolEst(doublereal TKelvin, doublereal& presGu //printf (" RedlichKwongMFTP::liquidVolEst %g %g converged in %d its\n", TKelvin, pres, i); return v; } -//==================================================================================================================== -// Calculates the density given the temperature and the pressure and a guess at the density. -/* - * Note, below T_c, this is a multivalued function. We do not cross the vapor dome in this. - * This is protected because it is called during setState_TP() routines. Infinite loops would result - * if it were not protected. - * - * -> why is this not const? - * - * parameters: - * @param TKelvin Temperature in Kelvin - * @param pressure Pressure in Pascals (Newton/m**2) - * @param phaseReqested int representing the phase whose density we are requesting. If we put - * a gas or liquid phase here, we will attempt to find a volume in that - * part of the volume space, only, in this routine. A value of FLUID_UNDEFINED - * means that we will accept anything. - * - * @param rhoguess Guessed density of the fluid. A value of -1.0 indicates that there - * is no guessed density - * - * - * @return We return the density of the fluid at the requested phase. If we have not found any - * acceptable density we return a -1. If we have found an acceptable density at a - * different phase, we return a -2. - */ + doublereal RedlichKwongMFTP::densityCalc(doublereal TKelvin, doublereal presPa, int phaseRequested, doublereal rhoguess) { @@ -1382,12 +1180,7 @@ doublereal RedlichKwongMFTP::densityCalc(doublereal TKelvin, doublereal presPa, densBase = mmw / molarVolLast; return densBase; } -//==================================================================================================================== -// Return the value of the density at the liquid spinodal point (on the liquid side) -// for the current temperature. -/* - * @return returns the density with units of kg m-3 - */ + doublereal RedlichKwongMFTP::densSpinodalLiquid() const { if (NSolns_ != 3) { @@ -1412,12 +1205,7 @@ doublereal RedlichKwongMFTP::densSpinodalLiquid() const doublereal rho = mmw / vbest; return rho; } -//==================================================================================================================== -// Return the value of the density at the gas spinodal point (on the gas side) -// for the current temperature. -/* - * @return returns the density with units of kg m-3 - */ + doublereal RedlichKwongMFTP::densSpinodalGas() const { if (NSolns_ != 3) { @@ -1442,16 +1230,7 @@ doublereal RedlichKwongMFTP::densSpinodalGas() const doublereal rho = mmw / vbest; return rho; } -//==================================================================================================================== -// Calculate the pressure given the temperature and the molar volume -/* - * Calculate the pressure given the temperature and the molar volume - * - * @param TKelvin temperature in kelvin - * @param molarVol molar volume ( m3/kmol) - * - * @return Returns the pressure. - */ + doublereal RedlichKwongMFTP::pressureCalc(doublereal TKelvin, doublereal molarVol) const { doublereal sqt = sqrt(TKelvin); @@ -1459,18 +1238,7 @@ doublereal RedlichKwongMFTP::pressureCalc(doublereal TKelvin, doublereal molarVo - m_a_current / (sqt * molarVol * (molarVol + m_b_current)); return pres; } -//==================================================================================================================== -// Calculate the pressure and the pressure derivative given the temperature and the molar volume -/* - * Temperature and mole number are held constant - * - * @param TKelvin temperature in kelvin - * @param molarVol molar volume ( m3/kmol) - * - * @param presCalc Returns the pressure. - * - * @return Returns the derivative of the pressure wrt the molar volume - */ + doublereal RedlichKwongMFTP::dpdVCalc(doublereal TKelvin, doublereal molarVol, doublereal& presCalc) const { doublereal sqt = sqrt(TKelvin); @@ -1483,7 +1251,6 @@ doublereal RedlichKwongMFTP::dpdVCalc(doublereal TKelvin, doublereal molarVol, + m_a_current * (2 * molarVol + m_b_current) / (sqt * molarVol * molarVol * vpb * vpb)); return dpdv; } -//==================================================================================================================== void RedlichKwongMFTP::pressureDerivatives() const { @@ -1501,12 +1268,12 @@ void RedlichKwongMFTP::pressureDerivatives() const dpdT_ = (GasConstant / (vmb) - fac / (sqt * mv * vpb)); } -//==================================================================================================================== + void RedlichKwongMFTP::updateMixingExpressions() { updateAB(); } -//==================================================================================================================== + void RedlichKwongMFTP::updateAB() { double temp = temperature(); @@ -1528,7 +1295,7 @@ void RedlichKwongMFTP::updateAB() } } } -//==================================================================================================================== + void RedlichKwongMFTP::calculateAB(doublereal temp, doublereal& aCalc, doublereal& bCalc) const { bCalc = 0.0; @@ -1553,7 +1320,7 @@ void RedlichKwongMFTP::calculateAB(doublereal temp, doublereal& aCalc, doublerea } } } -//==================================================================================================================== + doublereal RedlichKwongMFTP::da_dt() const { @@ -1568,7 +1335,7 @@ doublereal RedlichKwongMFTP::da_dt() const } return dadT; } -//==================================================================================================================== + void RedlichKwongMFTP::calcCriticalConditions(doublereal a, doublereal b, doublereal a0_coeff, doublereal aT_coeff, doublereal& pc, doublereal& tc, doublereal& vc) const { @@ -1612,18 +1379,6 @@ void RedlichKwongMFTP::calcCriticalConditions(doublereal a, doublereal b, double vc = omega_vc * GasConstant * tc / pc; } -//==================================================================================================================== -// Solve the cubic equation of state -/* - * The R-K equation of state may be solved via the following formula - * - * V**3 - V**2(RT/P) - V(RTb/P - a/(P T**.5) + b*b) - (a b / (P T**.5)) = 0 - * - - * Returns the number of solutions found. If it only finds the liquid branch solution, it will return a -1 or a -2 - * instead of 1 or 2. If it returns 0, then there is an error. - * - */ int RedlichKwongMFTP::NicholsSolve(double TKelvin, double pres, doublereal a, doublereal b, doublereal Vroot[3]) const {