使用 AVX 但不使用 AVX2 分别计算多个 64 位位掩码上的每个位位置

2024-07-03

（有关的：如何在 Sandy Bridge 上将一系列整数中的位快速计数到单独的容器中？ https://stackoverflow.com/questions/7793997/how-to-quickly-count-bits-into-separate-bins-in-a-series-of-ints-on-sandy-bridge是这个的早期副本，有一些不同的答案。编者注：这里的答案可能更好。

此外，类似问题的 AVX2 版本，有一整行位的许多箱比一个宽得多uint64_t: 改进列总体计数算法 https://stackoverflow.com/questions/58486138/improve-column-population-count-algorithm)

我正在做一个 C 项目，我需要遍历数千万个掩码（ulong 类型（64 位））并更新一个数组（称为target）基于简单规则的 64 个短整数 (uint16)：

// for any given mask, do the following loop
for (i = 0; i < 64; i++) {
    if (mask & (1ull << i)) {
        target[i]++
    }
}

问题是我需要在数千万个掩模上执行上述循环，并且我需要在不到一秒的时间内完成。想知道是否有任何方法可以加快速度，例如使用某种代表上述循环的特殊汇编指令。

目前我在ubuntu 14.04（i7-2670QM，支持AVX，不支持AVX2）上使用gcc 4.8.4编译并运行以下代码，大约花费了2秒。希望使其运行时间低于 200 毫秒。

#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <sys/time.h>
#include <sys/stat.h>

double getTS() {
    struct timeval tv;
    gettimeofday(&tv, NULL);
    return tv.tv_sec + tv.tv_usec / 1000000.0;
}
unsigned int target[64];

int main(int argc, char *argv[]) {
    int i, j;
    unsigned long x = 123;
    unsigned long m = 1;
    char *p = malloc(8 * 10000000);
    if (!p) {
        printf("failed to allocate\n");
        exit(0);
    }
    memset(p, 0xff, 80000000);
    printf("p=%p\n", p);
    unsigned long *pLong = (unsigned long*)p;
    double start = getTS();
    for (j = 0; j < 10000000; j++) {
        m = 1;
        for (i = 0; i < 64; i++) {
            if ((pLong[j] & m) == m) {
                target[i]++;
            }
            m = (m << 1);
        }
    }
    printf("took %f secs\n", getTS() - start);
    return 0;
}

提前致谢！

有关的：

早期的副本有一些替代想法：如何在 Sandy Bridge 上将一系列整数中的位快速计数到单独的容器中？ https://stackoverflow.com/questions/7793997/how-to-quickly-count-bits-into-separate-bins-in-a-series-of-ints-on-sandy-bridge.
哈罗德的回答AVX2 列总体计数算法分别针对每个位列 https://stackoverflow.com/questions/58486138/improve-column-population-count-algorithm.
矩阵转置和总体计数 https://stackoverflow.com/questions/51475704/matrix-transpose-and-population-count有一些关于 AVX2 的有用答案，包括基准测试。它使用 32 位块而不是 64 位。

Also: https://github.com/mklarqvist/positional-popcount https://github.com/mklarqvist/positional-popcount具有 SSE 混合、各种 AVX2、各种 AVX512（包括非常适合大型数组的 Harley-Seal）以及用于位置弹出计数的各种其他算法。可能只为了uint16_t，但大多数都可以适应其他字宽。我认为我下面提出的算法就是他们所说的adder_forest.

最好的选择是 SIMD，在 Sandybridge CPU 上使用 AVX1。编译器不够智能，无法为您自动矢量化循环位，即使您无分支地编写它以给它们更好的机会。

不幸的是，它不够智能，无法自动矢量化逐渐加宽和添加的快速版本。

See intel avx2 中是否有 movemask 指令的逆指令？ https://stackoverflow.com/q/36488675了解不同大小的位图 -> 矢量解包方法的摘要。 Ext3h 在另一个答案中的建议很好：将位解压缩为比最终计数数组更窄的内容，可以为每条指令提供更多元素。使用 SIMD 字节效率很高，然后您最多可以执行 255 个垂直paddb没有溢出，在解包之前累加到 32 位计数器数组中。

只需要 4x 16 字节__m128i保存所有 64 个向量uint8_t元素，因此这些累加器可以保留在寄存器中，仅在外循环中扩展到 32 位计数器时才添加到内存中。

解包不必按顺序进行: 你可以随时随机播放target[]在最后一次，在累积所有结果之后。

内部循环可以展开以从 64 或 128 位向量加载开始，并使用 4 或 8 种不同的方式解包pshufb (_mm_shuffle_epi8).

更好的策略是逐步扩大

从 2 位累加器开始，然后进行掩码/移位以将其扩展为 4 位。因此，在最内部的循环中，大多数操作都在处理“密集”数据，而不是立即“稀释”太多数据。更高的信息/熵密度意味着每条指令做了更多有用的工作。

Using SWAR https://en.wikipedia.org/wiki/SWAR在标量或 SIMD 寄存器内进行 32x 2 位加法的技术很简单/便宜，因为无论如何我们都需要避免执行元素顶部的可能性。使用正确的 SIMD，我们会丢失这些计数，而使用 SWAR，我们会破坏下一个元素。

uint64_t x = *(input++);        // load a new bitmask
const uint64_t even_1bits = 0x5555555555555555;  // 0b...01010101;

uint64_t lo = x & even_1bits;
uint64_t hi = (x>>1) & even_1bits;            // or use ANDN before shifting to avoid a MOV copy

accum2_lo += lo;   // can do up to 3 iterations of this without overflow
accum2_hi += hi;   // because a 2-bit integer overflows at 4

然后你重复最多 4 个 4 位元素的向量，然后是 8 个 8 位元素的向量，然后你应该一直扩大到 32 并累积到内存中的数组中，因为无论如何你都会用完寄存器，这外层循环工作并不频繁，因此我们不需要费心去使用 16 位。（特别是如果我们手动矢量化）。

最大的缺点：这个doesn't自动矢量化，与@njuffa 的版本不同。但与gcc -O3 -march=sandybridge对于AVX1（然后在Skylake上运行代码），这个运行标量64位实际上仍然稍微faster比 @njuffa 代码中的 128 位 AVX 自动矢量化 asm 更好。

但这是 Skylake 上的时序，它有 4 个标量 ALU 端口（和 mov-elimination），而 Sandybridge 缺乏 mov-elimination 并且只有 3 个 ALU 端口，因此标量代码可能会遇到后端执行端口瓶颈。（但是 SIMD 代码可能几乎一样快，因为有大量的 AND / ADD 与移位混合在一起，并且 SnB 的所有 3 个端口上确实都有 SIMD 执行单元，这些端口上有 ALU。Haswell 刚刚添加了端口 6，用于标量-仅包括轮班和分行。）

通过良好的手动矢量化，速度应该提高近 2 或 4 倍。

但是，如果您必须在这个标量或带有 AVX2 自动矢量化的 @njuffa 之间进行选择，@njuffa 在 Skylake 上速度更快-march=native

如果可以/需要在 32 位目标上构建，这会受到很大影响（由于在 32 位寄存器中使用 uint64_t，因此没有矢量化），而矢量化代码几乎不会受到任何影响（因为所有工作都发生在相同的向量寄存器中）宽度）。

// TODO: put the target[] re-ordering somewhere
// TODO: cleanup for N not a multiple of 3*4*21 = 252
// TODO: manual vectorize with __m128i, __m256i, and/or __m512i

void sum_gradual_widen (const uint64_t *restrict input, unsigned int *restrict target, size_t length)
{
    const uint64_t *endp = input + length - 3*4*21;     // 252 masks per outer iteration
    while(input <= endp) {
        uint64_t accum8[8] = {0};     // 8-bit accumulators
        for (int k=0 ; k<21 ; k++) {
            uint64_t accum4[4] = {0};  // 4-bit accumulators can hold counts up to 15.  We use 4*3=12
            for(int j=0 ; j<4 ; j++){
                uint64_t accum2_lo=0, accum2_hi=0;
                for(int i=0 ; i<3 ; i++) {  // the compiler should fully unroll this
                    uint64_t x = *input++;    // load a new bitmask
                    const uint64_t even_1bits = 0x5555555555555555;
                    uint64_t lo = x & even_1bits; // 0b...01010101;
                    uint64_t hi = (x>>1) & even_1bits;  // or use ANDN before shifting to avoid a MOV copy
                    accum2_lo += lo;
                    accum2_hi += hi;   // can do up to 3 iterations of this without overflow
                }

                const uint64_t even_2bits = 0x3333333333333333;
                accum4[0] +=  accum2_lo       & even_2bits;  // 0b...001100110011;   // same constant 4 times, because we shift *first*
                accum4[1] += (accum2_lo >> 2) & even_2bits;
                accum4[2] +=  accum2_hi       & even_2bits;
                accum4[3] += (accum2_hi >> 2) & even_2bits;
            }
            for (int i = 0 ; i<4 ; i++) {
                accum8[i*2 + 0] +=   accum4[i] & 0x0f0f0f0f0f0f0f0f;
                accum8[i*2 + 1] +=  (accum4[i] >> 4) & 0x0f0f0f0f0f0f0f0f;
            }
        }

        // char* can safely alias anything.
        unsigned char *narrow = (uint8_t*) accum8;
        for (int i=0 ; i<64 ; i++){
            target[i] += narrow[i];
        }
    }
    /* target[0] = bit 0
     * target[1] = bit 8
     * ...
     * target[8] = bit 1
     * target[9] = bit 9
     * ...
     */
    // TODO: 8x8 transpose
}

我们不关心顺序，所以accum4[0]例如，每第 4 位就有 4 位累加器。最后需要的最终修复（但尚未实现）是 8x8 转置uint32_t target[64] array,这可以使用 unpck 和vshufps仅使用 AVX1。（使用 AVX/AVX2 转置 8x8 浮点 https://stackoverflow.com/questions/25622745/transpose-an-8x8-float-using-avx-avx2）。还有最后 251 个掩模的清理循环。

我们可以使用任何 SIMD 元素宽度来实现这些移位；无论如何，我们必须屏蔽低于 16 位的宽度（SSE/AVX 没有字节粒度移位，最小只有 16 位。）

Arch Linux i7-6700k 的基准测试结果来自 @njuffa 的测试工具，添加了这个。（Godbolt https://godbolt.org/#z:OYLghAFBqd5QCxAYwPYBMCmBRdBLAF1QCcAaPECAKxAEZSAbAQwDtRkBSAJgCFufSAZ1QBXYskwgA5AHoZAagQECAB0Eg5ggk2QBrVADdMxAGYNUAdwB0aALYyAjiMxa8qFoJkBWLwAYAHLReXF4yaCIsBAC0mDoIUQBGhFEqqIKEbixRgpgqTMRMBJgMAJ5RhsZRtqxlAGwALInJSQTVgrqCURaE8UwGAB6JItEsqNFM3n6BtVwA7NwAzD4BtDPzvgCCi3gsyAwiWPIcCwDCWvioVgjH2Byb27v7h8dnBOgMeAlXN3dbXAs7PYHTBHU7nHYEb4LW73f6Ap4gl5aYg7YBQmFbWEAkzyLAmHaYdAQAD6AHUAJIAOQWXAAlL9tjiwGA8QSiRTqVxiQAZbAbSnE/kAEWJAFk%2BZT6VjWSwQRyaTyJYLKSLxfyGf9MCx8CYNQDHsDQSdutrLIJ0b8dh9Zb90KIEgwQTk0Nr5BADKg8OgpVtZnxNvJA/JuRsAEoAcWwxKpABVsJHQ/ICMd/Rsg/ItIU8MhcfbHfJUKgTMRMA4U7905mCNn5BD5MgEJg9ISAGIkAASeGACFDLhjeFsxnLAaDVZrPAA8hPuYomIJO93e4J%2B4PiMPMWmg3gcRBmQ2m7pWx2uz2%2BwPjLSjn6K%2Bmgwg5wvT8vz8RQUL5ABFZzEEoABWMJgkNUuyYC2JZOFqyAlG63C1AQ9ILKmt5BoWxalm%2B8i0FYvjyAoEB2iIDqYLSkJfkw6C/vkyaITeyH7s26BtsQj5LiuxgYbQ67phwsxCrRtY7ve84nqxL6XjxSHIZ%2B35/gBQGsBIJyiJE7EQLB8FcVJJYEGILBugRREkVYZEUVR8gAFQFkW4GaYGPHvsUORXpJyHabp%2Bl5sR4aYAQ/Z6EpEQEBAl6WVhviYFECy2VefH3LxerFNuuKYPispEsSxLWiI/QZeJswnDxJzJalhIkoKv6/ryuV6vChpIm8eIWlitXPGCJSeNWg5NRsVoErankZk27joG6Hpej6En8ciIjIAQSbngYTAMEmBjRcAPmdZgRboEwJRqVwcEGKQ8iUgAqty3IIS58hucQen4Z5JEGFYBAGMSzpHLwHmEY6T0vW9IgfaF2ERbU672QlDA5AlxDECQ8gRIIIgqKkxBFCNKjMAQgHELYCXatuDKbDI5k3mBqC2CA8jragRrANUxDmlYVgAAL00weAMDYFM/COlKYBYgjAHDyPqBmyB4FY1QEAgpB2CoVjMGwNhyxTCtK2iOMEAULA3gAyoRVBNgQVMNE0c0ANLknrevyKGlLhoIN5CoUkjyHrhTHVw/jyC2mAJPIXC%2BL4ACc8i%2BPUIALL4ID1D7UQBMHbq/nrMYTSOMYIHggjyGbLTyFbNt2w7s45zLJYgvL7hagQgjHbEDbyLK%2BSlDewCFiNWqiN2Sa0zkxBGPIS3V1Y6aZ5XavV5EOf5JIN6iiIDDVpjZSkj0UQnPkP5uqKpInLSx0qMYbgjRAXAcF42C0FwnG8Ofl%2B1AsUS0D6gYABokIIWcmHNEBv3r7YD7yCPiiDAgcL7YDNpxEcSk2DEGcJEPAS03QnAdkAkBJ9wEP3qL8cyMhfhjhzCICEDRiRzV0NnQQpD%2BgcX%2BPULwtRZj%2BBDqHfwsxaheHqDSWg51uTrkIQjEh9QyHyAoYIKhwiczHHfNfBY9DaiyPkYohhvD%2BHaGrEQoRIixESOJNBaR8gFgzC4Y/YxRijFcBMaomimwBHEMiKQ8hlDSEAC8OK%2BFqAo%2BoIcb6zHYbQQIvhrGpjsVopx4jHHg3%2BDKEEu8TgNEDBAHRjiMJJOccImhLwXjyC8P4cSX1kmSOOhfAqI4pJSUKcSKRCx3xpIiRk0EtxoS5yAZUzJvADGVIIMUrwpTNzlOQpUqRHSalujaUaUR6SyGtKmf0Ca0SUoEnkP/dsCSxlTP0ZfAxdTdH6NONk2gCwZn1L0UcLZoydmkL2U07AmFZhAJKfxAZQZKmbKabU15RpslcNpPMhYMT5CoPDGsy5wi3EGNqCwyFkKPEh14RZSZJzwVfWUYw3hfyAWFz1gkiAcSEn8GWQA/FX0gUNHTn8f5izZQnWeRAWgwcGVJzwgscy9RzI3xCoY1l7KX6BjkPIEOgqQ7%2BAaHqAFPByQxltshCAZKxVUpBDwbkE4TgW2JHrckAAtW5Z8fDkpJvWZGO0ihJkbPWZSc0iyYXNjnHY8gG4IHkPnVI6RqzuHkDjBGKhe6Bx8E3EQtgEjGBzng34Y0RpI1sMSB0qA9BuhdFoQRDjhFzXMiobk7hgDHURl2NKtZIgWW0MQda3T81zXMMdOsWcJrXjKWW%2BQVAqDHV0M23Q0V7EEBSQkEoRRiTdByO9ANZyeASqlTIfwF93wGIkr4CGNj%2BkevhtQKgGFzApgbSu7JWd12Nv4PwPK110wdpSXNAx6bM0Xx4LurwfF53lM9UkjCs7EKiImaOvW4711to6bwA9TypLdt7f2zAg68ZeB4G2m9n0%2BCjNPQdcO/RfD0uQ0h1DKHOJ3ueUmRpPxRkYcPUGCGdaiMLoNToZAAbF6uydT2p0AbbWRFpkWktFl8F1ofboJ9X632So/ROl936%2BC/ucv%2BwMHHOMGOfTwUREnTjyHfV%2BiTIz3x5JE3W8pzGfKXt0Ep6TkGp3KbdIBkDwHQPacnTh5pOn8m1AQyYXUmHkIke4vFOKsUNzhozAG4kwtyIiCWn2r0Wp43uETcelNFkSzImzHNHYKhhjZo8Lmwk9bzJRe1jFpM%2BQS3HXSC4kD5atTABljW1MdwQ4JrmuFkR5ktToG9QYuLwxoPyEdGwGW8gohcrZTfYcIdyv8pCFweQbQOjAPYqIIor5CDGCzO4ImIcLBZ0dBAJrp79mjLqyoP99xQ7leq3NcjAaJ3gZO1O0ZM653SfTPy/w5sh7IAo7YKjRBGYLfKw%2Busbaak4S/S8G%2BoI9N7uE5NXb/W9uhwOw9p7ODwOw/O1Ov0s74ovtwgoRo%2BcjvPaxp/esrBFCoAYCNcI08vU%2BqCKPeQpIQSAxBGyqKNTr7vYhyHHGq2C1UGkb9l9nPTj1EBw24HPB6S1sxOD8X%2B2wnQ4DVyNdP3jpY65Nun74MWflY%2ByQdnsWucC7wC8BYuuhcHr5QoGWE9bAqE5uxL%2BogicesXqUBGLA4YMGWjLbOzOJcs6h5k0Z5k1tC/XUGfl5hyJDybgLJ1hBRue/VyzyrSbO0RcwEYFgxJaAtBzpJ/oPhc95/zz4VXXvi9Q/MBhdptmU9anT5n9d/LfBfGZmhtDRe48l6l1nVJ/RcNXxs/a1PNfCBO1R/y%2BGtP5DCkpE6lKJAnTf2rGwH1fRPQjSYPIUUE4ABq5qVAlFj2r0OivMq0z3aMuXZWD/F6P530/75lfXbR3jvSdoydEEMbWKbc2PAFhxO7nO3QZZJsCwjBTBzALB98JdXMNw2828E8ocq808uBM8n1%2BgFg0D0CMDMCGcL9i9yssdYceBZ0oNb9Awj8y8nN4MEDiQkCh8g969G9sIkN6UmDg5aB8Ng8FBBAmBBxzUPBtAC1%2BdNo65p9kAmBx8LA59txU18RGYCBSYwcYDD9HsA0CDOJiDDMIAyCT9oQbhA4%2B8qCaDa5W9L98DL174DMYN3xpdbAlc8ApJYJ%2B9q9DDh8cDFCQ5TDwMop1DLC3Rr87Ce89DPpK8B9nDjDID3M3DPsC1/DRludpM9c%2BdDcf1hc1Mxc3CscTseA8B2UWsiCLCDFSDlDbACC9d1DbNfBEMTBfAqiajqi6iaiwi28MjL1sjhsCU1D8jRlfCiiSiLMAj6g%2B8Kj6jhjaiGj503CSMYCoDPd%2BUGx8hLJRC9IuCTBihoIlokEZ4WAe0s5lZPcc1gA805jXxzIWAt5LBUkO1/AyFzJLxmjxiWcojtcftDdTh8UX09dkiRdXCWdNNkxwNSjOj3xTjYZLAWjJ1Gjpi9tISQ4DVfjL08iMJ85Z0wcEU4TwMOjETCB5AJ0UTLJmYrAFtUTsstNTs%2BjRl85oExciTi0SSeB%2BtiDySsTwcqS8TmZCTQ03MiYNw6xqg7V3QV9St%2BJ9i806xfjiQSwcRL0FMGTEdCCrshSksDiUtRTiTT1wNpSEdnJkdb1rphTlSC1fi2j1TeMyTZTtTVcF0oc00M02Bopz1F8DFqhXdY03Q8stoTAIBrTM1OVJRookpdwwB7TgAdsF10wVAURIgPTuAuATAOZHQRo38lpzBRCigSldYuA6Rop0xbo9JsA35JViQWwNhyRuRTpQwmlrpnMgxwyIQoyMyVBpFuAvAGzel0yuBD4bTgyLSbtLIdgsTjU18OSF1HiG0uMecJlKQd0jdUjykgzL1OcZSC5rZsUcFHMYouSezzVAoo9a4MwwDHchz0wDJ8xMw0ZxSUoMJnRhpgpooRzN1Yid0JypzPjUjGRMITd7V%2BhtYdBy1zj858ASxZpW51Mgw9SRpzBF8ILgAkwLcMI5zwMFydTRNF1Xwtdawdd10Yiip30kihMUiJIOD5AlsrcIBWhtsty7dbBs5wzCQYtjoEhmt24wFPVBA8hGYQQBzkL0wxSJSwTvCDEyLYJ8MuKgwyLLN3lMIszCNpjNQPgcQcJ%2BUEhCxloohdC6wv4pDURw8Y8QLAwodbAOIpKxMl06wYj3wpNawJlRV3jpzQdQypIeKUo%2BK%2BIND4Kr0%2Bj4NbBLxmRJMjLbwDLtkAqNsioX4/LJjZLoZdLE8UkgqZE/KRyzLw5MKeMpVMLbLRcsMBI3QIA3LEKgiRtvKwBJMQzMruLVTzzdR/jwTki/KXN3NMrYq3QgqCo5NQq1y6r8YdR%2BIqzAxjyQQJSyFzxLyhptRgpOsMwi0CAKqLT%2BUgw%2BqPU5wprNoudjpfN0B/MGBAssA09lqVc1zZiLUdyc5YytB%2BJ5rTz1t3wrzRqrp%2BI6xedroH17zzLHzskp8FAlUVU1UNVtVny8KSqpJI1o1kzOMcrOzjpfjjoV1LJPrVV1UtVsBjpqB%2BBeUYblU4afrsBbriMoDKxvMY040wbM0IbVSoaEVYbvqEbjpfS1z5qTqlqhqDFrqiRLwusLruzCL6JOMlsfJGxXwotF5Ysc40BYZjZ%2BIEqxz4iUq9Y0qXy7Lyl/S0SsiLMfLYNyreKqqb0Qy28TAazIzXS3hjAyBPoYzFqqYjb4YWZPovB0AqYotGyQgRAbonKakmyRA0zozK0SaaS/jlab1vaS1prNahRsb7K7JcbpKIi60LrhqXQWbooga1qNqtqtQia2AA6fJ2yTpQ6jyBok6AtuhtrBqeCmaRqWbxr2b7jSMFjGw40eazd%2BaXBBbaxhaSBALU02NhyTLojJbLLskFMbK5aMqpJFbVT74/b3xVb3xHLKqJ7tbQ50xdaIzsYDasBYZikMz86GBzaQTXwramzbbnbh8XKvAnbeLXbHaPaMyvassfbx6ASM6pqNa56wqI7w6o7q7gFl7p9dgEAGZOMBal4Q1O6wzl66yuAJSqYmycRnRhD6aoGQgYGmxhCt6EGvAkHkAnZWzPbnaTBi7MBjp6b8HVqCh1qC6gsdrzwc6gwcz5A8yCy9ZToTgThsAbZwYoCpAD5t6pAvApBSAWBpBfA%2BHUBpACpkiMxRBxBER/haA%2BGTYpAhHflSBdAQBLErAfFQ4jE2DH5/Ag4Fh/BGBpB6g%2BGBGFHhHpA%2BH1BfBSB5HFHSA4BYAkB5YrcyAKAIBnHHRiAQBsxkBaAQ56B8Ql5g1KAEhBG%2BGkhgSShpBZHVZbBBxIgJwWBSgwnSAsBgJgBHQUmALjY8AjB1AzHSBMB%2BgmxhhJApAYmIRigUmPgEgCgfwlIsBom5GURKZynOHGBWB2AhdGBPh1BIAD5UAVA3UPBpAoh%2Bgcwqh8gGxpEuDtQShamvR1pOsJwDcoh7NFrJmZZLHJGJA6AuHpBeH%2BGUmRGpB%2BgRUogEkfHMIQ5sIrBaA3RcBCB4ZFh6BAU1YXHPoFheUxG8K5GwmlHGxyJjBKAD4VGQhsJZgoV9H/AtHIXPFDGpBjGjmCmTnLGQBrHbH2mHHEAUB3nPHyBKAPHgXgAYXSBAmpt1AIBQmCmIn8gom2nYn4mCBEnkmCm0nOnMm2W8B27cmXAUmimSmigmnyAVJuGhGenam6WGmymYmMtWnZGuHOmUBumam%2BngpSBBnhn8mxmJmGZpmalZn0B5mUR0AlmogVnOt1mtBNmEBtmxBdmX4EXDnTHxWTmznagLn%2BdgBHtsSrADcIBHnXtPnXmlILcPnFg6RAUhc/mzGlGwX6h1GhUk3k2Q5agEWkWXXzGpA0WMX/n7GYAcWiXXHCW8XgW9hOnZhg4yXOYKWQmUnaWfxhW7AmWWX6XxX2W2BOX23uWcm8n%2BXimKMhWGXKmxW%2BGam6mShpXhW5WmnFW2BlXkiJW1WBmhnMhtXxnOs9XrgDXWAjWFnTWQRzXVmrXohpZbWhAdnJBHXuHnXjnpB3XPX6wNZ5BZhbmcIA38Ag2XnjpQ3LdPHPmjko3kiY3FHQXVGvB02TG73s2hB0WbH/n9mpAuAoOUWLH4PY2D4QD0h3BY4gA) N = (10000000 / (3*4*21) * 3*4*21) = 9999864（即 10000000 向下舍入为 252 次迭代“展开”因子的倍数，因此我的简单实现是执行相同的工作量，不计算重新排序target[]它不这样做，所以它确实打印不匹配的结果。但打印的计数与参考数组的另一个位置匹配。）

我连续运行该程序 4 次（以确保 CPU 预热到最大睿频），其中一次运行看起来不错（3 次都没有异常高）。

参考：最佳位循环（下一节）
快：@njuffa 的代码。（使用 128 位 AVX 整数指令自动矢量化）。
渐进：我的版本（不是由 gcc 或 clang 自动矢量化，至少不是在内部循环中。）gcc 和 clang 完全展开内部 12 次迭代。

gcc8.2 -O3 -march=sandybridge -fpie -no-pie
参考：0.331373 秒，快速：0.011387 秒，渐进：0.009966 秒
gcc8.2 -O3 -march=sandybridge -fno-pie -no-pie
参考：0.397175 秒，快速：0.011255 秒，渐进：0.010018 秒
clang7.0 -O3 -march=sandybridge -fpie -no-pie
参考：0.352381 秒，快速：0.011926 秒，渐进：0.009269 秒（端口 7 uop 的计数非常低，clang 使用索引寻址进行存储）
clang7.0 -O3 -march=sandybridge -fno-pie -no-pie
参考：0.293014 秒，快：0.011777 秒，渐变：0.009235 秒

-三月=天湖（允许 AVX2 用于 256 位整数向量）对两者都有帮助，但 @njuffa 的帮助最大，因为它的更多向量化（包括其最内部的循环）：

gcc8.2 -O3 -march=skylake -fpie -no-pie
参考：0.328725秒，快速：0.007621秒，渐进：0.010054秒（gcc显示“渐进”没有增益，只有“快速”）
gcc8.2 -O3 -march=skylake -fno-pie -no-pie
参考：0.333922 秒，快速：0.007620 秒，渐进：0.009866 秒
clang7.0 -O3 -march=skylake -fpie -no-pie
参考：0.260616 秒，快速：0.007521 秒，渐进：0.008535 秒（不知道为什么渐进比 -march=sandybridge 更快；它没有使用 BMI1andn。我猜是因为它使用 256 位 AVX2 作为 k=0..20 外循环vpaddq)
clang7.0 -O3 -march=skylake -fno-pie -no-pie
参考：0.259159 秒, 快：0.007496 秒，渐变：0.008671 秒

没有 AVX，只需 SSE4.2： (-march=nehalem），奇怪的是 clang 的渐变速度比 AVX/tune=sandybridge 更快。 “快”仅比 AVX 慢一点点。

gcc8.2 -O3 -march=skylake -fno-pie -no-pie
参考：0.337178 秒，快速：0.011983 秒，渐进：0.010587 秒
clang7.0 -O3 -march=skylake -fno-pie -no-pie
参考：0.293555 秒，快：0.012549 秒，渐变：0.008697 秒

-fprofile-generate / -fprofile-use为 GCC 提供一些帮助，特别是对于默认情况下根本不展开的“ref”版本。

我强调了最好的，但它们通常都在彼此的测量噪声容限内。这并不奇怪-fno-pie -no-pie有时更快：索引静态数组[disp32 + reg] is not索引寻址模式，只需基址 + disp32，因此它不会在 Sandybridge 系列 CPU 上分层。

但有时用 gcc-fpie更快；我没有检查，但我认为当 32 位绝对寻址成为可能时，gcc 只是搬起了石头砸自己的脚。或者只是代码生成中看似无辜的差异恰好导致了对齐或 uop 缓存问题；我没有详细检查。

对于 SIMD，我们可以简单地执行 2 或 4xuint64_t并行地，仅在最后一步中水平累加，将字节扩展为 32 位元素。（也许通过在车道上洗牌然后使用pmaddubsw乘数为_mm256_set1_epi8(1)将水平字节对添加到 16 位元素中。）

TODO：手动矢量化__m128i and __m256i (and __m512i）的版本。应比上述“渐进”时间快接近 2 倍、4 倍甚至 8 倍。硬件预取可能仍然可以跟上它，除了数据来自 DRAM 的 AVX512 版本，特别是在存在来自其他线程的争用的情况下。我们阅读的每个 qword 都做了大量的工作。

过时的代码：对位循环的改进

您的便携式标量版本也可以改进，将其速度从约 1.92 秒（总体分支误预测率为 34%，快速循环被注释掉！）到〜0.35秒（clang7.0 -O3 -march=sandybridge）在 3.9GHz Skylake 上具有适当的随机输入。或者 1.83 秒的分支版本!= 0代替== m，因为编译器无法证明m始终恰好设置 1 位并/或进行相应优化。

（相对于 @njuffa 或我上面的快速版本为 0.01 秒，所以这在绝对意义上是相当无用的，但作为何时使用无分支代码的一般优化示例值得一提。）

如果您期望零和一的随机组合，那么您需要一些不会错误预测的无分支的东西。正在做+= 0for 零元素避免了这种情况，并且也意味着 C 抽象机肯定会触及该内存，而不管数据如何。

编译器不允许发明写入，所以如果他们想自动矢量化你的if() target[i]++版本，他们必须使用像 x86 这样的屏蔽存储vmaskmovps避免对未修改的元素进行非原子读取/重写target。因此，一些假设的未来编译器可以自动向量化普通标量代码，这样会更容易。

无论如何，一种写法是target[i] += (pLong[j] & m != 0);，使用 bool->int 转换得到一个 0 / 1 整数。

但是，如果我们只是移动数据并隔离低位，那么对于 x86（可能对于大多数其他架构）我们会得到更好的 asm&1。编译器有点愚蠢，似乎没有发现这种优化。他们很好地优化了额外的循环计数器，并把m <<= 1 into add same,same有效地左移，但他们仍然使用异或零/test / setne创建一个 0 / 1 整数。

像这样的内部循环编译效率稍高（但仍然很多很多比我们使用 SSE2 或 AVX 所做的更糟糕，甚至比使用 @chrqlie 的查找表的标量更糟糕，当像这样重复使用时，它会在 L1d 中保持热状态，从而允许 SWARuint64_t):

    for (int j = 0; j < 10000000; j++) {
#if 1  // extract low bit directly
        unsigned long long tmp = pLong[j];
        for (int i=0 ; i<64 ; i++) {   // while(tmp) could mispredict, but good for sparse data
            target[i] += tmp&1;
            tmp >>= 1;
        }
#else // bool -> int shifting a mask
        unsigned long m = 1;
        for (i = 0; i < 64; i++) {
            target[i]+= (pLong[j] & m) != 0;
            m = (m << 1);
        }
#endif

注意unsigned long不保证是 64 位类型，并且不在 x86-64 System V x32（64 位模式下的 ILP32）和 Windows x64 中。或者在 32 位 ABI 中，例如 i386 System V。

已编译在 gcc、clang 和 ICC 的 Godbolt 编译器浏览器上 https://godbolt.org/#z:OYLghAFBqd5QCxAYwPYBMCmBRdBLAF1QCcAaPECAKxAEZSAbAQwDtRkBSAJgCFufSAZ1QBXYskwgA5NwDMeFsgYisAag6yAwoIL5UAOgQbsHAAwBBOQqUrM6rTvQM8AI0PGzlrvMXK1G7QJiBWB3WRMLK19be20AT0EAegI8AFtMMIivHxt/BwTEnSYCTM9PdFEXBjtgTAIAFQBlCABKdQB2PgtVHtUdYhFkAlUU9IA3JgYRsY0u817VWoJRzFQAM3QmOIhuADYCMdJVADkAVQAZc5bZzwXiOrEWaf0DgH1BTGR1Xme3kQ%2BvolVLRTKCwfpTDdIu0ACKeEQsQR4YAsTDoVQKYYEJjEJYcACsPF2ABYCXDZHNPJjVKkmAoINSccBkEdkAgcaoAFRMmaEsltDidW69BFIlFo1QMVBsVQAD3sMOB3ih816bI5nIADgqaZMpcgIAAOLnAsFm64U4U9PBrVQQMBgTUCoXdBa9TXBFgENY7LhcNZ06rooiqPWoZDFTAEzQsbhcC1zN29TCywgQUwJq0dOGunrpVIfAgQTVHUyytZrI6Gs1gzO51QezE%2BuOajRwrj41v4mNxo5OlULUXI1HoqUyrXnaXAHUQIfi0dTzktVuW%2BsVERVOxFYjDNuLOpNVoD3prEh26lUHWQimqS8BU01688W/8fjOyneG3AnqJIEpoJMEMkqoAA7qoLiEKo%2BD3EMDBxFmg6IsOEpjtOqEjKk2p7pqk5sASPBUGSx5JqexDnl6GJtqY9jPngAQkjRlG8G%2BHSJkmbrYridT4XR%2BKKq%2BsiKgQmF7LQxHsT0wlYeEHiCcC4kLIKObZJgDAfKov7gagqBTAAtMYGIUYICA2ikMpMLqggANYISKSHzsBMqpDqYmrqqJFngyV6zBisSqCSPl0cxvDvrZSacXihK8e2fBycWuHAPhhF8d8uw0m0Dp7k%2BYVus5e4QHlWj3rQdbuW6SllN4mAsPgaxZhV9aNl6zZ%2BkQqBWd8%2BK2gCgjRrGfpHEsh5tLpfScQQpV3A8xBPNl0IwlILSMNI%2BJSKQLDSKYa2oNImivj8whiBI3yyLQa0EJti1LVZIBcMS%2BgAJxcKYD2yLstC0LssiGs933LVIxJrRtUhbaQO1SGtgggKYpAXSDi2kHAsBIGgmF4NUZAUBAqOaujmDECAeDIMgtAPfQazowQ%2BNQxALiXaQEEsDicTSGdpCo%2BkXoAPIsHB9NYLSbDVPzeAwSkYyYFD8OkCmnwiFTrNrZiqn084LjEMzmgYJIUhs0EaSK0tzBsCg%2B0CGrUOQEtqCaik0pS7pspfLptLiEYgmCKw6BxOreDoLUqi6VzsiBxWTA6C7xQIJDojiJIJX/at630%2BDsqGrsukMUTXykxC%2Bi0HauCEGecj0KoWtoxjJ3xuXZvnZdLRLQgmBMFgBOtKQN0dhC7QPaYX2Goab20L3uy7P9gPJ9L4OQ9DsMN4jMCICgqCV/j5CUDjeME8AQ%2BkBTDBU8QNN09LjPM4r7Or5zBA83z0sC6wwDCw/oufOLkv07LyDyzrbPKwwVWrgNbEDiFrLAl99apENowJ%2Bptgrm1cJbDuNs7aImkI7Z2rs2Rtk9jVH2wR/Z2CDiHXSYcI60gINHIQscJB0CNitIGKdpBpwzgxZkXxDT6BDhAIuRAyKlyOBXXGVdS5tD2gg%2Bu8NG6d1uvdB6CjFFKKUePQBAMmHT2kLPGGcMtpLSRsvLeGMN7Y1XiI/GKBjbAHaKCfelNqaUFPqDc%2BoDL4c2qrfXmLMH6YEFs/HWoNoLvzwBLKWoNv6/0vgAoB6tNba0gcEaBusEZWPgXwXgjAkHwGtrbPA9sMFO0Dtg92MI8He19kQwOwdQ4BgoVHGOR144MKkEnYGoNU7p0zsSVQSgn6qHaHnaivD8D8JOmXYR28TqyHEXXee0jrq3XxBPDR7StFCDnroq6/0uArO2mszZMiJbHzyRtYkQA%3D%3D%3D，gcc 循环中的微指令少了 1 个。但所有这些都只是简单的标量，其中 clang 和 ICC 展开为 2。

# clang7.0 -O3 -march=sandybridge
.LBB1_2:                            # =>This Loop Header: Depth=1
   # outer loop loads a uint64 from the src
    mov     rdx, qword ptr [r14 + 8*rbx]
    mov     rsi, -256
.LBB1_3:                            #   Parent Loop BB1_2 Depth=1
                                    # do {
    mov     edi, edx
    and     edi, 1                              # isolate the low bit
    add     dword ptr [rsi + target+256], edi   # and += into target

    mov     edi, edx
    shr     edi
    and     edi, 1                              # isolate the 2nd bit
    add     dword ptr [rsi + target+260], edi

    shr     rdx, 2                              # tmp >>= 2;

    add     rsi, 8
    jne     .LBB1_3                       # } while(offset += 8 != 0);

这比我们得到的要好一些test / setnz。无需展开，bt / setc可能是相等的，但编译器不擅长使用bt实施bool (x & (1ULL << n)), or bts实施x |= 1ULL << n.

如果许多字的最高设置位远低于位 63，则循环while(tmp)可能是一场胜利。如果大多数情况下它只节省约 0 到 4 次迭代，那么分支错误预测就变得不值得，但如果它经常节省 32 次迭代，那确实是值得的。也许在源代码中展开，因此循环仅测试tmp每 2 次迭代（因为编译器不会为您执行该转换），但是循环分支可以是shr rdx, 2 / jnz.

在 Sandybridge 系列上，每 2 位输入前端有 11 个融合域微指令。（add [mem], reg使用非索引寻址模式微融合加载+ALU和存储地址+存储数据，其他一切都是单微指令。添加/jcc 宏熔断器。请参阅 Agner Fog 的指南，以及https://stackoverflow.com/tags/x86/info https://stackoverflow.com/tags/x86/info）。因此，它应该以每 2 位 3 个周期 = 每 96 个周期 1 个 uint64_t 的速度运行。（Sandybridge 不会在其循环缓冲区内部“展开”，因此非 4 倍微指令计数基本上会向上舍入，这与 Haswell 及更高版本不同）。

与 gcc 的未展开版本相比，每 1 位 7 个微指令 = 每位 2 个周期。如果你编译的是gcc -O3 -march=native -fprofile-generate/ 测试运行 /gcc -O3 -march=native -fprofile-use，配置文件引导的优化将启用循环展开。

对于完全可预测的数据，这可能比分支版本慢，就像您从memset具有任何重复字节模式。我建议用快速 PRNG（如 SSE2 xorshift+）随机生成的数据填充数组，或者如果您只是对计数循环进行计时，则使用您想要的任何数据，例如rand().

本文内容由网友自发贡献，版权归原作者所有，本站不承担相应法律责任。如您发现有涉嫌抄袭侵权的内容，请联系:hwhale#tublm.com(使用前将#替换为@)