Floating-point arithmetic is very different from real math. In particular, floating point arithmetic is not associative or distributive. Therefore, mathematically equivalent expressions are not necessarily equivalent when evaluated using finite point precision floating point arithmetic. This happens even when individual operations, such as addition, subtraction, multiplication, division and square root, give correctly rounded results, as required by IEEE-754, the standard standard .
In this case, g() will be much more accurate than f() on average, as well as for the specific case x = 500 . The reason is that f() suffers from subtractive annulment. This happens when you effectively subtract two floating point numbers that are almost the same in magnitude. Leading digits are canceled during subtraction, leaving only a few remaining lagging digits that are included in subsequent calculations. In addition, any rounding error accumulated in the final digits of the subtracted source operands can be increased by subsequent calculation. An extended explanation with an example can be found in this Wikipedia article .
In this case, sqrt(x+1) and sqrt(x) have almost the same value, in particular, as the value of x increases. Using the example x = 500 and using IEEE-754 precision arithmetic, we find:
x = 500 f(x) = 0x1.659ae0p+3 (11.175156) reference = 0x1.659798p+3 (11.174755) x = 500 g(x) = 0x1.659798p+3 (11.174755) reference = 0x1.659798p+3 (11.174755)
The error in f(500) is 420 ulps , and g(500) provides the correct rounded result with one precision.
source share