如果没有 AVX2 指令中的快速收集和分散,你会怎么做?

What do you do without fast gather and scatter in AVX2 instructions?

我正在编写一个程序来检测素数。一部分是筛选出可能的候选人。我写了一个相当快的程序,但我想我会看看是否有人有更好的想法。我的程序可以使用一些快速收集和分散指令,但我仅限于 x86 架构的 AVX2 硬件(我知道 AVX-512 有这些,但我不确定它们有多快)。

#include <stdint.h>
#include <immintrin.h>

#define USE_AVX2

// Sieve the bits in array sieveX for later use
void sieveFactors(uint64_t *sieveX)
{
    const uint64_t totalX = 5000000;
#ifdef USE_AVX2
    uint64_t indx[4], bits[4];

    const __m256i sieveX2 = _mm256_set1_epi64x((uint64_t)(sieveX));
    const __m256i total = _mm256_set1_epi64x(totalX - 1);
    const __m256i mask = _mm256_set1_epi64x(0x3f);

    // Just filling with some typical values (not really constant)
    __m256i ans = _mm256_set_epi64x(58, 52, 154, 1);
    __m256i ans2 = _mm256_set_epi64x(142, 70, 136, 100);

    __m256i sum = _mm256_set_epi64x(201, 213, 219, 237);    // 3x primes
    __m256i sum2 = _mm256_set_epi64x(201, 213, 219, 237);   // This aren't always the same

    // Actually algorithm can changes these
    __m256i mod1 = _mm256_set1_epi64x(1);
    __m256i mod3 = _mm256_set1_epi64x(1);

    __m256i mod2, mod4, sum3;

    // Sieve until all factors (start under 32-bit threshold) exceed the limit
    do {
        // Sieve until one of the factors exceeds the limit
        do {
            // Compiler does a nice job converting these into extracts
            *(__m256i *)(&indx[0]) = _mm256_add_epi64(_mm256_srli_epi64(_mm256_andnot_si256(mask, ans), 3), sieveX2);
            *(__m256i *)(&bits[0]) = _mm256_sllv_epi64(mod1, _mm256_and_si256(mask, ans));

            ans = _mm256_add_epi64(ans, sum);

            // Early on these locations can overlap
            *(uint64_t *)(indx[0]) |= bits[0];
            *(uint64_t *)(indx[1]) |= bits[1];
            *(uint64_t *)(indx[2]) |= bits[2];
            *(uint64_t *)(indx[3]) |= bits[3];

            mod2 = _mm256_sub_epi64(total, ans);

            *(__m256i *)(&indx[0]) = _mm256_add_epi64(_mm256_srli_epi64(_mm256_andnot_si256(mask, ans2), 3), sieveX2);
            *(__m256i *)(&bits[0]) = _mm256_sllv_epi64(mod3, _mm256_and_si256(mask, ans2));

            ans2 = _mm256_add_epi64(ans2, sum2);

            // Two types of candidates are being performed at once
            *(uint64_t *)(indx[0]) |= bits[0];
            *(uint64_t *)(indx[1]) |= bits[1];
            *(uint64_t *)(indx[2]) |= bits[2];
            *(uint64_t *)(indx[3]) |= bits[3];

            mod4 = _mm256_sub_epi64(total, ans2);
        } while (!_mm256_movemask_pd(_mm256_castsi256_pd(_mm256_or_si256(mod2, mod4))));

        // Remove one factor
        mod2 = _mm256_castpd_si256(_mm256_blendv_pd(_mm256_setzero_pd(), _mm256_castsi256_pd(sum), _mm256_castsi256_pd(mod2)));
        mod4 = _mm256_castpd_si256(_mm256_blendv_pd(_mm256_setzero_pd(), _mm256_castsi256_pd(sum2), _mm256_castsi256_pd(mod4)));
        ans = _mm256_sub_epi64(ans, mod2);
        ans2 = _mm256_sub_epi64(ans2, mod4);
        sum = _mm256_sub_epi64(sum, mod2);
        sum2 = _mm256_sub_epi64(sum2, mod4);
        sum3 = _mm256_or_si256(sum, sum2);
     } while (!_mm256_testz_si256(sum3, sum3));
#else
     // Just some example values (not really constant - compiler will optimize away code incorrectly)
     uint64_t cur = 58;
     uint64_t cur2 = 142;
     uint64_t factor = 67;

     if (cur < cur2) {
        std::swap(cur, cur2);
    }
    while (cur < totalX) {
        sieveX[cur >> 6] |= (1ULL << (cur & 0x3f));
        sieveX[cur2 >> 6] |= (1ULL << (cur2 & 0x3f));
        cur += factor;
        cur2 += factor;
    }
    while (cur2 < totalX) {
        sieveX[cur2 >> 6] |= (1ULL << (cur2 & 0x3f));
        cur2 += factor;
    }
#endif
}

请注意,这些位置最初可能会重叠。在循环中片刻之后,情况并非如此。如果可能的话,我很乐意使用不同的方法。这部分算法中大约 82% 的时间都在这个循环中。希望这与其他已发布的问题不太接近。

IDK 为什么使用同一个 cur[8] 数组的不同部分作为索引和值;弄清楚只有一个真正的数组,这让源代码更难理解。另一个只是将矢量反弹为标量。

看起来您只是在进行矢量 -> 标量,而不是将标量插入回矢量。而且循环内的任何内容都不依赖于 sieveX[] 中的任何数据;我不熟悉您的筛选算法,但我想这样做的目的是在内存中创建数据供以后使用。


AVX2 有收集(不是分散),但它们仅在 Skylake 和更新的 上速度快。他们在 Broadwell 上还好,在 Haswell 上速度慢,在 AMD 上慢。 (就像 Ryzen 的 vpgatherqq 每 12 个时钟一个)。参见 http://agner.org/optimize/ and other performance links in the x86 tag wiki

Intel 的优化手册中有一小部分是关于手动收集/分散(使用 insert/extract 或 movhps)与硬件指令,可能值得一读。在这种情况下,索引是 运行 时间变量(不是恒定步长或其他东西),我认为 Skylake 可以从 AVX2 此处收集指令中受益。

请参阅 Intel 的内在函数指南以查找像 movhps 这样的 asm 指令的内在函数。我只是在谈论你想让你的编译器发出什么,因为这很重要,而且 asm 助记符更短,不需要转换。你必须知道 asm 助记符才能在 Agner Fog 的指令表中查找它们,或者从 auto-vectorization 读取编译器输出,所以我通常在 asm 中思考然后将其转换为内在函数。


使用 AVX,您有 3 个主要选项:

  • 做一切标量。注册压力可能是个问题,但根据需要生成索引(而不是一次执行所有 4 个添加或子生成 curr[4..7])可能会有所帮助。除非那些 mask 向量在不同的元素中具有不同的值。

(将内存源用于标量常量可能还不错,但是,如果它们不适合 32 位立即数,并且如果您不在每个时钟的 2 个内存操作上遇到瓶颈。memory-destination or 指令将使用索引寻址模式,因此无法使用 Haswell 端口 7 上的专用 store-AGU 和更高版本。因此 AGU 吞吐量可能成为瓶颈。)

提取一个向量的所有 4 个元素作为标量比 4 倍标量 add 或移位指令更昂贵,但您所做的工作比这更多。不过,对于 1 微秒 variable-count 的 BMI2 偏移(而不是英特尔的 3 微秒),它可能并不可怕。不过,我认为我们可以使用 SIMD 做得更好,尤其是在仔细调整的情况下。

  • 将索引和值提取为标量,就像您现在所做的那样,因此 sieveX[] 中的 OR 是纯标量 。即使两个或多个索引相同也有效。

    每个 ymm 向量花费大约 7 微指令 -> 使用提取 ALU 指令的 4x 标量寄存器,或使用 store/reload 的 5 微指令(编译器值得考虑,可能是 4 个向量提取中的一个或两个,因为此代码可能无法解决加载/存储端口吞吐量的瓶颈。)但是,如果编译器将 C 源代码中的 store/reload 转换为 shuffle/extract 指令,则您无法轻易覆盖其策略除了可能使用 volatile。顺便说一句,您需要使用 alignas(32) cur[8] 来确保实际的矢量存储不会跨越 cache-line 边界。

or [rdi + rax*8], rdx (with an indexed addressing mode preventing full micro-fusion) 在现代英特尔 CPU(Haswell 及更高版本)上为 3 微指令。 我们可以通过缩放+使用 SIMD 添加到数组基地址来避免索引寻址模式(使 front-end 为 2 微指令):例如srli 乘以 3 而不是 6,屏蔽掉低 3 位 (vpand),并用 set1_epi64(sieveX) 屏蔽掉 vpaddq。因此,每个索引向量需要 2 条额外的 SIMD 指令才能在 SnB-family 上节省 4 微指令。 (您将提取 uint64_t* 指针元素而不是 uint64_t 索引。或者如果 sieveX 可以是 32 位绝对地址 1,您可以跳过 vpaddq 并提取 already-scaled 索引以获得相同的增益。)

它还会在端口 7(Haswell 和更高版本)上启用 store-address 微指令到 运行;端口 7 上的简单 AGU 只能处理 non-indexed 寻址模式。 (这使得使用 store+reload 将值提取为标量更具吸引力。您希望提取索引的延迟更低,因为在 memory-dst or 的加载部分完成之前不需要这些值。)对于调度程序/执行单元来说确实意味着更多 unfused-domain 微指令,但很值得权衡。

这不是其他 AVX2 CPU(Excavator/Ryzen 或 Xeon Phi)的胜利;只有 SnB-family 具有 front-end 成本和 execution-port 索引寻址模式的限制。

  • 提取索引,手动将 vmovq / vmovhps 收集到一个向量中用于 SIMD vpor,然后使用 [= 分散回来34=] / vmovhps.

    就像硬件 gather/scatter 一样,正确性要求所有索引都是唯一的,因此您需要使用上述选项之一,直到达到那个点t 在你的算法中。 (矢量冲突检测 + 回退与总是提取标量相比不值得付出代价:)。

    请参阅 了解内在函数版本。 (我知道我最近写了一个手动收集/分散的答案,但我花了一段时间才找到它!)在那种情况下,我只使用了 128 位向量,因为没有任何额外的 SIMD 工作来证明额外的vinserti128 / vextracti128.

实际上,我想在这里你想要提取 _mm256_sllv_epi64 结果的高半部分,这样你就可以将 cur[4..5]cur[6..7] (数据)分成两个单独的部分__m128i 个变量。你会有 vextracti128 / 2x vpor xmm 而不是 vinserti128 / vpor ymm / vextracti128.

前者端口 5 压力较小,instruction-level 并行性更好:两个 128 位的一半是独立的依赖链,不会相互耦合,因此 store/reload 瓶颈(和缓存未命中)影响较少的依赖微指令,允许 out-of-order 执行在等待时继续处理更多内容。

在 256b 向量中进行地址计算并提取指针而不是索引将使 vmovhps 负载在 Intel 上更便宜(索引负载不能保持 micro-fused 到 vmovhps2)。请参阅上一个要点。但是 vmovq loads/stores 始终是单个 uop,并且 vmovhps 索引存储可以在 Haswell 及更高版本上保留 micro-fused,因此 break-even 对于 front-end在 AMD 或 KNL 上吞吐量更差。这也意味着调度程序/执行单元有更多 unfused-domain 微指令,这看起来比 port2/3 AGU 压力更像是一个潜在瓶颈。唯一的好处是 store-address 微指令可以 运行 在端口 7 上,减轻一些压力。

AVX2 给了我们一个新的选择:

  • AVX2 vpgatherqq 用于聚集(_mm256_i64gather_epi64(sieveX, srli_result, 8)),然后提取索引并手动分散。所以它与手动聚集/手动分散完全一样,除非您用 AVX2 硬件收集替换手动收集。 (两次 128 位收集的成本高于一次 256 位收集,因此您可能希望采用 instruction-level 并行度并收集到单个 256 位寄存器中)。

可能在 Skylake 上获胜(其中 vpgatherqq ymm 是 4 uops / 4c 吞吐量,加上 1 uop 设置),但甚至不是 Broadwell(9 uops,每 6c 吞吐量一个)而且绝对不是 Haswell(22 uops / 9c 吞吐量)。无论如何,您确实需要标量寄存器中的索引,因此您 节省了 manual-gather 部分工作。真便宜。


Skylake 上每个策略的总成本

看起来这不会在任何一个端口上造成严重的瓶颈。 GP reg->xmm 需要端口 5,但 xmm->int 需要 SnB-family CPU 上的端口 0,因此当与提取所需的洗牌混合时,它不太可能在端口 5 上出现瓶颈。 (例如 vpextrq rax, xmm0, 1 是一条 2 uop 指令,一个端口 5 shuffle uop 用于获取高 qword,一个端口 0 uop 用于将该数据从 SIMD 发送到整数域。)

所以你的 SIMD 计算需要经常提取一个向量到标量 比需要频繁地将标量计算结果插入向量中要好。另请参阅 ,但这是在谈论从 GP regs 开始的数据,而不是内存。

  • 提取两者/标量或:总计 = 24 uops = 6 个周期的 front-end 吞吐量。

  • vpaddq + vpand 地址计算(Skylake 上的端口 0/1/5 为 2 微指令)

  • 2x vextracti128(端口 5 为 2 微指令)

  • 4x vmovq (4 p0)

  • 4x vpextrq (8: 4p0 4p5)

  • 4x or [r], r (4x2 = 8 front-end 每个微指令。后端: 4p0156 4p23 (load) 4p237 (store-addres) 4p4 (store-data )). Non-indexed寻址模式。

p5 总计 = 6 微指令,勉强适合。 Store/reload 对于数据提取看起来很合理,如果你能让你的编译器这样做的话。 (但编译器通常不会对管道进行足够详细的建模,以在同一循环中使用混合策略来平衡端口压力。)

  • 手动 gather/scatter:20 微指令,5 个周期 front-end 吞吐量 (Haswell / BDW / Skylake)。在 Ryzen 上也不错。

  • (可选,可能不值得):vpaddq + vpand 地址计算(Skylake 上的端口 0/1/5 为 2 微指令)如果可以使用 non-VEX 则跳过这些 movhps 用于 1-uop micro-fused 索引负载。 (但随后 p237 商店变成了 p23)。

  • vextracti128 指针(端口 5 1 uop)

  • 2x vmovq 提取 (2p0)

  • 2x vpextrq (4 = 2p0 2p5)

  • 2x vmovq 加载 (2p23)

  • 2x vmovhps xmm, xmm, [r] non-indexed 负载 (2 front-end uops micro-fused: 2p23 + 2p5)

  • vextracti128 拆分数据 (p5)

  • 2x vpor xmm (2p015)

  • 2x vmovq 存储(2x 1 micro-fused uop,2p237 + 2p4)

  • 2x vmovhps 存储(2x 1 micro-fused uop,2p237 + 2p)

端口瓶颈:4 p0 和 4 p5 适合 5 个周期,特别是当你将它与你的循环混合时,它可以 运行 端口 1 上的几个微指令。在 Haswell paddq 上只有 p15(不是 p015),班次只有 p0(不是 p01)。 AVX2 _mm256_sllv_epi64 在 Skylake 上是 1 uop (p01),但在 Haswell 上是 3 uops = 2p0 + p5。因此 Haswell 可能更接近此循环的 p0 或 p5 瓶颈,在这种情况下,您可能需要查看一个索引向量的 store/reload 提取策略。

跳过 SIMD 地址计算可能很好,因为除非您使用 store/reload 提取,否则 AGU 压力看起来不是问题。这意味着 uop 缓存中的指令更少/更小 code-size 和更少的 uops。 (un-lamination 直到解码器/uop 缓存之后才会发生,因此您仍然可以在 front-end 的早期部分受益于 micro-fusion,只是不在问题瓶颈处。)

  • Skylake AVX2 收集/手动分散:总计 = 18 uops,front-end 吞吐量的 4.5 个周期。(更早的 uarch 或AMD).

  • vextracti128 索引(端口 5 为 1 uop)

  • 2x vmovq 提取 (2p0)

  • 2x vpextrq (4 = 2p0 2p5)

  • vpcmpeqd ymm0,ymm0,ymm0vpgatherqq 创建一个 all-ones 掩码 (p015)

  • vpgatherqq ymm1, [rdi + ymm2*8], ymm0 某些端口为 4 微码。

  • vpor ymm (p015)

  • vextracti128 或结果 (p5)

  • 2x vmovq 存储(2x 1 micro-fused uop,2p23 + 2p4)。注意没有端口 7,我们正在使用索引存储。

  • 2x vmovhps 存储(2x 1 micro-fused uop,2p23 + 2p4)。

因此,即使选择了最佳吞吐量,我们仍然每 4.5 个周期仅管理 4 次加载/4 次存储,而且还没有考虑循环中的 SIMD 工作,这会消耗一些 front-end 吞吐量。 所以我们还没有接近 AGU 吞吐量的瓶颈并且不必担心使用端口 7。

我们可能会考虑 store/reload 对于其中一个摘录(如果我们是编译器),将 7 uop 5 指令 vextracti128 / 2x vmovq / 2x vpextrq 序列替换为 5 uops 存储 / 4x 加载.


总体:一个循环直到我们处理完冲突,然后是一个 SIMD 收集循环

你说在某个点之后,你没有像cur[0] == cur[2]这样的索引之间的冲突(重叠)。

您肯定需要一个完全不检查冲突的单独循环来利用这一点。即使您有 AVX512,Skylake 的 vpconflictq 也是 micro-code 并且速度不快。 (KNL有single-uopvpconflictq但完全避开它还是更快)。

我会留给你(或一个单独的问题)如何确定你何时处理完冲突并离开解决这种可能性的循环。

您可能需要提取索引 + 数据策略,但可能会发生冲突。 SIMD 冲突检查是可能的,但并不便宜,32 位元素需要 11 微指令:。一个 qword 版本显然比 dword 便宜得多(更少的洗牌和比较来获得所有对所有),但你可能仍然只想在你的提取循环中每 10 次迭代左右做一次。

从最好的 scalar-or 版本到最好的收集版本并没有很大的加速(6 周期对比 4.5 没有考虑循环中的其他工作,所以比率甚至比那个还小)。尽快离开稍慢的版本不值得让它慢很多。

因此,如果您可以可靠地检测到您何时处理完冲突,请使用类似

的方法
int conflictcheck = 10;

do {

    if (--conflictcheck == 0) {
       vector stuff to check for conflicts
       if (no conflicts now or in the future)
           break;

       conflictcheck = 10;  // reset the down-counter
    }

    main loop body,  extract -> scalar OR strategy

} while(blah);


// then fall into the gather/scatter loop.
do {
    main loop body, gather + manual scatter strategy
} while();

应该编译成 dec / je,在 not-taken 情况下只需要 1 uop。

对 slightly-slower 循环进行总共 9 次额外迭代 比进行数千次额外昂贵的冲突检查要好得多


脚注 1:

如果 sieveX 是静态的并且您正在 Linux(不是 MacOS)上构建 non-PIC 代码,那么它的地址将适合 disp32 作为 [reg+disp32] 寻址模式的一部分。在这种情况下,您可以省略 vpaddq。但是让编译器将 uint64_t 视为 already-scaled 数组索引(清除其低位)会很丑陋。可能必须将 sieveX 转换为 uintptr_t 并添加,然后再转换回来。

这在 PIE 可执行文件或共享库(不允许使用 32 位绝对地址)或 OS X 上根本不可能(其中静态地址始终高于 2^32) ).我不确定 Windows 允许什么。请注意 [disp32 + reg*8] 只有 1 个寄存器,但仍然是索引寻址模式,因此所有 SnB-family 惩罚都适用。但是如果你不需要缩放,reg + disp32就是base + disp32.

脚注 2:有趣的事实:non-VEX movhps 负载可以留在 micro-fused Haswell 上。它不会在 Skylake 上导致 SSE/AVX 停顿,但是 你不会让编译器在 m 中发出 non-VEX 版本AVX2 函数的 ddle.

不过,IACA(英特尔的静态分析工具)弄错了。 :( What is IACA and how do I use it?.

这基本上是 -mtune=skylake 的 missed-optimization,但它 在 Haswell 上停滞:

Skylake 上的“penalty A”(用脏上层执行 SSE)仅仅是对那个寄存器的错误依赖。 (还有一个用于指令的合并 uop,否则 write-only,但 movhps 已经是其目的地的 read-modify-write。)我在 Skylake 上用 Linux [=84 测试了这个=] 计算 uops,用这个循环:

    mov     r15d, 100000000

.loop:
    vpaddq  ymm0, ymm1, ymm2      ; dirty the upper part
    vpaddq  ymm3, ymm1, ymm2      ; dirty another register for good measure

    vmovq  xmm0, [rdi+rbx*8]       ; zero the full register, breaking dependencies
    movhps xmm0, [rdi+rbx*8+8]     ; RMW the low 128 bits
                          ; fast on Skylake, will stall on Haswell

    dec r15d
    jnz .loop

循环 运行s 在 Skylake (i7-6700k) 上每次迭代约 1.25 个周期,最大化每个时钟 4 微指令的 front-end 吞吐量。总共 5 个 fused-domain 微指令 (uops_issued.any),6 个 unfused-domain 微指令 (uops_executed.thread)。所以 micro-fusion 确实发生了 movhps,没有任何 SSE/AVX 问题。

将其更改为 vmovhps xmm0, xmm0, [rdi+rbx*8+8] 将其减慢到每次迭代 1.50 个周期,现在是 6 fused-domain,但仍然是 6 unfused-domain 微指令。

如果 ymm0 的上半部分在 movhps xmm0, [mem] 运行 时变脏,则没有额外的 uop。我通过注释掉 vmovq 进行测试。但是将 vmovq 更改为 movq 确实 会导致额外的 uop:movq 变成 micro-fused 加载+合并,取代低 64位(并且仍然将 xmm0 的高 64 位归零,因此它不完全 movlps)。


另请注意,即使没有 VEX,pinsrq xmm0, [mem], 1 也无法微熔断。但是对于 VEX,出于 code-size 的原因,您应该更喜欢 vmovhps

您的编译器可能希望将整数数据上 movhps 的内在函数“优化”为 vpinsrq,不过,我没有检查。

我刚刚仔细查看了您在这里所做的事情:对于 mod1 = mod3 = _mm256_set1_epi64x(1); 的情况,您只是在位图中设置单个位,并将 ans 的元素作为索引。

它由两个展开,ans 和 ans2 运行ning 并行,使用 mod1 << ansmod3 << ans2。评论您的代码并使用英文文本解释大局中发生的事情!这只是普通埃拉托色尼筛法的位设置循环的一个非常复杂的实现。 (所以如果问题一开始就这么说就好了。)

并行展开多个 start/strides 是一个非常好的优化,因此您通常在高速缓存行中设置多个位,而它在 L1d 中仍然很热。 一次针对较少不同因素的缓存阻塞具有类似的优势。在移动到下一个之前,针对多个因素(步幅)重复迭代相同的 8kiB 或 16kiB 内存块。为 2 个不同的步幅分别展开 4 个偏移量可能是创建更多 ILP 的好方法。

您 运行 的并行步幅越大,您在第一次接触新缓存行时的速度就越慢。 (提供高速缓存/TLB 预取空间以避免初始停顿)。所以缓存阻塞不会消除多步的所有好处。


步幅 < 256 的可能特例

单个 256 位向量 load/VPOR/store 可以设置多个位。诀窍是创建一个矢量常量或一组矢量常量,其中位位于正确的位置。不过,重复模式的长度类似于 LCM(256, bit_stride) 位。例如,stride=3 将以 3 个向量长的模式重复。除非有更聪明的东西,否则这很快就无法用于奇数/主要步幅:(

64 位标量很有趣,因为按位旋转可用于创建一系列模式,但 SnB 系列 CPU 上的可变计数旋转成本为 2 微指令。

您可以用它做更多的事情;也许未对齐的负载可以以某种方式提供帮助。

位掩码的重复模式甚至对于大步幅的情况也是有用的,例如每次旋转stride % 8。但是,如果您正在 JIT 一个将模式硬编码为 or byte [mem], imm8 的循环,并且选择与重复长度一致的展开因子,那将更有用。


减少冲突 loads/stores

当您只设置一个位时,您不必 load/modify/store 64 位块。您的 RMW 操作越窄,您的位索引就越接近而不会冲突。

(但是您在同一位置没有长循环携带的 dep 链;您将在 OoO exec 停止等待长链末端的重新加载之前继续前进。所以如果冲突不是一个正确性问题,它不太可能在这里产生很大的性能差异。不像位图直方图或其他很容易在附近位上重复命中的长链。)

32 位元素将是一个显而易见的选择。 x86 可以有效地 load/store 双字 to/from SIMD 寄存器和标量。 (标量字节操作也很有效,但是来自 SIMD regs 的字节存储总是需要 pextrb 的多个微指令。)

如果您不收集到 SIMD 寄存器中,ans / ans2 的 SIMD 元素宽度不必与 RMW 宽度匹配。如果您想将位索引拆分为标量中的地址/位偏移,则 32 位 RMW 比 8 位有优势,使用移位或 bts 将移位计数隐式屏蔽为 32 位(或 64 位用于64 位移位)。但是 8 位 shlxbts 不存在。

使用 64 位 SIMD 元素的主要优点是,如果您计算的是指针而不仅仅是索引。如果您可以将 sieveX 限制为 32 位,您仍然可以这样做。例如在 Linux 上分配 mmap(..., MAP_32BIT|MAP_ANONYMOUS, ...)这是假设您不需要超过 2^32 位 (512MiB) 的筛选 space,因此您的位索引始终适合 32 位元素。 如果不是这种情况,到那时您仍然可以使用 32 位元素向量,然后将当前循环用于高部分。

如果您使用 32 位 SIMD 元素而不将 sieveX 限制为 32 位点指针,您将不得不放弃使用 SIMD 指针计算并仅提取位索引,或者仍然在 SIMD 中拆分进入 idx/bit 并提取两者。

(对于 32 位元素,基于 store/reload 的 SIMD -> 标量策略看起来更有吸引力,但在 C 中这主要取决于您的编译器。)

如果您手动收集到 32 位元素,则不能再使用 movhps。您必须对高 3 元素使用 pinsrd / pextrd,而那些从不使用微型保险丝/在 SnB 系列上始终需要 port5 uop。 (不像 movhps 这是一个纯粹的商店)。但这意味着 vpinsrd 在索引寻址模式下仍然是 2 微指令。您仍然可以对元素 2 使用 vmovhps(然后用 vpinsrd 覆盖向量的顶部双字);未对齐的负载很便宜,可以重叠下一个元素。但是你不能做 movhps 商店,这就是它非常好的地方。


您当前的策略有两个性能问题:

有时您将此与 mod1mod3 的某些元素一起使用 0,导致完全无用的浪费工作,做 [mem] |= 0 对于这些进步。

我认为一旦ansans2中的元素达到total,你就会跳出内部循环并执行ans -= sum 1 次通过内部循环。 您不一定要将其重置回 ans = sum(对于该元素)以重做 ORing(设置已经设置的位),因为该内存在缓存中将是冷的。我们真正想要的是将剩余的仍在使用的元素打包到已知位置,并输入只执行 7 个、然后 6 个、然后 5 个元素的循环的其他版本。然后我们只剩下 1 个向量。

这看起来真的很笨拙。对于一个元素结束的更好策略可能是用标量完成该向量中的其余三个,一次一个,然后 运行 剩余的单个 __m256i 向量。如果步幅都在附近,您可能会获得良好的缓存位置。


更便宜的标量,或者可能仍然是 SIMD 但只提取一个位索引

将位索引拆分为一个 qword 索引和一个带有 SIMD 的位掩码,然后分别提取两者对于标量或情况会花费大量微指令:如此之多以至于您不会在每时钟 1 上遇到瓶颈存储吞吐量,即使在我的 scatter/gather 答案中进行了所有优化。 (缓存未命中有时可能会减慢速度,但更少的前端 uops 意味着更大的乱序 window 以找到并行性并保持更多内存操作在运行中。)

如果我们能让编译器生成好的标量代码来拆分位索引,我们可以考虑使用纯标量。或者至少只提取位索引并跳过 SIMD shift/mask 东西。

太糟糕了,标量内存目标 bts 速度不快。 bts [rdi], rax 会在位串中设置该位,即使它在 [rdi] 选择的双字之外。 (这种疯狂的 CISC 行为是 为什么 它并不快,但是!就像 Skylake 上的 10 微码。)

不过,纯标量可能并不理想。我在玩 with this on Godbolt:

#include <immintrin.h>
#include <stdint.h>
#include <algorithm>

// Sieve the bits in array sieveX for later use
void sieveFactors(uint64_t *sieveX64, unsigned cur1, unsigned cur2, unsigned factor1, unsigned factor2)
{
    const uint64_t totalX = 5000000;
#ifdef USE_AVX2
//...
#else
     //uint64_t cur = 58;
     //uint64_t cur2 = 142;
     //uint64_t factor = 67;
     uint32_t *sieveX = (uint32_t*)sieveX64;

    if (cur1 > cur2) {
        // TODO: if factors can be different, properly check which will end first
        std::swap(cur1, cur2);
        std::swap(factor1, factor2);
    }
    // factor1 = factor2;  // is this always true?

    while (cur2 < totalX) {
         sieveX[cur1 >> 5] |= (1U << (cur1 & 0x1f));
         sieveX[cur2 >> 5] |= (1U << (cur2 & 0x1f));
         cur1 += factor1;
         cur2 += factor2;
    }
    while (cur1 < totalX) {
         sieveX[cur1 >> 5] |= (1U << (cur1 & 0x1f));
         cur1 += factor1;
    }
#endif
}

请注意我是如何用排序 cur1、cur2 替换你的外部 if() 以在循环之间进行选择的。

GCC 和 clang 将 1 放在循环外的寄存器中,并在循环内使用 shlx r9d, ecx, esi 在单个 uop 中执行 1U << (cur1 & 0x1f) 而不会破坏 1. (MSVC 使用 load / BTS / store,但是它有很多 mov 指令很笨重。我不知道如何告诉 MSVC 它允许使用 BMI2。)

如果 or [mem], reg 的索引寻址模式不需要额外的 uop,那就太好了。

问题是你需要一个 shr reg, 5 在那里的某个地方,这是破坏性的。将 5 放入寄存器并使用它来复制+移位位索引将是加载/BTS/存储的理想设置,但编译器似乎不知道这种优化。

最佳(?)标量分割和位索引的使用

   mov   ecx, 5    ; outside the loop

.loop:
    ; ESI is the bit-index.
    ; Could be pure scalar, or could come from an extract of ans directly

    shrx  edx, esi, ecx           ; EDX = ESI>>5 = dword index
    mov   eax, [rdi + rdx*4]
    bts   eax, esi                ; set esi % 32 in EAX
    mov   [rdi + rdx*4]


    more unrolled iterations

    ; add   esi, r10d               ; ans += factor if we're doing scalar

    ...
    cmp/jb .loop

因此在 GP 寄存器中给定一个位索引,即 4 微指令来设置内存中的位。请注意,加载和存储都使用 mov,因此索引寻址模式对 Haswell 和更高版本没有影响。

但我认为,使用 shlx / shr / or [mem], reg,我能让编译器生成的最好结果是 5。 (使用索引寻址模式,or 是 3 微指令而不是 2。)

我认为如果您愿意使用手写 asm,您可以使用这个标量并完全放弃 SIMD 来更快。冲突从来都不是正确性问题。

也许你甚至可以让编译器发出类似的东西,但即使是每个展开的 RMW 一个额外的 uop 也是一个大问题。