我正在学习如何使用SIMD内在函数和自动向量化。幸运的是,我有一个有用的项目,我正在努力,似乎非常适合SIMD,但对于像我这样的新手来说仍然很棘手。
我正在为图像编写一个过滤器来计算2x2像素的平均值。我正在通过将两个像素的总和累加到单个像素中来进行计算的一部分。
template <typename T, typename U>
inline void accumulate_2x2_x_pass(
T* channel, U* accum,
const size_t sx, const size_t sy,
const size_t osx, const size_t osy,
const size_t yoff, const size_t oyoff
) {
const bool odd_x = (sx & 0x01);
size_t i_idx, o_idx;
// Should be vectorizable somehow...
for (size_t x = 0, ox = 0; x < sx - (size_t)odd_x; x += 2, ox++) {
i_idx = x + yoff;
o_idx = ox + oyoff;
accum[o_idx] += channel[i_idx];
accum[o_idx] += channel[i_idx + 1];
}
if (odd_x) {
// << 1 bc we need to multiply by two on the edge
// to avoid darkening during render
accum[(osx - 1) + oyoff] += (U)(channel[(sx - 1) + yoff]) * 2;
}
}
然而,godbolt表明我的循环不能自动发展。 (https://godbolt.org/z/qZxvof)我如何构建SIMD内在函数来解决这个问题?我控制了accum的对齐,但没有控制通道。
(我知道有一个平均内在的,但这里不合适,因为我需要生成多个mip级别,该命令会导致下一级别的精度损失。)
感谢大家。 :)
狭窄类型T
= uint8_t
或uint16_t
的扩大情况可能最好用SSSE3 pmaddubsw
或SSE2 pmaddwd
以及1
的乘数来实现。 (Intrinsics guide)那些指令单一uop并完全做水平扩展,你需要比改组更有效率。
如果可以在不损失精度的情况下执行此操作,请先在行之间添加垂直添加,然后再加宽水平添加。 (例如,[u]int16_t
中的10,12或14位像素分量不能溢出)。在大多数CPU上,负载和垂直添加(每个时钟吞吐量至少为2),而pmadd*
每个时钟只有1个时钟,Skylake及更高版本的吞吐量只有2个。这意味着你只需要1x加+ 1x pmadd对2x pmadd + 1x加,所以即使在Skylake上它也是一个重要的胜利。 (对于第二种方式,如果你有AVX,两个加载都可以折叠到pmadd的内存操作数中。对于pmadd方式之前的加法,你首先需要一个纯加载,然后将第二个加载折叠成add,这样你就可能无法保存前端uops,除非你使用索引寻址模式并且它们没有层压。)
理想情况下,您不需要将+=
转换为累加器数组,而是可以并行读取2行而累加器是只写的,因此您的循环只有2个输入流和1个输出流。
// SSSE3
__m128i hadd_widen8_to_16(__m128i a) {
// uint8_t, int8_t (doesn't matter when multiplier is +1)
return _mm_maddubs_epi16(a, _mm_set_epi8(1));
}
// SSE2
__m128i hadd_widen16_to_32(__m128i a) {
// int16_t, int16_t
return _mm_madd_epi16(a, _mm_set_epi16(1));
}
这些端口直接连接到256位AVX2,因为输入和输出宽度是相同的。修理车道内包装不需要洗牌。
是的,他们都是_epi16
。英特尔可能与内在名称大相径庭。 asm助记符更一致,更容易记住什么是什么。 (ubsw
=无符号字节到有符号字,除了其中一个输入是有符号字节.pmaddwd
是打包乘法加字到dword,与punpcklwd
相同的命名方案等)
使用uint16_t
或uint32_t
的T = U情况是SSSE3 _mm_hadd_epi16
或_mm_hadd_epi32
的用例。它的成本与2 shuffle +一个垂直添加相同,但无论如何你需要将2个输入打包为1。
如果你想在Haswell及其后的某个洗牌端瓶颈上工作,你可以考虑在输入上使用qword shift,然后用shufps
(_mm_shuffle_ps
+ some casting)将结果混合在一起。这可能是Skylake的胜利(每个时钟移位吞吐量为2),即使它总共花费5个uop而不是3个。它可以在每个矢量输出最多运行5/3个周期,而不是每个矢量2个周期,如果有的话没有前端瓶颈
// UNTESTED
//Only any good with AVX, otherwise the extra movdqa instructions kill this
//Only worth considering for Skylake, not Haswell (1/c shifts) or Sandybridge (2/c shuffle)
__m128i hadd32_emulated(__m128i a, __m128i b) {
__m128i a_shift = _mm_srli_epi64(a, 32);
__m128i b_shift = _mm_srli_epi64(b, 32);
a = _mm_add_epi32(a, a_shift);
b = _mm_add_epi32(b, b_shift);
__m128 combined = _mm_shuffle_ps(_mm_castsi128_ps(a), _mm_castsi128_ps(b), _MM_SHUFFLE(2,0,2,0));
return _mm_castps_si128(combined);
}
对于AVX2版本,你需要一个交叉的shuffle来修复vphadd
结果。因此,效率变化可能会带来更大的胜利。
// 3x shuffle 1x add uops
__m256i hadd32_avx2(__m256i a, __m256i b) {
__m256i hadd = _mm256_hadd_epi32(a, b); // 2x in-lane hadd
return _mm256_permutex_epi64( hadd, _MM_SHUFFLE(3,1,2,0) );
}
// UNTESTED
// 2x shift, 2x add, 1x blend-immediate (any ALU port), 1x shuffle
__m256i hadd32_emulated_avx2(__m256i a, __m256i b)
{
__m256i a_shift = _mm256_srli_epi64(a, 32); // useful result in the low half of each qword
__m256i b_shift = _mm256_slli_epi64(b, 32); // ... high half of each qword
a = _mm256_add_epi32(a, a_shift);
b = _mm256_add_epi32(b, b_shift);
__m256i blended = _mm256_blend_epi32(a,b, 0b10101010); // alternating low/high results
return _mm256_permutexvar_epi32(_mm256_set_epi32(7,5,3,1, 6,4,2,0), blended);
}
在Haswell和Skylake上,hadd32_emulated_avx2
可以每2个时钟运行一次(使所有矢量ALU端口饱和)。将add_epi32
加到accum[]
中的额外hadd32_avx2
将减慢到每256位结果向量的最佳7/3周期,并且您需要展开(或使用展开的编译器)而不仅仅是前端的瓶颈。
add_epi32
可以每3个时钟运行一次(在端口5上用于洗牌)。加载+存储+额外的https://agner.org/optimize/ uops来实现你的循环可以轻松地在阴影中运行。
(https://stackoverflow.com/tags/x86/info,见qazxswpoi)