From fb2fca90be39783f76ba05b521360afeeda265f2 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Damiano=20Franz=C3=B2?= Date: Thu, 1 May 2025 23:17:21 +0000 Subject: [PATCH] Avoid unnecessary matrix copy in BDCSVD and JacobiSVD --- Eigen/src/SVD/BDCSVD.h | 21 ++++++++++++++++----- Eigen/src/SVD/JacobiSVD.h | 32 ++++++++++++++++++++++++-------- 2 files changed, 40 insertions(+), 13 deletions(-) diff --git a/Eigen/src/SVD/BDCSVD.h b/Eigen/src/SVD/BDCSVD.h index 64fec7590..85661a30d 100644 --- a/Eigen/src/SVD/BDCSVD.h +++ b/Eigen/src/SVD/BDCSVD.h @@ -165,7 +165,8 @@ class BDCSVD : public SVDBase > { * * \param matrix the matrix to decompose */ - BDCSVD(const MatrixType& matrix) : m_algoswap(16), m_numIters(0) { + template + BDCSVD(const MatrixBase& matrix) : m_algoswap(16), m_numIters(0) { compute_impl(matrix, internal::get_computation_options(Options)); } @@ -193,7 +194,10 @@ class BDCSVD : public SVDBase > { * * \param matrix the matrix to decompose */ - BDCSVD& compute(const MatrixType& matrix) { return compute_impl(matrix, m_computationOptions); } + template + BDCSVD& compute(const MatrixBase& matrix) { + return compute_impl(matrix, m_computationOptions); + } /** \brief Method performing the decomposition of given matrix, as specified by * the `computationOptions` parameter. @@ -204,7 +208,8 @@ class BDCSVD : public SVDBase > { * \deprecated Will be removed in the next major Eigen version. Options should * be specified in the \a Options template parameter. */ - EIGEN_DEPRECATED BDCSVD& compute(const MatrixType& matrix, unsigned int computationOptions) { + template + EIGEN_DEPRECATED BDCSVD& compute(const MatrixBase& matrix, unsigned int computationOptions) { internal::check_svd_options_assertions(computationOptions, matrix.rows(), matrix.cols()); return compute_impl(matrix, computationOptions); } @@ -215,7 +220,8 @@ class BDCSVD : public SVDBase > { } private: - BDCSVD& compute_impl(const MatrixType& matrix, unsigned int computationOptions); + template + BDCSVD& compute_impl(const MatrixBase& matrix, unsigned int computationOptions); void divide(Index firstCol, Index lastCol, Index firstRowW, Index firstColW, Index shift); void computeSVDofM(Index firstCol, Index n, MatrixXr& U, VectorType& singVals, MatrixXr& V); void computeSingVals(const ArrayRef& col0, const ArrayRef& diag, const IndicesRef& perm, VectorType& singVals, @@ -307,8 +313,13 @@ void BDCSVD::allocate(Index rows, Index cols, unsigned int } // end allocate template -BDCSVD& BDCSVD::compute_impl(const MatrixType& matrix, +template +BDCSVD& BDCSVD::compute_impl(const MatrixBase& matrix, unsigned int computationOptions) { + EIGEN_STATIC_ASSERT_SAME_MATRIX_SIZE(Derived, MatrixType); + EIGEN_STATIC_ASSERT((std::is_same::value), + Input matrix must have the same Scalar type as the BDCSVD object.); + #ifdef EIGEN_BDCSVD_DEBUG_VERBOSE std::cout << "\n\n\n=================================================================================================" "=====================\n\n\n"; diff --git a/Eigen/src/SVD/JacobiSVD.h b/Eigen/src/SVD/JacobiSVD.h index ecbc3638a..1abde17fd 100644 --- a/Eigen/src/SVD/JacobiSVD.h +++ b/Eigen/src/SVD/JacobiSVD.h @@ -565,7 +565,10 @@ class JacobiSVD : public SVDBase > { * * \param matrix the matrix to decompose */ - explicit JacobiSVD(const MatrixType& matrix) { compute_impl(matrix, internal::get_computation_options(Options)); } + template + explicit JacobiSVD(const MatrixBase& matrix) { + compute_impl(matrix, internal::get_computation_options(Options)); + } /** \brief Constructor performing the decomposition of given matrix using specified options * for computing unitaries. @@ -580,8 +583,10 @@ class JacobiSVD : public SVDBase > { * be specified in the \a Options template parameter. */ // EIGEN_DEPRECATED // TODO(cantonios): re-enable after fixing a few 3p libraries that error on deprecation warnings. - JacobiSVD(const MatrixType& matrix, unsigned int computationOptions) { - internal::check_svd_options_assertions(computationOptions, matrix.rows(), matrix.cols()); + template + JacobiSVD(const MatrixBase& matrix, unsigned int computationOptions) { + internal::check_svd_options_assertions, Options>(computationOptions, matrix.rows(), + matrix.cols()); compute_impl(matrix, computationOptions); } @@ -590,7 +595,10 @@ class JacobiSVD : public SVDBase > { * * \param matrix the matrix to decompose */ - JacobiSVD& compute(const MatrixType& matrix) { return compute_impl(matrix, m_computationOptions); } + template + JacobiSVD& compute(const MatrixBase& matrix) { + return compute_impl(matrix, m_computationOptions); + } /** \brief Method performing the decomposition of given matrix, as specified by * the `computationOptions` parameter. @@ -601,8 +609,10 @@ class JacobiSVD : public SVDBase > { * \deprecated Will be removed in the next major Eigen version. Options should * be specified in the \a Options template parameter. */ - EIGEN_DEPRECATED JacobiSVD& compute(const MatrixType& matrix, unsigned int computationOptions) { - internal::check_svd_options_assertions(m_computationOptions, matrix.rows(), matrix.cols()); + template + EIGEN_DEPRECATED JacobiSVD& compute(const MatrixBase& matrix, unsigned int computationOptions) { + internal::check_svd_options_assertions, Options>(m_computationOptions, matrix.rows(), + matrix.cols()); return compute_impl(matrix, computationOptions); } @@ -626,7 +636,8 @@ class JacobiSVD : public SVDBase > { } private: - JacobiSVD& compute_impl(const MatrixType& matrix, unsigned int computationOptions); + template + JacobiSVD& compute_impl(const MatrixBase& matrix, unsigned int computationOptions); protected: using Base::m_computationOptions; @@ -664,8 +675,13 @@ class JacobiSVD : public SVDBase > { }; template -JacobiSVD& JacobiSVD::compute_impl(const MatrixType& matrix, +template +JacobiSVD& JacobiSVD::compute_impl(const MatrixBase& matrix, unsigned int computationOptions) { + EIGEN_STATIC_ASSERT_SAME_MATRIX_SIZE(Derived, MatrixType); + EIGEN_STATIC_ASSERT((std::is_same::value), + Input matrix must have the same Scalar type as the BDCSVD object.); + using std::abs; allocate(matrix.rows(), matrix.cols(), computationOptions);