Cantera  3.1.0
Loading...
Searching...
No Matches
Derivative Calculations
Collaboration diagram for Derivative Calculations:

Detailed Description

Methods for calculating analytical and/or numerical derivatives.

Classes

class  MultiJac
 Class MultiJac evaluates the Jacobian of a system of equations defined by a residual function supplied by an instance of class OneDim. More...
 

Functions

Eigen::SparseMatrix< double > jacobian ()
 Return semi-analytical Jacobian of an AdaptivePreconditioner object.
 
MultiJacjacobian ()
 Return a reference to the Jacobian evaluator of an OneDim object.
 
virtual Eigen::SparseMatrix< double > jacobian ()
 Calculate the Jacobian of a specific Reactor specialization.
 

Routines to Calculate Kinetics Derivatives (Jacobians)

Kinetics derivatives are calculated with respect to temperature, pressure, molar concentrations and species mole fractions for forward/reverse/net rates of progress as well as creation/destruction and net production of species.

The following suffixes are used to indicate derivatives:

  • _ddT: derivative with respect to temperature (a vector)
  • _ddP: derivative with respect to pressure (a vector)
  • _ddC: derivative with respect to molar concentration (a vector)
  • _ddX: derivative with respect to species mole fractions (a matrix)
  • _ddCi: derivative with respect to species concentrations (a matrix)
Since
New in Cantera 2.6
Warning
The calculation of kinetics derivatives is an experimental part of the Cantera API and may be changed or removed without notice.

Source term derivatives are based on a generic rate-of-progress expression for the \( i \)-th reaction \( R_i \), which is a function of temperature \( T \), pressure \( P \) and molar concentrations \( C_j \):

\[ R_i = k_{f,i} C_M^{\nu_{M,i}} \prod_j C_j^{\nu_{ji}^\prime} - k_{r,i} C_M^{\nu_{M,i}} \prod_j C_j^{\nu_{ji}^{\prime\prime}} \]

Forward/reverse rate expressions \( k_{f,i} \) and \( k_{r,i} \) are implemented by ReactionRate specializations; forward/reverse stoichiometric coefficients are \( \nu_{ji}^\prime \) and \( \nu_{ji}^{\prime\prime} \). Unless the reaction involves third-body colliders, \( \nu_{M,i} = 0 \). For three-body reactions, effective ThirdBody collider concentrations \( C_M \) are considered with \( \nu_{M,i} = 1 \). For more detailed information on relevant theory, see, for example, Perini, et al. [34] or Niemeyer, et al. [31], although specifics of Cantera's implementation may differ.

Partial derivatives are obtained from the product rule, where resulting terms consider reaction rate derivatives, derivatives of the concentration product term, and, if applicable, third-body term derivatives. ReactionRate specializations may implement exact derivatives (example: ArrheniusRate::ddTScaledFromStruct) or approximate them numerically (examples: ReactionData::perturbTemperature, PlogData::perturbPressure, FalloffData::perturbThirdBodies). Derivatives of concentration and third-body terms are based on analytic expressions.

Species creation and destruction rates are obtained by multiplying rate-of-progress vectors by stoichiometric coefficient matrices. As this is a linear operation, it is possible to calculate derivatives the same way.

All derivatives are calculated for source terms while holding other properties constant, independent of whether equation of state or \( \sum X_k = 1 \) constraints are satisfied. Thus, derivatives deviate from Jacobians and numerical derivatives that implicitly enforce these constraints. Depending on application and equation of state, derivatives can nevertheless be used to obtain Jacobians, for example:

  • The Jacobian of net production rates \( \dot{\omega}_{k,\mathrm{net}} \) with respect to temperature at constant pressure needs to consider changes of molar density \( C \) due to temperature

    \[ \left. \frac{\partial \dot{\omega}_{k,\mathrm{net}}}{\partial T} \right|_{P=\mathrm{const}} = \frac{\partial \dot{\omega}_{k,\mathrm{net}}}{\partial T} + \frac{\partial \dot{\omega}_{k,\mathrm{net}}}{\partial C} \left. \frac{\partial C}{\partial T} \right|_{P=\mathrm{const}} \]

    where for an ideal gas \( \partial C / \partial T = - C / T \). The remaining partial derivatives are obtained from getNetProductionRates_ddT() and getNetProductionRates_ddC(), respectively.
  • The Jacobian of \( \dot{\omega}_{k,\mathrm{net}} \) with respect to temperature at constant volume needs to consider pressure changes due to temperature

    \[ \left. \frac{\partial \dot{\omega}_{k,\mathrm{net}}}{\partial T} \right|_{V=\mathrm{const}} = \frac{\partial \dot{\omega}_{k,\mathrm{net}}}{\partial T} + \frac{\partial \dot{\omega}_{k,\mathrm{net}}}{\partial P} \left. \frac{\partial P}{\partial T} \right|_{V=\mathrm{const}} \]

    where for an ideal gas \( \partial P / \partial T = P / T \). The remaining partial derivatives are obtained from getNetProductionRates_ddT() and getNetProductionRates_ddP(), respectively.
  • Similar expressions can be derived for other derivatives and source terms.

While some applications require exact derivatives, others can tolerate approximate derivatives that neglect terms to increase computational speed and/or improve Jacobian sparsity (example: AdaptivePreconditioner). Derivative evaluations settings are accessible by keyword/value pairs using the methods getDerivativeSettings() and setDerivativeSettings().

For BulkKinetics, the following keyword/value pairs are supported:

  • skip-third-bodies (boolean): if false (default), third body concentrations are considered for the evaluation of Jacobians
  • skip-falloff (boolean): if false (default), third-body effects on rate constants are considered for the evaluation of derivatives.
  • rtol-delta (double): relative tolerance used to perturb properties when calculating numerical derivatives. The default value is 1e-8.

For InterfaceKinetics, the following keyword/value pairs are supported:

  • skip-coverage-dependence (boolean): if false (default), rate constant coverage dependence is not considered when evaluating derivatives.
  • skip-electrochemistry (boolean): if false (default), electrical charge is not considered in evaluating the derivatives and these reactions are treated as normal surface reactions.
  • rtol-delta (double): relative tolerance used to perturb properties when calculating numerical derivatives. The default value is 1e-8.
virtual void getDerivativeSettings (AnyMap &settings) const
 Retrieve derivative settings.
 
virtual void setDerivativeSettings (const AnyMap &settings)
 Set/modify derivative settings.
 
virtual void getFwdRateConstants_ddT (double *dkfwd)
 Calculate derivatives for forward rate constants with respect to temperature at constant pressure, molar concentration and mole fractions.
 
virtual void getFwdRateConstants_ddP (double *dkfwd)
 Calculate derivatives for forward rate constants with respect to pressure at constant temperature, molar concentration and mole fractions.
 
virtual void getFwdRateConstants_ddC (double *dkfwd)
 Calculate derivatives for forward rate constants with respect to molar concentration at constant temperature, pressure and mole fractions.
 
virtual void getFwdRatesOfProgress_ddT (double *drop)
 Calculate derivatives for forward rates-of-progress with respect to temperature at constant pressure, molar concentration and mole fractions.
 
virtual void getFwdRatesOfProgress_ddP (double *drop)
 Calculate derivatives for forward rates-of-progress with respect to pressure at constant temperature, molar concentration and mole fractions.
 
virtual void getFwdRatesOfProgress_ddC (double *drop)
 Calculate derivatives for forward rates-of-progress with respect to molar concentration at constant temperature, pressure and mole fractions.
 
virtual Eigen::SparseMatrix< double > fwdRatesOfProgress_ddX ()
 Calculate derivatives for forward rates-of-progress with respect to species mole fractions at constant temperature, pressure and molar concentration.
 
virtual Eigen::SparseMatrix< double > fwdRatesOfProgress_ddCi ()
 Calculate derivatives for forward rates-of-progress with respect to species concentration at constant temperature, pressure and remaining species concentrations.
 
virtual void getRevRatesOfProgress_ddT (double *drop)
 Calculate derivatives for reverse rates-of-progress with respect to temperature at constant pressure, molar concentration and mole fractions.
 
virtual void getRevRatesOfProgress_ddP (double *drop)
 Calculate derivatives for reverse rates-of-progress with respect to pressure at constant temperature, molar concentration and mole fractions.
 
virtual void getRevRatesOfProgress_ddC (double *drop)
 Calculate derivatives for reverse rates-of-progress with respect to molar concentration at constant temperature, pressure and mole fractions.
 
virtual Eigen::SparseMatrix< double > revRatesOfProgress_ddX ()
 Calculate derivatives for reverse rates-of-progress with respect to species mole fractions at constant temperature, pressure and molar concentration.
 
virtual Eigen::SparseMatrix< double > revRatesOfProgress_ddCi ()
 Calculate derivatives for forward rates-of-progress with respect to species concentration at constant temperature, pressure and remaining species concentrations.
 
virtual void getNetRatesOfProgress_ddT (double *drop)
 Calculate derivatives for net rates-of-progress with respect to temperature at constant pressure, molar concentration and mole fractions.
 
virtual void getNetRatesOfProgress_ddP (double *drop)
 Calculate derivatives for net rates-of-progress with respect to pressure at constant temperature, molar concentration and mole fractions.
 
virtual void getNetRatesOfProgress_ddC (double *drop)
 Calculate derivatives for net rates-of-progress with respect to molar concentration at constant temperature, pressure and mole fractions.
 
virtual Eigen::SparseMatrix< double > netRatesOfProgress_ddX ()
 Calculate derivatives for net rates-of-progress with respect to species mole fractions at constant temperature, pressure and molar concentration.
 
virtual Eigen::SparseMatrix< double > netRatesOfProgress_ddCi ()
 Calculate derivatives for net rates-of-progress with respect to species concentration at constant temperature, pressure, and remaining species concentrations.
 
void getCreationRates_ddT (double *dwdot)
 Calculate derivatives for species creation rates with respect to temperature at constant pressure, molar concentration and mole fractions.
 
void getCreationRates_ddP (double *dwdot)
 Calculate derivatives for species creation rates with respect to pressure at constant temperature, molar concentration and mole fractions.
 
void getCreationRates_ddC (double *dwdot)
 Calculate derivatives for species creation rates with respect to molar concentration at constant temperature, pressure and mole fractions.
 
Eigen::SparseMatrix< double > creationRates_ddX ()
 Calculate derivatives for species creation rates with respect to species mole fractions at constant temperature, pressure and molar concentration.
 
Eigen::SparseMatrix< double > creationRates_ddCi ()
 Calculate derivatives for species creation rates with respect to species concentration at constant temperature, pressure, and concentration of all other species.
 
void getDestructionRates_ddT (double *dwdot)
 Calculate derivatives for species destruction rates with respect to temperature at constant pressure, molar concentration and mole fractions.
 
void getDestructionRates_ddP (double *dwdot)
 Calculate derivatives for species destruction rates with respect to pressure at constant temperature, molar concentration and mole fractions.
 
void getDestructionRates_ddC (double *dwdot)
 Calculate derivatives for species destruction rates with respect to molar concentration at constant temperature, pressure and mole fractions.
 
Eigen::SparseMatrix< double > destructionRates_ddX ()
 Calculate derivatives for species destruction rates with respect to species mole fractions at constant temperature, pressure and molar concentration.
 
Eigen::SparseMatrix< double > destructionRates_ddCi ()
 Calculate derivatives for species destruction rates with respect to species concentration at constant temperature, pressure, and concentration of all other species.
 
void getNetProductionRates_ddT (double *dwdot)
 Calculate derivatives for species net production rates with respect to temperature at constant pressure, molar concentration and mole fractions.
 
void getNetProductionRates_ddP (double *dwdot)
 Calculate derivatives for species net production rates with respect to pressure at constant temperature, molar concentration and mole fractions.
 
void getNetProductionRates_ddC (double *dwdot)
 Calculate derivatives for species net production rates with respect to molar concentration at constant temperature, pressure and mole fractions.
 
Eigen::SparseMatrix< double > netProductionRates_ddX ()
 Calculate derivatives for species net production rates with respect to species mole fractions at constant temperature, pressure and molar concentration.
 
Eigen::SparseMatrix< double > netProductionRates_ddCi ()
 Calculate derivatives for species net production rates with respect to species concentration at constant temperature, pressure, and concentration of all other species.
 

Function Documentation

◆ getDerivativeSettings()

virtual void getDerivativeSettings ( AnyMap settings) const
inlinevirtual

Retrieve derivative settings.

Parameters
settingsAnyMap containing settings determining derivative evaluation.

Reimplemented in BulkKinetics, and InterfaceKinetics.

Definition at line 692 of file Kinetics.h.

◆ setDerivativeSettings()

virtual void setDerivativeSettings ( const AnyMap settings)
inlinevirtual

Set/modify derivative settings.

Parameters
settingsAnyMap containing settings determining derivative evaluation.

Reimplemented in BulkKinetics, and InterfaceKinetics.

Definition at line 703 of file Kinetics.h.

◆ getFwdRateConstants_ddT()

virtual void getFwdRateConstants_ddT ( double *  dkfwd)
inlinevirtual

Calculate derivatives for forward rate constants with respect to temperature at constant pressure, molar concentration and mole fractions.

Parameters
[out]dkfwdOutput vector of derivatives. Length: nReactions().

Reimplemented in BulkKinetics.

Definition at line 714 of file Kinetics.h.

◆ getFwdRateConstants_ddP()

virtual void getFwdRateConstants_ddP ( double *  dkfwd)
inlinevirtual

Calculate derivatives for forward rate constants with respect to pressure at constant temperature, molar concentration and mole fractions.

Parameters
[out]dkfwdOutput vector of derivatives. Length: nReactions().

Reimplemented in BulkKinetics.

Definition at line 725 of file Kinetics.h.

◆ getFwdRateConstants_ddC()

virtual void getFwdRateConstants_ddC ( double *  dkfwd)
inlinevirtual

Calculate derivatives for forward rate constants with respect to molar concentration at constant temperature, pressure and mole fractions.

Parameters
[out]dkfwdOutput vector of derivatives. Length: nReactions().
Warning
This method is an experimental part of the Cantera API and may be changed or removed without notice.

Reimplemented in BulkKinetics.

Definition at line 739 of file Kinetics.h.

◆ getFwdRatesOfProgress_ddT()

virtual void getFwdRatesOfProgress_ddT ( double *  drop)
inlinevirtual

Calculate derivatives for forward rates-of-progress with respect to temperature at constant pressure, molar concentration and mole fractions.

Parameters
[out]dropOutput vector of derivatives. Length: nReactions().

Reimplemented in BulkKinetics.

Definition at line 750 of file Kinetics.h.

◆ getFwdRatesOfProgress_ddP()

virtual void getFwdRatesOfProgress_ddP ( double *  drop)
inlinevirtual

Calculate derivatives for forward rates-of-progress with respect to pressure at constant temperature, molar concentration and mole fractions.

Parameters
[out]dropOutput vector of derivatives. Length: nReactions().

Reimplemented in BulkKinetics.

Definition at line 761 of file Kinetics.h.

◆ getFwdRatesOfProgress_ddC()

virtual void getFwdRatesOfProgress_ddC ( double *  drop)
inlinevirtual

Calculate derivatives for forward rates-of-progress with respect to molar concentration at constant temperature, pressure and mole fractions.

Parameters
[out]dropOutput vector of derivatives. Length: nReactions().
Warning
This method is an experimental part of the Cantera API and may be changed or removed without notice.

Reimplemented in BulkKinetics.

Definition at line 775 of file Kinetics.h.

◆ fwdRatesOfProgress_ddX()

virtual Eigen::SparseMatrix< double > fwdRatesOfProgress_ddX ( )
inlinevirtual

Calculate derivatives for forward rates-of-progress with respect to species mole fractions at constant temperature, pressure and molar concentration.

The method returns a matrix with nReactions() rows and nTotalSpecies() columns. For a derivative with respect to \( X_i \), all other \( X_j \) are held constant, rather than enforcing \( \sum X_j = 1 \).

Warning
This method is an experimental part of the Cantera API and may be changed or removed without notice.

Reimplemented in BulkKinetics.

Definition at line 792 of file Kinetics.h.

◆ fwdRatesOfProgress_ddCi()

virtual Eigen::SparseMatrix< double > fwdRatesOfProgress_ddCi ( )
inlinevirtual

Calculate derivatives for forward rates-of-progress with respect to species concentration at constant temperature, pressure and remaining species concentrations.

The method returns a matrix with nReactions() rows and nTotalSpecies() columns. For a derivative with respect to \( c_i \), all other \( c_j \) are held constant.

Warning
This method is an experimental part of the Cantera API and may be changed or removed without notice.
Since
New in Cantera 3.0.

Reimplemented in BulkKinetics, and InterfaceKinetics.

Definition at line 812 of file Kinetics.h.

◆ getRevRatesOfProgress_ddT()

virtual void getRevRatesOfProgress_ddT ( double *  drop)
inlinevirtual

Calculate derivatives for reverse rates-of-progress with respect to temperature at constant pressure, molar concentration and mole fractions.

Parameters
[out]dropOutput vector of derivatives. Length: nReactions().

Reimplemented in BulkKinetics.

Definition at line 823 of file Kinetics.h.

◆ getRevRatesOfProgress_ddP()

virtual void getRevRatesOfProgress_ddP ( double *  drop)
inlinevirtual

Calculate derivatives for reverse rates-of-progress with respect to pressure at constant temperature, molar concentration and mole fractions.

Parameters
[out]dropOutput vector of derivatives. Length: nReactions().

Reimplemented in BulkKinetics.

Definition at line 834 of file Kinetics.h.

◆ getRevRatesOfProgress_ddC()

virtual void getRevRatesOfProgress_ddC ( double *  drop)
inlinevirtual

Calculate derivatives for reverse rates-of-progress with respect to molar concentration at constant temperature, pressure and mole fractions.

Parameters
[out]dropOutput vector of derivatives. Length: nReactions().
Warning
This method is an experimental part of the Cantera API and may be changed or removed without notice.

Reimplemented in BulkKinetics.

Definition at line 848 of file Kinetics.h.

◆ revRatesOfProgress_ddX()

virtual Eigen::SparseMatrix< double > revRatesOfProgress_ddX ( )
inlinevirtual

Calculate derivatives for reverse rates-of-progress with respect to species mole fractions at constant temperature, pressure and molar concentration.

The method returns a matrix with nReactions() rows and nTotalSpecies() columns. For a derivative with respect to \( X_i \), all other \( X_j \) are held constant, rather than enforcing \( \sum X_j = 1 \).

Warning
This method is an experimental part of the Cantera API and may be changed or removed without notice.

Reimplemented in BulkKinetics.

Definition at line 865 of file Kinetics.h.

◆ revRatesOfProgress_ddCi()

virtual Eigen::SparseMatrix< double > revRatesOfProgress_ddCi ( )
inlinevirtual

Calculate derivatives for forward rates-of-progress with respect to species concentration at constant temperature, pressure and remaining species concentrations.

The method returns a matrix with nReactions() rows and nTotalSpecies() columns. For a derivative with respect to \( c_i \), all other \( c_j \) are held constant.

Warning
This method is an experimental part of the Cantera API and may be changed or removed without notice.
Since
New in Cantera 3.0.

Reimplemented in BulkKinetics, and InterfaceKinetics.

Definition at line 885 of file Kinetics.h.

◆ getNetRatesOfProgress_ddT()

virtual void getNetRatesOfProgress_ddT ( double *  drop)
inlinevirtual

Calculate derivatives for net rates-of-progress with respect to temperature at constant pressure, molar concentration and mole fractions.

Parameters
[out]dropOutput vector of derivatives. Length: nReactions().

Reimplemented in BulkKinetics.

Definition at line 896 of file Kinetics.h.

◆ getNetRatesOfProgress_ddP()

virtual void getNetRatesOfProgress_ddP ( double *  drop)
inlinevirtual

Calculate derivatives for net rates-of-progress with respect to pressure at constant temperature, molar concentration and mole fractions.

Parameters
[out]dropOutput vector of derivatives. Length: nReactions().

Reimplemented in BulkKinetics.

Definition at line 907 of file Kinetics.h.

◆ getNetRatesOfProgress_ddC()

virtual void getNetRatesOfProgress_ddC ( double *  drop)
inlinevirtual

Calculate derivatives for net rates-of-progress with respect to molar concentration at constant temperature, pressure and mole fractions.

Parameters
[out]dropOutput vector of derivatives. Length: nReactions().
Warning
This method is an experimental part of the Cantera API and may be changed or removed without notice.

Reimplemented in BulkKinetics.

Definition at line 921 of file Kinetics.h.

◆ netRatesOfProgress_ddX()

virtual Eigen::SparseMatrix< double > netRatesOfProgress_ddX ( )
inlinevirtual

Calculate derivatives for net rates-of-progress with respect to species mole fractions at constant temperature, pressure and molar concentration.

The method returns a matrix with nReactions() rows and nTotalSpecies() columns. For a derivative with respect to \( X_i \), all other \( X_j \) are held constant, rather than enforcing \( \sum X_j = 1 \).

Warning
This method is an experimental part of the Cantera API and may be changed or removed without notice.

Reimplemented in BulkKinetics.

Definition at line 938 of file Kinetics.h.

◆ netRatesOfProgress_ddCi()

virtual Eigen::SparseMatrix< double > netRatesOfProgress_ddCi ( )
inlinevirtual

Calculate derivatives for net rates-of-progress with respect to species concentration at constant temperature, pressure, and remaining species concentrations.

The method returns a matrix with nReactions() rows and nTotalSpecies() columns. For a derivative with respect to \( c_i \), all other \( c_j \) are held constant.

Warning
This method is an experimental part of the Cantera API and may be changed or removed without notice.
Since
New in Cantera 3.0.

Reimplemented in BulkKinetics, and InterfaceKinetics.

Definition at line 958 of file Kinetics.h.

◆ getCreationRates_ddT()

void getCreationRates_ddT ( double *  dwdot)

Calculate derivatives for species creation rates with respect to temperature at constant pressure, molar concentration and mole fractions.

Parameters
[out]dwdotOutput vector of derivatives. Length: m_kk.

Definition at line 424 of file Kinetics.cpp.

◆ getCreationRates_ddP()

void getCreationRates_ddP ( double *  dwdot)

Calculate derivatives for species creation rates with respect to pressure at constant temperature, molar concentration and mole fractions.

Parameters
[out]dwdotOutput vector of derivatives. Length: m_kk.

Definition at line 436 of file Kinetics.cpp.

◆ getCreationRates_ddC()

void getCreationRates_ddC ( double *  dwdot)

Calculate derivatives for species creation rates with respect to molar concentration at constant temperature, pressure and mole fractions.

Parameters
[out]dwdotOutput vector of derivatives. Length: m_kk.
Warning
This method is an experimental part of the Cantera API and may be changed or removed without notice.

Definition at line 448 of file Kinetics.cpp.

◆ creationRates_ddX()

Eigen::SparseMatrix< double > creationRates_ddX ( )

Calculate derivatives for species creation rates with respect to species mole fractions at constant temperature, pressure and molar concentration.

The method returns a square matrix with nTotalSpecies() rows and columns. For a derivative with respect to \( X_i \), all other \( X_j \) are held constant, rather than enforcing \( \sum X_j = 1 \).

Warning
This method is an experimental part of the Cantera API and may be changed or removed without notice.

Definition at line 460 of file Kinetics.cpp.

◆ creationRates_ddCi()

Eigen::SparseMatrix< double > creationRates_ddCi ( )

Calculate derivatives for species creation rates with respect to species concentration at constant temperature, pressure, and concentration of all other species.

The method returns a square matrix with nTotalSpecies() rows and columns. For a derivative with respect to \( c_i \), all other \( c_j \) are held constant.

Warning
This method is an experimental part of the Cantera API and may be changed or removed without notice.
Since
New in Cantera 3.0.

Definition at line 470 of file Kinetics.cpp.

◆ getDestructionRates_ddT()

void getDestructionRates_ddT ( double *  dwdot)

Calculate derivatives for species destruction rates with respect to temperature at constant pressure, molar concentration and mole fractions.

Parameters
[out]dwdotOutput vector of derivatives. Length: m_kk.

Definition at line 480 of file Kinetics.cpp.

◆ getDestructionRates_ddP()

void getDestructionRates_ddP ( double *  dwdot)

Calculate derivatives for species destruction rates with respect to pressure at constant temperature, molar concentration and mole fractions.

Parameters
[out]dwdotOutput vector of derivatives. Length: m_kk.

Definition at line 492 of file Kinetics.cpp.

◆ getDestructionRates_ddC()

void getDestructionRates_ddC ( double *  dwdot)

Calculate derivatives for species destruction rates with respect to molar concentration at constant temperature, pressure and mole fractions.

Parameters
[out]dwdotOutput vector of derivatives. Length: m_kk.
Warning
This method is an experimental part of the Cantera API and may be changed or removed without notice.

Definition at line 504 of file Kinetics.cpp.

◆ destructionRates_ddX()

Eigen::SparseMatrix< double > destructionRates_ddX ( )

Calculate derivatives for species destruction rates with respect to species mole fractions at constant temperature, pressure and molar concentration.

The method returns a square matrix with nTotalSpecies() rows and columns. For a derivative with respect to \( X_i \), all other \( X_j \) are held constant, rather than enforcing \( \sum X_j = 1 \).

Warning
This method is an experimental part of the Cantera API and may be changed or removed without notice.

Definition at line 516 of file Kinetics.cpp.

◆ destructionRates_ddCi()

Eigen::SparseMatrix< double > destructionRates_ddCi ( )

Calculate derivatives for species destruction rates with respect to species concentration at constant temperature, pressure, and concentration of all other species.

The method returns a square matrix with nTotalSpecies() rows and columns. For a derivative with respect to \( c_i \), all other \( c_j \) are held constant.

Warning
This method is an experimental part of the Cantera API and may be changed or removed without notice.
Since
New in Cantera 3.0.

Definition at line 526 of file Kinetics.cpp.

◆ getNetProductionRates_ddT()

void getNetProductionRates_ddT ( double *  dwdot)

Calculate derivatives for species net production rates with respect to temperature at constant pressure, molar concentration and mole fractions.

Parameters
[out]dwdotOutput vector of derivatives. Length: m_kk.

Definition at line 536 of file Kinetics.cpp.

◆ getNetProductionRates_ddP()

void getNetProductionRates_ddP ( double *  dwdot)

Calculate derivatives for species net production rates with respect to pressure at constant temperature, molar concentration and mole fractions.

Parameters
[out]dwdotOutput vector of derivatives. Length: m_kk.

Definition at line 544 of file Kinetics.cpp.

◆ getNetProductionRates_ddC()

void getNetProductionRates_ddC ( double *  dwdot)

Calculate derivatives for species net production rates with respect to molar concentration at constant temperature, pressure and mole fractions.

Parameters
[out]dwdotOutput vector of derivatives. Length: m_kk.
Warning
This method is an experimental part of the Cantera API and may be changed or removed without notice.

Definition at line 552 of file Kinetics.cpp.

◆ netProductionRates_ddX()

Eigen::SparseMatrix< double > netProductionRates_ddX ( )

Calculate derivatives for species net production rates with respect to species mole fractions at constant temperature, pressure and molar concentration.

The method returns a square matrix with nTotalSpecies() rows and columns. For a derivative with respect to \( X_i \), all other \( X_j \) are held constant, rather than enforcing \( \sum X_j = 1 \).

Warning
This method is an experimental part of the Cantera API and may be changed or removed without notice.

Definition at line 560 of file Kinetics.cpp.

◆ netProductionRates_ddCi()

Eigen::SparseMatrix< double > netProductionRates_ddCi ( )

Calculate derivatives for species net production rates with respect to species concentration at constant temperature, pressure, and concentration of all other species.

The method returns a square matrix with nTotalSpecies() rows and columns. For a derivative with respect to \( c_i \), all other \( c_j \) are held constant.

Warning
This method is an experimental part of the Cantera API and may be changed or removed without notice.
Since
New in Cantera 3.0.

Definition at line 565 of file Kinetics.cpp.

◆ jacobian() [1/3]

Eigen::SparseMatrix< double > jacobian ( )
inline

Return semi-analytical Jacobian of an AdaptivePreconditioner object.

Definition at line 53 of file AdaptivePreconditioner.h.

◆ jacobian() [2/3]

MultiJac & jacobian ( )

Return a reference to the Jacobian evaluator of an OneDim object.

Definition at line 87 of file OneDim.cpp.

◆ jacobian() [3/3]

virtual Eigen::SparseMatrix< double > jacobian ( )
inlinevirtual

Calculate the Jacobian of a specific Reactor specialization.

Warning
Depending on the particular implementation, this may return an approximate Jacobian intended only for use in forming a preconditioner for iterative solvers.
This method is an experimental part of the Cantera API and may be changed or removed without notice.

Reimplemented in IdealGasConstPressureMoleReactor, and IdealGasMoleReactor.

Definition at line 212 of file Reactor.h.