From f90580815fdd333312dab5c716d71268c930d085 Mon Sep 17 00:00:00 2001 From: Harry Moffat Date: Thu, 29 May 2008 16:22:15 +0000 Subject: [PATCH] Fixed an error with the calculation of the total gibbs free energy and the damping algorithm when there are unknowns representing interface voltages at electrodes. --- Cantera/src/equil/vcs_Gibbs.cpp | 10 +++++++--- Cantera/src/equil/vcs_solve_TP.cpp | 14 ++++++++++---- 2 files changed, 17 insertions(+), 7 deletions(-) diff --git a/Cantera/src/equil/vcs_Gibbs.cpp b/Cantera/src/equil/vcs_Gibbs.cpp index 49ad78f12..494cb2741 100644 --- a/Cantera/src/equil/vcs_Gibbs.cpp +++ b/Cantera/src/equil/vcs_Gibbs.cpp @@ -51,7 +51,9 @@ namespace VCSnonideal { } for (int kspec = 0; kspec < m_numSpeciesRdc; ++kspec) { - g += molesSp[kspec] * chemPot[kspec]; + if (m_speciesUnknownType[kspec] != VCS_SPECIES_TYPE_INTERFACIALVOLTAGE) { + g += molesSp[kspec] * chemPot[kspec]; + } } return g; @@ -72,8 +74,10 @@ namespace VCSnonideal { double phaseMols = 0.0; for (int kspec = 0; kspec < m_numSpeciesRdc; ++kspec) { if (m_phaseID[kspec] == iphase) { - g += w[kspec] * fe[kspec]; - phaseMols += w[kspec]; + if (m_speciesUnknownType[kspec] != VCS_SPECIES_TYPE_INTERFACIALVOLTAGE) { + g += w[kspec] * fe[kspec]; + phaseMols += w[kspec]; + } } } diff --git a/Cantera/src/equil/vcs_solve_TP.cpp b/Cantera/src/equil/vcs_solve_TP.cpp index e4122a837..dad961379 100644 --- a/Cantera/src/equil/vcs_solve_TP.cpp +++ b/Cantera/src/equil/vcs_solve_TP.cpp @@ -2115,7 +2115,7 @@ namespace VCSnonideal { else { /* * Voltage calculation - * HKM -> Need to check the sign + * Need to check the sign -> This is good for electrons */ dx = m_deltaGRxn_new[irxn]/ m_Faraday_dim; #ifdef DEBUG_MODE @@ -2823,7 +2823,9 @@ namespace VCSnonideal { s2 = 0.0; for (irxn = 0; irxn < m_numRxnRdc; ++irxn) { kspec = irxn + m_numComponents; - s2 += dptr[irxn] * m_deltaMolNumSpecies[kspec]; + if (m_speciesUnknownType[kspec] != VCS_SPECIES_TYPE_INTERFACIALVOLTAGE) { + s2 += dptr[irxn] * m_deltaMolNumSpecies[kspec]; + } } @@ -2834,7 +2836,9 @@ namespace VCSnonideal { dptr = VCS_DATA_PTR(m_deltaGRxn_old); for (irxn = 0; irxn < m_numRxnRdc; ++irxn) { kspec = irxn + m_numComponents; - s1 += dptr[irxn] * m_deltaMolNumSpecies[kspec]; + if (m_speciesUnknownType[kspec] != VCS_SPECIES_TYPE_INTERFACIALVOLTAGE) { + s1 += dptr[irxn] * m_deltaMolNumSpecies[kspec]; + } } #ifdef DEBUG_MODE @@ -2934,7 +2938,9 @@ namespace VCSnonideal { s2 = 0.0; for (irxn = 0; irxn < m_numRxnRdc; ++irxn) { kspec = irxn + m_numComponents; - s2 += dptr[irxn] * m_deltaMolNumSpecies[kspec]; + if (m_speciesUnknownType[kspec] != VCS_SPECIES_TYPE_INTERFACIALVOLTAGE) { + s2 += dptr[irxn] * m_deltaMolNumSpecies[kspec]; + } }