17

有谁知道为什么 GCC/Clang 不会在下面的代码示例中优化函数test1以在使用快速数学选项时仅使用 RCPPS 指令？是否有另一个编译器标志会生成此代码？

typedef float float4 __attribute__((vector_size(16)));

float4 test1(float4 v)
{
    return 1.0f / v;
}

您可以在此处查看编译后的输出：https ://goo.gl/jXsqat

4

2 回答 2

19

因为精度RCPPS比除法低很多。float

启用该优化的选项不适合作为-ffast-math.

gcc 手册的x86 目标选项说实际上有一个选项 (with -ffast-math) 确实让 gcc 使用它们（使用 Newton-Raphson 迭代 -快速矢量化 rsqrt 并与 SSE/AVX 互惠，具体取决于精度/ Newton Raphson 与 SSE2 - 有人可以向我解释这 3 行- SIMD 和标量每条指令的性能基本相同，牛顿迭代数学是相同的）：

  • -mrecip 此选项允许使用 RCPSS 和 RSQRTSS 指令（及其矢量化变体 RCPPS 和 RSQRTPS）以及额外的 Newton-Raphson 步骤来提高单精度浮点参数的精度，而不是 DIVSS 和 SQRTSS（及其矢量化变体）。仅当 -funsafe-math-optimizations 与 -finite-math-only 和 -fno-trapping-math 一起启用时，才会生成这些指令。请注意，虽然序列的吞吐量高于单向指令的吞吐量，但序列的精度最多可降低 2 ulp（即 1.0 的倒数等于 0.99999994）。

请注意，GCC 已经使用 -ffast-math（或上述选项组合）根据 RSQRTSS（或 RSQRTPS）实现 1.0f/sqrtf(x)，并且不需要 -mrecip。

另请注意，GCC 已使用 -ffast-math（或上述选项组合）通过额外的 Newton-Raphson 步骤发出上述序列，用于矢量化单浮点除法和矢量化 sqrtf(x)，并且不需要 -mrecip。

  • -mrecip=opt

此选项控制可以使用哪些倒数估计指令。opt 是以逗号分隔的选项列表，前面可能有一个“！” 反转选项：

’all’
      Enable all estimate instructions.
‘default’
    Enable the default instructions, equivalent to -mrecip.
‘none’
    Disable all estimate instructions, equivalent to -mno-recip.
‘div’
    Enable the approximation for scalar division.
‘vec-div’
    Enable the approximation for vectorized division.
‘sqrt’
    Enable the approximation for scalar square root.
‘vec-sqrt’
    Enable the approximation for vectorized square root.

因此，例如，-mrecip=all,!sqrt 启用除平方根之外的所有倒数近似。

请注意，英特尔全新的 Skylake 设计进一步提高了 FP 划分性能，达到 8-11c 延迟，1/3c 吞吐量。（或 256b 向量每 5c 吞吐量一个，但延迟相同vdivps）。他们扩大了分隔线，因此 AVXvdivps ymm现在的延迟与 128b 向量的延迟相同。

（SnB 到 Haswell 的 256b div 和 sqrt 具有大约两倍的延迟/接收吞吐量，因此它们显然只有 128b 宽的分隔线。）Skylake 还对这两个操作进行了更多的管道化，因此大约 4 个 div 操作可以在运行中。sqrt 也更快。

所以几年后，一旦 Skylake 普及，只有在rcpps需要多次除以同一事物时才值得这样做。 rcpps一对夫妇fma的吞吐量可能略高，但延迟更差。此外，vdivps只有一个 uop；因此，将有更多的执行资源可用于与部门同时发生的事情。

AVX512 的初始实施会是什么样子还有待观察。rcpps如果 FP 除法性能是瓶颈，那么 Newton-Raphson 迭代的几个 FMA 可能会是一个胜利。如果 uop 吞吐量是一个瓶颈，并且在除法运行时还有很多其他工作要做，那么vdivps zmm可能仍然很好（当然，除非重复使用相同的除数）。

浮点除法与浮点乘法更多地是关于 FP 吞吐量与延迟。

于 2015-08-14T04:31:54.930 回答
2

我正在我的一个应用程序中尝试使用浮点数学重的热路径，并发现了类似的东西。我通常不看编译器发出的指令，所以我有点惊讶并深入研究了数学细节。

这是 gcc 生成的一组指令，由我用携带的计算进行注释：

test1(float __vector): ; xmm0               = a
    rcpps   xmm1, xmm0 ; xmm1 = 1 / xmm0    = 1/a
    mulps   xmm0, xmm1 ; xmm0 = xmm0 * xmm1 = a * 1/a
    mulps   xmm0, xmm1 ; xmm0 = xmm0 * xmm1 = a * (1/a)^2
    addps   xmm1, xmm1 ; xmm1 = xmm1 + xmm1 = 2 * (1/a)
    subps   xmm1, xmm0 ; xmm1 = xmm1 - xmm0 = 2 * (1/a) - a * (1/a)^2
    movaps  xmm0, xmm1 ; xmm0 = xmm1        = 2 * (1/a) - a * (1/a)^2
    ret

那么这里发生了什么？为什么要浪费额外的 4 条指令来计算在数学上等于倒数的表达式？

好吧，这些rcpps指令只是一个近似的倒数。其他算术指令 ( mulps, addps, subps) 精确到单精度。让我们写出r(x)近似倒数函数。最后变成y = 2*r(a) - a*r(a)^2。如果我们用相对误差代替r(x) = (1 + eps) * (1/x)，我们得到：eps

y = 2 * (1 + eps) * (1/a) - a * (1 + eps)^2 * (1/a)^2
  = (2 + 2*eps - (1 + eps)^2) * (1/a)
  = (2 + 2*eps - (1 + 2*eps + eps^2)) * (1/a)
  = (1 - eps^2) * (1/a)

的相对误差rcpps 小于 1.5 * 2^-12，所以eps <= 1.5 * 2^-12，所以：

eps^2 <= 2.25 * 2^-24
      <  1.5  * 2^-23

因此，通过执行这些额外的指令，我们从 12 位精度提高到 23 位精度。请注意，单精度浮点数具有 24 位精度，因此我们在这里几乎可以得到全精度。

那么这只是一些神奇的指令序列，恰好可以让我们获得额外的精度吗？不完全的。它基于牛顿的方法（我认为经常使用装配的人将其称为 Newton-Raphson）。

牛顿法是一种寻根法。给定某个函数f(x)，它f(x) = 0通过从一个近似解开始x_0并对其进行迭代改进来找到 的近似解。牛顿迭代由下式给出：

x_n+1 = x_n - f(x_n) / f'(x_n)

在我们的例子中，我们可以用导数将求函数的倒数重新表述1/a为a求函数的根。将其代入牛顿迭代方程，我们得到：f(x) = a*x - 1f'(x) = a

x_n+1 = x_n - (a*x_n - 1) / a

两个观察：

  1. 在这种情况下，牛顿迭代实际上给了我们一个精确的结果，而不仅仅是一个更好的近似值。这是有道理的，因为牛顿的方法是通过对f周围进行线性近似来工作的x_n。在这种情况下f是完全线性的，所以近似是完美的。然而...

  2. 计算牛顿迭代需要我们除以a，这是我们试图近似的精确计算。这就产生了一个循环问题。我们通过修改牛顿迭代来打破循环，使用我们的近似倒数x_n除以a：

    x_n+1 =  x_n - (a*x_n - 1) * x_n
      ~= x_n - (a*x_n - 1) / a

这个迭代可以正常工作，但从向量数学的角度来看并不是很好：它需要减去1. 要使用向量数学来做到这一点，需要准备一个带有 1 序列的向量寄存器。这需要额外的指令和额外的寄存器。

我们可以重写迭代来避免这种情况：

x_n+1 = x_n - (a*x_n - 1) * x_n
      = x_n - (a*x_n^2 - x_n)
      = 2*x_n - a*x_n^2

现在替换x_0 = r(a)，我们从上面恢复我们的表达式：

y = x_1 = 2*r(a) - a*r(a)^2
于 2020-11-29T01:33:21.687 回答