[Bug-gsl] problem in cubic solver

[Lorenzo Moneta]

Hello,

 I have found a problem with solver for cubic equation (both complex  
and real one)
 In some particular case a NaN is returned, as shown in the attached  
text example.
The problem is observed mainly on 64 bit architectures (for example  
Linux 64-bit  with gcc 4.3)  and not on 32 bit architectures.

This is due to a problem in a sqrt. A patch is attached fixing the  
problem for the complex routine (gsl_poly_complex_solve_cubic ).
A similar patch can be probably done also for the real routine

 This patch fails the current test for polynomial in gsl, however in  
my opinion, this is acceptable, because the test condition is too  
strict. With a slight change of the test coefficient, you will have a  
similar failure also with the current version.

 Best Regards

 Lorenzo Moneta
 ROOT project  (http://root.cern.ch )
 CERN

test_cubic_gsl.c
Description: Binary data

cubic.patch
Description: Binary data