用于二元相关的SSE向量的Popcount？

Question

我有这个简单的二进制相关方法,它比g3的__builtin_popcount(我认为在启用SSE4时映射到popcnt指令)更好地比x3-4和%25更好的表查找和Hakmem位twiddling方法.

这是简化的代码:

int correlation(uint64_t *v1, uint64_t *v2, int size64) {
  __m128i* a = reinterpret_cast<__m128i*>(v1);
  __m128i* b = reinterpret_cast<__m128i*>(v2);
  int count = 0;
  for (int j = 0; j < size64 / 2; ++j, ++a, ++b) {
    union { __m128i s; uint64_t b[2]} x;
    x.s = _mm_xor_si128(*a, *b);
    count += _mm_popcnt_u64(x.b[0]) +_mm_popcnt_u64(x.b[1]);
  }
  return count;
}

我尝试展开循环,但我认为GCC已经自动执行此操作,因此我最终获得了相同的性能.您是否认为性能进一步提高而不会使代码过于复杂？假设v1和v2具有相同的尺寸和大小是均匀的.

我对它目前的表现感到满意,但我很好奇,看看它是否可以进一步改进.

谢谢.

编辑:修复了联合中的错误,事实证明这个错误使得这个版本比内置__builtin_popcount更快,无论如何我再次修改了代码,它现在再次比内置稍快(15%),但我认为不值得投资值得时间在此.感谢所有意见和建议.

for (int j = 0; j < size64 / 4; ++j, a+=2, b+=2) {
  __m128i x0 = _mm_xor_si128(_mm_load_si128(a), _mm_load_si128(b));
  count += _mm_popcnt_u64(_mm_extract_epi64(x0, 0))
        +_mm_popcnt_u64(_mm_extract_epi64(x0, 1));
  __m128i x1 = _mm_xor_si128(_mm_load_si128(a + 1), _mm_load_si128(b + 1));
  count += _mm_popcnt_u64(_mm_extract_epi64(x1, 0))
        +_mm_popcnt_u64(_mm_extract_epi64(x1, 1));
}

第二编辑:事实证明,内置是最快的,叹息.特别是使用-funroll-loops和-fprefetch-loop-arrays args.像这样的东西:

for (int j = 0; j < size64; ++j) {
  count += __builtin_popcountll(a[j] ^ b[j]);
}

第三编辑:

这是一个有趣的SSE3并行4位查找算法.想法来自WojciechMuła,实施是来自Marat Dukhan的回答.感谢@Apriori提醒我这个算法.下面是算法的核心,它非常聪明,基本上使用SSE寄存器作为16路查找表来计算字节的位数,使用较低的半字节作为选择表格单元格的索引.然后总结计数.

static inline __m128i hamming128(__m128i a, __m128i b) {
  static const __m128i popcount_mask = _mm_set1_epi8(0x0F);
  static const __m128i popcount_table = _mm_setr_epi8(0, 1, 1, 2, 1, 2, 2, 3, 1, 2, 2, 3, 2, 3, 3, 4);
  const __m128i x = _mm_xor_si128(a, b);
  const __m128i pcnt0 = _mm_shuffle_epi8(popcount_table, _mm_and_si128(x, popcount_mask));
  const __m128i pcnt1 = _mm_shuffle_epi8(popcount_table, _mm_and_si128(_mm_srli_epi16(x, 4), popcount_mask));
  return _mm_add_epi8(pcnt0, pcnt1);
}

在我的测试中,这个版本是平价的; 较小的输入稍微快一些,比较大的输入稍微慢一点,使用hw popcount.我认为如果它是在AVX中实现的话,它应该真的很闪耀.但是我没时间做这件事,如果有人愿意,我很乐意听到他们的结果.

Answer 1

问题是popcnt(这是__builtin_popcnt在intel CPU上编译的)对整数寄存器进行操作.这会导致编译器发出指令以在SSE和整数寄存器之间移动数据.我并不感到惊讶,因为在向量和整数寄存器之间移动数据的能力非常有限/慢,所以非sse版本更快.

uint64_t count_set_bits(const uint64_t *a, const uint64_t *b, size_t count)
{
    uint64_t sum = 0;
    for(size_t i = 0; i < count; i++) {
        sum += popcnt(a[i] ^ b[i]);
    }
    return sum;
}

这运行大约.小数据集上的每个循环2.36个时钟(适合缓存).我认为它运行缓慢是因为'long'依赖链sum限制了CPU处理更多乱序的能力.我们可以通过手动流水线循环来改进它:

uint64_t count_set_bits_2(const uint64_t *a, const uint64_t *b, size_t count)
{
    uint64_t sum = 0, sum2 = 0;

    for(size_t i = 0; i < count; i+=2) {
        sum  += popcnt(a[i  ] ^ b[i  ]);
        sum2 += popcnt(a[i+1] ^ b[i+1]);
    }
    return sum + sum2;
}

每项运行1.75个时钟.我的CPU是Sandy Bridge型号(i7-2820QM固定为@ 2.4Ghz).

四向流水线怎么样？这是每件1.65个时钟.8路怎么样？每件1.57个时钟.我们可以推导出每个项目的运行时间,(1.5n + 0.5) / n其中n是循环中管道的数量.我应该注意到,由于某种原因,当数据集增长时,8路流水线的性能比其他流量更差,我不知道为什么.生成的代码看起来没问题.

现在,如果仔细观察,每个项目都有一条xor,一条add,一条popcnt和一条mov指令.每个循环还有一个lea指令(和一个分支和减量,我忽略它因为它们几乎是免费的).

$LL3@count_set_:
; Line 50
    mov rcx, QWORD PTR [r10+rax-8]
    lea rax, QWORD PTR [rax+32]
    xor rcx, QWORD PTR [rax-40]
    popcnt  rcx, rcx
    add r9, rcx
; Line 51
    mov rcx, QWORD PTR [r10+rax-32]
    xor rcx, QWORD PTR [rax-32]
    popcnt  rcx, rcx
    add r11, rcx
; Line 52
    mov rcx, QWORD PTR [r10+rax-24]
    xor rcx, QWORD PTR [rax-24]
    popcnt  rcx, rcx
    add rbx, rcx
; Line 53
    mov rcx, QWORD PTR [r10+rax-16]
    xor rcx, QWORD PTR [rax-16]
    popcnt  rcx, rcx
    add rdi, rcx
    dec rdx
    jne SHORT $LL3@count_set_

您可以查看Agner Fog的优化手册,其中lea整个时钟周期为半个,而mov/ xor/ popcnt/ addcombo显然是1.5个时钟周期,但我并不完全理解为什么.

不幸的是,我认为我们被困在这里.该PEXTRQ指令通常用于将数据从向量寄存器移动到整数寄存器,我们可以在一个时钟周期内巧妙地拟合该指令和一条popcnt指令.添加一个整数add指令,我们的管道至少需要1.33个周期,我们仍然需要在某处添加向量加载和xor ...如果intel提供了在向量和整数寄存器之间移动多个寄存器的指令,它将是一个不同的故事.

我手头没有AVX2 cpu(xor on 256-bit vector register是AVX2的一个特性),但我的矢量化加载实现在低数据大小时执行得非常差,并且每个项目的时钟周期最少达到1.97个.

作为参考,这些是我的基准:

"管道2","管道4"和"管道8"是上面所示代码的2路,4路和8路流水线版本."sse load"的糟糕表现似乎是lzcnt/tzcnt/popcnt false依赖性bug的一种表现,gcc通过使用相同的寄存器来避免输入和输出."sse load 2"如下:

uint64_t count_set_bits_4sse_load(const uint64_t *a, const uint64_t *b, size_t count)
{
    uint64_t sum1 = 0, sum2 = 0;

    for(size_t i = 0; i < count; i+=4) {
        __m128i tmp = _mm_xor_si128(
                    _mm_load_si128(reinterpret_cast<const __m128i*>(a + i)),
                    _mm_load_si128(reinterpret_cast<const __m128i*>(b + i)));
        sum1 += popcnt(_mm_extract_epi64(tmp, 0));
        sum2 += popcnt(_mm_extract_epi64(tmp, 1));
        tmp = _mm_xor_si128(
                    _mm_load_si128(reinterpret_cast<const __m128i*>(a + i+2)),
                    _mm_load_si128(reinterpret_cast<const __m128i*>(b + i+2)));
        sum1 += popcnt(_mm_extract_epi64(tmp, 0));
        sum2 += popcnt(_mm_extract_epi64(tmp, 1));
    }

    return sum1 + sum2;
}