Cantera  3.1.0a1
Loading...
Searching...
No Matches
Kinetics.cpp
Go to the documentation of this file.
1/**
2 * @file Kinetics.cpp Declarations for the base class for kinetics managers
3 * (see @ref kineticsmgr and class @link Cantera::Kinetics Kinetics @endlink).
4 *
5 * Kinetics managers calculate rates of progress of species due to
6 * homogeneous or heterogeneous kinetics.
7 */
8
9// This file is part of Cantera. See License.txt in the top-level directory or
10// at https://cantera.org/license.txt for license and copyright information.
11
18#include "cantera/base/global.h"
19#include <unordered_set>
20#include <boost/algorithm/string.hpp>
21
22using namespace std;
23
24namespace Cantera
25{
26
27void Kinetics::checkReactionIndex(size_t i) const
28{
29 if (i >= nReactions()) {
30 throw IndexError("Kinetics::checkReactionIndex", "reactions", i,
31 nReactions()-1);
32 }
33}
34
36{
37 size_t nRxn = nReactions();
38
39 // Stoichiometry managers
43
44 m_rbuf.resize(nRxn);
45
46 // products are created for positive net rate of progress
48 // reactants are destroyed for positive net rate of progress
50
51 m_ready = true;
52}
53
55{
56 if (nReactions() > ii) {
57 throw ArraySizeError("Kinetics::checkReactionArraySize", ii,
58 nReactions());
59 }
60}
61
62void Kinetics::checkPhaseIndex(size_t m) const
63{
64 if (m >= nPhases()) {
65 throw IndexError("Kinetics::checkPhaseIndex", "phase", m, nPhases()-1);
66 }
67}
68
69void Kinetics::checkPhaseArraySize(size_t mm) const
70{
71 if (nPhases() > mm) {
72 throw ArraySizeError("Kinetics::checkPhaseArraySize", mm, nPhases());
73 }
74}
75
77{
78 warn_deprecated("Kinetics::reactionPhaseIndex", "The reacting phase is always "
79 "the first phase. To be removed after Cantera 3.1.");
80 return 0;
81}
82
83shared_ptr<ThermoPhase> Kinetics::reactionPhase() const
84{
85 return m_thermo[0];
86}
87
88void Kinetics::checkSpeciesIndex(size_t k) const
89{
90 if (k >= m_kk) {
91 throw IndexError("Kinetics::checkSpeciesIndex", "species", k, m_kk-1);
92 }
93}
94
96{
97 if (m_kk > kk) {
98 throw ArraySizeError("Kinetics::checkSpeciesArraySize", kk, m_kk);
99 }
100}
101
102pair<size_t, size_t> Kinetics::checkDuplicates(bool throw_err) const
103{
104 //! Map of (key indicating participating species) to reaction numbers
105 map<size_t, vector<size_t>> participants;
106 vector<map<int, double>> net_stoich;
107 std::unordered_set<size_t> unmatched_duplicates;
108 for (size_t i = 0; i < m_reactions.size(); i++) {
109 if (m_reactions[i]->duplicate) {
110 unmatched_duplicates.insert(i);
111 }
112 }
113
114 for (size_t i = 0; i < m_reactions.size(); i++) {
115 // Get data about this reaction
116 unsigned long int key = 0;
117 Reaction& R = *m_reactions[i];
118 net_stoich.emplace_back();
119 map<int, double>& net = net_stoich.back();
120 for (const auto& [name, stoich] : R.reactants) {
121 int k = static_cast<int>(kineticsSpeciesIndex(name));
122 key += k*(k+1);
123 net[-1 -k] -= stoich;
124 }
125 for (const auto& [name, stoich] : R.products) {
126 int k = static_cast<int>(kineticsSpeciesIndex(name));
127 key += k*(k+1);
128 net[1+k] += stoich;
129 }
130
131 // Compare this reaction to others with similar participants
132 vector<size_t>& related = participants[key];
133 for (size_t m : related) {
134 Reaction& other = *m_reactions[m];
135 if (R.duplicate && other.duplicate) {
136 // marked duplicates
137 unmatched_duplicates.erase(i);
138 unmatched_duplicates.erase(m);
139 continue;
140 } else if (R.type() != other.type()) {
141 continue; // different reaction types
142 } else if (R.rate() && other.rate()
143 && R.rate()->type() != other.rate()->type())
144 {
145 continue; // different rate parameterizations
146 }
147 double c = checkDuplicateStoich(net_stoich[i], net_stoich[m]);
148 if (c == 0) {
149 continue; // stoichiometries differ (not by a multiple)
150 } else if (c < 0.0 && !R.reversible && !other.reversible) {
151 continue; // irreversible reactions in opposite directions
152 } else if (R.usesThirdBody() && other.usesThirdBody()) {
153 ThirdBody& tb1 = *(R.thirdBody());
154 ThirdBody& tb2 = *(other.thirdBody());
155 bool thirdBodyOk = true;
156 for (size_t k = 0; k < nTotalSpecies(); k++) {
157 string s = kineticsSpeciesName(k);
158 if (tb1.efficiency(s) * tb2.efficiency(s) != 0.0) {
159 // non-zero third body efficiencies for species `s` in
160 // both reactions
161 thirdBodyOk = false;
162 break;
163 }
164 }
165 if (thirdBodyOk) {
166 continue; // No overlap in third body efficiencies
167 }
168 }
169 if (throw_err) {
170 throw InputFileError("Kinetics::checkDuplicates",
171 R.input, other.input,
172 "Undeclared duplicate reactions detected:\n"
173 "Reaction {}: {}\nReaction {}: {}\n",
174 i+1, other.equation(), m+1, R.equation());
175 } else {
176 return {i,m};
177 }
178 }
179 participants[key].push_back(i);
180 }
181 if (unmatched_duplicates.size()) {
182 size_t i = *unmatched_duplicates.begin();
183 if (throw_err) {
184 throw InputFileError("Kinetics::checkDuplicates",
185 m_reactions[i]->input,
186 "No duplicate found for declared duplicate reaction number {}"
187 " ({})", i, m_reactions[i]->equation());
188 } else {
189 return {i, i};
190 }
191 }
192 return {npos, npos};
193}
194
195double Kinetics::checkDuplicateStoich(map<int, double>& r1, map<int, double>& r2) const
196{
197 std::unordered_set<int> keys; // species keys (k+1 or -k-1)
198 for (auto& [speciesKey, stoich] : r1) {
199 keys.insert(speciesKey);
200 }
201 for (auto& [speciesKey, stoich] : r2) {
202 keys.insert(speciesKey);
203 }
204 int k1 = r1.begin()->first;
205 // check for duplicate written in the same direction
206 double ratio = 0.0;
207 if (r1[k1] && r2[k1]) {
208 ratio = r2[k1]/r1[k1];
209 bool different = false;
210 for (int k : keys) {
211 if ((r1[k] && !r2[k]) ||
212 (!r1[k] && r2[k]) ||
213 (r1[k] && fabs(r2[k]/r1[k] - ratio) > 1.e-8)) {
214 different = true;
215 break;
216 }
217 }
218 if (!different) {
219 return ratio;
220 }
221 }
222
223 // check for duplicate written in the reverse direction
224 if (r1[k1] == 0.0 || r2[-k1] == 0.0) {
225 return 0.0;
226 }
227 ratio = r2[-k1]/r1[k1];
228 for (int k : keys) {
229 if ((r1[k] && !r2[-k]) ||
230 (!r1[k] && r2[-k]) ||
231 (r1[k] && fabs(r2[-k]/r1[k] - ratio) > 1.e-8)) {
232 return 0.0;
233 }
234 }
235 return ratio;
236}
237
238string Kinetics::kineticsSpeciesName(size_t k) const
239{
240 for (size_t n = m_start.size()-1; n != npos; n--) {
241 if (k >= m_start[n]) {
242 return thermo(n).speciesName(k - m_start[n]);
243 }
244 }
245 return "<unknown>";
246}
247
248size_t Kinetics::kineticsSpeciesIndex(const string& nm) const
249{
250 for (size_t n = 0; n < m_thermo.size(); n++) {
251 // Check the ThermoPhase object for a match
252 size_t k = thermo(n).speciesIndex(nm);
253 if (k != npos) {
254 return k + m_start[n];
255 }
256 }
257 return npos;
258}
259
261{
262 for (size_t n = 0; n < m_thermo.size(); n++) {
263 size_t k = thermo(n).speciesIndex(nm);
264 if (k != npos) {
265 return thermo(n);
266 }
267 }
268 throw CanteraError("Kinetics::speciesPhase", "unknown species '{}'", nm);
269}
270
271const ThermoPhase& Kinetics::speciesPhase(const string& nm) const
272{
273 for (const auto& thermo : m_thermo) {
274 if (thermo->speciesIndex(nm) != npos) {
275 return *thermo;
276 }
277 }
278 throw CanteraError("Kinetics::speciesPhase", "unknown species '{}'", nm);
279}
280
281size_t Kinetics::speciesPhaseIndex(size_t k) const
282{
283 for (size_t n = m_start.size()-1; n != npos; n--) {
284 if (k >= m_start[n]) {
285 return n;
286 }
287 }
288 throw CanteraError("Kinetics::speciesPhaseIndex",
289 "illegal species index: {}", k);
290}
291
292double Kinetics::reactantStoichCoeff(size_t kSpec, size_t irxn) const
293{
294 return m_reactantStoich.stoichCoeffs().coeff(kSpec, irxn);
295}
296
297double Kinetics::productStoichCoeff(size_t kSpec, size_t irxn) const
298{
299 return m_productStoich.stoichCoeffs().coeff(kSpec, irxn);
300}
301
303{
304 updateROP();
305 std::copy(m_ropf.begin(), m_ropf.end(), fwdROP);
306}
307
309{
310 updateROP();
311 std::copy(m_ropr.begin(), m_ropr.end(), revROP);
312}
313
315{
316 updateROP();
317 std::copy(m_ropnet.begin(), m_ropnet.end(), netROP);
318}
319
320void Kinetics::getReactionDelta(const double* prop, double* deltaProp) const
321{
322 fill(deltaProp, deltaProp + nReactions(), 0.0);
323 // products add
324 m_productStoich.incrementReactions(prop, deltaProp);
325 // reactants subtract
326 m_reactantStoich.decrementReactions(prop, deltaProp);
327}
328
329void Kinetics::getRevReactionDelta(const double* prop, double* deltaProp) const
330{
331 fill(deltaProp, deltaProp + nReactions(), 0.0);
332 // products add
333 m_revProductStoich.incrementReactions(prop, deltaProp);
334 // reactants subtract
335 m_reactantStoich.decrementReactions(prop, deltaProp);
336}
337
339{
340 updateROP();
341
342 // zero out the output array
343 fill(cdot, cdot + m_kk, 0.0);
344
345 // the forward direction creates product species
346 m_productStoich.incrementSpecies(m_ropf.data(), cdot);
347
348 // the reverse direction creates reactant species
349 m_reactantStoich.incrementSpecies(m_ropr.data(), cdot);
350}
351
353{
354 updateROP();
355
356 fill(ddot, ddot + m_kk, 0.0);
357 // the reverse direction destroys products in reversible reactions
358 m_revProductStoich.incrementSpecies(m_ropr.data(), ddot);
359 // the forward direction destroys reactants
360 m_reactantStoich.incrementSpecies(m_ropf.data(), ddot);
361}
362
364{
365 updateROP();
366
367 fill(net, net + m_kk, 0.0);
368 // products are created for positive net rate of progress
369 m_productStoich.incrementSpecies(m_ropnet.data(), net);
370 // reactants are destroyed for positive net rate of progress
371 m_reactantStoich.decrementSpecies(m_ropnet.data(), net);
372}
373
375{
376 Eigen::Map<Eigen::VectorXd> out(dwdot, m_kk);
377 Eigen::Map<Eigen::VectorXd> buf(m_rbuf.data(), nReactions());
378 // the forward direction creates product species
379 getFwdRatesOfProgress_ddT(buf.data());
380 out = m_productStoich.stoichCoeffs() * buf;
381 // the reverse direction creates reactant species
382 getRevRatesOfProgress_ddT(buf.data());
383 out += m_reactantStoich.stoichCoeffs() * buf;
384}
385
387{
388 Eigen::Map<Eigen::VectorXd> out(dwdot, m_kk);
389 Eigen::Map<Eigen::VectorXd> buf(m_rbuf.data(), nReactions());
390 // the forward direction creates product species
391 getFwdRatesOfProgress_ddP(buf.data());
392 out = m_productStoich.stoichCoeffs() * buf;
393 // the reverse direction creates reactant species
394 getRevRatesOfProgress_ddP(buf.data());
395 out += m_reactantStoich.stoichCoeffs() * buf;
396}
397
399{
400 Eigen::Map<Eigen::VectorXd> out(dwdot, m_kk);
401 Eigen::Map<Eigen::VectorXd> buf(m_rbuf.data(), nReactions());
402 // the forward direction creates product species
403 getFwdRatesOfProgress_ddC(buf.data());
404 out = m_productStoich.stoichCoeffs() * buf;
405 // the reverse direction creates reactant species
406 getRevRatesOfProgress_ddC(buf.data());
407 out += m_reactantStoich.stoichCoeffs() * buf;
408}
409
410Eigen::SparseMatrix<double> Kinetics::creationRates_ddX()
411{
412 Eigen::SparseMatrix<double> jac;
413 // the forward direction creates product species
415 // the reverse direction creates reactant species
417 return jac;
418}
419
420Eigen::SparseMatrix<double> Kinetics::creationRates_ddCi()
421{
422 Eigen::SparseMatrix<double> jac;
423 // the forward direction creates product species
425 // the reverse direction creates reactant species
427 return jac;
428}
429
431{
432 Eigen::Map<Eigen::VectorXd> out(dwdot, m_kk);
433 Eigen::Map<Eigen::VectorXd> buf(m_rbuf.data(), nReactions());
434 // the reverse direction destroys products in reversible reactions
435 getRevRatesOfProgress_ddT(buf.data());
436 out = m_revProductStoich.stoichCoeffs() * buf;
437 // the forward direction destroys reactants
438 getFwdRatesOfProgress_ddT(buf.data());
439 out += m_reactantStoich.stoichCoeffs() * buf;
440}
441
443{
444 Eigen::Map<Eigen::VectorXd> out(dwdot, m_kk);
445 Eigen::Map<Eigen::VectorXd> buf(m_rbuf.data(), nReactions());
446 // the reverse direction destroys products in reversible reactions
447 getRevRatesOfProgress_ddP(buf.data());
448 out = m_revProductStoich.stoichCoeffs() * buf;
449 // the forward direction destroys reactants
450 getFwdRatesOfProgress_ddP(buf.data());
451 out += m_reactantStoich.stoichCoeffs() * buf;
452}
453
455{
456 Eigen::Map<Eigen::VectorXd> out(dwdot, m_kk);
457 Eigen::Map<Eigen::VectorXd> buf(m_rbuf.data(), nReactions());
458 // the reverse direction destroys products in reversible reactions
459 getRevRatesOfProgress_ddC(buf.data());
460 out = m_revProductStoich.stoichCoeffs() * buf;
461 // the forward direction destroys reactants
462 getFwdRatesOfProgress_ddC(buf.data());
463 out += m_reactantStoich.stoichCoeffs() * buf;
464}
465
466Eigen::SparseMatrix<double> Kinetics::destructionRates_ddX()
467{
468 Eigen::SparseMatrix<double> jac;
469 // the reverse direction destroys products in reversible reactions
471 // the forward direction destroys reactants
473 return jac;
474}
475
476Eigen::SparseMatrix<double> Kinetics::destructionRates_ddCi()
477{
478 Eigen::SparseMatrix<double> jac;
479 // the reverse direction destroys products in reversible reactions
481 // the forward direction destroys reactants
483 return jac;
484}
485
487{
488 Eigen::Map<Eigen::VectorXd> out(dwdot, m_kk);
489 Eigen::Map<Eigen::VectorXd> buf(m_rbuf.data(), nReactions());
490 getNetRatesOfProgress_ddT(buf.data());
491 out = m_stoichMatrix * buf;
492}
493
495{
496 Eigen::Map<Eigen::VectorXd> out(dwdot, m_kk);
497 Eigen::Map<Eigen::VectorXd> buf(m_rbuf.data(), nReactions());
498 getNetRatesOfProgress_ddP(buf.data());
499 out = m_stoichMatrix * buf;
500}
501
503{
504 Eigen::Map<Eigen::VectorXd> out(dwdot, m_kk);
505 Eigen::Map<Eigen::VectorXd> buf(m_rbuf.data(), nReactions());
506 getNetRatesOfProgress_ddC(buf.data());
507 out = m_stoichMatrix * buf;
508}
509
510Eigen::SparseMatrix<double> Kinetics::netProductionRates_ddX()
511{
513}
514
515Eigen::SparseMatrix<double> Kinetics::netProductionRates_ddCi()
516{
518}
519
520void Kinetics::addThermo(shared_ptr<ThermoPhase> thermo)
521{
522 // the phase with lowest dimensionality is assumed to be the
523 // phase/interface at which reactions take place
524 if (thermo->nDim() <= m_mindim) {
525 if (!m_thermo.empty()) {
526 throw CanteraError("Kinetics::addThermo",
527 "The reacting (lowest dimensional) phase must be added first.");
528 }
529 m_mindim = thermo->nDim();
530 }
531
532 m_thermo.push_back(thermo);
533 m_phaseindex[m_thermo.back()->name()] = nPhases();
535}
536
538{
539 AnyMap out;
540 string name = KineticsFactory::factory()->canonicalize(kineticsType());
541 if (name != "none") {
542 out["kinetics"] = name;
543 if (nReactions() == 0) {
544 out["reactions"] = "none";
545 }
547 out["skip-undeclared-third-bodies"] = true;
548 }
549 }
550 return out;
551}
552
554{
555 m_kk = 0;
556 m_start.resize(nPhases());
557
558 for (size_t i = 0; i < m_thermo.size(); i++) {
559 m_start[i] = m_kk; // global index of first species of phase i
560 m_kk += m_thermo[i]->nSpecies();
561 }
562 invalidateCache();
563}
564
565bool Kinetics::addReaction(shared_ptr<Reaction> r, bool resize)
566{
567 r->check();
568 r->validate(*this);
569
570 if (m_kk == 0) {
571 init();
572 }
574
575 // Check validity of reaction within the context of the Kinetics object
576 if (!r->checkSpecies(*this)) {
577 // Do not add reaction
578 return false;
579 }
580
581 // For reactions created outside the context of a Kinetics object, the units
582 // of the rate coefficient can't be determined in advance. Do that here.
583 if (r->rate_units.factor() == 0) {
584 r->rate()->setRateUnits(r->calculateRateCoeffUnits(*this));
585 }
586
587 size_t irxn = nReactions(); // index of the new reaction
588
589 // indices of reactant and product species within this Kinetics object
590 vector<size_t> rk, pk;
591
592 // Reactant and product stoichiometric coefficients, such that rstoich[i] is
593 // the coefficient for species rk[i]
594 vector<double> rstoich, pstoich;
595
596 for (const auto& [name, stoich] : r->reactants) {
597 rk.push_back(kineticsSpeciesIndex(name));
598 rstoich.push_back(stoich);
599 }
600
601 for (const auto& [name, stoich] : r->products) {
602 pk.push_back(kineticsSpeciesIndex(name));
603 pstoich.push_back(stoich);
604 }
605
606 // The default order for each reactant is its stoichiometric coefficient,
607 // which can be overridden by entries in the Reaction.orders map. rorder[i]
608 // is the order for species rk[i].
609 vector<double> rorder = rstoich;
610 for (const auto& [name, order] : r->orders) {
611 size_t k = kineticsSpeciesIndex(name);
612 // Find the index of species k within rk
613 auto rloc = std::find(rk.begin(), rk.end(), k);
614 if (rloc != rk.end()) {
615 rorder[rloc - rk.begin()] = order;
616 } else {
617 // If the reaction order involves a non-reactant species, add an
618 // extra term to the reactants with zero stoichiometry so that the
619 // stoichiometry manager can be used to compute the global forward
620 // reaction rate.
621 rk.push_back(k);
622 rstoich.push_back(0.0);
623 rorder.push_back(order);
624 }
625 }
626
627 m_reactantStoich.add(irxn, rk, rorder, rstoich);
628 // product orders = product stoichiometric coefficients
629 m_productStoich.add(irxn, pk, pstoich, pstoich);
630 if (r->reversible) {
631 m_revProductStoich.add(irxn, pk, pstoich, pstoich);
632 }
633
634 m_reactions.push_back(r);
635 m_rfn.push_back(0.0);
636 m_delta_gibbs0.push_back(0.0);
637 m_rkcn.push_back(0.0);
638 m_ropf.push_back(0.0);
639 m_ropr.push_back(0.0);
640 m_ropnet.push_back(0.0);
641 m_perturb.push_back(1.0);
642 m_dH.push_back(0.0);
643
644 if (resize) {
646 } else {
647 m_ready = false;
648 }
649
650 return true;
651}
652
653void Kinetics::modifyReaction(size_t i, shared_ptr<Reaction> rNew)
654{
656 shared_ptr<Reaction>& rOld = m_reactions[i];
657 if (rNew->type() != rOld->type()) {
658 throw CanteraError("Kinetics::modifyReaction",
659 "Reaction types are different: {} != {}.",
660 rOld->type(), rNew->type());
661 }
662
663 if (rNew->rate()->type() != rOld->rate()->type()) {
664 throw CanteraError("Kinetics::modifyReaction",
665 "ReactionRate types are different: {} != {}.",
666 rOld->rate()->type(), rNew->rate()->type());
667 }
668
669 if (rNew->reactants != rOld->reactants) {
670 throw CanteraError("Kinetics::modifyReaction",
671 "Reactants are different: '{}' != '{}'.",
672 rOld->reactantString(), rNew->reactantString());
673 }
674
675 if (rNew->products != rOld->products) {
676 throw CanteraError("Kinetics::modifyReaction",
677 "Products are different: '{}' != '{}'.",
678 rOld->productString(), rNew->productString());
679 }
680 m_reactions[i] = rNew;
681 invalidateCache();
682}
683
684shared_ptr<Reaction> Kinetics::reaction(size_t i)
685{
687 return m_reactions[i];
688}
689
690shared_ptr<const Reaction> Kinetics::reaction(size_t i) const
691{
693 return m_reactions[i];
694}
695
696}
Base class for kinetics managers and also contains the kineticsmgr module documentation (see Kinetics...
Header file for class ThermoPhase, the base class for phases with thermodynamic properties,...
A map of string keys to values whose type can vary at runtime.
Definition AnyMap.h:427
Array size error.
Base class for exceptions thrown by Cantera classes.
An array index is out of range.
Error thrown for problems processing information contained in an AnyMap or AnyValue.
Definition AnyMap.h:738
size_t reactionPhaseIndex() const
Phase where the reactions occur.
Definition Kinetics.cpp:76
virtual void resizeReactions()
Finalize Kinetics object and associated objects.
Definition Kinetics.cpp:35
void checkPhaseIndex(size_t m) const
Check that the specified phase index is in range Throws an exception if m is greater than nPhases()
Definition Kinetics.cpp:62
void checkSpeciesArraySize(size_t mm) const
Check that an array size is at least nSpecies() Throws an exception if kk is less than nSpecies().
Definition Kinetics.cpp:95
virtual void getFwdRatesOfProgress(double *fwdROP)
Return the forward rates of progress of the reactions.
Definition Kinetics.cpp:302
ThermoPhase & thermo(size_t n=0)
This method returns a reference to the nth ThermoPhase object defined in this kinetics mechanism.
Definition Kinetics.h:242
vector< shared_ptr< Reaction > > m_reactions
Vector of Reaction objects represented by this Kinetics manager.
Definition Kinetics.h:1455
void checkSpeciesIndex(size_t k) const
Check that the specified species index is in range Throws an exception if k is greater than nSpecies(...
Definition Kinetics.cpp:88
double checkDuplicateStoich(map< int, double > &r1, map< int, double > &r2) const
Check whether r1 and r2 represent duplicate stoichiometries This function returns a ratio if two reac...
Definition Kinetics.cpp:195
vector< double > m_perturb
Vector of perturbation factors for each reaction's rate of progress vector.
Definition Kinetics.h:1452
vector< double > m_ropf
Forward rate-of-progress for each reaction.
Definition Kinetics.h:1496
shared_ptr< Reaction > reaction(size_t i)
Return the Reaction object for reaction i.
Definition Kinetics.cpp:684
bool m_ready
Boolean indicating whether Kinetics object is fully configured.
Definition Kinetics.h:1444
virtual void getRevReactionDelta(const double *g, double *dg) const
Given an array of species properties 'g', return in array 'dg' the change in this quantity in the rev...
Definition Kinetics.cpp:329
virtual void getNetRatesOfProgress(double *netROP)
Net rates of progress.
Definition Kinetics.cpp:314
vector< size_t > m_start
m_start is a vector of integers specifying the beginning position for the species vector for the n'th...
Definition Kinetics.h:1473
virtual string kineticsType() const
Identifies the Kinetics manager type.
Definition Kinetics.h:144
vector< double > m_rkcn
Reciprocal of the equilibrium constant in concentration units.
Definition Kinetics.h:1493
shared_ptr< ThermoPhase > reactionPhase() const
Return pointer to phase where the reactions occur.
Definition Kinetics.cpp:83
size_t m_kk
The number of species in all of the phases that participate in this kinetics mechanism.
Definition Kinetics.h:1448
virtual pair< size_t, size_t > checkDuplicates(bool throw_err=true) const
Check for unmarked duplicate reactions and unmatched marked duplicates.
Definition Kinetics.cpp:102
virtual void getReactionDelta(const double *property, double *deltaProperty) const
Change in species properties.
Definition Kinetics.cpp:320
void checkPhaseArraySize(size_t mm) const
Check that an array size is at least nPhases() Throws an exception if mm is less than nPhases().
Definition Kinetics.cpp:69
vector< double > m_dH
The enthalpy change for each reaction to calculate Blowers-Masel rates.
Definition Kinetics.h:1505
vector< double > m_ropr
Reverse rate-of-progress for each reaction.
Definition Kinetics.h:1499
size_t nPhases() const
The number of phases participating in the reaction mechanism.
Definition Kinetics.h:184
vector< shared_ptr< ThermoPhase > > m_thermo
m_thermo is a vector of pointers to ThermoPhase objects that are involved with this kinetics operator
Definition Kinetics.h:1467
virtual bool addReaction(shared_ptr< Reaction > r, bool resize=true)
Add a single reaction to the mechanism.
Definition Kinetics.cpp:565
string kineticsSpeciesName(size_t k) const
Return the name of the kth species in the kinetics manager.
Definition Kinetics.cpp:238
virtual void getDestructionRates(double *ddot)
Species destruction rates [kmol/m^3/s or kmol/m^2/s].
Definition Kinetics.cpp:352
AnyMap parameters()
Return the parameters for a phase definition which are needed to reconstruct an identical object usin...
Definition Kinetics.cpp:537
virtual void init()
Prepare the class for the addition of reactions, after all phases have been added.
Definition Kinetics.h:1280
virtual void modifyReaction(size_t i, shared_ptr< Reaction > rNew)
Modify the rate expression associated with a reaction.
Definition Kinetics.cpp:653
vector< double > m_ropnet
Net rate-of-progress for each reaction.
Definition Kinetics.h:1502
virtual double productStoichCoeff(size_t k, size_t i) const
Stoichiometric coefficient of species k as a product in reaction i.
Definition Kinetics.cpp:297
virtual void addThermo(shared_ptr< ThermoPhase > thermo)
Add a phase to the kinetics manager object.
Definition Kinetics.cpp:520
vector< double > m_rbuf
Buffer used for storage of intermediate reaction-specific results.
Definition Kinetics.h:1508
Eigen::SparseMatrix< double > m_stoichMatrix
Net stoichiometry (products - reactants)
Definition Kinetics.h:1440
map< string, size_t > m_phaseindex
Mapping of the phase name to the position of the phase within the kinetics object.
Definition Kinetics.h:1481
size_t m_mindim
number of spatial dimensions of lowest-dimensional phase.
Definition Kinetics.h:1484
StoichManagerN m_productStoich
Stoichiometry manager for the products for each reaction.
Definition Kinetics.h:1434
StoichManagerN m_revProductStoich
Stoichiometry manager for the products of reversible reactions.
Definition Kinetics.h:1437
size_t nReactions() const
Number of reactions in the reaction mechanism.
Definition Kinetics.h:152
virtual void getRevRatesOfProgress(double *revROP)
Return the Reverse rates of progress of the reactions.
Definition Kinetics.cpp:308
bool m_hasUndeclaredThirdBodies
Flag indicating whether reactions include undeclared third bodies.
Definition Kinetics.h:1517
size_t kineticsSpeciesIndex(size_t k, size_t n) const
The location of species k of phase n in species arrays.
Definition Kinetics.h:276
vector< double > m_rfn
Forward rate constant for each reaction.
Definition Kinetics.h:1487
StoichManagerN m_reactantStoich
Stoichiometry manager for the reactants for each reaction.
Definition Kinetics.h:1431
virtual void resizeSpecies()
Resize arrays with sizes that depend on the total number of species.
Definition Kinetics.cpp:553
void checkReactionArraySize(size_t ii) const
Check that an array size is at least nReactions() Throws an exception if ii is less than nReactions()...
Definition Kinetics.cpp:54
size_t nTotalSpecies() const
The total number of species in all phases participating in the kinetics mechanism.
Definition Kinetics.h:254
virtual double reactantStoichCoeff(size_t k, size_t i) const
Stoichiometric coefficient of species k as a reactant in reaction i.
Definition Kinetics.cpp:292
void checkReactionIndex(size_t m) const
Check that the specified reaction index is in range Throws an exception if i is greater than nReactio...
Definition Kinetics.cpp:27
virtual void getNetProductionRates(double *wdot)
Species net production rates [kmol/m^3/s or kmol/m^2/s].
Definition Kinetics.cpp:363
ThermoPhase & speciesPhase(const string &nm)
This function looks up the name of a species and returns a reference to the ThermoPhase object of the...
Definition Kinetics.cpp:260
vector< double > m_delta_gibbs0
Delta G^0 for all reactions.
Definition Kinetics.h:1490
virtual void getCreationRates(double *cdot)
Species creation rates [kmol/m^3/s or kmol/m^2/s].
Definition Kinetics.cpp:338
size_t speciesPhaseIndex(size_t k) const
This function takes as an argument the kineticsSpecies index (that is, the list index in the list of ...
Definition Kinetics.cpp:281
size_t nDim() const
Returns the number of spatial dimensions (1, 2, or 3)
Definition Phase.h:546
string speciesName(size_t k) const
Name of the species with index k.
Definition Phase.cpp:142
size_t speciesIndex(const string &name) const
Returns the index of a species named 'name' within the Phase object.
Definition Phase.cpp:129
Abstract base class which stores data about a reaction and its rate parameterization so that it can b...
Definition Reaction.h:25
shared_ptr< ReactionRate > rate()
Get reaction rate pointer.
Definition Reaction.h:147
bool usesThirdBody() const
Check whether reaction involves third body collider.
Definition Reaction.h:161
shared_ptr< ThirdBody > thirdBody()
Get pointer to third-body handler.
Definition Reaction.h:155
bool reversible
True if the current reaction is reversible. False otherwise.
Definition Reaction.h:126
string type() const
The type of reaction, including reaction rate information.
Definition Reaction.cpp:504
string equation() const
The chemical equation for this reaction.
Definition Reaction.cpp:345
Composition products
Product species and stoichiometric coefficients.
Definition Reaction.h:114
Composition reactants
Reactant species and stoichiometric coefficients.
Definition Reaction.h:111
AnyMap input
Input data used for specific models.
Definition Reaction.h:139
bool duplicate
True if the current reaction is marked as duplicate.
Definition Reaction.h:129
void resizeCoeffs(size_t nSpc, size_t nRxn)
Resize the sparse coefficient matrix.
void add(size_t rxn, const vector< size_t > &k)
Add a single reaction to the list of reactions that this stoichiometric manager object handles.
const Eigen::SparseMatrix< double > & stoichCoeffs() const
Return matrix containing stoichiometric coefficients.
Base class for a phase with thermodynamic properties.
A class for managing third-body efficiencies, including default values.
Definition Reaction.h:192
double efficiency(const string &k) const
Get the third-body efficiency for species k
Definition Reaction.cpp:792
This file contains definitions for utility functions and text for modules, inputfiles and logging,...
virtual Eigen::SparseMatrix< double > fwdRatesOfProgress_ddCi()
Calculate derivatives for forward rates-of-progress with respect to species concentration at constant...
Definition Kinetics.h:812
Eigen::SparseMatrix< double > creationRates_ddCi()
Calculate derivatives for species creation rates with respect to species concentration at constant te...
Definition Kinetics.cpp:420
virtual Eigen::SparseMatrix< double > netRatesOfProgress_ddCi()
Calculate derivatives for net rates-of-progress with respect to species concentration at constant tem...
Definition Kinetics.h:958
void getCreationRates_ddT(double *dwdot)
Calculate derivatives for species creation rates with respect to temperature at constant pressure,...
Definition Kinetics.cpp:374
virtual Eigen::SparseMatrix< double > netRatesOfProgress_ddX()
Calculate derivatives for net rates-of-progress with respect to species mole fractions at constant te...
Definition Kinetics.h:938
Eigen::SparseMatrix< double > destructionRates_ddX()
Calculate derivatives for species destruction rates with respect to species mole fractions at constan...
Definition Kinetics.cpp:466
void getCreationRates_ddC(double *dwdot)
Calculate derivatives for species creation rates with respect to molar concentration at constant temp...
Definition Kinetics.cpp:398
void getDestructionRates_ddP(double *dwdot)
Calculate derivatives for species destruction rates with respect to pressure at constant temperature,...
Definition Kinetics.cpp:442
void getNetProductionRates_ddC(double *dwdot)
Calculate derivatives for species net production rates with respect to molar concentration at constan...
Definition Kinetics.cpp:502
void getDestructionRates_ddT(double *dwdot)
Calculate derivatives for species destruction rates with respect to temperature at constant pressure,...
Definition Kinetics.cpp:430
virtual void getFwdRatesOfProgress_ddP(double *drop)
Calculate derivatives for forward rates-of-progress with respect to pressure at constant temperature,...
Definition Kinetics.h:761
Eigen::SparseMatrix< double > netProductionRates_ddX()
Calculate derivatives for species net production rates with respect to species mole fractions at cons...
Definition Kinetics.cpp:510
Eigen::SparseMatrix< double > creationRates_ddX()
Calculate derivatives for species creation rates with respect to species mole fractions at constant t...
Definition Kinetics.cpp:410
Eigen::SparseMatrix< double > destructionRates_ddCi()
Calculate derivatives for species destruction rates with respect to species concentration at constant...
Definition Kinetics.cpp:476
void getNetProductionRates_ddT(double *dwdot)
Calculate derivatives for species net production rates with respect to temperature at constant pressu...
Definition Kinetics.cpp:486
virtual Eigen::SparseMatrix< double > fwdRatesOfProgress_ddX()
Calculate derivatives for forward rates-of-progress with respect to species mole fractions at constan...
Definition Kinetics.h:792
virtual void getRevRatesOfProgress_ddT(double *drop)
Calculate derivatives for reverse rates-of-progress with respect to temperature at constant pressure,...
Definition Kinetics.h:823
void getCreationRates_ddP(double *dwdot)
Calculate derivatives for species creation rates with respect to pressure at constant temperature,...
Definition Kinetics.cpp:386
virtual void getNetRatesOfProgress_ddT(double *drop)
Calculate derivatives for net rates-of-progress with respect to temperature at constant pressure,...
Definition Kinetics.h:896
virtual void getRevRatesOfProgress_ddP(double *drop)
Calculate derivatives for reverse rates-of-progress with respect to pressure at constant temperature,...
Definition Kinetics.h:834
Eigen::SparseMatrix< double > netProductionRates_ddCi()
Calculate derivatives for species net production rates with respect to species concentration at const...
Definition Kinetics.cpp:515
virtual void getRevRatesOfProgress_ddC(double *drop)
Calculate derivatives for reverse rates-of-progress with respect to molar concentration at constant t...
Definition Kinetics.h:848
void getNetProductionRates_ddP(double *dwdot)
Calculate derivatives for species net production rates with respect to pressure at constant temperatu...
Definition Kinetics.cpp:494
virtual void getNetRatesOfProgress_ddP(double *drop)
Calculate derivatives for net rates-of-progress with respect to pressure at constant temperature,...
Definition Kinetics.h:907
virtual Eigen::SparseMatrix< double > revRatesOfProgress_ddCi()
Calculate derivatives for forward rates-of-progress with respect to species concentration at constant...
Definition Kinetics.h:885
virtual void getFwdRatesOfProgress_ddT(double *drop)
Calculate derivatives for forward rates-of-progress with respect to temperature at constant pressure,...
Definition Kinetics.h:750
void getDestructionRates_ddC(double *dwdot)
Calculate derivatives for species destruction rates with respect to molar concentration at constant t...
Definition Kinetics.cpp:454
virtual void getNetRatesOfProgress_ddC(double *drop)
Calculate derivatives for net rates-of-progress with respect to molar concentration at constant tempe...
Definition Kinetics.h:921
virtual Eigen::SparseMatrix< double > revRatesOfProgress_ddX()
Calculate derivatives for reverse rates-of-progress with respect to species mole fractions at constan...
Definition Kinetics.h:865
virtual void getFwdRatesOfProgress_ddC(double *drop)
Calculate derivatives for forward rates-of-progress with respect to molar concentration at constant t...
Definition Kinetics.h:775
Namespace for the Cantera kernel.
Definition AnyMap.cpp:564
const size_t npos
index returned by functions to indicate "no position"
Definition ct_defs.h:180
void warn_deprecated(const string &source, const AnyBase &node, const string &message)
A deprecation warning for syntax in an input file.
Definition AnyMap.cpp:1926
Contains declarations for string manipulation functions within Cantera.
Various templated functions that carry out common vector and polynomial operations (see Templated Arr...