1c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath// This file is part of Eigen, a lightweight C++ template library 2c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath// for linear algebra. 3c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath// 4c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath// Copyright (C) 2008-2009 Gael Guennebaud <gael.guennebaud@inria.fr> 5c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath// Copyright (C) 2009 Benoit Jacob <jacob.benoit.1@gmail.com> 6c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath// 7c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath// This Source Code Form is subject to the terms of the Mozilla 8c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath// Public License v. 2.0. If a copy of the MPL was not distributed 9c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath// with this file, You can obtain one at http://mozilla.org/MPL/2.0/. 10c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 11c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath#ifndef EIGEN_FULLPIVOTINGHOUSEHOLDERQR_H 12c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath#define EIGEN_FULLPIVOTINGHOUSEHOLDERQR_H 13c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 14c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamathnamespace Eigen { 15c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 16c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamathnamespace internal { 17c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 18c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamathtemplate<typename MatrixType> struct FullPivHouseholderQRMatrixQReturnType; 19c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 20c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamathtemplate<typename MatrixType> 21c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamathstruct traits<FullPivHouseholderQRMatrixQReturnType<MatrixType> > 22c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath{ 23c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typedef typename MatrixType::PlainObject ReturnType; 24c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath}; 25c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 26c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath} 27c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 28c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath/** \ingroup QR_Module 29c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * 30c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \class FullPivHouseholderQR 31c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * 32c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \brief Householder rank-revealing QR decomposition of a matrix with full pivoting 33c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * 34c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \param MatrixType the type of the matrix of which we are computing the QR decomposition 35c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * 36c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * This class performs a rank-revealing QR decomposition of a matrix \b A into matrices \b P, \b Q and \b R 37c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * such that 38c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \f[ 39c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \mathbf{A} \, \mathbf{P} = \mathbf{Q} \, \mathbf{R} 40c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \f] 41c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * by using Householder transformations. Here, \b P is a permutation matrix, \b Q a unitary matrix and \b R an 42c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * upper triangular matrix. 43c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * 44c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * This decomposition performs a very prudent full pivoting in order to be rank-revealing and achieve optimal 45c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * numerical stability. The trade-off is that it is slower than HouseholderQR and ColPivHouseholderQR. 46c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * 47c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \sa MatrixBase::fullPivHouseholderQr() 48c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath */ 49c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamathtemplate<typename _MatrixType> class FullPivHouseholderQR 50c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath{ 51c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath public: 52c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 53c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typedef _MatrixType MatrixType; 54c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath enum { 55c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath RowsAtCompileTime = MatrixType::RowsAtCompileTime, 56c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath ColsAtCompileTime = MatrixType::ColsAtCompileTime, 57c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath Options = MatrixType::Options, 58c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath MaxRowsAtCompileTime = MatrixType::MaxRowsAtCompileTime, 59c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath MaxColsAtCompileTime = MatrixType::MaxColsAtCompileTime 60c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath }; 61c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typedef typename MatrixType::Scalar Scalar; 62c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typedef typename MatrixType::RealScalar RealScalar; 63c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typedef typename MatrixType::Index Index; 64c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typedef internal::FullPivHouseholderQRMatrixQReturnType<MatrixType> MatrixQReturnType; 65c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typedef typename internal::plain_diag_type<MatrixType>::type HCoeffsType; 667faaa9f3f0df9d23790277834d426c3d992ac3baCarlos Hernandez typedef Matrix<Index, 1, 677faaa9f3f0df9d23790277834d426c3d992ac3baCarlos Hernandez EIGEN_SIZE_MIN_PREFER_DYNAMIC(ColsAtCompileTime,RowsAtCompileTime), RowMajor, 1, 687faaa9f3f0df9d23790277834d426c3d992ac3baCarlos Hernandez EIGEN_SIZE_MIN_PREFER_FIXED(MaxColsAtCompileTime,MaxRowsAtCompileTime)> IntDiagSizeVectorType; 69c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typedef PermutationMatrix<ColsAtCompileTime, MaxColsAtCompileTime> PermutationType; 70c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typedef typename internal::plain_row_type<MatrixType>::type RowVectorType; 71c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typedef typename internal::plain_col_type<MatrixType>::type ColVectorType; 72c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 73c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath /** \brief Default Constructor. 74c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * 75c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * The default constructor is useful in cases in which the user intends to 76c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * perform decompositions via FullPivHouseholderQR::compute(const MatrixType&). 77c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath */ 78c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath FullPivHouseholderQR() 79c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath : m_qr(), 80c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath m_hCoeffs(), 81c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath m_rows_transpositions(), 82c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath m_cols_transpositions(), 83c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath m_cols_permutation(), 84c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath m_temp(), 85c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath m_isInitialized(false), 86c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath m_usePrescribedThreshold(false) {} 87c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 88c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath /** \brief Default Constructor with memory preallocation 89c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * 90c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * Like the default constructor but with preallocation of the internal data 91c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * according to the specified problem \a size. 92c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \sa FullPivHouseholderQR() 93c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath */ 94c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath FullPivHouseholderQR(Index rows, Index cols) 95c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath : m_qr(rows, cols), 96c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath m_hCoeffs((std::min)(rows,cols)), 977faaa9f3f0df9d23790277834d426c3d992ac3baCarlos Hernandez m_rows_transpositions((std::min)(rows,cols)), 987faaa9f3f0df9d23790277834d426c3d992ac3baCarlos Hernandez m_cols_transpositions((std::min)(rows,cols)), 99c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath m_cols_permutation(cols), 1007faaa9f3f0df9d23790277834d426c3d992ac3baCarlos Hernandez m_temp(cols), 101c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath m_isInitialized(false), 102c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath m_usePrescribedThreshold(false) {} 103c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 1047faaa9f3f0df9d23790277834d426c3d992ac3baCarlos Hernandez /** \brief Constructs a QR factorization from a given matrix 1057faaa9f3f0df9d23790277834d426c3d992ac3baCarlos Hernandez * 1067faaa9f3f0df9d23790277834d426c3d992ac3baCarlos Hernandez * This constructor computes the QR factorization of the matrix \a matrix by calling 1077faaa9f3f0df9d23790277834d426c3d992ac3baCarlos Hernandez * the method compute(). It is a short cut for: 1087faaa9f3f0df9d23790277834d426c3d992ac3baCarlos Hernandez * 1097faaa9f3f0df9d23790277834d426c3d992ac3baCarlos Hernandez * \code 1107faaa9f3f0df9d23790277834d426c3d992ac3baCarlos Hernandez * FullPivHouseholderQR<MatrixType> qr(matrix.rows(), matrix.cols()); 1117faaa9f3f0df9d23790277834d426c3d992ac3baCarlos Hernandez * qr.compute(matrix); 1127faaa9f3f0df9d23790277834d426c3d992ac3baCarlos Hernandez * \endcode 1137faaa9f3f0df9d23790277834d426c3d992ac3baCarlos Hernandez * 1147faaa9f3f0df9d23790277834d426c3d992ac3baCarlos Hernandez * \sa compute() 1157faaa9f3f0df9d23790277834d426c3d992ac3baCarlos Hernandez */ 116c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath FullPivHouseholderQR(const MatrixType& matrix) 117c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath : m_qr(matrix.rows(), matrix.cols()), 118c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath m_hCoeffs((std::min)(matrix.rows(), matrix.cols())), 1197faaa9f3f0df9d23790277834d426c3d992ac3baCarlos Hernandez m_rows_transpositions((std::min)(matrix.rows(), matrix.cols())), 1207faaa9f3f0df9d23790277834d426c3d992ac3baCarlos Hernandez m_cols_transpositions((std::min)(matrix.rows(), matrix.cols())), 121c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath m_cols_permutation(matrix.cols()), 1227faaa9f3f0df9d23790277834d426c3d992ac3baCarlos Hernandez m_temp(matrix.cols()), 123c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath m_isInitialized(false), 124c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath m_usePrescribedThreshold(false) 125c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 126c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath compute(matrix); 127c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 128c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 129c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath /** This method finds a solution x to the equation Ax=b, where A is the matrix of which 1307faaa9f3f0df9d23790277834d426c3d992ac3baCarlos Hernandez * \c *this is the QR decomposition. 131c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * 132c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \param b the right-hand-side of the equation to solve. 133c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * 1347faaa9f3f0df9d23790277834d426c3d992ac3baCarlos Hernandez * \returns the exact or least-square solution if the rank is greater or equal to the number of columns of A, 1357faaa9f3f0df9d23790277834d426c3d992ac3baCarlos Hernandez * and an arbitrary solution otherwise. 136c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * 137c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \note The case where b is a matrix is not yet implemented. Also, this 138c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * code is space inefficient. 139c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * 140c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \note_about_checking_solutions 141c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * 142c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \note_about_arbitrary_choice_of_solution 143c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * 144c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * Example: \include FullPivHouseholderQR_solve.cpp 145c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * Output: \verbinclude FullPivHouseholderQR_solve.out 146c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath */ 147c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath template<typename Rhs> 148c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath inline const internal::solve_retval<FullPivHouseholderQR, Rhs> 149c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath solve(const MatrixBase<Rhs>& b) const 150c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 151c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath eigen_assert(m_isInitialized && "FullPivHouseholderQR is not initialized."); 152c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath return internal::solve_retval<FullPivHouseholderQR, Rhs>(*this, b.derived()); 153c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 154c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 155c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath /** \returns Expression object representing the matrix Q 156c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath */ 157c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath MatrixQReturnType matrixQ(void) const; 158c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 159c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath /** \returns a reference to the matrix where the Householder QR decomposition is stored 160c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath */ 161c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath const MatrixType& matrixQR() const 162c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 163c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath eigen_assert(m_isInitialized && "FullPivHouseholderQR is not initialized."); 164c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath return m_qr; 165c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 166c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 167c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath FullPivHouseholderQR& compute(const MatrixType& matrix); 168c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 1697faaa9f3f0df9d23790277834d426c3d992ac3baCarlos Hernandez /** \returns a const reference to the column permutation matrix */ 170c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath const PermutationType& colsPermutation() const 171c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 172c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath eigen_assert(m_isInitialized && "FullPivHouseholderQR is not initialized."); 173c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath return m_cols_permutation; 174c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 175c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 1767faaa9f3f0df9d23790277834d426c3d992ac3baCarlos Hernandez /** \returns a const reference to the vector of indices representing the rows transpositions */ 1777faaa9f3f0df9d23790277834d426c3d992ac3baCarlos Hernandez const IntDiagSizeVectorType& rowsTranspositions() const 178c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 179c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath eigen_assert(m_isInitialized && "FullPivHouseholderQR is not initialized."); 180c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath return m_rows_transpositions; 181c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 182c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 183c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath /** \returns the absolute value of the determinant of the matrix of which 184c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * *this is the QR decomposition. It has only linear complexity 185c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * (that is, O(n) where n is the dimension of the square matrix) 186c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * as the QR decomposition has already been computed. 187c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * 188c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \note This is only for square matrices. 189c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * 190c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \warning a determinant can be very big or small, so for matrices 191c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * of large enough dimension, there is a risk of overflow/underflow. 192c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * One way to work around that is to use logAbsDeterminant() instead. 193c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * 194c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \sa logAbsDeterminant(), MatrixBase::determinant() 195c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath */ 196c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typename MatrixType::RealScalar absDeterminant() const; 197c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 198c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath /** \returns the natural log of the absolute value of the determinant of the matrix of which 199c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * *this is the QR decomposition. It has only linear complexity 200c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * (that is, O(n) where n is the dimension of the square matrix) 201c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * as the QR decomposition has already been computed. 202c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * 203c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \note This is only for square matrices. 204c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * 205c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \note This method is useful to work around the risk of overflow/underflow that's inherent 206c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * to determinant computation. 207c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * 208c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \sa absDeterminant(), MatrixBase::determinant() 209c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath */ 210c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typename MatrixType::RealScalar logAbsDeterminant() const; 211c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 212c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath /** \returns the rank of the matrix of which *this is the QR decomposition. 213c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * 214c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \note This method has to determine which pivots should be considered nonzero. 215c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * For that, it uses the threshold value that you can control by calling 216c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * setThreshold(const RealScalar&). 217c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath */ 218c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath inline Index rank() const 219c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 2207faaa9f3f0df9d23790277834d426c3d992ac3baCarlos Hernandez using std::abs; 221c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath eigen_assert(m_isInitialized && "FullPivHouseholderQR is not initialized."); 2227faaa9f3f0df9d23790277834d426c3d992ac3baCarlos Hernandez RealScalar premultiplied_threshold = abs(m_maxpivot) * threshold(); 223c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath Index result = 0; 224c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath for(Index i = 0; i < m_nonzero_pivots; ++i) 2257faaa9f3f0df9d23790277834d426c3d992ac3baCarlos Hernandez result += (abs(m_qr.coeff(i,i)) > premultiplied_threshold); 226c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath return result; 227c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 228c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 229c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath /** \returns the dimension of the kernel of the matrix of which *this is the QR decomposition. 230c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * 231c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \note This method has to determine which pivots should be considered nonzero. 232c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * For that, it uses the threshold value that you can control by calling 233c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * setThreshold(const RealScalar&). 234c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath */ 235c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath inline Index dimensionOfKernel() const 236c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 237c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath eigen_assert(m_isInitialized && "FullPivHouseholderQR is not initialized."); 238c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath return cols() - rank(); 239c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 240c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 241c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath /** \returns true if the matrix of which *this is the QR decomposition represents an injective 242c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * linear map, i.e. has trivial kernel; false otherwise. 243c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * 244c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \note This method has to determine which pivots should be considered nonzero. 245c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * For that, it uses the threshold value that you can control by calling 246c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * setThreshold(const RealScalar&). 247c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath */ 248c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath inline bool isInjective() const 249c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 250c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath eigen_assert(m_isInitialized && "FullPivHouseholderQR is not initialized."); 251c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath return rank() == cols(); 252c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 253c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 254c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath /** \returns true if the matrix of which *this is the QR decomposition represents a surjective 255c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * linear map; false otherwise. 256c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * 257c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \note This method has to determine which pivots should be considered nonzero. 258c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * For that, it uses the threshold value that you can control by calling 259c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * setThreshold(const RealScalar&). 260c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath */ 261c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath inline bool isSurjective() const 262c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 263c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath eigen_assert(m_isInitialized && "FullPivHouseholderQR is not initialized."); 264c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath return rank() == rows(); 265c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 266c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 267c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath /** \returns true if the matrix of which *this is the QR decomposition is invertible. 268c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * 269c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \note This method has to determine which pivots should be considered nonzero. 270c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * For that, it uses the threshold value that you can control by calling 271c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * setThreshold(const RealScalar&). 272c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath */ 273c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath inline bool isInvertible() const 274c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 275c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath eigen_assert(m_isInitialized && "FullPivHouseholderQR is not initialized."); 276c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath return isInjective() && isSurjective(); 277c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 278c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 279c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath /** \returns the inverse of the matrix of which *this is the QR decomposition. 280c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * 281c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \note If this matrix is not invertible, the returned matrix has undefined coefficients. 282c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * Use isInvertible() to first determine whether this matrix is invertible. 283c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath */ inline const 284c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath internal::solve_retval<FullPivHouseholderQR, typename MatrixType::IdentityReturnType> 285c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath inverse() const 286c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 287c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath eigen_assert(m_isInitialized && "FullPivHouseholderQR is not initialized."); 288c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath return internal::solve_retval<FullPivHouseholderQR,typename MatrixType::IdentityReturnType> 289c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath (*this, MatrixType::Identity(m_qr.rows(), m_qr.cols())); 290c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 291c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 292c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath inline Index rows() const { return m_qr.rows(); } 293c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath inline Index cols() const { return m_qr.cols(); } 2947faaa9f3f0df9d23790277834d426c3d992ac3baCarlos Hernandez 2957faaa9f3f0df9d23790277834d426c3d992ac3baCarlos Hernandez /** \returns a const reference to the vector of Householder coefficients used to represent the factor \c Q. 2967faaa9f3f0df9d23790277834d426c3d992ac3baCarlos Hernandez * 2977faaa9f3f0df9d23790277834d426c3d992ac3baCarlos Hernandez * For advanced uses only. 2987faaa9f3f0df9d23790277834d426c3d992ac3baCarlos Hernandez */ 299c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath const HCoeffsType& hCoeffs() const { return m_hCoeffs; } 300c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 301c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath /** Allows to prescribe a threshold to be used by certain methods, such as rank(), 302c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * who need to determine when pivots are to be considered nonzero. This is not used for the 303c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * QR decomposition itself. 304c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * 305c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * When it needs to get the threshold value, Eigen calls threshold(). By default, this 306c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * uses a formula to automatically determine a reasonable threshold. 307c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * Once you have called the present method setThreshold(const RealScalar&), 308c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * your value is used instead. 309c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * 310c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \param threshold The new value to use as the threshold. 311c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * 312c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * A pivot will be considered nonzero if its absolute value is strictly greater than 313c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \f$ \vert pivot \vert \leqslant threshold \times \vert maxpivot \vert \f$ 314c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * where maxpivot is the biggest pivot. 315c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * 316c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * If you want to come back to the default behavior, call setThreshold(Default_t) 317c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath */ 318c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath FullPivHouseholderQR& setThreshold(const RealScalar& threshold) 319c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 320c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath m_usePrescribedThreshold = true; 321c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath m_prescribedThreshold = threshold; 322c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath return *this; 323c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 324c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 325c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath /** Allows to come back to the default behavior, letting Eigen use its default formula for 326c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * determining the threshold. 327c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * 328c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * You should pass the special object Eigen::Default as parameter here. 329c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \code qr.setThreshold(Eigen::Default); \endcode 330c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * 331c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * See the documentation of setThreshold(const RealScalar&). 332c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath */ 333c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath FullPivHouseholderQR& setThreshold(Default_t) 334c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 335c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath m_usePrescribedThreshold = false; 336c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath return *this; 337c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 338c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 339c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath /** Returns the threshold that will be used by certain methods such as rank(). 340c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * 341c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * See the documentation of setThreshold(const RealScalar&). 342c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath */ 343c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath RealScalar threshold() const 344c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 345c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath eigen_assert(m_isInitialized || m_usePrescribedThreshold); 346c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath return m_usePrescribedThreshold ? m_prescribedThreshold 347c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath // this formula comes from experimenting (see "LU precision tuning" thread on the list) 348c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath // and turns out to be identical to Higham's formula used already in LDLt. 3497faaa9f3f0df9d23790277834d426c3d992ac3baCarlos Hernandez : NumTraits<Scalar>::epsilon() * RealScalar(m_qr.diagonalSize()); 350c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 351c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 352c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath /** \returns the number of nonzero pivots in the QR decomposition. 353c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * Here nonzero is meant in the exact sense, not in a fuzzy sense. 354c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * So that notion isn't really intrinsically interesting, but it is 355c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * still useful when implementing algorithms. 356c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * 357c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \sa rank() 358c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath */ 359c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath inline Index nonzeroPivots() const 360c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 361c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath eigen_assert(m_isInitialized && "LU is not initialized."); 362c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath return m_nonzero_pivots; 363c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 364c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 365c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath /** \returns the absolute value of the biggest pivot, i.e. the biggest 366c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * diagonal coefficient of U. 367c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath */ 368c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath RealScalar maxPivot() const { return m_maxpivot; } 369c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 370c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath protected: 371c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath MatrixType m_qr; 372c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath HCoeffsType m_hCoeffs; 3737faaa9f3f0df9d23790277834d426c3d992ac3baCarlos Hernandez IntDiagSizeVectorType m_rows_transpositions; 3747faaa9f3f0df9d23790277834d426c3d992ac3baCarlos Hernandez IntDiagSizeVectorType m_cols_transpositions; 375c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath PermutationType m_cols_permutation; 376c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath RowVectorType m_temp; 377c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath bool m_isInitialized, m_usePrescribedThreshold; 378c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath RealScalar m_prescribedThreshold, m_maxpivot; 379c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath Index m_nonzero_pivots; 380c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath RealScalar m_precision; 381c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath Index m_det_pq; 382c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath}; 383c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 384c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamathtemplate<typename MatrixType> 385c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamathtypename MatrixType::RealScalar FullPivHouseholderQR<MatrixType>::absDeterminant() const 386c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath{ 3877faaa9f3f0df9d23790277834d426c3d992ac3baCarlos Hernandez using std::abs; 388c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath eigen_assert(m_isInitialized && "FullPivHouseholderQR is not initialized."); 389c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath eigen_assert(m_qr.rows() == m_qr.cols() && "You can't take the determinant of a non-square matrix!"); 3907faaa9f3f0df9d23790277834d426c3d992ac3baCarlos Hernandez return abs(m_qr.diagonal().prod()); 391c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath} 392c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 393c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamathtemplate<typename MatrixType> 394c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamathtypename MatrixType::RealScalar FullPivHouseholderQR<MatrixType>::logAbsDeterminant() const 395c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath{ 396c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath eigen_assert(m_isInitialized && "FullPivHouseholderQR is not initialized."); 397c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath eigen_assert(m_qr.rows() == m_qr.cols() && "You can't take the determinant of a non-square matrix!"); 398c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath return m_qr.diagonal().cwiseAbs().array().log().sum(); 399c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath} 400c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 4017faaa9f3f0df9d23790277834d426c3d992ac3baCarlos Hernandez/** Performs the QR factorization of the given matrix \a matrix. The result of 4027faaa9f3f0df9d23790277834d426c3d992ac3baCarlos Hernandez * the factorization is stored into \c *this, and a reference to \c *this 4037faaa9f3f0df9d23790277834d426c3d992ac3baCarlos Hernandez * is returned. 4047faaa9f3f0df9d23790277834d426c3d992ac3baCarlos Hernandez * 4057faaa9f3f0df9d23790277834d426c3d992ac3baCarlos Hernandez * \sa class FullPivHouseholderQR, FullPivHouseholderQR(const MatrixType&) 4067faaa9f3f0df9d23790277834d426c3d992ac3baCarlos Hernandez */ 407c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamathtemplate<typename MatrixType> 408c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan KamathFullPivHouseholderQR<MatrixType>& FullPivHouseholderQR<MatrixType>::compute(const MatrixType& matrix) 409c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath{ 4107faaa9f3f0df9d23790277834d426c3d992ac3baCarlos Hernandez using std::abs; 411c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath Index rows = matrix.rows(); 412c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath Index cols = matrix.cols(); 413c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath Index size = (std::min)(rows,cols); 414c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 415c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath m_qr = matrix; 416c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath m_hCoeffs.resize(size); 417c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 418c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath m_temp.resize(cols); 419c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 4207faaa9f3f0df9d23790277834d426c3d992ac3baCarlos Hernandez m_precision = NumTraits<Scalar>::epsilon() * RealScalar(size); 421c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 4227faaa9f3f0df9d23790277834d426c3d992ac3baCarlos Hernandez m_rows_transpositions.resize(size); 4237faaa9f3f0df9d23790277834d426c3d992ac3baCarlos Hernandez m_cols_transpositions.resize(size); 424c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath Index number_of_transpositions = 0; 425c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 426c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath RealScalar biggest(0); 427c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 428c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath m_nonzero_pivots = size; // the generic case is that in which all pivots are nonzero (invertible case) 429c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath m_maxpivot = RealScalar(0); 430c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 431c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath for (Index k = 0; k < size; ++k) 432c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 433c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath Index row_of_biggest_in_corner, col_of_biggest_in_corner; 434c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath RealScalar biggest_in_corner; 435c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 436c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath biggest_in_corner = m_qr.bottomRightCorner(rows-k, cols-k) 437c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath .cwiseAbs() 438c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath .maxCoeff(&row_of_biggest_in_corner, &col_of_biggest_in_corner); 439c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath row_of_biggest_in_corner += k; 440c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath col_of_biggest_in_corner += k; 441c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath if(k==0) biggest = biggest_in_corner; 442c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 443c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath // if the corner is negligible, then we have less than full rank, and we can finish early 444c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath if(internal::isMuchSmallerThan(biggest_in_corner, biggest, m_precision)) 445c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 446c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath m_nonzero_pivots = k; 447c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath for(Index i = k; i < size; i++) 448c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 449c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath m_rows_transpositions.coeffRef(i) = i; 450c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath m_cols_transpositions.coeffRef(i) = i; 451c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath m_hCoeffs.coeffRef(i) = Scalar(0); 452c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 453c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath break; 454c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 455c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 456c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath m_rows_transpositions.coeffRef(k) = row_of_biggest_in_corner; 457c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath m_cols_transpositions.coeffRef(k) = col_of_biggest_in_corner; 458c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath if(k != row_of_biggest_in_corner) { 459c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath m_qr.row(k).tail(cols-k).swap(m_qr.row(row_of_biggest_in_corner).tail(cols-k)); 460c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath ++number_of_transpositions; 461c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 462c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath if(k != col_of_biggest_in_corner) { 463c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath m_qr.col(k).swap(m_qr.col(col_of_biggest_in_corner)); 464c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath ++number_of_transpositions; 465c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 466c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 467c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath RealScalar beta; 468c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath m_qr.col(k).tail(rows-k).makeHouseholderInPlace(m_hCoeffs.coeffRef(k), beta); 469c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath m_qr.coeffRef(k,k) = beta; 470c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 471c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath // remember the maximum absolute value of diagonal coefficients 4727faaa9f3f0df9d23790277834d426c3d992ac3baCarlos Hernandez if(abs(beta) > m_maxpivot) m_maxpivot = abs(beta); 473c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 474c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath m_qr.bottomRightCorner(rows-k, cols-k-1) 475c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath .applyHouseholderOnTheLeft(m_qr.col(k).tail(rows-k-1), m_hCoeffs.coeffRef(k), &m_temp.coeffRef(k+1)); 476c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 477c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 478c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath m_cols_permutation.setIdentity(cols); 479c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath for(Index k = 0; k < size; ++k) 480c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath m_cols_permutation.applyTranspositionOnTheRight(k, m_cols_transpositions.coeff(k)); 481c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 482c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath m_det_pq = (number_of_transpositions%2) ? -1 : 1; 483c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath m_isInitialized = true; 484c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 485c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath return *this; 486c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath} 487c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 488c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamathnamespace internal { 489c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 490c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamathtemplate<typename _MatrixType, typename Rhs> 491c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamathstruct solve_retval<FullPivHouseholderQR<_MatrixType>, Rhs> 492c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath : solve_retval_base<FullPivHouseholderQR<_MatrixType>, Rhs> 493c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath{ 494c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath EIGEN_MAKE_SOLVE_HELPERS(FullPivHouseholderQR<_MatrixType>,Rhs) 495c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 496c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath template<typename Dest> void evalTo(Dest& dst) const 497c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 498c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath const Index rows = dec().rows(), cols = dec().cols(); 499c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath eigen_assert(rhs().rows() == rows); 500c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 501c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath // FIXME introduce nonzeroPivots() and use it here. and more generally, 502c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath // make the same improvements in this dec as in FullPivLU. 503c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath if(dec().rank()==0) 504c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 505c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath dst.setZero(); 506c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath return; 507c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 508c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 509c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typename Rhs::PlainObject c(rhs()); 510c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 511c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath Matrix<Scalar,1,Rhs::ColsAtCompileTime> temp(rhs().cols()); 512c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath for (Index k = 0; k < dec().rank(); ++k) 513c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 514c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath Index remainingSize = rows-k; 515c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath c.row(k).swap(c.row(dec().rowsTranspositions().coeff(k))); 516c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath c.bottomRightCorner(remainingSize, rhs().cols()) 517c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath .applyHouseholderOnTheLeft(dec().matrixQR().col(k).tail(remainingSize-1), 518c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath dec().hCoeffs().coeff(k), &temp.coeffRef(0)); 519c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 520c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 521c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath dec().matrixQR() 522c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath .topLeftCorner(dec().rank(), dec().rank()) 523c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath .template triangularView<Upper>() 524c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath .solveInPlace(c.topRows(dec().rank())); 525c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 526c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath for(Index i = 0; i < dec().rank(); ++i) dst.row(dec().colsPermutation().indices().coeff(i)) = c.row(i); 527c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath for(Index i = dec().rank(); i < cols; ++i) dst.row(dec().colsPermutation().indices().coeff(i)).setZero(); 528c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 529c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath}; 530c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 531c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath/** \ingroup QR_Module 532c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * 533c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \brief Expression type for return value of FullPivHouseholderQR::matrixQ() 534c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * 535c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \tparam MatrixType type of underlying dense matrix 536c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath */ 537c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamathtemplate<typename MatrixType> struct FullPivHouseholderQRMatrixQReturnType 538c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath : public ReturnByValue<FullPivHouseholderQRMatrixQReturnType<MatrixType> > 539c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath{ 540c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamathpublic: 541c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typedef typename MatrixType::Index Index; 5427faaa9f3f0df9d23790277834d426c3d992ac3baCarlos Hernandez typedef typename FullPivHouseholderQR<MatrixType>::IntDiagSizeVectorType IntDiagSizeVectorType; 543c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typedef typename internal::plain_diag_type<MatrixType>::type HCoeffsType; 544c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typedef Matrix<typename MatrixType::Scalar, 1, MatrixType::RowsAtCompileTime, RowMajor, 1, 545c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath MatrixType::MaxRowsAtCompileTime> WorkVectorType; 546c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 547c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath FullPivHouseholderQRMatrixQReturnType(const MatrixType& qr, 548c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath const HCoeffsType& hCoeffs, 5497faaa9f3f0df9d23790277834d426c3d992ac3baCarlos Hernandez const IntDiagSizeVectorType& rowsTranspositions) 550c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath : m_qr(qr), 551c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath m_hCoeffs(hCoeffs), 552c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath m_rowsTranspositions(rowsTranspositions) 553c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath {} 554c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 555c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath template <typename ResultType> 556c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath void evalTo(ResultType& result) const 557c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 558c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath const Index rows = m_qr.rows(); 559c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath WorkVectorType workspace(rows); 560c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath evalTo(result, workspace); 561c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 562c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 563c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath template <typename ResultType> 564c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath void evalTo(ResultType& result, WorkVectorType& workspace) const 565c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 5667faaa9f3f0df9d23790277834d426c3d992ac3baCarlos Hernandez using numext::conj; 567c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath // compute the product H'_0 H'_1 ... H'_n-1, 568c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath // where H_k is the k-th Householder transformation I - h_k v_k v_k' 569c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath // and v_k is the k-th Householder vector [1,m_qr(k+1,k), m_qr(k+2,k), ...] 570c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath const Index rows = m_qr.rows(); 571c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath const Index cols = m_qr.cols(); 572c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath const Index size = (std::min)(rows, cols); 573c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath workspace.resize(rows); 574c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath result.setIdentity(rows, rows); 575c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath for (Index k = size-1; k >= 0; k--) 576c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath { 577c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath result.block(k, k, rows-k, rows-k) 5787faaa9f3f0df9d23790277834d426c3d992ac3baCarlos Hernandez .applyHouseholderOnTheLeft(m_qr.col(k).tail(rows-k-1), conj(m_hCoeffs.coeff(k)), &workspace.coeffRef(k)); 579c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath result.row(k).swap(result.row(m_rowsTranspositions.coeff(k))); 580c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 581c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath } 582c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 583c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath Index rows() const { return m_qr.rows(); } 584c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath Index cols() const { return m_qr.rows(); } 585c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 586c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamathprotected: 587c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typename MatrixType::Nested m_qr; 588c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath typename HCoeffsType::Nested m_hCoeffs; 5897faaa9f3f0df9d23790277834d426c3d992ac3baCarlos Hernandez typename IntDiagSizeVectorType::Nested m_rowsTranspositions; 590c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath}; 591c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 592c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath} // end namespace internal 593c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 594c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamathtemplate<typename MatrixType> 595c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamathinline typename FullPivHouseholderQR<MatrixType>::MatrixQReturnType FullPivHouseholderQR<MatrixType>::matrixQ() const 596c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath{ 597c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath eigen_assert(m_isInitialized && "FullPivHouseholderQR is not initialized."); 598c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath return MatrixQReturnType(m_qr, m_hCoeffs, m_rows_transpositions); 599c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath} 600c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 601c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath/** \return the full-pivoting Householder QR decomposition of \c *this. 602c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * 603c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath * \sa class FullPivHouseholderQR 604c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath */ 605c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamathtemplate<typename Derived> 606c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamathconst FullPivHouseholderQR<typename MatrixBase<Derived>::PlainObject> 607c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan KamathMatrixBase<Derived>::fullPivHouseholderQr() const 608c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath{ 609c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath return FullPivHouseholderQR<PlainObject>(eval()); 610c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath} 611c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 612c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath} // end namespace Eigen 613c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath 614c981c48f5bc9aefeffc0bcb0cc3934c2fae179ddNarayan Kamath#endif // EIGEN_FULLPIVOTINGHOUSEHOLDERQR_H 615