HessenbergDecomposition.h
Go to the documentation of this file.
1 // This file is part of Eigen, a lightweight C++ template library
2 // for linear algebra.
3 //
4 // Copyright (C) 2008-2009 Gael Guennebaud <gael.guennebaud@inria.fr>
5 // Copyright (C) 2010 Jitse Niesen <jitse@maths.leeds.ac.uk>
6 //
7 // This Source Code Form is subject to the terms of the Mozilla
8 // Public License v. 2.0. If a copy of the MPL was not distributed
9 // with this file, You can obtain one at http://mozilla.org/MPL/2.0/.
10 
11 #ifndef EIGEN_HESSENBERGDECOMPOSITION_H
12 #define EIGEN_HESSENBERGDECOMPOSITION_H
13 
14 namespace Eigen {
15 
16 namespace internal {
17 
18 template<typename MatrixType> struct HessenbergDecompositionMatrixHReturnType;
19 template<typename MatrixType>
21 {
22  typedef MatrixType ReturnType;
23 };
24 
25 }
26 
57 template<typename _MatrixType> class HessenbergDecomposition
58 {
59  public:
60 
62  typedef _MatrixType MatrixType;
63 
64  enum {
65  Size = MatrixType::RowsAtCompileTime,
66  SizeMinusOne = Size == Dynamic ? Dynamic : Size - 1,
67  Options = MatrixType::Options,
68  MaxSize = MatrixType::MaxRowsAtCompileTime,
69  MaxSizeMinusOne = MaxSize == Dynamic ? Dynamic : MaxSize - 1
70  };
71 
73  typedef typename MatrixType::Scalar Scalar;
74  typedef typename MatrixType::Index Index;
75 
83 
86 
88 
100  HessenbergDecomposition(Index size = Size==Dynamic ? 2 : Size)
101  : m_matrix(size,size),
102  m_temp(size),
103  m_isInitialized(false)
104  {
105  if(size>1)
106  m_hCoeffs.resize(size-1);
107  }
108 
118  HessenbergDecomposition(const MatrixType& matrix)
119  : m_matrix(matrix),
120  m_temp(matrix.rows()),
121  m_isInitialized(false)
122  {
123  if(matrix.rows()<2)
124  {
125  m_isInitialized = true;
126  return;
127  }
128  m_hCoeffs.resize(matrix.rows()-1,1);
129  _compute(m_matrix, m_hCoeffs, m_temp);
130  m_isInitialized = true;
131  }
132 
150  HessenbergDecomposition& compute(const MatrixType& matrix)
151  {
152  m_matrix = matrix;
153  if(matrix.rows()<2)
154  {
155  m_isInitialized = true;
156  return *this;
157  }
158  m_hCoeffs.resize(matrix.rows()-1,1);
159  _compute(m_matrix, m_hCoeffs, m_temp);
160  m_isInitialized = true;
161  return *this;
162  }
163 
177  const CoeffVectorType& householderCoefficients() const
178  {
179  eigen_assert(m_isInitialized && "HessenbergDecomposition is not initialized.");
180  return m_hCoeffs;
181  }
182 
212  const MatrixType& packedMatrix() const
213  {
214  eigen_assert(m_isInitialized && "HessenbergDecomposition is not initialized.");
215  return m_matrix;
216  }
217 
232  HouseholderSequenceType matrixQ() const
233  {
234  eigen_assert(m_isInitialized && "HessenbergDecomposition is not initialized.");
235  return HouseholderSequenceType(m_matrix, m_hCoeffs.conjugate())
236  .setLength(m_matrix.rows() - 1)
237  .setShift(1);
238  }
239 
260  MatrixHReturnType matrixH() const
261  {
262  eigen_assert(m_isInitialized && "HessenbergDecomposition is not initialized.");
263  return MatrixHReturnType(*this);
264  }
265 
266  private:
267 
270  static void _compute(MatrixType& matA, CoeffVectorType& hCoeffs, VectorType& temp);
271 
272  protected:
273  MatrixType m_matrix;
274  CoeffVectorType m_hCoeffs;
275  VectorType m_temp;
277 };
278 
291 template<typename MatrixType>
293 {
294  eigen_assert(matA.rows()==matA.cols());
295  Index n = matA.rows();
296  temp.resize(n);
297  for (Index i = 0; i<n-1; ++i)
298  {
299  // let's consider the vector v = i-th column starting at position i+1
300  Index remainingSize = n-i-1;
301  RealScalar beta;
302  Scalar h;
303  matA.col(i).tail(remainingSize).makeHouseholderInPlace(h, beta);
304  matA.col(i).coeffRef(i+1) = beta;
305  hCoeffs.coeffRef(i) = h;
306 
307  // Apply similarity transformation to remaining columns,
308  // i.e., compute A = H A H'
309 
310  // A = H A
311  matA.bottomRightCorner(remainingSize, remainingSize)
312  .applyHouseholderOnTheLeft(matA.col(i).tail(remainingSize-1), h, &temp.coeffRef(0));
313 
314  // A = A H'
315  matA.rightCols(remainingSize)
316  .applyHouseholderOnTheRight(matA.col(i).tail(remainingSize-1).conjugate(), numext::conj(h), &temp.coeffRef(0));
317  }
318 }
319 
320 namespace internal {
321 
337 template<typename MatrixType> struct HessenbergDecompositionMatrixHReturnType
338 : public ReturnByValue<HessenbergDecompositionMatrixHReturnType<MatrixType> >
339 {
340  typedef typename MatrixType::Index Index;
341  public:
347 
353  template <typename ResultType>
354  inline void evalTo(ResultType& result) const
355  {
356  result = m_hess.packedMatrix();
357  Index n = result.rows();
358  if (n>2)
359  result.bottomLeftCorner(n-2, n-2).template triangularView<Lower>().setZero();
360  }
361 
362  Index rows() const { return m_hess.packedMatrix().rows(); }
363  Index cols() const { return m_hess.packedMatrix().cols(); }
364 
365  protected:
367 };
368 
369 } // end namespace internal
370 
371 } // end namespace Eigen
372 
373 #endif // EIGEN_HESSENBERGDECOMPOSITION_H
const MatrixType & packedMatrix() const
Returns the internal representation of the decomposition.
HouseholderSequenceType matrixQ() const
Reconstructs the orthogonal matrix Q in the decomposition.
const AutoDiffScalar< DerType > & conj(const AutoDiffScalar< DerType > &x)
MatrixHReturnType matrixH() const
Constructs the Hessenberg matrix H in the decomposition.
static void _compute(MatrixType &matA, CoeffVectorType &hCoeffs, VectorType &temp)
iterative scaling algorithm to equilibrate rows and column norms in matrices
Definition: matrix.hpp:471
Holds information about the various numeric (i.e. scalar) types allowed by Eigen. ...
Definition: NumTraits.h:88
internal::HessenbergDecompositionMatrixHReturnType< MatrixType > MatrixHReturnType
_MatrixType MatrixType
Synonym for the template parameter _MatrixType.
MatrixType::Scalar Scalar
Scalar type for matrices of type MatrixType.
Sequence of Householder reflections acting on subspaces with decreasing size.
Matrix< Scalar, 1, Size, Options|RowMajor, 1, MaxSize > VectorType
void evalTo(ResultType &result) const
Hessenberg matrix in decomposition.
NumTraits< Scalar >::Real RealScalar
HessenbergDecomposition & compute(const MatrixType &matrix)
Computes Hessenberg decomposition of given matrix.
EIGEN_STRONG_INLINE Scalar & coeffRef(Index rowId, Index colId)
HessenbergDecompositionMatrixHReturnType(const HessenbergDecomposition< MatrixType > &hess)
Constructor.
EIGEN_STRONG_INLINE void resize(Index nbRows, Index nbCols)
Provides a generic way to set and pass user-specified options.
Definition: options.hpp:65
const HessenbergDecomposition< MatrixType > & m_hess
HessenbergDecomposition(const MatrixType &matrix)
Constructor; computes Hessenberg decomposition of given matrix.
Expression type for return value of HessenbergDecomposition::matrixH()
HessenbergDecomposition(Index size=Size==Dynamic?2:Size)
Default constructor; the decomposition will be computed later.
Reduces a square matrix to Hessenberg form by an orthogonal similarity transformation.
HouseholderSequence< MatrixType, typename internal::remove_all< typename CoeffVectorType::ConjugateReturnType >::type > HouseholderSequenceType
Return type of matrixQ()
#define eigen_assert(x)
Matrix< Scalar, SizeMinusOne, 1, Options &~RowMajor, MaxSizeMinusOne, 1 > CoeffVectorType
Type for vector of Householder coefficients.
const CoeffVectorType & householderCoefficients() const
Returns the Householder coefficients.


acado
Author(s): Milan Vukov, Rien Quirynen
autogenerated on Mon Jun 10 2019 12:34:40