diff --git a/include/cantera/kinetics/ThirdBodyCalc.h b/include/cantera/kinetics/ThirdBodyCalc.h index 8f793950f8..2f50286084 100644 --- a/include/cantera/kinetics/ThirdBodyCalc.h +++ b/include/cantera/kinetics/ThirdBodyCalc.h @@ -32,6 +32,22 @@ class ThirdBodyCalc m_eff.back().push_back(eff.second - dflt); } } + + void modifyEfficiencies(size_t rxnNumber, const std::map& enhanced, + double dflt=1.0) { + auto it = find(m_reaction_index.begin(), m_reaction_index.end(), rxnNumber); + if (it != m_reaction_index.end()) { + int index = it - m_reaction_index.begin(); + for (const auto& eff : enhanced) { + auto si = find(m_species[index].begin(), + m_species[index].end(), eff.first); + if (it != m_reaction_index.end()) { + int j = si - m_species[index].begin(); + m_eff[index][j] = (eff.second - m_default[index]); + } + } + } + } void update(const vector_fp& conc, double ctot, double* work) { for (size_t i = 0; i < m_species.size(); i++) { diff --git a/src/kinetics/GasKinetics.cpp b/src/kinetics/GasKinetics.cpp index 3d8042e299..f109959869 100644 --- a/src/kinetics/GasKinetics.cpp +++ b/src/kinetics/GasKinetics.cpp @@ -395,6 +395,16 @@ void GasKinetics::modifyReaction(size_t i, shared_ptr rNew) void GasKinetics::modifyThreeBodyReaction(size_t i, ThreeBodyReaction2& r) { m_rates.replace(i, r.rate); + + map efficiencies; + for (const auto& eff : r.third_body.efficiencies) { + size_t k = kineticsSpeciesIndex(eff.first); + if (k != npos) { + efficiencies[k] = eff.second; + } + } + m_3b_concm.modifyEfficiencies(i, efficiencies, + r.third_body.default_efficiency); } void GasKinetics::modifyFalloffReaction(size_t i, FalloffReaction& r) @@ -403,6 +413,16 @@ void GasKinetics::modifyFalloffReaction(size_t i, FalloffReaction& r) m_falloff_high_rates.replace(iFall, r.high_rate); m_falloff_low_rates.replace(iFall, r.low_rate); m_falloffn.replace(iFall, r.falloff); + + map efficiencies; + for (const auto& eff : r.third_body.efficiencies) { + size_t k = kineticsSpeciesIndex(eff.first); + if (k != npos) { + efficiencies[k] = eff.second; + } + } + m_falloff_concm.modifyEfficiencies(iFall, efficiencies, + r.third_body.default_efficiency); } void GasKinetics::modifyPlogReaction(size_t i, PlogReaction2& r) diff --git a/test/kinetics/kineticsFromScratch.cpp b/test/kinetics/kineticsFromScratch.cpp index 969628d154..e567415ddd 100644 --- a/test/kinetics/kineticsFromScratch.cpp +++ b/test/kinetics/kineticsFromScratch.cpp @@ -46,6 +46,24 @@ class KineticsFromScratch : public testing::Test kin_ref->getRevRateConstants(&k_ref[0]); EXPECT_DOUBLE_EQ(k_ref[iRef], k[0]); } + + void check_uneq_rates(int iRef) { + ASSERT_EQ((size_t) 1, kin.nReactions()); + + std::string X = "O:0.02 H2:0.2 O2:0.5 H:0.03 OH:0.05 H2O:0.1 HO2:0.01"; + p.setState_TPX(1200, 5*OneAtm, X); + p_ref.setState_TPX(1200, 5*OneAtm, X); + + vector_fp k(1), k_ref(kin_ref->nReactions()); + + kin.getFwdRateConstants(&k[0]); + kin_ref->getFwdRateConstants(&k_ref[0]); + ASSERT_NE(k_ref[iRef], k[0]); + + kin.getRevRateConstants(&k[0]); + kin_ref->getRevRateConstants(&k_ref[0]); + ASSERT_NE(k_ref[iRef], k[0]); + } }; TEST_F(KineticsFromScratch, add_elementary_reaction) @@ -118,10 +136,18 @@ TEST_F(KineticsFromScratch, add_falloff_reaction) Arrhenius low_rate(2.3e12, -0.9, -1700.0 / GasConst_cal_mol_K); vector_fp falloff_params { 0.7346, 94.0, 1756.0, 5182.0 }; ThirdBody tbody; - tbody.efficiencies = parseCompString("AR:0.7 H2:2.0 H2O:6.0"); + tbody.efficiencies = parseCompString("AR:0.7 H2:2.0 H2O:60.0"); auto R = make_shared(reac, prod, low_rate, high_rate, tbody); R->falloff = newFalloff("Troe", falloff_params); + + ThirdBody tbody_correct; + tbody_correct.efficiencies = parseCompString("AR:0.7 H2:2.0 H2O:6.0"); + auto R_correct = make_shared(reac, prod, low_rate, high_rate, tbody_correct); + R_correct->falloff = newFalloff("Troe", falloff_params); + kin.addReaction(R); + check_uneq_rates(2); + kin.modifyReaction(0, R_correct); check_rates(2); }