[NumPy-Tickets] [NumPy] #1454: numpy.linalg.eig computes wrong eigenvalues

[NumPy Trac]

#1454: numpy.linalg.eig computes wrong eigenvalues
-----------------------+----------------------------------------------------
  Reporter:  skavhaug  |       Owner:  somebody
      Type:  defect    |      Status:  closed  
  Priority:  normal    |   Milestone:          
 Component:  Other     |     Version:          
Resolution:  invalid   |    Keywords:          
-----------------------+----------------------------------------------------
Changes (by pv):

  * status:  new => closed
  * resolution:  => invalid


Comment:

 Not a bug in Numpy, I think.

 I get large errors also in Octave, and moreover, with Fortran. This
 probably means your Numpy and Octave installations are not linked against
 the same BLAS or LAPACK libraries, or, they call different LAPACK
 routines.

 {{{
 octave:46> a=randn(28,28);
 octave:47> [v,w]=eig(a);
 octave:48> norm(a*v - v*w)
 ans =  3.7966e-14
 octave:49> a=load('mat.txt');
 octave:50> [v,w]=eig(a);
 octave:51> norm(a*v - v*w)
 ans =  5.1489
 octave:57> for i=1:28; disp(norm(a*v(:,i) - v(:,i)*w(i,i))); end;
  1.4087e-14
  2.5419
  4.9980e-12
  1.8571
  1.8611
  2.9552
  2.9860
  6.0958e-17
  1.8616
  1.2179e-16
  0.13314
  1.8877e-16
  0.14456
  0.14048
  1.1777e-16
  0.14744
  0.16791
  0.16818
  0.16812
  0.16809
  0.17824
  0.17853
  0.17810
  0.17841
  0.18085
  0.18040
  0.18081
  0.18115
 }}}

 {{{
 ~/tmp/eig_bug$ python debug_eig.py
 1.40870336032e-14
 2.54188367754
 4.99798123543e-12
 1.85710182267
 1.86110743511
 2.95523600121
 2.98604486963
 6.09582527197e-17
 1.86161121926
 1.21789448793e-16
 0.133135020496
 1.88766292639e-16
 0.144558641178
 0.140477977597
 1.17769711153e-16
 0.147443420398
 0.167908344177
 0.168175494243
 0.168115056732
 0.168085500384
 0.178240497784
 0.178527005307
 0.178100327262
 0.178405422448
 0.180849858857
 0.180396150447
 0.180807521072
 0.181148599112
 }}}

 Note that Numpy gives exactly the same numbers as Octave here.

 Finally, Fortran:

 {{{
 ~/tmp/eig_bug$ gfortran -o test test.f90  -lblas -llapack
 ~/tmp/eig_bug$ ./test < mat.txt
   1.98444515737284859E-028
    6.4611726301492300
   2.49798164296772676E-023
    3.4488271797617207
    3.4637208850032497
    8.7334198228484805
    8.9164639634284466
   3.71590857463634918E-033
    3.4655963316918923
   1.48326698373198211E-032
   1.77249336825933100E-002
   3.56327132367019683E-032
   2.08972007390916328E-002
   1.97340621896311824E-002
   1.38697048651262032E-032
   2.17395622186454772E-002
   2.81932120443128577E-002
   2.82829968637567722E-002
   2.82626723001265789E-002
   2.82527354393586408E-002
   3.17696750502141287E-002
   3.18718916238386121E-002
   3.17197265706650952E-002
   3.18284947588092018E-002
   3.27066714486422411E-002
   3.25427710961019787E-002
   3.26913596763166309E-002
   3.28148149603390490E-002
 }}}

 So, seems like a LAPACK/BLAS bug, which could be submitted to LAPACK/BLAS
 authors.

-- 
Ticket URL: <http://projects.scipy.org/numpy/ticket/1454#comment:1>
NumPy <http://projects.scipy.org/numpy>
My example project