Problem with acos returning NaN, and I am sure parameter is in [-1;1]...

Tuesday, February 10 2004 @ 12:46 PM GMT

Contributed by: unix

Sometimes a program returns false results, the problem has been tracked down to acos returning NaN, but the math behind the algorithm proves the parameter should be correctly bounded in [-1; 1]. What's up ?

Even if you are sure of the theory behind an algorithm implementation, doesn't mean that some rounding may not break your assumptions, especially when using the double data type.
Here an example (code extracted from a color conversion algorithm), assuming A, B, C are the results of other calculations of no direct interest here:


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

int main(int argc, char **argv)
{
        double  A = 1.333333333333333   ;       //   4.0 / 3.0       ;
        double  B = 0.33333333333333343 ;       //      theorically 1.0 / 3.0 but due to calculation errors...     
        double  C = 1.3333333333333335  ;       //      theorically 4.0 / 3.0 but due to calculation errors...    

        double x = (A - 2*B + C) / 3.0  ;
        double y = (A - C) / 2.0        ;

        double r = sqrt( 9 * x * x - 3 * y * y) ;
	 // in the real code there is a if (r > .0) before the divide by r ... 
        // but we don't care here...
	
        double t = acos( (3 * x) / r )  ;	// supposed to be 1.0 BUT is 1.0000000000000002

        printf("A %f B %f C %f x %f y %f r %f t %fn",A,B,C,x,y,r,t);

        return 0        ;
}

A 1.333333 B 0.333333 C 1.333333 x 0.666667 y -0.000000 r 2.000000 t nan
The solution is application specific, and even in an application you could apply different methods according to the local context...
Here you could just "clamp" the value to be in the valid range: MIN(1.0, MAX(-1.0, (3 * x) / r )).
Note that similar precaution should be taken in the sqrt call against a potential "minus zero" value...

0 comments



http://www.macosxguru.net/article.php?story=20040210124637626