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 <[email protected]> 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 > [email protected] > http://lists.gnu.org/mailman/listinfo/bug-gsl >
