最快的水平SSE向量求和(或其他缩减)方法

83

给定一个包含三个(或四个)浮点数的向量,最快的求和方法是什么?

SSE(movaps、shuffle、add、movd)是否总是比x87更快?SSE3中的水平加法指令是否值得使用?

从FPU移动到faddp、faddp的成本是多少?最快的特定指令序列是什么?

“尝试安排事物,以便您可以同时对四个向量求和”不会被接受为答案。 :-) 比如对于一个数组的求和,您可以使用多个矢量累加器进行垂直求和(以隐藏addps的延迟),并在循环结束后将其缩减为一个,但然后需要对该最后一个向量进行水平求和。


2
如果你认为横向加法对性能至关重要,那么你可能并没有以最佳方式接近SIMD编码。请发一些代码以展示你需要在哪里进行优化。 - Paul R
1
向量之间的点积,主要用于计算它们之间的夹角。请注意最后一句话。 - FeepingCreature
我读了最后一句话,但我仍然认为可能有更好的方法。 - Paul R
我知道有一种更好的方法,那就是“每次执行四个元素的循环,这样你就可以并行化所有操作”。问题是,如果不采用这种复杂而难以理解的方式,我们还能做到最好的是什么? - FeepingCreature
3
“在x86上没有最快的方法……”不适用。不同的x86处理器具有不同的执行特性。你的目标处理器是什么?你的“三个浮点数向量”最初是在内存中,还是连续地存储在SSE寄存器中,或者其他地方? - Stephen Canon
显示剩余3条评论
5个回答

144
一般来说,对于任何类型的向量水平缩减,提取/洗牌高半部分以与低部分对齐,然后进行垂直加法(或min/max/or/and/xor/multiply/whatever); 重复此过程直到只剩下一个元素(其余向量中含有高垃圾)。
如果您从宽度大于128位的向量开始,将其缩小一半直到达到128位(然后可以在该向量上使用此答案中的一个函数)。但是,如果您需要在最后广播结果到所有元素,则可以考虑一直使用全宽度洗牌。
与更宽的向量、整数和 FP 相关的 Q&A。

整数


这个问题的主要答案:大多数是浮点数和__m128

以下是根据Agner Fog的微架构指南和指令表进行调整的一些版本。请参见标签wiki。它们应该在任何CPU上都很高效,没有主要瓶颈。(例如,我避免了一些可能对某个uarch有所帮助但在另一个uarch上速度较慢的东西)。代码大小也被最小化。

常见的SSE3 / SSSE3 2x hadd习惯用法只适用于代码大小,而不是任何现有CPU的速度。它有其用例(如转置和加法,请参见下文),但单个向量不是其中之一。

我还包括了一个AVX版本。任何使用AVX / AVX2进行水平约简的操作都应该从一个vextractf128和一个“垂直”操作开始,以将其减少到一个XMM(__m128)向量。通常对于宽向量,最好是重复将其缩小一半,直到缩小到128位向量,而不管元素类型如何。(除了8位整数,如果要进行hsum而不会溢出到更宽的元素,则应首先使用vpsadbw。)

从所有这些代码中查看汇编输出 在Godbolt Compiler Explorer上 还请查看我对Agner Fog的C++向量类库horizontal_add函数的改进(留言板线程,以及github上的代码)。我使用CPP宏来选择适用于SSE2、SSE4和AVX的最佳洗牌方式,并避免在AVX不可用时使用movdqa


需要考虑权衡:

  • 代码大小:出于L1 I-cache的原因,越小越好,并且对于从磁盘获取代码(更小的二进制文件)也是如此。总二进制文件大小在整个程序中反复进行编译器决策时非常重要。如果你正在费心手动编写具有内部函数的代码,则值得花费一些代码字节,如果它能为整个程序带来任何加速效果(要注意不要让展开看起来很好的微基准测试误导自己)
  • uop-cache大小:通常比L1 I$更珍贵。4个单一uop指令所占用的空间可以比2个haddps少,因此这在这里非常相关。
  • 延迟:有时相关
  • 吞吐量(后端端口):通常不相关,水平求和不应该是最内层循环中唯一的事情。端口压力只在包含此内容的整个循环中才有影响。
  • 吞吐量(总前端融合域uops):如果周围代码没有在与hsum使用相同的端口上瓶颈,那么这是hsum对整个事物吞吐量影响的代理。

当水平相加不频繁时:

如果CPU没有uop缓存,则很少使用的情况下,可能会偏爱2x haddps:当它运行时速度较慢,但这种情况并不经常发生。只有2个指令可以最小化对周围代码(I$大小)的影响。

如果CPU有uop缓存,则可能更喜欢需要较少uop的东西,即使它需要更多指令/更多的x86代码大小。我们要最小化使用的总uop缓存行,这并不像最小化总uop那么简单(采取分支和32B边界始终启动新的uop缓存行)。

无论如何,水平求和经常出现,因此我尝试精心制作了一些可以编译的版本。没有在任何真正的硬件上进行基准测试,甚至没有仔细测试。洗牌常数或其他方面可能存在错误。


如果您正在制作代码的回退/基线版本,请记住只有旧CPU才能运行它; 新er的CPU将运行您的AVX版本,或SSE4.1或其他版本。
旧CPU(如K8和Core2(merom)及更早版本)仅具有64位洗牌单元。Core2对于大多数指令具有128位执行单元,但对于洗牌而言不是这样。(Pentium M和K8将所有128b矢量指令视为两个64位半部分处理)。
像movhlps这样以64位块移动数据的洗牌(在64位半部分内没有洗牌)也很快。
相关:新CPU上的洗牌,以及避免Haswell及更高版本中1/clock shuffle吞吐量瓶颈的技巧:在AVX512中执行128位跨lane操作是否会提供更好的性能? 在旧CPU上,由于洗牌速度较慢:
  • movhlps(Merom:1uop)比shufps(Merom:3uops)快得多。在Pentium-M上,比movaps便宜。此外,在Core2上,它在FP域中运行,避免了来自其他洗牌的旁路延迟。
  • unpcklpdunpcklps快。
  • pshufd很慢,pshuflw/pshufhw很快(因为它们只混洗64位的一半)
  • pshufb mm0(MMX)很快,pshufb xmm0很慢。
  • haddps非常慢(Merom和Pentium M上的6uops)
  • movshdup(Merom:1uop)很有趣:它是唯一在64b元素内混洗的1uop指令。
shufps 在Core2(包括Penryn)上将数据带入整数域,导致绕过延迟以使其返回FP执行单元用于addps,但movhlps完全在FP域中。 shufpd也在浮点域中运行。 movshdup在整数域中运行,但只有一个uop。
AMD K10、Intel Core2(Penryn/Wolfdale)和所有更高版本的CPU都将所有xmm洗牌作为单个uop运行。(但请注意在Penryn上使用shufps时的绕行延迟,可通过movhlps避免)

没有AVX,避免浪费的movaps/movdqa指令需要仔细选择洗牌操作。只有少数几个洗牌操作可以作为复制和洗牌操作,而不是修改目标操作。将两个输入数据组合的洗牌操作(例如unpck*movhlps)可以与一个不再需要的tmp变量一起使用,而不是使用_mm_movehl_ps(same,same)

其中一些可以通过使用一个虚拟参数作为初始洗牌的目标来使其更快(节省MOVAPS),但会让代码更丑陋/不够"简洁"。例如:

// Use dummy = a recently-dead variable that vec depends on,
//  so it doesn't introduce a false dependency,
//  and the compiler probably still has it in a register
__m128d highhalf_pd(__m128d dummy, __m128d vec) {
#ifdef __AVX__
    // With 3-operand AVX instructions, don't create an extra dependency on something we don't need anymore.
    (void)dummy;
    return _mm_unpackhi_pd(vec, vec);
#else
    // Without AVX, we can save a MOVAPS with MOVHLPS into a dead register
    __m128 tmp = _mm_castpd_ps(dummy);
    __m128d high = _mm_castps_pd(_mm_movehl_ps(tmp, _mm_castpd_ps(vec)));
    return high;
#endif
}

SSE1(又称SSE):

float hsum_ps_sse1(__m128 v) {                                  // v = [ D C | B A ]
    __m128 shuf   = _mm_shuffle_ps(v, v, _MM_SHUFFLE(2, 3, 0, 1));  // [ C D | A B ]
    __m128 sums   = _mm_add_ps(v, shuf);      // sums = [ D+C C+D | B+A A+B ]
    shuf          = _mm_movehl_ps(shuf, sums);      //  [   C   D | D+C C+D ]  // let the compiler avoid a mov by reusing shuf
    sums          = _mm_add_ss(sums, shuf);
    return    _mm_cvtss_f32(sums);
}
    # gcc 5.3 -O3:  looks optimal
    movaps  xmm1, xmm0     # I think one movaps is unavoidable, unless we have a 2nd register with known-safe floats in the upper 2 elements
    shufps  xmm1, xmm0, 177
    addps   xmm0, xmm1
    movhlps xmm1, xmm0     # note the reuse of shuf, avoiding a movaps
    addss   xmm0, xmm1

    # clang 3.7.1 -O3:  
    movaps  xmm1, xmm0
    shufps  xmm1, xmm1, 177
    addps   xmm1, xmm0
    movaps  xmm0, xmm1
    shufpd  xmm0, xmm0, 1
    addss   xmm0, xmm1

我报告了一个有关 clang bug about pessimizing the shuffles 的问题。它有自己的内部表示来进行洗牌,并将其转换回洗牌。gcc更经常使用与您使用的内置函数直接匹配的指令。

通常情况下,当指令选择没有手动调整或常量传播可以简化非常量情况下的最优内置函数时,clang比gcc表现更好。总的来说,编译器像正常的编译器一样处理内置函数是件好事,而不仅仅是汇编程序。编译器通常可以从标量C生成良好的汇编代码,即使该代码不能像良好的汇编代码那样运行。最终,编译器将把内置函数视为输入优化器的另一个C运算符。


SSE3

float hsum_ps_sse3(__m128 v) {
    __m128 shuf = _mm_movehdup_ps(v);        // broadcast elements 3,1 to 2,0
    __m128 sums = _mm_add_ps(v, shuf);
    shuf        = _mm_movehl_ps(shuf, sums); // high half -> low half
    sums        = _mm_add_ss(sums, shuf);
    return        _mm_cvtss_f32(sums);
}

    # gcc 5.3 -O3: perfectly optimal code
    movshdup    xmm1, xmm0
    addps       xmm0, xmm1
    movhlps     xmm1, xmm0
    addss       xmm0, xmm1

这样做有几个优点:
  • 不需要任何movaps复制来解决破坏性洗牌的问题(没有AVX):movshdup xmm1,xmm2的目标是只写的,因此它为我们创建了一个死寄存器的tmp。这也是我使用movehl_ps(tmp, sums)而不是movehl_ps(sums, sums)的原因。

  • 代码大小小。混洗指令很小:movhlps为3字节,movshdup为4字节(与shufps相同)。不需要立即字节,所以使用AVX,vshufps为5字节,但vmovhlpsvmovshdup都是4字节。

我可以使用addps代替addss来节省一个字节。由于这不会在内部循环中使用,因此切换额外的晶体管所需的额外能量可能是可忽略的。来自上3个元素的FP异常并不是风险,因为所有元素都包含有效的FP数据。然而,clang/LLVM实际上“理解”向量混洗,并且如果它知道仅低元素重要,则会发出更好的代码。
与SSE1版本一样,将奇数元素加起来可能会导致FP异常(如溢出),否则不会发生,但这不应该是问题。Denormals很慢,但我IRC产生+Inf结果在大多数uarches上并不常见。

针对代码大小的SSE3优化

如果代码大小是您关注的主要问题,那么两个haddps (_mm_hadd_ps)指令可以解决问题(Paul R的答案)。这也是最容易输入和记忆的。然而,它不快。即使是Intel Skylake仍将每个haddps解码为3个uop,并具有6个周期的延迟。因此,即使它节省了机器码字节(L1 I-cache),它在更有价值的uop-cache中占用更多空间。 haddps的实际用例:一个转置和求和问题, 或在中间步骤执行一些缩放 SSE atoi() 实现中


AVX:

这个版本相比Marat对AVX问题的回答可以节省一个代码字节。

#ifdef __AVX__
float hsum256_ps_avx(__m256 v) {
    __m128 vlow  = _mm256_castps256_ps128(v);
    __m128 vhigh = _mm256_extractf128_ps(v, 1); // high 128
           vlow  = _mm_add_ps(vlow, vhigh);     // add the low 128
    return hsum_ps_sse3(vlow);         // and inline the sse3 version, which is optimal for AVX
    // (no wasted instructions, and all of them are the 4B minimum)
}
#endif

 vmovaps xmm1,xmm0               # huh, what the heck gcc?  Just extract to xmm1
 vextractf128 xmm0,ymm0,0x1
 vaddps xmm0,xmm1,xmm0
 vmovshdup xmm1,xmm0
 vaddps xmm0,xmm1,xmm0
 vmovhlps xmm1,xmm1,xmm0
 vaddss xmm0,xmm0,xmm1
 vzeroupper 
 ret

双精度:

double hsum_pd_sse2(__m128d vd) {                      // v = [ B | A ]
    __m128 undef  = _mm_undefined_ps();                       // don't worry, we only use addSD, never touching the garbage bits with an FP add
    __m128 shuftmp= _mm_movehl_ps(undef, _mm_castpd_ps(vd));  // there is no movhlpd
    __m128d shuf  = _mm_castps_pd(shuftmp);
    return  _mm_cvtsd_f64(_mm_add_sd(vd, shuf));
}

# gcc 5.3.0 -O3
    pxor    xmm1, xmm1          # hopefully when inlined, gcc could pick a register it knew wouldn't cause a false dep problem, and avoid the zeroing
    movhlps xmm1, xmm0
    addsd   xmm0, xmm1


# clang 3.7.1 -O3 again doesn't use movhlps:
    xorpd   xmm2, xmm2          # with  #define _mm_undefined_ps _mm_setzero_ps
    movapd  xmm1, xmm0
    unpckhpd        xmm1, xmm2
    addsd   xmm1, xmm0
    movapd  xmm0, xmm1    # another clang bug: wrong choice of operand order


// This doesn't compile the way it's written
double hsum_pd_scalar_sse2(__m128d vd) {
    double tmp;
    _mm_storeh_pd(&tmp, vd);       // store the high half
    double lo = _mm_cvtsd_f64(vd); // cast the low half
    return lo+tmp;
}

    # gcc 5.3 -O3
    haddpd  xmm0, xmm0   # Lower latency but less throughput than storing to memory

    # ICC13
    movhpd    QWORD PTR [-8+rsp], xmm0    # only needs the store port, not the shuffle unit
    addsd     xmm0, QWORD PTR [-8+rsp]

将数据存储到内存并再次读取可以避免 ALU uop。如果洗牌端口压力或 ALU uop 是瓶颈的话,这是非常好的。(请注意,它不需要 sub rsp, 8 或任何其他操作,因为 x86-64 SysV ABI 提供了一个红区,信号处理程序不会使用它。)
有些人将数据存储到数组中,并对所有元素求和,但编译器通常无法意识到数组的低位元素仍然保存在之前的寄存器中。

整数:

pshufd 是一个方便的复制和重排操作。位移和字节移位不幸地是原地操作,而 punpckhqdq 将目标的高半部分放入结果的低半部分,与 movhlps 把高半部分提取到另一个寄存器的方式相反。

在某些 CPU 上使用 movhlps 作为第一步可能是个好选择,但只有我们有一个临时寄存器时才行。pshufd 是一个安全的选择,在 Merom 之后的所有 CPU 上速度都很快。

int hsum_epi32_sse2(__m128i x) {
#ifdef __AVX__
    __m128i hi64  = _mm_unpackhi_epi64(x, x);           // 3-operand non-destructive AVX lets us save a byte without needing a mov
#else
    __m128i hi64  = _mm_shuffle_epi32(x, _MM_SHUFFLE(1, 0, 3, 2));
#endif
    __m128i sum64 = _mm_add_epi32(hi64, x);
    __m128i hi32  = _mm_shufflelo_epi16(sum64, _MM_SHUFFLE(1, 0, 3, 2));    // Swap the low two elements
    __m128i sum32 = _mm_add_epi32(sum64, hi32);
    return _mm_cvtsi128_si32(sum32);       // SSE2 movd
    //return _mm_extract_epi32(hl, 0);     // SSE4, even though it compiles to movd instead of a literal pextrd r32,xmm,0
}

    # gcc 5.3 -O3
    pshufd xmm1,xmm0,0x4e
    paddd  xmm0,xmm1
    pshuflw xmm1,xmm0,0x4e
    paddd  xmm0,xmm1
    movd   eax,xmm0

int hsum_epi32_ssse3_slow_smallcode(__m128i x){
    x = _mm_hadd_epi32(x, x);
    x = _mm_hadd_epi32(x, x);
    return _mm_cvtsi128_si32(x);
}

在一些CPU上,使用FP shuffle处理整数数据是安全的。但我没有这样做,因为在现代CPU上,这最多只能节省1或2个代码字节,并且没有速度提升(除了代码大小/对齐效果)。


2
@plasmacel:在许多CPU上,包括英特尔SnB系列,存在额外的旁路延迟来将FP指令的结果转发到整数洗牌器,以及从PSHUFD到ADDPS。如果您关心吞吐量和uop计数而不是延迟,则非常好。 (在SnB系列上,在整数指令之间进行SHUFPS没有惩罚(与Nehalem不同),但反过来则不成立。) - Peter Cordes
2
如果您有特定的微架构和编译器,您可以并且应该制作一个更适合它们的版本。这个答案试图在现代CPU(如Haswell)上优化(延迟、吞吐量和代码大小),同时尽可能少地影响旧的CPU。也就是说,我的SSE1 / SSE2版本不会做任何对Haswell更差的事情,只为了在旧的SlowShuffle CPU(如Merom)上运行得更快。对于Merom来说,PSHUFD可能是一个胜利,因为它和SHUFPS都在flt->int域中运行。 - Peter Cordes
2
@plasmacel:不会,除非你的向量一开始就在内存中,因为VPERMILPS可以加载+洗牌。使用旧指令的AVX版本可以获得更小的代码大小,因为您不需要立即数,它们只需要2字节的VEX前缀(C5..而不是C4....)。像VSHUFPS和VMOVHLPS这样的双源洗牌与像VPSHUFD或VPERMILPS这样的单源洗牌一样快。如果有能耗差异,那可能是微不足道的。 - Peter Cordes
2
@plasmacel:正如我的回答所指出的那样,我的SSE3版本在AVX下编译最优,但是clang将其悲观地转换为VPERMILPD: https://godbolt.org/g/ZH88wH。gcc的版本是四个4B指令(不包括RET)。clang的版本比较长2个字节,但速度相同。你认为VPERMILPS比SHUFPS更好吗?据我所知,clang错误地偏爱它用于源已经在寄存器中的立即洗牌。Agner Fog的表格没有显示任何差异。它对于加载+洗牌和变量洗牌很有用,也许对于编译器来说更容易,因为它是一个1输入指令,但不会更快。 - Peter Cordes
1
@jww:不,我不想把任何东西定下来,如果/当我意识到我的建议并不是最优的时候,我无法回来编辑。我曾经收到一封电子邮件,问我是否想参与编写汇编书,但我从未回复他们 >.< 无论如何,收集我和其他人编写的“食谱”式的更有用的SO答案链接将是一个好项目,如果我有时间的话。 - Peter Cordes
显示剩余18条评论

18

SSE2

全部四个:

const __m128 t = _mm_add_ps(v, _mm_movehl_ps(v, v));
const __m128 sum = _mm_add_ss(t, _mm_shuffle_ps(t, t, 1));

r1+r2+r3:

const __m128 t1 = _mm_movehl_ps(v, v);
const __m128 t2 = _mm_add_ps(v, t1);
const __m128 sum = _mm_add_ss(t1, _mm_shuffle_ps(t2, t2, 1));

我发现这些与双 HADDPS 大致相同的速度(但我没有太仔细地测量)。


15

你可以在SSE3中用两个HADDPS指令完成它:

v = _mm_hadd_ps(v, v);
v = _mm_hadd_ps(v, v);

这将所有元素中的总和放置其中。


2
总和难道不会出现在所有元素中吗? - Jens Björnhager
@Jens:是的,谢谢 - 我想你是对的 - 我会更新我的回答。 - Paul R
对于一个三维向量求和,我需要先将第四个分量设置为零。有什么最快的方法吗?我倾向于使用“加载掩码和ps” - 有没有一种快速屏蔽元素的方法? - FeepingCreature
@FeepingCreature __m128 vector3 = _mm_castps_si128(_mm_castsi128_ps(_mm_srli_si128(vector4, 4))); - 根据您的掩码是否已经从内存中加载,这可能比遮罩更快 - awdz9nld
1
@Royi:请查看Peter在他的答案中,标题为“SSE3优化代码大小”的评论。 - Paul R
显示剩余2条评论

4
我建议您一定要尝试使用SSE 4.2。如果您需要多次执行此操作(我假设您需要,因为性能是一个问题),那么您可以将一个寄存器预先加载为(1,1,1,1),然后对其进行几个dot4(my_vec(s),one_vec)。是的,它会进行一个多余的乘法运算,但这些在今天来说相对便宜,并且这样的操作可能被水平依赖所支配,而这些依赖可能在新的SSE点积函数中得到更优化。您应该测试是否比Paul R发布的双重水平加法更快。
我还建议将其与直接标量(或标量SSE)代码进行比较-奇怪的是,它通常更快(通常是因为内部序列化但使用寄存器旁路紧密管道化,在其中特殊的水平指令可能没有被快速路径化(但))。除非您正在运行类似SIMT的代码,否则它听起来并不是(否则您将执行四个点积)。

3
即使在Skylake处理器中,一个dpps的操作需要4个微操作指令,延迟为13个时钟周期(但吞吐量为每1.5个时钟周期执行一次)。而haddps操作需要3个微操作指令,延迟为6个时钟周期(但吞吐量为每2个时钟周期执行一次)。存储和标量操作的性能不算太差,因为它们的微操作指令数量不多,但与Kornel的回答相比,延迟较高。标量操作和向量操作具有相同的延迟。你关于“使用寄存器旁路严密流水线化”的猜测是不正确的。除了除法操作之外,所有操作都是完全流水线化的,但你说得对,水平指令并没有快速执行路径。它们被解码为内部洗牌微操作指令。 - Peter Cordes

1

通常,关于“最快的方式”的问题都预设了一个需要在时间紧迫的循环中多次完成的任务。

那么,最快的方法可能是一种成对迭代工作的迭代方法,可以在迭代之间分摊一些工作。

将向量分裂为低/高部分的减少总成本为O(log2(N)),而将向量分裂为偶数/奇数序列的分摊成本为O(1)。

inline vec update(vec context, vec data) {
    vec even = get_evens(context, data);
    vec odd = get_odds(context, data);
    return vertical_operation(even, odd);
}

void my_algo(vec *data, int N, vec_element_type *out) {

   vec4 context{0,0,0,0};
   context = update(context, data[0]);
   int i;
   for (int i = 0; i < N-1; i++) {
       context = update(context, data[i+1]);
       output[i] = extract_lane(context, 1);
   }
   context = update(context, anything);
   output[N-1] = extract_lane(context, 1);
}

所需的总和将从累加器的第二个元素(索引1)中找到(经过1次迭代后),而第一个元素将包含迄今为止所有元素的总减少量。
Reduct = [ -- ][ -- ][ -- ][ -- ]
New input = [i0 ][ i1 ][ i2 ][ i3 ]

evens = [ -- ][ -- ][ i0 ][ i2 ]
odds  = [ -- ][ -- ][ i1 ][ i3 ]
-------   vertical arithmetic reduction ----
Reduct = [ -- ][ -- ][ 01 ][ 23 ]


input = [ 4 ][ 5 ][ 6 ][ 7 ]

evens = [ -- ][ 01 ][ 4 ][ 6 ]
odds  = [ -- ][ 23 ][ 5 ][ 7 ]

Reduct = [ -- ][ 0123 ][ 45 ][ 67 ]

New input: [ 8 ] [ 9 ] [ a ] [ b ]
evens = [ -- ][ 45 ][ 8 ][ a ]
odds =  [0123][ 67 ][ 9 ][ b ]
------------------------------
Reduct = [0123][4567][ 89 ][ ab ]
        

我有疑问,如果向量长度为3或4,则此方法是否比Cordes先生介绍的更快,但对于16位或8位数据,此方法应该证明是值得的。当然,需要分别执行3或4轮才能获得结果。

如果水平操作恰好是求和 - 那么每次迭代实际上只需使用一个hadd


网页内容由stack overflow 提供, 点击上面的
可以查看英文原文,
原文链接