This is the mail archive of the
gsl-discuss@sources.redhat.com
mailing list for the GSL project.
qr method failed to converge
- To: gsl-discuss at sourceware dot cygnus dot com
- Subject: qr method failed to converge
- From: jean-michel dot rouet at philips dot com
- Date: Mon, 29 Oct 2001 11:54:32 +0100
Hi to all,
although I've been looking through the archive of this mailing list, I
found nothing related to my issue.
I have a 6th degree polynomial expression where I want to find real
roots. I found that the following one was not solved correctly:
-1 + 7450 x - -17929093.75 x^2 + 12320152500 x^3
+ 11774690828125 x^4 - 19748906100000000 x^5
+ 7150798462499999744 x^6
the following prog (see attach also) gives an example of the problem
encountered :
---- polytest.c
#include <gsl/gsl_math.h>
#include <gsl/gsl_poly.h>
#include <gsl/gsl_errno.h>
int main ()
{
double P[7] = { -1, 7450, -17929093.75, 12320152500,
11774690828125, -19748906100000000,
7150798462499999744};
double comp[12];
int status, i;
gsl_poly_complex_workspace * w
= gsl_poly_complex_workspace_alloc (7) ;
status = gsl_poly_complex_solve (P, 7, w, comp);
if (status) {
fprintf (stderr, "failed, gsl_errno=%d (%s)\n",
status, gsl_strerror(status));
}
gsl_poly_complex_workspace_free (w) ;
for (i = 0; i < 6 ; i++)
if (comp[2*i+1] == 0)
printf ("real root: %f\n", comp[2*i]);
exit (1);
}
----
The output is then:
gsl: zsolve.c:79: ERROR: root solving qr method failed to converge
zsh: 2334 IOT instruction (core dumped) /tmp/polytest
I've then tried to change qr.c as following:
--- qr.c Mon Oct 29 11:33:09 2001
+++ /tmp/qr.c Mon Oct 29 11:50:40 2001
@@ -101,14 +101,14 @@
/* No more roots found yet, do another iteration */
- if (iterations == 40)
+ if (iterations == 30)
{
/* too many iterations - give up! */
return GSL_FAILURE ;
}
- if (iterations == 10 || iterations == 20 || iterations == 30)
+ if (iterations == 10 || iterations == 20)
{
/* use an exceptional shift */
and then it happened to find my solutions:
rouet@neptune ~>/tmp/polytest
real root: -0.001000
real root: 0.000290
Can anybody explain why such a hardcoded value "30" is used in the qr method ?
and how it can we get rid of it (using an other hardcoded value such "40" as
I did is not a good solution!)
JM
--
Jean-Michel Rouet (PhD) | Tel : +33 147 28 36 13
Philips Research France (PRF) | Fax : +33 147 28 36 00
51, rue Carnot | email : Jean-Michel.Rouet@philips.com
BP 301 - 92156 Suresnes Cedex
polytest.c