[Bug-gsl] problem in cubic solver

Andrew W. Steiner awsteiner@gmail.com
Tue Apr 28 16:31:00 GMT 2009


On the topic of cubics, have you tried translating CERNLIB's rrteq3
into C for root? I have found it to give better results in general
than GSL's cubic routine, and thus it might be worth adding a
translation of the CERNLIB code to GSL.

My attempt at this translation is given at line 517 of
http://o2scl.svn.sourceforge.net/viewvc/o2scl/trunk/src/other/poly.cpp?revision=35&view=markup

Take care,
Andrew

On Tue, Apr 28, 2009 at 11:11 AM, Lorenzo Moneta <Lorenzo.Moneta@cern.ch> wrote:
>
> 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
>
>
>
>
>
>
> _______________________________________________
> Bug-gsl mailing list
> Bug-gsl@gnu.org
> http://lists.gnu.org/mailman/listinfo/bug-gsl
>



More information about the Gsl-discuss mailing list