Cantera
3.1.0a1
|
AdaptivePreconditioner a preconditioner designed for use with large mechanisms that leverages sparse solvers. More...
#include <AdaptivePreconditioner.h>
AdaptivePreconditioner a preconditioner designed for use with large mechanisms that leverages sparse solvers.
It does this by pruning the preconditioner by a threshold value. It also neglects pressure dependence and third body contributions in its formation and has a finite difference approximation for temperature.
Definition at line 26 of file AdaptivePreconditioner.h.
Public Member Functions | |
void | initialize (size_t networkSize) override |
Called during setup for any processes that need to be completed prior to setup functions used in sundials. More... | |
void | reset () override |
Reset preconditioner parameters as needed. More... | |
void | setup () override |
Perform preconditioner specific post-reactor setup operations such as factorize. More... | |
void | solve (const size_t stateSize, double *rhs_vector, double *output) override |
Solve a linear system Ax=b where A is the preconditioner. More... | |
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. More... | |
void | stateAdjustment (vector< double > &state) override |
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. More... | |
void | updatePreconditioner () override |
Transform Jacobian vector and write into preconditioner, P = (I - gamma * J) More... | |
void | prunePreconditioner () |
Prune preconditioner elements. More... | |
Eigen::SparseMatrix< double > | jacobian () |
Return semi-analytical Jacobian of an AdaptivePreconditioner object. More... | |
Eigen::SparseMatrix< double > | matrix () |
Return the internal preconditioner matrix. More... | |
double | threshold () |
Get the threshold value for setting elements. More... | |
double | ilutFillFactor () |
Get ILUT fill factor. More... | |
double | ilutDropTol () |
Get ILUT drop tolerance. More... | |
void | setThreshold (double threshold) |
Set the threshold value to compare elements against. More... | |
void | setIlutDropTol (double droptol) |
Set drop tolerance for ILUT. More... | |
void | setIlutFillFactor (int fillFactor) |
Set the fill factor for ILUT solver. More... | |
void | printPreconditioner () override |
Print preconditioner contents. More... | |
void | printJacobian () |
Print jacobian contents. More... | |
Public Member Functions inherited from PreconditionerBase | |
string | preconditionerSide () const |
Get preconditioner application side for CVODES. More... | |
virtual void | setPreconditionerSide (const string &preconSide) |
virtual void | setGamma (double gamma) |
Set gamma used in preconditioning. More... | |
virtual double | gamma () |
Get gamma used in preconditioning. More... | |
virtual void | setAbsoluteTolerance (double atol) |
Set the absolute tolerance in the solver outside of the network initialization. More... | |
Protected Attributes | |
double | m_fill_factor = 0 |
ILUT fill factor. More... | |
double | m_drop_tol = 0 |
ILUT drop tolerance. More... | |
vector< Eigen::Triplet< double > > | m_jac_trips |
Vector of triples representing the jacobian used in preconditioning. More... | |
Eigen::SparseMatrix< double > | m_identity |
Storage of appropriately sized identity matrix for making the preconditioner. More... | |
Eigen::SparseMatrix< double > | m_precon_matrix |
Container that is the sparse preconditioner. More... | |
Eigen::IncompleteLUT< double > | m_solver |
Solver used in solving the linear system. More... | |
double | m_threshold = 0.0 |
Minimum value a non-diagonal element must be to be included in the preconditioner. More... | |
double | m_prune_precon = true |
Bool set whether to prune the matrix or not. More... | |
Protected Attributes inherited from PreconditionerBase | |
size_t | m_dim |
Dimension of the preconditioner. More... | |
double | m_gamma = 1.0 |
gamma value used in M = I - gamma*J More... | |
bool | m_init = false |
bool saying whether or not the preconditioner is initialized More... | |
double | m_atol = 0 |
Absolute tolerance of the ODE solver. More... | |
string | m_precon_side = "none" |
|
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 PreconditionerBase.
Definition at line 29 of file AdaptivePreconditioner.cpp.
|
inlineoverridevirtual |
Reset preconditioner parameters as needed.
Reimplemented from PreconditionerBase.
Definition at line 33 of file AdaptivePreconditioner.h.
|
overridevirtual |
Perform preconditioner specific post-reactor setup operations such as factorize.
Reimplemented from PreconditionerBase.
Definition at line 57 of file AdaptivePreconditioner.cpp.
|
overridevirtual |
Solve a linear system Ax=b where A is the preconditioner.
[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 PreconditionerBase.
Definition at line 96 of file AdaptivePreconditioner.cpp.
|
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 PreconditionerBase.
Definition at line 17 of file AdaptivePreconditioner.cpp.
|
overridevirtual |
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.
state | a vector containing the state to be updated |
Reimplemented from PreconditionerBase.
Definition at line 22 of file AdaptivePreconditioner.cpp.
|
overridevirtual |
Transform Jacobian vector and write into preconditioner, P = (I - gamma * J)
Reimplemented from PreconditionerBase.
Definition at line 72 of file AdaptivePreconditioner.cpp.
void prunePreconditioner | ( | ) |
Prune preconditioner elements.
Definition at line 84 of file AdaptivePreconditioner.cpp.
|
inline |
Return the internal preconditioner matrix.
Definition at line 60 of file AdaptivePreconditioner.h.
|
inline |
Get the threshold value for setting elements.
Definition at line 66 of file AdaptivePreconditioner.h.
|
inline |
Get ILUT fill factor.
Definition at line 69 of file AdaptivePreconditioner.h.
|
inline |
Get ILUT drop tolerance.
Definition at line 72 of file AdaptivePreconditioner.h.
|
inline |
Set the threshold value to compare elements against.
threshold | double value used in setting by threshold |
Definition at line 76 of file AdaptivePreconditioner.h.
|
inline |
Set drop tolerance for ILUT.
droptol | double value used in setting solver drop tolerance |
Definition at line 83 of file AdaptivePreconditioner.h.
|
inline |
Set the fill factor for ILUT solver.
fillFactor | fill in factor for ILUT solver |
Definition at line 90 of file AdaptivePreconditioner.h.
|
overridevirtual |
Print preconditioner contents.
Reimplemented from PreconditionerBase.
Definition at line 110 of file AdaptivePreconditioner.cpp.
void printJacobian | ( | ) |
Print jacobian contents.
Definition at line 117 of file AdaptivePreconditioner.cpp.
|
protected |
ILUT fill factor.
Definition at line 103 of file AdaptivePreconditioner.h.
|
protected |
ILUT drop tolerance.
Definition at line 106 of file AdaptivePreconditioner.h.
|
protected |
Vector of triples representing the jacobian used in preconditioning.
Definition at line 109 of file AdaptivePreconditioner.h.
|
protected |
Storage of appropriately sized identity matrix for making the preconditioner.
Definition at line 112 of file AdaptivePreconditioner.h.
|
protected |
Container that is the sparse preconditioner.
Definition at line 115 of file AdaptivePreconditioner.h.
|
protected |
Solver used in solving the linear system.
Definition at line 118 of file AdaptivePreconditioner.h.
|
protected |
Minimum value a non-diagonal element must be to be included in the preconditioner.
Definition at line 122 of file AdaptivePreconditioner.h.
|
protected |
Bool set whether to prune the matrix or not.
Definition at line 125 of file AdaptivePreconditioner.h.