每个计算机科学家都应该知道的浮点算术(1991)

注-本附录是David Goldberg的论文“每个计算机科学家都应该知道的浮点算术”的编辑再版，该论文发表在1991年3月的“计算调查”上。版权所有1991，计算机械协会，Inc.，经许可转载。

浮点算术被许多人认为是一门深奥的学科。这相当令人惊讶，因为浮点在计算机系统中无处不在。几乎每种语言都有浮点数据类型；从PC到超级计算机的计算机都有浮点加速器；大多数编译器将不时被要求编译浮点算法；几乎每个操作系统都必须对浮点异常(如溢出)做出响应。这篇文章介绍了浮点技术对计算机系统设计者有直接影响的那些方面的教程。它从浮点表示和舍入误差的背景开始，接着讨论IEEE浮点标准，并以计算机构建者如何更好地支持浮点的大量示例结束。

类别和主题描述符：(主要)c.0[计算机系统组织]：一般--指令集设计；D.3.4[编程语言]：处理器--编译器，优化；G.1.0[数值分析]：一般--计算机算术、误差分析、数值算法(次要)。

D.2.1[软件工程]：需求/规范--语言；D.3.4编程语言]：形式定义和理论--语义；D.4.1操作系统]：进程管理--同步。

附加关键词和短语：非规范化数、异常、浮点、浮点标准、渐进式下溢、保护位、NaN、上溢、相对误差、舍入误差、舍入方式、ULP、下溢。

计算机系统的构建者通常需要有关浮点算术的信息。然而，关于它的详细信息来源非常少。Pat Sterbenz的“浮点计算”是为数不多的关于这一主题的书籍之一，它已经绝版很久了。本文是关于浮点算术(以下简称浮点)与系统构建有直接联系的那些方面的教程。它由三个松散连接的部分组成。第一部分，舍入误差，讨论了使用不同的舍入策略对基本的加、减、乘、除运算的影响。它还包含测量舍入误差(ULPS)和相对误差的两种方法的背景信息。第二部分讨论了正在迅速被商用硬件制造商接受的IEEE浮点标准。IEEE标准中包括基本运算的舍入方法。标准的讨论借鉴了截面舍入误差中的材料。第三部分论述了浮点运算与计算机系统各个方面的设计之间的联系。主题包括指令集设计、优化编译器和异常处理。

我试图避免发表关于浮点数的声明，也没有给出声明为真的原因，特别是因为理由只涉及初等微积分。那些不是主要论点核心的解释已被归入一个名为“细节”的部分，以便在需要时可以跳过。特别地，许多定理的证明都出现在这一节中。每个校样的末尾都标有z符号。如果没有包含证明，则z紧跟在定理语句之后。

将无穷多个实数压缩成有限个比特需要近似表示法。虽然有无限多的整数，但在大多数程序中，整数计算的结果可以存储在32位中。相反，在给定任何固定位数的情况下，大多数使用实数的计算都会产生无法使用那么多位精确表示的量。因此，浮点计算的结果通常必须进行四舍五入，以便重新适合其有限表示。这种舍入误差是浮点计算的特征。相对误差和ULPS一节描述了如何测量它。

由于大多数浮点计算都有舍入误差，如果基本算术运算引入的舍入误差比必要的多一点，这有什么关系吗？这个问题是本节的主题。保护数字一节讨论了保护数字，这是一种在减去两个相邻数字时减少误差的方法。IBM认为保护数字非常重要，以至于在1968年，它在System/360体系结构的双精度格式中增加了一个保护数字(单精度已经有了一个保护数字)，并对该领域现有的所有机器进行了改造。给出了两个例子来说明保护数字的用途。

IEEE标准不仅仅要求使用保护数字。它给出了加、减、乘、除和平方根的算法，并要求实现产生与该算法相同的结果。因此，当程序从一台机器转移到另一台机器时，如果两台机器都支持IEEE标准，则基本操作的结果在每一位上都是相同的。这大大简化了程序的移植。此精确规范的其他用法在精确四舍五入运算中给出。

已经提出了几种不同的实数表示法，但到目前为止最广泛使用的是浮点表示法。1浮点表示法具有基数(始终假定为偶数)和精度p。如果=10且p=3，则数字0.1表示为1.00×10-1。如果=2且p=0.24，则不能精确表示十进制数0.1，而是近似为1.1001100110011001101×2-4。

通常，浮点数将表示为±d.dd.。d×e，其中d.dd.。D称为有效位2，具有p位。更准确地说，±d0。D%1%D%2..。d p-1×e代表数字。

(一)。术语浮点数将用于表示能够以所讨论的格式精确表示的实数。与浮点表示相关联的另外两个参数是允许的最大和最小指数，即emax和鄂敏。由于有p个可能的有效位，并且emax-鄂敏+1个可能的指数，所以浮点数可以编码为。

位，其中最后的+1表示符号位。目前，精确的编码并不重要。

实数可能无法精确表示为浮点数的原因有两个。最常见的情况由十进制数字0.1说明。虽然它有一个有限的十进制表示法，但在二进制中它有一个无限重复的表示法。因此，当=2时，数字0.1严格位于两个浮点数之间，并且两者都不能精确表示。一种不太常见的情况是实数超出范围，即其绝对值大于×或小于1.0×。由于第一个原因，本文大部分时间都在讨论问题。但是，超出范围的数字将在“无穷大”和“非正规化数字”一节中讨论。

浮点表示形式不一定是唯一的。例如，0.01×10-1和1.00×10-1都代表0.1.。如果前导数字不是零(上面等式(1)中的d0)，则表示被认为是规范化的。浮点数1.00×10-1被归一化，而0.01×101不归一化。当t=2，p=3时，鄂敏x=-1，e max=2时，归一化浮点数有16个，如图D-1所示。粗体散列标记对应于有效位数为1.00的数字。要求对浮点表示进行标准化使得该表示是唯一的。不幸的是，此限制使得不可能表示零！表示0的一种自然方式是使用1.0×，因为这保留了这样一个事实，即非负实数的数字顺序对应于它们的浮点表示的字典顺序。3当指数存储在k位字段中时，这意味着只有2个k-1值可用作指数，因为必须保留一个值来表示0。

请注意，浮点数中的×是表示法的一部分，不同于浮点乘法运算。×符号的意思应该从上下文中看得很清楚。例如，表达式(2.5×10-3)×(4.0×102)仅涉及单个浮点乘法。

图D-1当=2，p=3，鄂敏=-1，emax=2时的归一化数由于四舍五入误差是浮点计算的固有误差，因此有一种测量该误差的方法非常重要。考虑一下浮点格式，其中的参数=10和p=3，这将在本节中使用。如果浮点计算的结果是3.12×10-2，而当计算到无限精度时的答案是.0314，则很明显这最后一个位置的误差是2个单位。类似地，如果实数.0314159表示为3.14×10-2，则它的误差为最后一位的0.159个单位。一般来说，如果浮点数D.D.。用d×e来表示z，则以d×e来表示是错误的。最后一位是d-(z/e)p-1个单位。4、5术语ULP将用作最后一位单位的速记。如果计算结果是最接近正确结果的浮点数，则仍可能误差高达0.5个ULP。另一种测量浮点数和它所逼近的实数之间的差异的方法是相对误差，即两个数字除以实数的差值。例如，将3.14159近似为3.14×100时的相对误差为0.00159/3.14159和0.0005。

要计算对应于0.5ULP的相对误差，请观察到，当实数由可能最接近的浮点数d.dd近似时…。dd×e，误差最大可达0.00.00×e，其中#2为数字/2，浮点数的有效位有p个单位，误差的有效位有p个单位为0。此错误为((/2)-p)×e。由于d.dd格式的数字.。Dd×e都具有相同的绝对误差，但是具有范围在e和×e之间的值，相对误差范围在((/2)-p)×e/e和((/2)-p)×e/e+1之间。即，

(2)特别地，对应于0.5ULP的相对误差可以按因数变化。这个因素称为摆动。将=(/2)-p设置为上面(2)中的最大界限，我们可以说，当实数四舍五入到最接近的浮点数时，相对误差始终以e为界，这称为机器ε。

在上例中，相对误差为.00159/3.14159.0005。为了避免这样的小数字，通常将相对误差写成因子倍，在这种情况下为=(/2)-p=5(10)-3=0.005。因此，相对误差将表示为(.00159/3.14159)/.005)0.1。

为了说明ULP和相对误差之间的差异，请考虑实数x=12.35。其近似表达式为=1.24×101，误差为0.5ULPS，相对误差为0。8.。接下来，考虑计算8。精确值为8x=98.8，计算值为8=9.92×101，误差为4.0ULPS，但相对误差仍为0。8.。即使相对误差相同，以ULPS测量的误差也要大8倍。通常，当基数为时，以ULPS表示的固定相对误差可以抖动高达1倍的抖动。相反，如上面的公式(2)所示，0.5ULPS的固定误差导致可抖动的相对误差。

测量舍入误差最自然的方法是使用ULPS。例如，四舍五入到最接近的浮点数对应于小于或等于0.5ulp的误差。然而，在分析各种公式引起的舍入误差时，相对误差是较好的衡量标准。定理9中的分析就是一个很好的例证。由于可以高估抖动因子舍入到最接近的浮点数的影响，因此公式的误差估计在具有较小误差的机器上会更紧密。

当只关注舍入误差的数量级时，ULP和可以互换使用，因为它们最多相差1倍。例如，当一个浮点数的误差为n个ULP时，意味着受污染的位数是logn。如果计算中的相对误差是n，那么

(3)污染数字logn。计算两个浮点数之差的一种方法是精确计算两个浮点数之间的差，然后四舍五入到最接近的浮点数。如果操作数的大小差别很大，这是非常昂贵的。假设p=3，2.15×1012-1.25×10-5将计算为。

x=2.15×1012y=0.000000000000000012x1012x-y=2.1499999999999999875×1012。

它舍入为2.15×1012。浮点硬件通常在固定位数上操作，而不是使用所有这些位。假设保留的位数是p，并且当较小的操作数向右移位时，数字将被简单地丢弃(与舍入相反)。那么2.15亿×1012亿-1.25亿×1010-5就变成了。

x=2.15×10 12 y=0.00×10 12 x-y=2.15×10 12。

答案完全相同，就好像差额经过精确计算，然后四舍五入一样。再举一个例子：10.1-9.93。这就变成了。

x=1.01×101 y=0.99×101 x-y=0.02×101。

正确答案是0.17，所以计算出的差值相差30ULP，而且每个数字都是错误的！错误会有多严重呢？

使用带参数和p的浮点格式，并使用p位计算差值，结果的相对误差可达-1。

当x=1.00.0且y=0时，表达式x-y中的相对误差为-1。.，其中=-1。这里y有p位(全部等于)。精确的差值是x-y=-p。但是，当仅使用p位计算答案时，y的最右边的数字被移位，因此计算的差值为-p+1。因此，误差为-p=-p+1=-p(-1)，相对误差为-p(-1)/-p=-1。z。

当=2时，相对误差可以和结果一样大，当=10时，相对误差可以是结果的9倍。或者换句话说，当=2时，公式(3)显示污染位数为log2(1/)=log2(2p)=p，也就是说，结果中的所有p位都是错误的！假设增加了一个额外的数字来防止这种情况(保护数字)。也就是说，将较小的数字截断为p+1位，然后将减法结果舍入为p位。使用保护数字，上一个示例变为。

x=1.010×101 y=0.993×101 x-y=0.017×101。

答案是准确的。对于单个保护数字，结果的相对误差可能大于，如110-8.59。

x=1.10×102 y=0.085×102 x-y=1.015×102。

与正确答案101.41相比，取整为102时的相对误差为0.006，大于=0.005。一般来说，结果的相对误差只能略大于。更准确地说，

如果x和y是具有参数和p的格式的浮点数，并且如果用p+1位(即一个保护位)进行减法，则结果中的相对舍入误差小于2。

这个定理将在舍入误差中得到证明。由于x和y可以是正的，也可以是负的，所以上述定理中包含了加法。

最后一节可以总结说，如果没有保护数字，当减去两个相邻的量时，所犯的相对误差可能非常大。换句话说，任何包含减法(或加法符号相反的量)的表达式的求值都可能导致较大的相对误差，以至于所有的数字都是没有意义的(定理1)。当减去附近的数量时，操作数中的最高有效位相互匹配并相互抵消。取消有两种：灾难性的和良性的。

当操作数受到舍入误差的影响时，会发生灾难性的取消。例如，在二次公式中，出现表达式b2-4ac。量B2和4Ac受到舍入误差的影响，因为它们是浮点乘法的结果。假设它们四舍五入为最接近的浮点数，因此精度在0.5ULP以内。当它们被减去时，抵消可能会导致许多准确的数字消失，留下的主要是被舍入误差污染的数字。因此，差值可能会有许多ULP的误差。例如，考虑b=3.34，a=1.22，c=2.28。b2+-b4 ac的精确值为0.0292.。但是b2舍入到11.2，4ac舍入到11.1，因此最终的答案是.1，这是70ULP的误差，即使11.2-11.1正好等于1.6.减法没有引入任何误差，而是暴露了较早的乘法中引入的误差。

良性抵消发生在精确减去已知量时。如果x和y没有舍入误差，那么根据定理2，如果减法是用保护位完成的，则x-y的差具有非常小的相对误差(小于2)。

有时可以重新排列显示灾难性取消的公式来消除问题。再次考虑二次公式。

(四)。但其中一个公式中的另一个加法(减法)将产生灾难性的取消。要避免这种情况，请将r 1的分子和分母乘以。

(5)如果与，则使用公式(4)计算r1将涉及取消。因此，使用公式(5)计算r1，使用公式(4)计算r2。另一方面，如果b<；0，则使用公式(4)计算r1，使用(5)计算r2。

表达式x2-y2是表现出灾难性抵消的另一个公式。将其计算为(x-y)(x+y)更为准确。7与二次公式不同的是，这种改进的形式仍然有减法，但它是没有舍入误差的良性数量抵消，而不是灾难性的。根据定理2，x≤y的相对误差至多为2。x+y也是如此。将两个量乘以较小的相对误差得到的乘积具有较小的相对误差(请参阅舍入误差一节)。

为了避免将精确值和计算值混淆，使用以下表示法。x-y表示x和y的精确差值，xy表示计算的差值(即，具有舍入误差)。类似地，和分别表示计算的加法、乘法和除法。所有大写表示函数的计算值，如ln(X)或sqrt(X)。小写函数和传统的数学符号表示它们的精确值，如ln(X)和。

虽然(x，y)(x，y)是x2-y 2的极佳近似值，但浮点数x和y本身可能是近似值。

..