Subject: [Boost-bugs] [Boost C++ Libraries] #10905: sinpx function in boost/math/special_functions/gamma.hpp returns incorrect results for positive z arguments
From: Boost C++ Libraries (noreply_at_[hidden])
Date: 2014-12-29 09:09:52
#10905: sinpx function in boost/math/special_functions/gamma.hpp returns incorrect
results for positive z arguments
------------------------------------+-------------------------
Reporter: Brent Rector <brent@…> | Owner: johnmaddock
Type: Bugs | Status: new
Milestone: To Be Determined | Component: math
Version: Boost 1.57.0 | Severity: Showstopper
Keywords: |
------------------------------------+-------------------------
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 . . .
}}}
-- Ticket URL: <https://svn.boost.org/trac/boost/ticket/10905> Boost C++ Libraries <http://www.boost.org/> Boost provides free peer-reviewed portable C++ source libraries.
This archive was generated by hypermail 2.1.7 : 2017-02-16 18:50:17 UTC