张银峰的编程课堂

浮点型比较

10个0.1相加不等于1.0

浮点数因为自身在计算机中的表达方式,会产生一些有趣的表现。比如示例程序计算10个0.1的和,再与1.0做相等比较,但程序给出了不相等的答案。

#include <stdio.h>

int main()
{
    float sum = 0.0f;
    float one = 1.0f;

    printf("index \t %%.1f \t %%.4f \t\t %%.6f \t\t %%.8f\n");
    printf("============================================================\n");

    for (int i = 0; i < 10; i++)
    {
        sum += 0.1f;
        printf("%3d \t %.1f \t %.4f \t %.6f \t %.8f \t\n", i, sum, sum, sum, sum);
    }

    printf("\n");
    printf(sum == one ? "Equal!" : "Not Equal!");

    return 0;
}

glimix.com

程序在格式化方面有些小技巧,格式符%3d表示以3个字符宽度输出,不足3个字符时,前面填充空白。因为表头文本"index"是5个字符的宽度,这使得索引值恰好落在了这一列的中间。在显示精度<=6位时,格式化输出使得sum表现出等于1.0的假象;但在执行如比较之类的运算时,实际值1.00000012显然不等于1.0,所以程序输出了"Not Equal!"。

glimix.com

探索误差来源

在找到答案之前,我们需要先掌握如何将小数转换为对应的二进制表示:方法是将小数乘以2,当结果>=1时取1,反之取0;再用结果的小数部分继续乘以2,直到结果为0,否则一直计算下去。对0.1的计算如下:

0.1                     结果
==================================
0.1 * 2 = 0.2           0
0.2 * 2 = 0.4           0
0.4 * 2 = 0.8           0
0.8 * 2 = 1.6           1
0.6 * 2 = 1.2           1
0.2 * 2 = 0.4           0
0.4 * 2 = 0.8           0
0.8 * 2 = 1.6           1
0.6 * 2 = 1.2           1
0.2 * 2 = 0.4           0
...
==================================
0.1 = .0001100110

让我们将0.1的二进制.0001100110反向验证一下。将二进制小数转换为十进制小数时,对于小数部分的二进制位,从左到右,每位值分别是1/2n,其中(n=1,2,3...),即小数点后第1位是1/2,第2位是1/4....,之后将这些被置位的值相加即可。

    1/2     1/4     1/8     1/16    1/32        1/64    1/128   1/256       1/512        1/1024
-----------------------------------------------------------------------------------------------
    0       0       0       1       1           0       0       1           1            0
                            0.0625  0.03125                     0.00390625  0.001953125
-----------------------------------------------------------------------------------------------
.0001100110 = 0.0625 + 0.03125 + 0.00390625 + 0.001953125 = 0.099609375

可以看到(0.099609375)相当接近原始值(0.1),接下来我们看一下0.5。

0.5                     结果
==================================
0.5 * 2 = 1.0           1
0.0 * 2 = 0.0           0
==================================
0.5 = .100000000

可以看到,0.5不会像0.1那样无限下去。我们将0.5的二进制.10反向验证一下。

    1/2     1/4     1/8
--------------------------
    1       0       0
    0.5
--------------------------
.100000000 = 0.5

通过这两个示例可以看出,有些小数(0.5)可以准确的表达;而一些小数(0.1)只能近似表达。对于浮点数的比较,这会带来一些有趣的结果。

#include <stdio.h>

int main()
{
    float f1 = 0.0f,
          f2 = 0.0f,
          f3 = 0.0f;

    printf("index \t f1 \t\t f2 \t\t f3 \n");
    printf("============================================================\n");

    for (int i = 0; i < 10; i++)
    {
        f1 += 0.2f;
        f2 += 0.5f;
        f3 += 0.3125f;

        printf("%3d \t %.8f \t %.8f \t %.8f \n", i, f1, f2, f3);
    }

    printf("============================================================\n");
    printf("   \t ");

    if (f1 == 2.0) printf("EQ(2.0)\t "); else printf("NEQ(2.0)\t ");
    if (f2 == 5.0) printf("EQ(5.0)\t "); else printf("NEQ(5.0)\t ");
    if (f3 == 3.125) printf("EQ(3.125)\n"); else printf("NEQ(3.125)\n");

    return 0;
}

glimix.com

有点出乎意料吧!

  • 看似温和的0.2由于只能近似表达,所以f1==2.0为假。
  • 在上面的分析中,0.5可以精确定表达,因此f2==5.0为真。
  • 感觉只能近似的0.3125,其实由1/4+1/16=0.25+0.0625组成,因为能精确表达,所以f3==3.125为真。

f1在理论与实际上存在精确值2.0与近似等于2.0两种状态,使得浮点数的相等比较操作,不应该使用==、!=运算符(但可以使用>=、<=、>、<运算符)。这时我们需要一个折中的表达:对于浮点数a、b和参考值delta,如果满足不等式:(a - delta) < b < (a + delta), 就可以认为a等于b;这种手法移除了等于条件,使得在浮点数比较中仅使用<和>运算符即可

==============a==============      // a
==============b==============      // b在同一位置
=============ab==============      // a在b的左侧
==============ba=============      // a在b的右则
=== a-delta   b   delta+a ===      // 满足不等式时 a == b

下面的示例程序通过不停的收缩delta的来测试相等性,可以看到delta范围的宽松粒度(精度)直接影响比较结果。

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

int main()
{
    float f1 = 1.002f;
    float f2 = 1.000f;
    float delta = 0.3f;

    printf("      F1(%.3f)    EQ    F2(%.3f)   ? \n\n", f1, f2);
    printf("===============================================\n");
    printf("delta \t f1-delta \t f1+delta \t result\n");
    printf("===============================================\n");

    for (int i = 0; i < 5; i++)
    {
        printf("%.5f", delta);
        printf("\t %.5f \t %.5f", (f1 - delta), (f1 + delta));

        if ((f1 - delta) < f2 && f2 < (f1 + delta))
            printf("\t EQ\n");
        else
            printf("\t NEQ\n");

        delta /= 10.0f;
    }

    return 0;
}

glimix.com

当delta=0.3时,可以认为对f1的测试有30%的上下浮动范围;当delta=0.03时,则有3%的浮动范围。随着范围的逐步缩小,代表对精度要求的提高,整个过程就是浮动范围从粗略到精细,匹配度由高到低。(理解delta的另一种方式是:你向外借出了100RMB,然后告诉对方还70就等同于全还完了,这就是30%的忽略范围;相应的,你要求对方偿还90,即delta更小,则偿还值更接近借出值。)程序中的if语句可以简化一下,因为两个值都大于0,所以它们的差如果小于(注意不是小于等于)delta,那说明位于可接受范围内,此时就认为两个值是相等的。

if ((f1 - f2) < delta)
    printf("\t EQ\n");

考虑一下f1-f2小于0的情况,假如差的绝对值也小于delta,也表示两个值相等,所以更准确的表达就是使用浮点数绝对值函数fabs()来改写程序。

if (fabs(f1 - f2) < delta)
    printf("\t EQ\n");

或许自己输入感觉一下会更切贴。下面的正确提示你输入pi的值,当你的答案与参考值3.14159之间的误差大于delta 0.0001时,将一直重复,直到输入满足精度。

#include <stdio.h>
#include <math.h>   // for fabs()

int main()
{
    const double ANSWER = 3.14159;
    const double epsilon = 0.0001;

    double response;
    double differ;

    printf("What is the value of pi?\n");

    while (1)
    {
        scanf("%lf", &response);
        differ = fabs(response - ANSWER);

        if (differ > epsilon)
        {
            printf("fabs(%lf - 3.14159) = %lf > 0.0001 => Try again!\n", response, differ);
        }
        else
        {
            printf("fabs(%lf - 3.14159) = %lf < 0.0001 => Close enough!\n", response, differ);
            break;
        }
    }

    return 0;
}

glimix.com

我们可以使用头文件float.h中的预定义常量值FLT_EPSILON作为相等比较的参考值。

#include <stdio.h>
#include <math.h>   // for fabs()
#include <float.h>  // for FLT_EPSILON

int main()
{
    printf("%.12f\n", FLT_EPSILON);

    double sum = 1.0000001;

    if (fabs(sum - 1.0) < FLT_EPSILON)
        printf("equal\n");
    else
        printf("not equal\n");

    return 0;
}

glimix.com