At Tue, 10 Nov 2009 08:02:38 +0100, Katrin Wolff wrote: > I sent the wrong version of the program before, which will not reproduce > the error. I attach the correctly faulty one to this mail.
I have been able to reproduce the problem now, which occurs because of the presence of very small non-zero values that prevent convergence. The patch below should avoid the problem. Thanks for the bug report. -- Brian Gough GNU Scientific Library - http://www.gnu.org/software/gsl/ diff --git a/eigen/qrstep.c b/eigen/qrstep.c index 26705a2..a8b51c4 100644 --- a/eigen/qrstep.c +++ b/eigen/qrstep.c @@ -104,6 +104,16 @@ qrstep (const size_t n, double d[], double sd[], double gc[], double gs[]) double mu = trailing_eigenvalue (n, d, sd); + /* If mu is large relative to d_0 and sd_0 then the Givens rotation + will have no effect, leading to an infinite loop. + + We set mu to zero in this case, which at least diagonalises the + submatrix [d_0, sd_0 ; sd_0, d_0] and allows further progress. */ + + if (GSL_DBL_EPSILON * fabs(mu) > (fabs(d[0]) + fabs(sd[0]))) { + mu = 0; + } + x = d[0] - mu; z = sd[0]; _______________________________________________ Bug-gsl mailing list [email protected] http://lists.gnu.org/mailman/listinfo/bug-gsl
