Re: [eigen] Issue in SelfAdjointEigenSolver for 3-by-3 matrices
• To: eigen@xxxxxxxxxxxxxxxxxxx
• Subject: Re: [eigen] Issue in SelfAdjointEigenSolver for 3-by-3 matrices
• From: Manoj Rajagopalan <rmanoj@xxxxxxxxx>
• Date: Thu, 15 Jul 2010 17:12:34 -0400
• Organization: EECS Dept., University of Michigan, Ann Arbor, MI, USA

```Just an observation: you probably meant A(2,0) = 0.001;
The documentation says only the lower triangular part is referenced.

Anyway, I can confirm what you have noticed with this correction.

-- Manoj

On Thursday 15 July 2010 03:59:09 pm Jitse Niesen wrote:
> Hi all,
>
> There seems to be an issue with SelfAdjointEigenSolver when applied to
> 3-by-3 matrices, because the unit test eigensolver_selfadjoint_1 fails
> occasionally on my computer.
>
> I think I boiled it down to this example:
>
> ------------------------------------------------
>
> #include <iostream>
> #include <Eigen/Core>
> #include <Eigen/Eigenvalues>
>
> using namespace Eigen;
> using namespace std;
>
> int main()
> {
>    Matrix3f A = Matrix3f::Identity(3,3);
>    A(0,2) = 0.001;
>    cout << eigOfA.eigenvalues() << endl;
> }
>
> ------------------------------------------------
>
> The output of this program is
>
>    1
>    1
>    1
>
> This program computes the eigenvalues of the matrix
>
>    1      0     0.001
>    0      1     0
>    0.001  0     1
>
> They are easily computed by hand; the exact eigenvalues are 0.999, 1, and
> 1.001. I think this is rather too big an error.
>
> The issue seems to be in
>
>    ei_tridiagonalization_inplace_selector<MatrixType,3>::run()
>
> in Tridiagonalization.h, where we have
>
>      RealScalar v1norm2 = ei_abs2(mat(2,0));
>      if (ei_isMuchSmallerThan(v1norm2, RealScalar(1)))
>      {
>  	...
>
> If the condition in the if statement is true, then mat(2,0) is simply
> ignored. This does not seem the right condition, and I'm not even sure why