计算两个 _m128i SIMD 向量之间的匹配字节数

Posted 2023-02-16

【中文标题】计算两个 _m128i SIMD 向量之间的匹配字节数

【英文标题】：Count number of matching bytes between two _m128i SIMD vectors

【发布时间】：2021-04-26 22:43:44

【问题描述】：

我正在开发一种生物信息学工具，并尝试使用 SIMD 来提高其速度。

给定两个长度为 16 的 char 数组，我需要快速计算字符串匹配的索引数。例如，以下两个字符串“TTTTTTTTTTTTTTTT”和“AAAAGGGGTTTTCCCC”从第 9 位到第 12 位（“TTTT”）匹配，因此输出应为 4。

如以下函数 foo 所示（工作正常但速度慢），我将 seq1 和 seq2 中的每个字符打包到 __m128i 变量 s1 和 s2 中，并使用 _mm_cmpeq_epi8 来同时比较每个位置。然后，使用 popcnt128（来自 Marat Dukhan 的 Fast counting the number of set bits in __m128i register）将匹配位数相加。

float foo(char* seq1, char* seq2) 
    __m128i s1, s2, ceq;
    int match;
    s1 =  _mm_load_si128((__m128i*)(seq1));
    s2 =  _mm_load_si128((__m128i*)(seq2));
    ceq = _mm_cmpeq_epi8(s1, s2);
    match = (popcnt128(ceq)/8);
    return match;

虽然 Marat Dukhan 的 popcnt128 比在 __m128i 中简单地累加每一位要快得多，但 __popcnt128() 是函数中最慢的瓶颈，占用了大约 80% 的计算速度。所以，我想提出一个 popcnt128 的替代方案。

---

我试图将__m128i ceq 解释为字符串，并将其用作预先计算的查找表的键，该查找表将字符串映射到总位数。如果 char 数组是可散列的，我可以做类似的事情

union__m128i ceq; char c_arr[16];
match = table[c_arr] // table = unordered map

如果我尝试对字符串（即union__m128i ceq; string s;;）执行类似的操作，我会收到以下错误消息“::()' 被隐式删除，因为默认定义格式不正确”。当我尝试其他事情时，我遇到了分段错误。

有什么方法可以告诉编译器将 __m128i 读取为字符串，以便我可以直接使用 __m128i 作为 unordered_map 的键？我不明白为什么它不应该工作，因为字符串是一个连续的字符数组，可以自然地用 __m128i 表示。但我无法让它工作，也无法在线找到任何解决方案。

【问题讨论】：

我认为使用表格不是一个好主意，因为我猜有 2^128 种可能性，不是吗？ （远远超过您的 RAM 可以容纳的容量）。此外，散列可能会比当前的 popcnt128 函数慢。除了访问 c_arr 之外，您生成格式错误的程序的方式会导致 C++ 中未定义的行为（如果需要，请改用 memcpy）。

popcnt128 所做的远远超出您的实际需要。只需__builtin_popcnt(_mm_movemask_epi8(cnt))。

@chtz: 或者psadbw 反对 0（比较结果的字节hsum），但这需要qword 一半的最后hsum 步骤，所以如果硬件 popcnt 可用，情况会更糟，但如果您需要仅使用 SSE2 的基线 x86-64，则值得考虑。 （另外你需要在 psadbw 之前取反，或者将总和缩小 1/255）。

@JérômeRichard：OP 提出的是std::unordered_map<char[16]> 或其他东西，而不是平面数组，所以它是可能的； “仅”需要 2^16 个条目。但与 movemask / popcnt 相比，当然慢得可怕。计算 16 字节 char 数组的哈希函数所需的工作量与获得答案相当。 :P

@chtz 我试过 __builtin_popcnt(_mm_movemask_epi8(cnt)) ，现在速度快了很多。

【参考方案1】：

您可能正在为更长的序列、多个 SIMD 数据向量执行此操作。在这种情况下，您可以在一个向量中累积计数，您只能在最后总结。 单独计算每个向量的效率要低得多。

参见How to count character occurrences using SIMD - 代替_mm256_set1_epi8(c); 来搜索特定字符，从另一个字符串加载。做其他所有事情，包括counts = _mm_sub_epi8(counts, _mm_cmpeq_epi8(s1, s2)); 在内部循环中，循环展开。 （比较结果是整数 0 / -1，因此减去它会将 0 或 1 添加到另一个向量。）这在 256 次迭代后有溢出的风险，因此最多 255 次。该链接问题使用 AVX2，但 @987654325 @ 这些内在函数的版本只需要 SSE2。 （当然，AVX2 可以让每条向量指令完成两倍的工作量。）

使用_mm_sad_epu8(v, _mm_setzero_si128()); 对外部循环中的字节计数器进行水平求和，然后累加到另一个计数向量中。 再一次，这全部都在链接问答中的代码中，所以只需复制/粘贴它并将另一个字符串的负载添加到内部循环中，而不是使用广播常量。

---

对于单个向量：

您不需要计算__m128i 中的所有位；通过将每个元素的 1 位提取为标量整数，利用每个字节中的所有 8 位都相同的事实。 （与其他一些 SIMD ISA 不同，x86 SIMD 可以有效地做到这一点）

    count = __builtin_popcnt(_mm_movemask_epi8(cmp_result));

另一个可能选项是psadbw反对0（比较结果的字节hsum），但这需要qword一半的最后一个hsum步骤，所以这将比硬件popcnt更糟糕。但是，如果您无法使用-mpopcnt 进行编译，那么您是否需要仅使用 SSE2 的基线 x86-64 是值得考虑的。 （另外你需要在 psadbw 之前取反，或者将总和缩小 1/255...）

（请注意，psadbw 策略基本上是我在答案的第一部分中描述的，但仅适用于单个向量，没有利用廉价地将多个计数添加到一个向量累加器中的能力。）

如果您确实需要 float 形式的结果，那么 psadbw 策略就不那么糟糕了：您可以始终将值保留在 SIMD 向量中，使用 _mm_cvtepi32_ps 对水平和进行打包转换结果（甚至比cvtsi2ss int->float 标量转换便宜）。 _mm_cvtps_f32 是免费的；标量浮点数只是 XMM 寄存器的低位元素。

但是说真的，你真的需要一个整数作为float现在吗？你不能至少等到你得到所有向量的总和，还是保持整数？

-mpopcnt 由gcc -msse4.2 或-march=native 隐含在小于10 岁的任何东西上。 Core 2 缺少硬件 popcnt，但 Nehalem 为 Intel 提供了它。

【讨论】：

谢谢！这很有帮助。我不需要整数计数作为浮点数。我有点草率，将其保留为整数确实提高了一点性能。

@JWO：好的。使用 sub(count, cmp(load1,load2)) 的内部循环应该有更多帮助，例如3 条指令 (movdqa / pcmpeqb xmm,mem / psubb 从 5 (movdqa / pcmeqb xmm, mem / pmovmskb / popcnt / add) 减少，加上循环开销。特别是如果你的字符串在 L1d 甚至 L2 缓存中很热；如果不是那么无论哪种方式都可能成为内存带宽的瓶颈，特别是如果您可以使用 AVX2。

新的第三代英特尔可扩展处理器（Ice Lake Xeon）支持 VPOPCNT 系列指令。这些返回为宽度为 128/256/512 位的向量中 BYTE/WORD/DWORD/QWORD 元素的所有组合的每个元素设置的位数。我还没有看到延迟/吞吐量数据，但这些肯定会比软件方法更快。

@JohnDMcCalpin：但是这个问题不需要甚至真的想要计算SIMD向量中的所有位；相反，他们想计算匹配的字节数。为此，pcmpeqb / psubb 是多向量的最佳选择；我看不出用 SIMD popcnt 做得更好的方法。也不是单个向量的标量； pmovmskb / popcnt 显然比 vpopcntq / shuffle / vpaddd / vmovd / shr eax,3 好，因为你在每个字节中计算了 8 位。

@JohnDMcCalpin：IceLake-client 也有这个扩展，uops.info 有性能数据。端口 5 的单个 uop，3 周期延迟，可以微熔断内存源操作数。 （或使用字节或字掩码的 5 周期延迟，否则所有元素大小的性能相同。）