|
| 1 | +/* |
| 2 | + * Copyright 2020 INRIA |
| 3 | + */ |
| 4 | + |
| 5 | +#ifndef __eigenpy_decomposition_ldlt_hpp__ |
| 6 | +#define __eigenpy_decomposition_ldlt_hpp__ |
| 7 | + |
| 8 | +#include <boost/python.hpp> |
| 9 | +#include <Eigen/Core> |
| 10 | +#include <Eigen/Cholesky> |
| 11 | + |
| 12 | +#include "eigenpy/utils/scalar-name.hpp" |
| 13 | + |
| 14 | +namespace eigenpy |
| 15 | +{ |
| 16 | + |
| 17 | + template<typename _MatrixType> |
| 18 | + struct LDLTSolverVisitor |
| 19 | + : public boost::python::def_visitor< LDLTSolverVisitor<_MatrixType> > |
| 20 | + { |
| 21 | + |
| 22 | + typedef _MatrixType MatrixType; |
| 23 | + typedef typename MatrixType::Scalar Scalar; |
| 24 | + typedef typename MatrixType::RealScalar RealScalar; |
| 25 | + typedef Eigen::Matrix<Scalar,Eigen::Dynamic,1,MatrixType::Options> VectorType; |
| 26 | + typedef Eigen::LDLT<MatrixType> Solver; |
| 27 | + |
| 28 | + template<class PyClass> |
| 29 | + void visit(PyClass& cl) const |
| 30 | + { |
| 31 | + namespace bp = boost::python; |
| 32 | + cl |
| 33 | + .def(bp::init<>("Default constructor")) |
| 34 | + .def(bp::init<Eigen::DenseIndex>(bp::arg("size"), |
| 35 | + "Default constructor with memory preallocation")) |
| 36 | + .def(bp::init<MatrixType>(bp::arg("matrix"), |
| 37 | + "Constructs a LDLT factorization from a given matrix.")) |
| 38 | + |
| 39 | + .def("isNegative",&Solver::isNegative,bp::arg("self"), |
| 40 | + "Returns true if the matrix is negative (semidefinite).") |
| 41 | + .def("isPositive",&Solver::isPositive,bp::arg("self"), |
| 42 | + "Returns true if the matrix is positive (semidefinite).") |
| 43 | + |
| 44 | + .def("matrixL",&matrixL,bp::arg("self"), |
| 45 | + "Returns the lower triangular matrix L.") |
| 46 | + .def("matrixU",&matrixU,bp::arg("self"), |
| 47 | + "Returns the upper triangular matrix U.") |
| 48 | + .def("vectorD",&vectorD,bp::arg("self"), |
| 49 | + "Returns the coefficients of the diagonal matrix D.") |
| 50 | + .def("transpositionsP",&transpositionsP,bp::arg("self"), |
| 51 | + "Returns the permutation matrix P.") |
| 52 | + |
| 53 | + .def("matrixLDLT",&Solver::matrixLDLT,bp::arg("self"), |
| 54 | + "Returns the LDLT decomposition matrix.", |
| 55 | + bp::return_value_policy<bp::return_by_value>()) |
| 56 | + |
| 57 | + .def("rankUpdate",(Solver & (Solver::*)(const Eigen::MatrixBase<VectorType> &, const RealScalar &))&Solver::template rankUpdate<VectorType>, |
| 58 | + bp::args("self","vector","sigma"), |
| 59 | + bp::return_value_policy<bp::reference_existing_object>() ) |
| 60 | + |
| 61 | +#if EIGEN_VERSION_AT_LEAST(3,3,0) |
| 62 | + .def("adjoint",&Solver::adjoint,bp::arg("self"), |
| 63 | + "Returns the adjoint, that is, a reference to the decomposition itself as if the underlying matrix is self-adjoint.", |
| 64 | + bp::return_value_policy<bp::reference_existing_object>()) |
| 65 | +#endif |
| 66 | + |
| 67 | + .def("compute",(Solver & (Solver::*)(const Eigen::EigenBase<MatrixType> & matrix))&Solver::compute, |
| 68 | + bp::args("self","matrix"), |
| 69 | + "Computes the LDLT of given matrix.", |
| 70 | + bp::return_value_policy<bp::reference_existing_object>()) |
| 71 | + |
| 72 | + .def("info",&Solver::info,bp::arg("self"), |
| 73 | + "NumericalIssue if the input contains INF or NaN values or overflow occured. Returns Success otherwise.") |
| 74 | +#if EIGEN_VERSION_AT_LEAST(3,3,0) |
| 75 | + .def("rcond",&Solver::rcond,bp::arg("self"), |
| 76 | + "Returns an estimate of the reciprocal condition number of the matrix.") |
| 77 | +#endif |
| 78 | + .def("reconstructedMatrix",&Solver::reconstructedMatrix,bp::arg("self"), |
| 79 | + "Returns the matrix represented by the decomposition, i.e., it returns the product: L L^*. This function is provided for debug purpose.") |
| 80 | + .def("solve",&solve<VectorType>,bp::args("self","b"), |
| 81 | + "Returns the solution x of A x = b using the current decomposition of A.") |
| 82 | + |
| 83 | + .def("setZero",&Solver::setZero,bp::arg("self"), |
| 84 | + "Clear any existing decomposition.") |
| 85 | + ; |
| 86 | + } |
| 87 | + |
| 88 | + static void expose() |
| 89 | + { |
| 90 | + static const std::string classname = "LDLT" + scalar_name<Scalar>::shortname(); |
| 91 | + expose(classname); |
| 92 | + } |
| 93 | + |
| 94 | + static void expose(const std::string & name) |
| 95 | + { |
| 96 | + namespace bp = boost::python; |
| 97 | + bp::class_<Solver>(name.c_str(), |
| 98 | + "Robust Cholesky decomposition of a matrix with pivoting.\n\n" |
| 99 | + "Perform a robust Cholesky decomposition of a positive semidefinite or negative semidefinite matrix $ A $ such that $ A = P^TLDL^*P $, where P is a permutation matrix, L is lower triangular with a unit diagonal and D is a diagonal matrix.\n\n" |
| 100 | + "The decomposition uses pivoting to ensure stability, so that L will have zeros in the bottom right rank(A) - n submatrix. Avoiding the square root on D also stabilizes the computation.", |
| 101 | + bp::no_init) |
| 102 | + .def(LDLTSolverVisitor()); |
| 103 | + } |
| 104 | + |
| 105 | + private: |
| 106 | + |
| 107 | + static MatrixType matrixL(const Solver & self) { return self.matrixL(); } |
| 108 | + static MatrixType matrixU(const Solver & self) { return self.matrixU(); } |
| 109 | + static VectorType vectorD(const Solver & self) { return self.vectorD(); } |
| 110 | + |
| 111 | + static MatrixType transpositionsP(const Solver & self) |
| 112 | + { |
| 113 | + return self.transpositionsP() * MatrixType::Identity(self.matrixL().rows(), |
| 114 | + self.matrixL().rows()); |
| 115 | + } |
| 116 | + |
| 117 | + template<typename VectorType> |
| 118 | + static VectorType solve(const Solver & self, const VectorType & vec) |
| 119 | + { |
| 120 | + return self.solve(vec); |
| 121 | + } |
| 122 | + }; |
| 123 | + |
| 124 | +} // namespace eigenpy |
| 125 | + |
| 126 | +#endif // ifndef __eigenpy_decomposition_ldlt_hpp__ |
0 commit comments