Cantera  2.3.0
StoichManager.h
Go to the documentation of this file.
1 /**
2  * @file StoichManager.h
3  */
4 
5 // This file is part of Cantera. See License.txt in the top-level directory or
6 // at http://www.cantera.org/license.txt for license and copyright information.
7 
8 #ifndef CT_STOICH_MGR_H
9 #define CT_STOICH_MGR_H
10 
13 
14 namespace Cantera
15 {
16 
17 /**
18  * @defgroup Stoichiometry Stoichiometry
19  *
20  * Note: these classes are designed for internal use in class
21  * ReactionStoichManager.
22  *
23  * The classes defined here implement simple operations that are used by class
24  * ReactionStoichManager to compute things like rates of progress, species
25  * production rates, etc. In general, a reaction mechanism may involve many
26  * species and many reactions, but any given reaction typically only involves
27  * a few species as reactants, and a few as products. Therefore, the matrix of
28  * stoichiometric coefficients is very sparse. Not only is it sparse, but the
29  * non-zero matrix elements often have the value 1, and in many cases no more
30  * than three coefficients are non-zero for the reactants and/or the products.
31  *
32  * For the present purposes, we will consider each direction of a reversible
33  * reaction to be a separate reaction. We often need to compute quantities
34  * that can formally be written as a matrix product of a stoichiometric
35  * coefficient matrix and a vector of reaction rates. For example, the species
36  * creation rates are given by
37  *
38  * \f[
39  * \dot C_k = \sum_k \nu^{(p)}_{k,i} R_i
40  * \f]
41  *
42  * where \f$ \nu^{(p)_{k,i}} \f$ is the product-side stoichiometric
43  * coefficient of species \a k in reaction \a i. This could be done be
44  * straightforward matrix multiplication, but would be inefficient, since most
45  * of the matrix elements of \f$ \nu^{(p)}_{k,i} \f$ are zero. We could do
46  * better by using sparse-matrix algorithms to compute this product.
47  *
48  * If the reactions are general ones, with non-integral stoichiometric
49  * coefficients, this is about as good as we can do. But we are particularly
50  * concerned here with the performance for very large reaction mechanisms,
51  * which are usually composed of elementary reactions, which have integral
52  * stoichiometric coefficients. Furthermore, very few elementary reactions
53  * involve more than 3 product or reactant molecules.
54  *
55  * But we can do even better if we take account of the special structure of
56  * this matrix for elementary reactions involving three or fewer product
57  * molecules (or reactant molecules).
58  *
59  * To take advantage of this structure, reactions are divided into four groups.
60  * These classes are designed to take advantage of this sparse structure when
61  * computing quantities that can be written as matrix multiplies.
62  *
63  * They are designed to explicitly unroll loops over species or reactions for
64  * Operations on reactions that require knowing the reaction stoichiometry.
65  *
66  * This module consists of class StoichManager, and classes C1, C2, and C3.
67  * Classes C1, C2, and C3 handle operations involving one, two, or three
68  * species, respectively, in a reaction. Instances are instantiated with a
69  * reaction number, and n species numbers (n = 1 for C1, etc.). All three
70  * classes have the same interface.
71  *
72  * These classes are designed for use by StoichManager, and the operations
73  * implemented are those needed to efficiently compute quantities such as
74  * rates of progress, species production rates, reaction thermochemistry, etc.
75  * The compiler will inline these methods into the body of the corresponding
76  * StoichManager method, and so there is no performance penalty (unless
77  * inlining is turned off).
78  *
79  * To describe the methods, consider class C3 and suppose an instance is
80  * created with reaction number irxn and species numbers k0, k1, and k2.
81  *
82  * - multiply(in, out) : out[irxn] is multiplied by
83  * in[k0] * in[k1] * in[k2]
84  *
85  * - power(in, out) : out[irxn] is multiplied by
86  * (in[k0]^order0) * (in[k1]^order1) * (in[k2]^order2)
87  *
88  * - incrementReaction(in, out) : out[irxn] is incremented by
89  * in[k0] + in[k1] + in[k2]
90  *
91  * - decrementReaction(in, out) : out[irxn] is decremented by
92  * in[k0] + in[k1] + in[k2]
93  *
94  * - incrementSpecies(in, out) : out[k0], out[k1], and out[k2]
95  * are all incremented by in[irxn]
96  *
97  * - decrementSpecies(in, out) : out[k0], out[k1], and out[k2]
98  * are all decremented by in[irxn]
99  *
100  * The function multiply() is usually used when evaluating the forward and
101  * reverse rates of progress of reactions. The rate constants are usually
102  * loaded into out[]. Then multiply() is called to add in the dependence of
103  * the species concentrations to yield a forward and reverse rop.
104  *
105  * The function incrementSpecies() and its cousin decrementSpecies() is used
106  * to translate from rates of progress to species production rates. The vector
107  * in[] is preloaded with the rates of progress of all reactions. Then
108  * incrementSpecies() is called to increment the species production vector,
109  * out[], with the rates of progress.
110  *
111  * The functions incrementReaction() and decrementReaction() are used to find
112  * the standard state equilibrium constant for a reaction. Here, output[] is a
113  * vector of length number of reactions, usually the standard Gibbs free
114  * energies of reaction, while input, usually the standard state Gibbs free
115  * energies of species, is a vector of length number of species.
116  *
117  * Note the stoichiometric coefficient for a species in a reaction is handled
118  * by always assuming it is equal to one and then treating reactants and
119  * products for a reaction separately. Bimolecular reactions involving the
120  * identical species are treated as involving separate species.
121  */
122 
123 /**
124  * Handles one species in a reaction.
125  * See @ref Stoichiometry
126  * @ingroup Stoichiometry
127  * @internal
128  */
129 class C1
130 {
131 public:
132  C1(size_t rxn = 0, size_t ic0 = 0) :
133  m_rxn(rxn),
134  m_ic0(ic0) {
135  }
136 
137  size_t data(std::vector<size_t>& ic) {
138  ic.resize(1);
139  ic[0] = m_ic0;
140  return m_rxn;
141  }
142 
143  void incrementSpecies(const doublereal* R, doublereal* S) const {
144  S[m_ic0] += R[m_rxn];
145  }
146 
147  void decrementSpecies(const doublereal* R, doublereal* S) const {
148  S[m_ic0] -= R[m_rxn];
149  }
150 
151  void multiply(const doublereal* S, doublereal* R) const {
152  R[m_rxn] *= S[m_ic0];
153  }
154 
155  void incrementReaction(const doublereal* S, doublereal* R) const {
156  R[m_rxn] += S[m_ic0];
157  }
158 
159  void decrementReaction(const doublereal* S, doublereal* R) const {
160  R[m_rxn] -= S[m_ic0];
161  }
162 
163  size_t rxnNumber() const {
164  return m_rxn;
165  }
166  size_t speciesIndex(size_t n) const {
167  return m_ic0;
168  }
169  size_t nSpecies() {
170  return 1;
171  }
172 
173 private:
174  //! Reaction number
175  size_t m_rxn;
176  //! Species number
177  size_t m_ic0;
178 };
179 
180 /**
181  * Handles two species in a single reaction.
182  * See @ref Stoichiometry
183  * @ingroup Stoichiometry
184  */
185 class C2
186 {
187 public:
188  C2(size_t rxn = 0, size_t ic0 = 0, size_t ic1 = 0)
189  : m_rxn(rxn), m_ic0(ic0), m_ic1(ic1) {}
190 
191  size_t data(std::vector<size_t>& ic) {
192  ic.resize(2);
193  ic[0] = m_ic0;
194  ic[1] = m_ic1;
195  return m_rxn;
196  }
197 
198  void incrementSpecies(const doublereal* R, doublereal* S) const {
199  S[m_ic0] += R[m_rxn];
200  S[m_ic1] += R[m_rxn];
201  }
202 
203  void decrementSpecies(const doublereal* R, doublereal* S) const {
204  S[m_ic0] -= R[m_rxn];
205  S[m_ic1] -= R[m_rxn];
206  }
207 
208  void multiply(const doublereal* S, doublereal* R) const {
209  if (S[m_ic0] < 0 && S[m_ic1] < 0) {
210  R[m_rxn] = 0;
211  } else {
212  R[m_rxn] *= S[m_ic0] * S[m_ic1];
213  }
214  }
215 
216  void incrementReaction(const doublereal* S, doublereal* R) const {
217  R[m_rxn] += S[m_ic0] + S[m_ic1];
218  }
219 
220  void decrementReaction(const doublereal* S, doublereal* R) const {
221  R[m_rxn] -= (S[m_ic0] + S[m_ic1]);
222  }
223 
224  size_t rxnNumber() const {
225  return m_rxn;
226  }
227  size_t speciesIndex(size_t n) const {
228  return (n == 0 ? m_ic0 : m_ic1);
229  }
230  size_t nSpecies() {
231  return 2;
232  }
233 
234 private:
235  //! Reaction index -> index into the ROP vector
236  size_t m_rxn;
237 
238  //! Species index -> index into the species vector for the two species.
239  size_t m_ic0, m_ic1;
240 };
241 
242 /**
243  * Handles three species in a reaction.
244  * See @ref Stoichiometry
245  * @ingroup Stoichiometry
246  */
247 class C3
248 {
249 public:
250  C3(size_t rxn = 0, size_t ic0 = 0, size_t ic1 = 0, size_t ic2 = 0)
251  : m_rxn(rxn), m_ic0(ic0), m_ic1(ic1), m_ic2(ic2) {}
252 
253  size_t data(std::vector<size_t>& ic) {
254  ic.resize(3);
255  ic[0] = m_ic0;
256  ic[1] = m_ic1;
257  ic[2] = m_ic2;
258  return m_rxn;
259  }
260 
261  void incrementSpecies(const doublereal* R, doublereal* S) const {
262  S[m_ic0] += R[m_rxn];
263  S[m_ic1] += R[m_rxn];
264  S[m_ic2] += R[m_rxn];
265  }
266 
267  void decrementSpecies(const doublereal* R, doublereal* S) const {
268  S[m_ic0] -= R[m_rxn];
269  S[m_ic1] -= R[m_rxn];
270  S[m_ic2] -= R[m_rxn];
271  }
272 
273  void multiply(const doublereal* S, doublereal* R) const {
274  if ((S[m_ic0] < 0 && (S[m_ic1] < 0 || S[m_ic2] < 0)) ||
275  (S[m_ic1] < 0 && S[m_ic2] < 0)) {
276  R[m_rxn] = 0;
277  } else {
278  R[m_rxn] *= S[m_ic0] * S[m_ic1] * S[m_ic2];
279  }
280  }
281 
282  void incrementReaction(const doublereal* S, doublereal* R) const {
283  R[m_rxn] += S[m_ic0] + S[m_ic1] + S[m_ic2];
284  }
285 
286  void decrementReaction(const doublereal* S, doublereal* R) const {
287  R[m_rxn] -= (S[m_ic0] + S[m_ic1] + S[m_ic2]);
288  }
289 
290  size_t rxnNumber() const {
291  return m_rxn;
292  }
293  size_t speciesIndex(size_t n) const {
294  return (n == 0 ? m_ic0 : (n == 1 ? m_ic1 : m_ic2));
295  }
296  size_t nSpecies() {
297  return 3;
298  }
299 
300 private:
301  size_t m_rxn;
302  size_t m_ic0;
303  size_t m_ic1;
304  size_t m_ic2;
305 };
306 
307 /**
308  * Handles any number of species in a reaction, including fractional
309  * stoichiometric coefficients, and arbitrary reaction orders.
310  * See @ref Stoichiometry
311  * @ingroup Stoichiometry
312  */
313 class C_AnyN
314 {
315 public:
316  C_AnyN() :
317  m_n(0),
318  m_rxn(npos) {
319  }
320 
321  C_AnyN(size_t rxn, const std::vector<size_t>& ic, const vector_fp& order_, const vector_fp& stoich_) :
322  m_n(0),
323  m_rxn(rxn) {
324  m_n = ic.size();
325  m_ic.resize(m_n);
326  m_order.resize(m_n);
327  m_stoich.resize(m_n);
328  for (size_t n = 0; n < m_n; n++) {
329  m_ic[n] = ic[n];
330  m_order[n] = order_[n];
331  m_stoich[n] = stoich_[n];
332  }
333  }
334 
335  size_t data(std::vector<size_t>& ic) {
336  ic.resize(m_n);
337  for (size_t n = 0; n < m_n; n++) {
338  ic[n] = m_ic[n];
339  }
340  return m_rxn;
341  }
342 
343  doublereal order(size_t n) const {
344  return m_order[n];
345  }
346  doublereal stoich(size_t n) const {
347  return m_stoich[n];
348  }
349  size_t speciesIndex(size_t n) const {
350  return m_ic[n];
351  }
352 
353  void multiply(const doublereal* input, doublereal* output) const {
354  for (size_t n = 0; n < m_n; n++) {
355  double order = m_order[n];
356  if (order != 0.0) {
357  double c = input[m_ic[n]];
358  if (c > 0.0) {
359  output[m_rxn] *= std::pow(c, order);
360  } else {
361  output[m_rxn] = 0.0;
362  }
363  }
364  }
365  }
366 
367  void incrementSpecies(const doublereal* input,
368  doublereal* output) const {
369  doublereal x = input[m_rxn];
370  for (size_t n = 0; n < m_n; n++) {
371  output[m_ic[n]] += m_stoich[n]*x;
372  }
373  }
374 
375  void decrementSpecies(const doublereal* input,
376  doublereal* output) const {
377  doublereal x = input[m_rxn];
378  for (size_t n = 0; n < m_n; n++) {
379  output[m_ic[n]] -= m_stoich[n]*x;
380  }
381  }
382 
383  void incrementReaction(const doublereal* input,
384  doublereal* output) const {
385  for (size_t n = 0; n < m_n; n++) {
386  output[m_rxn] += m_stoich[n]*input[m_ic[n]];
387  }
388  }
389 
390  void decrementReaction(const doublereal* input,
391  doublereal* output) const {
392  for (size_t n = 0; n < m_n; n++) {
393  output[m_rxn] -= m_stoich[n]*input[m_ic[n]];
394  }
395  }
396 
397 private:
398  //! Length of the m_ic vector
399  /*!
400  * This is the number of species which participate in the reaction order
401  * and stoichiometric coefficient vectors for the reactant or product
402  * description of the reaction.
403  */
404  size_t m_n;
405 
406  //! ID of the reaction corresponding to this stoichiometric manager
407  /*!
408  * This is used within the interface to select the array position to read
409  * and write to Normally this is associated with the reaction number in an
410  * array of quantities indexed by the reaction number, e.g., ROP[irxn].
411  */
412  size_t m_rxn;
413 
414  //! Vector of species which are involved with this stoichiometric manager
415  //! calculations
416  /*!
417  * This is an integer list of species which participate in either the
418  * reaction order matrix or the stoichiometric order matrix for this
419  * reaction, m_rxn.
420  */
421  std::vector<size_t> m_ic;
422 
423  //! Reaction orders for the reaction
424  /*!
425  * This is either for the reactants or products. Length = m_n. Species
426  * number, m_ic[n], has a reaction order of m_order[n].
427  */
429 
430  //! Stoichiometric coefficients for the reaction, reactant or product side.
431  /*!
432  * This is either for the reactants or products. Length = m_n. Species
433  * number m_ic[m], has a stoichiometric coefficient of m_stoich[n].
434  */
436 };
437 
438 template<class InputIter, class Vec1, class Vec2>
439 inline static void _multiply(InputIter begin, InputIter end,
440  const Vec1& input, Vec2& output)
441 {
442  for (; begin != end; ++begin) {
443  begin->multiply(input, output);
444  }
445 }
446 
447 template<class InputIter, class Vec1, class Vec2>
448 inline static void _incrementSpecies(InputIter begin,
449  InputIter end, const Vec1& input, Vec2& output)
450 {
451  for (; begin != end; ++begin) {
452  begin->incrementSpecies(input, output);
453  }
454 }
455 
456 template<class InputIter, class Vec1, class Vec2>
457 inline static void _decrementSpecies(InputIter begin,
458  InputIter end, const Vec1& input, Vec2& output)
459 {
460  for (; begin != end; ++begin) {
461  begin->decrementSpecies(input, output);
462  }
463 }
464 
465 template<class InputIter, class Vec1, class Vec2>
466 inline static void _incrementReactions(InputIter begin,
467  InputIter end, const Vec1& input, Vec2& output)
468 {
469  for (; begin != end; ++begin) {
470  begin->incrementReaction(input, output);
471  }
472 }
473 
474 template<class InputIter, class Vec1, class Vec2>
475 inline static void _decrementReactions(InputIter begin,
476  InputIter end, const Vec1& input, Vec2& output)
477 {
478  for (; begin != end; ++begin) {
479  begin->decrementReaction(input, output);
480  }
481 }
482 
483 /*
484  * This class handles operations involving the stoichiometric coefficients on
485  * one side of a reaction (reactant or product) for a set of reactions
486  * comprising a reaction mechanism. This class is used by class Kinetics, which
487  * contains three instances of this class (one to handle operations on the
488  * reactions, one for the products of reversible reactions, and one for the
489  * products of irreversible reactions).
490  *
491  * This class is designed for use with elementary reactions, or at least ones
492  * with integral stoichiometric coefficients. Let \f$ M(i) \f$ be the number of
493  * molecules on the product or reactant side of reaction number i.
494  * \f[
495  * r_i = \sum_m^{M_i} s_{k_{m,i}}
496  * \f]
497  * To understand the operations performed by this class, let \f$ N_{k,i}\f$
498  * denote the stoichiometric coefficient of species k on one side (reactant or
499  * product) in reaction i. Then \b N is a sparse K by I matrix of stoichiometric
500  * coefficients.
501  *
502  * The following matrix operations may be carried out with a vector S of length
503  * K, and a vector R of length I:
504  *
505  * - \f$ S = S + N R\f$ (incrementSpecies)
506  * - \f$ S = S - N R\f$ (decrementSpecies)
507  * - \f$ R = R + N^T S \f$ (incrementReaction)
508  * - \f$ R = R - N^T S \f$ (decrementReaction)
509  *
510  * The actual implementation, however, does not compute these quantities by
511  * matrix multiplication. A faster algorithm is used that makes use of the fact
512  * that the \b integer-valued N matrix is very sparse, and the non-zero terms
513  * are small positive integers.
514  * \f[
515  * S_k = R_{i1} + \dots + R_{iM}
516  * \f]
517  * where M is the number of molecules, and $\f i(m) \f$ is the
518  * See @ref Stoichiometry
519  * @ingroup Stoichiometry
520  */
521 class StoichManagerN
522 {
523 public:
524  /**
525  * Constructor for the StoichManagerN class.
526  *
527  * @internal Consider adding defaulted entries here that supply the total
528  * number of reactions in the mechanism and the total number of species
529  * in the species list. Then, we could use those numbers to provide
530  * error checks during the construction of the object. Those numbers
531  * would also provide some clarity to the purpose and utility of this
532  * class.
533  *
534  * DGG - the problem is that the number of reactions and species are not
535  * known initially.
536  */
537  StoichManagerN() {
538  }
539 
540  /**
541  * Add a single reaction to the list of reactions that this stoichiometric
542  * manager object handles.
543  *
544  * This function is the same as the add() function below. However, the order
545  * of each species in the power list expression is set to one automatically.
546  */
547  void add(size_t rxn, const std::vector<size_t>& k) {
548  vector_fp order(k.size(), 1.0);
549  vector_fp stoich(k.size(), 1.0);
550  add(rxn, k, order, stoich);
551  }
552 
553  void add(size_t rxn, const std::vector<size_t>& k, const vector_fp& order) {
554  vector_fp stoich(k.size(), 1.0);
555  add(rxn, k, order, stoich);
556  }
557 
558  //! Add a single reaction to the list of reactions that this
559  //! stoichiometric manager object handles.
560  /*!
561  * @param rxn Reaction index of the current reaction. This is used as an
562  * index into vectors which have length n_total_rxn.
563  * @param k This is a vector of integer values specifying the species
564  * indices. The length of this vector species the number of different
565  * species in the description. The value of the entries are the species
566  * indices. These are used as indexes into vectors which have length
567  * n_total_species.
568  * @param order This is a vector of the same length as vector k. The order
569  * is used for the routine power(), which produces a power law
570  * expression involving the species vector.
571  * @param stoich This is used to handle fractional stoichiometric
572  * coefficients on the product side of irreversible reactions.
573  */
574  void add(size_t rxn, const std::vector<size_t>& k, const vector_fp& order,
575  const vector_fp& stoich) {
576  if (order.size() != k.size()) {
577  throw CanteraError("StoichManagerN::add()", "size of order and species arrays differ");
578  }
579  if (stoich.size() != k.size()) {
580  throw CanteraError("StoichManagerN::add()", "size of stoich and species arrays differ");
581  }
582  bool frac = false;
583  for (size_t n = 0; n < stoich.size(); n++) {
584  if (fmod(stoich[n], 1.0) || stoich[n] != order[n]) {
585  frac = true;
586  break;
587  }
588  }
589  if (frac || k.size() > 3) {
590  m_cn_list.emplace_back(rxn, k, order, stoich);
591  } else {
592  // Try to express the reaction with unity stoichiometric
593  // coefficients (by repeating species when necessary) so that the
594  // simpler 'multiply' function can be used to compute the rate
595  // instead of 'power'.
596  std::vector<size_t> kRep;
597  for (size_t n = 0; n < k.size(); n++) {
598  for (size_t i = 0; i < stoich[n]; i++) {
599  kRep.push_back(k[n]);
600  }
601  }
602 
603  switch (kRep.size()) {
604  case 1:
605  m_c1_list.emplace_back(rxn, kRep[0]);
606  break;
607  case 2:
608  m_c2_list.emplace_back(rxn, kRep[0], kRep[1]);
609  break;
610  case 3:
611  m_c3_list.emplace_back(rxn, kRep[0], kRep[1], kRep[2]);
612  break;
613  default:
614  m_cn_list.emplace_back(rxn, k, order, stoich);
615  }
616  }
617  }
618 
619  void multiply(const doublereal* input, doublereal* output) const {
620  _multiply(m_c1_list.begin(), m_c1_list.end(), input, output);
621  _multiply(m_c2_list.begin(), m_c2_list.end(), input, output);
622  _multiply(m_c3_list.begin(), m_c3_list.end(), input, output);
623  _multiply(m_cn_list.begin(), m_cn_list.end(), input, output);
624  }
625 
626  void incrementSpecies(const doublereal* input, doublereal* output) const {
627  _incrementSpecies(m_c1_list.begin(), m_c1_list.end(), input, output);
628  _incrementSpecies(m_c2_list.begin(), m_c2_list.end(), input, output);
629  _incrementSpecies(m_c3_list.begin(), m_c3_list.end(), input, output);
630  _incrementSpecies(m_cn_list.begin(), m_cn_list.end(), input, output);
631  }
632 
633  void decrementSpecies(const doublereal* input, doublereal* output) const {
634  _decrementSpecies(m_c1_list.begin(), m_c1_list.end(), input, output);
635  _decrementSpecies(m_c2_list.begin(), m_c2_list.end(), input, output);
636  _decrementSpecies(m_c3_list.begin(), m_c3_list.end(), input, output);
637  _decrementSpecies(m_cn_list.begin(), m_cn_list.end(), input, output);
638  }
639 
640  void incrementReactions(const doublereal* input, doublereal* output) const {
641  _incrementReactions(m_c1_list.begin(), m_c1_list.end(), input, output);
642  _incrementReactions(m_c2_list.begin(), m_c2_list.end(), input, output);
643  _incrementReactions(m_c3_list.begin(), m_c3_list.end(), input, output);
644  _incrementReactions(m_cn_list.begin(), m_cn_list.end(), input, output);
645  }
646 
647  void decrementReactions(const doublereal* input, doublereal* output) const {
648  _decrementReactions(m_c1_list.begin(), m_c1_list.end(), input, output);
649  _decrementReactions(m_c2_list.begin(), m_c2_list.end(), input, output);
650  _decrementReactions(m_c3_list.begin(), m_c3_list.end(), input, output);
651  _decrementReactions(m_cn_list.begin(), m_cn_list.end(), input, output);
652  }
653 
654 private:
655  std::vector<C1> m_c1_list;
656  std::vector<C2> m_c2_list;
657  std::vector<C3> m_c3_list;
658  std::vector<C_AnyN> m_cn_list;
659 };
660 
661 }
662 
663 #endif
size_t m_n
Length of the m_ic vector.
size_t m_rxn
Reaction index -> index into the ROP vector.
const size_t npos
index returned by functions to indicate "no position"
Definition: ct_defs.h:165
size_t m_ic0
Species number.
size_t m_rxn
Reaction number.
size_t m_ic0
Species index -> index into the species vector for the two species.
Handles two species in a single reaction.
std::vector< size_t > m_ic
Vector of species which are involved with this stoichiometric manager calculations.
Handles any number of species in a reaction, including fractional stoichiometric coefficients, and arbitrary reaction orders.
void multiply(const DenseMatrix &A, const double *const b, double *const prod)
Multiply A*b and return the result in prod. Uses BLAS routine DGEMV.
Handles one species in a reaction.
std::vector< double > vector_fp
Turn on the use of stl vectors for the basic array type within cantera Vector of doubles.
Definition: ct_defs.h:157
Contains declarations for string manipulation functions within Cantera.
Handles three species in a reaction.
Namespace for the Cantera kernel.
Definition: application.cpp:29
size_t m_rxn
ID of the reaction corresponding to this stoichiometric manager.
vector_fp m_stoich
Stoichiometric coefficients for the reaction, reactant or product side.
vector_fp m_order
Reaction orders for the reaction.
Definitions for the classes that are thrown when Cantera experiences an error condition (also contain...