为什么在这种情况下功能层会给出不同的结果？答案

【问题标题】：Why is the function floor giving different results in this case?为什么在这种情况下功能层会给出不同的结果？
【发布时间】：2016-12-30 10:41:48
【问题描述】：

在这个例子中，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？

【问题讨论】：

当然是舍入问题，如果您需要精确计算，请使用任意精度数学库（例如 gmp）。
因为你正在处理二进制浮点数。
@cHao 请注意，即使使用 decimal floating point，当值接近格式的精度限制时，也会出现此类问题。这就是浮点算术的特性，在任何基础上。十进制文本和二进制 double 肯定更常见。

标签： c floating-point rounding

【解决方案1】：

C 中的double 并不代表可以在文本中表达的所有可能数字。

double 通常可以表示大约 2⁶⁴ 个不同的数字。当 double 编码为 binary floating point 数字时，34000000.535 和 33000000.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) 很少使用 十进制 浮点 double。所以代码的结果可能会有所不同。

【讨论】：

chux，你说有 2^64 个不同的数字表示是错误的。实际上，还有很多。因为在这 64 个里面有指数部分，与我在回答中分析的 float 的情况相同。编译器可以将数字 100 表示为 10e2 或 100e1。指数部分没有关于编译器应该如何设置它的property。请参阅我的答案，该答案打印了 double/long double/float 的表示方式。因此，当您编写 x=100 时，您既不能确定编译器将给出什么表示，也不能确定浮点硬件
@alinsoar 这个答案没有说有 2^64 个不同的数字。在您进行不必要的编辑之前，这个答案说“大约 2^64 个不同的数字”。实际上，有不少很多的值。 2^64 个不同的位模式中有一部分是非数字，它只会将值模式计数减少到 0.1% 以下。您将 100 打印为 10e2 或 100e1 的指数讨论（仍然只是 2^64 模式的计数之一）是无关紧要的。不要编辑此答案 - 编辑和评论意图会误导正确的意图。
你可以使用“小于”，因为它比“约”更正确。
@alinsoar “大约 2^64 个不同的数字”和“少于 2^64 个不同的数字”之间的差异与 OP 的问题无关。 “约”具有表示接近的优点，“小于”则没有。可以使用“近 2^64 个不同的数字”来传达这两种想法，但同样，这不是 OP 问题的核心。

【解决方案2】：

问题已经回答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比较小。

【讨论】：

【解决方案3】：

使用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

【讨论】：

我认为您在这里混淆了一些东西。在最初的问题中有两个不同的数字（34.... vs 33....），这是问题的根本原因。你有两次相同的数字，当然在应用floor() 时会给出相同的结果。这里不考虑 printf 实现。
@Ctx 谢谢。我没有注意到数字不同。已更正。
void main 调用未定义的行为。
@RolandIllig 这里是关于输出浮点数的表示，主要返回在这个问题中并不重要。此外，不是返回类型是 gnu-extension，但这与返回错误代码无关。我认为您很难理解问题所在。
为什么设置 scale=20 而不是匹配具有 24 位二进制精度的 typical float？使用 scale=24，bc 报告 scale=24 ( 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 3400000000.000000000000000000000000 与这些分数的数学总和匹配。