0.1 + 0.2 = 0.30000000000000004

Deep dive into IEEE 754

I: fractional2binary

首先，让我们考虑如何在计算机的二进制世界中表示小数。毫无疑问，我们可以使用的方式来表示小数/非整数。比如：可以把所有表示成。

接下来让我们把这个过程自动化。

# 仅能处理小于0的小数
import fractions

f = fractions.Fraction("0.2")

step = 0
limits = 24
# exponent = 0
# mantissa = ""
ans = ""

while(step < limits):
    ans += str(int((f * 2) > 1))
    f = f * 2
    if f > 1:
        f -= 1
    step += 1

print(ans)

因为Python和C/C++一样，小数在内存中都是采用IEEE 754来表示，故我们这里选择了Python的库fractions来进行纯数学意义上的运算。

II: Intros Exponent

上面的程序有一个很大的问题，就是它只能表示小于1的小数。当然这个问题很容易解决：我们将大于1的部分用普通的，表示整数的二进制形式来表示就行了。

让我们稍微修改一下上面这个程序。

import fractions
from fractional2binary import fractional2binary

f = fractions.Fraction("3062.4770")

def interger2binary(interger):
    ans = ""
    while interger:
        ans += str(interger % 2)
        interger = interger // 2
    return ans[::-1]

if __name__ == "__main__":
    print("{}.{}".format(interger2binary(f - f % 1), fractional2binary(f % 1)))

另外一个问题出现了。如果按照这种储存方式，储存的空间是不定的：虽然在fractional2binary中我们可以通过控制limits来限制位数。但是在整数部分————恰恰是因为整数部分可以完美地被转换为二进制，所以不好去限制它转换后的位数。我们还需要一段空间来表示小数点在哪里，因为此，我们不知道整数部分有多长，也就不知道小数点的位置会有多大（范围）。

这时我们就不得不放弃一些东西了：我们需要放弃整数的准确性————当然在这之前肯定得先放弃小数的准确性。让我们把整数+小数转换成的二进制字符串强行缩短到一个长度，不妨先定为24个比特吧。

import fractions
from fractional2binary import fractional2binary

f = fractions.Fraction("3062.4770")

def interger2binary(interger):
    ans = ""
    while interger:
        ans += str(interger % 2)
        interger = interger // 2
    return ans[::-1]

def add1(string: str) -> str:
    # if string == '1':
    #     return '10'
    # if string[-1] == '0':
    #     return string[0:-1] + '1'
    # else:
    #     return add1(string[0:-1]) + '0'
    return bin(int(string, 2) + 1)[2:]

def limitsbinary(binary: str, limits: int):
    binary = binary.replace(".", "")
    ans = binary[0:limits]
    if binary[limits] == '1':
        ans = add1(ans)
        if len(ans) != limits:
            ans = ans[0:-1]
    assert(len(ans) == limits)
    return ans

if __name__ == "__main__":
    ans = "{}.{}".format(interger2binary(f - f % 1), fractional2binary(f % 1))
    print(ans)
    print(limitsbinary(ans, 24))

output:

101111110110.011110100001110010101100
101111110110011110100010

很不错！我们现在有了一个长度固定的二进制串，那么现在我们的问题就是如何表示出小数点的位置。

且慢，我们能用更数学的方式描述一下小数点的意义吗？

那样的话我们是不是就可以把这个数表示为这样了？

此时。

根据观察，应该说，exponent使用8个比特来表达就绰绰有余了。完全按照binary2interger的方法去计算，exponent的取值范围应该是。但是我们还需要负数exponent。解决负数的一个方法是使用补码来表达数字，但是在IEEE 754中我们选择了更加粗暴的方法，也就是直接对二进制表示的直接数减去127以达到可以表达的正数和负数数量基本相等。

这么做之后我们exponent的取值范围就变成了。

最后我们还需要一个标志位来标志我们的数的正负，一个比特足矣。

进行完这些后，我们可以说，我们发现/发明了自己的小数表达方式。总结一下就是：

• 一个float_tiger1218所占用的内存是的是比特，和int类型一样。
• 第一个比特表示数正负，其中0表示正数，1表示负数。
• 后八个比特表示指数。将后八个比特看成一个unsigned char后再减去127作为真正的指数。
• 最后二十三个比特表示需要乘的整数。将其看成一个无符号整数后再将其乘上就得到了无符号的整数。再乘上就是我们最终得到的浮点数。

我们不妨把这个标准称作Pre IEEE 754 / float tiger1218，作为后面重构的基础。

这是float2ieee754的Python脚本，可以看出在很多地方已经和IEEE 754一样了。

from intros_exponent import interger2binary
from fractional2binary import fractional2binary
import fractions

f = fractions.Fraction("3062.4770")

def add1(string: str) -> str:
    if string == '1':
        return '10'
    if string[-1] == '0':
        return string[0:-1] + '1'
    else:
        return add1(string[0:-1]) + '0'

def limitsbinary(string: str, bits: int) -> str:
    # print(string)
    ans = string[:bits + 1]
    if string[bits + 1] == '1':
        ans = add1(ans)
    ans = ans[1:bits + 1]
    return ans
    # ans = string[1:]
    

def sign(f):
    return str(int(f < 0))

def exponent(f):
    interger = interger2binary(f - f % 1)
    original_exponent = len(interger) - 1
    exponent = interger2binary(original_exponent + 127)
    return exponent

def fraction(f):
    interger = interger2binary(f - f % 1)
    fraction = fractional2binary(f % 1)
    tmp = interger + fraction
    # print(len(tmp))
    fraction = limitsbinary(tmp, 23)
    return fraction
    # print(len(fraction))

def pre_ieee_754(f):
    return sign(f) + exponent(f) + fraction(f)


if __name__ == "__main__":
    print(pre_ieee_754(f))

III: Reorganization

我们注意到上面这个Pre IEEE 754存在两个小问题：

1. 对于绝大多数（还有一小部分在后面会详细阐述）小数的表示总是会损失掉一个比特的信息熵。这个比特就是最前面一定是1。
2. exponent并不对称：也就是说，exponent大概率是负数————也是浪费储存空间和信息的地方。

解决第一个问题的方式很简单，就是不存储第一个1继而往后存一位。而解决第二个问题的方式，就是默认一串字符串（例如11100101010011010）的小数点打在第一位后面（就是1.1100101010011010），那样就能基本保证exponent对称。

将这两个方式结合起来，我们就得到了精简版的Pre IEEE 754.

我们不会把有着23比特的第三段串看成一个完整的整数，取而代之的，我们将其看作。其它两段的意义与作用不变。这样我们就完美解决了上面这两个问题。

最后我们得到的公式就是：

IV: Denormalized

考虑下面两个二进制浮点数：

将它们相减得到：

显然我们的Pre IEEE 754不能表示这个数，因为exponent的最小值是-127.

这就是在III中所说的「一小部分」。

这个值显然不能近似为零。不然下面这个语句就会出现问题：

if(b != a) z = 1 / (b - a)

因此我们需要引入Denormalized数来解决这个问题。考虑到exponent位的最小值-127用到时就已经很接近这种情况了，我们就单独挑出exponent=-127，也就是original_exponent=0（在减去127之前的原始exponent）这个情况作为特殊情况处理。

我们定义当exponent=-127 / original_exponent=0时，我们把有着23比特的第三段串看作。并且此时（即使）仍然作乘计算。

V: Particular Value

首先我们来讨论在Pre IEEE 754中0应该如何表示。

运行IV中给出的脚本，0应该是 0111111000000000000000000000000，也就是。然而float中的0和int中的0不一样本来就让人迷惑，更别提在memset的时候将内存清空为0甚至不是0了。因此我们需要定义一下：

• 在original_exponent=0时，若第三段也全为0,令表示的数为0。

仔细一看这与IV并不冲突，因为。

然后，为了方便，我们单独剖出original_exponent=255时的情况，此时exponent=128. 定义在此时，第三段若为0,则该数表示无穷大。若不为零，则表示NaN。

• 注：QNaN指Quiet NaN, SNaN指Signaling NaN. 两者的差别见此处。

VI: Epilogue

把这些都加进Pre IEEE 754后，我们就能得到IEEE 754关于float32的实现了。当然，double的实现也基本差不多。双精度浮点数（也就是double）的Sign, Exponent和Fraction的分配是：1+11+52。

Reference

IEEE 754浮点表示法详解

IEEE 754 - Wikipedia

计算机中的浮点数在数轴上分布均匀吗？ - 知乎