Fix jacobi svd for TriangularBase

This commit is contained in:
Damiano Franzò
2025-10-08 00:41:26 +02:00
committed by Rasmus Munk Larsen
parent dbe9e6961e
commit 5bc944a3ef
2 changed files with 37 additions and 0 deletions

View File

@@ -571,6 +571,11 @@ class JacobiSVD : public SVDBase<JacobiSVD<MatrixType_, Options_> > {
compute_impl(matrix, internal::get_computation_options(Options));
}
template <typename Derived>
explicit JacobiSVD(const TriangularBase<Derived>& matrix) {
compute_impl(matrix, internal::get_computation_options(Options));
}
/** \brief Constructor performing the decomposition of given matrix using specified options
* for computing unitaries.
*
@@ -601,6 +606,11 @@ class JacobiSVD : public SVDBase<JacobiSVD<MatrixType_, Options_> > {
return compute_impl(matrix, m_computationOptions);
}
template <typename Derived>
JacobiSVD& compute(const TriangularBase<Derived>& matrix) {
return compute_impl(matrix, m_computationOptions);
}
/** \brief Method performing the decomposition of given matrix, as specified by
* the `computationOptions` parameter.
*
@@ -638,6 +648,8 @@ class JacobiSVD : public SVDBase<JacobiSVD<MatrixType_, Options_> > {
}
private:
template <typename Derived>
JacobiSVD& compute_impl(const TriangularBase<Derived>& matrix, unsigned int computationOptions);
template <typename Derived>
JacobiSVD& compute_impl(const MatrixBase<Derived>& matrix, unsigned int computationOptions);
@@ -676,6 +688,13 @@ class JacobiSVD : public SVDBase<JacobiSVD<MatrixType_, Options_> > {
WorkMatrixType m_workMatrix;
};
template <typename MatrixType, int Options>
template <typename Derived>
JacobiSVD<MatrixType, Options>& JacobiSVD<MatrixType, Options>::compute_impl(const TriangularBase<Derived>& matrix,
unsigned int computationOptions) {
return compute_impl(matrix.toDenseMatrix(), computationOptions);
}
template <typename MatrixType, int Options>
template <typename Derived>
JacobiSVD<MatrixType, Options>& JacobiSVD<MatrixType, Options>::compute_impl(const MatrixBase<Derived>& matrix,

View File

@@ -94,6 +94,19 @@ void jacobisvd_verify_inputs(const MatrixType& input = MatrixType()) {
(int)ColPivHouseholderQRPreconditioner));
}
template <typename MatrixType>
void svd_triangular_matrix(const MatrixType& input = MatrixType()) {
MatrixType matrix(input.rows(), input.cols());
svd_fill_random(matrix);
// Make sure that we only consider the 'Lower' part of the matrix.
MatrixType matrix_self_adj = matrix.template selfadjointView<Lower>().toDenseMatrix();
JacobiSVD<MatrixType, ComputeFullV> svd_triangular(matrix.template selfadjointView<Lower>());
JacobiSVD<MatrixType, ComputeFullV> svd_full(matrix_self_adj);
VERIFY_IS_APPROX(svd_triangular.singularValues(), svd_full.singularValues());
}
namespace Foo {
// older compiler require a default constructor for Bar
// cf: https://stackoverflow.com/questions/7411515/
@@ -211,5 +224,10 @@ EIGEN_DECLARE_TEST(jacobisvd) {
CALL_SUBTEST_55(svd_underoverflow<void>());
// Check that the TriangularBase constructor works
CALL_SUBTEST_56((svd_triangular_matrix<Matrix3d>()));
CALL_SUBTEST_57((svd_triangular_matrix<Matrix4f>()));
CALL_SUBTEST_58((svd_triangular_matrix<Matrix<double, 10, 10>>()));
msvc_workaround();
}