Cantera  3.2.0a1
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules Pages
AdaptivePreconditioner.cpp
Go to the documentation of this file.
1//! @file AdaptivePreconditioner.cpp
2
3// This file is part of Cantera. See License.txt in the top-level directory or
4// at https://cantera.org/license.txt for license and copyright information.
5
8
9namespace Cantera
10{
11
12AdaptivePreconditioner::AdaptivePreconditioner()
13{
14 setPreconditionerSide("right");
15}
16
17void AdaptivePreconditioner::stateAdjustment(vector<double>& state) {
18 // Only keep positive composition based on given tol
19 for (size_t i = 0; i < state.size(); i++) {
20 state[i] = std::max(state[i], m_atol);
21 }
22}
23
24void AdaptivePreconditioner::initialize(size_t networkSize)
25{
27 // don't use legacy rate constants
29 // setting default ILUT parameters
30 if (m_drop_tol == 0) {
31 setIlutDropTol(1e-12);
32 }
33 if (m_drop_tol == 0) {
34 setIlutFillFactor(static_cast<int>(m_dim) / 2);
35 }
36 // update initialized status
37 m_init = true;
38}
39
41{
42 warn_deprecated("AdaptivePreconditioner::setup",
43 "To be removed after Cantera 3.2. Use updatePreconditioner() instead.");
45 factorize();
46}
47
49{
50 if (m_prune_precon) {
52 }
53 // compress sparse matrix structure
54 m_matrix.makeCompressed();
55 // analyze and factorize
56 m_solver.compute(m_matrix);
57 // check for errors
58 if (m_solver.info() != Eigen::Success) {
59 throw CanteraError("AdaptivePreconditioner::factorize",
60 "error code: {}", static_cast<int>(m_solver.info()));
61 }
62}
63
65{
66 for (int k=0; k<m_matrix.outerSize(); ++k) {
67 for (Eigen::SparseMatrix<double>::InnerIterator it(m_matrix, k); it;
68 ++it) {
69 if (std::abs(it.value()) < m_threshold && it.row() != it.col()) {
70 it.valueRef() = 0;
71 }
72 }
73 }
74}
75
76void AdaptivePreconditioner::solve(const size_t stateSize, double* rhs_vector, double*
77 output)
78{
79 // creating vectors in the form of Ax=b
80 Eigen::Map<Eigen::VectorXd> bVector(rhs_vector, stateSize);
81 Eigen::Map<Eigen::VectorXd> xVector(output, stateSize);
82 // solve for xVector
83 xVector = m_solver.solve(bVector);
84 if (m_solver.info() != Eigen::Success) {
85 throw CanteraError("AdaptivePreconditioner::solve",
86 "error code: {}", static_cast<int>(m_solver.info()));
87 }
88}
89
90}
Declarations for the class AdaptivePreconditioner which is a child class of SystemJacobian for precon...
Eigen::IncompleteLUT< double > m_solver
Solver used in solving the linear system.
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...
void prunePreconditioner()
Prune preconditioner elements.
Base class for exceptions thrown by Cantera classes.
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 updatePreconditioner() override
Transform Jacobian vector and write into preconditioner, P = (I - gamma * J)
Eigen::SparseMatrix< double > m_matrix
Container that is the sparse preconditioner.
size_t m_dim
Dimension of the system.
virtual void setPreconditionerSide(const string &preconSide)
For iterative solvers, set the side where the preconditioner is applied.
bool m_init
bool saying whether or not the system is initialized
double m_atol
Absolute tolerance of the ODE solver.
This file contains definitions for utility functions and text for modules, inputfiles and logging,...
void use_legacy_rate_constants(bool legacy)
Set definition used for rate constant calculation.
Definition global.cpp:102
Namespace for the Cantera kernel.
Definition AnyMap.cpp:595
void warn_deprecated(const string &source, const AnyBase &node, const string &message)
A deprecation warning for syntax in an input file.
Definition AnyMap.cpp:1997