Loading [MathJax]/extensions/tex2jax.js
Cantera  3.2.0a1
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules Pages
AdaptivePreconditioner.h
Go to the documentation of this file.
1/**
2 * @file AdaptivePreconditioner.h Declarations for the class
3 * AdaptivePreconditioner which is a child class of SystemJacobian
4 * for preconditioners used by sundials
5 */
6
7// This file is part of Cantera. See License.txt in the top-level directory or
8// at https://cantera.org/license.txt for license and copyright information.
9
10#ifndef ADAPTIVEPRECONDITIONER_H
11#define ADAPTIVEPRECONDITIONER_H
12
14
15namespace Cantera
16{
17
18//! AdaptivePreconditioner a preconditioner designed for use with large
19//! mechanisms that leverages sparse solvers. It does this by pruning
20//! the preconditioner by a threshold value. It also neglects pressure
21//! dependence and third body contributions in its formation and has a
22//! finite difference approximation for temperature.
24{
25public:
27
28 void initialize(size_t networkSize) override;
29 const string type() const override { return "Adaptive"; }
30 void setup() override; // deprecated
31 void factorize() override;
32 void solve(const size_t stateSize, double* rhs_vector, double* output) override;
33 void stateAdjustment(vector<double>& state) override;
34
35 int info() const override {
36 return static_cast<int>(m_solver.info());
37 }
38
39 //! Prune preconditioner elements
41
42 //! Get the threshold value for setting elements
43 double threshold() { return m_threshold; }
44
45 //! Get ILUT fill factor
46 double ilutFillFactor() { return m_fill_factor; }
47
48 //! Get ILUT drop tolerance
49 double ilutDropTol() { return m_drop_tol; }
50
51 //! Set the threshold value to compare elements against
52 //! @param threshold double value used in setting by threshold
53 void setThreshold(double threshold) {
55 m_prune_precon = (threshold <= 0) ? false : true;
56 }
57
58 //! Set drop tolerance for ILUT
59 //! @param droptol double value used in setting solver drop tolerance
60 void setIlutDropTol(double droptol) {
61 m_drop_tol = droptol;
62 m_solver.setDroptol(droptol);
63 }
64
65 //! Set the fill factor for ILUT solver
66 //! @param fillFactor fill in factor for ILUT solver
67 void setIlutFillFactor(int fillFactor) {
68 m_fill_factor = fillFactor;
69 m_solver.setFillfactor(fillFactor);
70 }
71
72protected:
73 //! ILUT fill factor
74 double m_fill_factor = 0;
75
76 //! ILUT drop tolerance
77 double m_drop_tol = 0;
78
79 //! Solver used in solving the linear system
80 Eigen::IncompleteLUT<double> m_solver;
81
82 //! Minimum value a non-diagonal element must be to be included in
83 //! the preconditioner
84 double m_threshold = 0.0;
85
86 //! Bool set whether to prune the matrix or not
87 double m_prune_precon = true;
88};
89
90}
91
92#endif
AdaptivePreconditioner a preconditioner designed for use with large mechanisms that leverages sparse ...
double threshold()
Get the threshold value for setting elements.
Eigen::IncompleteLUT< double > m_solver
Solver used in solving the linear system.
double ilutDropTol()
Get ILUT drop tolerance.
void setIlutDropTol(double droptol)
Set drop tolerance for ILUT.
double m_drop_tol
ILUT drop tolerance.
void solve(const size_t stateSize, double *rhs_vector, double *output) override
Solve a linear system using the system matrix M
double m_prune_precon
Bool set whether to prune the matrix or not.
double m_threshold
Minimum value a non-diagonal element must be to be included in the preconditioner.
void factorize() override
Factorize the system matrix.
void initialize(size_t networkSize) override
Called during setup for any processes that need to be completed prior to setup functions used in sund...
void setIlutFillFactor(int fillFactor)
Set the fill factor for ILUT solver.
void setup() override
Perform preconditioner specific post-reactor setup operations such as factorize.
void stateAdjustment(vector< double > &state) override
Adjust the state vector based on the preconditioner, e.g., Adaptive preconditioning uses a strictly p...
const string type() const override
Derived type, corresponding to names registered with SystemJacobianFactory.
void prunePreconditioner()
Prune preconditioner elements.
void setThreshold(double threshold)
Set the threshold value to compare elements against.
int info() const override
Get latest status of linear solver.
double ilutFillFactor()
Get ILUT fill factor.
System Jacobians that use Eigen sparse matrices for storage.
Namespace for the Cantera kernel.
Definition AnyMap.cpp:595