快速计算两个数组之间相等的字节数[重复]

Question

我编写了函数int compare_16bytes(__m128i lhs, __m128i rhs)，以便使用 SSE 指令比较两个 16 字节数：该函数在执行比较后返回有多少字节相等。

现在我想使用上面的函数来比较任意长度的两个字节数组：长度可能不是 16 字节的倍数，所以我需要处理这个问题。我怎样才能完成下面功能的实现？如何改进下面的功能？

int fast_compare(const char* s, const char* t, int length)
{
    int result = 0;

    const char* sPtr = s;
    const char* tPtr = t;

    while(...)
    {
        const __m128i* lhs = (const __m128i*)sPtr;
        const __m128i* rhs = (const __m128i*)tPtr;

        // compare the next 16 bytes of s and t
        result += compare_16bytes(*lhs,*rhs);

        sPtr += 16;
        tPtr += 16;
    }

    return result;
}

Answer 1

SSE 中的整数比较产生全零或全一的字节。如果要计数，首先需要将比较结果右移（不是算术）7，然后添加到结果向量中。 最后，您仍然需要通过对其元素求和来减少结果向量。这种减少必须在标量代码中完成，或者通过一系列添加/移位来完成。通常这部分不值得麻烦。

Answer 2

正如@Mysticial 在上面的 cmets 中所说，进行比较和垂直求和，然后在主循环结束时水平求和：

#include <stdio.h>
#include <stdlib.h>
#include <time.h>
#include <emmintrin.h>

// reference implementation
int fast_compare_ref(const char *s, const char *t, int length)
{
    int result = 0;
    int i;

    for (i = 0; i < length; ++i)
    {
        if (s[i] == t[i])
            result++;
    }
    return result;
}

// optimised implementation
int fast_compare(const char *s, const char *t, int length)
{
    int result = 0;
    int i;

    __m128i vsum = _mm_set1_epi32(0);
    for (i = 0; i < length - 15; i += 16)
    {
        __m128i vs, vt, v, vh, vl, vtemp;

        vs = _mm_loadu_si128((__m128i *)&s[i]); // load 16 chars from input
        vt = _mm_loadu_si128((__m128i *)&t[i]);
        v = _mm_cmpeq_epi8(vs, vt);             // compare
        vh = _mm_unpackhi_epi8(v, v);           // unpack compare result into 2 x 8 x 16 bit vectors
        vl = _mm_unpacklo_epi8(v, v);
        vtemp = _mm_madd_epi16(vh, vh);         // accumulate 16 bit vectors into 4 x 32 bit partial sums
        vsum = _mm_add_epi32(vsum, vtemp);
        vtemp = _mm_madd_epi16(vl, vl);
        vsum = _mm_add_epi32(vsum, vtemp);
    }

    // get sum of 4 x 32 bit partial sums
    vsum = _mm_add_epi32(vsum, _mm_srli_si128(vsum, 8));
    vsum = _mm_add_epi32(vsum, _mm_srli_si128(vsum, 4));
    result = _mm_cvtsi128_si32(vsum);

    // handle any residual bytes ( < 16)
    if (i < length)
    {
        result += fast_compare_ref(&s[i], &t[i], length - i);
    }

    return result;
}

// test harness
int main(void)
{
    const int n = 1000000;
    char *s = malloc(n);
    char *t = malloc(n);
    int i, result_ref, result;

    srand(time(NULL));

    for (i = 0; i < n; ++i)
    {
        s[i] = rand();
        t[i] = rand();
    }

    result_ref = fast_compare_ref(s, t, n);
    result = fast_compare(s, t, n);

    printf("result_ref = %d, result = %d\n", result_ref, result);;

    return 0;
}

编译并运行上述测试工具：

$ gcc -Wall -O3 -msse3 fast_compare.c -o fast_compare
$ ./fast_compare
result_ref = 3955, result = 3955
$ ./fast_compare
result_ref = 3947, result = 3947
$ ./fast_compare
result_ref = 3945, result = 3945

---

请注意，在上述 SSE 代码中可能存在一个不明显的技巧，我们使用 _mm_madd_epi16 将 16 位 0/-1 值解包并累加为 32 位部分和。我们利用-1*-1 = 1（当然还有0*0 = 0）这一事实——我们并没有真正在这里进行乘法运算，只是在一条指令中解包和求和。

---

更新：如下面的 cmets 所述，此解决方案不是最佳解决方案 - 我只是采用了一个相当最佳的 16 位解决方案，并将 8 位添加到 16 位解包以使其适用于 8 位数据。但是对于 8 位数据，有更有效的方法，例如使用psadbw/_mm_sad_epu8。我将把这个答案留给后人，以及任何可能想用 16 位数据做这种事情的人，但实际上不需要解压输入数据的其他答案之一应该是公认的答案。

Answer 3

在 16 x uint8 元素中使用部分和可能会提供更好的性能。
我将循环分为内循环和外循环。
内部循环求和 uint8 元素（每个 uint8 元素最多可以求和 255 个“1”）。
小技巧：_mm_cmpeq_epi8 将相等的元素设置为 0xFF，并且 (char)0xFF = -1，因此您可以从总和中减去结果（减去 -1 为加 1）。

这是我对 fast_compare 的优化版本：

int fast_compare2(const char *s, const char *t, int length)
{
    int result = 0;
    int inner_length = length;
    int i;
    int j = 0;

    //Points beginning of 4080 elements block.
    const char *s0 = s;
    const char *t0 = t;


    __m128i vsum = _mm_setzero_si128();

    //Outer loop sum result of 4080 sums.
    for (i = 0; i < length; i += 4080)
    {
        __m128i vsum_uint8 = _mm_setzero_si128(); //16 uint8 sum elements (each uint8 element can sum up to 255).
        __m128i vh, vl, vhl, vhl_lo, vhl_hi;

        //Points beginning of 4080 elements block.
        s0 = s + i;
        t0 = t + i;

        if (i + 4080 <= length)
        {
            inner_length = 4080;
        }
        else
        {
            inner_length = length - i;
        }

        //Inner loop - sum up to 4080 (compared) results.
        //Each uint8 element can sum up to 255. 16 uint8 elements can sum up to 255*16 = 4080 (compared) results.
        //////////////////////////////////////////////////////////////////////////
        for (j = 0; j < inner_length-15; j += 16)
        {
              __m128i vs, vt, v;

              vs = _mm_loadu_si128((__m128i *)&s0[j]); // load 16 chars from input
              vt = _mm_loadu_si128((__m128i *)&t0[j]);
              v = _mm_cmpeq_epi8(vs, vt);             // compare - set to 0xFF where equal, and 0 otherwise.

              //Consider this: (char)0xFF = (-1)
              vsum_uint8 = _mm_sub_epi8(vsum_uint8, v); //Subtract the comparison result - subtract (-1) where equal.
        }
        //////////////////////////////////////////////////////////////////////////

        vh = _mm_unpackhi_epi8(vsum_uint8, _mm_setzero_si128());        // unpack result into 2 x 8 x 16 bit vectors
        vl = _mm_unpacklo_epi8(vsum_uint8, _mm_setzero_si128());
        vhl = _mm_add_epi16(vh, vl);    //Sum high and low as uint16 elements.

        vhl_hi = _mm_unpackhi_epi16(vhl, _mm_setzero_si128());   //unpack sum of vh an vl into 2 x 4 x 32 bit vectors
        vhl_lo = _mm_unpacklo_epi16(vhl, _mm_setzero_si128());   //unpack sum of vh an vl into 2 x 4 x 32 bit vectors

        vsum = _mm_add_epi32(vsum, vhl_hi);
        vsum = _mm_add_epi32(vsum, vhl_lo);
    }

    // get sum of 4 x 32 bit partial sums
    vsum = _mm_add_epi32(vsum, _mm_srli_si128(vsum, 8));
    vsum = _mm_add_epi32(vsum, _mm_srli_si128(vsum, 4));
    result = _mm_cvtsi128_si32(vsum);

    // handle any residual bytes ( < 16)
    if (j < inner_length)
    {
        result += fast_compare_ref(&s0[j], &t0[j], inner_length - j);
    }

    return result;
}

Answer 4

大输入的最快方法是 Rotem 的答案，其中内部循环是 pcmpeqb / psubb，在向量累加器的任何字节元素溢出之前突破水平求和。使用psadbw 对全零向量进行无符号字节的hsum。

另请参阅How to count character occurrences using SIMD，您可以在其中将 C++ 与 AVX2 的内在函数一起使用，以使用从另一个数组加载的向量而不是该问题的 _mm_set1_epi8(char_to_count) 来计算匹配项。有效地将比较结果相加是相同的，使用 psadbw 进行水平总和。

---

没有展开/嵌套循环，最好的选择可能是

pcmpeqb   -> vector of  0  or  0xFF  elements
psadbw    -> two 64bit sums of  (0*no_matches + 0xFF*matches)
paddq     -> accumulate the psadbw result in a vector accumulator

#outside the loop:
horizontal sum
divide the result by 255

如果您的循环中没有太大的寄存器压力，请psadbw 反对0x7f 的向量而不是全零。

• psadbw(0x00, set1(0x7f)) => sum += 0x7f
• psadbw(0xff, set1(0x7f)) => sum += 0x80

因此，您不必除以 255（编译器应该在没有实际 div 的情况下有效地执行此操作），而只需减去 n * 0x7f，其中 n 是元素的数量。

另请注意，paddq 在 Nehalem 之前和 Atom 上运行缓慢，因此如果您不希望 128 * 计数会溢出 32 位整数，则可以使用 paddd (_mm_add_epi32)。

这与 Paul R 的 pcmpeqb / 2x punpck / 2x pmaddwd / 2x paddw 相比非常好。

---

但通过小幅展开，您可以在 psadbw / paddq 之前累积 4 或 8 个与 psubb 的比较结果。