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

Implement MaskellSolidSolnPhase::standardConcentrations() and getActivityConcentrations()

Browse source
This commit is contained in:
Victor Brunini 2014-03-06 22:30:19 +00:00
parent 5898bf025f
commit 7abb459bbe
3 changed files with 103 additions and 27 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

76
include/cantera/thermo/MaskellSolidSolnPhase.h
View file

@ -34,26 +34,8 @@ namespace Cantera
class MaskellSolidSolnPhase : public VPStandardStateTP
{
public:
/**
* The generalized concentrations can have three different forms
* depending on the value of the member attribute #m_formGC, which
* is supplied in the constructor or read from the xml data file.
*
* @param formCG This parameter initializes the #m_formGC variable.
*/
MaskellSolidSolnPhase();
//! Construct and initialize an MaskellSolidSolnPhase ThermoPhase object
//! directly from an XML database
/*!
* @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.
*/
//MaskellSolidSolnPhase(XML_Node& root, const std::string& id="");
//! Copy Constructor
MaskellSolidSolnPhase(const MaskellSolidSolnPhase&);
@ -67,6 +49,61 @@ public:
*/
virtual ThermoPhase* duplMyselfAsThermoPhase() const;
/**
* 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<SUP>-3</SUP>. 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.
*
* 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.
*/
virtual void getActivityConcentrations(doublereal* c) const;
//! Return the standard concentration for the kth species
/*!
* The standard concentration \f$ C^0_k \f$ used to normalize the
* 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<SUP>-3</SUP>.
*
* @param k Species number: this is a require parameter,
* a change from the ThermoPhase base class, where it was
* an optional parameter.
*/
virtual doublereal standardConcentration(size_t k) const;
//! @name Molar Thermodynamic Properties of the Solution
//! @{
/**
@ -310,6 +347,9 @@ protected:
//! Value of the enthalpy change on mixing due to protons changing from type B to type A configurations.
doublereal h_mixing;
//! Index of the species whose mole fraction defines the extent of reduction r
int product_species_index;
private:
// Functions to calculate some of the pieces of the mixing terms.
doublereal s() const;

52
src/thermo/MaskellSolidSolnPhase.cpp
View file

@ -28,7 +28,8 @@ MaskellSolidSolnPhase::MaskellSolidSolnPhase() :
m_cp0_R(2),
m_g0_RT(2),
m_s0_R(2),
h_mixing(0.0)
h_mixing(0.0),
product_species_index(0)
{
}
//=====================================================================================================
@ -40,7 +41,8 @@ MaskellSolidSolnPhase::MaskellSolidSolnPhase(const MaskellSolidSolnPhase& b) :
m_cp0_R(2),
m_g0_RT(2),
m_s0_R(2),
h_mixing(0.0)
h_mixing(0.0),
product_species_index(0)
{
*this = b;
}
@ -59,6 +61,26 @@ ThermoPhase* MaskellSolidSolnPhase::duplMyselfAsThermoPhase() const
return new MaskellSolidSolnPhase(*this);
}
//=====================================================================================================
void MaskellSolidSolnPhase::
getActivityConcentrations(doublereal* c) const
{
std::vector<doublereal> pmv(m_kk);
getPartialMolarVolumes(&pmv[0]);
const doublereal* const dtmp = moleFractdivMMW();
const double mmw = meanMolecularWeight();
for (size_t k = 0; k < m_kk; k++) {
c[k] = dtmp[k] * mmw / pmv[k];
}
}
doublereal MaskellSolidSolnPhase::standardConcentration(size_t k) const
{
std::vector<doublereal> pmv(m_kk);
getPartialMolarVolumes(&pmv[0]);
doublereal result = 1.0 / pmv[k];
return result;
}
/********************************************************************
* Molar Thermodynamic Properties of the Solution
********************************************************************/
@ -68,21 +90,21 @@ enthalpy_mole() const
{
_updateThermo();
const doublereal h0 = GasConstant * temperature() * mean_X(&m_h0_RT[0]);
const doublereal r = moleFraction(0);
const doublereal r = moleFraction(product_species_index);
const doublereal fmval = fm(r);
return h0 + r * fmval * h_mixing;
}
//=====================================================================================================
doublereal xlogx(doublereal x)
{
return x * std::log(x);
return x * std::log(x);
}
doublereal MaskellSolidSolnPhase::entropy_mole() const
{
_updateThermo();
const doublereal s0 = GasConstant * mean_X(&m_s0_R[0]);
const doublereal r = moleFraction(0);
const doublereal r = moleFraction(product_species_index);
const doublereal fmval = fm(r);
const doublereal rfm = r * fmval;
return s0 + GasConstant * (xlogx(1-rfm) - xlogx(rfm) - xlogx(1-r-rfm) - xlogx((1-fmval)*r) - xlogx(1-r) - xlogx(r));
@ -149,7 +171,7 @@ void MaskellSolidSolnPhase::
getChemPotentials(doublereal* mu) const
{
_updateThermo();
const doublereal r = moleFraction(0);
const doublereal r = moleFraction(product_species_index);
const doublereal pval = p(r);
const doublereal fmval = fm(r);
const doublereal rfm = r * fmval;
@ -157,8 +179,9 @@ getChemPotentials(doublereal* mu) const
const doublereal muDelta = -pval * h_mixing - GasConstant * temperature()
* std::log( (std::pow(1 - rfm, pval) * std::pow(rfm, pval) * std::pow(r - rfm, 1 - pval) * r) /
(std::pow(1 - r - rfm, 1 + pval) * (1 - r)) );
mu[0] = RT * m_g0_RT[0] + muDelta;
mu[1] = RT * m_g0_RT[1] - muDelta;
const int sign = (product_species_index == 0) ? 1 : -1;
mu[0] = RT * m_g0_RT[0] + sign * muDelta;
mu[1] = RT * m_g0_RT[1] - sign * muDelta;
}
void MaskellSolidSolnPhase::
@ -252,6 +275,19 @@ void MaskellSolidSolnPhase::initThermoXML(XML_Node& phaseNode, const std::string
throw CanteraError(subname.c_str(),
"Mixing enthalpy parameter not specified.");
}
if (thNode.hasChild("product_species")) {
XML_Node& scNode = thNode.child("product_species");
std::string product_species_name = scNode.value();
product_species_index = speciesIndex(product_species_name);
if( product_species_index == npos )
{
throw CanteraError(subname.c_str(),
"Species " + product_species_name + " not found.");
}
std::cout << "parsed product_species_index = " << product_species_index << std::endl;
}
} else {
throw CanteraError(subname.c_str(),
"Unspecified thermo model");

2
test/data/MaskellSolidSolnPhase_valid.xml
View file

@ -12,7 +12,7 @@
</speciesArray>
<thermo model="MaskellSolidSolnPhase">
<h_mix>-1000.</h_mix>
<density units="kg/m3"> 1.0 </density>
<product_species>H(s)</product_species>
</thermo>
<kinetics model="none"/>
<state>