Re: [Bug-gsl] Bug in polynomial root finder

[Brian Gough]

At Sun, 6 Nov 2011 08:38:22 -0800 (PST),
Munagala Ramanath wrote:
> 
> Thanks for the response.
> 
> I have run it on more than 1.6 million polynomials of order 15 and this is
> the only one on which it has problems.
> 
> I've also run it on some polynomials of order 18 and again it had no problems.
> 
> So it seems the order is not the problem, just some corner case with this
> particular polynomial.
> 

Ah, thanks. That is useful information.

Can you confirm what platform you are one and try the program as a
standalone test case in C -- it works fine for me.

#include <stdio.h>
#include <gsl/gsl_poly.h>

int main ()
{
    double a[16] = { 1, -4, 2, 10, -17, 10, 6, -16, 12, -16, 16, -8, 28, -8, 
-48, 32 };
    double z[16*2];
    int i;

    gsl_poly_complex_workspace *w = gsl_poly_complex_workspace_alloc (16);

    int status = gsl_poly_complex_solve (a, 16, w, z);

    for (i = 0; i<8; i++)
      printf("%.18e %.18e\n", z[2*i], z[2*i+1]);

    gsl_poly_complex_workspace_free (w);
}



$ ./a.out 
-1.000000000000002442e+00 0.000000000000000000e+00
-5.827430765394314705e-01 7.380899912666477602e-01
-5.827430765394314705e-01 -7.380899912666477602e-01
-7.652876256274119271e-01 0.000000000000000000e+00
-6.219637508703932394e-01 0.000000000000000000e+00
-5.472734432737202948e-02 8.769244425642294116e-01
-5.472734432737202948e-02 -8.769244425642294116e-01
3.245172121152936628e-01 6.791075284395361455e-01
$ $ uname -a
Linux nc 2.6.35-30-generic #61-Ubuntu SMP Tue Oct 11 15:29:15 UTC 2011 i686 
GNU/Linux