与-O3相比，gcc -Ofast的汇编代码中计算不精确的根源在哪里？ [关闭]答案

【问题标题】：Where is the source of imprecise calculation in the assembler code of gcc -Ofast compared with -O3? [closed]与-O3相比，gcc -Ofast的汇编代码中计算不精确的根源在哪里？ [关闭]
【发布时间】：2021-07-31 16:46:54
【问题描述】：

以下 3 行使用 "gcc -Ofast -march=skylake" 给出了不精确的结果：

int32_t  i = -5;
const double  sqr_N_min_1 = (double)i * i;
1. - ((double)i * i) / sqr_N_min_1

显然，sqr_N_min_1 得到 25.，并且在第三行中，(-5 * -5) / 25 应该变为 1.，因此第三行的总体结果正好是 0.。事实上，编译器选项 "gcc -O3 -march=skylake" 也是如此。

但是使用 "-Ofast" 最后一行产生 -2.081668e-17 而不是 0. 和其他 i 而不是 -5（例如 6 或 7）它得到与0. 的其他非常小的正或负随机偏差。我的问题是：这种不精确的根源在哪里？

为了调查这个问题，我用 C 编写了一个小测试程序：

#include <stdint.h>      /* int32_t */
#include <stdio.h>
#define MAX_SIZE 10

double W[MAX_SIZE];

int main( int argc, char *argv[] )
{
  volatile int32_t n = 6; /* try 6 7 or argv[1][0]-'0' */
  double           *w = W;
  int32_t          i = 1 - n;
  const int32_t    end = n - 1;
  const double     sqr_N_min_1 = (double)i * i;

  /* Here is the crucial part. The loop avoids the compiler replacing it with constants: */
  do {
    *w++ = 1. - ((double)i * i) / sqr_N_min_1;
  } while ( (i+=2) <= end );

  /* Then, show the results (only the 1st and last output line matters): */
  w = W;
  i = 1 - n;
  do {
    fprintf( stderr, "%e\n", *w++ );
  } while ( (i+=2) <= end );

  return( 0 );
}

Godbolt 向我展示了由 "x86-64 gcc9.3" 生成的程序集，带有选项 "-Ofast -march=skylake" 与 " -O3 -march=skylake"。请检查网站的五个栏目（1. 源代码，2. "-Ofast" 汇编，3. "-O3" 汇编，4. 输出第 1 次装配，5. 第 2 次装配的输出）：

Godbolt site with five columns

正如您所见，程序集的差异很明显，但我无法弄清楚不精确的确切来源。那么，问题是，哪些汇编指令对此负责？

后续问题是：是否有可能通过重新编写 C 程序来避免使用“-Ofast -march=skylake”的这种不精确性？

【问题讨论】：

什么不精确？它比“正常”的 FP 不精确性更糟吗？请有问题的详细信息，而不是链接。
-Ofast 不是很稳定，可能会偏离标准 C。它可能会在某处偷工减料以降低速度的准确性。
您是否尝试过其他编译器，例如 Clang 或 CompCert ？您可能还对CADNA 软件项目或FLUCTUAT 感兴趣（请通过电子邮件联系我了解更多信息）
一眼看去，-Ofast 使用vfnmadd132sd 计算1- i*i/sqr_N_min_1 作为1 - y*z 其中y 确实是i*i 但z 是1 / sqr_N_min_1（计算在循环之前）。另一个版本使用普通的vmulsd/vmulsb/vsubsd。取倒数会影响精度，以及 FMA 比等效的三个指令序列具有更高的精度这一事实。
你知道-Ofast 是-O3 -ffast-math 的同义词，对吧？ -ffast-math 的一部分是 -funsafe-math-optimizations。这正是您通过使用该选项所要求的那种速度超过精度的优化。如果您不想这样做，请不要启用所有 -ffast-math 子选项。

标签： c assembly gcc x86-64 fast-math

【解决方案1】：

从 cmets 看来，对于 -O3，编译器会计算 1. - ((double)i * i) / sqr_N_min_1：

将i 转换为双倍并平方。
除以sqr_N_min_1。
从 1 中减去。

对于-Ofast，计算它：

在循环之前，计算sqr_N_min_1 的倒数。
将i 转换为加倍并平方。
计算融合乘减法 1 减去平方乘以倒数。

后者提高了速度，因为它只计算一次除法，并且在目标处理器中乘法比除法快得多。最重要的是，融合运算比单独的乘法和减法运算更快。

发生错误是因为倒数运算引入了原始表达式中不存在的舍入误差（1/25 不能以二进制格式精确表示，而 25/25 当然可以）。这就是为什么编译器在尝试提供严格的浮点语义时不进行这种优化的原因。

此外，只需将倒数乘以 25 即可消除错误。（这有点“偶然”，因为舍入误差以复杂的方式变化。1./25*25 产生 1，但 1./49*49 没有。）但是融合运算产生了更准确的结果（它产生的结果就像产品是精确计算，仅在减法之后进行舍入），因此它保留了错误。

【讨论】：

【解决方案2】：

评论和另一个答案指出了在您的案例中发生的具体转变，使用倒数和 FMA 而不是除法。

是否有可能通过重新编写 C 程序来避免使用“-Ofast -march=skylake”的这种不精确性？

一般不。

-Ofast（目前）是-O3 -ffast-math 的同义词。
参见https://gcc.gnu.org/wiki/FloatingPointMath

-ffast-math 的一部分是-funsafe-math-optimizations，顾名思义，可以改变数值结果。（目标是允许更多优化，例如将 FP 数学视为关联，以允许使用 SIMD 自动矢量化数组的总和，和/或展开多个累加器，或者甚至只是在一个表达式中重新排列一系列操作以组合两个单独的常量。）

这正是您通过使用该选项所要求的那种速度超过精度的优化。如果您不想这样，请不要启用所有 -ffast-math 子选项，只启用像 -fno-math-errno / -fno-trapping-math 这样的安全选项。（见How to force GCC to assume that a floating-point expression is non-negative?）

没有办法制定您的来源来避免所有可能的问题。

您可能可以在各处使用volatile tmp 变量来破坏语句之间的优化，但这会使您的代码比使用默认-fno-fast-math 的常规-O3 慢。即便如此，由于-ffinite-math-only，对sin 或log 等库函数的调用可能会解析为假定args 是有限的，而不是NaN 或无穷大的版本。

GCC issue with -Ofast? 指出另一个效果：isnan() 被优化为编译时0。

【讨论】：