aboutsummaryrefslogtreecommitdiffhomepage
diff options
context:
space:
mode:
-rw-r--r--Eigen/src/Core/util/ForwardDeclarations.h2
-rw-r--r--Eigen/src/Householder/Householder.h2
-rw-r--r--Eigen/src/Householder/HouseholderSequence.h80
-rw-r--r--test/householder.cpp22
4 files changed, 69 insertions, 37 deletions
diff --git a/Eigen/src/Core/util/ForwardDeclarations.h b/Eigen/src/Core/util/ForwardDeclarations.h
index 6959a1c1b..b6bba04e6 100644
--- a/Eigen/src/Core/util/ForwardDeclarations.h
+++ b/Eigen/src/Core/util/ForwardDeclarations.h
@@ -134,7 +134,7 @@ template<typename MatrixType> class SVD;
template<typename MatrixType, unsigned int Options = 0> class JacobiSVD;
template<typename MatrixType, int UpLo = Lower> class LLT;
template<typename MatrixType> class LDLT;
-template<typename VectorsType, typename CoeffsType> class HouseholderSequence;
+template<typename VectorsType, typename CoeffsType, int Side=OnTheLeft> class HouseholderSequence;
template<typename Scalar> class PlanarRotation;
// Geometry module:
diff --git a/Eigen/src/Householder/Householder.h b/Eigen/src/Householder/Householder.h
index 1e549633a..cfd3935fc 100644
--- a/Eigen/src/Householder/Householder.h
+++ b/Eigen/src/Householder/Householder.h
@@ -1,7 +1,7 @@
// This file is part of Eigen, a lightweight C++ template library
// for linear algebra.
//
-// Copyright (C) 2009 Benoit Jacob <jacob.benoit.1@gmail.com>
+// Copyright (C) 2010 Benoit Jacob <jacob.benoit.1@gmail.com>
// Copyright (C) 2009 Gael Guennebaud <g.gael@free.fr>
//
// Eigen is free software; you can redistribute it and/or
diff --git a/Eigen/src/Householder/HouseholderSequence.h b/Eigen/src/Householder/HouseholderSequence.h
index 7bd403373..61a8e8c1d 100644
--- a/Eigen/src/Householder/HouseholderSequence.h
+++ b/Eigen/src/Householder/HouseholderSequence.h
@@ -2,6 +2,7 @@
// for linear algebra.
//
// Copyright (C) 2009 Gael Guennebaud <g.gael@free.fr>
+// Copyright (C) 2010 Benoit Jacob <jacob.benoit.1@gmail.com>
//
// Eigen is free software; you can redistribute it and/or
// modify it under the terms of the GNU Lesser General Public
@@ -48,31 +49,61 @@
* \sa MatrixBase::applyOnTheLeft(), MatrixBase::applyOnTheRight()
*/
-template<typename VectorsType, typename CoeffsType>
-struct ei_traits<HouseholderSequence<VectorsType,CoeffsType> >
+template<typename VectorsType, typename CoeffsType, int Side>
+struct ei_traits<HouseholderSequence<VectorsType,CoeffsType,Side> >
{
typedef typename VectorsType::Scalar Scalar;
enum {
- RowsAtCompileTime = ei_traits<VectorsType>::RowsAtCompileTime,
- ColsAtCompileTime = ei_traits<VectorsType>::RowsAtCompileTime,
- MaxRowsAtCompileTime = ei_traits<VectorsType>::MaxRowsAtCompileTime,
- MaxColsAtCompileTime = ei_traits<VectorsType>::MaxRowsAtCompileTime,
+ RowsAtCompileTime = Side==OnTheLeft ? ei_traits<VectorsType>::RowsAtCompileTime
+ : ei_traits<VectorsType>::ColsAtCompileTime,
+ ColsAtCompileTime = RowsAtCompileTime,
+ MaxRowsAtCompileTime = Side==OnTheLeft ? ei_traits<VectorsType>::MaxRowsAtCompileTime
+ : ei_traits<VectorsType>::MaxColsAtCompileTime,
+ MaxColsAtCompileTime = MaxRowsAtCompileTime,
Flags = 0
};
};
-template<typename VectorsType, typename CoeffsType> class HouseholderSequence
- : public AnyMatrixBase<HouseholderSequence<VectorsType,CoeffsType> >
+template<typename VectorsType, typename CoeffsType, int Side>
+struct ei_hseq_side_dependent_impl
+{
+ typedef Block<VectorsType, Dynamic, 1> EssentialVectorType;
+ typedef HouseholderSequence<VectorsType, CoeffsType, OnTheLeft> HouseholderSequenceType;
+ static inline const EssentialVectorType essentialVector(const HouseholderSequenceType& h, int k)
+ {
+ const int start = k+1+h.m_shift;
+ return Block<VectorsType,Dynamic,1>(h.m_vectors, start, k, h.rows()-start, 1);
+ }
+};
+
+template<typename VectorsType, typename CoeffsType>
+struct ei_hseq_side_dependent_impl<VectorsType, CoeffsType, OnTheRight>
+{
+ typedef Transpose<Block<VectorsType, 1, Dynamic> > EssentialVectorType;
+ typedef HouseholderSequence<VectorsType, CoeffsType, OnTheRight> HouseholderSequenceType;
+ static inline const EssentialVectorType essentialVector(const HouseholderSequenceType& h, int k)
+ {
+ const int start = k+1+h.m_shift;
+ return Block<VectorsType,1,Dynamic>(h.m_vectors, k, start, 1, h.rows()-start).transpose();
+ }
+};
+
+template<typename VectorsType, typename CoeffsType, int Side> class HouseholderSequence
+ : public AnyMatrixBase<HouseholderSequence<VectorsType,CoeffsType,Side> >
{
typedef typename VectorsType::Scalar Scalar;
- typedef Block<VectorsType, Dynamic, 1> EssentialVectorType;
+ typedef typename ei_hseq_side_dependent_impl<VectorsType,CoeffsType,Side>::EssentialVectorType
+ EssentialVectorType;
public:
- typedef HouseholderSequence<VectorsType,
+ typedef HouseholderSequence<
+ VectorsType,
typename ei_meta_if<NumTraits<Scalar>::IsComplex,
typename ei_cleantype<typename CoeffsType::ConjugateReturnType>::type,
- CoeffsType>::ret> ConjugateReturnType;
+ CoeffsType>::ret,
+ 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()),
@@ -85,14 +116,13 @@ template<typename VectorsType, typename CoeffsType> class HouseholderSequence
{
}
- int rows() const { return m_vectors.rows(); }
- int cols() const { return m_vectors.rows(); }
+ int rows() const { return Side==OnTheLeft ? m_vectors.rows() : m_vectors.cols(); }
+ int cols() const { return rows(); }
const EssentialVectorType essentialVector(int k) const
{
- ei_assert(k>= 0 && k < m_actualVectors);
- const int start = k+1+m_shift;
- return Block<VectorsType,Dynamic,1>(m_vectors, start, k, rows()-start, 1);
+ ei_assert(k >= 0 && k < m_actualVectors);
+ return ei_hseq_side_dependent_impl<VectorsType,CoeffsType,Side>::essentialVector(*this, k);
}
HouseholderSequence transpose() const
@@ -164,6 +194,8 @@ template<typename VectorsType, typename CoeffsType> class HouseholderSequence
return res;
}
+ template<typename _VectorsType, typename _CoeffsType, int _Side> friend struct ei_hseq_side_dependent_impl;
+
protected:
typename VectorsType::Nested m_vectors;
typename CoeffsType::Nested m_coeffs;
@@ -175,13 +207,25 @@ template<typename VectorsType, typename CoeffsType> class HouseholderSequence
template<typename VectorsType, typename CoeffsType>
HouseholderSequence<VectorsType,CoeffsType> householderSequence(const VectorsType& v, const CoeffsType& h, bool trans=false)
{
- return HouseholderSequence<VectorsType,CoeffsType>(v, h, trans);
+ return HouseholderSequence<VectorsType,CoeffsType,OnTheLeft>(v, h, trans);
}
template<typename VectorsType, typename CoeffsType>
HouseholderSequence<VectorsType,CoeffsType> householderSequence(const VectorsType& v, const CoeffsType& h, bool trans, int actualVectors, int shift)
{
- return HouseholderSequence<VectorsType,CoeffsType>(v, h, trans, actualVectors, shift);
+ return HouseholderSequence<VectorsType,CoeffsType,OnTheLeft>(v, h, trans, actualVectors, shift);
+}
+
+template<typename VectorsType, typename CoeffsType>
+HouseholderSequence<VectorsType,CoeffsType> rightHouseholderSequence(const VectorsType& v, const CoeffsType& h, bool trans=false)
+{
+ return HouseholderSequence<VectorsType,CoeffsType,OnTheRight>(v, h, trans);
+}
+
+template<typename VectorsType, typename CoeffsType>
+HouseholderSequence<VectorsType,CoeffsType> rightHouseholderSequence(const VectorsType& v, const CoeffsType& h, bool trans, int actualVectors, int shift)
+{
+ return HouseholderSequence<VectorsType,CoeffsType,OnTheRight>(v, h, trans, actualVectors, shift);
}
#endif // EIGEN_HOUSEHOLDER_SEQUENCE_H
diff --git a/test/householder.cpp b/test/householder.cpp
index c4c95e648..d98780872 100644
--- a/test/householder.cpp
+++ b/test/householder.cpp
@@ -106,27 +106,15 @@ template<typename MatrixType> void householder(const MatrixType& m)
m5.block(shift,0,brows,cols).template triangularView<StrictlyLower>().setZero();
VERIFY_IS_APPROX(hseq * m5, m1); // test applying hseq directly
m3 = hseq;
- VERIFY_IS_APPROX(m3*m5, m1); // test evaluating hseq to a dense matrix, then applying
+ VERIFY_IS_APPROX(m3 * m5, m1); // test evaluating hseq to a dense matrix, then applying
-#if 0
// test householder sequence on the right with a shift
- TMatrixType tm1 = m1.transpose();
TMatrixType tm2 = m2.transpose();
-
- int bcols = cols - shift;
- VBlockMatrixType vbm =
- HouseholderQR<HBlockMatrixType> qr(hbm);
- m2 = m1;
- m2.block(shift,0,brows,cols) = qr.matrixQR();
- HCoeffsVectorType hc = qr.hCoeffs().conjugate();
- HouseholderSequence<MatrixType, HCoeffsVectorType> hseq(m2, hc, false, hc.size(), shift);
- MatrixType m5 = m2;
- m5.block(shift,0,brows,cols).template triangularView<StrictlyLower>().setZero();
- VERIFY_IS_APPROX(hseq * m5, m1); // test applying hseq directly
- m3 = hseq;
- VERIFY_IS_APPROX(m3*m5, m1); // test evaluating hseq to a dense matrix, then applying
-#endif
+ HouseholderSequence<TMatrixType, HCoeffsVectorType, OnTheRight> rhseq(tm2, hc, false, hc.size(), shift);
+ VERIFY_IS_APPROX(rhseq * m5, m1); // test applying rhseq directly
+ m3 = rhseq;
+ VERIFY_IS_APPROX(m3 * m5, m1); // test evaluating rhseq to a dense matrix, then applying
}
void test_householder()