This is the mail archive of the
gsl-discuss@sourceware.org
mailing list for the GSL project.
Re: [Bug-gsl] problem in cubic solver
- From: "Andrew W. Steiner" <awsteiner at gmail dot com>
- To: Lorenzo Moneta <Lorenzo dot Moneta at cern dot ch>, GSL Discuss Mailing List <gsl-discuss at sourceware dot org>
- Date: Tue, 28 Apr 2009 12:31:24 -0400
- Subject: Re: [Bug-gsl] problem in cubic solver
- References: <254EE974-5E07-4EE0-87F2-B2A4FAF11472@cern.ch>
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
>