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