浮点数除法的软件实现,舍入问题

9
作为一个学习项目,我正在使用c++在软件中实现浮点数操作(加、减、乘、除)。目标是更加熟悉浮点数行为的基本细节。
我尝试将处理器操作与精确的比特匹配,即IEEE 754标准。到目前为止,它运行得很好,加、减和乘法表现完美,我在大约1.1亿个随机操作中进行了测试,并得到了与处理器在硬件上执行的完全相同的结果。(虽然没有考虑边缘情况、溢出等)。
之后,我开始转向最后一个操作,即除法。它可以正常工作并达到预期结果,但是有时候我会发现最后一位尾数错误,未对其进行四舍五入。我有点难以理解为什么。
我主要参考的是John Farrier 的优秀演讲(时间戳显示如何进行四舍五入)。

https://youtu.be/k12BJGSc2Nc?t=1153

那个四舍五入对于所有运算都很有效,但对于除法却给我带来了麻烦。 让我举个具体的例子。 我正试图将 645.68011474609375 除以 493.20962524414063。
我得到的最终结果是:
我的:0-01111111-01001111001000111100000
C++:0-01111111-01001111001000111100001
可以看到除了最后一位之外,其他都匹配。 我计算除法的方式基于这个视频: https://www.youtube.com/watch?v=fi8A4zz1d-s

接下来,我计算了28位精度为24的尾数(隐含1 + 23个尾数)和3位保护位、四舍五入位再加上一个可能的移位的额外位。 使用视频中的算法,我最多可以获得1个归一化移位,这就是为什么我在结尾处有一个额外的位,以防它在归一化中被移动,因此在舍入中可用。现在,以下是我从除法算法中得到的结果:

 010100111100100011110000 0100
 ------------------------ ----
 ^                        grs^
 |__ to be normalized        |____ extra bit

正如您所看到的,我在第24个位置得到了0,因此我需要向左移动一个位置以获得正确的归一化。这意味着我将得到:

10100111100100011110000 100

根据John Farrier的视频,在100 grs位的情况下,如果尾数的最低有效位是1,我只会进行规格化。在我的情况下,它是0,这就是为什么我不会将结果四舍五入的原因。
我有点迷茫的原因是我确信我的算法正在计算正确的尾数,我已经用在线计算器仔细检查过了,舍入策略对于所有其他操作都有效。此外,以这种方式计算会触发规范化,最终得到正确的指数。
我错过了什么吗?哪里有小细节?
有一件事让我感到奇怪的是粘性位,在加法和乘法中,您会获得不同程度的移位,这导致更高的粘性位触发几率,在这种情况下,我最多只能移动一个位数,这使得粘性位不是真正的粘性。
我希望我提供足够的细节来理解我的问题。在此处,您可以在底部找到我的除法实现,其中填充了我用于调试的打印内容,但应该可以了解我正在做什么,代码从第374行开始:

https://gist.github.com/giordi91/1388504fadcf94b3f6f42103dfd1f938

PS:同时,我正在学习“科学家应该了解的浮点数知识”,以便查看是否有遗漏的内容。


2
“应该能给你个我在做什么的概念。”- 我们强烈建议将确切的代码整理后放在问题本身。现在我能明白为什么你不愿意把它添加到问题中了,因为问题已经非常长了。但那是因为它用词太多了。我们相信你可以实现加法,无需告诉我们。 - MSalters
没问题,我从代码片段中删除了与除法无直接关系的代码。 - Marco Giordano
@RichardCritten 感谢您的回复,我知道那篇博客文章,正在慢慢地阅读。在我的特定情况下,我已经计算了额外的4位尾数,在grs舍入方案中,它是唯一需要的位数。因此,即使我将尾数的额外位四舍五入,也不会发生进位,从而导致粘性位翻转并进行舍入。这就是我目前遇到的问题,额外的精度位对我没有帮助。同时,我将尝试从那篇博客中挖掘有用的信息。如果感兴趣,我可以打印出更多尾数的精度位。 - Marco Giordano
1
也许您想在进行定点除法之前在分子上添加/移动几个比特,从而得到一个大约4倍的结果,然后您可以在将它们移回去之前进行检查/舍入(x/y = (4x/4y) = (4x/y)/4)。可以保留许多位,如5或8或16位... - old_timer
2
你为什么要使用逐位相除?你可以使用整数除法,这样你的代码会更快。并使用余数来决定舍入方向。 - geza
显示剩余6条评论
1个回答

7
从除法算法得到的结果是不充分的。你展示了如下内容:
 010100111100100011110000 0100
 ------------------------ ----
 ^                        grs^
 |__ to be normalized        |____ extra bit

数学上精确的商继续如下:
 010100111100100011110000 0100 110000111100100100011110

因此,在你进行四舍五入的位置处的余数超过了½ ULP,所以应该向上舍入。我没有详细研究你的代码,但看起来你可能只是计算了一个或两个有效数字的额外位数。实际上,你需要知道余数是否非零,而不仅仅是它的下一位或两位是否为零。如果精确的数学结果中在该位置或之后的任何位都是非零的,则最终的粘性位应为1。

脚注

1“Significand”是首选术语。“Mantissa”是对于对数的分数部分的传统术语。浮点值的有效数字是线性的。幂指数是对数的。


感谢您的回复,也感谢您在mantissa与significand方面的更正。我不知道其“遗留”的含义。您提出的解决方案让在grs之后添加“任何”额外位来影响粘滞,这个想法非常有可行性!我会尝试一下并告知您结果。我想知道还需要多少额外的精度位?有人知道现代CPU在硬件方面是怎样处理的吗? - Marco Giordano
非常感谢您的帮助,先生。看起来您的解决方案奏效了。现在它是脆弱的分割,所有的测试都通过了,已经完成了8600万个,而且还在继续。非常感谢! - Marco Giordano
1
@MarcoGiordano,您不应使用任何有限数量的额外位。 "在精确的数学结果中,该位置或之后的任何位都将是非零"与“当前余数为非零”完全等效。 - Patricia Shanahan
@PatriciaShanahan 很有趣,所以您是告诉我如果在计算grs位后我收到一个提醒,我只需设置粘性位? - Marco Giordano

网页内容由stack overflow 提供, 点击上面的
可以查看英文原文,
原文链接