diff --git a/Cantera/src/equil/MultiPhase.cpp b/Cantera/src/equil/MultiPhase.cpp index f139083c5..81dcb4865 100644 --- a/Cantera/src/equil/MultiPhase.cpp +++ b/Cantera/src/equil/MultiPhase.cpp @@ -230,6 +230,18 @@ namespace Cantera { return sum; } + int MultiPhase::speciesIndex(std::string speciesName, std::string phaseName) { + int p = phaseIndex(phaseName); + if (p < 0) { + throw CanteraError("MultiPhase::speciesIndex", "phase not found: " + phaseName); + } + int k = m_phase[p]->speciesIndex(speciesName); + if (k < 0) { + throw CanteraError("MultiPhase::speciesIndex", "species not found: " + speciesName); + } + return m_spstart[p] + k; + } + /// Net charge of one phase (Coulombs). The net charge is computed as /// \f[ Q_p = N_p \sum_k F z_k X_k \f] /// where the sum runs only over species in phase \a p. @@ -459,6 +471,16 @@ namespace Cantera { } } + void MultiPhase::addSpeciesMoles(const int indexS, const doublereal addedMoles) { + vector_fp tmpMoles(m_nsp, 0.0); + getMoles(DATA_PTR(tmpMoles)); + tmpMoles[indexS] += addedMoles; + if (tmpMoles[indexS] < 0.0) { + tmpMoles[indexS] = 0.0; + } + setMoles(DATA_PTR(tmpMoles)); + } + void MultiPhase::setState_TP(const doublereal T, const doublereal Pres) { if (!m_init) init(); m_temp = T; diff --git a/Cantera/src/equil/MultiPhase.h b/Cantera/src/equil/MultiPhase.h index 224493a1e..69c9812aa 100644 --- a/Cantera/src/equil/MultiPhase.h +++ b/Cantera/src/equil/MultiPhase.h @@ -231,6 +231,21 @@ namespace Cantera { return m_spstart[p] + k; } + //! Return the global index of the species belonging to phase name \c phaseName + //! with species name \c speciesName + /*! + * Returns the index of the global species + * + * @param speciesName Species Name + * @param phaseName Phase Name + * + * @return returns the global index + * + * If the species or phase name is not recognized, this routine throws + * a CanteraError. + */ + int speciesIndex(std::string speciesName, std::string phaseName); + /// Minimum temperature for which all solution phases have /// valid thermo data. Stoichiometric phases are not /// considered, since they may have thermo data only valid for @@ -478,6 +493,13 @@ namespace Cantera { */ void setMoles(const doublereal* n); + + //! Adds moles of a certain species to the mixture + /*! + * + */ + void addSpeciesMoles(const int indexS, const doublereal addedMoles); + //! Retrieves a vector of element abundances /*! * @param elemAbundances Vector of element abundances diff --git a/Cantera/src/equil/vcs_solve.h b/Cantera/src/equil/vcs_solve.h index 4b28f146a..67663950b 100644 --- a/Cantera/src/equil/vcs_solve.h +++ b/Cantera/src/equil/vcs_solve.h @@ -1278,8 +1278,8 @@ private: //! energy by making sure the slope of the following functional stays //! negative: /*! - * The slope of the following functional is equivalent to the slope of the total - * Gibbs free energy of the system: + * The slope of the following functional is equivalent to the slope + * of the total Gibbs free energy of the system: * * d_Gibbs/ds = sum_k( m_deltaGRxn * m_deltaMolNumSpecies[k] ) * @@ -1768,17 +1768,51 @@ public: * to is * kspec = irxn + m_numComponents * - * kspec : 1 -> Major species VCS_SPECIES_MAJOR + * kspec : 2 -> Component species VCS_SPECIES_COMPONENT + * -> deprecated, want to assign -2 to some + * component species. We can already determine + * whether the species is a component from + * its position in the species vector. + * 1 -> Major species VCS_SPECIES_MAJOR * 0 -> Minor species VCS_SPECIES_MINOR - * -1 -> Mole number is zero - * in inactive phase VCS_SPECIES_ZEROEDPHASE - * -2 -> Deleted species in an - * active multicom phase VCS_SPECIES_ZEROEDMS - * -3 -> Mole number is zero - * in a stoich phase - VCS_SPECIES_ZEREODSS - * -4 -> Species is deleted + * -1 -> The species lies in a multicomponent phase + * that exists. Its concentration is currently + * very low, necessitating a different method + * of calculation. + * - VCS_SPECIES_ZEROEDPHASE + * -2 -> The species lies in a multicomponent phase + * which currently doesn't exist. + * Its concentration is currently zero. + * - VCS_SPECIES_ZEROEDMS + * -3 -> Species lies in a single-species phase which + * is currently zereod out. + * - VCS_SPECIES_ZEREODSS + * -4 -> Species has such a small mole fraction it is + * deleted even though its phase may possibly exist. + * The species is believed to have such a small + * mole fraction that it best to throw the + * calculation of it out. It will be added back in + * at the end of the calculation. * - VCS_SPECIES_DELETED - * -> Length equal to number of species + * -5 -> Species refers to an electron in the metal + * The unknown is equal to the interfacial voltage + * drop across the interface on the SHE (standard + * hydroogen electrode) scale (volts). + * - VCS_SPECIES_INTERFACIALVOLTAGE + * -6 -> Species lies in a multicomponent phase that + * is zeroed atm and will stay deleted due to a + * choice from a higher level. + * These species will formally always have zero + * mole numbers in the solution vector. + * - VCS_SPECIES_ZEROEDPHASE + * -7 -> The species lies in a multicomponent phase which + * currently does exist. Its concentration is currently + * identically zero, though the phase exists. Note, this + * is a temporary condition that exists at the start + * of an equilibrium problem. + * The species is soon "birthed" or "deleted". + * - VCS_SPECIES_ACTIVEBUTZERO + * */ std::vector m_speciesStatus; diff --git a/Cantera/src/equil/vcs_solve_TP.cpp b/Cantera/src/equil/vcs_solve_TP.cpp index 3cc1c9e28..bdc9bd20a 100644 --- a/Cantera/src/equil/vcs_solve_TP.cpp +++ b/Cantera/src/equil/vcs_solve_TP.cpp @@ -292,24 +292,25 @@ namespace VCSnonideal { MajorSpeciesHaveConverged = false; /*************************************************************************/ - /************** EVALUATE INITIAL MAJOR-MINOR VECTOR **********************/ + /************** EVALUATE INITIAL SPECIES STATUS VECTOR *******************/ /*************************************************************************/ m_numRxnMinorZeroed = 0; #ifdef DEBUG_MODE if (m_debug_print_lvl >= 2) { - plogf(" --- MAJOR-MINOR decision is reavaluated: All species are minor except for:\n"); + plogf(" --- Species Status decision is reavaluated: All species are minor except for:\n"); } else if (m_debug_print_lvl >= 5) { - plogf(" --- MAJOR-MINOR decision is reavaluated"); + plogf(" --- Species Status decision is reavaluated"); plogendl(); } #endif - for (irxn = 0; irxn < m_numRxnRdc; ++irxn) { - kspec = m_indexRxnToSpecies[irxn]; + for (kspec = 0; kspec < m_numSpeciesTot; ++kspec) { m_speciesStatus[kspec] = vcs_species_type(kspec); #ifdef DEBUG_MODE if (m_debug_print_lvl >= 2) { if (m_speciesStatus[kspec] != VCS_SPECIES_MINOR) { switch (m_speciesStatus[kspec]) { + case VCS_SPECIES_COMPONENT: + break; case VCS_SPECIES_MAJOR: plogf(" --- Major Species : %-s\n", m_speciesName[kspec].c_str()); break; @@ -327,7 +328,8 @@ namespace VCSnonideal { plogf(" --- Deleted-Small Species : %-s\n", m_speciesName[kspec].c_str()); break; case VCS_SPECIES_ACTIVEBUTZERO: - plogf(" --- Zeroed Species in an active MS phase (tmp): %-s\n", m_speciesName[kspec].c_str()); + plogf(" --- Zeroed Species in an active MS phase (tmp): %-s\n", + m_speciesName[kspec].c_str()); break; case VCS_SPECIES_INTERFACIALVOLTAGE: plogf(" --- InterfaceVoltage Species: %-s\n", m_speciesName[kspec].c_str()); @@ -340,8 +342,10 @@ namespace VCSnonideal { } } #endif - if (m_speciesStatus[kspec] != VCS_SPECIES_MAJOR) { - ++m_numRxnMinorZeroed; + if (kspec >= m_numComponents) { + if (m_speciesStatus[kspec] != VCS_SPECIES_MAJOR) { + ++m_numRxnMinorZeroed; + } } } #ifdef DEBUG_MODE @@ -1619,7 +1623,7 @@ namespace VCSnonideal { if (m_speciesStatus[kspec] != VCS_SPECIES_MINOR) { if (m_speciesStatus[kspec] == VCS_SPECIES_MAJOR) { plogf(" --- Noncomponent turned from major to minor: "); - } else if (m_speciesStatus[kspec] == VCS_SPECIES_COMPONENT) { + } else if (kspec < m_numComponents) { plogf(" --- Component turned into a minor species: "); } else { plogf(" --- Zeroed Species turned into a " @@ -1636,7 +1640,7 @@ namespace VCSnonideal { if (m_debug_print_lvl >= 2) { if (m_speciesStatus[kspec] == VCS_SPECIES_MINOR) { plogf(" --- Noncomponent turned from minor to major: "); - } else if (m_speciesStatus[kspec] == VCS_SPECIES_COMPONENT) { + } else if (kspec < m_numComponents) { plogf(" --- Component turned into a major: "); } else { plogf(" --- Noncomponent turned from zeroed to major: ");