为什么 function floor 在这种情况下会给出不同的结果?

Why is the function floor giving different results in this case?

在此示例中,floor 的行为不同,我不明白为什么:

printf("floor(34000000.535 * 100 + 0.5) : %lf \n", floor(34000000.535 * 100 + 0.5));
printf("floor(33000000.535 * 100 + 0.5) : %lf \n", floor(33000000.535 * 100 + 0.5));

这段代码的输出是:

floor(34000000.535 * 100 + 0.5) : 3400000053.000000
floor(33000000.535 * 100 + 0.5) : 3300000054.000000

为什么第一个结果不等于我们预期的 3400000054.0?

问题已经回答here

基本的浮点数只是近似值。如果我们有这样的程序:

float a = 0.2 + 0.3;
float b = 0.25 + 0.25;

if (a == b) {
    //might happen
}
if (a != b) {
    // also might happen
}

唯一可以保证的是a-b比较小

C 中的

double 并不表示所有可以用文本表达的可能数字。

double 通常可以表示大约 264 个不同的数字。当 double 被编码为 binary floating point 数字时,34000000.53533000000.535 都不在该集合中。而是使用最接近的可表示数字。

Text             34000000.535
closest double   34000000.534999996423...
Text             33000000.535
closest double   33000000.535000000149...

double 作为二进制浮点数,乘以非 2 的幂,如 100.0,会引入额外的舍入差异。然而在这些情况下,它仍然会导致产品,一个刚好高于 xxx.5,另一个低于 xxx.5。

添加 0.5,一个简单的 2 次幂,不会导致舍入问题,因为与 3x00000053.5 相比,该值并不极端。

看到更高打印精度的中间结果很好地展示了典型的逐步过程。

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

 void fma_test(double a, double b, double c) {
   int n = DBL_DIG + 3;
   printf("a b c      %.*e %.*e %.*e\n", n, a, n, b, n, c);
   printf("a*b        %.*e\n", n, a*b);
   printf("a*b+c      %.*e\n", n, a*b+c);
   printf("a*b+c      %.*e\n", n, floor(a*b+c));
   puts("");
 }

int main(void) {
  fma_test(34000000.535, 100, 0.5);
  fma_test(33000000.535, 100, 0.5);
}

输出

a b c      3.400000053499999642e+07 1.000000000000000000e+02 5.000000000000000000e-01
a*b        3.400000053499999523e+09
a*b+c      3.400000053999999523e+09
a*b+c      3.400000053000000000e+09

a b c      3.300000053500000015e+07 1.000000000000000000e+02 5.000000000000000000e-01
a*b        3.300000053500000000e+09
a*b+c      3.300000054000000000e+09
a*b+c      3.300000054000000000e+09

问题比这个简单的答案更复杂,因为各种平台可以 1) 使用更高精度的数学,如 long double 或 2) 很少使用 decimal 浮点数double。所以代码的结果可能会有所不同。

使用 code that shows the representation of floats in memory as sum of terms,我们得到:

main()
{
    float x=floor(34000000.535 * 100 + 0.5);
    float y=floor(33000000.535 * 100 + 0.5);
    xx(&x);
    xx(&y);
    yy(x);
    yy(y);
}

此代码将在两种情况下输出 floor 返回值在内存中的表示形式。

使用 bc 计算器,我们可以看到近似值确实不错,但由于楼层表示背后的数学原因,存在一些扰动。

注意:我在bc中设置了scale=20,也就是说,每次中间计算保留点后20位。

./a.out
1ST NUMBER=>    sign:0 exponent:1 0 0 1 1 1 1 fraction:0 1 0 0 1 0 1 0 1 0 1 0 0 1 1 1 1 1 1 0 0 0 1 0
2ND NUMBER=>    sign:0 exponent:1 0 0 1 1 1 1 fraction:0 1 0 0 0 1 0 0 1 0 1 1 0 0 1 0 0 0 0 0 0 0 0 1
1ST NUMBER=>    positive ( 1+1/(2) +1/(16) +1/(64) +1/(256) +1/(1024) +1/(8192) +1/(16384) +1/(32768) +1/(65536) +1/(131072) +1/(262144) +1/(4194304) )*2^31
2ND NUMBER=>    positive ( 1+1/(2) +1/(32) +1/(256) +1/(1024) +1/(2048) +1/(16384) +1/(8388608) )*2^31
@ bc
scale=20
( 1+1/(2) +1/(16) +1/(64) +1/(256) +1/(1024) +1/(8192) +1/(16384) +1/(32768) +1/(65536) +1/(131072) +1/(262144) +1/(4194304) )*2^31
3399999999.99999999999463129088
 ( 1+1/(2) +1/(32) +1/(256) +1/(1024) +1/(2048) +1/(16384) +1/(8388608) )*2^31
3299999999.99999999999731564544