Skip to content

Commit 9133ed0

Browse files
committed
sparse: add sparse Cholesky decompositions with Cholmod
1 parent 1451c72 commit 9133ed0

File tree

9 files changed

+386
-4
lines changed

9 files changed

+386
-4
lines changed

CMakeLists.txt

Lines changed: 28 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -109,6 +109,7 @@ if(BUILD_WITH_CHOLMOD_SUPPORT)
109109
STATUS
110110
"Build with CHOLDOD support (LGPL). See CHOLMOD/Doc/License.txt for further details."
111111
)
112+
add_definitions(-DEIGENPY_WITH_CHOLMOD_SUPPORT)
112113
endif(BUILD_WITH_CHOLMOD_SUPPORT)
113114

114115
# ----------------------------------------------------
@@ -131,11 +132,26 @@ set(${PROJECT_NAME}_SOLVERS_HEADERS
131132

132133
set(${PROJECT_NAME}_EIGEN_HEADERS include/eigenpy/eigen/EigenBase.hpp)
133134

134-
set(${PROJECT_NAME}_DECOMPOSITIONS_HEADERS
135-
include/eigenpy/decompositions/decompositions.hpp
135+
set(${PROJECT_NAME}_DECOMPOSITIONS_SPARSE_CHOLMOD_HEADERS
136+
include/eigenpy/decompositions/sparse/cholmod/CholmodBase.hpp
137+
include/eigenpy/decompositions/sparse/cholmod/CholmodDecomposition.hpp
138+
include/eigenpy/decompositions/sparse/cholmod/CholmodSimplicialLDLT.hpp
139+
include/eigenpy/decompositions/sparse/cholmod/CholmodSimplicialLLT.hpp
140+
include/eigenpy/decompositions/sparse/cholmod/CholmodSupernodalLLT.hpp)
141+
142+
set(${PROJECT_NAME}_DECOMPOSITIONS_SPARSE_HEADERS
136143
include/eigenpy/decompositions/sparse/LLT.hpp
137144
include/eigenpy/decompositions/sparse/LDLT.hpp
138-
include/eigenpy/decompositions/sparse/SimplicialCholesky.hpp
145+
include/eigenpy/decompositions/sparse/SimplicialCholesky.hpp)
146+
147+
if(BUILD_WITH_CHOLMOD_SUPPORT)
148+
list(APPEND ${PROJECT_NAME}_DECOMPOSITIONS_SPARSE_HEADERS
149+
${${PROJECT_NAME}_DECOMPOSITIONS_SPARSE_CHOLMOD_HEADERS})
150+
endif(BUILD_WITH_CHOLMOD_SUPPORT)
151+
152+
set(${PROJECT_NAME}_DECOMPOSITIONS_HEADERS
153+
${${PROJECT_NAME}_DECOMPOSITIONS_SPARSE_HEADERS}
154+
include/eigenpy/decompositions/decompositions.hpp
139155
include/eigenpy/decompositions/EigenSolver.hpp
140156
include/eigenpy/decompositions/PermutationMatrix.hpp
141157
include/eigenpy/decompositions/LDLT.hpp
@@ -207,6 +223,11 @@ set(${PROJECT_NAME}_SOLVERS_SOURCES src/solvers/preconditioners.cpp
207223
set(${PROJECT_NAME}_DECOMPOSITIONS_SOURCES
208224
src/decompositions/decompositions.cpp)
209225

226+
if(BUILD_WITH_CHOLMOD_SUPPORT)
227+
list(APPEND ${PROJECT_NAME}_DECOMPOSITIONS_SOURCES
228+
src/decompositions/cholmod.cpp)
229+
endif(BUILD_WITH_CHOLMOD_SUPPORT)
230+
210231
set(${PROJECT_NAME}_SOURCES
211232
${${PROJECT_NAME}_SOLVERS_SOURCES}
212233
${${PROJECT_NAME}_DECOMPOSITIONS_SOURCES}
@@ -259,8 +280,11 @@ modernize_target_link_libraries(
259280
${PYTHON_INCLUDE_DIR})
260281

261282
if(BUILD_WITH_CHOLMOD_SUPPORT)
283+
if(TARGET CHOLMOD::CHOLMOD)
284+
message(STATUS "CHOLMOD_FOUND: ${CHOLMOD_FOUND}")
285+
endif()
262286
modernize_target_link_libraries(${PROJECT_NAME} SCOPE PUBLIC TARGETS
263-
SuiteSparse::CHOLDMOD)
287+
CHOLMOD::CHOLMOD)
264288
endif(BUILD_WITH_CHOLMOD_SUPPORT)
265289
if(SUFFIX_SO_VERSION)
266290
set_target_properties(${PROJECT_NAME} PROPERTIES SOVERSION ${PROJECT_VERSION})

include/eigenpy/decompositions/decompositions.hpp

Lines changed: 5 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -9,6 +9,11 @@
99

1010
namespace eigenpy {
1111
void EIGENPY_DLLAPI exposeDecompositions();
12+
13+
#ifdef EIGENPY_WITH_CHOLMOD_SUPPORT
14+
void EIGENPY_DLLAPI exposeCholmod();
15+
#endif
16+
1217
} // namespace eigenpy
1318

1419
#endif // define __eigenpy_decompositions_decompositions_hpp__
Lines changed: 96 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,96 @@
1+
/*
2+
* Copyright 2024 INRIA
3+
*/
4+
5+
#ifndef __eigenpy_decomposition_sparse_cholmod_cholmod_base_hpp__
6+
#define __eigenpy_decomposition_sparse_cholmod_cholmod_base_hpp__
7+
8+
#include "eigenpy/eigenpy.hpp"
9+
#include "eigenpy/eigen/EigenBase.hpp"
10+
11+
#include <Eigen/CholmodSupport>
12+
13+
namespace eigenpy {
14+
15+
template <typename CholdmodDerived>
16+
struct CholmodBaseVisitor
17+
: public boost::python::def_visitor<CholmodBaseVisitor<CholdmodDerived> > {
18+
typedef CholdmodDerived Solver;
19+
20+
typedef typename CholdmodDerived::MatrixType MatrixType;
21+
typedef typename MatrixType::Scalar Scalar;
22+
typedef typename MatrixType::RealScalar RealScalar;
23+
typedef MatrixType CholMatrixType;
24+
typedef typename MatrixType::StorageIndex StorageIndex;
25+
26+
typedef Eigen::Matrix<Scalar, Eigen::Dynamic, 1, MatrixType::Options>
27+
DenseVectorXs;
28+
typedef Eigen::Matrix<Scalar, Eigen::Dynamic, Eigen::Dynamic,
29+
MatrixType::Options>
30+
DenseMatrixXs;
31+
32+
template <class PyClass>
33+
void visit(PyClass &cl) const {
34+
cl.def("analyzePattern", &Solver::analyzePattern,
35+
bp::args("self", "matrix"),
36+
"Performs a symbolic decomposition on the sparcity of matrix.\n"
37+
"This function is particularly useful when solving for several "
38+
"problems having the same structure.")
39+
40+
.def(EigenBaseVisitor<Solver>())
41+
42+
.def("compute",
43+
(Solver & (Solver::*)(const MatrixType &matrix)) & Solver::compute,
44+
bp::args("self", "matrix"),
45+
"Computes the sparse Cholesky decomposition of a given matrix.",
46+
bp::return_self<>())
47+
48+
.def("determinant", &Solver::determinant, bp::arg("self"),
49+
"Returns the determinant of the underlying matrix from the "
50+
"current factorization.")
51+
52+
.def("factorize", &Solver::factorize, bp::args("self", "matrix"),
53+
"Performs a numeric decomposition of a given matrix.\n"
54+
"The given matrix must has the same sparcity than the matrix on "
55+
"which the symbolic decomposition has been performed.\n"
56+
"See also analyzePattern().")
57+
58+
.def("info", &Solver::info, bp::arg("self"),
59+
"NumericalIssue if the input contains INF or NaN values or "
60+
"overflow occured. Returns Success otherwise.")
61+
62+
.def("logDeterminant", &Solver::logDeterminant, bp::arg("self"),
63+
"Returns the log determinant of the underlying matrix from the "
64+
"current factorization.")
65+
66+
.def("setShift", &Solver::setShift, (bp::args("self", "offset")),
67+
"Sets the shift parameters that will be used to adjust the "
68+
"diagonal coefficients during the numerical factorization.\n"
69+
"During the numerical factorization, the diagonal coefficients "
70+
"are transformed by the following linear model: d_ii = offset + "
71+
"d_ii.\n"
72+
"The default is the identity transformation with offset=0.",
73+
bp::return_self<>())
74+
75+
.def("solve", &solve<DenseVectorXs>, bp::args("self", "b"),
76+
"Returns the solution x of A x = b using the current "
77+
"decomposition of A.")
78+
.def("solve", &solve<DenseMatrixXs>, bp::args("self", "B"),
79+
"Returns the solution X of A X = B using the current "
80+
"decomposition of A where B is a right hand side matrix.")
81+
82+
.def("solve", &solve<MatrixType>, bp::args("self", "B"),
83+
"Returns the solution X of A X = B using the current "
84+
"decomposition of A where B is a right hand side matrix.");
85+
}
86+
87+
private:
88+
template <typename MatrixOrVector>
89+
static MatrixOrVector solve(const Solver &self, const MatrixOrVector &vec) {
90+
return self.solve(vec);
91+
}
92+
};
93+
94+
} // namespace eigenpy
95+
96+
#endif // ifndef __eigenpy_decomposition_sparse_cholmod_cholmod_base_hpp__
Lines changed: 32 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,32 @@
1+
/*
2+
* Copyright 2024 INRIA
3+
*/
4+
5+
#ifndef __eigenpy_decomposition_sparse_cholmod_cholmod_decomposition_hpp__
6+
#define __eigenpy_decomposition_sparse_cholmod_cholmod_decomposition_hpp__
7+
8+
#include "eigenpy/eigenpy.hpp"
9+
#include "eigenpy/decompositions/sparse/cholmod/CholmodBase.hpp"
10+
11+
namespace eigenpy {
12+
13+
template <typename CholdmodDerived>
14+
struct CholmodDecompositionVisitor
15+
: public boost::python::def_visitor<
16+
CholmodDecompositionVisitor<CholdmodDerived> > {
17+
typedef CholdmodDerived Solver;
18+
19+
template <class PyClass>
20+
void visit(PyClass &cl) const {
21+
cl
22+
23+
.def(CholmodBaseVisitor<Solver>())
24+
.def("setMode", &Solver::setMode, bp::args("self", "mode"),
25+
"Set the mode for the Cholesky decomposition.");
26+
}
27+
};
28+
29+
} // namespace eigenpy
30+
31+
#endif // ifndef
32+
// __eigenpy_decomposition_sparse_cholmod_cholmod_decomposition_hpp__
Lines changed: 63 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,63 @@
1+
/*
2+
* Copyright 2024 INRIA
3+
*/
4+
5+
#ifndef __eigenpy_decomposition_sparse_cholmod_cholmod_simplicial_ldlt_hpp__
6+
#define __eigenpy_decomposition_sparse_cholmod_cholmod_simplicial_ldlt_hpp__
7+
8+
#include "eigenpy/eigenpy.hpp"
9+
#include "eigenpy/decompositions/sparse/cholmod/CholmodDecomposition.hpp"
10+
#include "eigenpy/utils/scalar-name.hpp"
11+
12+
namespace eigenpy {
13+
14+
template <typename MatrixType_, int UpLo_ = Eigen::Lower>
15+
struct CholmodSimplicialLDLTVisitor
16+
: public boost::python::def_visitor<
17+
CholmodSimplicialLDLTVisitor<MatrixType_, UpLo_> > {
18+
typedef MatrixType_ MatrixType;
19+
typedef typename MatrixType::Scalar Scalar;
20+
typedef typename MatrixType::RealScalar RealScalar;
21+
22+
typedef Eigen::CholmodSimplicialLDLT<MatrixType_, UpLo_> Solver;
23+
24+
template <class PyClass>
25+
void visit(PyClass &cl) const {
26+
cl
27+
28+
.def(CholmodBaseVisitor<Solver>())
29+
.def(bp::init<>(bp::arg("self"), "Default constructor"))
30+
.def(bp::init<MatrixType>(bp::args("self", "matrix"),
31+
"Constructs and performs the LDLT "
32+
"factorization from a given matrix."))
33+
34+
;
35+
}
36+
37+
static void expose() {
38+
static const std::string classname =
39+
"CholmodSimplicialLDLT_" + scalar_name<Scalar>::shortname();
40+
expose(classname);
41+
}
42+
43+
static void expose(const std::string &name) {
44+
bp::class_<Solver, boost::noncopyable>(
45+
name.c_str(),
46+
"A simplicial direct Cholesky (LDLT) factorization and solver based on "
47+
"Cholmod.\n\n"
48+
"This class allows to solve for A.X = B sparse linear problems via a "
49+
"simplicial LL^T Cholesky factorization using the Cholmod library."
50+
"This simplicial variant is equivalent to Eigen's built-in "
51+
"SimplicialLDLT class."
52+
"Therefore, it has little practical interest. The sparse matrix A must "
53+
"be selfadjoint and positive definite."
54+
"The vectors or matrices X and B can be either dense or sparse.",
55+
bp::no_init)
56+
.def(CholmodSimplicialLDLTVisitor());
57+
}
58+
};
59+
60+
} // namespace eigenpy
61+
62+
#endif // ifndef
63+
// __eigenpy_decomposition_sparse_cholmod_cholmod_simplicial_ldlt_hpp__
Lines changed: 63 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,63 @@
1+
/*
2+
* Copyright 2024 INRIA
3+
*/
4+
5+
#ifndef __eigenpy_decomposition_sparse_cholmod_cholmod_simplicial_llt_hpp__
6+
#define __eigenpy_decomposition_sparse_cholmod_cholmod_simplicial_llt_hpp__
7+
8+
#include "eigenpy/eigenpy.hpp"
9+
#include "eigenpy/decompositions/sparse/cholmod/CholmodDecomposition.hpp"
10+
#include "eigenpy/utils/scalar-name.hpp"
11+
12+
namespace eigenpy {
13+
14+
template <typename MatrixType_, int UpLo_ = Eigen::Lower>
15+
struct CholmodSimplicialLLTVisitor
16+
: public boost::python::def_visitor<
17+
CholmodSimplicialLLTVisitor<MatrixType_, UpLo_> > {
18+
typedef MatrixType_ MatrixType;
19+
typedef typename MatrixType::Scalar Scalar;
20+
typedef typename MatrixType::RealScalar RealScalar;
21+
22+
typedef Eigen::CholmodSimplicialLLT<MatrixType_, UpLo_> Solver;
23+
24+
template <class PyClass>
25+
void visit(PyClass &cl) const {
26+
cl
27+
28+
.def(CholmodBaseVisitor<Solver>())
29+
.def(bp::init<>(bp::arg("self"), "Default constructor"))
30+
.def(bp::init<MatrixType>(bp::args("self", "matrix"),
31+
"Constructs and performs the LLT "
32+
"factorization from a given matrix."))
33+
34+
;
35+
}
36+
37+
static void expose() {
38+
static const std::string classname =
39+
"CholmodSimplicialLLT_" + scalar_name<Scalar>::shortname();
40+
expose(classname);
41+
}
42+
43+
static void expose(const std::string &name) {
44+
bp::class_<Solver, boost::noncopyable>(
45+
name.c_str(),
46+
"A simplicial direct Cholesky (LLT) factorization and solver based on "
47+
"Cholmod.\n\n"
48+
"This class allows to solve for A.X = B sparse linear problems via a "
49+
"simplicial LL^T Cholesky factorization using the Cholmod library."
50+
"This simplicial variant is equivalent to Eigen's built-in "
51+
"SimplicialLLT class."
52+
"Therefore, it has little practical interest. The sparse matrix A must "
53+
"be selfadjoint and positive definite."
54+
"The vectors or matrices X and B can be either dense or sparse.",
55+
bp::no_init)
56+
.def(CholmodSimplicialLLTVisitor());
57+
}
58+
};
59+
60+
} // namespace eigenpy
61+
62+
#endif // ifndef
63+
// __eigenpy_decomposition_sparse_cholmod_cholmod_simplicial_llt_hpp__
Lines changed: 62 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,62 @@
1+
/*
2+
* Copyright 2024 INRIA
3+
*/
4+
5+
#ifndef __eigenpy_decomposition_sparse_cholmod_cholmod_supernodal_llt_hpp__
6+
#define __eigenpy_decomposition_sparse_cholmod_cholmod_supernodal_llt_hpp__
7+
8+
#include "eigenpy/eigenpy.hpp"
9+
#include "eigenpy/decompositions/sparse/cholmod/CholmodDecomposition.hpp"
10+
#include "eigenpy/utils/scalar-name.hpp"
11+
12+
namespace eigenpy {
13+
14+
template <typename MatrixType_, int UpLo_ = Eigen::Lower>
15+
struct CholmodSupernodalLLTVisitor
16+
: public boost::python::def_visitor<
17+
CholmodSupernodalLLTVisitor<MatrixType_, UpLo_> > {
18+
typedef MatrixType_ MatrixType;
19+
typedef typename MatrixType::Scalar Scalar;
20+
typedef typename MatrixType::RealScalar RealScalar;
21+
22+
typedef Eigen::CholmodSupernodalLLT<MatrixType_, UpLo_> Solver;
23+
24+
template <class PyClass>
25+
void visit(PyClass &cl) const {
26+
cl
27+
28+
.def(CholmodBaseVisitor<Solver>())
29+
.def(bp::init<>(bp::arg("self"), "Default constructor"))
30+
.def(bp::init<MatrixType>(bp::args("self", "matrix"),
31+
"Constructs and performs the LLT "
32+
"factorization from a given matrix."))
33+
34+
;
35+
}
36+
37+
static void expose() {
38+
static const std::string classname =
39+
"CholmodSupernodalLLT_" + scalar_name<Scalar>::shortname();
40+
expose(classname);
41+
}
42+
43+
static void expose(const std::string &name) {
44+
bp::class_<Solver, boost::noncopyable>(
45+
name.c_str(),
46+
"A supernodal direct Cholesky (LLT) factorization and solver based on "
47+
"Cholmod.\n\n"
48+
"This class allows to solve for A.X = B sparse linear problems via a "
49+
"supernodal LL^T Cholesky factorization using the Cholmod library."
50+
"This supernodal variant performs best on dense enough problems, e.g., "
51+
"3D FEM, or very high order 2D FEM."
52+
"The sparse matrix A must be selfadjoint and positive definite. The "
53+
"vectors or matrices X and B can be either dense or sparse.",
54+
bp::no_init)
55+
.def(CholmodSupernodalLLTVisitor());
56+
}
57+
};
58+
59+
} // namespace eigenpy
60+
61+
#endif // ifndef
62+
// __eigenpy_decomposition_sparse_cholmod_cholmod_supernodal_llt_hpp__

0 commit comments

Comments
 (0)