IEEE-754浮点精度:允许多少误差？

Question

我正在努力将sqrt函数(64位双倍)从fdlibm移植到我目前使用的模型检查工具(cbmc)。

作为我工作的一部分，我读了很多关于ieee-754标准的文章，但我想我不理解基本操作的精确性保证(包括。(平方米)

通过测试我的fdlibm的sqrt端口，我得到了在64位双倍上使用sqrt的计算结果：

sqrt(1977061516825203605555216616167125005658976571589721139027150498657494589171970335387417823661417383745964289845929120708819092392090053015474001800648403714048.0) = 44464159913633855548904943164666890000299422761159637702558734139742800916250624.0

(在我关于精度的测试中，这个案例打破了一个简单的后条件；我不再确定这种后状态是否可以在IEEE-754中实现)

作为比较，几种多精度工具计算的结果如下：

sqrt(1977061516825203605555216616167125005658976571589721139027150498657494589171970335387417823661417383745964289845929120708819092392090053015474001800648403714048.0) =44464159913633852501611468455197640079591886932526256694498106717014555047373210.truncated

可以看出，左边的第17位数字是不同的，这意味着一个错误，如：

3047293474709469249920707535828633381008060627422728245868877413.0

问题1:允许这么大数量的错误吗？

标准是说每个基本操作(+，-，*，/，sqrt)都应该在0.5ULP以内，这意味着它应该等于一个数学上精确的结果四舍五入到最接近的fp表示(wiki是说有些库只保证一个ulp，但这目前并不那么重要)。

问题2:这是否意味着，每个基本操作都应该有一个错误< 2.220446e-16和64位双(机器-epsilon)？。

我在x86-32 linux系统(glibc / eglibc)中进行了同样的计算，得到了与fdlibm相同的结果，这让我认为：

• 答:我做错了什么(但如何：printf会成为候选人，但我不知道这是否是原因)
• b:错误/精度在这些库中很常见

Answer 1

IEEE-754标准要求所谓的“基本操作”(包括加法、乘法、除法和平方根)是正确四舍五入的。这意味着有一个唯一允许的答案，而且它是与操作的所谓“无限精确”结果最接近的可表示的浮点数字。

在双精度中，数字有53个二进制数字的精度，所以正确的答案是精确的答案四舍五入为53个重要数字。正如 Regan在他的回答中所显示的，这正是你得到的结果。

你的问题的答案是：

问题1:允许这么大数量的错误吗？

是的，但是把这个错误称为“巨大的”是很误导的。事实上，没有一个可以返回的双精度值会有较小的错误。

问题2:这是否意味着，每个基本操作都应该有一个错误< 2.220446e-16和64位双(机器-epsilon)？。

不是的。这意味着根据当前的舍入模式，每个基本操作都应该四舍五入到(唯一的)最近可表示的浮点数。这与说相对误差被机器epsilon限制是不完全一样的。

问题3:您的x86硬件和gcc + libc?的结果是什么？

同样的答案也是如此，因为sqrt在任何合理的平台上都是正确的。

Answer 2

在二进制中，任意精确答案的前58位是1011111111111111111111110101010101111111111111111011010001...。

这个双值的53位意义是

10111111111111111111111101010101011111111111111110111

这意味着双值被正确地舍入到53位有效位，并且在1/2的ULP内。(错误“大”只是因为数字本身很大)。