diff --git a/Cantera/src/thermo/GibbsExcessVPSSTP.cpp b/Cantera/src/thermo/GibbsExcessVPSSTP.cpp index 5ee5c16e8..94b366c79 100644 --- a/Cantera/src/thermo/GibbsExcessVPSSTP.cpp +++ b/Cantera/src/thermo/GibbsExcessVPSSTP.cpp @@ -40,7 +40,7 @@ namespace Cantera { d2lnActCoeffdT2_Scaled_(0), dlnActCoeffdlnN_diag_(0), dlnActCoeffdlnX_diag_(0), - dlnActCoeffdN_Scaled_(0,0), + dlnActCoeffdlnN_(0,0), m_pp(0) { } @@ -59,7 +59,7 @@ namespace Cantera { d2lnActCoeffdT2_Scaled_(0), dlnActCoeffdlnN_diag_(0), dlnActCoeffdlnX_diag_(0), - dlnActCoeffdN_Scaled_(0,0), + dlnActCoeffdlnN_(0,0), m_pp(0) { GibbsExcessVPSSTP::operator=(b); @@ -85,7 +85,7 @@ namespace Cantera { d2lnActCoeffdT2_Scaled_ = b.d2lnActCoeffdT2_Scaled_; dlnActCoeffdlnX_diag_ = b.dlnActCoeffdlnX_diag_; dlnActCoeffdlnN_diag_ = b.dlnActCoeffdlnN_diag_; - dlnActCoeffdN_Scaled_ = b.dlnActCoeffdN_Scaled_; + dlnActCoeffdlnN_ = b.dlnActCoeffdlnN_; m_pp = b.m_pp; return *this; @@ -345,7 +345,7 @@ namespace Cantera { d2lnActCoeffdT2_Scaled_.resize(m_kk); dlnActCoeffdlnX_diag_.resize(m_kk); dlnActCoeffdlnN_diag_.resize(m_kk); - dlnActCoeffdN_Scaled_.resize(m_kk, m_kk); + dlnActCoeffdlnN_.resize(m_kk, m_kk); m_pp.resize(m_kk); } diff --git a/Cantera/src/thermo/GibbsExcessVPSSTP.h b/Cantera/src/thermo/GibbsExcessVPSSTP.h index 63b387da6..d98b95241 100644 --- a/Cantera/src/thermo/GibbsExcessVPSSTP.h +++ b/Cantera/src/thermo/GibbsExcessVPSSTP.h @@ -304,6 +304,28 @@ namespace Cantera { err("getdlnActCoeffdT"); } + //! Get the array of derivatives of the log activity coefficients with respect to the log of the species mole numbers + /*! + * Implementations should take the derivative of the logarithm of the activity coefficient with respect to a + * species log mole number (with all other species mole numbers held constant). The default treatment in the + * %ThermoPhase object is to set this vector to zero. + * + * units = 1 / kmol + * + * dlnActCoeffdlnN[ ld * k + m] will contain the derivative of log act_coeff for the mth + * species with respect to the number of moles of the kth species. + * + * \f[ + * \frac{d \ln(\gamma_m) }{d \ln( n_k ) }\Bigg|_{n_i} + * \f] + * + * @param ld Number of rows in the matrix + * @param dlnActCoeffdlnN Output vector of derivatives of the + * log Activity Coefficients. length = m_kk * m_kk + */ + virtual void getdlnActCoeffdlnN(const int ld, doublereal * const dlnActCoeffdlnN) const { + err(" getdlnActCoeffdlnN: nonzero and nonimplemented"); + } //@} @@ -539,9 +561,9 @@ namespace Cantera { //! 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 + * dlnActCoeffdlnN_(k, m) is the derivative of ln(gamma_k) wrt ln mole number of species m */ - mutable Array2D dlnActCoeffdN_Scaled_; + mutable Array2D dlnActCoeffdlnN_; //! Temporary storage space that is fair game mutable std::vector m_pp; diff --git a/Cantera/src/thermo/IonsFromNeutralVPSSTP.cpp b/Cantera/src/thermo/IonsFromNeutralVPSSTP.cpp index 8aaeb95da..57acfd7ed 100644 --- a/Cantera/src/thermo/IonsFromNeutralVPSSTP.cpp +++ b/Cantera/src/thermo/IonsFromNeutralVPSSTP.cpp @@ -1221,7 +1221,7 @@ namespace Cantera { } return fMax; } - + //==================================================================================================================== /* * initThermoXML() (virtual from ThermoPhase) * Import and initialize a ThermoPhase object @@ -1386,7 +1386,7 @@ namespace Cantera { * have charge conservation. */ } - + //==================================================================================================================== // Update the activity coefficients /* * This function will be called to update the internally storred @@ -1577,7 +1577,7 @@ namespace Cantera { } } - + //==================================================================================================================== /* * This function will be called to update the internally storred * temperature derivative of the natural logarithm of the activity coefficients @@ -1635,7 +1635,7 @@ namespace Cantera { } } - + //==================================================================================================================== /* * This function will be called to update the internally storred * temperature derivative of the natural logarithm of the activity coefficients @@ -1693,7 +1693,70 @@ namespace Cantera { } } + //==================================================================================================================== + // Update the derivative of the log of the activity coefficients + // wrt log(number of moles) - diagonal components + /* + * This function will be called to update the internally storred + * derivative of the natural logarithm of the activity coefficients + * wrt logarithm of the number of moles of given species. + */ + void IonsFromNeutralVPSSTP::s_update_dlnActCoeff_dlnN() const { + int k, icat, jNeut; + doublereal fmij; + dlnActCoeffdlnN_.zero(); + /* + * Get the activity coefficients of the neutral molecules + */ + GibbsExcessVPSSTP *geThermo = dynamic_cast(neutralMoleculePhase_); + + if (!geThermo) { + + return; + } + int nsp_ge = geThermo->nSpecies(); + geThermo->getdlnActCoeffdlnN(nsp_ge, &(dlnActCoeffdlnN_NeutralMolecule_(0,0))); + switch (ionSolnType_) { + case cIonSolnType_PASSTHROUGH: + break; + case cIonSolnType_SINGLEANION: + + // Do the cation list + for (k = 0; k < (int) cationList_.size(); k++) { + //! Get the id for the next cation + icat = cationList_[k]; + jNeut = fm_invert_ionForNeutral[icat]; + fmij = fm_neutralMolec_ions_[icat + jNeut * m_kk]; + dlnActCoeffdlnN_diag_[icat] = dlnActCoeffdlnN_diag_NeutralMolecule_[jNeut]/fmij; + } + + // Do the anion list + icat = anionList_[0]; + jNeut = fm_invert_ionForNeutral[icat]; + dlnActCoeffdlnN_diag_[icat]= 0.0; + + // Do the list of neutral molecules + for (k = 0; k < numPassThroughSpecies_; k++) { + icat = passThroughList_[k]; + jNeut = fm_invert_ionForNeutral[icat]; + dlnActCoeffdlnN_diag_[icat] = dlnActCoeffdlnN_diag_NeutralMolecule_[jNeut]; + } + break; -} + case cIonSolnType_SINGLECATION: + throw CanteraError("IonsFromNeutralVPSSTP::s_update_lnActCoeff", "Unimplemented type"); + break; + case cIonSolnType_MULTICATIONANION: + throw CanteraError("IonsFromNeutralVPSSTP::s_update_lnActCoeff", "Unimplemented type"); + break; + default: + throw CanteraError("IonsFromNeutralVPSSTP::s_update_lnActCoeff", "Unimplemented type"); + break; + } + + } + //==================================================================================================================== +} +//====================================================================================================================== diff --git a/Cantera/src/thermo/IonsFromNeutralVPSSTP.h b/Cantera/src/thermo/IonsFromNeutralVPSSTP.h index f73105263..964891acd 100644 --- a/Cantera/src/thermo/IonsFromNeutralVPSSTP.h +++ b/Cantera/src/thermo/IonsFromNeutralVPSSTP.h @@ -751,6 +751,15 @@ namespace Cantera { */ void s_update_dlnActCoeff_dlnN_diag() const; + //! Update the derivative of the log of the activity coefficients + //! wrt log(number of moles) - diagonal components + /*! + * This function will be called to update the internally storred + * derivative of the natural logarithm of the activity coefficients + * wrt logarithm of the number of moles of given species. + */ + void s_update_dlnActCoeff_dlnN() const; + private: //! Error function @@ -890,6 +899,8 @@ namespace Cantera { mutable std::vector dlnActCoeffdlnX_diag_NeutralMolecule_; mutable std::vector dlnActCoeffdlnN_diag_NeutralMolecule_; + mutable Array2D dlnActCoeffdlnN_NeutralMolecule_; + }; diff --git a/Cantera/src/thermo/MargulesVPSSTP.cpp b/Cantera/src/thermo/MargulesVPSSTP.cpp index bb50b31db..bc2a11b7f 100644 --- a/Cantera/src/thermo/MargulesVPSSTP.cpp +++ b/Cantera/src/thermo/MargulesVPSSTP.cpp @@ -667,7 +667,7 @@ namespace Cantera { // been identified. void MargulesVPSSTP::initLengths() { m_kk = nSpecies(); - dlnActCoeffdN_Scaled_.resize(m_kk, m_kk); + dlnActCoeffdlnN_.resize(m_kk, m_kk); } /* @@ -916,9 +916,25 @@ namespace Cantera { 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; - dlnActCoeffdlnN_diag_[iK] += 2*(delBK-XB)*(g0*(delAK-XA)+g1*(2*(delAK-XA)*XB+XA*(delBK-XB))); + // dlnActCoeffdlnN_diag_[iK] += 2*(delBK-XB)*(g0*(delAK-XA)+g1*(2*(delAK-XA)*XB+XA*(delBK-XB))); + + double gfac = g0 + g1 * XB; + double gggg = (delBK - XB) * g1; + + + dlnActCoeffdlnN_diag_[iK] += gfac * delAK * ( - XB + delBK); + + dlnActCoeffdlnN_diag_[iK] += gfac * delBK * ( - XA + delAK); + + dlnActCoeffdlnN_diag_[iK] += gfac * (2.0 * XA * XB - delAK * XB - XA * delBK); + + dlnActCoeffdlnN_diag_[iK] += (delAK * XB + XA * delBK - XA * XB) * g1 * (-XB + delBK); + + dlnActCoeffdlnN_diag_[iK] += gggg * ( - 2.0 * XA * XB + delAK * XB + XA * delBK); + + dlnActCoeffdlnN_diag_[iK] += - g1 * XA * XB * (- XB + delBK); } - dlnActCoeffdlnN_diag_[iK] = XK*dlnActCoeffdlnN_diag_[iK]-XK; + // dlnActCoeffdlnN_diag_[iK] = XK*dlnActCoeffdlnN_diag_[iK]-XK; } } @@ -929,9 +945,8 @@ namespace Cantera { * 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 { + void MargulesVPSSTP::s_update_dlnActCoeff_dlnN() const { int iA, iB; doublereal delAK, delBK; double XA, XB, g0 , g1; @@ -940,7 +955,7 @@ namespace Cantera { doublereal delAM, delBM; - dlnActCoeffdN_Scaled_.zero(); + dlnActCoeffdlnN_.zero(); /* * Loop over the activity coefficient gamma_k @@ -971,19 +986,18 @@ namespace Cantera { 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); + dlnActCoeffdlnN_(iK, iM) += gfac * delAK * ( - XB + delBM); - dlnActCoeffdN_Scaled_(iK, iM) += gfac * delBK * ( - XA + delAM); + dlnActCoeffdlnN_(iK, iM) += gfac * delBK * ( - XA + delAM); - dlnActCoeffdN_Scaled_(iK, iM) += gfac * (2.0 * XA * XB - delAM * XB - XA * delBM); + dlnActCoeffdlnN_(iK, iM) += gfac * (2.0 * XA * XB - delAM * XB - XA * delBM); - dlnActCoeffdN_Scaled_(iK, iM) += (delAK * XB + XA * delBK - XA * XB) * g1 * (-XB + delBM); + dlnActCoeffdlnN_(iK, iM) += (delAK * XB + XA * delBK - XA * XB) * g1 * (-XB + delBM); - dlnActCoeffdN_Scaled_(iK, iM) += gggg * ( - 2.0 * XA * XB + delAM * XB + XA * delBM); + dlnActCoeffdlnN_(iK, iM) += gggg * ( - 2.0 * XA * XB + delAM * XB + XA * delBM); - dlnActCoeffdN_Scaled_(iK, iM) += - g1 * XA * XB * (- XB + delBM); + dlnActCoeffdlnN_(iK, iM) += - g1 * XA * XB * (- XB + delBM); } } } @@ -1031,12 +1045,12 @@ namespace Cantera { } } //==================================================================================================================== - void MargulesVPSSTP::getdlnActCoeffdN(const int ld, doublereal *dlnActCoeffdN) const { - s_update_dlnActCoeff_dN(); - double *data = & dlnActCoeffdN_Scaled_(0,0); + void MargulesVPSSTP::getdlnActCoeffdlnN(const int ld, doublereal *dlnActCoeffdlnN) const { + s_update_dlnActCoeff_dlnN(); + double *data = & dlnActCoeffdlnN_(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]; + dlnActCoeffdlnN[ld * k + m] = data[m_kk * k + m]; } } } diff --git a/Cantera/src/thermo/MargulesVPSSTP.h b/Cantera/src/thermo/MargulesVPSSTP.h index 8e3f4b198..af2edca83 100644 --- a/Cantera/src/thermo/MargulesVPSSTP.h +++ b/Cantera/src/thermo/MargulesVPSSTP.h @@ -828,25 +828,25 @@ namespace Cantera { virtual void getdlnActCoeffdlnN_diag(doublereal *dlnActCoeffdlnN_diag) const; - //! Get the array of derivatives of the log activity coefficients with respect to the species mole numbers + //! Get the array of derivatives of the log activity coefficients with respect to the ln species mole numbers /*! * Implementations should take the derivative of the logarithm of the activity coefficient with respect to a - * species mole number (with all other species mole numbers held constant) + * log of a species mole number (with all other species mole numbers held constant) * * units = 1 / kmol * - * dlnActCoeffdN[ ld * k + m] will contain the derivative of log act_coeff for the mth - * species with respect to the number of moles of the kth species. + * dlnActCoeffdlnN[ ld * k + m] will contain the derivative of log act_coeff for the mth + * species with respect to the number of moles of the kth species. * * \f[ - * \frac{d \ln(\gamma_m) }{d n_k }\Bigg|_{n_i} + * \frac{d \ln(\gamma_m) }{d \ln( n_k ) }\Bigg|_{n_i} * \f] * * @param ld Number of rows in the matrix - * @param dlnActCoeffdN Output vector of derivatives of the + * @param dlnActCoeffdlnN Output vector of derivatives of the * log Activity Coefficients. length = m_kk * m_kk */ - virtual void getdlnActCoeffdN(const int ld, doublereal * const dlnActCoeffdN) const; + virtual void getdlnActCoeffdlnN(const int ld, doublereal * const dlnActCoeffdlnN) const; //@} @@ -915,7 +915,7 @@ namespace Cantera { * derivative of the natural logarithm of the activity coefficients * wrt logarithm of the mole number of species */ - void s_update_dlnActCoeff_dN() const; + void s_update_dlnActCoeff_dlnN() const; private: diff --git a/Cantera/src/thermo/MolalityVPSSTP.h b/Cantera/src/thermo/MolalityVPSSTP.h index 60ea18feb..4bb286c64 100644 --- a/Cantera/src/thermo/MolalityVPSSTP.h +++ b/Cantera/src/thermo/MolalityVPSSTP.h @@ -793,6 +793,29 @@ namespace Cantera { */ void setState_TPM(doublereal t, doublereal p, const std::string& m); + //! Get the array of derivatives of the log activity coefficients with respect to the log of the species mole numbers + /*! + * Implementations should take the derivative of the logarithm of the activity coefficient with respect to a + * species log mole number (with all other species mole numbers held constant). The default treatment in the + * %ThermoPhase object is to set this vector to zero. + * + * units = 1 / kmol + * + * dlnActCoeffdlnN[ ld * k + m] will contain the derivative of log act_coeff for the mth + * species with respect to the number of moles of the kth species. + * + * \f[ + * \frac{d \ln(\gamma_m) }{d \ln( n_k ) }\Bigg|_{n_i} + * \f] + * + * @param ld Number of rows in the matrix + * @param dlnActCoeffdlnN Output vector of derivatives of the + * log Activity Coefficients. length = m_kk * m_kk + */ + virtual void getdlnActCoeffdlnN(const int ld, doublereal * const dlnActCoeffdlnN) const { + err(" getdlnActCoeffdlnN: nonzero and nonimplemented"); + } + //! returns a summary of the state of the phase as a string /*! * @param show_thermo If true, extra information is printed out diff --git a/Cantera/src/thermo/ThermoPhase.cpp b/Cantera/src/thermo/ThermoPhase.cpp index 20e5deff7..d7952e724 100644 --- a/Cantera/src/thermo/ThermoPhase.cpp +++ b/Cantera/src/thermo/ThermoPhase.cpp @@ -81,7 +81,7 @@ namespace Cantera { */ *this = operator=(right); } - + //==================================================================================================================== /* * operator=() * @@ -135,7 +135,7 @@ namespace Cantera { m_ssConvention = right.m_ssConvention; return *this; } - +//==================================================================================================================== /* * Duplication routine for objects which inherit from * ThermoPhase. @@ -151,7 +151,7 @@ namespace Cantera { ThermoPhase* tp = new ThermoPhase(*this); return tp; } - +//==================================================================================================================== int ThermoPhase::activityConvention() const { return cAC_CONVENTION_MOLAR; } @@ -1031,7 +1031,33 @@ namespace Cantera { } return (m_hasElementPotentials); } - + //==================================================================================================================== + // Get the array of derivatives of the log activity coefficients with respect to the species mole numbers + /* + * Implementations should take the derivative of the logarithm of the activity coefficient with respect to a + * species mole number (with all other species mole numbers held constant) + * + * units = 1 / kmol + * + * dlnActCoeffdN[ ld * k + m] will contain the derivative of log act_coeff for the mth + * species with respect to the number of moles of the kth species. + * + * \f[ + * \frac{d \ln(\gamma_m) }{d n_k }\Bigg|_{n_i} + * \f] + * + * @param ld Number of rows in the matrix + * @param dlnActCoeffdN Output vector of derivatives of the + * log Activity Coefficients. length = m_kk * m_kk + */ + void ThermoPhase::getdlnActCoeffdlnN(const int ld, doublereal * const dlnActCoeffdlnN) const { + for (int m = 0; m < m_kk; m++) { + for (int k = 0; k < m_kk; k++) { + dlnActCoeffdlnN[ld * k + m] = 0.0; + } + } + } + //==================================================================================================================== /* * Format a summary of the mixture state for output. */ @@ -1140,7 +1166,7 @@ namespace Cantera { } return s; } - +//==================================================================================================================== /* * Format a summary of the mixture state for output. */ diff --git a/Cantera/src/thermo/ThermoPhase.h b/Cantera/src/thermo/ThermoPhase.h index fc5f868c5..274047eed 100644 --- a/Cantera/src/thermo/ThermoPhase.h +++ b/Cantera/src/thermo/ThermoPhase.h @@ -2115,27 +2115,26 @@ namespace Cantera { err("getdlnActCoeffdlnN_diag"); } - //! Get the array of derivatives of the log activity coefficients with respect to the species mole numbers + //! Get the array of derivatives of the log activity coefficients with respect to the log of the species mole numbers /*! * Implementations should take the derivative of the logarithm of the activity coefficient with respect to a - * species mole number (with all other species mole numbers held constant) + * species log mole number (with all other species mole numbers held constant). The default treatment in the + * %ThermoPhase object is to set this vector to zero. * * units = 1 / kmol * - * dlnActCoeffdN[ ld * k + m] will contain the derivative of log act_coeff for the mth + * dlnActCoeffdlnN[ ld * k + m] will contain the derivative of log act_coeff for the mth * species with respect to the number of moles of the kth species. * * \f[ - * \frac{d \ln(\gamma_m) }{d n_k }\Bigg|_{n_i} + * \frac{d \ln(\gamma_m) }{d \ln( n_k ) }\Bigg|_{n_i} * \f] * * @param ld Number of rows in the matrix - * @param dlnActCoeffdN Output vector of derivatives of the - * log Activity Coefficients. length = m_kk * m_kk + * @param dlnActCoeffdlnN 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 getdlnActCoeffdlnN(const int ld, doublereal * const dlnActCoeffdlnN) const; /** * @}