Cantera  3.2.0a1
Loading...
Searching...
No Matches

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

#include <MultiJac.h>

Inheritance diagram for MultiJac:
[legend]

Detailed Description

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

The residual function may consist of several linked 1D domains, with different variables in each domain.

Definition at line 24 of file MultiJac.h.

Public Member Functions

 MultiJac (OneDim &r)
 Constructor.
 
void eval (double *x0, double *resid0, double rdt)
 Evaluates the Jacobian at x0 using finite differences.
 
void reset () override
 Reset parameters as needed.
 
const string type () const override
 Derived type, corresponding to names registered with SystemJacobianFactory.
 
void setValue (size_t row, size_t col, double value) override
 Set a value at the specified row and column of the jacobian triplet vector.
 
void initialize (size_t nVars) override
 Called during setup for any processes that need to be completed prior to setup functions used in sundials.
 
void setBandwidth (size_t bw) override
 Used to provide system bandwidth for implementations that use banded matrix storage.
 
void updateTransient (double rdt, integer *mask) override
 Update the diagonal terms in the Jacobian by using the transient mask \( \alpha \).
 
void factorize () override
 Factorize the system matrix.
 
double & value (size_t i, size_t j)
 
double value (size_t i, size_t j) const
 
void solve (const double *const b, double *const x)
 
void solve (const size_t stateSize, double *b, double *x) override
 Solve a linear system using the system matrix M
 
int info () const override
 Get latest status of linear solver.
 
vector< int > & transientMask ()
 Return the transient mask.
 
- Public Member Functions inherited from SystemJacobian
virtual const string type () const =0
 Derived type, corresponding to names registered with SystemJacobianFactory.
 
virtual void setValue (size_t row, size_t col, double value)
 Set a value at the specified row and column of the jacobian triplet vector.
 
virtual void stateAdjustment (vector< double > &state)
 Adjust the state vector based on the preconditioner, e.g., Adaptive preconditioning uses a strictly positive composition when preconditioning which is handled by this function.
 
string preconditionerSide () const
 Get preconditioner application side for CVODES.
 
virtual void setPreconditionerSide (const string &preconSide)
 For iterative solvers, set the side where the preconditioner is applied.
 
virtual void updatePreconditioner ()
 Transform Jacobian vector and write into preconditioner, P = (I - gamma * J)
 
virtual void updateTransient (double rdt, int *mask)
 Update the diagonal terms in the Jacobian by using the transient mask \( \alpha \).
 
virtual void solve (const size_t stateSize, double *rhs_vector, double *output)
 Solve a linear system using the system matrix M
 
virtual void setup ()
 Perform preconditioner specific post-reactor setup operations such as factorize.
 
virtual void reset ()
 Reset parameters as needed.
 
virtual void initialize (size_t networkSize)
 Called during setup for any processes that need to be completed prior to setup functions used in sundials.
 
virtual void setBandwidth (size_t bw)
 Used to provide system bandwidth for implementations that use banded matrix storage.
 
virtual void printPreconditioner ()
 Print preconditioner contents.
 
virtual void setGamma (double gamma)
 Set gamma used in preconditioning.
 
virtual double gamma ()
 Get gamma used in preconditioning.
 
virtual void setAbsoluteTolerance (double atol)
 Set the absolute tolerance in the solver outside of the network initialization.
 
virtual int info () const
 Get latest status of linear solver.
 
double elapsedTime () const
 Elapsed CPU time spent computing the Jacobian elements.
 
void updateElapsed (double evalTime)
 Increase the elapsed time by the specified time.
 
int nEvals () const
 Number of Jacobian evaluations.
 
void incrementEvals ()
 Increment the number of times the Jacobian has been evaluated.
 
int age () const
 Number of times 'incrementAge' has been called since the last evaluation.
 
void incrementAge ()
 Increment the Jacobian age.
 
void setAge (int age)
 Set the Jacobian age.
 
void clearStats ()
 Clear collected stats about number of Jacobian evaluations, CPU time spent on Jacobian updates, and the number of Newton steps that have been taken since updating the Jacobian.
 

Protected Attributes

OneDimm_resid = nullptr
 Residual evaluator for this Jacobian.
 
BandMatrix m_mat
 Underlying matrix storage.
 
vector< double > m_ssdiag
 Diagonal of the steady-state Jacobian.
 
vector< int > m_mask
 Transient mask for transient terms, 1 if transient, 0 if steady-state.
 
- Protected Attributes inherited from SystemJacobian
size_t m_dim
 Dimension of the system.
 
double m_gamma = 1.0
 gamma value used in M = I - gamma*J
 
bool m_init = false
 bool saying whether or not the system is initialized
 
double m_atol = 0
 Absolute tolerance of the ODE solver.
 
double m_elapsed = 0.0
 Elapsed CPU time taken to compute the Jacobian.
 
int m_nevals = 0
 Number of Jacobian evaluations.
 
int m_age = 100000
 Age of the Jacobian (times incrementAge() has been called)
 
string m_precon_side = "none"
 For iterative solvers, side of the system to apply the preconditioner on.
 

Additional Inherited Members

virtual void factorize ()
 Factorize the system matrix.
 

Constructor & Destructor Documentation

◆ MultiJac()

MultiJac ( OneDim r)

Constructor.

Parameters
rThe nonlinear system for which to compute the Jacobian.
Deprecated:
To be removed after Cantera 3.2. Use default constructor instead.

Definition at line 12 of file MultiJac.cpp.

Member Function Documentation

◆ eval()

void eval ( double *  x0,
double *  resid0,
double  rdt 
)

Evaluates the Jacobian at x0 using finite differences.

The unperturbed residual function is resid0, which must be supplied on input. The parameter 'rdt' is the reciprocal of the time step. If zero, the steady-state Jacobian is evaluated, otherwise the transient Jacobian is evaluated.

Each variable in the input vector x0 is perturbed to compute the corresponding column of the Jacobian matrix. The Jacobian is computed by perturbing each variable, evaluating the residual function, and then estimating the partial derivatives numerically using finite differences.

Parameters
x0Solution vector at which to evaluate the Jacobian
resid0Residual vector at x0
rdtReciprocal of the time step
Deprecated:
To be removed after Cantera 3.2. Jacobian evaluation moved to OneDim::evalJacobian().

Definition at line 56 of file MultiJac.cpp.

◆ reset()

void reset ( )
overridevirtual

Reset parameters as needed.

Reimplemented from SystemJacobian.

Definition at line 21 of file MultiJac.cpp.

◆ type()

const string type ( ) const
inlineoverridevirtual

Derived type, corresponding to names registered with SystemJacobianFactory.

Implements SystemJacobian.

Definition at line 54 of file MultiJac.h.

◆ setValue()

void setValue ( size_t  row,
size_t  col,
double  value 
)
overridevirtual

Set a value at the specified row and column of the jacobian triplet vector.

Parameters
rowrow in the jacobian matrix
colcolumn in the jacobian matrix
valuevalue of the element at row and col

Reimplemented from SystemJacobian.

Definition at line 40 of file MultiJac.cpp.

◆ initialize()

void initialize ( size_t  networkSize)
overridevirtual

Called during setup for any processes that need to be completed prior to setup functions used in sundials.

Parameters
networkSizethe number of variables in the associated reactor network

Reimplemented from SystemJacobian.

Definition at line 27 of file MultiJac.cpp.

◆ setBandwidth()

void setBandwidth ( size_t  bw)
overridevirtual

Used to provide system bandwidth for implementations that use banded matrix storage.

Ignored if not needed.

Reimplemented from SystemJacobian.

Definition at line 35 of file MultiJac.cpp.

◆ updateTransient()

void updateTransient ( double  rdt,
integer *  mask 
)
overridevirtual

Update the diagonal terms in the Jacobian by using the transient mask \( \alpha \).

Parameters
rdtReciprocal of the time step [1/s]
maskMask for transient terms: 1 if transient, 0 if algebraic.

Reimplemented from SystemJacobian.

Definition at line 48 of file MultiJac.cpp.

◆ factorize()

void factorize ( )
inlineoverridevirtual

Factorize the system matrix.

This method should be called at the end of implementations of updatePreconditioner() and updateTransient().

Reimplemented from SystemJacobian.

Definition at line 59 of file MultiJac.h.

◆ value() [1/2]

double & value ( size_t  i,
size_t  j 
)
inline

Definition at line 63 of file MultiJac.h.

◆ value() [2/2]

double value ( size_t  i,
size_t  j 
) const
inline

Definition at line 67 of file MultiJac.h.

◆ solve() [1/2]

void solve ( const double *const  b,
double *const  x 
)
inline

Definition at line 71 of file MultiJac.h.

◆ solve() [2/2]

void solve ( const size_t  stateSize,
double *  rhs_vector,
double *  output 
)
inlineoverridevirtual

Solve a linear system using the system matrix M

The matrix M can either be the preconditioner or the transient Jacobian matrix, depending on whether updateTransient() or updatePreconditioner() was called last.

Parameters
[in]stateSizelength of the rhs and output vectors
[in]rhs_vectorright hand side vector used in linear system
[out]outputoutput vector for solution

Reimplemented from SystemJacobian.

Definition at line 75 of file MultiJac.h.

◆ info()

int info ( ) const
inlineoverridevirtual

Get latest status of linear solver.

Zero indicates success. Meaning of non-zero values is solver dependent.

Reimplemented from SystemJacobian.

Definition at line 79 of file MultiJac.h.

◆ transientMask()

vector< int > & transientMask ( )
inline

Return the transient mask.

Deprecated:
Unused. To be removed after Cantera 3.2.

Definition at line 85 of file MultiJac.h.

Member Data Documentation

◆ m_resid

OneDim* m_resid = nullptr
protected

Residual evaluator for this Jacobian.

This is a pointer to the residual evaluator. This object isn't owned by this Jacobian object.

Deprecated:
Unused. To be removed after Cantera 3.2.

Definition at line 98 of file MultiJac.h.

◆ m_mat

BandMatrix m_mat
protected

Underlying matrix storage.

Definition at line 100 of file MultiJac.h.

◆ m_ssdiag

vector<double> m_ssdiag
protected

Diagonal of the steady-state Jacobian.

Definition at line 101 of file MultiJac.h.

◆ m_mask

vector<int> m_mask
protected

Transient mask for transient terms, 1 if transient, 0 if steady-state.

Deprecated:
Unused. To be removed after Cantera 3.2.

Definition at line 105 of file MultiJac.h.


The documentation for this class was generated from the following files: