1993年，英特尔发布了强大的奔腾处理器，确立了其高性能处理器品牌的长期运营。

奔腾处理器包含一个浮点单元，可以快速计算正弦、余弦、对数和指数等函数。

但奔腾是如何计算这些函数的呢？

早期的英特尔芯片使用称为CORDIC的二进制算法，但奔腾切换到多项式来逼近这些超越函数，速度快得多。

这些多项式具有经过精心优化的系数，存储在芯片浮点单元内的专用ROM中。

尽管奔腾是一个拥有310万个晶体管的复杂芯片，但仍然可以通过显微镜观察这些晶体管并读取这些常数。

本文的第一部分讨论了浮点常数ROM在硬件中的实现方式。

第二部分解释了奔腾如何使用这些常数来计算正弦、对数和其他函数。

下图显示了在显微镜下观察到的奔腾的拇指大小的硅芯片。

我标记了主要的功能模块；浮点单元位于右下方。

常数ROM（突出显示）位于浮点单元的底部。

在浮点单元上方，微代码ROM保存微指令，即复杂指令的各个步骤。为了执行诸如正弦之类的指令，微代码ROM引导浮点单元执行数十个步骤，以使用来自常数ROM的常数计算近似多项式。

**在常数ROM中查找圆周率**

在二进制中，圆周率为11.00100100001111110...，但这意味着什么？

要解释这一点，二进制点左侧的11在二进制中只是3。（“二进制点”与小数点相同，只是针对二进制。）

二进制点右侧的数字的值为1/2、1/4、1/8等。

因此，二进制值`11.001001000011...对应于3 + 1/8 + 1/64 + 1/4096 + 1/8192 + ...，这与圆周率的小数值相匹配。

由于圆周率是无理数，因此位序列是无限且不重复的；ROM中的值被截断为67位并存储为浮点数。

浮点数由两部分表示：指数和尾数。

浮点数包括非常大的数字，例如6.02×1023，以及非常小的数字，例如1.055×10−34。

在十进制中，6.02×1023的尾数（或尾数）为6.02，乘以10的23次幂。

在二进制中，浮点数以类似的方式表示，具有尾数和指数，只是尾数乘以2的幂而不是10的幂。

例如，圆周率以浮点数表示为1.1001001...×21。

下图显示了圆周率在奔腾芯片中的编码方式。放大

显示常数ROM。

放大ROM的一小部分显示存储常数的晶体管行。

箭头指向表示位序列11001001的晶体管，其中0位由晶体管（垂直白线）表示，1位由没有晶体管（实心深色硅）表示。

底部每个放大的黑色矩形有两个潜在的晶体管，存储两位。

关键点是，通过查看条纹图案，我们可以确定晶体管的图案，从而确定每个常数的值，在本例中为圆周率。

**常数ROM的实现**

ROM

由MOS（金属氧化物半导体）晶体管构成，这是所有现代计算机中使用的晶体管。

下图显示了MOS晶体管的结构。

集成电路由硅衬底构成。

硅的区域掺杂有杂质以创建具有所需电性能的“扩散”区域。

晶体管可以看作是一个开关，允许电流在称为源极和漏极的两个扩散区域之间流动。

晶体管由栅极控制，栅极由一种称为多晶硅的特殊类型的硅制成。

在栅极上施加电压可使电流在源极和漏极之间流动，否则会被阻挡。

大多数计算机使用两种类型的MOS晶体管：NMOS和PMOS。这两种类型具有相似的结构，但反转了

掺杂；NMOS使用n型扩散区域，如下所示，而PMOS使用p型扩散区域。

由于这两种类型是互补的（C），

使用这两种类型的晶体管构建的电路称为CMOS。

**浮点单元的结构**

浮点单元的结构是数据垂直流过水平功能单元，如下所示。

功能单元——加法器、移位器、寄存器和比较器——按行排列。

这个具有数据流过它的功能单元集合称为数据通路。

每个功能单元由单元构成，每个单元一位，高位在左侧，低位在右侧。

每个单元具有相同的宽度——38.5 µm——因此功能单元可以像乐高积木一样卡在一起，最大限度地减少布线。

功能单元的高度根据需要而变化，具体取决于电路的复杂性。

功能单元通常有69位，但有些更宽，因此数据通路电路的边缘参差不齐。

这种基于单元的构建解释了为什么ROM每行有八个常数。

ROM位需要单个晶体管，这比加法器窄得多。因此，在每个38.5 µm单元中放置一位会浪费大部分空间。

将ROM位压缩成一个狭窄的块效率也很低，需要斜线连接每个ROM位到相应的数据通路位。

通过在每个单元中放置八个常数的八位，ROM单元的宽度与数据通路的其余部分匹配，并且位对齐

得到保留。

因此，ROM在硅中的布局密集、高效，并且与浮点单元其余部分的宽度相匹配。

**多项式逼近：不要使用泰勒级数**

现在我将从硬件转向常数。

如果您查看附录中的常数ROM内容，您可能会注意到许多常数接近倒数或倒数阶乘，但并不完全

匹配。例如，一个常数是0.1111111089，接近1/9，但明显错误。

另一个常数几乎是1/13！（阶乘），但误差为0.1%。这是怎么回事？

奔腾使用多项式来逼近超越函数（正弦、余弦、正切、反正切以及以2为底的幂和对数）。

英特尔的早期浮点单元，从8087到486，使用了一种称为CORDIC的算法，该算法一次生成一个结果

一位。

但是，奔腾利用其快速的乘法器和更大的ROM，并使用多项式，

计算结果的速度是486算法的两到三倍。

您可能还记得微积分中，泰勒级数多项式逼近一个点（通常为0）附近的函数。

例如，下面的公式给出了正弦的泰勒级数。

使用上面显示的五个项生成一个在下面的图形中看起来与正弦无法区分的函数。

但是，事实证明，这种近似值的误差太大而无法使用。

问题在于泰勒级数在0附近非常准确，但在参数范围的边缘误差急剧上升，如下面的左侧图形所示。

在实现函数时，我们希望函数在任何地方都准确，而不仅仅是在0附近，因此泰勒

级数不够好。

一种改进称为范围缩减：将参数缩减到更小的范围，以便您处于准确的平坦部分。

右图查看了较小范围[-1/32, 1/32]上的泰勒级数。

这显著降低了误差，大约降低了22个数量级（注意比例尺的变化）。

但是，误差仍然以完全相同的方式在范围的边缘上升。

无论你多少

缩减范围，中间几乎没有误差，但边缘有很多误差。

我们如何消除边缘附近的误差？诀窍是以一种特殊的方式调整泰勒级数的系数，这将

增加中间的误差，但更大幅度地减少边缘的误差。

由于我们希望最大限度地减少整个范围内的最大误差（称为极小化极大），因此这种权衡是有益的。

具体来说，可以通过称为Remez算法的过程优化系数。

如下所示，将系数更改不到1%即可显着提高精度。

优化后的函数（蓝色）在整个范围内具有更低的误差，因此比泰勒级数（橙色）的近似值要好得多。

总之，泰勒级数在微积分中很有用，但不应用于逼近函数。通过使用Remez算法对系数进行非常小的修改，您可以获得

更好的逼近。这解释了

为什么ROM中的系数几乎但不完全匹配泰勒级数。

**反正切**

我现在将研究奔腾用于不同超越函数的常数。

常数ROM包含两个反正切多项式的系数，一个用于单精度，另一个

用于双精度。

这些多项式几乎与泰勒级数匹配，但经过修改以提高精度。

ROM还保存了反正切(1/32)到反正切(32/32)的值；

范围缩减过程使用这些常数与三角恒等式一起将参数范围缩减到

[-1/64, 1/64]。

您可以在附录中看到反正切常数。

下图显示了奔腾的反正切多项式（蓝色）与相同长度的泰勒级数（橙色）的误差。

**三角函数**

正弦和余弦各有两种多项式实现，一种在ROM中包含4项，另一种在ROM中包含6项。

（请注意，系数1未存储在ROM中。）

常数表还保存了16个常数，例如sin(36/64)和cos(18/64)，用于参数范围缩减。

奔腾通过将正弦除以余弦来计算正切。

我没有显示图形，因为奔腾的误差结果比泰勒级数差，所以要么我的系数有误，要么我的操作方式错误。

**指数**

奔腾有一个指令来计算2的幂。

ROM中有两组指数多项式系数，一组包含6项，另一组包含11项。

奇怪的是，ROM中的多项式计算ex，而不是2x。

因此，奔腾必须将参数按ln(2)缩放，这是一个存储在ROM中的常数。

下图显示了奔腾的多项式相对于泰勒级数多项式的优势。

多项式处理狭窄的参数范围[-1/128, 1/128]。

观察到，当以二进制计算2的幂时，对参数的整数部分求幂是微不足道的，因为它成为结果的指数。

因此，该函数只需要处理范围[1, 2]。

对于范围缩减，常数ROM保存64个形式为2n/128-1的常数值。

要将范围从[1, 2]缩减到[-1/128, 1/128]，将最接近的n/128从参数中减去，然后将结果乘以ROM中相应的常数。

常数间隔不规则，大概是为了提高精度；有些以4/128为步长，另一些以2/128为步长。

**对数**

奔腾可以计算以2为底的对数。

常数ROM有9个系数，大概用于对数多项式（或多项式），但我无法

从中形成有用的多项式。

与其他多项式不同，此多项式与相应的泰勒级数不相似。

ROM还包含用于范围缩减的64个常数：log2(1+n/64)，对于从1到63的奇数n。

这些常数的独特之处在于每个常数被分成两部分以提高精度：

顶部部分具有40位的精度，底部部分具有67位的精度，总共提供107位的常数。

需要额外的位，因为对数难以精确计算。

**其他常数**

X87浮点指令集可以直接访问少数几个常数——0、1、圆周率、

log2(10)、log2(e)、log10(2)和loge(2)——因此这些常数

存储在ROM中。

（这些对数用于更改对数和指数的底数。）

ROM保存其他用于浮点单元内部使用的常数，例如-1、2、7/8、9/8、圆周率/2、圆周率/4和2log2(e)。

ROM还保存用于提取字的一部分的位掩码，例如访问字中的4位BCD数字。

虽然我可以解释大多数值，但有一些谜团，例如具有难以理解的值的掩码

0x5c3bd5191b525a249。

ROM在末尾有34个未使用的条目；这些条目保存包含描述性十六进制值0xbad的字。

**我如何检查ROM**

为了分析奔腾，我使用各种化学物质（硫酸、磷酸、Whink）去除了金属和氧化物层。

（我后来发现，简单地打磨芯片效果出奇地好。）

接下来，我用显微镜拍摄了ROM的许多照片。

这个奔腾的特征尺寸是800 nm，略大于可见光（380-700 nm）。

因此，可以在光学显微镜下检查芯片，但它已经接近极限。

为了确定ROM内容，我辛辛苦苦地浏览了ROM图像，检查了26144位中的每一个位并标记了每个晶体管。

在弄清楚ROM格式后，

我编写了程序，将简单函数以多种不同的组合组合起来，以

确定数学表达式，例如反正切(19/32)或log2(10)。

由于多项式常数是优化的，并且我的ROM数据存在位错误，因此我的程序需要

检查不精确的匹配，包括数值和按位匹配。

最后，我必须确定如何将这些常数用于算法。

**结论**

通过在显微镜下检查奔腾的浮点ROM，可以提取存储在ROM中的304个常数。

我能够确定这些常数的大部分含义并推断出奔腾使用的一些浮点算法。

这些常数说明了多项式如何有效地计算超越函数。

虽然泰勒级数多项式众所周知，但它们令人惊讶地不准确，应避免使用。

但是，通过Remez算法对系数进行微小的更改会产生更好的多项式。

在之前的一篇文章中，我检查了存储在8087协处理器中的浮点常数。

奔腾在奔腾中拥有304个常数，而8087中只有42个，支持更有效的算法。

此外，8087是一个外部浮点单元，而奔腾的浮点单元是处理器的一部分。

8087（1980年，65,000个晶体管）和奔腾（1993年，310万个晶体管）之间的变化是由于

晶体管数量呈指数级增长，如摩尔定律所述。

我计划撰写更多关于奔腾的文章，因此请在Bluesky (@righto.com)或

RSS 上关注我以获取更新。（我不再使用Twitter了。）

我还写过关于奔腾除法错误和奔腾纳瓦霍地毯的文章。

感谢CuriousMarc提供的显微镜帮助。