GeneralizedSelfAdjointEigenSolver.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-2010 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_GENERALIZEDSELFADJOINTEIGENSOLVER_H
12 #define EIGEN_GENERALIZEDSELFADJOINTEIGENSOLVER_H
13 
14 #include "./Tridiagonalization.h"
15 
16 namespace Eigen {
17 
47 template<typename _MatrixType>
49 {
51  public:
52 
53  typedef typename Base::Index Index;
54  typedef _MatrixType MatrixType;
55 
64 
78  : Base(size)
79  {}
80 
107  GeneralizedSelfAdjointEigenSolver(const MatrixType& matA, const MatrixType& matB,
109  : Base(matA.cols())
110  {
111  compute(matA, matB, options);
112  }
113 
154  GeneralizedSelfAdjointEigenSolver& compute(const MatrixType& matA, const MatrixType& matB,
156 
157  protected:
158 
159 };
160 
161 
162 template<typename MatrixType>
164 compute(const MatrixType& matA, const MatrixType& matB, int options)
165 {
166  eigen_assert(matA.cols()==matA.rows() && matB.rows()==matA.rows() && matB.cols()==matB.rows());
167  eigen_assert((options&~(EigVecMask|GenEigMask))==0
168  && (options&EigVecMask)!=EigVecMask
169  && ((options&GenEigMask)==0 || (options&GenEigMask)==Ax_lBx
170  || (options&GenEigMask)==ABx_lx || (options&GenEigMask)==BAx_lx)
171  && "invalid option parameter");
172 
173  bool computeEigVecs = ((options&EigVecMask)==0) || ((options&EigVecMask)==ComputeEigenvectors);
174 
175  // Compute the cholesky decomposition of matB = L L' = U'U
176  LLT<MatrixType> cholB(matB);
177 
178  int type = (options&GenEigMask);
179  if(type==0)
180  type = Ax_lBx;
181 
182  if(type==Ax_lBx)
183  {
184  // compute C = inv(L) A inv(L')
185  MatrixType matC = matA.template selfadjointView<Lower>();
186  cholB.matrixL().template solveInPlace<OnTheLeft>(matC);
187  cholB.matrixU().template solveInPlace<OnTheRight>(matC);
188 
189  Base::compute(matC, computeEigVecs ? ComputeEigenvectors : EigenvaluesOnly );
190 
191  // transform back the eigen vectors: evecs = inv(U) * evecs
192  if(computeEigVecs)
193  cholB.matrixU().solveInPlace(Base::m_eivec);
194  }
195  else if(type==ABx_lx)
196  {
197  // compute C = L' A L
198  MatrixType matC = matA.template selfadjointView<Lower>();
199  matC = matC * cholB.matrixL();
200  matC = cholB.matrixU() * matC;
201 
202  Base::compute(matC, computeEigVecs ? ComputeEigenvectors : EigenvaluesOnly);
203 
204  // transform back the eigen vectors: evecs = inv(U) * evecs
205  if(computeEigVecs)
206  cholB.matrixU().solveInPlace(Base::m_eivec);
207  }
208  else if(type==BAx_lx)
209  {
210  // compute C = L' A L
211  MatrixType matC = matA.template selfadjointView<Lower>();
212  matC = matC * cholB.matrixL();
213  matC = cholB.matrixU() * matC;
214 
215  Base::compute(matC, computeEigVecs ? ComputeEigenvectors : EigenvaluesOnly);
216 
217  // transform back the eigen vectors: evecs = L * evecs
218  if(computeEigVecs)
219  Base::m_eivec = cholB.matrixL() * Base::m_eivec;
220  }
221 
222  return *this;
223 }
224 
225 } // end namespace Eigen
226 
227 #endif // EIGEN_GENERALIZEDSELFADJOINTEIGENSOLVER_H
GeneralizedSelfAdjointEigenSolver(Index size)
Constructor, pre-allocates memory for dynamic-size matrices.
Traits::MatrixU matrixU() const
Definition: LLT.h:97
Computes eigenvalues and eigenvectors of selfadjoint matrices.
iterative scaling algorithm to equilibrate rows and column norms in matrices
Definition: matrix.hpp:471
Standard Cholesky decomposition (LL^T) of a matrix and associated features.
Definition: LLT.h:50
Computes eigenvalues and eigenvectors of the generalized selfadjoint eigen problem.
SelfAdjointEigenSolver & compute(const MatrixType &matrix, int options=ComputeEigenvectors)
Computes eigendecomposition of given matrix.
GeneralizedSelfAdjointEigenSolver & compute(const MatrixType &matA, const MatrixType &matB, int options=ComputeEigenvectors|Ax_lBx)
Computes generalized eigendecomposition of given matrix pencil.
#define eigen_assert(x)
Traits::MatrixL matrixL() const
Definition: LLT.h:104
GeneralizedSelfAdjointEigenSolver(const MatrixType &matA, const MatrixType &matB, int options=ComputeEigenvectors|Ax_lBx)
Constructor; computes generalized eigendecomposition of given matrix pencil.
GeneralizedSelfAdjointEigenSolver()
Default constructor for fixed-size matrices.


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