在 X87 和 SSE FPU 中用户定义的点之后截断浮点数和双精度数

Question

我已经创建了一个函数g，它能够在一定程度上逼近一个函数，这个函数给出了精确到小数点后 5 位的结果（1,23456xxxxxxxxxxxx，其中 x 位置只是舍入错误/垃圾）。

为避免将错误传播到将使用g 结果的其他计算，我想将所有 x 位置设置为零，更好的是，只需将小数点后 5 位的所有内容设置为 0。

我在 X87 和 SSE 文献中没有找到任何可以让我以我想要的方式使用 IEEE 754 位或其表示的内容。

对 X87 的 FISTP 指令有一个旧的引用，它显然在 SSE 世界中与 FISTTP 相映成趣，其好处是 FISTTP 不需要修改控制字，因此速度更快。

我注意到FISTTP 被称为“切碎模式”，但现在在更现代的文献中只是“向零舍入”或“截断”，这让我感到困惑，因为“切碎”意味着完全删除“向零舍入”的东西零”对我来说不一定意味着同样的事情。

我不需要四舍五入到零，我只需要保留最多 5 位小数作为函数的最后一步，然后将结果存储到内存中；如何在 X87（标量 FPU）和 SSE（矢量 FPU）中执行此操作？

Answer 1

正如几个人评论的那样，more early rounding doesn't help the final result be more accurate。如果您想阅读更多关于浮点比较和怪异/陷阱的信息，我强烈推荐 Bruce Dawson 的浮点系列文章。这是the one with the index的引述

我们终于在这个系列中达到了我一直在等待的点 为了。在this post我将分享最关键的部分 我拥有的浮点数学知识。这里是：

[浮点] 数学很难。

你简直不会相信它是多么巨大、巨大、令人难以置信的困难。 我的意思是，你可能认为很难计算火车从 芝加哥和洛杉矶会发生碰撞，但这只是小菜一碟 浮点数学。

（如果您认识到最后一段是对有关空间的著名行的释义，则可以加分。）

---

如何实际实施你的坏主意：

没有任何机器指令或 C 标准库函数可以截断或舍入到整数以外的任何值。

请注意，有些机器指令（和 C 函数）会将 double 舍入为最接近的（可表示的）整数，而不会将其转换为 intmax_t 或任何东西，只需 double->double。所以不会通过固定宽度的 2 的补码整数往返。

所以要使用它们，您可以将浮点数放大某个因子，四舍五入到最接近的整数，然后再缩小。 （就像 chux 的基于 round() 的函数，但我建议使用 C99 double rint(double) 而不是 round()。round 具有奇怪的舍入语义，与 x86 上的任何可用舍入模式都不匹配，因此它编译得更糟代码。

---

你一直提到的 x86 asm 指令没什么特别的，不要做任何你不能要求编译器用纯 C 做的事情。

FISTP (Float Integer STore (and Pop the x87 stack) 是编译器或 asm 程序员实现long lrint(double) or (int)nearbyint(double). 的一种方法@ 一些编译器为其中一个或另一个生成更好的代码。它使用当前的 x87 舍入模式进行舍入（默认：四舍五入），即same semantics as those ISO C standard functions。

FISTTP (Float Integer STore with Truncation (and Pop the x87 stack) is part of SSE3, 即使它在 x87 堆栈上运行。它让编译器避免将舍入模式设置为截断（向零舍入）实现(long)x的C截断语义，然后恢复旧的舍入模式。

这就是“不修改控制字”的内容。两条指令都没有这样做，但要在没有 FISTTP 的情况下实现 (int)x，编译器必须使用 other 指令来修改和恢复围绕 FIST 指令的舍入模式。或者只是使用SSE2 CVTTSD2SI 将 xmm 寄存器中的双精度值转换为截断，而不是传统 x87 堆栈上的 FP 值。

由于 FISTTP 仅适用于 SSE3，因此您只能将其用于 long double，或者用于在 x87 寄存器中具有 FP 值的 32 位代码，因为旧的 32 位 ABI 返回 FP x87 堆栈上的值。

---

PS。如果你不认识 Bruce 的 HHGTG 参考，the original is：

空间很大。真的很大。你不会相信有多么巨大 大得令人难以置信。我的意思是你可能认为这是一个很长的路要走 通往化学家之路，但那只是通往太空的路。

Answer 2

如何在 X87（标量 FPU）和 SSE（向量 FPU）中执行此操作？

以下不使用 X87 和 SSE。我已将其作为通用代码的社区参考。如果有的话，它可以用来测试 X87 解决方案。

1. g() 的结果的任何“斩波”至少会略微增加错误，希望可以容忍，因为 OP 说“为了避免将错误传播到其他计算......”

2. 不清楚 OP 是否希望“精确到小数点后 5 位的结果”来反映绝对精度 (+/- 0.000005) 或相对精度 (+/- 0.000005 * 结果)。将假定“绝对精度”。

3. 由于 float、double 通常是 二进制 浮点数，因此任何“印章”都会将 FP 数 最接近 的倍数0.00001.

4. 文字方法：

   //       - x    xxx...xxx    . xxxxx \0
   char buf[1+1+ DBL_MAX_10_EXP+1  +5   +1];
   sprintf(buf, "%.5f", x);
   x = atof(buf);
   

5. round() rint() 方法：

   #define SCALE 100000.0
   if (fabs(x) < DBL_MAX/SCALE) {
     x = x*SCALE;
     x = rint(x)/SCALE;
   }
   

6. x 的直接位操作。只需将有效数字中的选择位置零即可。

   TBD code.