Golang 笔记之深刻浮点数

引言

下面的一段简单程序 0.3 + 0.6 结果是什么？

1 var f1 float64 = 0.32 var f2 float64 = 0.63 fmt.Println(f1 + f2)

有人会天真的认为是0.9，但实际输出倒是0.8999999999999999（go 1.13.5）

问题在于大多数小数表示成二进制以后是近似且无限的。
以0.1为例。它多是你能想到的最简单的十进制之一，可是二进制看起来却很是复杂：0.0001100110011001100…
其是一串连续循环无限的数字（涉及到10进制转换为2进制，暂不介绍）。
结果的荒诞性告诉咱们，必须深刻理解浮点数在计算机中的存储方式及其性质，才能正确处理数字的计算。
golang 与其余不少语言（C、C++、Python…）同样，使用了IEEE-754标准存储浮点数。

IEEE-754 如何存储浮点数

IEEE-754规范使用特殊的以2为基数的科学表示法表示浮点数。

32位的单精度浮点数 与 64位的双精度浮点数的差别

符号位：1 为 负数， 0 为正数。
指数位：存储 指数加上偏移量，偏移量是为了表达负数而设计的。
小数位：存储系数的小数位的准确或者最接近的值。

以 数字 0.085 为例。

小位数的表达方式

以0.36 为例:
010 1110 0001 0100 0111 1011 = 0.36 (第一位数字表明1/2,第二位数字是1/4…，0.36 是全部位相加)
分解后的计算步骤为:

go语言显示浮点数 - 验证以前的理论

接下来用一个案例有助于咱们理解并验证IEEE-754 浮点数的表示方式。

math.Float32bits 能够为咱们打印出32位数据的二进制表示。(注：math.Float64bits能够打印64位数据的二进制)

下面的go代码将输出0.085的浮点数二进制表达，而且为了验证以前理论的正确性，根据二进制表示反向推导出其所表示的原始十进制0.085

 

输出：代表咱们对于浮点数的理解正确。

1 Starting Number: 0.0850002 Bit Pattern: 0 | 0111 1011 | 010 1110 0001 0100 0111 10113 Sign: 0 Exponent: 123 (-4) Mantissa: 0.360000 Value: 0.085000

经典问题：如何判断一个浮点数其实存储的是整数

下面是一个有趣的问题，如何判断一个浮点数其实存储的是整数？

思考10秒钟…

下面是一段判断浮点数是否为整数的go代码实现，咱们接下来逐行分析函数。
它能够加深对于浮点数的理解

 

一、要保证是整数，一个重要的条件是必需要指数位大于127，若是指数位为127，表明指数为0. 指数位大于127，表明指数大于0， 反之小于0.

下面咱们以数字234523为例子：

第一步,计算指数。因为 多减去了23，因此在第一个判断中 判断条件为 exponent < -23
exponent := int(bits >> 23) - bias - 23

第二步，
(bits & ((1 << 23) - 1)) 计算小数位。

| (1 << 23) 表明 将1加在前方。
1 + 小数 = 系数。

 第三步，计算intTest 只有当指数的倍数能够弥补最小的小数位的时候，才是一个整数。
以下，指数是17位，其不可以弥补最后6位的小数。即不能弥补1/2^18 的小数。
因为2^18位以后为0.因此是整数。

golang decimal 包详解

要理解decimal包，首先须要知道两个重要的概念，Normal number、denormal (or subnormal) number 以及精度。

1.概念：Normal number and denormal (or subnormal) number

wiki的解释是：

 

什么意思呢？在IEEE-754中指数位有一个偏移量，偏移量是为了表达负数而设计的。好比单精度中的0.085，实际的指数是 -3， 存储到指数位是123。
因此表达的负数就是有上限的。这个上限就是2^-126。若是比这个负数还要小，例如2^-127,这个时候应该表达为0.1 * 2 ^ -126. 这时系数变为了避免是1为前导的数，这个数就叫作denormal (or subnormal) number。
正常的系数是以1为前导的数就叫作Normal number。

2.概念：精度

精度是一个很是复杂的概念，在这里笔者讨论的是2进制浮点数的10进制精度。
精度为d表示的是在一个范围内，若是咱们将d位10进制（按照科学计数法表达）转换为二进制。再将二进制转换为d位10进制。数据不损失意味着在此范围内是有d精度的。
精度的缘由在于，数据在进制之间相互转换时，是不可以精准匹配的，而是匹配到一个最近的数。如图所示：

                                          精度转换

 

在这里暂时不深刻探讨，而是给出结论：（注：精度是动态变化的，不一样的范围可能有不一样的精度。这是因为 2的幂 与 10的幂之间的交错是不一样的。）

float32的精度为6-8位，

float64的精度为15-17位

目前使用比较多的精准操做浮点数的decimal包是shopspring/decimal。连接:https://github.com/shopspring/decimal

decimal包使用math/big包存储大整数并进行大整数的计算。

好比对于字符串 “123.45” 咱们能够将其转换为12345这个大整数，以及-2表明指数。参考decimal结构体：

在本文中，笔者不会探讨math/big是如何进行大整数运算的，而是探讨decimal包一个很是重要的函数：

NewFromFloat(value float64) Decimal

其主要调用了下面的函数：

 

此函数会将浮点数转换为Decimal结构。
读者想象一下这个问题：若是存储到浮点数中的值（例如0.1）自己就是一个近似值，为何decimal包可以解决计算的准确性？
缘由在于，deciimal包能够精准的将一个浮点数转换为10进制。这就是NewFromFloat为咱们作的事情。
下面我将对此函数作逐行分析。

 

第5行：剥离出IEEE浮点数的指数位
exp := int(bits>>flt.mantbits) & (1<<flt.expbits - 1)

第6行：剥离出浮点数的系数的小数位
mant := bits & (uint64(1)<<flt.mantbits - 1)

第7行：若是是指数位为0，表明浮点数是denormal (or subnormal) number；
默认状况下会在mant以前加上1，由于mant只是系数的小数，在前面加上1后，表明真正的小数位。
如今 mant = IEEE浮点数系数 * 2^53

第13行：加上偏移量，exp如今表明真正的指数。
第14行：引入了一个中间结构decimal

第15行：调用d.Assign(mant) , 将mant做为10进制数，存起来。
10进制数的每一位都做为一个字符存储到 decimal的byte数组中

第16行：调用shift函数，这个函数很是难理解。

 

此函数的功能是为了获取此浮点数表明的10进制数据的整数位个数以及小数位个数，此函数的完整证实附后。（注1)
exp是真实的指数，其也是可以覆盖小数部分2进制位的个数。（参考前面如何判断浮点数是整数）
exp - int(flt.mantbits)表明不能被exp覆盖的2进制位的个数
若是exp - int(flt.mantbits) > 0 表明exp可以彻底覆盖小数位 所以 浮点数是一个很是大的整数，这时会调用leftShift(a, uint(k))。不然将调用rightShift(a, uint(-k)), 常规rightShift会调用得更多。所以咱们来看看rightShift函数的实现。

第5行：此for循环将计算浮点数10进制表示的小数部分的有效位为 r-1 。
n >> k 是一个重要的衡量指标，表明了小数部分与整数部分的分割。此函数的完整证实附后。（注1)

第21行：此时整数部分所占的有效位数为a.dp -=（r-1）
第24行：这两个循环作了2件事情：
一、计算10进制表示的有效位数
二、将10进制表示存入bytes数组中。例如对于浮点数64.125，如今byte数组存储的前5位就是64125

 

继续回到newFromFloat函数，第18行，调用了roundShortest函数，
此函数很是关键。其会将浮点数转换为离其最近的十进制数。
这是为何decimal.NewFromFloat(0.1)可以精准表达0.1的缘由。

参考上面的精度，此函数主要考察了2的幂与10的幂之间的交错关系。四舍五入到最接近的10进制值。
此函数实质实现的是Grisu3 算法,有想深刻了解的能够去看看论文。笔者在这里提示几点：
一、2^exp <= d < 10^dp。
二、10进制数之间至少相聚10^(dp-nd)
三、2的幂之间的最小间距至少为2^(exp-mantbits)
四、何时d就是最接近2进制的10进制数？
若是10^(dp-nd) > 2^(exp-mantbits)，代表 当十进制降低一个最小位数时，匹配到的是更小的数字value - 2^(exp-mantbits)，因此d就是最接近浮点数的10进制数。

 

 

继续回到newFromFloat函数，第19行 若是精度小于19，是位于int64范围内的，可使用快速路径，不然使用math/big包进行赋值操做，效率稍微要慢一些。
第36行，正常状况几乎不会发生。若是setstring在异常的状况下会调用NewFromFloatWithExponent 指定精度进行四舍五入截断。

注一：快速的获取一个浮点数表明的十进制

以典型的数字64.125 为例 ， 它能够被浮点数二进制精准表达为：
Bit Patterns: 0 | 10000000101 | 0000000010000000000000000000000000000000000000000000
Sign: 0 | Exponent: 1029 (6) | Mantissa: 0.001953

即 64.125 = 1.001953125 * 2^6
注意观察浮点数的小数位在第九位有1, 表明2^-9 即 0.001953125.

咱们在浮点数的小数位前 附上数字1，10000000010000000000000000000000000000000000000000000 表明其为1 / 2^0 .

此时咱们能够认为这个数表明的是1.001953125. 那么这样长的二进制数变为10进制又是多少呢:4512395720392704。

即 1.001953125 = 4512395720392704 * 2^(-52)

因此64.125 = 4512395720392704 * 2^(-52) * 2^6 = 4512395720392704 * 2^(-46)
在这里，有一种重要的等式。即 (2 ^ -46) 等价于向左移动了46位。而且移动后剩下的部分即为64,而舍弃的部分实际上是小数部分0.125。
这个等式看似复杂其实很好证实，即第46位其实表明的是2^45。其除以2^46后是一个小数。依次类推…

所以对于数字 4512395720392704 ， 咱们能够用4，45，451，4512 … 依次除以 2 ^ 46. 一直到找到数451239572039270 其除以2^46不为0。这个不为0的数必定为6。
接着咱们保留后46位，实际上是保留了小数位。

假设 4512395720392704 / 2^46 = (6 + num)
64.125 =(6 + num) * 10 + C = 60 + 10* num + C

当咱们将经过位运算保留后46位，设为A, 则 A / 2^46 = num
因此 (A * 10 + C) / 2 ^46 =(num * 10 +C) = 4.125
此咱们又能够把4提取出来。实在精彩。
10进制小数位的提取是同样的，留给读者本身探索。

总结

一、本文介绍了go语言使用的IEEE-754标准存储浮点数的具体存储方式。

二、本文经过实际代码片断和一个脑筋急转弯帮助读者理解浮点数的存储方式。

三、本文介绍了normal number 以及精度这两个重要概念。

四、本文详细介绍了shopspring/decimal的实现方式，即借助了big.int，以及进制的巧妙精准转换。

五、shopspring/decimal其实在精度的巧妙转换方面参考了go源码ftoa函数的实现。读者能够参考go源码

六、shopspring/decimal目前roundShortest函数有一个bug，笔者已经提交了pr，此bug已在go源码中获得了修复。

七、big.int计算存在效率问题，若是遇到特殊的快速大量计算的场景可能不太适合。

八、还有一些decimal的实现，例如tibd/decimal,代码实在不忍淬读。

九、浮点数计算，除了要解决进制的转换外，还须要解决重要的溢出问题，例如相乘经常要超过int64的范围，这就是为何shopspring/decimal使用了big.int,而tibd/decimal将数据转换为了不少的word（int32），致使其计算很是复杂。

参考资料

1.Why 0.1 Does Not Exist In Floating-Point

2.Normal number

3.7-bits-are-not-enough-for-2-digit-accuracy

4.Decimal Precision of Binary Floating-Point Numbers

5.Introduction To Numeric Constants In Go