eigensolver_generalized_real.cpp
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) 2012-2016 Gael Guennebaud <gael.guennebaud@inria.fr>
5 //
6 // This Source Code Form is subject to the terms of the Mozilla
7 // Public License v. 2.0. If a copy of the MPL was not distributed
8 // with this file, You can obtain one at http://mozilla.org/MPL/2.0/.
9 
10 #define EIGEN_RUNTIME_NO_MALLOC
11 #include "main.h"
12 #include <limits>
13 #include <Eigen/Eigenvalues>
14 #include <Eigen/LU>
15 
16 template<typename MatrixType> void generalized_eigensolver_real(const MatrixType& m)
17 {
18  /* this test covers the following files:
19  GeneralizedEigenSolver.h
20  */
21  Index rows = m.rows();
22  Index cols = m.cols();
23 
24  typedef typename MatrixType::Scalar Scalar;
25  typedef std::complex<Scalar> ComplexScalar;
27 
28  MatrixType a = MatrixType::Random(rows,cols);
29  MatrixType b = MatrixType::Random(rows,cols);
30  MatrixType a1 = MatrixType::Random(rows,cols);
31  MatrixType b1 = MatrixType::Random(rows,cols);
32  MatrixType spdA = a.adjoint() * a + a1.adjoint() * a1;
33  MatrixType spdB = b.adjoint() * b + b1.adjoint() * b1;
34 
35  // lets compare to GeneralizedSelfAdjointEigenSolver
36  {
39 
40  VERIFY_IS_EQUAL(eig.eigenvalues().imag().cwiseAbs().maxCoeff(), 0);
41 
42  VectorType realEigenvalues = eig.eigenvalues().real();
43  std::sort(realEigenvalues.data(), realEigenvalues.data()+realEigenvalues.size());
44  VERIFY_IS_APPROX(realEigenvalues, symmEig.eigenvalues());
45 
46  // check eigenvectors
47  typename GeneralizedEigenSolver<MatrixType>::EigenvectorsType D = eig.eigenvalues().asDiagonal();
49  VERIFY_IS_APPROX(spdA*V, spdB*V*D);
50  }
51 
52  // non symmetric case:
53  {
55  // TODO enable full-prealocation of required memory, this probably requires an in-place mode for HessenbergDecomposition
56  //Eigen::internal::set_is_malloc_allowed(false);
57  eig.compute(a,b);
58  //Eigen::internal::set_is_malloc_allowed(true);
59  for(Index k=0; k<cols; ++k)
60  {
61  Matrix<ComplexScalar,Dynamic,Dynamic> tmp = (eig.betas()(k)*a).template cast<ComplexScalar>() - eig.alphas()(k)*b;
62  if(tmp.size()>1 && tmp.norm()>(std::numeric_limits<Scalar>::min)())
63  tmp /= tmp.norm();
64  VERIFY_IS_MUCH_SMALLER_THAN( std::abs(tmp.determinant()), Scalar(1) );
65  }
66  // check eigenvectors
67  typename GeneralizedEigenSolver<MatrixType>::EigenvectorsType D = eig.eigenvalues().asDiagonal();
69  VERIFY_IS_APPROX(a*V, b*V*D);
70  }
71 
72  // regression test for bug 1098
73  {
74  GeneralizedSelfAdjointEigenSolver<MatrixType> eig1(a.adjoint() * a,b.adjoint() * b);
75  eig1.compute(a.adjoint() * a,b.adjoint() * b);
76  GeneralizedEigenSolver<MatrixType> eig2(a.adjoint() * a,b.adjoint() * b);
77  eig2.compute(a.adjoint() * a,b.adjoint() * b);
78  }
79 
80  // check without eigenvectors
81  {
82  GeneralizedEigenSolver<MatrixType> eig1(spdA, spdB, true);
83  GeneralizedEigenSolver<MatrixType> eig2(spdA, spdB, false);
85  }
86 }
87 
88 EIGEN_DECLARE_TEST(eigensolver_generalized_real)
89 {
90  for(int i = 0; i < g_repeat; i++) {
91  int s = 0;
93  s = internal::random<int>(1,EIGEN_TEST_MAX_SIZE/4);
95 
96  // some trivial but implementation-wise special cases
102  }
103 }
VERIFY_IS_MUCH_SMALLER_THAN
#define VERIFY_IS_MUCH_SMALLER_THAN(a, b)
Definition: main.h:390
simple_graph::b1
Vector2 b1(2, -1)
D
MatrixXcd D
Definition: EigenSolver_EigenSolver_MatrixType.cpp:14
s
RealScalar s
Definition: level1_cplx_impl.h:126
Eigen::GeneralizedEigenSolver
Computes the generalized eigenvalues and eigenvectors of a pair of general matrices.
Definition: GeneralizedEigenSolver.h:58
MatrixType
MatrixXf MatrixType
Definition: benchmark-blocking-sizes.cpp:52
VERIFY_IS_EQUAL
#define VERIFY_IS_EQUAL(a, b)
Definition: main.h:386
b
Scalar * b
Definition: benchVecAdd.cpp:17
Eigen::GeneralizedEigenSolver::compute
GeneralizedEigenSolver & compute(const MatrixType &A, const MatrixType &B, bool computeEigenvectors=true)
Computes generalized eigendecomposition of given matrix.
Definition: GeneralizedEigenSolver.h:287
Eigen::SelfAdjointEigenSolver::eigenvalues
const EIGEN_DEVICE_FUNC RealVectorType & eigenvalues() const
Returns the eigenvalues of given matrix.
Definition: SelfAdjointEigenSolver.h:300
EIGEN_DECLARE_TEST
EIGEN_DECLARE_TEST(eigensolver_generalized_real)
Definition: eigensolver_generalized_real.cpp:88
rows
int rows
Definition: Tutorial_commainit_02.cpp:1
CALL_SUBTEST_4
#define CALL_SUBTEST_4(FUNC)
Definition: split_test_helper.h:22
align_3::a1
Point2 a1
Definition: testPose2.cpp:769
CALL_SUBTEST_3
#define CALL_SUBTEST_3(FUNC)
Definition: split_test_helper.h:16
CALL_SUBTEST_1
#define CALL_SUBTEST_1(FUNC)
Definition: split_test_helper.h:4
generalized_eigensolver_real
void generalized_eigensolver_real(const MatrixType &m)
Definition: eigensolver_generalized_real.cpp:16
Eigen::g_repeat
static int g_repeat
Definition: main.h:169
m
Matrix3f m
Definition: AngleAxis_mimic_euler.cpp:1
eig
SelfAdjointEigenSolver< PlainMatrixType > eig(mat, computeVectors?ComputeEigenvectors:EigenvaluesOnly)
CALL_SUBTEST_2
#define CALL_SUBTEST_2(FUNC)
Definition: split_test_helper.h:10
VERIFY_IS_APPROX
#define VERIFY_IS_APPROX(a, b)
Definition: integer_types.cpp:15
a
ArrayXXi a
Definition: Array_initializer_list_23_cxx11.cpp:1
main.h
Eigen::GeneralizedSelfAdjointEigenSolver
Computes eigenvalues and eigenvectors of the generalized selfadjoint eigen problem.
Definition: GeneralizedSelfAdjointEigenSolver.h:48
EIGEN_TEST_MAX_SIZE
#define EIGEN_TEST_MAX_SIZE
Definition: boostmultiprec.cpp:16
min
#define min(a, b)
Definition: datatypes.h:19
V
MatrixXcd V
Definition: EigenSolver_EigenSolver_MatrixType.cpp:15
Eigen::Matrix
The matrix class, also used for vectors and row-vectors.
Definition: 3rdparty/Eigen/Eigen/src/Core/Matrix.h:178
abs
#define abs(x)
Definition: datatypes.h:17
VectorType
Definition: FFTW.cpp:65
Eigen::GeneralizedEigenSolver::eigenvalues
EigenvalueType eigenvalues() const
Returns an expression of the computed generalized eigenvalues.
Definition: GeneralizedEigenSolver.h:202
Eigen::GeneralizedSelfAdjointEigenSolver::compute
GeneralizedSelfAdjointEigenSolver & compute(const MatrixType &matA, const MatrixType &matB, int options=ComputeEigenvectors|Ax_lBx)
Computes generalized eigendecomposition of given matrix pencil.
Definition: GeneralizedSelfAdjointEigenSolver.h:163
cols
int cols
Definition: Tutorial_commainit_02.cpp:1
i
int i
Definition: BiCGSTAB_step_by_step.cpp:9
TEST_SET_BUT_UNUSED_VARIABLE
#define TEST_SET_BUT_UNUSED_VARIABLE(X)
Definition: main.h:121
Scalar
SCALAR Scalar
Definition: bench_gemm.cpp:46
Eigen::Index
EIGEN_DEFAULT_DENSE_INDEX_TYPE Index
The Index type as used for the API.
Definition: Meta.h:74


gtsam
Author(s):
autogenerated on Fri Nov 1 2024 03:32:29