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>
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. | |
![]() | |
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 | |
OneDim * | m_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. | |
![]() | |
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.
r | The nonlinear system for which to compute the Jacobian. |
Definition at line 12 of file MultiJac.cpp.
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.
x0 | Solution vector at which to evaluate the Jacobian |
resid0 | Residual vector at x0 |
rdt | Reciprocal of the time step |
Definition at line 56 of file MultiJac.cpp.
|
overridevirtual |
Reset parameters as needed.
Reimplemented from SystemJacobian.
Definition at line 21 of file MultiJac.cpp.
|
inlineoverridevirtual |
Derived type, corresponding to names registered with SystemJacobianFactory.
Implements SystemJacobian.
Definition at line 54 of file MultiJac.h.
|
overridevirtual |
Set a value at the specified row and column of the jacobian triplet vector.
row | row in the jacobian matrix |
col | column in the jacobian matrix |
value | value of the element at row and col |
Reimplemented from SystemJacobian.
Definition at line 40 of file MultiJac.cpp.
|
overridevirtual |
Called during setup for any processes that need to be completed prior to setup functions used in sundials.
networkSize | the number of variables in the associated reactor network |
Reimplemented from SystemJacobian.
Definition at line 27 of file MultiJac.cpp.
|
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.
|
overridevirtual |
Update the diagonal terms in the Jacobian by using the transient mask \( \alpha \).
rdt | Reciprocal of the time step [1/s] |
mask | Mask for transient terms: 1 if transient, 0 if algebraic. |
Reimplemented from SystemJacobian.
Definition at line 48 of file MultiJac.cpp.
|
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.
|
inline |
Definition at line 63 of file MultiJac.h.
|
inline |
Definition at line 67 of file MultiJac.h.
|
inline |
Definition at line 71 of file MultiJac.h.
|
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.
[in] | stateSize | length of the rhs and output vectors |
[in] | rhs_vector | right hand side vector used in linear system |
[out] | output | output vector for solution |
Reimplemented from SystemJacobian.
Definition at line 75 of file MultiJac.h.
|
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.
|
inline |
Return the transient mask.
Definition at line 85 of file MultiJac.h.
|
protected |
Residual evaluator for this Jacobian.
This is a pointer to the residual evaluator. This object isn't owned by this Jacobian object.
Definition at line 98 of file MultiJac.h.
|
protected |
Underlying matrix storage.
Definition at line 100 of file MultiJac.h.
|
protected |
Diagonal of the steady-state Jacobian.
Definition at line 101 of file MultiJac.h.
|
protected |
Transient mask for transient terms, 1 if transient, 0 if steady-state.
Definition at line 105 of file MultiJac.h.