qr method failed to converge

[jean-michel dot rouet at philips dot com]

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