From: Lorenzo Moneta <Lorenzo.Moneta@cern.ch>
To: Andrew W.Steiner <awsteiner@gmail.com>
Cc: GSL Discuss Mailing List <gsl-discuss@sourceware.org>
Subject: Re: [Bug-gsl] problem in cubic solver
Date: Tue, 28 Apr 2009 20:05:00 -0000 [thread overview]
Message-ID: <99A23E5B-6A5F-476F-B98E-516EA7DB80C4@cern.ch> (raw)
In-Reply-To: <63c059b10904280931o6600dacg80f040ff32a1ac71@mail.gmail.com>
Hi Andrew,
Thank you for the link, we have not in Root performed this
translation. We have currently an implementation very similar to the
GSL one
(but only finding real roots).
Thank you for the link, I will try your code and check as well if
it gives better numerical stability. In that case it wold be good also
to have it in GSL.
I remember also, you were working on importing the GSL quartic in GSL.
Is this still planned ?
We are using presently in ROOT this code, it would be good if we
could use instead GSL directly.
Best Regards
Lorenzo
On Apr 28, 2009, at 6:31 PM, Andrew W. Steiner wrote:
> 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
>>
next prev parent reply other threads:[~2009-04-28 20:05 UTC|newest]
Thread overview: 3+ messages / expand[flat|nested] mbox.gz Atom feed top
[not found] <254EE974-5E07-4EE0-87F2-B2A4FAF11472@cern.ch>
2009-04-28 16:31 ` Andrew W. Steiner
2009-04-28 20:05 ` Lorenzo Moneta [this message]
2009-04-28 20:36 ` Andrew W. Steiner
Reply instructions:
You may reply publicly to this message via plain-text email
using any one of the following methods:
* Save the following mbox file, import it into your mail client,
and reply-to-all from there: mbox
Avoid top-posting and favor interleaved quoting:
https://en.wikipedia.org/wiki/Posting_style#Interleaved_style
* Reply using the --to, --cc, and --in-reply-to
switches of git-send-email(1):
git send-email \
--in-reply-to=99A23E5B-6A5F-476F-B98E-516EA7DB80C4@cern.ch \
--to=lorenzo.moneta@cern.ch \
--cc=awsteiner@gmail.com \
--cc=gsl-discuss@sourceware.org \
/path/to/YOUR_REPLY
https://kernel.org/pub/software/scm/git/docs/git-send-email.html
* If your mail client supports setting the In-Reply-To header
via mailto: links, try the mailto: link
Be sure your reply has a Subject: header at the top and a blank line
before the message body.
This is a public inbox, see mirroring instructions
for how to clone and mirror all data and code used for this inbox;
as well as URLs for read-only IMAP folder(s) and NNTP newsgroup(s).