From 21c0a2ce0c317fa258380c92bea4e9e16e0840f2 Mon Sep 17 00:00:00 2001 From: Gael Guennebaud Date: Wed, 29 Oct 2014 11:29:33 +0100 Subject: [PATCH] Move D&C SVD to official SVD module. --- Eigen/SVD | 9 ++++-- Eigen/src/Core/MatrixBase.h | 1 + Eigen/src/Core/util/ForwardDeclarations.h | 1 + .../src/BDCSVD => Eigen/src/SVD}/BDCSVD.h | 23 +++++++------- test/CMakeLists.txt | 1 + {unsupported/test => test}/bdcsvd.cpp | 30 +++++++++---------- unsupported/Eigen/BDCSVD | 26 ---------------- unsupported/Eigen/src/BDCSVD/CMakeLists.txt | 6 ---- unsupported/Eigen/src/BDCSVD/TODOBdcsvd.txt | 13 -------- unsupported/Eigen/src/BDCSVD/doneInBDCSVD.txt | 8 ----- unsupported/Eigen/src/CMakeLists.txt | 1 - unsupported/test/CMakeLists.txt | 1 - 12 files changed, 37 insertions(+), 83 deletions(-) rename {unsupported/Eigen/src/BDCSVD => Eigen/src/SVD}/BDCSVD.h (99%) rename {unsupported/test => test}/bdcsvd.cpp (84%) delete mode 100644 unsupported/Eigen/BDCSVD delete mode 100644 unsupported/Eigen/src/BDCSVD/CMakeLists.txt delete mode 100644 unsupported/Eigen/src/BDCSVD/TODOBdcsvd.txt delete mode 100644 unsupported/Eigen/src/BDCSVD/doneInBDCSVD.txt diff --git a/Eigen/SVD b/Eigen/SVD index c13472e82..dbd37b17a 100644 --- a/Eigen/SVD +++ b/Eigen/SVD @@ -12,20 +12,25 @@ * * * This module provides SVD decomposition for matrices (both real and complex). - * This decomposition is accessible via the following MatrixBase method: + * Two decomposition algorithms are provided: + * - JacobiSVD implementing two-sided Jacobi iterations is numerically very accurate, fast for small matrices, but very slow for larger ones. + * - BDCSVD implementing a recursive divide & conquer strategy on top of an upper-bidiagonalization which remains fast for large problems. + * These decompositions are accessible via the respective classes and following MatrixBase methods: * - MatrixBase::jacobiSvd() + * - MatrixBase::bdcSvd() * * \code * #include * \endcode */ +#include "src/SVD/UpperBidiagonalization.h" #include "src/SVD/SVDBase.h" #include "src/SVD/JacobiSVD.h" +#include "src/SVD/BDCSVD.h" #if defined(EIGEN_USE_LAPACKE) && !defined(EIGEN_USE_LAPACKE_STRICT) #include "src/SVD/JacobiSVD_MKL.h" #endif -#include "src/SVD/UpperBidiagonalization.h" #include "src/Core/util/ReenableStupidWarnings.h" diff --git a/Eigen/src/Core/MatrixBase.h b/Eigen/src/Core/MatrixBase.h index 001513187..135d8dfc8 100644 --- a/Eigen/src/Core/MatrixBase.h +++ b/Eigen/src/Core/MatrixBase.h @@ -390,6 +390,7 @@ template class MatrixBase /////////// SVD module /////////// JacobiSVD jacobiSvd(unsigned int computationOptions = 0) const; + BDCSVD bdcSvd(unsigned int computationOptions = 0) const; /////////// Geometry module /////////// diff --git a/Eigen/src/Core/util/ForwardDeclarations.h b/Eigen/src/Core/util/ForwardDeclarations.h index be156a44a..1f7503dfa 100644 --- a/Eigen/src/Core/util/ForwardDeclarations.h +++ b/Eigen/src/Core/util/ForwardDeclarations.h @@ -251,6 +251,7 @@ template class HouseholderQR; template class ColPivHouseholderQR; template class FullPivHouseholderQR; template class JacobiSVD; +template class BDCSVD; template class LLT; template class LDLT; template class HouseholderSequence; diff --git a/unsupported/Eigen/src/BDCSVD/BDCSVD.h b/Eigen/src/SVD/BDCSVD.h similarity index 99% rename from unsupported/Eigen/src/BDCSVD/BDCSVD.h rename to Eigen/src/SVD/BDCSVD.h index 94bc03aac..498541050 100644 --- a/unsupported/Eigen/src/BDCSVD/BDCSVD.h +++ b/Eigen/src/SVD/BDCSVD.h @@ -514,6 +514,7 @@ void BDCSVD::divide (Index firstCol, Index lastCol, Index firstRowW, // // TODO Opportunities for optimization: better root finding algo, better stopping criterion, better // handling of round-off errors, be consistent in ordering +// For instance, to solve the secular equation using FMM, see http://www.stat.uchicago.edu/~lekheng/courses/302/classics/greengard-rokhlin.pdf template void BDCSVD::computeSVDofM(Index firstCol, Index n, MatrixXr& U, VectorType& singVals, MatrixXr& V) { @@ -528,7 +529,10 @@ void BDCSVD::computeSVDofM(Index firstCol, Index n, MatrixXr& U, Vec U.resize(n+1, n+1); if (m_compV) V.resize(n, n); - if (col0.hasNaN() || diag.hasNaN()) { std::cout << "\n\nHAS NAN\n\n"; return; } +#ifdef EIGEN_BDCSVD_DEBUG_VERBOSE + if (col0.hasNaN() || diag.hasNaN()) + std::cout << "\n\nHAS NAN\n\n"; +#endif // Many singular values might have been deflated, the zero ones have been moved to the end, // but others are interleaved and we must ignore them at this stage. @@ -544,7 +548,7 @@ void BDCSVD::computeSVDofM(Index firstCol, Index n, MatrixXr& U, Vec ArrayXr shifts(n), mus(n), zhat(n); -#ifdef EIGEN_BDCSVD_DEBUG_VERBOSE +#ifdef EIGEN_BDCSVD_DEBUG_VERBOSE std::cout << "computeSVDofM using:\n"; std::cout << " z: " << col0.transpose() << "\n"; std::cout << " d: " << diag.transpose() << "\n"; @@ -1149,21 +1153,18 @@ void BDCSVD::deflation(Index firstCol, Index lastCol, Index k, Index }//end deflation - /** \svd_module - * - * \return the singular value decomposition of \c *this computed by - * BDC Algorithm - * - * \sa class BDCSVD - */ -/* +/** \svd_module + * + * \return the singular value decomposition of \c *this computed by Divide & Conquer algorithm + * + * \sa class BDCSVD + */ template BDCSVD::PlainObject> MatrixBase::bdcSvd(unsigned int computationOptions) const { return BDCSVD(*this, computationOptions); } -*/ } // end namespace Eigen diff --git a/test/CMakeLists.txt b/test/CMakeLists.txt index 63ed2c7a4..f57d8ce36 100644 --- a/test/CMakeLists.txt +++ b/test/CMakeLists.txt @@ -211,6 +211,7 @@ ei_add_test(real_qz) ei_add_test(eigensolver_generalized_real) ei_add_test(jacobi) ei_add_test(jacobisvd) +ei_add_test(bdcsvd) ei_add_test(householder) ei_add_test(geo_orthomethods) ei_add_test(geo_quaternion) diff --git a/unsupported/test/bdcsvd.cpp b/test/bdcsvd.cpp similarity index 84% rename from unsupported/test/bdcsvd.cpp rename to test/bdcsvd.cpp index 97d7880bb..d58d259d6 100644 --- a/unsupported/test/bdcsvd.cpp +++ b/test/bdcsvd.cpp @@ -15,7 +15,7 @@ #define EIGEN_RUNTIME_NO_MALLOC #include "main.h" -#include +#include #include #include @@ -35,18 +35,18 @@ void bdcsvd(const MatrixType& a = MatrixType(), bool pickrandom = true) CALL_SUBTEST(( svd_test_all_computation_options >(m, false) )); } -// template -// void bdcsvd_method() -// { -// enum { Size = MatrixType::RowsAtCompileTime }; -// typedef typename MatrixType::RealScalar RealScalar; -// typedef Matrix RealVecType; -// MatrixType m = MatrixType::Identity(); -// VERIFY_IS_APPROX(m.bdcSvd().singularValues(), RealVecType::Ones()); -// VERIFY_RAISES_ASSERT(m.bdcSvd().matrixU()); -// VERIFY_RAISES_ASSERT(m.bdcSvd().matrixV()); -// VERIFY_IS_APPROX(m.bdcSvd(ComputeFullU|ComputeFullV).solve(m), m); -// } +template +void bdcsvd_method() +{ + enum { Size = MatrixType::RowsAtCompileTime }; + typedef typename MatrixType::RealScalar RealScalar; + typedef Matrix RealVecType; + MatrixType m = MatrixType::Identity(); + VERIFY_IS_APPROX(m.bdcSvd().singularValues(), RealVecType::Ones()); + VERIFY_RAISES_ASSERT(m.bdcSvd().matrixU()); + VERIFY_RAISES_ASSERT(m.bdcSvd().matrixV()); + VERIFY_IS_APPROX(m.bdcSvd(ComputeFullU|ComputeFullV).solve(m), m); +} // compare the Singular values returned with Jacobi and Bdc template @@ -97,8 +97,8 @@ void test_bdcsvd() } // test matrixbase method -// CALL_SUBTEST_1(( bdcsvd_method() )); -// CALL_SUBTEST_3(( bdcsvd_method() )); + CALL_SUBTEST_1(( bdcsvd_method() )); + CALL_SUBTEST_3(( bdcsvd_method() )); // Test problem size constructors CALL_SUBTEST_7( BDCSVD(10,10) ); diff --git a/unsupported/Eigen/BDCSVD b/unsupported/Eigen/BDCSVD deleted file mode 100644 index 44649dbd0..000000000 --- a/unsupported/Eigen/BDCSVD +++ /dev/null @@ -1,26 +0,0 @@ -#ifndef EIGEN_BDCSVD_MODULE_H -#define EIGEN_BDCSVD_MODULE_H - -#include - -#include "../../Eigen/src/Core/util/DisableStupidWarnings.h" - -/** \defgroup BDCSVD_Module BDCSVD module - * - * - * - * This module provides Divide & Conquer SVD decomposition for matrices (both real and complex). - * This decomposition is accessible via the following MatrixBase method: - * - MatrixBase::bdcSvd() - * - * \code - * #include - * \endcode - */ - -#include "src/BDCSVD/BDCSVD.h" - -#include "../../Eigen/src/Core/util/ReenableStupidWarnings.h" - -#endif // EIGEN_BDCSVD_MODULE_H -/* vim: set filetype=cpp et sw=2 ts=2 ai: */ diff --git a/unsupported/Eigen/src/BDCSVD/CMakeLists.txt b/unsupported/Eigen/src/BDCSVD/CMakeLists.txt deleted file mode 100644 index 1045512f9..000000000 --- a/unsupported/Eigen/src/BDCSVD/CMakeLists.txt +++ /dev/null @@ -1,6 +0,0 @@ -FILE(GLOB Eigen_BDCSVD_SRCS "*.h") - -INSTALL(FILES - ${Eigen_BDCSVD_SRCS} - DESTINATION ${INCLUDE_INSTALL_DIR}/unsupported/Eigen/src/BDCSVD COMPONENT Devel - ) diff --git a/unsupported/Eigen/src/BDCSVD/TODOBdcsvd.txt b/unsupported/Eigen/src/BDCSVD/TODOBdcsvd.txt deleted file mode 100644 index 9b7bf9314..000000000 --- a/unsupported/Eigen/src/BDCSVD/TODOBdcsvd.txt +++ /dev/null @@ -1,13 +0,0 @@ -TO DO LIST - -- check more carefully single precision -- check with duplicated singularvalues -- no-malloc mode - -(optional optimization) -- do all the allocations in the allocate part -- support static matrices -- return a error at compilation time when using integer matrices (int, long, std::complex, ...) -- To solve the secular equation using FMM: -http://www.stat.uchicago.edu/~lekheng/courses/302/classics/greengard-rokhlin.pdf - diff --git a/unsupported/Eigen/src/BDCSVD/doneInBDCSVD.txt b/unsupported/Eigen/src/BDCSVD/doneInBDCSVD.txt deleted file mode 100644 index 29ab9cd40..000000000 --- a/unsupported/Eigen/src/BDCSVD/doneInBDCSVD.txt +++ /dev/null @@ -1,8 +0,0 @@ -This unsupported package is about a divide and conquer algorithm to compute SVD. - -The implementation follows as closely as possible the following reference paper : -http://www.cs.yale.edu/publications/techreports/tr933.pdf - -To solve the secular equation using FMM: -http://www.stat.uchicago.edu/~lekheng/courses/302/classics/greengard-rokhlin.pdf - diff --git a/unsupported/Eigen/src/CMakeLists.txt b/unsupported/Eigen/src/CMakeLists.txt index 654a2327f..8eb2808e3 100644 --- a/unsupported/Eigen/src/CMakeLists.txt +++ b/unsupported/Eigen/src/CMakeLists.txt @@ -12,4 +12,3 @@ ADD_SUBDIRECTORY(Skyline) ADD_SUBDIRECTORY(SparseExtra) ADD_SUBDIRECTORY(KroneckerProduct) ADD_SUBDIRECTORY(Splines) -ADD_SUBDIRECTORY(BDCSVD) diff --git a/unsupported/test/CMakeLists.txt b/unsupported/test/CMakeLists.txt index 48b61cde0..97849a25a 100644 --- a/unsupported/test/CMakeLists.txt +++ b/unsupported/test/CMakeLists.txt @@ -92,7 +92,6 @@ ei_add_test(splines) ei_add_test(gmres) ei_add_test(minres) ei_add_test(levenberg_marquardt) -ei_add_test(bdcsvd) ei_add_test(kronecker_product) option(EIGEN_TEST_CXX11 "Enable testing of C++11 features (e.g. Tensor module)." OFF)