Opened 8 years ago
Closed 8 years ago
#10905 closed Bugs (fixed)
sinpx function in boost/math/special_functions/gamma.hpp returns incorrect results for positive z arguments
Reported by: | Owned by: | John Maddock | |
---|---|---|---|
Milestone: | To Be Determined | Component: | math |
Version: | Boost 1.57.0 | Severity: | Showstopper |
Keywords: | Cc: |
Description
The sinpx function in gamma.hpp produces incorrect results for all positive arguments. The sign of the result is exactly the opposite of the correct value. As the sinpx function is called by all the gamma functions, their results are affected as well.
Here is the 1.57 version of sinpx. Removing the "else { sign = -sign; }" block fixes the problem.
Attached is a source file that calculates z * sin (z * PI) using the sinpx function and the C++ math.h sin() function. Note that the sinpx function produces the right absolute value, but with the wrong sign, for all positive z arguments.
template <class T> T sinpx(T z) { // Ad hoc function calculates x * sin(pi * x), // taking extra care near when x is near a whole number. BOOST_MATH_STD_USING int sign = 1; if(z < 0) { z = -z; } else { sign = -sign; } T fl = floor(z); T dist; if(is_odd(fl)) { fl += 1; dist = fl - z; sign = -sign; } else { dist = z - fl; } BOOST_ASSERT(fl >= 0); if(dist > 0.5) dist = 1 - dist; T result = sin(dist*boost::math::constants::pi<T>()); return sign*z*result; }
The following program demonstrates the bug:
#include <boost/math/special_functions/gamma.hpp> #include <boost/format.hpp> #include <iostream> #include <math.h> int main() { using namespace boost::math::detail; for (double z = -10.0; z < 10.1; z += 0.1) { double result = sinpx(z); double result2 = z * sin(z * 3.14159265); bool error = abs(result - result2) > 0.001; auto errormsg = error ? "\t*** Error" : ""; std::cout << boost::format ("%1$4.1f * sinpx (%1$4.1f * PI) is %2$11.8f\tC++ sin() is %3$11.8f%4%\n") % z % result % result2 % errormsg; } }
It produces the following output:
-10.0 * sinpx (-10.0 * PI) is 0.00000000 C++ sin() is -0.00000036 -9.9 * sinpx (-9.9 * PI) is -3.05926824 C++ sin() is -3.05926858 -9.8 * sinpx (-9.8 * PI) is -5.76029547 C++ sin() is -5.76029575 -9.7 * sinpx (-9.7 * PI) is -7.84746485 C++ sin() is -7.84746504 -9.6 * sinpx (-9.6 * PI) is -9.13014256 C++ sin() is -9.13014266 -9.5 * sinpx (-9.5 * PI) is -9.50000000 C++ sin() is -9.50000000 -9.4 * sinpx (-9.4 * PI) is -8.93993125 C++ sin() is -8.93993116 -9.3 * sinpx (-9.3 * PI) is -7.52385805 C++ sin() is -7.52385787 -9.2 * sinpx (-9.2 * PI) is -5.40762432 C++ sin() is -5.40762408 -9.1 * sinpx (-9.1 * PI) is -2.81205465 C++ sin() is -2.81205437 -9.0 * sinpx (-9.0 * PI) is -0.00000000 C++ sin() is 0.00000029 -8.9 * sinpx (-8.9 * PI) is 2.75025125 C++ sin() is 2.75025152 -8.8 * sinpx (-8.8 * PI) is 5.17251022 C++ sin() is 5.17251045 -8.7 * sinpx (-8.7 * PI) is 7.03844785 C++ sin() is 7.03844801 -8.6 * sinpx (-8.6 * PI) is 8.17908604 C++ sin() is 8.17908612 -8.5 * sinpx (-8.5 * PI) is 8.50000000 C++ sin() is 8.50000000 -8.4 * sinpx (-8.4 * PI) is 7.98887474 C++ sin() is 7.98887466 -8.3 * sinpx (-8.3 * PI) is 6.71484105 C++ sin() is 6.71484091 -8.2 * sinpx (-8.2 * PI) is 4.81983907 C++ sin() is 4.81983887 -8.1 * sinpx (-8.1 * PI) is 2.50303765 C++ sin() is 2.50303743 -8.0 * sinpx (-8.0 * PI) is 0.00000000 C++ sin() is -0.00000023 -7.9 * sinpx (-7.9 * PI) is -2.44123426 C++ sin() is -2.44123447 -7.8 * sinpx (-7.8 * PI) is -4.58472497 C++ sin() is -4.58472514 -7.7 * sinpx (-7.7 * PI) is -6.22943086 C++ sin() is -6.22943098 -7.6 * sinpx (-7.6 * PI) is -7.22802952 C++ sin() is -7.22802959 -7.5 * sinpx (-7.5 * PI) is -7.50000000 C++ sin() is -7.50000000 -7.4 * sinpx (-7.4 * PI) is -7.03781822 C++ sin() is -7.03781816 -7.3 * sinpx (-7.3 * PI) is -5.90582406 C++ sin() is -5.90582395 -7.2 * sinpx (-7.2 * PI) is -4.23205382 C++ sin() is -4.23205367 -7.1 * sinpx (-7.1 * PI) is -2.19402066 C++ sin() is -2.19402049 -7.0 * sinpx (-7.0 * PI) is -0.00000000 C++ sin() is 0.00000018 -6.9 * sinpx (-6.9 * PI) is 2.13221726 C++ sin() is 2.13221742 -6.8 * sinpx (-6.8 * PI) is 3.99693972 C++ sin() is 3.99693985 -6.7 * sinpx (-6.7 * PI) is 5.42041386 C++ sin() is 5.42041396 -6.6 * sinpx (-6.6 * PI) is 6.27697301 C++ sin() is 6.27697306 -6.5 * sinpx (-6.5 * PI) is 6.50000000 C++ sin() is 6.50000000 -6.4 * sinpx (-6.4 * PI) is 6.08676170 C++ sin() is 6.08676166 -6.3 * sinpx (-6.3 * PI) is 5.09680706 C++ sin() is 5.09680698 -6.2 * sinpx (-6.2 * PI) is 3.64426856 C++ sin() is 3.64426845 -6.1 * sinpx (-6.1 * PI) is 1.88500367 C++ sin() is 1.88500354 -6.0 * sinpx (-6.0 * PI) is 0.00000000 C++ sin() is -0.00000013 -5.9 * sinpx (-5.9 * PI) is -1.82320027 C++ sin() is -1.82320039 -5.8 * sinpx (-5.8 * PI) is -3.40915446 C++ sin() is -3.40915456 -5.7 * sinpx (-5.7 * PI) is -4.61139687 C++ sin() is -4.61139694 -5.6 * sinpx (-5.6 * PI) is -5.32591649 C++ sin() is -5.32591653 -5.5 * sinpx (-5.5 * PI) is -5.50000000 C++ sin() is -5.50000000 -5.4 * sinpx (-5.4 * PI) is -5.13570519 C++ sin() is -5.13570516 -5.3 * sinpx (-5.3 * PI) is -4.28779007 C++ sin() is -4.28779001 -5.2 * sinpx (-5.2 * PI) is -3.05648331 C++ sin() is -3.05648323 -5.1 * sinpx (-5.1 * PI) is -1.57598667 C++ sin() is -1.57598658 -5.0 * sinpx (-5.0 * PI) is -0.00000000 C++ sin() is 0.00000009 -4.9 * sinpx (-4.9 * PI) is 1.51418327 C++ sin() is 1.51418335 -4.8 * sinpx (-4.8 * PI) is 2.82136921 C++ sin() is 2.82136928 -4.7 * sinpx (-4.7 * PI) is 3.80237987 C++ sin() is 3.80237992 -4.6 * sinpx (-4.6 * PI) is 4.37485997 C++ sin() is 4.37486000 -4.5 * sinpx (-4.5 * PI) is 4.50000000 C++ sin() is 4.50000000 -4.4 * sinpx (-4.4 * PI) is 4.18464867 C++ sin() is 4.18464865 -4.3 * sinpx (-4.3 * PI) is 3.47877308 C++ sin() is 3.47877304 -4.2 * sinpx (-4.2 * PI) is 2.46869806 C++ sin() is 2.46869801 -4.1 * sinpx (-4.1 * PI) is 1.26696968 C++ sin() is 1.26696962 -4.0 * sinpx (-4.0 * PI) is 0.00000000 C++ sin() is -0.00000006 -3.9 * sinpx (-3.9 * PI) is -1.20516628 C++ sin() is -1.20516633 -3.8 * sinpx (-3.8 * PI) is -2.23358396 C++ sin() is -2.23358400 -3.7 * sinpx (-3.7 * PI) is -2.99336288 C++ sin() is -2.99336291 -3.6 * sinpx (-3.6 * PI) is -3.42380346 C++ sin() is -3.42380347 -3.5 * sinpx (-3.5 * PI) is -3.50000000 C++ sin() is -3.50000000 -3.4 * sinpx (-3.4 * PI) is -3.23359216 C++ sin() is -3.23359214 -3.3 * sinpx (-3.3 * PI) is -2.66975608 C++ sin() is -2.66975606 -3.2 * sinpx (-3.2 * PI) is -1.88091281 C++ sin() is -1.88091278 -3.1 * sinpx (-3.1 * PI) is -0.95795268 C++ sin() is -0.95795265 -3.0 * sinpx (-3.0 * PI) is -0.00000000 C++ sin() is 0.00000003 -2.9 * sinpx (-2.9 * PI) is 0.89614928 C++ sin() is 0.89614931 -2.8 * sinpx (-2.8 * PI) is 1.64579871 C++ sin() is 1.64579873 -2.7 * sinpx (-2.7 * PI) is 2.18434588 C++ sin() is 2.18434590 -2.6 * sinpx (-2.6 * PI) is 2.47274694 C++ sin() is 2.47274695 -2.5 * sinpx (-2.5 * PI) is 2.50000000 C++ sin() is 2.50000000 -2.4 * sinpx (-2.4 * PI) is 2.28253564 C++ sin() is 2.28253563 -2.3 * sinpx (-2.3 * PI) is 1.86073909 C++ sin() is 1.86073908 -2.2 * sinpx (-2.2 * PI) is 1.29312756 C++ sin() is 1.29312754 -2.1 * sinpx (-2.1 * PI) is 0.64893569 C++ sin() is 0.64893567 -2.0 * sinpx (-2.0 * PI) is 0.00000000 C++ sin() is -0.00000001 -1.9 * sinpx (-1.9 * PI) is -0.58713229 C++ sin() is -0.58713230 -1.8 * sinpx (-1.8 * PI) is -1.05801345 C++ sin() is -1.05801346 -1.7 * sinpx (-1.7 * PI) is -1.37532889 C++ sin() is -1.37532890 -1.6 * sinpx (-1.6 * PI) is -1.52169043 C++ sin() is -1.52169043 -1.5 * sinpx (-1.5 * PI) is -1.50000000 C++ sin() is -1.50000000 -1.4 * sinpx (-1.4 * PI) is -1.33147912 C++ sin() is -1.33147912 -1.3 * sinpx (-1.3 * PI) is -1.05172209 C++ sin() is -1.05172209 -1.2 * sinpx (-1.2 * PI) is -0.70534230 C++ sin() is -0.70534230 -1.1 * sinpx (-1.1 * PI) is -0.33991869 C++ sin() is -0.33991869 -1.0 * sinpx (-1.0 * PI) is -0.00000000 C++ sin() is 0.00000000 -0.9 * sinpx (-0.9 * PI) is 0.27811529 C++ sin() is 0.27811530 -0.8 * sinpx (-0.8 * PI) is 0.47022820 C++ sin() is 0.47022820 -0.7 * sinpx (-0.7 * PI) is 0.56631190 C++ sin() is 0.56631190 -0.6 * sinpx (-0.6 * PI) is 0.57063391 C++ sin() is 0.57063391 -0.5 * sinpx (-0.5 * PI) is 0.50000000 C++ sin() is 0.50000000 -0.4 * sinpx (-0.4 * PI) is 0.38042261 C++ sin() is 0.38042261 -0.3 * sinpx (-0.3 * PI) is 0.24270510 C++ sin() is 0.24270510 -0.2 * sinpx (-0.2 * PI) is 0.11755705 C++ sin() is 0.11755705 -0.1 * sinpx (-0.1 * PI) is 0.03090170 C++ sin() is 0.03090170 -0.0 * sinpx (-0.0 * PI) is 0.00000000 C++ sin() is 0.00000000 0.1 * sinpx ( 0.1 * PI) is -0.03090170 C++ sin() is 0.03090170 *** Error 0.2 * sinpx ( 0.2 * PI) is -0.11755705 C++ sin() is 0.11755705 *** Error 0.3 * sinpx ( 0.3 * PI) is -0.24270510 C++ sin() is 0.24270510 *** Error 0.4 * sinpx ( 0.4 * PI) is -0.38042261 C++ sin() is 0.38042261 *** Error 0.5 * sinpx ( 0.5 * PI) is -0.50000000 C++ sin() is 0.50000000 *** Error 0.6 * sinpx ( 0.6 * PI) is -0.57063391 C++ sin() is 0.57063391 *** Error 0.7 * sinpx ( 0.7 * PI) is -0.56631190 C++ sin() is 0.56631190 *** Error 0.8 * sinpx ( 0.8 * PI) is -0.47022820 C++ sin() is 0.47022820 *** Error 0.9 * sinpx ( 0.9 * PI) is -0.27811529 C++ sin() is 0.27811530 *** Error 1.0 * sinpx ( 1.0 * PI) is -0.00000000 C++ sin() is 0.00000000 1.1 * sinpx ( 1.1 * PI) is 0.33991869 C++ sin() is -0.33991869 *** Error 1.2 * sinpx ( 1.2 * PI) is 0.70534230 C++ sin() is -0.70534230 *** Error 1.3 * sinpx ( 1.3 * PI) is 1.05172209 C++ sin() is -1.05172209 *** Error 1.4 * sinpx ( 1.4 * PI) is 1.33147912 C++ sin() is -1.33147912 *** Error 1.5 * sinpx ( 1.5 * PI) is 1.50000000 C++ sin() is -1.50000000 *** Error 1.6 * sinpx ( 1.6 * PI) is 1.52169043 C++ sin() is -1.52169043 *** Error 1.7 * sinpx ( 1.7 * PI) is 1.37532889 C++ sin() is -1.37532890 *** Error 1.8 * sinpx ( 1.8 * PI) is 1.05801345 C++ sin() is -1.05801346 *** Error 1.9 * sinpx ( 1.9 * PI) is 0.58713229 C++ sin() is -0.58713230 *** Error 2.0 * sinpx ( 2.0 * PI) is 0.00000000 C++ sin() is -0.00000001 2.1 * sinpx ( 2.1 * PI) is -0.64893569 C++ sin() is 0.64893567 *** Error 2.2 * sinpx ( 2.2 * PI) is -1.29312756 C++ sin() is 1.29312754 *** Error 2.3 * sinpx ( 2.3 * PI) is -1.86073909 C++ sin() is 1.86073908 *** Error 2.4 * sinpx ( 2.4 * PI) is -2.28253564 C++ sin() is 2.28253563 *** Error 2.5 * sinpx ( 2.5 * PI) is -2.50000000 C++ sin() is 2.50000000 *** Error 2.6 * sinpx ( 2.6 * PI) is -2.47274694 C++ sin() is 2.47274695 *** Error 2.7 * sinpx ( 2.7 * PI) is -2.18434588 C++ sin() is 2.18434590 *** Error 2.8 * sinpx ( 2.8 * PI) is -1.64579871 C++ sin() is 1.64579873 *** Error 2.9 * sinpx ( 2.9 * PI) is -0.89614928 C++ sin() is 0.89614931 *** Error 3.0 * sinpx ( 3.0 * PI) is -0.00000000 C++ sin() is 0.00000003 3.1 * sinpx ( 3.1 * PI) is 0.95795268 C++ sin() is -0.95795265 *** Error 3.2 * sinpx ( 3.2 * PI) is 1.88091281 C++ sin() is -1.88091278 *** Error 3.3 * sinpx ( 3.3 * PI) is 2.66975608 C++ sin() is -2.66975606 *** Error 3.4 * sinpx ( 3.4 * PI) is 3.23359216 C++ sin() is -3.23359214 *** Error 3.5 * sinpx ( 3.5 * PI) is 3.50000000 C++ sin() is -3.50000000 *** Error 3.6 * sinpx ( 3.6 * PI) is 3.42380346 C++ sin() is -3.42380347 *** Error 3.7 * sinpx ( 3.7 * PI) is 2.99336288 C++ sin() is -2.99336291 *** Error 3.8 * sinpx ( 3.8 * PI) is 2.23358396 C++ sin() is -2.23358400 *** Error 3.9 * sinpx ( 3.9 * PI) is 1.20516628 C++ sin() is -1.20516633 *** Error 4.0 * sinpx ( 4.0 * PI) is 0.00000000 C++ sin() is -0.00000006 4.1 * sinpx ( 4.1 * PI) is -1.26696968 C++ sin() is 1.26696962 *** Error 4.2 * sinpx ( 4.2 * PI) is -2.46869806 C++ sin() is 2.46869801 *** Error 4.3 * sinpx ( 4.3 * PI) is -3.47877308 C++ sin() is 3.47877304 *** Error 4.4 * sinpx ( 4.4 * PI) is -4.18464867 C++ sin() is 4.18464865 *** Error 4.5 * sinpx ( 4.5 * PI) is -4.50000000 C++ sin() is 4.50000000 *** Error 4.6 * sinpx ( 4.6 * PI) is -4.37485997 C++ sin() is 4.37486000 *** Error 4.7 * sinpx ( 4.7 * PI) is -3.80237987 C++ sin() is 3.80237992 *** Error 4.8 * sinpx ( 4.8 * PI) is -2.82136921 C++ sin() is 2.82136928 *** Error 4.9 * sinpx ( 4.9 * PI) is -1.51418327 C++ sin() is 1.51418335 *** Error 5.0 * sinpx ( 5.0 * PI) is -0.00000000 C++ sin() is 0.00000009 5.1 * sinpx ( 5.1 * PI) is 1.57598667 C++ sin() is -1.57598658 *** Error 5.2 * sinpx ( 5.2 * PI) is 3.05648331 C++ sin() is -3.05648323 *** Error 5.3 * sinpx ( 5.3 * PI) is 4.28779007 C++ sin() is -4.28779001 *** Error 5.4 * sinpx ( 5.4 * PI) is 5.13570519 C++ sin() is -5.13570516 *** Error 5.5 * sinpx ( 5.5 * PI) is 5.50000000 C++ sin() is -5.50000000 *** Error 5.6 * sinpx ( 5.6 * PI) is 5.32591649 C++ sin() is -5.32591653 *** Error 5.7 * sinpx ( 5.7 * PI) is 4.61139687 C++ sin() is -4.61139694 *** Error 5.8 * sinpx ( 5.8 * PI) is 3.40915446 C++ sin() is -3.40915456 *** Error 5.9 * sinpx ( 5.9 * PI) is 1.82320027 C++ sin() is -1.82320039 *** Error 6.0 * sinpx ( 6.0 * PI) is 0.00000000 C++ sin() is -0.00000013 6.1 * sinpx ( 6.1 * PI) is -1.88500367 C++ sin() is 1.88500354 *** Error 6.2 * sinpx ( 6.2 * PI) is -3.64426856 C++ sin() is 3.64426845 *** Error 6.3 * sinpx ( 6.3 * PI) is -5.09680706 C++ sin() is 5.09680698 *** Error 6.4 * sinpx ( 6.4 * PI) is -6.08676170 C++ sin() is 6.08676166 *** Error 6.5 * sinpx ( 6.5 * PI) is -6.50000000 C++ sin() is 6.50000000 *** Error 6.6 * sinpx ( 6.6 * PI) is -6.27697301 C++ sin() is 6.27697306 *** Error 6.7 * sinpx ( 6.7 * PI) is -5.42041386 C++ sin() is 5.42041396 *** Error 6.8 * sinpx ( 6.8 * PI) is -3.99693972 C++ sin() is 3.99693985 *** Error 6.9 * sinpx ( 6.9 * PI) is -2.13221726 C++ sin() is 2.13221742 *** Error 7.0 * sinpx ( 7.0 * PI) is -0.00000000 C++ sin() is 0.00000018 7.1 * sinpx ( 7.1 * PI) is 2.19402066 C++ sin() is -2.19402049 *** Error 7.2 * sinpx ( 7.2 * PI) is 4.23205382 C++ sin() is -4.23205367 *** Error 7.3 * sinpx ( 7.3 * PI) is 5.90582406 C++ sin() is -5.90582395 *** Error 7.4 * sinpx ( 7.4 * PI) is 7.03781822 C++ sin() is -7.03781816 *** Error 7.5 * sinpx ( 7.5 * PI) is 7.50000000 C++ sin() is -7.50000000 *** Error 7.6 * sinpx ( 7.6 * PI) is 7.22802952 C++ sin() is -7.22802959 *** Error 7.7 * sinpx ( 7.7 * PI) is 6.22943086 C++ sin() is -6.22943098 *** Error 7.8 * sinpx ( 7.8 * PI) is 4.58472497 C++ sin() is -4.58472514 *** Error 7.9 * sinpx ( 7.9 * PI) is 2.44123426 C++ sin() is -2.44123447 *** Error 8.0 * sinpx ( 8.0 * PI) is 0.00000000 C++ sin() is -0.00000023 8.1 * sinpx ( 8.1 * PI) is -2.50303765 C++ sin() is 2.50303743 *** Error 8.2 * sinpx ( 8.2 * PI) is -4.81983907 C++ sin() is 4.81983887 *** Error 8.3 * sinpx ( 8.3 * PI) is -6.71484105 C++ sin() is 6.71484091 *** Error 8.4 * sinpx ( 8.4 * PI) is -7.98887474 C++ sin() is 7.98887466 *** Error 8.5 * sinpx ( 8.5 * PI) is -8.50000000 C++ sin() is 8.50000000 *** Error 8.6 * sinpx ( 8.6 * PI) is -8.17908604 C++ sin() is 8.17908612 *** Error 8.7 * sinpx ( 8.7 * PI) is -7.03844785 C++ sin() is 7.03844801 *** Error 8.8 * sinpx ( 8.8 * PI) is -5.17251022 C++ sin() is 5.17251045 *** Error 8.9 * sinpx ( 8.9 * PI) is -2.75025125 C++ sin() is 2.75025152 *** Error 9.0 * sinpx ( 9.0 * PI) is -0.00000000 C++ sin() is 0.00000029 9.1 * sinpx ( 9.1 * PI) is 2.81205465 C++ sin() is -2.81205437 *** Error 9.2 * sinpx ( 9.2 * PI) is 5.40762432 C++ sin() is -5.40762408 *** Error 9.3 * sinpx ( 9.3 * PI) is 7.52385805 C++ sin() is -7.52385787 *** Error 9.4 * sinpx ( 9.4 * PI) is 8.93993125 C++ sin() is -8.93993116 *** Error 9.5 * sinpx ( 9.5 * PI) is 9.50000000 C++ sin() is -9.50000000 *** Error 9.6 * sinpx ( 9.6 * PI) is 9.13014256 C++ sin() is -9.13014266 *** Error 9.7 * sinpx ( 9.7 * PI) is 7.84746485 C++ sin() is -7.84746504 *** Error 9.8 * sinpx ( 9.8 * PI) is 5.76029547 C++ sin() is -5.76029575 *** Error 9.9 * sinpx ( 9.9 * PI) is 3.05926824 C++ sin() is -3.05926858 *** Error 10.0 * sinpx (10.0 * PI) is 0.00000000 C++ sin() is -0.00000036 10.1 * sinpx (10.1 * PI) is -3.12107164 C++ sin() is 3.12107129 *** Error Press any key to continue . . .
Attachments (1)
Change History (2)
by , 8 years ago
Attachment: | Booster.cpp added |
---|
comment:1 by , 8 years ago
Resolution: | → fixed |
---|---|
Status: | new → closed |
Yes you're correct, fixed in Git develop.
However, note that:
If you want a similar/documented function, try sin_pi and cos_pi: http://www.boost.org/doc/libs/1_57_0/libs/math/doc/html/math_toolkit/powers/sin_pi.html