Re: [Bug-gsl] Re: infinite loop in gsl_eigen_symm

[Brian Gough]

At Tue, 28 Aug 2007 08:55:56 +0200,
Andries E. Brouwer wrote:
> Note that my above "all is fine" only refers to the looping behaviour
> that is avoided by this test. The results are bad, in the sense that
> the computed eigenvalues already have errors in the second decimal
> (e.g. 2.292 instead of 2.303).
> 
> I have not seen much documentation for these functions - don't know whether
> there are any claims about a guaranteed accuracy. Is this basically a
> function that you give a matrix and it returns some numbers and nobody
> knows how far the real eigenvalues are from the returned numbers?

Thanks for the suggestion, I will add a short discussion of this to
the manual.  The algorithm is symmetric bidiagonalization followed by
QR reduction (see Golub & Van Loan, "Matrix Computations" (3rd ed),
Section 8.3 for details).  The eigenvalues are accurate to an absolute
error of eps*||A||_2, where eps = double precision accuracy.

> (It is easy to compute the exact eigenvalues in the example I gave.
> They are 0 with multiplicity 21, +-3sqrt{3}, +-(1+-sqrt{13})/2.)

For the matrix you sent, I get a different characteristic polynomial
(via pari/gp's matdet function)  p(u)=u^21 (u^6 -25*u^4 + 134*u^2 - 159) 

Solving this gives the same values as gsl_eigen_symm to about 15
digits precision, as expected from the error bound.

-- 
Brian Gough

Network Theory Ltd,
Publishing the GSL Manual - http://www.network-theory.co.uk/gsl/manual/