SparseCholesky.h
Go to the documentation of this file.
1 // Copyright (C) 2016-2019 Yixuan Qiu <yixuan.qiu@cos.name>
2 //
3 // This Source Code Form is subject to the terms of the Mozilla
4 // Public License v. 2.0. If a copy of the MPL was not distributed
5 // with this file, You can obtain one at https://mozilla.org/MPL/2.0/.
6 
7 #ifndef SPARSE_CHOLESKY_H
8 #define SPARSE_CHOLESKY_H
9 
10 #include <Eigen/Core>
11 #include <Eigen/SparseCore>
12 #include <Eigen/SparseCholesky>
13 #include <stdexcept>
14 #include "../Util/CompInfo.h"
15 
16 namespace Spectra {
17 
26 template <typename Scalar, int Uplo = Eigen::Lower, int Flags = 0, typename StorageIndex = int>
28 {
29 private:
36 
37  const Index m_n;
39  int m_info; // status of the decomposition
40 
41 public:
49  SparseCholesky(ConstGenericSparseMatrix& mat) :
50  m_n(mat.rows())
51  {
52  if (mat.rows() != mat.cols())
53  throw std::invalid_argument("SparseCholesky: matrix must be square");
54 
55  m_decomp.compute(mat);
56  m_info = (m_decomp.info() == Eigen::Success) ?
57  SUCCESSFUL :
59  }
60 
64  Index rows() const { return m_n; }
68  Index cols() const { return m_n; }
69 
74  int info() const { return m_info; }
75 
82  // y_out = inv(L) * x_in
83  void lower_triangular_solve(const Scalar* x_in, Scalar* y_out) const
84  {
85  MapConstVec x(x_in, m_n);
86  MapVec y(y_out, m_n);
87  y.noalias() = m_decomp.permutationP() * x;
88  m_decomp.matrixL().solveInPlace(y);
89  }
90 
97  // y_out = inv(L') * x_in
98  void upper_triangular_solve(const Scalar* x_in, Scalar* y_out) const
99  {
100  MapConstVec x(x_in, m_n);
101  MapVec y(y_out, m_n);
102  y.noalias() = m_decomp.matrixU().solve(x);
103  y = m_decomp.permutationPinv() * y;
104  }
105 };
106 
107 } // namespace Spectra
108 
109 #endif // SPARSE_CHOLESKY_H
SCALAR Scalar
Definition: bench_gemm.cpp:46
const MatrixU matrixU() const
Scalar * y
Eigen::SparseMatrix< Scalar, Flags, StorageIndex > SparseMatrix
void lower_triangular_solve(const Scalar *x_in, Scalar *y_out) const
A versatible sparse matrix representation.
Definition: SparseMatrix.h:96
SimplicialLLT & compute(const MatrixType &matrix)
A matrix or vector expression mapping an existing array of data.
Definition: Map.h:94
void upper_triangular_solve(const Scalar *x_in, Scalar *y_out) const
const Eigen::Ref< const SparseMatrix > ConstGenericSparseMatrix
Computation was successful.
Definition: CompInfo.h:19
Eigen::Map< const Vector > MapConstVec
EIGEN_DEFAULT_DENSE_INDEX_TYPE Index
The Index type as used for the API.
Definition: Meta.h:74
SparseCholesky(ConstGenericSparseMatrix &mat)
Eigen::SimplicialLLT< SparseMatrix, Uplo > m_decomp
Eigen::Matrix< Scalar, Eigen::Dynamic, 1 > Vector
ComputationInfo info() const
Reports whether previous computation was successful.
const PermutationMatrix< Dynamic, Dynamic, StorageIndex > & permutationPinv() const
const PermutationMatrix< Dynamic, Dynamic, StorageIndex > & permutationP() const
Eigen::Map< Vector > MapVec
A direct sparse LLT Cholesky factorizations.
set noclip points set clip one set noclip two set bar set border lt lw set xdata set ydata set zdata set x2data set y2data set boxwidth set dummy x
const MatrixL matrixL() const


gtsam
Author(s):
autogenerated on Tue Jul 4 2023 02:35:55