Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Minor enhancements to custom LDLT implementations #489

Merged
merged 2 commits into from
Aug 2, 2024
+101 −17
Merged

Minor enhancements to custom LDLT implementations #489

Changes from 1 commit
Commits
Show all changes
2 commits
Select commit Hold shift + click to select a range
2b48438
Sparse matrix support in block-diagonal LDLT
peddie Jun 5, 2024
e0ba89c
Condition number estimation for custom LDLTs
peddie Jun 5, 2024
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
Next Next commit
Sparse matrix support in block-diagonal LDLT
@peddie
peddie committed Aug 2, 2024
commit 2b484380456117992954982a31e79480a1a8e4db
40 changes: 40 additions & 0 deletions include/albatross/src/linalg/block_diagonal.hpp
View file
Original file line number Diff line number Diff line change
@@ -32,6 +32,14 @@ struct BlockDiagonalLDLT {
Eigen::Matrix<_Scalar, _Rows, _Cols>
sqrt_solve(const Eigen::Matrix<_Scalar, _Rows, _Cols> &rhs) const;

template <class _Scalar, int _Options, typename _StorageIndex>
Eigen::Matrix<_Scalar, Eigen::Dynamic, Eigen::Dynamic>
solve(const Eigen::SparseMatrix<_Scalar, _Options, _StorageIndex> &rhs) const;

template <class _Scalar, int _Options, typename _StorageIndex>
Eigen::Matrix<_Scalar, Eigen::Dynamic, Eigen::Dynamic> sqrt_solve(
const Eigen::SparseMatrix<_Scalar, _Options, _StorageIndex> &rhs) const;

template <class _Scalar, int _Rows, int _Cols>
Eigen::Matrix<_Scalar, _Rows, _Cols>
solve(const Eigen::Matrix<_Scalar, _Rows, _Cols> &rhs,
@@ -102,6 +110,38 @@ inline Eigen::Matrix<_Scalar, _Rows, _Cols> BlockDiagonalLDLT::sqrt_solve(
return output;
}

template <class _Scalar, int _Options, typename _StorageIndex>
inline Eigen::Matrix<_Scalar, Eigen::Dynamic, Eigen::Dynamic>
BlockDiagonalLDLT::solve(
const Eigen::SparseMatrix<_Scalar, _Options, _StorageIndex> &rhs) const {
ALBATROSS_ASSERT(cols() == rhs.rows());
Eigen::Index i = 0;
Eigen::Matrix<_Scalar, Eigen::Dynamic, Eigen::Dynamic> output(rows(),
rhs.cols());
for (const auto &b : blocks) {
const auto rhs_chunk = rhs.block(i, 0, b.rows(), rhs.cols());
output.block(i, 0, b.rows(), rhs.cols()) = b.solve(rhs_chunk);
i += b.rows();
}
return output;
}

template <class _Scalar, int _Options, typename _StorageIndex>
inline Eigen::Matrix<_Scalar, Eigen::Dynamic, Eigen::Dynamic>
BlockDiagonalLDLT::sqrt_solve(
const Eigen::SparseMatrix<_Scalar, _Options, _StorageIndex> &rhs) const {
ALBATROSS_ASSERT(cols() == rhs.rows());
Eigen::Index i = 0;
Eigen::Matrix<_Scalar, Eigen::Dynamic, Eigen::Dynamic> output(rows(),
rhs.cols());
for (const auto &b : blocks) {
const auto rhs_chunk = rhs.block(i, 0, b.rows(), rhs.cols());
output.block(i, 0, b.rows(), rhs.cols()) = b.sqrt_solve(rhs_chunk);
i += b.rows();
}
return output;
}

inline std::map<size_t, Eigen::Index>
BlockDiagonalLDLT::block_to_row_map() const {
Eigen::Index row = 0;