Fleshed out the gradient of actcoeff wrt mole number.

This commit is contained in:
Harry Moffat 2010-08-03 21:07:28 +00:00
parent 673cf207ed
commit b6cf6b3139
5 changed files with 144 additions and 21 deletions

View file

@ -223,6 +223,17 @@ namespace Cantera {
for (; b != end(); ++b, ++xb, ++yb) *b = a*(*xb) + *yb;
}
//! Set all of the entries to zero
inline void zero() {
int nn = m_nrows * m_ncols;
if (nn > 0) {
/*
* Using memset is the fastest way to zero a contiguous
* section of memory.
*/
(void) memset((void *) &m_data[0], 0, nn * sizeof(doublereal));
}
}
//! Allows setting elements using the syntax A(i,j) = x.
/*!

View file

@ -33,7 +33,15 @@ namespace Cantera {
*
*/
GibbsExcessVPSSTP::GibbsExcessVPSSTP() :
VPStandardStateTP()
VPStandardStateTP(),
moleFractions_(0),
lnActCoeff_Scaled_(0),
dlnActCoeffdT_Scaled_(0),
d2lnActCoeffdT2_Scaled_(0),
dlnActCoeffdlnN_Scaled_(0),
dlnActCoeffdlnX_Scaled_(0),
dlnActCoeffdN_Scaled_(0,0),
m_pp(0)
{
}
@ -44,7 +52,15 @@ namespace Cantera {
* has a working copy constructor
*/
GibbsExcessVPSSTP::GibbsExcessVPSSTP(const GibbsExcessVPSSTP &b) :
VPStandardStateTP()
VPStandardStateTP(),
moleFractions_(0),
lnActCoeff_Scaled_(0),
dlnActCoeffdT_Scaled_(0),
d2lnActCoeffdT2_Scaled_(0),
dlnActCoeffdlnN_Scaled_(0),
dlnActCoeffdlnX_Scaled_(0),
dlnActCoeffdN_Scaled_(0,0),
m_pp(0)
{
GibbsExcessVPSSTP::operator=(b);
}
@ -69,6 +85,7 @@ namespace Cantera {
d2lnActCoeffdT2_Scaled_ = b.d2lnActCoeffdT2_Scaled_;
dlnActCoeffdlnX_Scaled_ = b.dlnActCoeffdlnX_Scaled_;
dlnActCoeffdlnN_Scaled_ = b.dlnActCoeffdlnN_Scaled_;
dlnActCoeffdN_Scaled_ = b.dlnActCoeffdN_Scaled_;
m_pp = b.m_pp;
return *this;
@ -328,6 +345,7 @@ namespace Cantera {
d2lnActCoeffdT2_Scaled_.resize(m_kk);
dlnActCoeffdlnX_Scaled_.resize(m_kk);
dlnActCoeffdlnN_Scaled_.resize(m_kk);
dlnActCoeffdN_Scaled_.resize(m_kk, m_kk);
m_pp.resize(m_kk);
}

View file

@ -542,6 +542,8 @@ namespace Cantera {
protected:
// HKM get rid of _Scaled_ prefix
//! Storage for the current values of the mole fractions of the species
/*!
* This vector is kept up-to-date when the setState functions are called.
@ -553,29 +555,36 @@ namespace Cantera {
mutable std::vector<doublereal> moleFractions_;
//! Storage for the current values of the activity coefficients of the
//! species, divided by RT
//! species
mutable std::vector<doublereal> lnActCoeff_Scaled_;
//! Storage for the current derivative values of the
//! gradients with respect to temperature of the
//! log of theactivity coefficients of the species
//! log of the activity coefficients of the species
mutable std::vector<doublereal> dlnActCoeffdT_Scaled_;
//! Storage for the current derivative values of the
//! gradients with respect to temperature of the
//! log of theactivity coefficients of the species
//! log of the activity coefficients of the species
mutable std::vector<doublereal> d2lnActCoeffdT2_Scaled_;
//! Storage for the current derivative values of the
//! gradients with respect to logarithm of the mole fraction of the
//! log of theactivity coefficients of the species
//! log of the activity coefficients of the species @deprecated
mutable std::vector<doublereal> dlnActCoeffdlnN_Scaled_;
//! Storage for the current derivative values of the
//! gradients with respect to logarithm of the mole fraction of the
//! log of theactivity coefficients of the species
//! log of theactivity coefficients of the species @deprecated
mutable std::vector<doublereal> dlnActCoeffdlnX_Scaled_;
//! Storage for the current derivative values of the gradients with respect to logarithm of the species mole number of the
//! log of the activity coefficients of the species
/*!
* dlnActCoeffdN_Scaled_(k, m) is the derivative of ln(gamma_k) wrt ln mole number of species m
*/
mutable Array2D dlnActCoeffdN_Scaled_;
//! Temporary storage space that is fair game
mutable std::vector<doublereal> m_pp;

View file

@ -375,7 +375,7 @@ namespace Cantera {
//return -1.0;
return 0.0;
}
//====================================================================================================================
// Get the array of non-dimensional molar-based activity coefficients at
// the current solution temperature, pressure, and solution concentration.
/*
@ -667,7 +667,7 @@ namespace Cantera {
// been identified.
void MargulesVPSSTP::initLengths() {
m_kk = nSpecies();
dlnActCoeffdN_Scaled_.resize(m_kk, m_kk);
}
/*
@ -759,7 +759,7 @@ namespace Cantera {
double T = temperature();
double RT = GasConstant*T;
fvo_zero_dbl_1(lnActCoeff_Scaled_, m_kk);
for ( iK = 0; iK < m_kk; iK++ ){
for (iK = 0; iK < m_kk; iK++) {
XK = moleFractions_[iK];
for (int i = 0; i < numBinaryInteractions_; i++) {
iA = m_pSpecies_A_ij[i];
@ -786,13 +786,12 @@ namespace Cantera {
*/
void MargulesVPSSTP::s_update_dlnActCoeff_dT() const {
int iA, iB, iK, delAK, delBK;
doublereal XA, XB, XK, g0, g1;
doublereal XA, XB, g0, g1;
doublereal T = temperature();
doublereal RTT = GasConstant*T*T;
fvo_zero_dbl_1(dlnActCoeffdT_Scaled_, m_kk);
fvo_zero_dbl_1(d2lnActCoeffdT2_Scaled_, m_kk);
for ( iK = 0; iK < m_kk; iK++ ){
XK = moleFractions_[iK];
for (iK = 0; iK < m_kk; iK++) {
for (int i = 0; i < numBinaryInteractions_; i++) {
iA = m_pSpecies_A_ij[i];
iB = m_pSpecies_B_ij[i];
@ -923,6 +922,73 @@ namespace Cantera {
}
}
//====================================================================================================================
// Update the derivative of the log of the activity coefficients wrt dlnN
/*
* This function will be called to update the internally stored gradients of the
* logarithm of the activity coefficients. These are used in the determination
* of the diffusion coefficients.
*
* he = X_A X_B(B + C X_B)
*/
void MargulesVPSSTP::s_update_dlnActCoeff_dN() const {
int iA, iB;
doublereal delAK, delBK;
double XA, XB, g0 , g1;
double T = temperature();
double RT = GasConstant*T;
doublereal delAM, delBM;
dlnActCoeffdN_Scaled_.zero();
/*
* Loop over the activity coefficient gamma_k
*/
for (int iK = 0; iK < m_kk; iK++) {
for (int iM = 0; iM < m_kk; iM++) {
for (int i = 0; i < numBinaryInteractions_; i++) {
iA = m_pSpecies_A_ij[i];
iB = m_pSpecies_B_ij[i];
delAK = 0.0;
delBK = 0.0;
delAM = 0.0;
delBM = 0.0;
if (iA==iK) delAK = 1.0;
else if (iB==iK) delBK = 1.0;
if (iA==iM) delAM = 1.0;
else if (iB==iM) delBM = 1.0;
XA = moleFractions_[iA];
XB = moleFractions_[iB];
g0 = (m_HE_b_ij[i] - T * m_SE_b_ij[i]) / RT;
g1 = (m_HE_c_ij[i] - T * m_SE_c_ij[i]) / RT;
double gfac = g0 + g1 * XB;
double gggg = (delBK - XB) * g1;
// all values of dlnActCoeffdN_Scaled_(iK, iM) hare an additional divisor of n_total
dlnActCoeffdN_Scaled_(iK, iM) += gfac * delAK * ( - XB + delBM);
dlnActCoeffdN_Scaled_(iK, iM) += gfac * delBK * ( - XA + delAM);
dlnActCoeffdN_Scaled_(iK, iM) += gfac * (2.0 * XA * XB - delAM * XB - XA * delBM);
dlnActCoeffdN_Scaled_(iK, iM) += (delAK * XB + XA * delBK - XA * XB) * g1 * (-XB + delBM);
dlnActCoeffdN_Scaled_(iK, iM) += gggg * ( - 2.0 * XA * XB + delAM * XB + XA * delBM);
dlnActCoeffdN_Scaled_(iK, iM) += - g1 * XA * XB * (- XB + delBM);
}
}
}
}
//====================================================================================================================
void MargulesVPSSTP::s_update_dlnActCoeff_dlnX() const {
int iA, iB;
@ -950,21 +1016,31 @@ namespace Cantera {
}
}
//====================================================================================================================
void MargulesVPSSTP::getdlnActCoeffdlnN(doublereal *dlnActCoeffdlnN) const {
s_update_dlnActCoeff_dlnN();
for (int k = 0; k < m_kk; k++) {
dlnActCoeffdlnN[k] = dlnActCoeffdlnN_Scaled_[k];
}
}
//====================================================================================================================
void MargulesVPSSTP::getdlnActCoeffdlnX(doublereal *dlnActCoeffdlnX) const {
s_update_dlnActCoeff_dlnX();
for (int k = 0; k < m_kk; k++) {
dlnActCoeffdlnX[k] = dlnActCoeffdlnX_Scaled_[k];
}
}
//====================================================================================================================
void MargulesVPSSTP::getdlnActCoeffdN(const int ld, doublereal *dlnActCoeffdN) const {
s_update_dlnActCoeff_dN();
double *data = & dlnActCoeffdN_Scaled_(0,0);
for (int k = 0; k < m_kk; k++) {
for (int m = 0; m < m_kk; m++) {
dlnActCoeffdN[ld * k + m] = data[m_kk * k + m];
}
}
}
//====================================================================================================================
void MargulesVPSSTP::resizeNumInteractions(const int num) {
numBinaryInteractions_ = num;
m_HE_b_ij.resize(num, 0.0);
@ -984,7 +1060,7 @@ namespace Cantera {
m_pSpecies_B_ij.resize(num, -1);
}
//====================================================================================================================
/*
* Process an XML node called "binaryNeutralSpeciesParameters"

View file

@ -113,7 +113,7 @@ namespace Cantera {
* \f[
* R T \ln( \gamma_k )= \sum_i \left( \left( \delta_{Ai,k} X_{Bi} + \delta_{Bi,k} X_{Ai} - X_{Ai} X_{Bi} \right)
* \left( g^E_{o,i} + g^E_{1,i} X_{Bi} \right) +
* \left( \delta_{Ai,k} - X_{Bi} \right) X_{Ai} X_{Bi} g^E_{1,i} \right)
* \left( \delta_{Bi,k} - X_{Bi} \right) X_{Ai} X_{Bi} g^E_{1,i} \right)
* \f]
* where
* \f$ g^E_{o,i} = h_{o,i} - T s_{o,i} \f$ and \f$ g^E_{1,i} = h_{1,i} - T s_{1,i} \f$
@ -834,9 +834,8 @@ namespace Cantera {
* @param dlnActCoeffdN Output vector of derivatives of the
* log Activity Coefficients. length = m_kk * m_kk
*/
virtual void getdlnActCoeffdN(const int ld, doublereal * const dlnActCoeffdN) const {
err("getdlnActCoeffdN");
}
virtual void getdlnActCoeffdN(const int ld, doublereal * const dlnActCoeffdN) const;
//@}
private:
@ -898,6 +897,14 @@ namespace Cantera {
*/
void s_update_dlnActCoeff_dlnN() const;
//! Update the derivative of the log of the activity coefficients wrt log(moles_m)
/*!
* This function will be called to update the internally storred
* derivative of the natural logarithm of the activity coefficients
* wrt logarithm of the mole number of species
*/
void s_update_dlnActCoeff_dN() const;
private:
//! Error function
@ -961,6 +968,8 @@ namespace Cantera {
//! Entropy term for the quaternary mole fraction interaction of the
//! excess gibbs free energy expression
mutable vector_fp m_VSE_d_ij;
//! vector of species indices representing species A in the interaction
/*!