This bugzilla service is closed. All entries have been migrated to https://gitlab.com/libeigen/eigen

Bug 1557

Summary: Fail to compute eigenvalues for a simple 3x3 companion matrix for root finding
Product: Eigen Reporter: Gary Tan <tgntanguangning>
Component: EigenvaluesAssignee: Nobody <eigen.nobody>
Status: RESOLVED FIXED    
Severity: Wrong Result CC: chtz, gael.guennebaud, jacob.benoit.1, jitseniesen, pcv22, tgntanguangning
Priority: Normal Keywords: test-needed
Version: 3.3 (current stable)   
Hardware: All   
OS: All   
Whiteboard:
Bug Depends on:    
Bug Blocks: 814    

Description Gary Tan 2018-06-15 07:16:59 UTC
The following

#include <Eigen/Eigenvalues>
#include <iostream>
using namespace Eigen;
using namespace std;

int main() {
  Matrix3d A;
  A << 0, 0, 0, 1, 0, 0.5887907064808635127, 0, 1, 0;
  cout << "Here is a random 3x3 matrix, A:" << endl << A << endl << endl;
  EigenSolver<MatrixXd> es(A);
  cout << "The eigenvalues of A are:" << endl << es.eigenvalues() << endl;
  cout << "The matrix of eigenvectors, V, is:" << endl
       << es.eigenvectors() << endl
       << endl;
  complex<double> lambda = es.eigenvalues()[0];
  cout << "Consider the first eigenvalue, lambda = " << lambda << endl;
  VectorXcd v = es.eigenvectors().col(0);
  cout << "If v is the corresponding eigenvector, then lambda * v = " << endl
       << lambda * v << endl;
  cout << "... and A * v = " << endl
       << A.cast<complex<double> >() * v << endl
       << endl;
  MatrixXcd D = es.eigenvalues().asDiagonal();
  MatrixXcd V = es.eigenvectors();
  cout << "Finally, V * D * V^(-1) = " << endl << V * D * V.inverse() << endl;

  return 0;
}

produces

Here is a random 3x3 matrix, A:
       0        0        0
       1        0 0.588791
       0        1        0

The eigenvalues of A are:
(0,0)
(0,0)
(0,0)
The matrix of eigenvectors, V, is:
(1,0) (0,0) (0,0)
(0,0) (0,0) (0,0)
(0,0) (0,0) (0,0)

Consider the first eigenvalue, lambda = (0,0)
If v is the corresponding eigenvector, then lambda * v = 
(0,0)
(0,0)
(0,0)
... and A * v = 
(0,0)
(1,0)
(0,0)

Finally, V * D * V^(-1) = 
(-nan,-nan) (-nan,-nan) (-nan,-nan)
(-nan,-nan) (-nan,-nan) (-nan,-nan)
(-nan,-nan) (-nan,-nan) (-nan,-nan)

However, it's expected to obtain
sqrt(0.5887907064808635127) = 0.76733
and its opposite.

Thanks,
Gary
Comment 1 Gary Tan 2018-06-15 07:58:21 UTC
Matrix3d A;
A << 0, 0, 0, 1, 0, 0, 0, 1, 0;
double a = double(0.5887907064808635127);
A(1,2) = a;

This does not work either. If the number "0.5887907..." were "shorter", then it'd work fine. However, in general that number is assigned from other variable, say a. I tried to debug, but could not figure out why. Sorry for creating this weird problem.
Comment 2 Christoph Hertzberg 2018-06-15 08:06:12 UTC
If you check `es.info()` you will notice that it did not converge -- thus the computed values will not be useful.

It is strange that it does not converge, however. It may be failing, because all Eigenvalues are indeed real, e.g., you get a valid result if you set `A(0,2)=1e-30;`
Comment 3 Gary Tan 2018-06-15 08:19:40 UTC
(In reply to Christoph Hertzberg from comment #2)
> If you check `es.info()` you will notice that it did not converge -- thus
> the computed values will not be useful.
> 
> It is strange that it does not converge, however. It may be failing, because
> all Eigenvalues are indeed real, e.g., you get a valid result if you set
> `A(0,2)=1e-30;`

Hi, Chris!
If I use
0.588790706480863
or
0.588790706480864
then it computes fine and produces expected values. This strange problem starts with
0.5887907064808635
which is just a few epsilons away. It seems internally something wrongly handles this number. I ran the same code on macOS and got the same results. 

I am using g++ (Debian 6.3.0) and Apple LLVM v10.0.0 (clang-1000.10.25.5 on macOS 10.14).
Comment 4 Christoph Hertzberg 2018-12-04 11:32:55 UTC
Bisecting leads to this commit breaking your example:
https://bitbucket.org/eigen/eigen/commits/c90f5b999226

There is hopefully a way to fix this without breaking Bug 933.
Comment 5 Gael Guennebaud 2018-12-12 16:35:51 UTC
Thank you both for the findings.

Fixed with regression tests for both 933 and this one:

https://bitbucket.org/eigen/eigen/commits/4d387242247e/
https://bitbucket.org/eigen/eigen/commits/7e7fef7bac37/ (3.3)
Comment 6 Gael Guennebaud 2018-12-12 17:05:07 UTC
*** Bug 1174 has been marked as a duplicate of this bug. ***
Comment 7 Nobody 2019-12-04 17:42:58 UTC
-- GitLab Migration Automatic Message --

This bug has been migrated to gitlab.com's GitLab instance and has been closed from further activity.

You can subscribe and participate further through the new bug through this link to our GitLab instance: https://gitlab.com/libeigen/eigen/issues/1557.