From dbd9c5fd50cde5d5beaae44147eca3ba11934721 Mon Sep 17 00:00:00 2001 From: Benoit Jacob Date: Thu, 30 Dec 2010 04:18:40 -0500 Subject: fix HouseholderSequence API, bug #50: * remove ctors taking more than 2 ints * rename actualVectors to length * add length/shift/trans accessors/mutators --- Eigen/src/Eigenvalues/HessenbergDecomposition.h | 4 +- Eigen/src/Eigenvalues/Tridiagonalization.h | 8 ++- Eigen/src/Householder/HouseholderSequence.h | 87 +++++++++++++++---------- Eigen/src/QR/ColPivHouseholderQR.h | 13 ++-- Eigen/src/SVD/UpperBidiagonalization.h | 5 +- 5 files changed, 70 insertions(+), 47 deletions(-) (limited to 'Eigen') diff --git a/Eigen/src/Eigenvalues/HessenbergDecomposition.h b/Eigen/src/Eigenvalues/HessenbergDecomposition.h index 9333d81a2..c17f155a5 100644 --- a/Eigen/src/Eigenvalues/HessenbergDecomposition.h +++ b/Eigen/src/Eigenvalues/HessenbergDecomposition.h @@ -245,7 +245,9 @@ template class HessenbergDecomposition HouseholderSequenceType matrixQ() const { eigen_assert(m_isInitialized && "HessenbergDecomposition is not initialized."); - return HouseholderSequenceType(m_matrix, m_hCoeffs.conjugate(), false, m_matrix.rows() - 1, 1); + return HouseholderSequenceType(m_matrix, m_hCoeffs.conjugate()) + .setLength(m_matrix.rows() - 1) + .setShift(1); } /** \brief Constructs the Hessenberg matrix H in the decomposition diff --git a/Eigen/src/Eigenvalues/Tridiagonalization.h b/Eigen/src/Eigenvalues/Tridiagonalization.h index 2296f0a9f..755bca1aa 100644 --- a/Eigen/src/Eigenvalues/Tridiagonalization.h +++ b/Eigen/src/Eigenvalues/Tridiagonalization.h @@ -251,7 +251,9 @@ template class Tridiagonalization HouseholderSequenceType matrixQ() const { eigen_assert(m_isInitialized && "Tridiagonalization is not initialized."); - return HouseholderSequenceType(m_matrix, m_hCoeffs.conjugate(), false, m_matrix.rows() - 1, 1); + return HouseholderSequenceType(m_matrix, m_hCoeffs.conjugate()) + .setLength(m_matrix.rows() - 1) + .setShift(1); } /** \brief Returns an expression of the tridiagonal matrix T in the decomposition @@ -459,7 +461,9 @@ struct tridiagonalization_inplace_selector diag = mat.diagonal().real(); subdiag = mat.template diagonal<-1>().real(); if(extractQ) - mat = HouseholderSequenceType(mat, hCoeffs.conjugate(), false, mat.rows() - 1, 1); + mat = HouseholderSequenceType(mat, hCoeffs.conjugate()) + .setLength(mat.rows() - 1) + .setShift(1); } }; diff --git a/Eigen/src/Householder/HouseholderSequence.h b/Eigen/src/Householder/HouseholderSequence.h index d6260cf77..d3616ed70 100644 --- a/Eigen/src/Householder/HouseholderSequence.h +++ b/Eigen/src/Householder/HouseholderSequence.h @@ -129,14 +129,18 @@ template class HouseholderS Side > ConjugateReturnType; - HouseholderSequence(const VectorsType& v, const CoeffsType& h, bool trans = false) - : m_vectors(v), m_coeffs(h), m_trans(trans), m_actualVectors(v.diagonalSize()), + HouseholderSequence(const VectorsType& v, const CoeffsType& h) + : m_vectors(v), m_coeffs(h), m_trans(false), m_length(v.diagonalSize()), m_shift(0) { } - HouseholderSequence(const VectorsType& v, const CoeffsType& h, bool trans, Index actualVectors, Index shift) - : m_vectors(v), m_coeffs(h), m_trans(trans), m_actualVectors(actualVectors), m_shift(shift) + HouseholderSequence(const HouseholderSequence& other) + : m_vectors(other.m_vectors), + m_coeffs(other.m_coeffs), + m_trans(other.m_trans), + m_length(other.m_length), + m_shift(other.m_shift) { } @@ -145,25 +149,34 @@ template class HouseholderS const EssentialVectorType essentialVector(Index k) const { - eigen_assert(k >= 0 && k < m_actualVectors); + eigen_assert(k >= 0 && k < m_length); return internal::hseq_side_dependent_impl::essentialVector(*this, k); } HouseholderSequence transpose() const - { return HouseholderSequence(m_vectors, m_coeffs, !m_trans, m_actualVectors, m_shift); } + { + return HouseholderSequence(*this).setTrans(!m_trans); + } ConjugateReturnType conjugate() const - { return ConjugateReturnType(m_vectors, m_coeffs.conjugate(), m_trans, m_actualVectors, m_shift); } + { + return ConjugateReturnType(m_vectors, m_coeffs.conjugate()) + .setTrans(m_trans) + .setLength(m_length) + .setShift(m_shift); + } ConjugateReturnType adjoint() const - { return ConjugateReturnType(m_vectors, m_coeffs.conjugate(), !m_trans, m_actualVectors, m_shift); } + { + return conjugate().setTrans(!m_trans); + } ConjugateReturnType inverse() const { return adjoint(); } /** \internal */ template void evalTo(DestType& dst) const { - Index vecs = m_actualVectors; + Index vecs = m_length; // FIXME find a way to pass this temporary if the user wants to Matrix temp(rows()); @@ -210,9 +223,9 @@ template class HouseholderS template inline void applyThisOnTheRight(Dest& dst) const { Matrix temp(dst.rows()); - for(Index k = 0; k < m_actualVectors; ++k) + for(Index k = 0; k < m_length; ++k) { - Index actual_k = m_trans ? m_actualVectors-k-1 : k; + Index actual_k = m_trans ? m_length-k-1 : k; dst.rightCols(rows()-m_shift-actual_k) .applyHouseholderOnTheRight(essentialVector(actual_k), m_coeffs.coeff(actual_k), &temp.coeffRef(0)); } @@ -222,9 +235,9 @@ template class HouseholderS template inline void applyThisOnTheLeft(Dest& dst) const { Matrix temp(dst.cols()); - for(Index k = 0; k < m_actualVectors; ++k) + for(Index k = 0; k < m_length; ++k) { - Index actual_k = m_trans ? k : m_actualVectors-k-1; + Index actual_k = m_trans ? k : m_length-k-1; dst.bottomRows(rows()-m_shift-actual_k) .applyHouseholderOnTheLeft(essentialVector(actual_k), m_coeffs.coeff(actual_k), &temp.coeffRef(0)); } @@ -250,40 +263,46 @@ template class HouseholderS template friend struct internal::hseq_side_dependent_impl; + HouseholderSequence& setTrans(bool trans) + { + m_trans = trans; + return *this; + } + + HouseholderSequence& setLength(Index length) + { + m_length = length; + return *this; + } + + HouseholderSequence& setShift(Index shift) + { + m_shift = shift; + return *this; + } + + bool trans() const { return m_trans; } + Index length() const { return m_length; } + Index shift() const { return m_shift; } + protected: typename VectorsType::Nested m_vectors; typename CoeffsType::Nested m_coeffs; bool m_trans; - Index m_actualVectors; + Index m_length; Index m_shift; }; template -HouseholderSequence householderSequence(const VectorsType& v, const CoeffsType& h, bool trans=false) -{ - return HouseholderSequence(v, h, trans); -} - -template -HouseholderSequence householderSequence - (const VectorsType& v, const CoeffsType& h, - bool trans, typename VectorsType::Index actualVectors, typename VectorsType::Index shift) -{ - return HouseholderSequence(v, h, trans, actualVectors, shift); -} - -template -HouseholderSequence rightHouseholderSequence(const VectorsType& v, const CoeffsType& h, bool trans=false) +HouseholderSequence householderSequence(const VectorsType& v, const CoeffsType& h) { - return HouseholderSequence(v, h, trans); + return HouseholderSequence(v, h); } template -HouseholderSequence rightHouseholderSequence - (const VectorsType& v, const CoeffsType& h, bool trans, - typename VectorsType::Index actualVectors, typename VectorsType::Index shift) +HouseholderSequence rightHouseholderSequence(const VectorsType& v, const CoeffsType& h) { - return HouseholderSequence(v, h, trans, actualVectors, shift); + return HouseholderSequence(v, h); } #endif // EIGEN_HOUSEHOLDER_SEQUENCE_H diff --git a/Eigen/src/QR/ColPivHouseholderQR.h b/Eigen/src/QR/ColPivHouseholderQR.h index 21ad0febe..c8ecf2c43 100644 --- a/Eigen/src/QR/ColPivHouseholderQR.h +++ b/Eigen/src/QR/ColPivHouseholderQR.h @@ -483,13 +483,10 @@ struct solve_retval, Rhs> typename Rhs::PlainObject c(rhs()); // Note that the matrix Q = H_0^* H_1^*... so its inverse is Q^* = (H_0 H_1 ...)^T - c.applyOnTheLeft(householderSequence( - dec().matrixQR(), - dec().hCoeffs(), - true, - dec().nonzeroPivots(), - 0 - )); + c.applyOnTheLeft(householderSequence(dec().matrixQR(), dec().hCoeffs()) + .setTrans(true) + .setLength(dec().nonzeroPivots()) + ); dec().matrixQR() .topLeftCorner(nonzero_pivots, nonzero_pivots) @@ -517,7 +514,7 @@ typename ColPivHouseholderQR::HouseholderSequenceType ColPivHousehol ::householderQ() const { eigen_assert(m_isInitialized && "ColPivHouseholderQR is not initialized."); - return HouseholderSequenceType(m_qr, m_hCoeffs.conjugate(), false, m_nonzero_pivots, 0); + return HouseholderSequenceType(m_qr, m_hCoeffs.conjugate()).setLength(m_nonzero_pivots); } /** \return the column-pivoting Householder QR decomposition of \c *this. diff --git a/Eigen/src/SVD/UpperBidiagonalization.h b/Eigen/src/SVD/UpperBidiagonalization.h index 5503356d3..eef92fcbe 100644 --- a/Eigen/src/SVD/UpperBidiagonalization.h +++ b/Eigen/src/SVD/UpperBidiagonalization.h @@ -87,8 +87,9 @@ template class UpperBidiagonalization const HouseholderVSequenceType householderV() // const here gives nasty errors and i'm lazy { eigen_assert(m_isInitialized && "UpperBidiagonalization is not initialized."); - return HouseholderVSequenceType(m_householder, m_householder.const_derived().template diagonal<1>(), - false, m_householder.cols()-1, 1); + return HouseholderVSequenceType(m_householder, m_householder.const_derived().template diagonal<1>()) + .setLength(m_householder.cols()-1) + .setShift(1); } protected: -- cgit v1.2.3