AVX2 中排序数组的高效稳定和

Efficient stable sum of a sorted array in AVX2

本文关键字:高效 排序 数组 AVX2      更新时间:2023-10-16

考虑一个double数字的排序(升序)数组。为了数值稳定性,数组应该被总结为从头到尾迭代它,在某个变量中累积总和。

如何使用 AVX2 有效地对此进行矢量化?

我已经研究了这种方法 使用 AVX 指令进行水平向量求和的最快方法,但将其扩展到数组似乎很棘手(可能需要一些除法和征服方法),同时通过确保在将它们添加到较大数字之前对小数字求和来保持浮点精度。

澄清1:我认为例如将前4个项目相加,然后将它们添加到接下来4个项目的总和中,等等。我愿意用一些稳定性来换取性能。但我更喜欢一种不会完全破坏稳定性的方法。

澄清 2:内存不应该成为瓶颈,因为数组位于 L3 缓存中(但不在 L1/L2 缓存中,因为数组的各个部分是从不同的线程填充的)。我不想求助于 Kahan 求和,因为我认为重要的是操作的数量,而 Kahan 求和会增加大约 4 倍。

如果您需要精度并行性,请使用 Kahan 求和或其他误差补偿技术来重新排序总和(使用多个累加器进入 SIMD 矢量元素步幅)。

正如 Twodouble Fast sumation - Evgeny Latkin 指出的那样,如果你在内存带宽上遇到瓶颈,误差补偿的总和并不比最大性能总和慢多少,因为 CPU 有大量的计算吞吐量在简单并行的总和中未使用,这在内存带宽上成为瓶颈。

另请参阅(谷歌搜索结果kahan summation avx)

  • https://github.com/rreusser/summation-algorithms

  • https://scicomp.stackexchange.com/questions/10869/which-algorithm-is-more-accurate-for-computing-the-sum-of-a-sorted-array-of-numb

  • 这种使用 SSE
  • 过度杀伤处理数组尾部的方式吗? 有一个 SSE Kahan 的示例实现,未展开,并比较了实际误差(无误差)与顺序总和(坏)与简单 SIMD 总和(总误差少得多),表明仅使用多个累加器进行矢量化(和/或展开)往往有助于准确性。


Re:您的想法:按顺序对 4 个数字的组求和将让您隐藏 FP 添加延迟和标量添加吞吐量的瓶颈。

在向量中进行水平求和需要大量的洗牌,因此这是一个潜在的瓶颈。 您可以考虑加载a0 a1 a2 a3,然后洗牌以获得a0+a1 x a2+a3 x,然后(a0+a1) + (a2+a3)。 你有锐龙,对吧? 最后一步只是vextractf128到128b。 这仍然是总共 3 个 ADD uops 和 3 个随机 uops,但指令比标量或 128b 向量少。


你的想法与成对求和非常相似

你总是会得到一些舍入误差,但添加类似量级的数字可以最大限度地减少它。

另请参阅Simd matmul 程序给出了不同的数值结果,以获取有关成对求和和简单高效 SIMD 的一些注释。

添加 4 个连续数字与垂直添加 4 个 SIMD 向量之间的差异可能可以忽略不计。 SIMD 矢量可在阵列中为您提供小步幅(SIMD 矢量宽度)。 除非阵列增长得非常快,否则它们仍将具有基本相似的幅度。

你不需要水平求和,直到最后仍然获得大部分好处。 您可以维护 1 个或 2 个 SIMD 矢量累加器,同时使用更多 SIMD 寄存器对短期运行(可能是 4 或 8 个 SIMD 矢量)求和,然后再添加到主累加器中。

事实上,以更多方式(跨 SIMD 矢量元素)拆分意味着它不会变得那么大。 因此,它有助于解决您试图避免的问题,并且水平求和到单个标量累加器实际上会使事情变得更糟,尤其是对于严格排序的数组。

通过无序执行,您不需要太多的 tmp 累加器来完成这项工作,并隐藏累加到主累加器的 FP 添加延迟。 您可以执行几组累加到一个新的tmp = _mm_load_ps()向量中并将其添加到总数中,OoO exec 将重叠这些执行。 因此,您的主循环不需要巨大的展开系数。

但它不应该太小,你不想在主累加器中添加延迟的瓶颈,等待上一个添加产生结果,然后才能开始下一个。 您希望成为 FP 添加吞吐量的瓶颈。 (或者,如果您关心 Broadwell/Haswell 并且您没有完全成为内存带宽的瓶颈,请混合一些带有1.0乘数的 FMA 以利用该吞吐量。

例如,Skylake SIMD FP 加法具有 4 个周期延迟,0.5 个周期吞吐量,因此您需要至少执行 7 次加法,这些加法是单个累加器中每次添加的短 dep 链的一部分。 最好更多,和/或最好有2个长期累加器,以更好地吸收资源冲突调度中的气泡。

请参阅_mm256_fmadd_ps比_mm256_mul_ps + _mm256_add_ps慢?,了解有关多个累加器的更多信息。

这是我到目前为止的解决方案:

double SumVects(const __m256d* pv, size_t n) {
if(n == 0) return 0.0;
__m256d sum = pv[0];
if(n == 1) {
sum = _mm256_permute4x64_pd(sum, _MM_SHUFFLE(3, 1, 2, 0));
} else {
for(size_t i=1; i+1 < n; i++) {
sum = _mm256_hadd_pd(sum, pv[i]);
sum = _mm256_permute4x64_pd(sum, _MM_SHUFFLE(3, 1, 2, 0));
}
sum = _mm256_hadd_pd(sum, pv[n-1]);
}
const __m128d laneSums = _mm_hadd_pd(_mm256_extractf128_pd(sum, 1),
_mm256_castpd256_pd128(sum));
return laneSums.m128d_f64[0] + laneSums.m128d_f64[1];
}

一些解释:它首先添加相邻double数组项,例如a[0]+a[1]a[2]+a[3]等。然后,它将相邻项的总和相加。

你想玩的游戏可能会适得其反。 尝试通过从您喜欢的分布中生成一堆 iid 样本进行实验,对它们进行排序,并将"按递增顺序求和"按递增顺序对每个泳道求和,然后对车道总和求和"进行比较。

对于 4 条车道和 16 个数据,逐巷求和在大约 28% 的时间内给我带来较小的误差,而按递增顺序求和给我的误差大约 17% ;对于 4 条车道和 256 个数据,巷道求和给我的误差大约 68% 的时间较小,而按递增顺序求和给我的误差大约 12% 的时间较小。 巷道求和也击败了您在自我回答中给出的算法。 为此,我在 [0,1] 上使用了均匀分布。