LLT.hpp
Go to the documentation of this file.
1 /*
2  * Copyright 2020 INRIA
3  */
4 
5 #ifndef __eigenpy_decomposition_llt_hpp__
6 #define __eigenpy_decomposition_llt_hpp__
7 
8 #include "eigenpy/eigenpy.hpp"
9 
10 #include <Eigen/Core>
11 #include <Eigen/Cholesky>
12 
14 
15 namespace eigenpy
16 {
17 
18  template<typename _MatrixType>
20  : public boost::python::def_visitor< LLTSolverVisitor<_MatrixType> >
21  {
22 
23  typedef _MatrixType MatrixType;
24  typedef typename MatrixType::Scalar Scalar;
25  typedef typename MatrixType::RealScalar RealScalar;
26  typedef Eigen::Matrix<Scalar,Eigen::Dynamic,1,MatrixType::Options> VectorType;
27  typedef Eigen::LLT<MatrixType> Solver;
28 
29  template<class PyClass>
30  void visit(PyClass& cl) const
31  {
32  namespace bp = boost::python;
33  cl
34  .def(bp::init<>("Default constructor"))
35  .def(bp::init<Eigen::DenseIndex>(bp::arg("size"),
36  "Default constructor with memory preallocation"))
37  .def(bp::init<MatrixType>(bp::arg("matrix"),
38  "Constructs a LLT factorization from a given matrix."))
39 
40  .def("matrixL",&matrixL,bp::arg("self"),
41  "Returns the lower triangular matrix L.")
42  .def("matrixU",&matrixU,bp::arg("self"),
43  "Returns the upper triangular matrix U.")
44  .def("matrixLLT",&Solver::matrixLLT,bp::arg("self"),
45  "Returns the LLT decomposition matrix.",
46  bp::return_internal_reference<>())
47 
48 #if EIGEN_VERSION_AT_LEAST(3,3,90)
49  .def("rankUpdate",(Solver& (Solver::*)(const VectorType &, const RealScalar &))&Solver::template rankUpdate<VectorType>,
50  bp::args("self","vector","sigma"), bp::return_self<>())
51 #else
52  .def("rankUpdate",(Solver (Solver::*)(const VectorType &, const RealScalar &))&Solver::template rankUpdate<VectorType>,
53  bp::args("self","vector","sigma"))
54 #endif
55 
56 #if EIGEN_VERSION_AT_LEAST(3,3,0)
57  .def("adjoint",&Solver::adjoint,bp::arg("self"),
58  "Returns the adjoint, that is, a reference to the decomposition itself as if the underlying matrix is self-adjoint.",
59  bp::return_self<>())
60 #endif
61 
62  .def("compute",(Solver & (Solver::*)(const Eigen::EigenBase<MatrixType> & matrix))&Solver::compute,
63  bp::args("self","matrix"),
64  "Computes the LLT of given matrix.",
65  bp::return_self<>())
66 
67  .def("info",&Solver::info,bp::arg("self"),
68  "NumericalIssue if the input contains INF or NaN values or overflow occured. Returns Success otherwise.")
69 #if EIGEN_VERSION_AT_LEAST(3,3,0)
70  .def("rcond",&Solver::rcond,bp::arg("self"),
71  "Returns an estimate of the reciprocal condition number of the matrix.")
72 #endif
73  .def("reconstructedMatrix",&Solver::reconstructedMatrix,bp::arg("self"),
74  "Returns the matrix represented by the decomposition, i.e., it returns the product: L L^*. This function is provided for debug purpose.")
75  .def("solve",&solve<VectorType>,bp::args("self","b"),
76  "Returns the solution x of A x = b using the current decomposition of A.")
77  ;
78  }
79 
80  static void expose()
81  {
82  static const std::string classname = "LLT" + scalar_name<Scalar>::shortname();
83  expose(classname);
84  }
85 
86  static void expose(const std::string & name)
87  {
88  namespace bp = boost::python;
89  bp::class_<Solver>(name.c_str(),
90  "Standard Cholesky decomposition (LL^T) of a matrix and associated features.\n\n"
91  "This class performs a LL^T Cholesky decomposition of a symmetric, positive definite matrix A such that A = LL^* = U^*U, where L is lower triangular.\n\n"
92  "While the Cholesky decomposition is particularly useful to solve selfadjoint problems like D^*D x = b, for that purpose, we recommend the Cholesky decomposition without square root which is more stable and even faster. Nevertheless, this standard Cholesky decomposition remains useful in many other situations like generalised eigen problems with hermitian matrices.",
93  bp::no_init)
94  .def(LLTSolverVisitor());
95  }
96 
97  private:
98 
99  static MatrixType matrixL(const Solver & self) { return self.matrixL(); }
100  static MatrixType matrixU(const Solver & self) { return self.matrixU(); }
101 
102  template<typename VectorType>
103  static VectorType solve(const Solver & self, const VectorType & vec)
104  {
105  return self.solve(vec);
106  }
107  };
108 
109 } // namespace eigenpy
110 
111 #endif // ifndef __eigenpy_decomposition_llt_hpp__
boost::python::object matrix()
Definition: bnpy.cpp:20
static void expose(const std::string &name)
Definition: LLT.hpp:86
void visit(PyClass &cl) const
Definition: LLT.hpp:30
MatrixType::RealScalar RealScalar
Definition: LLT.hpp:25
_MatrixType MatrixType
Definition: LLT.hpp:23
static std::string shortname()
Eigen::LLT< MatrixType > Solver
Definition: LLT.hpp:27
static void expose()
Definition: LLT.hpp:80
static MatrixType matrixL(const Solver &self)
Definition: LLT.hpp:99
static MatrixType matrixU(const Solver &self)
Definition: LLT.hpp:100
Eigen::Matrix< Scalar, Eigen::Dynamic, 1, MatrixType::Options > VectorType
Definition: LLT.hpp:26
static VectorType solve(const Solver &self, const VectorType &vec)
Definition: LLT.hpp:103
MatrixType::Scalar Scalar
Definition: LLT.hpp:24


eigenpy
Author(s): Justin Carpentier, Nicolas Mansard
autogenerated on Sat Apr 17 2021 02:37:59