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