math problem for experts

[miltmobley]

While testing python code to rotate 3d points, I tried:

cos_ = np.cos((90.0/360.0) * (2.0*np.pi))

cos(90) was reported as 6.123234e-17, instead of correct value of zero.

tried both math.pi and numpy.pi with same result

os = windows 10, cpu = Intel i5-7200U
python = 3.8.10 provided by anaconda spyder 5.4.3

The same behavior occurred on Linux Mint os, using python 3.8.10 provided by anaconda spyder 4.5.3, cpu = AMD whatever

A simple c program compiled by gcc and executed on the Linux machine reported cos(90) as zero.

What causes this problem? What can be done about it?

[jefsummers]

Classic question regarding floating point math. It is not exact, though 10 to the -17 is about as close to zero as it gets. Suggest rounding to the nearest 10 to the minus 14 and being happy with the result.

[deanhystad]

I would like to see the C program. Should have the exact same problem.

I ran this program on Windows

#define _USE_MATH_DEFINES

#include <iostream> 
#include <cmath>

using namespace std;

int main()
{
    cout << cos(M_PI / 2);;

    return 0;
}

And the output

Output:
6.12323e-17

[miltmobley]

I had trouble making this an attachment, so here is the Linux code:

#include <stdio.h>
#include <stdlib.h>
#include <math.h>

int main(void) {
double c = cos((90.0/360.0) * (2.0 * M_PI))
double s = sin(0)
printf("cos(90) = %lf, sin(0) = %lf\n", c, s)
return EXIT_SUCCESS
}

/*
To compile and run:
gcc trig_test.c -o trig_test
chmod +x trig_test
./trig_test

Output is:
cos(90) = 0.00000, sin(0) = 0.000000

So your cpp compiler is wrong too?

My concern is how many other angles might I have to treat specially?

If I rotate the point (100.0, 0.0, 0.0) by 90 degrees, the resulting x value is 6.<15 bits>e-15 away from the correct value.
So an error band of 10**(-14) might not be appropriate in such a case.

One would think that:

(1) Python developers would save themselves some effort by just using the gcc math library functions (at least for non bignums),
and that the gcc functions would save processing time by using available floating point hardware, which on my machines is 64 bits.

(2) As it appears they are not, the Python developers could save processing time by special casing multiples of pi/2 (90 degrees),
which are the well known maxima, minima and zero crossing points of the sin and cos functions,
instead of blindly evaluating some supposedly converging power series.

[deanhystad]

My cpp compiler is just fine. As already mentioned, there are limitations when performing floating point operations. 64 bits puts a limit on how many different float numbers can be represented, but there are an infinite number of floats between -1 and 1.

I ran your program and it prints:

Output:
cos(90) = 0.000000, sin(0) = 0.000000

But when I run it in the debugger,

Output:
c=6.1232339957367660e-17

The print statement is rounding down the very small float number. c is not zero, it only prints that way.

I modified the program to make the small number bigger.

#include <stdio.h>
#include <stdlib.h>
#include <math.h>

int main()
{
		double c = cos((90.0 / 360.0) * (2.0 * M_PI)) * 1e17;
		double s = sin(0) * 1e17;
		printf("cos(90) = %f, sin(0) = %f\n", c, s);
		return EXIT_SUCCESS;
}

If c == 0, multiplying c by a really large number should still be zero. But this is the output.

Output:
cos(90) = 6.123234, sin(0) = 0.000000

Notice the sin(0) is still zero.

As far as identifying special cases like pi/2, that suffers from the same precision limitations as any other float operation. To start, there is no way to get pi. Pi is an irrational number. It cannot be represented as a numerical constant, nor as a fraction. If we cannot identify pi, how can we identify pi/2? Even if we could identify pi pr pi/2 or "important" fractions of pi, how do we get those exact values? If I do some divisions and multiplications that should result in 0.5, there is no guarantee the result will be 0.5 because some of the intermediate results cannot be accurately represented as a double. Something really close to 0.5 * pi != pi / 2

[Gribouillis]

The floating point numbers limitations are not related to the Python language. They are an intrinsinc limitation of the fact that computations are made with a finite precision. All the other usual programming languages have the same limitations and most of them use the same IEEE754 representation of floating point numbers. This includes C++ as well as Python. The only way to get more precise results is to use a multiprecision library such as gmp, but these computations have the same limitations, only at a different order of magnitude.