Cantera  3.0.0
Loading...
Searching...
No Matches
Kinetics.h
Go to the documentation of this file.
1/**
2 * @file Kinetics.h
3 * Base class for kinetics managers and also contains the kineticsmgr
4 * module documentation (see @ref kineticsmgr and class
5 * @link Cantera::Kinetics Kinetics@endlink).
6 */
7
8// This file is part of Cantera. See License.txt in the top-level directory or
9// at https://cantera.org/license.txt for license and copyright information.
10
11#ifndef CT_KINETICS_H
12#define CT_KINETICS_H
13
14#include "StoichManager.h"
16
17namespace Cantera
18{
19
20class ThermoPhase;
21class Reaction;
22class Solution;
23class AnyMap;
24
25//! @defgroup derivGroup Derivative Calculations
26//! @details Methods for calculating analytical and/or numerical derivatives.
27
28/**
29 * @defgroup chemkinetics Chemical Kinetics
30 */
31
32//! @defgroup reactionGroup Reactions and Reaction Rates
33//! Classes for handling reactions and reaction rates.
34//! @ingroup chemkinetics
35
36//! @defgroup kineticsmgr Kinetics Managers
37//! Classes implementing models for chemical kinetics.
38//! @section kinmodman Models and Managers
39//!
40//! A kinetics manager is a C++ class that implements a kinetics model; a
41//! kinetics model is a set of mathematical equation describing how various
42//! kinetic quantities are to be computed -- reaction rates, species production
43//! rates, etc. Many different kinetics models might be defined to handle
44//! different types of kinetic processes. For example, one kinetics model might
45//! use expressions valid for elementary reactions in ideal gas mixtures. It
46//! might, for example, require the reaction orders to be integral and equal to
47//! the forward stoichiometric coefficients, require that each reaction be
48//! reversible with a reverse rate satisfying detailed balance, include
49//! pressure-dependent unimolecular reactions, etc. Another kinetics model might
50//! be designed for heterogeneous chemistry at interfaces, and might allow
51//! empirical reaction orders, coverage-dependent activation energies,
52//! irreversible reactions, and include effects of potential differences across
53//! the interface on reaction rates.
54//!
55//! A kinetics manager implements a kinetics model. Since the model equations
56//! may be complex and expensive to evaluate, a kinetics manager may adopt
57//! various strategies to 'manage' the computation and evaluate the expressions
58//! efficiently. For example, if there are rate coefficients or other quantities
59//! that depend only on temperature, a manager class may choose to store these
60//! quantities internally, and re-evaluate them only when the temperature has
61//! actually changed. Or a manager designed for use with reaction mechanisms
62//! with a few repeated activation energies might precompute the terms @f$
63//! \exp(-E/RT) @f$, instead of evaluating the exponential repeatedly for each
64//! reaction. There are many other possible 'management styles', each of which
65//! might be better suited to some reaction mechanisms than others.
66//!
67//! But however a manager structures the internal computation, the tasks the
68//! manager class must perform are, for the most part, the same. It must be able
69//! to compute reaction rates, species production rates, equilibrium constants,
70//! etc. Therefore, all kinetics manager classes should have a common set of
71//! public methods, but differ in how they implement these methods.
72//!
73//! A kinetics manager computes reaction rates of progress, species production
74//! rates, equilibrium constants, and similar quantities for a reaction
75//! mechanism. All kinetics manager classes derive from class Kinetics, which
76//! defines a common public interface for all kinetics managers. Each derived
77//! class overloads the virtual methods of Kinetics to implement a particular
78//! kinetics model.
79//!
80//! For example, class GasKinetics implements reaction rate expressions
81//! appropriate for homogeneous reactions in ideal gas mixtures, and class
82//! InterfaceKinetics implements expressions appropriate for heterogeneous
83//! mechanisms at interfaces, including how to handle reactions involving
84//! charged species of phases with different electric potentials --- something
85//! that class GasKinetics doesn't deal with at all.
86//!
87//! Many of the methods of class Kinetics write into arrays the values of some
88//! quantity for each species, for example the net production rate. These
89//! methods always write the results into flat arrays, ordered by phase in the
90//! order the phase was added, and within a phase in the order the species were
91//! added to the phase (which is the same ordering as in the input file).
92//! Example: suppose a heterogeneous mechanism involves three phases -- a bulk
93//! phase 'a', another bulk phase 'b', and the surface phase 'a:b' at the a/b
94//! interface. Phase 'a' contains 12 species, phase 'b' contains 3, and at the
95//! interface there are 5 adsorbed species defined in phase 'a:b'. Then methods
96//! like getNetProductionRates(double* net) will write and output array of
97//! length 20, beginning at the location pointed to by 'net'. The first 12
98//! values will be the net production rates for all 12 species of phase 'a'
99//! (even if some do not participate in the reactions), the next 3 will be for
100//! phase 'b', and finally the net production rates for the surface species will
101//! occupy the last 5 locations.
102//! @ingroup chemkinetics
103
104//! @defgroup rateEvaluators Rate Evaluators
105//! These classes are used to evaluate the rates of reactions.
106//! @ingroup chemkinetics
107
108
109//! Public interface for kinetics managers.
110/*!
111 * This class serves as a base class to derive 'kinetics managers', which are
112 * classes that manage homogeneous chemistry within one phase, or heterogeneous
113 * chemistry at one interface. The virtual methods of this class are meant to be
114 * overloaded in subclasses. The non-virtual methods perform generic functions
115 * and are implemented in Kinetics. They should not be overloaded. Only those
116 * methods required by a subclass need to be overloaded; the rest will throw
117 * exceptions if called.
118 *
119 * When the nomenclature "kinetics species index" is used below, this means that
120 * the species index ranges over all species in all phases handled by the
121 * kinetics manager.
122 *
123 * @ingroup kineticsmgr
124 */
126{
127public:
128 //! @name Constructors and General Information about Mechanism
129 //! @{
130
131 //! Default constructor.
132 Kinetics() = default;
133
134 virtual ~Kinetics() = default;
135
136 //! Kinetics objects are not copyable or assignable
137 Kinetics(const Kinetics&) = delete;
138 Kinetics& operator=(const Kinetics&)= delete;
139
140 //! Identifies the Kinetics manager type.
141 //! Each class derived from Kinetics should override this method to return
142 //! a meaningful identifier.
143 //! @since Starting in %Cantera 3.0, the name returned by this method corresponds
144 //! to the canonical name used in the YAML input format.
145 virtual string kineticsType() const {
146 return "none";
147 }
148
149 //! Finalize Kinetics object and associated objects
150 virtual void resizeReactions();
151
152 //! Number of reactions in the reaction mechanism.
153 size_t nReactions() const {
154 return m_reactions.size();
155 }
156
157 //! Check that the specified reaction index is in range
158 //! Throws an exception if i is greater than nReactions()
159 void checkReactionIndex(size_t m) const;
160
161 //! Check that an array size is at least nReactions()
162 //! Throws an exception if ii is less than nReactions(). Used before calls
163 //! which take an array pointer.
164 void checkReactionArraySize(size_t ii) const;
165
166 //! Check that the specified species index is in range
167 //! Throws an exception if k is greater than nSpecies()-1
168 void checkSpeciesIndex(size_t k) const;
169
170 //! Check that an array size is at least nSpecies()
171 //! Throws an exception if kk is less than nSpecies(). Used before calls
172 //! which take an array pointer.
173 void checkSpeciesArraySize(size_t mm) const;
174
175 //! @}
176 //! @name Information/Lookup Functions about Phases and Species
177 //! @{
178
179 /**
180 * The number of phases participating in the reaction mechanism. For a
181 * homogeneous reaction mechanism, this will always return 1, but for a
182 * heterogeneous mechanism it will return the total number of phases in the
183 * mechanism.
184 */
185 size_t nPhases() const {
186 return m_thermo.size();
187 }
188
189 //! Check that the specified phase index is in range
190 //! Throws an exception if m is greater than nPhases()
191 void checkPhaseIndex(size_t m) const;
192
193 //! Check that an array size is at least nPhases()
194 //! Throws an exception if mm is less than nPhases(). Used before calls
195 //! which take an array pointer.
196 void checkPhaseArraySize(size_t mm) const;
197
198 /**
199 * Return the phase index of a phase in the list of phases defined within
200 * the object.
201 *
202 * @param ph string name of the phase
203 *
204 * If a -1 is returned, then the phase is not defined in the Kinetics
205 * object.
206 */
207 size_t phaseIndex(const string& ph) const {
208 if (m_phaseindex.find(ph) == m_phaseindex.end()) {
209 return npos;
210 } else {
211 return m_phaseindex.at(ph) - 1;
212 }
213 }
214
215 /**
216 * This returns the integer index of the phase which has ThermoPhase type
217 * cSurf. For heterogeneous mechanisms, this identifies the one surface
218 * phase. For homogeneous mechanisms, this returns -1.
219 *
220 * @deprecated To be removed after %Cantera 3.0. Use reactionPhaseIndex instead.
221 */
222 size_t surfacePhaseIndex() const;
223
224 /**
225 * Phase where the reactions occur. For heterogeneous mechanisms, one of
226 * the phases in the list of phases represents the 2D interface or 1D edge
227 * at which the reactions take place. This method returns the index of the
228 * phase with the smallest spatial dimension (1, 2, or 3) among the list
229 * of phases. If there is more than one, the index of the first one is
230 * returned. For homogeneous mechanisms, the value 0 is returned.
231 *
232 * @deprecated Starting in %Cantera 3.0, the reacting phase will always be the
233 * first phase in the InterfaceKinetics object. To be removed after %Cantera 3.1.
234 */
235 size_t reactionPhaseIndex() const {
236 return m_rxnphase;
237 }
238
239 /**
240 * Return pointer to phase where the reactions occur.
241 * @since New in %Cantera 3.0
242 */
243 shared_ptr<ThermoPhase> reactionPhase() const;
244
245 /**
246 * This method returns a reference to the nth ThermoPhase object defined
247 * in this kinetics mechanism. It is typically used so that member
248 * functions of the ThermoPhase object may be called. For homogeneous
249 * mechanisms, there is only one object, and this method can be called
250 * without an argument to access it.
251 *
252 * @param n Index of the ThermoPhase being sought.
253 */
254 ThermoPhase& thermo(size_t n=0) {
255 return *m_thermo[n];
256 }
257 const ThermoPhase& thermo(size_t n=0) const {
258 return *m_thermo[n];
259 }
260
261 /**
262 * The total number of species in all phases participating in the kinetics
263 * mechanism. This is useful to dimension arrays for use in calls to
264 * methods that return the species production rates, for example.
265 */
266 size_t nTotalSpecies() const {
267 return m_kk;
268 }
269
270 /**
271 * The location of species k of phase n in species arrays. Kinetics manager
272 * classes return species production rates in flat arrays, with the species
273 * of each phases following one another, in the order the phases were added.
274 * This method is useful to find the value for a particular species of a
275 * particular phase in arrays returned from methods like getCreationRates
276 * that return an array of species-specific quantities.
277 *
278 * Example: suppose a heterogeneous mechanism involves three phases. The
279 * first contains 12 species, the second 26, and the third 3. Then species
280 * arrays must have size at least 41, and positions 0 - 11 are the values
281 * for the species in the first phase, positions 12 - 37 are the values for
282 * the species in the second phase, etc. Then kineticsSpeciesIndex(7, 0) =
283 * 7, kineticsSpeciesIndex(4, 1) = 16, and kineticsSpeciesIndex(2, 2) = 40.
284 *
285 * @param k species index
286 * @param n phase index for the species
287 */
288 size_t kineticsSpeciesIndex(size_t k, size_t n) const {
289 return m_start[n] + k;
290 }
291
292 //! Return the name of the kth species in the kinetics manager.
293 /*!
294 * k is an integer from 0 to ktot - 1, where ktot is the number of
295 * species in the kinetics manager, which is the sum of the number of
296 * species in all phases participating in the kinetics manager. If k is
297 * out of bounds, the string "<unknown>" is returned.
298 *
299 * @param k species index
300 */
301 string kineticsSpeciesName(size_t k) const;
302
303 /**
304 * This routine will look up a species number based on the input
305 * string nm. The lookup of species will occur for all phases
306 * listed in the kinetics object.
307 *
308 * return
309 * - If a match is found, the position in the species list is returned.
310 * - If no match is found, the value -1 is returned.
311 *
312 * @param nm Input string name of the species
313 */
314 size_t kineticsSpeciesIndex(const string& nm) const;
315
316 /**
317 * This routine will look up a species number based on the input
318 * string nm. The lookup of species will occur in the specified
319 * phase of the object, or all phases if ph is "<any>".
320 *
321 * return
322 * - If a match is found, the position in the species list is returned.
323 * - If no match is found, the value npos (-1) is returned.
324 *
325 * @param nm Input string name of the species
326 * @param ph Input string name of the phase.
327 * @deprecated To be removed after %Cantera 3.0. Species names should be unique
328 * across all phases.
329 */
330 size_t kineticsSpeciesIndex(const string& nm,
331 const string& ph) const;
332
333 /**
334 * This function looks up the name of a species and returns a
335 * reference to the ThermoPhase object of the phase where the species
336 * resides. Will throw an error if the species doesn't match.
337 *
338 * @param nm String containing the name of the species.
339 */
340 ThermoPhase& speciesPhase(const string& nm);
341 const ThermoPhase& speciesPhase(const string& nm) const;
342
343 /**
344 * This function takes as an argument the kineticsSpecies index
345 * (that is, the list index in the list of species in the kinetics
346 * manager) and returns the species' owning ThermoPhase object.
347 *
348 * @param k Species index
349 */
351 return thermo(speciesPhaseIndex(k));
352 }
353
354 /**
355 * This function takes as an argument the kineticsSpecies index (that is, the
356 * list index in the list of species in the kinetics manager) and returns
357 * the index of the phase owning the species.
358 *
359 * @param k Species index
360 */
361 size_t speciesPhaseIndex(size_t k) const;
362
363 //! @}
364 //! @name Reaction Rates Of Progress
365 //! @{
366
367 //! Return the forward rates of progress of the reactions
368 /*!
369 * Forward rates of progress. Return the forward rates of
370 * progress in array fwdROP, which must be dimensioned at
371 * least as large as the total number of reactions.
372 *
373 * @param fwdROP Output vector containing forward rates
374 * of progress of the reactions. Length: nReactions().
375 */
376 virtual void getFwdRatesOfProgress(double* fwdROP);
377
378 //! Return the Reverse rates of progress of the reactions
379 /*!
380 * Return the reverse rates of progress in array revROP, which must be
381 * dimensioned at least as large as the total number of reactions.
382 *
383 * @param revROP Output vector containing reverse rates
384 * of progress of the reactions. Length: nReactions().
385 */
386 virtual void getRevRatesOfProgress(double* revROP);
387
388 /**
389 * Net rates of progress. Return the net (forward - reverse) rates of
390 * progress in array netROP, which must be dimensioned at least as large
391 * as the total number of reactions.
392 *
393 * @param netROP Output vector of the net ROP. Length: nReactions().
394 */
395 virtual void getNetRatesOfProgress(double* netROP);
396
397 //! Return a vector of Equilibrium constants.
398 /*!
399 * Return the equilibrium constants of the reactions in concentration
400 * units in array kc, which must be dimensioned at least as large as the
401 * total number of reactions.
402 *
403 * @f[
404 * Kc_i = \exp [ \Delta G_{ss,i} ] \prod(Cs_k) \exp(\sum_k \nu_{k,i} F \phi_n)
405 * @f]
406 *
407 * @param kc Output vector containing the equilibrium constants.
408 * Length: nReactions().
409 */
410 virtual void getEquilibriumConstants(double* kc) {
411 throw NotImplementedError("Kinetics::getEquilibriumConstants");
412 }
413
414 /**
415 * Change in species properties. Given an array of molar species property
416 * values @f$ z_k, k = 1, \dots, K @f$, return the array of reaction values
417 * @f[
418 * \Delta Z_i = \sum_k \nu_{k,i} z_k, i = 1, \dots, I.
419 * @f]
420 * For example, if this method is called with the array of standard-state
421 * molar Gibbs free energies for the species, then the values returned in
422 * array @c deltaProperty would be the standard-state Gibbs free energies of
423 * reaction for each reaction.
424 *
425 * @param property Input vector of property value. Length: #m_kk.
426 * @param deltaProperty Output vector of deltaRxn. Length: nReactions().
427 */
428 virtual void getReactionDelta(const double* property, double* deltaProperty) const;
429
430 /**
431 * Given an array of species properties 'g', return in array 'dg' the
432 * change in this quantity in the reversible reactions. Array 'g' must
433 * have a length at least as great as the number of species, and array
434 * 'dg' must have a length as great as the total number of reactions.
435 * This method only computes 'dg' for the reversible reactions, and the
436 * entries of 'dg' for the irreversible reactions are unaltered. This is
437 * primarily designed for use in calculating reverse rate coefficients
438 * from thermochemistry for reversible reactions.
439 */
440 virtual void getRevReactionDelta(const double* g, double* dg) const;
441
442 //! Return the vector of values for the reaction Gibbs free energy change.
443 /*!
444 * (virtual from Kinetics.h)
445 * These values depend upon the concentration of the solution.
446 *
447 * units = J kmol-1
448 *
449 * @param deltaG Output vector of deltaG's for reactions Length:
450 * nReactions().
451 */
452 virtual void getDeltaGibbs(double* deltaG) {
453 throw NotImplementedError("Kinetics::getDeltaGibbs");
454 }
455
456 //! Return the vector of values for the reaction electrochemical free
457 //! energy change.
458 /*!
459 * These values depend upon the concentration of the solution and the
460 * voltage of the phases
461 *
462 * units = J kmol-1
463 *
464 * @param deltaM Output vector of deltaM's for reactions Length:
465 * nReactions().
466 */
467 virtual void getDeltaElectrochemPotentials(double* deltaM) {
468 throw NotImplementedError("Kinetics::getDeltaElectrochemPotentials");
469 }
470
471 /**
472 * Return the vector of values for the reactions change in enthalpy.
473 * These values depend upon the concentration of the solution.
474 *
475 * units = J kmol-1
476 *
477 * @param deltaH Output vector of deltaH's for reactions Length:
478 * nReactions().
479 */
480 virtual void getDeltaEnthalpy(double* deltaH) {
481 throw NotImplementedError("Kinetics::getDeltaEnthalpy");
482 }
483
484 /**
485 * Return the vector of values for the reactions change in entropy. These
486 * values depend upon the concentration of the solution.
487 *
488 * units = J kmol-1 Kelvin-1
489 *
490 * @param deltaS Output vector of deltaS's for reactions Length:
491 * nReactions().
492 */
493 virtual void getDeltaEntropy(double* deltaS) {
494 throw NotImplementedError("Kinetics::getDeltaEntropy");
495 }
496
497 /**
498 * Return the vector of values for the reaction standard state Gibbs free
499 * energy change. These values don't depend upon the concentration of the
500 * solution.
501 *
502 * units = J kmol-1
503 *
504 * @param deltaG Output vector of ss deltaG's for reactions Length:
505 * nReactions().
506 */
507 virtual void getDeltaSSGibbs(double* deltaG) {
508 throw NotImplementedError("Kinetics::getDeltaSSGibbs");
509 }
510
511 /**
512 * Return the vector of values for the change in the standard state
513 * enthalpies of reaction. These values don't depend upon the concentration
514 * of the solution.
515 *
516 * units = J kmol-1
517 *
518 * @param deltaH Output vector of ss deltaH's for reactions Length:
519 * nReactions().
520 */
521 virtual void getDeltaSSEnthalpy(double* deltaH) {
522 throw NotImplementedError("Kinetics::getDeltaSSEnthalpy");
523 }
524
525 /**
526 * Return the vector of values for the change in the standard state
527 * entropies for each reaction. These values don't depend upon the
528 * concentration of the solution.
529 *
530 * units = J kmol-1 Kelvin-1
531 *
532 * @param deltaS Output vector of ss deltaS's for reactions Length:
533 * nReactions().
534 */
535 virtual void getDeltaSSEntropy(double* deltaS) {
536 throw NotImplementedError("Kinetics::getDeltaSSEntropy");
537 }
538
539 /**
540 * Return a vector of values of effective concentrations of third-body
541 * collision partners of any reaction. Entries for reactions not involving
542 * third-body collision partners are not defined and set to not-a-number.
543 *
544 * @param concm Output vector of effective third-body concentrations.
545 * Length: nReactions().
546 */
547 virtual void getThirdBodyConcentrations(double* concm) {
548 throw NotImplementedError("Kinetics::getThirdBodyConcentrations",
549 "Not applicable/implemented for Kinetics object of type '{}'",
550 kineticsType());
551 }
552
553 /**
554 * Provide direct access to current third-body concentration values.
555 * @see getThirdBodyConcentrations.
556 */
557 virtual const vector<double>& thirdBodyConcentrations() const {
558 throw NotImplementedError("Kinetics::thirdBodyConcentrations",
559 "Not applicable/implemented for Kinetics object of type '{}'",
560 kineticsType());
561 }
562
563 //! @}
564 //! @name Species Production Rates
565 //! @{
566
567 /**
568 * Species creation rates [kmol/m^3/s or kmol/m^2/s]. Return the species
569 * creation rates in array cdot, which must be dimensioned at least as
570 * large as the total number of species in all phases. @see nTotalSpecies.
571 *
572 * @param cdot Output vector of creation rates. Length: #m_kk.
573 */
574 virtual void getCreationRates(double* cdot);
575
576 /**
577 * Species destruction rates [kmol/m^3/s or kmol/m^2/s]. Return the species
578 * destruction rates in array ddot, which must be dimensioned at least as
579 * large as the total number of species. @see nTotalSpecies.
580 *
581 * @param ddot Output vector of destruction rates. Length: #m_kk.
582 */
583 virtual void getDestructionRates(double* ddot);
584
585 /**
586 * Species net production rates [kmol/m^3/s or kmol/m^2/s]. Return the
587 * species net production rates (creation - destruction) in array wdot,
588 * which must be dimensioned at least as large as the total number of
589 * species. @see nTotalSpecies.
590 *
591 * @param wdot Output vector of net production rates. Length: #m_kk.
592 */
593 virtual void getNetProductionRates(double* wdot);
594
595 //! @}
596
597 //! @addtogroup derivGroup
598 //! @{
599
600 /**
601 * @anchor kinDerivs
602 * @par Routines to Calculate Kinetics Derivatives (Jacobians)
603 * @name
604 *
605 * Kinetics derivatives are calculated with respect to temperature, pressure,
606 * molar concentrations and species mole fractions for forward/reverse/net rates
607 * of progress as well as creation/destruction and net production of species.
608 *
609 * The following suffixes are used to indicate derivatives:
610 * - `_ddT`: derivative with respect to temperature (a vector)
611 * - `_ddP`: derivative with respect to pressure (a vector)
612 * - `_ddC`: derivative with respect to molar concentration (a vector)
613 * - `_ddX`: derivative with respect to species mole fractions (a matrix)
614 * - `_ddCi`: derivative with respect to species concentrations (a matrix)
615 *
616 * @since New in Cantera 2.6
617 *
618 * @warning The calculation of kinetics derivatives is an experimental part of the
619 * %Cantera API and may be changed or removed without notice.
620 *
621 * Source term derivatives are based on a generic rate-of-progress expression
622 * for the @f$ i @f$-th reaction @f$ R_i @f$, which is a function of temperature
623 * @f$ T @f$, pressure @f$ P @f$ and molar concentrations @f$ C_j @f$:
624 * @f[
625 * R_i = k_{f,i} C_M^{\nu_{M,i}} \prod_j C_j^{\nu_{ji}^\prime} -
626 * k_{r,i} C_M^{\nu_{M,i}} \prod_j C_j^{\nu_{ji}^{\prime\prime}}
627 * @f]
628 * Forward/reverse rate expressions @f$ k_{f,i} @f$ and @f$ k_{r,i} @f$ are
629 * implemented by ReactionRate specializations; forward/reverse stoichiometric
630 * coefficients are @f$ \nu_{ji}^\prime @f$ and @f$ \nu_{ji}^{\prime\prime} @f$.
631 * Unless the reaction involves third-body colliders, @f$ \nu_{M,i} = 0 @f$.
632 * For three-body reactions, effective ThirdBody collider concentrations @f$ C_M @f$
633 * are considered with @f$ \nu_{M,i} = 1 @f$. For more detailed information on
634 * relevant theory, see, for example, Perini, et al. @cite perini2012 or Niemeyer,
635 * et al. @cite niemeyer2017, although specifics of %Cantera's implementation may
636 * differ.
637 *
638 * Partial derivatives are obtained from the product rule, where resulting terms
639 * consider reaction rate derivatives, derivatives of the concentration product
640 * term, and, if applicable, third-body term derivatives. ReactionRate
641 * specializations may implement exact derivatives (example:
642 * ArrheniusRate::ddTScaledFromStruct) or approximate them numerically (examples:
643 * ReactionData::perturbTemperature, PlogData::perturbPressure,
644 * FalloffData::perturbThirdBodies). Derivatives of concentration and third-body
645 * terms are based on analytic expressions.
646 *
647 * %Species creation and destruction rates are obtained by multiplying
648 * rate-of-progress vectors by stoichiometric coefficient matrices. As this is a
649 * linear operation, it is possible to calculate derivatives the same way.
650 *
651 * All derivatives are calculated for source terms while holding other properties
652 * constant, independent of whether equation of state or @f$ \sum X_k = 1 @f$
653 * constraints are satisfied. Thus, derivatives deviate from Jacobians and
654 * numerical derivatives that implicitly enforce these constraints. Depending
655 * on application and equation of state, derivatives can nevertheless be used to
656 * obtain Jacobians, for example:
657 *
658 * - The Jacobian of net production rates @f$ \dot{\omega}_{k,\mathrm{net}} @f$
659 * with respect to temperature at constant pressure needs to consider changes
660 * of molar density @f$ C @f$ due to temperature
661 * @f[
662 * \left.
663 * \frac{\partial \dot{\omega}_{k,\mathrm{net}}}{\partial T}
664 * \right|_{P=\mathrm{const}} =
665 * \frac{\partial \dot{\omega}_{k,\mathrm{net}}}{\partial T} +
666 * \frac{\partial \dot{\omega}_{k,\mathrm{net}}}{\partial C}
667 * \left. \frac{\partial C}{\partial T} \right|_{P=\mathrm{const}}
668 * @f]
669 * where for an ideal gas @f$ \partial C / \partial T = - C / T @f$. The
670 * remaining partial derivatives are obtained from getNetProductionRates_ddT()
671 * and getNetProductionRates_ddC(), respectively.
672 *
673 * - The Jacobian of @f$ \dot{\omega}_{k,\mathrm{net}} @f$ with respect to
674 * temperature at constant volume needs to consider pressure changes due to
675 * temperature
676 * @f[
677 * \left.
678 * \frac{\partial \dot{\omega}_{k,\mathrm{net}}}{\partial T}
679 * \right|_{V=\mathrm{const}} =
680 * \frac{\partial \dot{\omega}_{k,\mathrm{net}}}{\partial T} +
681 * \frac{\partial \dot{\omega}_{k,\mathrm{net}}}{\partial P}
682 * \left. \frac{\partial P}{\partial T} \right|_{V=\mathrm{const}}
683 * @f]
684 * where for an ideal gas @f$ \partial P / \partial T = P / T @f$. The
685 * remaining partial derivatives are obtained from getNetProductionRates_ddT()
686 * and getNetProductionRates_ddP(), respectively.
687 *
688 * - Similar expressions can be derived for other derivatives and source terms.
689 *
690 * While some applications require exact derivatives, others can tolerate
691 * approximate derivatives that neglect terms to increase computational speed
692 * and/or improve Jacobian sparsity (example: AdaptivePreconditioner).
693 * Derivative evaluations settings are accessible by keyword/value pairs
694 * using the methods getDerivativeSettings() and setDerivativeSettings().
695 *
696 * For BulkKinetics, the following keyword/value pairs are supported:
697 * - `skip-third-bodies` (boolean): if `false` (default), third body
698 * concentrations are considered for the evaluation of Jacobians
699 * - `skip-falloff` (boolean): if `false` (default), third-body effects
700 * on rate constants are considered for the evaluation of derivatives.
701 * - `rtol-delta` (double): relative tolerance used to perturb properties
702 * when calculating numerical derivatives. The default value is 1e-8.
703 *
704 * For InterfaceKinetics, the following keyword/value pairs are supported:
705 * - `skip-coverage-dependence` (boolean): if `false` (default), rate constant
706 * coverage dependence is not considered when evaluating derivatives.
707 * - `skip-electrochemistry` (boolean): if `false` (default), electrical charge
708 * is not considered in evaluating the derivatives and these reactions are
709 * treated as normal surface reactions.
710 * - `rtol-delta` (double): relative tolerance used to perturb properties
711 * when calculating numerical derivatives. The default value is 1e-8.
712 *
713 * @{
714 */
715
716 /**
717 * Retrieve derivative settings.
718 *
719 * @param settings AnyMap containing settings determining derivative evaluation.
720 */
721 virtual void getDerivativeSettings(AnyMap& settings) const
722 {
723 throw NotImplementedError("Kinetics::getDerivativeSettings",
724 "Not implemented for kinetics type '{}'.", kineticsType());
725 }
726
727 /**
728 * Set/modify derivative settings.
729 *
730 * @param settings AnyMap containing settings determining derivative evaluation.
731 */
732 virtual void setDerivativeSettings(const AnyMap& settings)
733 {
734 throw NotImplementedError("Kinetics::setDerivativeSettings",
735 "Not implemented for kinetics type '{}'.", kineticsType());
736 }
737
738 /**
739 * Calculate derivatives for forward rate constants with respect to temperature
740 * at constant pressure, molar concentration and mole fractions
741 * @param[out] dkfwd Output vector of derivatives. Length: nReactions().
742 */
743 virtual void getFwdRateConstants_ddT(double* dkfwd)
744 {
745 throw NotImplementedError("Kinetics::getFwdRateConstants_ddT",
746 "Not implemented for kinetics type '{}'.", kineticsType());
747 }
748
749 /**
750 * Calculate derivatives for forward rate constants with respect to pressure
751 * at constant temperature, molar concentration and mole fractions.
752 * @param[out] dkfwd Output vector of derivatives. Length: nReactions().
753 */
754 virtual void getFwdRateConstants_ddP(double* dkfwd)
755 {
756 throw NotImplementedError("Kinetics::getFwdRateConstants_ddP",
757 "Not implemented for kinetics type '{}'.", kineticsType());
758 }
759
760 /**
761 * Calculate derivatives for forward rate constants with respect to molar
762 * concentration at constant temperature, pressure and mole fractions.
763 * @param[out] dkfwd Output vector of derivatives. Length: nReactions().
764 *
765 * @warning This method is an experimental part of the %Cantera API and
766 * may be changed or removed without notice.
767 */
768 virtual void getFwdRateConstants_ddC(double* dkfwd)
769 {
770 throw NotImplementedError("Kinetics::getFwdRateConstants_ddC",
771 "Not implemented for kinetics type '{}'.", kineticsType());
772 }
773
774 /**
775 * Calculate derivatives for forward rates-of-progress with respect to temperature
776 * at constant pressure, molar concentration and mole fractions.
777 * @param[out] drop Output vector of derivatives. Length: nReactions().
778 */
779 virtual void getFwdRatesOfProgress_ddT(double* drop)
780 {
781 throw NotImplementedError("Kinetics::getFwdRatesOfProgress_ddT",
782 "Not implemented for kinetics type '{}'.", kineticsType());
783 }
784
785 /**
786 * Calculate derivatives for forward rates-of-progress with respect to pressure
787 * at constant temperature, molar concentration and mole fractions.
788 * @param[out] drop Output vector of derivatives. Length: nReactions().
789 */
790 virtual void getFwdRatesOfProgress_ddP(double* drop)
791 {
792 throw NotImplementedError("Kinetics::getFwdRatesOfProgress_ddP",
793 "Not implemented for kinetics type '{}'.", kineticsType());
794 }
795
796 /**
797 * Calculate derivatives for forward rates-of-progress with respect to molar
798 * concentration at constant temperature, pressure and mole fractions.
799 * @param[out] drop Output vector of derivatives. Length: nReactions().
800 *
801 * @warning This method is an experimental part of the %Cantera API and
802 * may be changed or removed without notice.
803 */
804 virtual void getFwdRatesOfProgress_ddC(double* drop)
805 {
806 throw NotImplementedError("Kinetics::getFwdRatesOfProgress_ddC",
807 "Not implemented for kinetics type '{}'.", kineticsType());
808 }
809
810 /**
811 * Calculate derivatives for forward rates-of-progress with respect to species
812 * mole fractions at constant temperature, pressure and molar concentration.
813 *
814 * The method returns a matrix with nReactions() rows and nTotalSpecies() columns.
815 * For a derivative with respect to @f$ X_i @f$, all other @f$ X_j @f$ are held
816 * constant, rather than enforcing @f$ \sum X_j = 1 @f$.
817 *
818 * @warning This method is an experimental part of the %Cantera API and
819 * may be changed or removed without notice.
820 */
821 virtual Eigen::SparseMatrix<double> fwdRatesOfProgress_ddX()
822 {
823 throw NotImplementedError("Kinetics::fwdRatesOfProgress_ddX",
824 "Not implemented for kinetics type '{}'.", kineticsType());
825 }
826
827 /**
828 * Calculate derivatives for forward rates-of-progress with respect to species
829 * concentration at constant temperature, pressure and remaining species
830 * concentrations.
831 *
832 * The method returns a matrix with nReactions() rows and nTotalSpecies() columns.
833 * For a derivative with respect to @f$ c_i @f$, all other @f$ c_j @f$ are held
834 * constant.
835 *
836 * @warning This method is an experimental part of the %Cantera API and
837 * may be changed or removed without notice.
838 *
839 * @since New in %Cantera 3.0.
840 */
841 virtual Eigen::SparseMatrix<double> fwdRatesOfProgress_ddCi()
842 {
843 throw NotImplementedError("Kinetics::fwdRatesOfProgress_ddCi",
844 "Not implemented for kinetics type '{}'.", kineticsType());
845 }
846
847 /**
848 * Calculate derivatives for reverse rates-of-progress with respect to temperature
849 * at constant pressure, molar concentration and mole fractions.
850 * @param[out] drop Output vector of derivatives. Length: nReactions().
851 */
852 virtual void getRevRatesOfProgress_ddT(double* drop)
853 {
854 throw NotImplementedError("Kinetics::getRevRatesOfProgress_ddT",
855 "Not implemented for kinetics type '{}'.", kineticsType());
856 }
857
858 /**
859 * Calculate derivatives for reverse rates-of-progress with respect to pressure
860 * at constant temperature, molar concentration and mole fractions.
861 * @param[out] drop Output vector of derivatives. Length: nReactions().
862 */
863 virtual void getRevRatesOfProgress_ddP(double* drop)
864 {
865 throw NotImplementedError("Kinetics::getRevRatesOfProgress_ddP",
866 "Not implemented for kinetics type '{}'.", kineticsType());
867 }
868
869 /**
870 * Calculate derivatives for reverse rates-of-progress with respect to molar
871 * concentration at constant temperature, pressure and mole fractions.
872 * @param[out] drop Output vector of derivatives. Length: nReactions().
873 *
874 * @warning This method is an experimental part of the %Cantera API and
875 * may be changed or removed without notice.
876 */
877 virtual void getRevRatesOfProgress_ddC(double* drop)
878 {
879 throw NotImplementedError("Kinetics::getRevRatesOfProgress_ddC",
880 "Not implemented for kinetics type '{}'.", kineticsType());
881 }
882
883 /**
884 * Calculate derivatives for reverse rates-of-progress with respect to species
885 * mole fractions at constant temperature, pressure and molar concentration.
886 *
887 * The method returns a matrix with nReactions() rows and nTotalSpecies() columns.
888 * For a derivative with respect to @f$ X_i @f$, all other @f$ X_j @f$ are held
889 * constant, rather than enforcing @f$ \sum X_j = 1 @f$.
890 *
891 * @warning This method is an experimental part of the %Cantera API and
892 * may be changed or removed without notice.
893 */
894 virtual Eigen::SparseMatrix<double> revRatesOfProgress_ddX()
895 {
896 throw NotImplementedError("Kinetics::revRatesOfProgress_ddX",
897 "Not implemented for kinetics type '{}'.", kineticsType());
898 }
899
900 /**
901 * Calculate derivatives for forward rates-of-progress with respect to species
902 * concentration at constant temperature, pressure and remaining species
903 * concentrations.
904 *
905 * The method returns a matrix with nReactions() rows and nTotalSpecies() columns.
906 * For a derivative with respect to @f$ c_i @f$, all other @f$ c_j @f$ are held
907 * constant.
908 *
909 * @warning This method is an experimental part of the %Cantera API and
910 * may be changed or removed without notice.
911 *
912 * @since New in %Cantera 3.0.
913 */
914 virtual Eigen::SparseMatrix<double> revRatesOfProgress_ddCi()
915 {
916 throw NotImplementedError("Kinetics::revRatesOfProgress_ddCi",
917 "Not implemented for kinetics type '{}'.", kineticsType());
918 }
919
920 /**
921 * Calculate derivatives for net rates-of-progress with respect to temperature
922 * at constant pressure, molar concentration and mole fractions.
923 * @param[out] drop Output vector of derivatives. Length: nReactions().
924 */
925 virtual void getNetRatesOfProgress_ddT(double* drop)
926 {
927 throw NotImplementedError("Kinetics::getNetRatesOfProgress_ddT",
928 "Not implemented for kinetics type '{}'.", kineticsType());
929 }
930
931 /**
932 * Calculate derivatives for net rates-of-progress with respect to pressure
933 * at constant temperature, molar concentration and mole fractions.
934 * @param[out] drop Output vector of derivatives. Length: nReactions().
935 */
936 virtual void getNetRatesOfProgress_ddP(double* drop)
937 {
938 throw NotImplementedError("Kinetics::getNetRatesOfProgress_ddP",
939 "Not implemented for kinetics type '{}'.", kineticsType());
940 }
941
942 /**
943 * Calculate derivatives for net rates-of-progress with respect to molar
944 * concentration at constant temperature, pressure and mole fractions.
945 * @param[out] drop Output vector of derivatives. Length: nReactions().
946 *
947 * @warning This method is an experimental part of the %Cantera API and
948 * may be changed or removed without notice.
949 */
950 virtual void getNetRatesOfProgress_ddC(double* drop)
951 {
952 throw NotImplementedError("Kinetics::getNetRatesOfProgress_ddC",
953 "Not implemented for kinetics type '{}'.", kineticsType());
954 }
955
956 /**
957 * Calculate derivatives for net rates-of-progress with respect to species
958 * mole fractions at constant temperature, pressure and molar concentration.
959 *
960 * The method returns a matrix with nReactions() rows and nTotalSpecies() columns.
961 * For a derivative with respect to @f$ X_i @f$, all other @f$ X_j @f$ are held
962 * constant, rather than enforcing @f$ \sum X_j = 1 @f$.
963 *
964 * @warning This method is an experimental part of the %Cantera API and
965 * may be changed or removed without notice.
966 */
967 virtual Eigen::SparseMatrix<double> netRatesOfProgress_ddX()
968 {
969 throw NotImplementedError("Kinetics::netRatesOfProgress_ddX",
970 "Not implemented for kinetics type '{}'.", kineticsType());
971 }
972
973 /**
974 * Calculate derivatives for net rates-of-progress with respect to species
975 * concentration at constant temperature, pressure, and remaining species
976 * concentrations.
977 *
978 * The method returns a matrix with nReactions() rows and nTotalSpecies() columns.
979 * For a derivative with respect to @f$ c_i @f$, all other @f$ c_j @f$ are held
980 * constant.
981 *
982 * @warning This method is an experimental part of the %Cantera API and
983 * may be changed or removed without notice.
984 *
985 * @since New in %Cantera 3.0.
986 */
987 virtual Eigen::SparseMatrix<double> netRatesOfProgress_ddCi()
988 {
989 throw NotImplementedError("Kinetics::netRatesOfProgress_ddCi",
990 "Not implemented for kinetics type '{}'.", kineticsType());
991 }
992
993 /**
994 * Calculate derivatives for species creation rates with respect to temperature
995 * at constant pressure, molar concentration and mole fractions.
996 * @param[out] dwdot Output vector of derivatives. Length: #m_kk.
997 */
998 void getCreationRates_ddT(double* dwdot);
999
1000 /**
1001 * Calculate derivatives for species creation rates with respect to pressure
1002 * at constant temperature, molar concentration and mole fractions.
1003 * @param[out] dwdot Output vector of derivatives. Length: #m_kk.
1004 */
1005 void getCreationRates_ddP(double* dwdot);
1006
1007 /**
1008 * Calculate derivatives for species creation rates with respect to molar
1009 * concentration at constant temperature, pressure and mole fractions.
1010 * @param[out] dwdot Output vector of derivatives. Length: #m_kk.
1011 *
1012 * @warning This method is an experimental part of the %Cantera API and
1013 * may be changed or removed without notice.
1014 */
1015 void getCreationRates_ddC(double* dwdot);
1016
1017 /**
1018 * Calculate derivatives for species creation rates with respect to species
1019 * mole fractions at constant temperature, pressure and molar concentration.
1020 *
1021 * The method returns a square matrix with nTotalSpecies() rows and columns.
1022 * For a derivative with respect to @f$ X_i @f$, all other @f$ X_j @f$ are held
1023 * constant, rather than enforcing @f$ \sum X_j = 1 @f$.
1024 *
1025 * @warning This method is an experimental part of the %Cantera API and
1026 * may be changed or removed without notice.
1027 */
1028 Eigen::SparseMatrix<double> creationRates_ddX();
1029
1030 /**
1031 * Calculate derivatives for species creation rates with respect to species
1032 * concentration at constant temperature, pressure, and concentration of all other
1033 * species.
1034 *
1035 * The method returns a square matrix with nTotalSpecies() rows and columns.
1036 * For a derivative with respect to @f$ c_i @f$, all other @f$ c_j @f$ are held
1037 * constant.
1038 *
1039 * @warning This method is an experimental part of the %Cantera API and
1040 * may be changed or removed without notice.
1041 *
1042 * @since New in %Cantera 3.0.
1043 */
1044 Eigen::SparseMatrix<double> creationRates_ddCi();
1045
1046 /**
1047 * Calculate derivatives for species destruction rates with respect to temperature
1048 * at constant pressure, molar concentration and mole fractions.
1049 * @param[out] dwdot Output vector of derivatives. Length: #m_kk.
1050 */
1051 void getDestructionRates_ddT(double* dwdot);
1052
1053 /**
1054 * Calculate derivatives for species destruction rates with respect to pressure
1055 * at constant temperature, molar concentration and mole fractions.
1056 * @param[out] dwdot Output vector of derivatives. Length: #m_kk.
1057 */
1058 void getDestructionRates_ddP(double* dwdot);
1059
1060 /**
1061 * Calculate derivatives for species destruction rates with respect to molar
1062 * concentration at constant temperature, pressure and mole fractions.
1063 * @param[out] dwdot Output vector of derivatives. Length: #m_kk.
1064 *
1065 * @warning This method is an experimental part of the %Cantera API and
1066 * may be changed or removed without notice.
1067 */
1068 void getDestructionRates_ddC(double* dwdot);
1069
1070 /**
1071 * Calculate derivatives for species destruction rates with respect to species
1072 * mole fractions at constant temperature, pressure and molar concentration.
1073 *
1074 * The method returns a square matrix with nTotalSpecies() rows and columns.
1075 * For a derivative with respect to @f$ X_i @f$, all other @f$ X_j @f$ are held
1076 * constant, rather than enforcing @f$ \sum X_j = 1 @f$.
1077 *
1078 * @warning This method is an experimental part of the %Cantera API and
1079 * may be changed or removed without notice.
1080 */
1081 Eigen::SparseMatrix<double> destructionRates_ddX();
1082
1083 /**
1084 * Calculate derivatives for species destruction rates with respect to species
1085 * concentration at constant temperature, pressure, and concentration of all other
1086 * species.
1087 *
1088 * The method returns a square matrix with nTotalSpecies() rows and columns.
1089 * For a derivative with respect to @f$ c_i @f$, all other @f$ c_j @f$ are held
1090 * constant.
1091 *
1092 * @warning This method is an experimental part of the %Cantera API and
1093 * may be changed or removed without notice.
1094 *
1095 * @since New in %Cantera 3.0.
1096 */
1097 Eigen::SparseMatrix<double> destructionRates_ddCi();
1098
1099 /**
1100 * Calculate derivatives for species net production rates with respect to
1101 * temperature at constant pressure, molar concentration and mole fractions.
1102 * @param[out] dwdot Output vector of derivatives. Length: #m_kk.
1103 */
1104 void getNetProductionRates_ddT(double* dwdot);
1105
1106 /**
1107 * Calculate derivatives for species net production rates with respect to pressure
1108 * at constant temperature, molar concentration and mole fractions.
1109 * @param[out] dwdot Output vector of derivatives. Length: #m_kk.
1110 */
1111 void getNetProductionRates_ddP(double* dwdot);
1112
1113 /**
1114 * Calculate derivatives for species net production rates with respect to molar
1115 * concentration at constant temperature, pressure and mole fractions.
1116 * @param[out] dwdot Output vector of derivatives. Length: #m_kk.
1117 *
1118 * @warning This method is an experimental part of the %Cantera API and
1119 * may be changed or removed without notice.
1120 */
1121 void getNetProductionRates_ddC(double* dwdot);
1122
1123 /**
1124 * Calculate derivatives for species net production rates with respect to species
1125 * mole fractions at constant temperature, pressure and molar concentration.
1126 *
1127 * The method returns a square matrix with nTotalSpecies() rows and columns.
1128 * For a derivative with respect to @f$ X_i @f$, all other @f$ X_j @f$ are held
1129 * constant, rather than enforcing @f$ \sum X_j = 1 @f$.
1130 *
1131 * @warning This method is an experimental part of the %Cantera API and
1132 * may be changed or removed without notice.
1133 */
1134 Eigen::SparseMatrix<double> netProductionRates_ddX();
1135
1136 /**
1137 * Calculate derivatives for species net production rates with respect to species
1138 * concentration at constant temperature, pressure, and concentration of all other
1139 * species.
1140 *
1141 * The method returns a square matrix with nTotalSpecies() rows and columns.
1142 * For a derivative with respect to @f$ c_i @f$, all other @f$ c_j @f$ are held
1143 * constant.
1144 *
1145 * @warning This method is an experimental part of the %Cantera API and
1146 * may be changed or removed without notice.
1147 *
1148 * @since New in %Cantera 3.0.
1149 */
1150 Eigen::SparseMatrix<double> netProductionRates_ddCi();
1151
1152 /** @} End of Kinetics Derivatives */
1153 //! @} End of addtogroup derivGroup
1154
1155 //! @name Reaction Mechanism Informational Query Routines
1156 //! @{
1157
1158 /**
1159 * Stoichiometric coefficient of species k as a reactant in reaction i.
1160 *
1161 * @param k kinetic species index
1162 * @param i reaction index
1163 */
1164 virtual double reactantStoichCoeff(size_t k, size_t i) const;
1165
1166 /**
1167 * Stoichiometric coefficient matrix for reactants.
1168 */
1169 Eigen::SparseMatrix<double> reactantStoichCoeffs() const {
1171 }
1172
1173 /**
1174 * Stoichiometric coefficient of species k as a product in reaction i.
1175 *
1176 * @param k kinetic species index
1177 * @param i reaction index
1178 */
1179 virtual double productStoichCoeff(size_t k, size_t i) const;
1180
1181 /**
1182 * Stoichiometric coefficient matrix for products.
1183 */
1184 Eigen::SparseMatrix<double> productStoichCoeffs() const {
1186 }
1187
1188 /**
1189 * Stoichiometric coefficient matrix for products of reversible reactions.
1190 */
1191 Eigen::SparseMatrix<double> revProductStoichCoeffs() const {
1193 }
1194
1195 //! Reactant order of species k in reaction i.
1196 /*!
1197 * This is the nominal order of the activity concentration in
1198 * determining the forward rate of progress of the reaction
1199 *
1200 * @param k kinetic species index
1201 * @param i reaction index
1202 */
1203 virtual double reactantOrder(size_t k, size_t i) const {
1204 throw NotImplementedError("Kinetics::reactantOrder");
1205 }
1206
1207 //! product Order of species k in reaction i.
1208 /*!
1209 * This is the nominal order of the activity concentration of species k in
1210 * determining the reverse rate of progress of the reaction i
1211 *
1212 * For irreversible reactions, this will all be zero.
1213 *
1214 * @param k kinetic species index
1215 * @param i reaction index
1216 */
1217 virtual double productOrder(int k, int i) const {
1218 throw NotImplementedError("Kinetics::productOrder");
1219 }
1220
1221 //! Get the vector of activity concentrations used in the kinetics object
1222 /*!
1223 * @param[out] conc Vector of activity concentrations. Length is equal
1224 * to the number of species in the kinetics object
1225 */
1226 virtual void getActivityConcentrations(double* const conc) {
1227 throw NotImplementedError("Kinetics::getActivityConcentrations");
1228 }
1229
1230 /**
1231 * String specifying the type of reaction.
1232 *
1233 * @param i reaction index
1234 * @since Method returned magic number prior to %Cantera 3.0.
1235 * @deprecated To be removed after %Cantera 3.0. Replace with
1236 * `kin->reaction(i)->type()`
1237 */
1238 virtual string reactionType(size_t i) const;
1239
1240 /**
1241 * String specifying the type of reaction.
1242 *
1243 * @param i reaction index
1244 * @deprecated To be removed after %Cantera 3.0.
1245 */
1246 virtual string reactionTypeStr(size_t i) const;
1247
1248 /**
1249 * True if reaction i has been declared to be reversible. If isReversible(i)
1250 * is false, then the reverse rate of progress for reaction i is always
1251 * zero.
1252 *
1253 * @param i reaction index
1254 */
1255 virtual bool isReversible(size_t i) {
1256 throw NotImplementedError("Kinetics::isReversible");
1257 }
1258
1259 /**
1260 * Return a string representing the reaction.
1261 *
1262 * @param i reaction index
1263 * @deprecated To be removed after %Cantera 3.0. Replace with
1264 * `kin->reaction(i)->equation()`.
1265 */
1266 string reactionString(size_t i) const;
1267
1268 //! Returns a string containing the reactants side of the reaction equation.
1269 //! @deprecated To be removed after %Cantera 3.0. Replace with
1270 //! `kin->reaction(i)->reactantString()`
1271 string reactantString(size_t i) const;
1272
1273 //! Returns a string containing the products side of the reaction equation.
1274 //! @deprecated To be removed after %Cantera 3.0. Replace with
1275 //! `kin->reaction(i)->productString()`
1276 string productString(size_t i) const;
1277
1278 /**
1279 * Return the forward rate constants
1280 *
1281 * The computed values include all temperature-dependent and pressure-dependent
1282 * contributions. By default, third-body concentrations are only considered if
1283 * they are part of the reaction rate definition; for a legacy implementation that
1284 * includes third-body concentrations see Cantera::use_legacy_rate_constants().
1285 * Length is the number of reactions. Units are a combination of kmol, m^3 and s,
1286 * that depend on the rate expression for the reaction.
1287 *
1288 * @param kfwd Output vector containing the forward reaction rate
1289 * constants. Length: nReactions().
1290 */
1291 virtual void getFwdRateConstants(double* kfwd) {
1292 throw NotImplementedError("Kinetics::getFwdRateConstants");
1293 }
1294
1295 /**
1296 * Return the reverse rate constants.
1297 *
1298 * The computed values include all temperature-dependent and pressure-dependent
1299 * contributions. By default, third-body concentrations are only considered if
1300 * they are part of the reaction rate definition; for a legacy implementation that
1301 * includes third-body concentrations see Cantera::use_legacy_rate_constants().
1302 * Length is the number of reactions. Units are a combination of kmol, m^3 and s,
1303 * that depend on the rate expression for the reaction. Note, this routine will
1304 * return rate constants for irreversible reactions if the default for
1305 * `doIrreversible` is overridden.
1306 *
1307 * @param krev Output vector of reverse rate constants
1308 * @param doIrreversible boolean indicating whether irreversible reactions
1309 * should be included.
1310 */
1311 virtual void getRevRateConstants(double* krev,
1312 bool doIrreversible = false) {
1313 throw NotImplementedError("Kinetics::getRevRateConstants");
1314 }
1315
1316 //! @}
1317 //! @name Reaction Mechanism Construction
1318 //! @{
1319
1320 //! Add a phase to the kinetics manager object.
1321 /*!
1322 * This must be done before the function init() is called or before any
1323 * reactions are input. The following fields are updated:
1324 *
1325 * - #m_start -> vector of integers, containing the starting position of
1326 * the species for each phase in the kinetics mechanism.
1327 * - #m_rxnphase -> index of the phase where reactions occur, which is the lowest-
1328 * dimensional phase in the system, for example the surface in a surface
1329 * mechanism.
1330 * - #m_thermo -> vector of pointers to ThermoPhase phases that
1331 * participate in the kinetics mechanism.
1332 * - #m_phaseindex -> map containing the string id of each
1333 * ThermoPhase phase as a key and the index of the phase within the
1334 * kinetics manager object as the value.
1335 *
1336 * @param thermo Reference to the ThermoPhase to be added.
1337 * @since New in %Cantera 3.0. Replaces addPhase.
1338 */
1339 virtual void addThermo(shared_ptr<ThermoPhase> thermo);
1340
1341 //! @see Kinetics::addThermo(shared_ptr<ThermoPhase>)
1342 //! @deprecated To be removed after %Cantera 3.0. Replaced by addThermo
1343 //! pointer.
1344 virtual void addPhase(ThermoPhase& thermo);
1345
1346 /**
1347 * Prepare the class for the addition of reactions, after all phases have
1348 * been added. This method is called automatically when the first reaction
1349 * is added. It needs to be called directly only in the degenerate case
1350 * where there are no reactions. The base class method does nothing, but
1351 * derived classes may use this to perform any initialization (allocating
1352 * arrays, etc.) that requires knowing the phases.
1353 */
1354 virtual void init() {}
1355
1356 //! Return the parameters for a phase definition which are needed to
1357 //! reconstruct an identical object using the newKinetics function. This
1358 //! excludes the reaction definitions, which are handled separately.
1360
1361 /**
1362 * Resize arrays with sizes that depend on the total number of species.
1363 * Automatically called before adding each Reaction and Phase.
1364 */
1365 virtual void resizeSpecies();
1366
1367 /**
1368 * Add a single reaction to the mechanism. Derived classes should call the
1369 * base class method in addition to handling their own specialized behavior.
1370 *
1371 * @param r Pointer to the Reaction object to be added.
1372 * @param resize If `true`, resizeReactions is called after reaction is added.
1373 * @return `true` if the reaction is added or `false` if it was skipped
1374 */
1375 virtual bool addReaction(shared_ptr<Reaction> r, bool resize=true);
1376
1377 /**
1378 * Modify the rate expression associated with a reaction. The
1379 * stoichiometric equation, type of the reaction, reaction orders, third
1380 * body efficiencies, reversibility, etc. must be unchanged.
1381 *
1382 * @param i Index of the reaction to be modified
1383 * @param rNew Reaction with the new rate expressions
1384 */
1385 virtual void modifyReaction(size_t i, shared_ptr<Reaction> rNew);
1386
1387 /**
1388 * Return the Reaction object for reaction *i*. Changes to this object do
1389 * not affect the Kinetics object until the #modifyReaction function is
1390 * called.
1391 */
1392 shared_ptr<Reaction> reaction(size_t i);
1393
1394 shared_ptr<const Reaction> reaction(size_t i) const;
1395
1396 //! Determine behavior when adding a new reaction that contains species not
1397 //! defined in any of the phases associated with this kinetics manager. If
1398 //! set to true, the reaction will silently be ignored. If false, (the
1399 //! default) an exception will be raised.
1400 void skipUndeclaredSpecies(bool skip) {
1402 }
1403 bool skipUndeclaredSpecies() const {
1405 }
1406
1407 //! Determine behavior when adding a new reaction that contains third-body
1408 //! efficiencies for species not defined in any of the phases associated
1409 //! with this kinetics manager. If set to true, the given third-body
1410 //! efficiency will be ignored. If false, (the default) an exception will be
1411 //! raised.
1414 }
1415 bool skipUndeclaredThirdBodies() const {
1417 }
1418
1419 //! @}
1420 //! @name Altering Reaction Rates
1421 //!
1422 //! These methods alter reaction rates. They are designed primarily for
1423 //! carrying out sensitivity analysis, but may be used for any purpose
1424 //! requiring dynamic alteration of rate constants. For each reaction, a
1425 //! real-valued multiplier may be defined that multiplies the reaction rate
1426 //! coefficient. The multiplier may be set to zero to completely remove a
1427 //! reaction from the mechanism.
1428 //! @{
1429
1430 //! The current value of the multiplier for reaction i.
1431 /*!
1432 * @param i index of the reaction
1433 */
1434 double multiplier(size_t i) const {
1435 return m_perturb[i];
1436 }
1437
1438 //! Set the multiplier for reaction i to f.
1439 /*!
1440 * @param i index of the reaction
1441 * @param f value of the multiplier.
1442 */
1443 virtual void setMultiplier(size_t i, double f) {
1444 m_perturb[i] = f;
1445 }
1446
1447 virtual void invalidateCache() {
1448 m_cache.clear();
1449 };
1450
1451 //! @}
1452 //! Check for unmarked duplicate reactions and unmatched marked duplicates
1453 /**
1454 * If `throw_err` is true, then an exception will be thrown if either any
1455 * unmarked duplicate reactions are found, or if any marked duplicate
1456 * reactions do not have a matching duplicate reaction. If `throw_err` is
1457 * false, the indices of the first pair of duplicate reactions found will be
1458 * returned, or the index of the unmatched duplicate will be returned as
1459 * both elements of the pair. If no unmarked duplicates or unmatched marked
1460 * duplicate reactions are found, returns `(npos, npos)`.
1461 */
1462 virtual pair<size_t, size_t> checkDuplicates(bool throw_err=true) const;
1463
1464 /**
1465 * Takes as input an array of properties for all species in the mechanism
1466 * and copies those values belonging to a particular phase to the output
1467 * array.
1468 * @param data Input data array.
1469 * @param phase Pointer to one of the phase objects participating in this
1470 * reaction mechanism
1471 * @param phase_data Output array where the values for the the specified
1472 * phase are to be written.
1473 * @deprecated Unused. To be removed after %Cantera 3.0.
1474 */
1475 void selectPhase(const double* data, const ThermoPhase* phase,
1476 double* phase_data);
1477
1478 //! Set root Solution holding all phase information
1479 virtual void setRoot(shared_ptr<Solution> root) {
1480 m_root = root;
1481 }
1482
1483 //! Get the Solution object containing this Kinetics object and associated
1484 //! ThermoPhase objects
1485 shared_ptr<Solution> root() const {
1486 return m_root.lock();
1487 }
1488
1489 //! Calculate the reaction enthalpy of a reaction which
1490 //! has not necessarily been added into the Kinetics object
1491 //! @deprecated To be removed after %Cantera 3.0
1492 virtual double reactionEnthalpy(const Composition& reactants,
1493 const Composition& products);
1494
1495protected:
1496 //! Cache for saved calculations within each Kinetics object.
1498
1499 // Update internal rate-of-progress variables #m_ropf and #m_ropr.
1500 virtual void updateROP() {
1501 throw NotImplementedError("Kinetics::updateROP");
1502 }
1503
1504 //! Check whether `r1` and `r2` represent duplicate stoichiometries
1505 //! This function returns a ratio if two reactions are duplicates of
1506 //! one another, and 0.0 otherwise.
1507 /*!
1508 * `r1` and `r2` are maps of species key to stoichiometric coefficient, one
1509 * for each reaction, where the species key is `1+k` for reactants and
1510 * `-1-k` for products and `k` is the species index.
1511 *
1512 * @return 0.0 if the stoichiometries are not multiples of one another
1513 * Otherwise, it returns the ratio of the stoichiometric coefficients.
1514 */
1515 double checkDuplicateStoich(map<int, double>& r1, map<int, double>& r2) const;
1516
1517 //! @name Stoichiometry management
1518 //!
1519 //! These objects and functions handle turning reaction extents into species
1520 //! production rates and also handle turning thermo properties into reaction
1521 //! thermo properties.
1522 //! @{
1523
1524 //! Stoichiometry manager for the reactants for each reaction
1526
1527 //! Stoichiometry manager for the products for each reaction
1529
1530 //! Stoichiometry manager for the products of reversible reactions
1532
1533 //! Net stoichiometry (products - reactants)
1534 Eigen::SparseMatrix<double> m_stoichMatrix;
1535 //! @}
1536
1537 //! Boolean indicating whether Kinetics object is fully configured
1538 bool m_ready = false;
1539
1540 //! The number of species in all of the phases
1541 //! that participate in this kinetics mechanism.
1542 size_t m_kk = 0;
1543
1544 //! Vector of perturbation factors for each reaction's rate of
1545 //! progress vector. It is initialized to one.
1546 vector<double> m_perturb;
1547
1548 //! Vector of Reaction objects represented by this Kinetics manager
1549 vector<shared_ptr<Reaction>> m_reactions;
1550
1551 //! m_thermo is a vector of pointers to ThermoPhase objects that are
1552 //! involved with this kinetics operator
1553 /*!
1554 * For homogeneous kinetics applications, this vector will only have one
1555 * entry. For interfacial reactions, this vector will consist of multiple
1556 * entries; some of them will be surface phases, and the other ones will be
1557 * bulk phases. The order that the objects are listed determines the order
1558 * in which the species comprising each phase are listed in the source term
1559 * vector, originating from the reaction mechanism.
1560 *
1561 * Note that this kinetics object doesn't own these ThermoPhase objects
1562 * and is not responsible for creating or deleting them.
1563 */
1564 vector<ThermoPhase*> m_thermo;
1565
1566 //! vector of shared pointers, @see m_thermo
1567 //! @todo replace m_thermo with shared version after %Cantera 3.0
1568 vector<shared_ptr<ThermoPhase>> m_sharedThermo;
1569
1570 /**
1571 * m_start is a vector of integers specifying the beginning position for the
1572 * species vector for the n'th phase in the kinetics class.
1573 */
1574 vector<size_t> m_start;
1575
1576 /**
1577 * Mapping of the phase name to the position of the phase within the
1578 * kinetics object. Positions start with the value of 1. The member
1579 * function, phaseIndex() decrements by one before returning the index
1580 * value, so that missing phases return -1.
1581 */
1582 map<string, size_t> m_phaseindex;
1583
1584 //! Index in the list of phases of the one surface phase.
1585 //! @deprecated To be removed after %Cantera 3.0.
1587
1588 //! Phase Index where reactions are assumed to be taking place
1589 /*!
1590 * We calculate this by assuming that the phase with the lowest
1591 * dimensionality is the phase where reactions are taking place.
1592 */
1594
1595 //! number of spatial dimensions of lowest-dimensional phase.
1596 size_t m_mindim = 4;
1597
1598 //! Forward rate constant for each reaction
1599 vector<double> m_rfn;
1600
1601 //! Delta G^0 for all reactions
1602 vector<double> m_delta_gibbs0;
1603
1604 //! Reciprocal of the equilibrium constant in concentration units
1605 vector<double> m_rkcn;
1606
1607 //! Forward rate-of-progress for each reaction
1608 vector<double> m_ropf;
1609
1610 //! Reverse rate-of-progress for each reaction
1611 vector<double> m_ropr;
1612
1613 //! Net rate-of-progress for each reaction
1614 vector<double> m_ropnet;
1615
1616 //! The enthalpy change for each reaction to calculate Blowers-Masel rates
1617 vector<double> m_dH;
1618
1619 //! Buffer used for storage of intermediate reaction-specific results
1620 vector<double> m_rbuf;
1621
1622 //! See skipUndeclaredSpecies()
1624
1625 //! See skipUndeclaredThirdBodies()
1627
1628 //! Flag indicating whether reactions include undeclared third bodies
1630
1631 //! reference to Solution
1632 std::weak_ptr<Solution> m_root;
1633};
1634
1635}
1636
1637#endif
A map of string keys to values whose type can vary at runtime.
Definition AnyMap.h:427
Public interface for kinetics managers.
Definition Kinetics.h:126
size_t reactionPhaseIndex() const
Phase where the reactions occur.
Definition Kinetics.h:235
virtual void resizeReactions()
Finalize Kinetics object and associated objects.
Definition Kinetics.cpp:35
Eigen::SparseMatrix< double > reactantStoichCoeffs() const
Stoichiometric coefficient matrix for reactants.
Definition Kinetics.h:1169
double multiplier(size_t i) const
The current value of the multiplier for reaction i.
Definition Kinetics.h:1434
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:100
virtual void getFwdRatesOfProgress(double *fwdROP)
Return the forward rates of progress of the reactions.
Definition Kinetics.cpp:374
virtual void getEquilibriumConstants(double *kc)
Return a vector of Equilibrium constants.
Definition Kinetics.h:410
Kinetics(const Kinetics &)=delete
Kinetics objects are not copyable or assignable.
ThermoPhase & thermo(size_t n=0)
This method returns a reference to the nth ThermoPhase object defined in this kinetics mechanism.
Definition Kinetics.h:254
vector< shared_ptr< Reaction > > m_reactions
Vector of Reaction objects represented by this Kinetics manager.
Definition Kinetics.h:1549
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:93
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:200
ValueCache m_cache
Cache for saved calculations within each Kinetics object.
Definition Kinetics.h:1497
virtual bool isReversible(size_t i)
True if reaction i has been declared to be reversible.
Definition Kinetics.h:1255
vector< double > m_perturb
Vector of perturbation factors for each reaction's rate of progress vector.
Definition Kinetics.h:1546
vector< double > m_ropf
Forward rate-of-progress for each reaction.
Definition Kinetics.h:1608
bool m_ready
Boolean indicating whether Kinetics object is fully configured.
Definition Kinetics.h:1538
virtual const vector< double > & thirdBodyConcentrations() const
Provide direct access to current third-body concentration values.
Definition Kinetics.h:557
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:401
virtual void getNetRatesOfProgress(double *netROP)
Net rates of progress.
Definition Kinetics.cpp:386
size_t phaseIndex(const string &ph) const
Return the phase index of a phase in the list of phases defined within the object.
Definition Kinetics.h:207
virtual void getDeltaSSEntropy(double *deltaS)
Return the vector of values for the change in the standard state entropies for each reaction.
Definition Kinetics.h:535
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:1574
virtual string kineticsType() const
Identifies the Kinetics manager type.
Definition Kinetics.h:145
string reactionString(size_t i) const
Return a string representing the reaction.
Definition Kinetics.cpp:354
virtual void getDeltaElectrochemPotentials(double *deltaM)
Return the vector of values for the reaction electrochemical free energy change.
Definition Kinetics.h:467
vector< double > m_rkcn
Reciprocal of the equilibrium constant in concentration units.
Definition Kinetics.h:1605
void selectPhase(const double *data, const ThermoPhase *phase, double *phase_data)
Takes as input an array of properties for all species in the mechanism and copies those values belong...
Definition Kinetics.cpp:243
virtual string reactionTypeStr(size_t i) const
String specifying the type of reaction.
Definition Kinetics.cpp:349
vector< ThermoPhase * > m_thermo
m_thermo is a vector of pointers to ThermoPhase objects that are involved with this kinetics operator
Definition Kinetics.h:1564
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:1542
size_t m_rxnphase
Phase Index where reactions are assumed to be taking place.
Definition Kinetics.h:1593
virtual pair< size_t, size_t > checkDuplicates(bool throw_err=true) const
Check for unmarked duplicate reactions and unmatched marked duplicates.
Definition Kinetics.cpp:107
virtual void getReactionDelta(const double *property, double *deltaProperty) const
Change in species properties.
Definition Kinetics.cpp:392
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:1617
vector< double > m_ropr
Reverse rate-of-progress for each reaction.
Definition Kinetics.h:1611
size_t nPhases() const
The number of phases participating in the reaction mechanism.
Definition Kinetics.h:185
virtual void addPhase(ThermoPhase &thermo)
Definition Kinetics.cpp:616
virtual bool addReaction(shared_ptr< Reaction > r, bool resize=true)
Add a single reaction to the mechanism.
Definition Kinetics.cpp:665
virtual void getThirdBodyConcentrations(double *concm)
Return a vector of values of effective concentrations of third-body collision partners of any reactio...
Definition Kinetics.h:547
virtual string reactionType(size_t i) const
String specifying the type of reaction.
Definition Kinetics.cpp:344
void skipUndeclaredSpecies(bool skip)
Determine behavior when adding a new reaction that contains species not defined in any of the phases ...
Definition Kinetics.h:1400
string kineticsSpeciesName(size_t k) const
Return the name of the kth species in the kinetics manager.
Definition Kinetics.cpp:258
virtual void getDestructionRates(double *ddot)
Species destruction rates [kmol/m^3/s or kmol/m^2/s].
Definition Kinetics.cpp:424
AnyMap parameters()
Return the parameters for a phase definition which are needed to reconstruct an identical object usin...
Definition Kinetics.cpp:637
Kinetics()=default
Default constructor.
virtual void getDeltaEntropy(double *deltaS)
Return the vector of values for the reactions change in entropy.
Definition Kinetics.h:493
virtual void init()
Prepare the class for the addition of reactions, after all phases have been added.
Definition Kinetics.h:1354
string reactantString(size_t i) const
Returns a string containing the reactants side of the reaction equation.
Definition Kinetics.cpp:361
virtual void getFwdRateConstants(double *kfwd)
Return the forward rate constants.
Definition Kinetics.h:1291
shared_ptr< Solution > root() const
Get the Solution object containing this Kinetics object and associated ThermoPhase objects.
Definition Kinetics.h:1485
bool m_skipUndeclaredThirdBodies
See skipUndeclaredThirdBodies()
Definition Kinetics.h:1626
virtual double reactantOrder(size_t k, size_t i) const
Reactant order of species k in reaction i.
Definition Kinetics.h:1203
virtual void modifyReaction(size_t i, shared_ptr< Reaction > rNew)
Modify the rate expression associated with a reaction.
Definition Kinetics.cpp:753
vector< double > m_ropnet
Net rate-of-progress for each reaction.
Definition Kinetics.h:1614
void skipUndeclaredThirdBodies(bool skip)
Determine behavior when adding a new reaction that contains third-body efficiencies for species not d...
Definition Kinetics.h:1412
virtual double productStoichCoeff(size_t k, size_t i) const
Stoichiometric coefficient of species k as a product in reaction i.
Definition Kinetics.cpp:339
virtual void addThermo(shared_ptr< ThermoPhase > thermo)
Add a phase to the kinetics manager object.
Definition Kinetics.cpp:592
Eigen::SparseMatrix< double > productStoichCoeffs() const
Stoichiometric coefficient matrix for products.
Definition Kinetics.h:1184
vector< double > m_rbuf
Buffer used for storage of intermediate reaction-specific results.
Definition Kinetics.h:1620
Eigen::SparseMatrix< double > m_stoichMatrix
Net stoichiometry (products - reactants)
Definition Kinetics.h:1534
virtual void getActivityConcentrations(double *const conc)
Get the vector of activity concentrations used in the kinetics object.
Definition Kinetics.h:1226
bool m_skipUndeclaredSpecies
See skipUndeclaredSpecies()
Definition Kinetics.h:1623
map< string, size_t > m_phaseindex
Mapping of the phase name to the position of the phase within the kinetics object.
Definition Kinetics.h:1582
size_t m_mindim
number of spatial dimensions of lowest-dimensional phase.
Definition Kinetics.h:1596
StoichManagerN m_productStoich
Stoichiometry manager for the products for each reaction.
Definition Kinetics.h:1528
virtual void getDeltaSSGibbs(double *deltaG)
Return the vector of values for the reaction standard state Gibbs free energy change.
Definition Kinetics.h:507
std::weak_ptr< Solution > m_root
reference to Solution
Definition Kinetics.h:1632
virtual void getDeltaSSEnthalpy(double *deltaH)
Return the vector of values for the change in the standard state enthalpies of reaction.
Definition Kinetics.h:521
StoichManagerN m_revProductStoich
Stoichiometry manager for the products of reversible reactions.
Definition Kinetics.h:1531
Eigen::SparseMatrix< double > revProductStoichCoeffs() const
Stoichiometric coefficient matrix for products of reversible reactions.
Definition Kinetics.h:1191
virtual double reactionEnthalpy(const Composition &reactants, const Composition &products)
Calculate the reaction enthalpy of a reaction which has not necessarily been added into the Kinetics ...
Definition Kinetics.cpp:796
size_t nReactions() const
Number of reactions in the reaction mechanism.
Definition Kinetics.h:153
virtual void getRevRatesOfProgress(double *revROP)
Return the Reverse rates of progress of the reactions.
Definition Kinetics.cpp:380
virtual void setRoot(shared_ptr< Solution > root)
Set root Solution holding all phase information.
Definition Kinetics.h:1479
bool m_hasUndeclaredThirdBodies
Flag indicating whether reactions include undeclared third bodies.
Definition Kinetics.h:1629
size_t kineticsSpeciesIndex(size_t k, size_t n) const
The location of species k of phase n in species arrays.
Definition Kinetics.h:288
virtual double productOrder(int k, int i) const
product Order of species k in reaction i.
Definition Kinetics.h:1217
size_t m_surfphase
Index in the list of phases of the one surface phase.
Definition Kinetics.h:1586
virtual void setMultiplier(size_t i, double f)
Set the multiplier for reaction i to f.
Definition Kinetics.h:1443
vector< double > m_rfn
Forward rate constant for each reaction.
Definition Kinetics.h:1599
virtual void getDeltaEnthalpy(double *deltaH)
Return the vector of values for the reactions change in enthalpy.
Definition Kinetics.h:480
StoichManagerN m_reactantStoich
Stoichiometry manager for the reactants for each reaction.
Definition Kinetics.h:1525
size_t surfacePhaseIndex() const
This returns the integer index of the phase which has ThermoPhase type cSurf.
Definition Kinetics.cpp:76
virtual void getRevRateConstants(double *krev, bool doIrreversible=false)
Return the reverse rate constants.
Definition Kinetics.h:1311
virtual void resizeSpecies()
Resize arrays with sizes that depend on the total number of species.
Definition Kinetics.cpp:653
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
vector< shared_ptr< ThermoPhase > > m_sharedThermo
vector of shared pointers,
Definition Kinetics.h:1568
size_t nTotalSpecies() const
The total number of species in all phases participating in the kinetics mechanism.
Definition Kinetics.h:266
virtual double reactantStoichCoeff(size_t k, size_t i) const
Stoichiometric coefficient of species k as a reactant in reaction i.
Definition Kinetics.cpp:334
string productString(size_t i) const
Returns a string containing the products side of the reaction equation.
Definition Kinetics.cpp:368
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 getDeltaGibbs(double *deltaG)
Return the vector of values for the reaction Gibbs free energy change.
Definition Kinetics.h:452
virtual void getNetProductionRates(double *wdot)
Species net production rates [kmol/m^3/s or kmol/m^2/s].
Definition Kinetics.cpp:435
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:302
vector< double > m_delta_gibbs0
Delta G^0 for all reactions.
Definition Kinetics.h:1602
virtual void getCreationRates(double *cdot)
Species creation rates [kmol/m^3/s or kmol/m^2/s].
Definition Kinetics.cpp:410
ThermoPhase & speciesPhase(size_t k)
This function takes as an argument the kineticsSpecies index (that is, the list index in the list of ...
Definition Kinetics.h:350
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:323
An error indicating that an unimplemented function has been called.
This class handles operations involving the stoichiometric coefficients on one side of a reaction (re...
const Eigen::SparseMatrix< double > & stoichCoeffs() const
Return matrix containing stoichiometric coefficients.
Base class for a phase with thermodynamic properties.
Storage for cached values.
Definition ValueCache.h:153
void clear()
Clear all cached values.
virtual void setDerivativeSettings(const AnyMap &settings)
Set/modify derivative settings.
Definition Kinetics.h:732
virtual Eigen::SparseMatrix< double > fwdRatesOfProgress_ddCi()
Calculate derivatives for forward rates-of-progress with respect to species concentration at constant...
Definition Kinetics.h:841
Eigen::SparseMatrix< double > creationRates_ddCi()
Calculate derivatives for species creation rates with respect to species concentration at constant te...
Definition Kinetics.cpp:492
virtual Eigen::SparseMatrix< double > netRatesOfProgress_ddCi()
Calculate derivatives for net rates-of-progress with respect to species concentration at constant tem...
Definition Kinetics.h:987
void getCreationRates_ddT(double *dwdot)
Calculate derivatives for species creation rates with respect to temperature at constant pressure,...
Definition Kinetics.cpp:446
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:967
Eigen::SparseMatrix< double > destructionRates_ddX()
Calculate derivatives for species destruction rates with respect to species mole fractions at constan...
Definition Kinetics.cpp:538
void getCreationRates_ddC(double *dwdot)
Calculate derivatives for species creation rates with respect to molar concentration at constant temp...
Definition Kinetics.cpp:470
void getDestructionRates_ddP(double *dwdot)
Calculate derivatives for species destruction rates with respect to pressure at constant temperature,...
Definition Kinetics.cpp:514
void getNetProductionRates_ddC(double *dwdot)
Calculate derivatives for species net production rates with respect to molar concentration at constan...
Definition Kinetics.cpp:574
void getDestructionRates_ddT(double *dwdot)
Calculate derivatives for species destruction rates with respect to temperature at constant pressure,...
Definition Kinetics.cpp:502
virtual void getFwdRatesOfProgress_ddP(double *drop)
Calculate derivatives for forward rates-of-progress with respect to pressure at constant temperature,...
Definition Kinetics.h:790
Eigen::SparseMatrix< double > netProductionRates_ddX()
Calculate derivatives for species net production rates with respect to species mole fractions at cons...
Definition Kinetics.cpp:582
Eigen::SparseMatrix< double > creationRates_ddX()
Calculate derivatives for species creation rates with respect to species mole fractions at constant t...
Definition Kinetics.cpp:482
Eigen::SparseMatrix< double > destructionRates_ddCi()
Calculate derivatives for species destruction rates with respect to species concentration at constant...
Definition Kinetics.cpp:548
void getNetProductionRates_ddT(double *dwdot)
Calculate derivatives for species net production rates with respect to temperature at constant pressu...
Definition Kinetics.cpp:558
virtual Eigen::SparseMatrix< double > fwdRatesOfProgress_ddX()
Calculate derivatives for forward rates-of-progress with respect to species mole fractions at constan...
Definition Kinetics.h:821
virtual void getRevRatesOfProgress_ddT(double *drop)
Calculate derivatives for reverse rates-of-progress with respect to temperature at constant pressure,...
Definition Kinetics.h:852
virtual void getFwdRateConstants_ddP(double *dkfwd)
Calculate derivatives for forward rate constants with respect to pressure at constant temperature,...
Definition Kinetics.h:754
virtual void getFwdRateConstants_ddT(double *dkfwd)
Calculate derivatives for forward rate constants with respect to temperature at constant pressure,...
Definition Kinetics.h:743
void getCreationRates_ddP(double *dwdot)
Calculate derivatives for species creation rates with respect to pressure at constant temperature,...
Definition Kinetics.cpp:458
virtual void getNetRatesOfProgress_ddT(double *drop)
Calculate derivatives for net rates-of-progress with respect to temperature at constant pressure,...
Definition Kinetics.h:925
virtual void getRevRatesOfProgress_ddP(double *drop)
Calculate derivatives for reverse rates-of-progress with respect to pressure at constant temperature,...
Definition Kinetics.h:863
Eigen::SparseMatrix< double > netProductionRates_ddCi()
Calculate derivatives for species net production rates with respect to species concentration at const...
Definition Kinetics.cpp:587
virtual void getRevRatesOfProgress_ddC(double *drop)
Calculate derivatives for reverse rates-of-progress with respect to molar concentration at constant t...
Definition Kinetics.h:877
virtual void getDerivativeSettings(AnyMap &settings) const
Retrieve derivative settings.
Definition Kinetics.h:721
void getNetProductionRates_ddP(double *dwdot)
Calculate derivatives for species net production rates with respect to pressure at constant temperatu...
Definition Kinetics.cpp:566
virtual void getNetRatesOfProgress_ddP(double *drop)
Calculate derivatives for net rates-of-progress with respect to pressure at constant temperature,...
Definition Kinetics.h:936
virtual Eigen::SparseMatrix< double > revRatesOfProgress_ddCi()
Calculate derivatives for forward rates-of-progress with respect to species concentration at constant...
Definition Kinetics.h:914
virtual void getFwdRatesOfProgress_ddT(double *drop)
Calculate derivatives for forward rates-of-progress with respect to temperature at constant pressure,...
Definition Kinetics.h:779
void getDestructionRates_ddC(double *dwdot)
Calculate derivatives for species destruction rates with respect to molar concentration at constant t...
Definition Kinetics.cpp:526
virtual void getNetRatesOfProgress_ddC(double *drop)
Calculate derivatives for net rates-of-progress with respect to molar concentration at constant tempe...
Definition Kinetics.h:950
virtual Eigen::SparseMatrix< double > revRatesOfProgress_ddX()
Calculate derivatives for reverse rates-of-progress with respect to species mole fractions at constan...
Definition Kinetics.h:894
virtual void getFwdRateConstants_ddC(double *dkfwd)
Calculate derivatives for forward rate constants with respect to molar concentration at constant temp...
Definition Kinetics.h:768
virtual void getFwdRatesOfProgress_ddC(double *drop)
Calculate derivatives for forward rates-of-progress with respect to molar concentration at constant t...
Definition Kinetics.h:804
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:195
map< string, double > Composition
Map from string names to doubles.
Definition ct_defs.h:184