Logo Questions Linux Laravel Mysql Ubuntu Git Menu
 

Can cos(a) ever equal 0 in floating point

Given that PI/2 can never be represented precisely in floating point, is it safe to assume that cos(a) can never return exact zero?

If this is the case, then the following pseudo-code will never enter the block (and it could be safely removed):

...
y = h/cos(a);
if (!isfinite(a))
{
  // handle infinite y
}
like image 963
Chris Leishman Avatar asked Jul 31 '11 10:07

Chris Leishman


2 Answers

Other than zero, the double precision value that comes closest to an exact multiple of π/2 is 6381956970095103 * 2^797, which is equal to:

(an odd integer) * π/2 + 2.983942503748063...e−19

Thus, for all double-precision values x, we have the bound:

|cos(x)| >= cos(2.983942503748063...e−19)

Note that this is a bound on the mathematically exact value, not on the value returned by the library function cos. On a platform with a good-quality math library, this bound is sufficiently good that we can say that cos(x) is not zero for any double-precision x. In fact, it turns out that this is not unique to double; this property holds for all IEEE-754 basic types, if cos is faithfully rounded.

However, that's not to say that this could never occur on a platform that had a spectacularly poor implementation of trigonometric argument reduction.

Even more importantly, it's critical to note that in your example y can be infinite without cos(a) being zero:

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

int main(int argc, char *argv[]) {
    double a = 0x1.6ac5b262ca1ffp+849;
    double h = 0x1.0p1022;
    printf("cos(a) = %g\n", cos(a));
    printf("h/cos(a) = %g\n", h/cos(a));
    return 0;
}

compile and run:

scanon$ clang example.c && ./a.out
cos(a) = -4.68717e-19
h/cos(a) = -inf
like image 192
Stephen Canon Avatar answered Sep 23 '22 11:09

Stephen Canon


Zero is one of several values that can be represented exactly. Many systems have a lookup table for common values of sin and cos, so it's not inconceivable that exactly zero could be returned.

But you are safer using a delta compare, before performing the division:

if (Abs(cos(a)) < 0.0000001)
{

}
like image 41
Mitch Wheat Avatar answered Sep 23 '22 11:09

Mitch Wheat