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

Question

考虑double数字的排序(升序)数组。为了数值稳定性，数组应该被总结，就像从开始到结束迭代它，将和积累在某个变量中一样。

如何用AVX2有效地矢量化这一点？

我已经研究过这个方法用AVX指令实现水平向量和的最快方法，但是将其缩放到一个数组(可能需要一些分而治之的方法)似乎相当困难，同时通过确保在将小数字加到更大的数字之前，保持浮点精度。

澄清1：我认为应该可以把前4项和起来，然后把它们加到后面4项的和中，等等。我愿意用一些稳定性来换取性能。但我更喜欢一种不会完全破坏稳定性的方法。

澄清2：内存不应该成为瓶颈，因为数组位于L3缓存中(而不是在L1/L2缓存中，因为数组的片段是从不同的线程填充的)。我不想求助于Kahan求和，因为我认为真正重要的是操作的数量，而Kahan求和会使它增加4倍。

Answer 1

如果您需要精度和并行性，可以使用Kahan求和或其他错误补偿技术来重新排序和(使用多个累加器进入SIMD向量元素)。

正如两次快速求和- Evgeny Latkin所指出的，如果您在内存带宽方面遇到瓶颈，错误补偿和不会比最大性能和慢多少，因为CPU有大量的计算吞吐量，而这些计算吞吐量在内存带宽瓶颈的简单并行和中没有使用。

另见(kahan summation avx的google结果)

• 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的示例实现(而不是展开)，并比较了实际错误(无错误)与顺序和(糟糕)和简单SIMD和(更不用说总错误)，表明使用多个累加器进行矢量化(和/或展开)有助于提高精度。

您的想法:按顺序对4个数字进行求和将使您隐藏FP添加延迟，以及标量添加吞吐量上的瓶颈。

在向量中做水平和需要大量的调整，因此这是一个潜在的瓶颈。您可以考虑加载a0 a1 a2 a3，然后调整以获得a0+a1 x a2+a3 x，然后是(a0+a1) + (a2+a3)。你有一个瑞森，对吧？最后一步是将vextractf128降到128 B。这仍然是3个总添加uop，和3个洗牌uops，但指令比标量或128 B向量少。

你的想法非常类似于成对求和。

你总是会得到一些四舍五入的错误，但是如果加上类似的数值，它就会最小化。

关于成对求和和简单有效的SIMD的一些评论，请参见 Simd matmul程序给出了不同的数值结果。

添加4个连续数与垂直添加4个SIMD向量之间的差异可能可以忽略不计。SIMD矢量在数组中给你很小的步幅( SIMD矢量宽度)。除非阵列增长得非常快，否则它们的震级仍将基本相似。

，您不需要水平求和，直到结束仍然得到大部分的利益。您可以维护1或2个SIMD矢量累加器，而在添加到主累加器之前，可以使用更多SIMD寄存器来和短运行(可能是4或8个SIMD向量)。

事实上，拥有更多的方法(跨SIMD向量元素)意味着它不会增长那么大。因此，它确实有助于解决您想要避免的问题，而水平求和到单个标量累加器实际上会使情况变得更糟，特别是对于严格排序的数组。

有了无序的执行，您就不需要太多的tmp累加器来完成这项工作，并隐藏将累积到主累加器中的FP添加延迟。您可以在一个新的tmp = _mm_load_ps()向量中进行几组积累，并将其添加到总数中，OoO执行部分将与这些执行重叠。因此，您不需要一个巨大的展开因子为您的主回路。

但是它不应该太小，您不希望将延迟添加到主累加器上，等待前面的添加在下一个累加器开始之前生成一个结果。您希望限制FP-添加吞吐量。(或者，如果您关心布罗德韦尔/哈斯韦尔，并且您不完全限制内存带宽，那么可以将一些FMA与1.0乘数混合起来，以利用这种吞吐量。)

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

有关多个累加器的更多信息，请参见ps？。

Answer 2

到目前为止，我的解决方案如下：

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]等，然后添加相邻项的和。

Answer 3

你想玩的游戏可能适得其反。尝试从您最喜欢的分布中生成一组iid样本，对它们进行排序，并将“和按递增顺序”与“每条车道按递增顺序求和，然后将车道和和”进行比较。

对于4条车道和16条数据，求和法给我的误差较小--约28% --按递增顺序求和给我的误差较小--约17%；对于4条车道和256条数据，求和法给我的误差较小--约68% --而累加顺序给我的误差约为12%。求和算法也击败了你在自己的答案中给出的算法。为此，我在0,1上使用了统一的分布。