Deep dive into IEEE 754

本文共2009字,阅读完需要约8分钟。
版权声明: 知识共享-版权归属-相同方式共享 3.0 授权协议 | CC BY-SA 3.0 CN
展开

0.1 + 0.2 = 0.30000000000000004

Deep dive into IEEE 754

I: fractional2binary

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

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

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
# 仅能处理小于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的部分用普通的,表示整数的二进制形式来表示就行了。

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

1
2
3
4
5
6
7
8
9
10
11
12
13
14
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个比特吧。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
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:

1
2
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一样了。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
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中所说的「一小部分」。

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

1
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。

1676184124061

  • 注: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

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