Eigen  3.4.90 (git rev a4098ac676528a83cfb73d4d26ce1b42ec05f47c)
GeneralizedSelfAdjointEigenSolver.h
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#include "./InternalHeaderCheck.h"
17
18namespace Eigen {
19
49template<typename MatrixType_>
51{
53 public:
54
55 typedef MatrixType_ MatrixType;
56
65
79 : Base(size)
80 {}
81
108 GeneralizedSelfAdjointEigenSolver(const MatrixType& matA, const MatrixType& matB,
109 int options = ComputeEigenvectors|Ax_lBx)
110 : Base(matA.cols())
111 {
112 compute(matA, matB, options);
113 }
114
155 GeneralizedSelfAdjointEigenSolver& compute(const MatrixType& matA, const MatrixType& matB,
156 int options = ComputeEigenvectors|Ax_lBx);
157
158 protected:
159
160};
161
162
163template<typename MatrixType>
165compute(const MatrixType& matA, const MatrixType& matB, int options)
166{
167 eigen_assert(matA.cols()==matA.rows() && matB.rows()==matA.rows() && matB.cols()==matB.rows());
168 eigen_assert((options&~(EigVecMask|GenEigMask))==0
169 && (options&EigVecMask)!=EigVecMask
170 && ((options&GenEigMask)==0 || (options&GenEigMask)==Ax_lBx
171 || (options&GenEigMask)==ABx_lx || (options&GenEigMask)==BAx_lx)
172 && "invalid option parameter");
173
174 bool computeEigVecs = ((options&EigVecMask)==0) || ((options&EigVecMask)==ComputeEigenvectors);
175
176 // Compute the cholesky decomposition of matB = L L' = U'U
177 LLT<MatrixType> cholB(matB);
178
179 int type = (options&GenEigMask);
180 if(type==0)
181 type = Ax_lBx;
182
183 if(type==Ax_lBx)
184 {
185 // compute C = inv(L) A inv(L')
186 MatrixType matC = matA.template selfadjointView<Lower>();
187 cholB.matrixL().template solveInPlace<OnTheLeft>(matC);
188 cholB.matrixU().template solveInPlace<OnTheRight>(matC);
189
190 Base::compute(matC, computeEigVecs ? ComputeEigenvectors : EigenvaluesOnly );
191
192 // transform back the eigen vectors: evecs = inv(U) * evecs
193 if(computeEigVecs)
194 cholB.matrixU().solveInPlace(Base::m_eivec);
195 }
196 else if(type==ABx_lx)
197 {
198 // compute C = L' A L
199 MatrixType matC = matA.template selfadjointView<Lower>();
200 matC = matC * cholB.matrixL();
201 matC = cholB.matrixU() * matC;
202
203 Base::compute(matC, computeEigVecs ? ComputeEigenvectors : EigenvaluesOnly);
204
205 // transform back the eigen vectors: evecs = inv(U) * evecs
206 if(computeEigVecs)
207 cholB.matrixU().solveInPlace(Base::m_eivec);
208 }
209 else if(type==BAx_lx)
210 {
211 // compute C = L' A L
212 MatrixType matC = matA.template selfadjointView<Lower>();
213 matC = matC * cholB.matrixL();
214 matC = cholB.matrixU() * matC;
215
216 Base::compute(matC, computeEigVecs ? ComputeEigenvectors : EigenvaluesOnly);
217
218 // transform back the eigen vectors: evecs = L * evecs
219 if(computeEigVecs)
220 Base::m_eivec = cholB.matrixL() * Base::m_eivec;
221 }
222
223 return *this;
224}
225
226} // end namespace Eigen
227
228#endif // EIGEN_GENERALIZEDSELFADJOINTEIGENSOLVER_H
Computes eigenvalues and eigenvectors of the generalized selfadjoint eigen problem.
Definition: GeneralizedSelfAdjointEigenSolver.h:51
GeneralizedSelfAdjointEigenSolver()
Default constructor for fixed-size matrices.
Definition: GeneralizedSelfAdjointEigenSolver.h:64
GeneralizedSelfAdjointEigenSolver(Index size)
Constructor, pre-allocates memory for dynamic-size matrices.
Definition: GeneralizedSelfAdjointEigenSolver.h:78
GeneralizedSelfAdjointEigenSolver & compute(const MatrixType &matA, const MatrixType &matB, int options=ComputeEigenvectors|Ax_lBx)
Computes generalized eigendecomposition of given matrix pencil.
Definition: GeneralizedSelfAdjointEigenSolver.h:165
GeneralizedSelfAdjointEigenSolver(const MatrixType &matA, const MatrixType &matB, int options=ComputeEigenvectors|Ax_lBx)
Constructor; computes generalized eigendecomposition of given matrix pencil.
Definition: GeneralizedSelfAdjointEigenSolver.h:108
Standard Cholesky decomposition (LL^T) of a matrix and associated features.
Definition: LLT.h:70
Traits::MatrixU matrixU() const
Definition: LLT.h:130
Traits::MatrixL matrixL() const
Definition: LLT.h:137
Computes eigenvalues and eigenvectors of selfadjoint matrices.
Definition: SelfAdjointEigenSolver.h:79
Eigen::Index Index
Definition: SelfAdjointEigenSolver.h:92
@ Ax_lBx
Definition: Constants.h:412
@ ComputeEigenvectors
Definition: Constants.h:407
@ BAx_lx
Definition: Constants.h:418
@ ABx_lx
Definition: Constants.h:415
@ EigenvaluesOnly
Definition: Constants.h:404
Namespace containing all symbols from the Eigen library.
Definition: B01_Experimental.dox:1