[Boost-bugs] [Boost C++ Libraries] #4812: Bug in boost::math::special_functions::sign showing up in special_functions::bessel_jy()

Subject: [Boost-bugs] [Boost C++ Libraries] #4812: Bug in boost::math::special_functions::sign showing up in special_functions::bessel_jy()
From: Boost C++ Libraries (noreply_at_[hidden])
Date: 2010-11-03 19:53:37

#4812: Bug in boost::math::special_functions::sign showing up in
special_functions::bessel_jy()
----------------------------------------------+-----------------------------
 Reporter:  peter zwart <PHZwart@…>           |       Owner:  johnmaddock
     Type:  Bugs                              |      Status:  new        
Milestone:  To Be Determined                  |   Component:  math       
  Version:  Boost 1.45.0                      |    Severity:  Showstopper
 Keywords:  sign sph_bessel special_function  |  
----------------------------------------------+-----------------------------
 Hi,

 When computing scattering curves, I encountered a show-stopping bug that I
 tracked back to the boost::math::special_functions::sign function .

 In math::special_functions::bessel_jy a small number is generated via

 sqrt(tools::min_value<T>())

 which subsequently ends up in this statement:  sign(current)

 I see this happening:

 variable current: -1.24911e-2468
 variable sign(current): 1

 which causes problems for math::special_functions::sph_bessel(1,x) for x
 between 7.845 and 8 or so on my machines.

 The fix that worked for me is the following

           T sign_current;
           sign_current = std::fabs(current);
           if (current*2.0 < sign_current){
             sign_current = -1.0;
           } else {
             sign_current = 1.0;
           }

           Ju = sign_current * sqrt(W / (q + gamma * (p - t)));
           // Ju = sign(current) * sqrt(W / (q + gamma * (p - t)));


 For now, I'll run with this local fix but it seems that on some platforms
 the current sign function causes problems when computing spherical Bessel
 functions of any order above 0.

-- 
Ticket URL: <https://svn.boost.org/trac/boost/ticket/4812>
Boost C++ Libraries <http://www.boost.org/>
Boost provides free peer-reviewed portable C++ source libraries.