Loading [MathJax]/extensions/tex2jax.js
Cantera  3.2.0a1
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules Pages
EigenSparseJacobian.cpp
Go to the documentation of this file.
1//! @file EigenSparseJacobian.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
12void EigenSparseJacobian::initialize(size_t networkSize)
13{
14 // reset arrays in case of re-initialization
15 m_jac_trips.clear();
16 // set dimensions of preconditioner from network
17 m_dim = networkSize;
18 // reserve some space for vectors making up SparseMatrix
19 m_jac_trips.reserve(3 * networkSize);
20 // reserve space for preconditioner
21 m_matrix.resize(m_dim, m_dim);
22 // creating sparse identity matrix
23 m_identity.resize(m_dim, m_dim);
24 m_identity.setIdentity();
25 m_identity.makeCompressed();
26 // update initialized status
27 m_init = true;
28}
29
31{
32 m_matrix.setZero();
33 m_jac_trips.clear();
34}
35
36void EigenSparseJacobian::setValue(size_t row, size_t col, double value)
37{
38 m_jac_trips.emplace_back(static_cast<int>(row), static_cast<int>(col), value);
39}
40
42{
43 // set precon to jacobian
44 m_matrix.setFromTriplets(m_jac_trips.begin(), m_jac_trips.end());
45 // make into preconditioner as P = (I - gamma * J_bar)
47 // prune by threshold if desired
48 factorize();
49}
50
51void EigenSparseJacobian::updateTransient(double rdt, int* mask)
52{
53 // set matrix to steady Jacobian
54 m_matrix.setFromTriplets(m_jac_trips.begin(), m_jac_trips.end());
55 // update transient diagonal terms
56 Eigen::VectorXd diag = Eigen::Map<Eigen::VectorXi>(mask, m_dim).cast<double>();
57 m_matrix -= rdt * diag.matrix().asDiagonal();
58 factorize();
59}
60
61Eigen::SparseMatrix<double> EigenSparseJacobian::jacobian()
62{
63 Eigen::SparseMatrix<double> jacobian_mat(m_dim, m_dim);
64 jacobian_mat.setFromTriplets(m_jac_trips.begin(), m_jac_trips.end());
65 return jacobian_mat;
66}
67
69 std::stringstream ss;
70 Eigen::IOFormat HeavyFmt(Eigen::FullPrecision, 0, ", ", ";\n", "[", "]", "[", "]");
71 ss << Eigen::MatrixXd(m_matrix).format(HeavyFmt);
72 writelog(ss.str());
73}
74
76 std::stringstream ss;
77 Eigen::SparseMatrix<double> jacobian(m_dim, m_dim);
78 jacobian.setFromTriplets(m_jac_trips.begin(), m_jac_trips.end());
79 ss << Eigen::MatrixXd(jacobian);
80 writelog(ss.str());
81}
82
83}
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 networkSize) override
Called during setup for any processes that need to be completed prior to setup functions used in sund...
vector< Eigen::Triplet< double > > m_jac_trips
Vector of triples representing the jacobian used in preconditioning.
Eigen::SparseMatrix< double > m_identity
Storage of appropriately sized identity matrix for making the preconditioner.
void printPreconditioner() override
Print preconditioner contents.
void updatePreconditioner() override
Transform Jacobian vector and write into preconditioner, P = (I - gamma * J)
void updateTransient(double rdt, int *mask) override
Update the diagonal terms in the Jacobian by using the transient mask .
void printJacobian()
Print jacobian contents.
void reset() override
Reset parameters as needed.
Eigen::SparseMatrix< double > m_matrix
Container that is the sparse preconditioner.
size_t m_dim
Dimension of the system.
double m_gamma
gamma value used in M = I - gamma*J
virtual void factorize()
Factorize the system matrix.
bool m_init
bool saying whether or not the system is initialized
This file contains definitions for utility functions and text for modules, inputfiles and logging,...
Eigen::SparseMatrix< double > jacobian()
Return underlying Jacobian matrix.
void writelog(const string &fmt, const Args &... args)
Write a formatted message to the screen.
Definition global.h:171
Namespace for the Cantera kernel.
Definition AnyMap.cpp:595