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

Chi*_*ipK 7 algorithm optimization performance simd avx2

我正在编写一个程序来检测素数.其中一部分是筛选可能的候选人.我写了一个相当快的程序,但我想我会看到是否有人有更好的想法.我的程序可以使用一些快速收集和分散指令,但我只限于AVX2硬件用于x86架构(我知道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
}
Run Code Online (Sandbox Code Playgroud)

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

Pet*_*des 7

IDK为什么使用相同cur[8]数组的不同部分作为索引和值; 它使得源更难以理解,因为只有一个真正的数组.另一种只是将矢量反射到标量.

看起来你只是去向量 - >标量,而不是将标量插入向量中.并且循环中的任何内容都不依赖于任何数据sieveX[]; 我不熟悉你的筛选算法,但我想这一点是在内存中创建数据供以后使用.


AVX2已经收集(不是分散),但它们只能在Skylake和更新版本上快速播出.他们对Broadwell很好,对Haswell来说很慢,而在AMD上也很慢.(就像Ryzen的每12个时钟一样vpgatherqq).请参阅http://agner.org/optimize/以及x86标记wiki其他性能链接.

英特尔的优化手册中有一小部分关于手动收集/分散(使用插入/提取或movhps)与硬件指令,可能值得一读.在这种情况下,索引是运行时变量(不是一个恒定的步幅或其他东西),我认为Skylake可以从AVX2收集指令中受益.

请参阅英特尔的内在函数指南,以查找asm指令的内在函数movhps.我只是在谈论你想让你的编译器发出什么,因为这是重要的,asm助记符的类型更短,不需要转换.您必须知道asm助记符才能在Agner Fog的指令表中查找它们,或者从自动向量化中读取编译器输出,因此我通常在asm中进行思考,然后将其转换为内在函数.


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

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

    (使用标量常量的内存源可能并不坏,如果它们不适合32位立即数,并且如果你没有每个时钟2个内存操作瓶颈.内存目标or指令将使用索引寻址模式,因此,Haswell及以后的端口7上的专用存储器AGU无法使用.因此AGU吞吐量可能成为瓶颈.)

    以标量的形式提取矢量的所有4个元素比4x标量add或移位指令更昂贵,但是你做的工作比这更多.仍然,BMI2为1 uops可变计数转换(而不是Intel上的3),它可能并不可怕.我认为我们可以用SIMD做得更好,特别是仔细调整.

  • 像你现在所做的那样将标量和值提取到标量,所以OR into sieveX[]是纯粹的标量.即使两个或多个索引相同也可以工作.

    这会花费你每ymm向量大约7 uops - >使用提取ALU指令的4x标量寄存器,或使用存储/重载的5 uops(值得考虑编译器,可能是4个向量提取中的一个或两个,因为这个代码可能不会' t管理加载/存储端口吞吐量的瓶颈.)但是,如果编译器将C源中的存储/重新加载变为shuffle/extract指令,除非可能使用,否则不能轻易地覆盖其策略volatile.而顺便说一句,您希望使用它alignas(32) cur[8]来确保实际的矢量存储不会跨越缓存行边界.

    or [rdi + rax*8], rdx(采用索引寻址模式防止完全微融合)在现代Intel CPU(Haswell及更高版本)上是3 uops. 我们可以通过缩放+使用SIMD添加到阵列基地址来避免索引寻址模式(使其为前端2 uop):例如,srli通过3而不是6,屏蔽低3位(vpand),并vpaddq使用set1_epi64(sieveX).因此,每个索引向量需要2个额外的SIMD指令才能在SnB系列上节省4个uop.(您将提取uint64_t*指针元素而不是uint64_t索引.或者如果sieveX可以是32位绝对地址1,您可以跳过vpaddq并提取已经缩放的索引以获得相同的增益.)

    它还可以使存储地址微操作在端口7上运行(Haswell和更高版本) ; port7上的简单AGU只能处理非索引寻址模式.(这使得使用存储+重新加载将标量值提取到标量更具吸引力.您希望提取索引的延迟更低,因为直到memory-dst的加载部分or完成之后才需要这些值.)它确实意味着更多未融合的域uops对于调度程序/执行单元,但很可能值得权衡.

    这不是其他AVX2 CPU(挖掘机/ Ryzen或Xeon Phi)的胜利; 只有SnB系列对索引寻址模式有前端成本和执行端口限制.

  • 提取索引,手动收集到带有vmovq/ vmovhps用于SIMD 的向量vpor,然后用vmovq/ 向后散射vmovhps.

    就像硬件收集/分散一样,正确性要求所有索引都是唯一的,所以你需要使用上述选项之一,直到你的算法达到这一点.(向量冲突检测+回退不值得花费,而只是总是提取到标量:AVX2中冲突检测的后备实现).

    选择性异或与AVX2指令的列表的元素为一个内在的版本.(我知道我最近用手动收集/散布写了一个答案,但花了我一段时间才发现它!)在这种情况下我只使用128位向量,因为没有任何额外的SIMD工作来证明额外的vinserti128/ vextracti128.

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

    前者具有较少的port5压力,并且具有更好的指令级并行性:两个128位半部分是不相互耦合的独立依赖链,因此存储/重载瓶颈(和缓存未命中)会影响较少的依赖性uop,允许乱序执行,以便在等待时继续处理更多内容.

    在256b向量中进行地址计算并提取指针而不是索引会使vmovhps英特尔的负载更便宜(索引负载不能保持微融合到vmovhps2).请参阅上一个要点.但是vmovq加载/存储总是单个uop,并且vmovhps索引存储可以在Haswell以及之后保持微融合,因此它对于前端吞吐量来说是收支平衡的,而在AMD或KNL上则更差.它还意味着调度程序/执行单元的更多未融合域uop,这看起来更像是一个潜在的瓶颈而不是port2/3 AGU的压力.唯一的优势是商店地址微博可以在端口7上运行,减轻了一些压力.

AVX2为我们提供了一个新选择:

  • AVX2 vpgatherqq用于gather(_mm256_i64gather_epi64(sieveX, srli_result, 8)),然后提取索引并手动分散. 因此,除了用AVX2硬件集合替换手动集合之外,它与手动集合/手动分散完全相同.(两个128位集合的成本超过一个256位,因此您需要将指令级并行性命中并收集到一个256位寄存器中).

    可能是Skylake的胜利(其中vpgatherqq ymm4 uops/4c吞吐量,加上1 uop的设置),但甚至没有Broadwell(9 uops,每6c吞吐量一个),绝对不是Haswell(22 uops/9c吞吐量).无论如何,您确实需要标量寄存器中的索引,因此您只需保存手动收集部分工作.那很便宜.


Skylake上每项策略的总成本

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

因此,您需要经常将矢量提取到标量的SIMD计算不如您需要经常将标量计算结果插入矢量中那么糟糕.另请参阅从GP regs加载xmm,但这是在讨论从GP regs开始的数据,而不是内存.

  • 提取/标量OR:总计= 24 uops = 6个前端吞吐量周期.

    • vpaddq + vpand地址计算(Skylake上的端口0/1/5为2 uops)
    • 2x vextracti128(端口5为2 uops)
    • 4x vmovq(4 p0)
    • 4x vpextrq(8:4p0 4p5)
    • 4x or [r], r(每个4x2 = 8个前端uop.后端:4p0156 4p23(加载)4p237(store-addres)4p4(存储数据)).非索引寻址模式.

    p5的总计= 6 uops,几乎不适合.如果您可以让编译器执行此操作,则存储/重新加载数据提取看起来很合理.(但编译器通常不会对管道进行足够详细的建模,以在同一循环中使用混合策略来平衡端口压力.)

  • 手动采集/散射:20 uops,5个周期的前端吞吐量(Haswell/BDW/Skylake).Ryzen也很好.

    • (可选,可能不值得):vpaddq + vpand地址计算(Skylake上端口0/1/5的2 uops)如果你可以使用非VEX movhps用于1-uop微融合索引负载,请跳过这些.(但后来p237商店变成了p23).
    • vextracti128指针(端口5为1 uop)
    • 2x vmovq提取物(2p0)
    • 2x vpextrq(4 = 2p0 2p5)
    • 2x vmovq load(2p23)
    • 2x vmovhps xmm, xmm, [r]非索引负载(2个前端uops微融合:2p23 + 2p5)

    • vextracti128拆分数据(p5)

    • 2x vpor xmm(2p015)
    • 2x vmovq商店(2x 1微融合uop,2p237 + 2p4)
    • 2x vmovhps商店(2x 1微融合uop,2p237 + 2p4)

    端口瓶颈:4 p0和4 p5在5个循环中很舒适,特别是当你将它与你的循环混合时,它可以在端口1上运行几个uop.在Haswell paddq上只有p15(不是p015),并且移位只有p0(不是P01)._mm256_sllv_epi64在Skylake上AVX2 是1 uop(p01),但在Haswell上它是3 uops = 2p0 + p5.因此,Haswell可能更接近此循环的p0或p5瓶颈,在这种情况下,您可能希望查看一个索引向量的存储/重新加载提取策略.

    跳过SIMD地址计算可能是好的,因为除非您使用存储/重新加载提取,否则AGU压力看起来不是问题.它意味着uop缓存中的指令/更小的代码大小和更少的uop.(在解码器/ uop缓存之后才会发生解除分层,所以你仍然可以在前端的早期部分受益于微融合,而不是在问题瓶颈.)

  • Skylake AVX2采集/手动分散:总计= 18 uops,4.5个周期的前端吞吐量. (更糟糕的是任何早期的uarch或AMD).

    • vextracti128指数(端口5为1 uop)
    • 2x vmovq提取物(2p0)
    • 2x vpextrq(4 = 2p0 2p5)

    • vpcmpeqd ymm0,ymm0,ymm0vpgatherqq(p015)创建一个全屏掩码

    • vpgatherqq ymm1, [rdi + ymm2*8], ymm0 某些端口有4个uop.

    • vpor ymm (P015)

    • OR结果的vextracti128(p5)
    • 2x vmovq商店(2x 1微熔融uop,2p23 + 2p4).注意注意port7,我们正在使用索引存储.
    • 2x vmovhps存储(2x 1微融合uop,2p23 + 2p4).

因此,即使有最佳的吞吐量选择,我们仍然只能在每4.5个周期管理4个负载/ 4个存储,而且没有考虑SIMD在环路中的工作,这会花费一些前端吞吐量. 所以我们并没有接近AGU吞吐量的瓶颈而不必担心使用端口7.

我们可以考虑其中一个提取(如果我们是编译器)的存储/重新加载,用5 uops存储/ 4x加载替换7 uop 5指令vextracti128/2x vmovq/2x vpextrq序列.


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

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

你肯定想要一个单独的循环,它根本不检查冲突以利用它.即使你有AVX512,Skylake vpconflictq也是微代码并且速度不快.(KNL有单用途,vpconflictq但完全避免它仍然更快).

我将把它留给你(或一个单独的问题)如何在你完成冲突时确定并且可以留下解释这种可能性的循环.

您可能需要提取索引+数据策略,但可能存在冲突.SIMD冲突检查是可能的,但它并不便宜,32位元素为11 uop: AVX2中冲突检测的后备实现.一个qword版本显然比dword便宜得多(更少的shuffle和所有对比所有对比),但你可能仍然只想每10次迭代左右你的提取循环.

从最佳标量或版本到最佳聚集版本没有巨大的加速(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();
Run Code Online (Sandbox Code Playgroud)

这应该编译成一个dec / je只在未采取的情况下花费1 uop.

总共进行9次额外迭代总计稍慢的循环要比进行数千次额外的昂贵冲突检查好得多.


脚注1:

如果sieveX是静态的并且您在Linux(而不是MacOS)上构建非PIC代码,则其地址将适合disp32作为[reg+disp32]寻址模式的一部分.在那种情况下你可以省略vpaddq.但是让一个编译器将一个处理uint64_t为已经缩放的数组索引(清除其低位)将是丑陋的.大概要投sieveXuintptr_t和补充,然后投退.

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

脚注2:有趣的事实:非VEX movhps负载可以在Haswell上保持微融合.它不会导致Skylake上的SSE/AVX停顿,但是你不会让编译器在AVX2功能的中间发出非VEX版本.

但IACA(英特尔的静态分析工具)却错了.:( 什么是IACA以及如何使用它?.

这基本上是一个错过优化-mtune=skylake,但它停止在Haswell:为什么这个SSE代码在没有Skylake上的VZEROUPPER时要慢6倍?.

Skylake上的"惩罚A"(执行带有脏上层的SSE)仅仅是对该一个寄存器的错误依赖.(以及用于指令的合并uop,否则这些指令只能写入,但movhps已经是对其目标的读取 - 修改 - 写入.)我在Skylake上使用Linux对此进行了测试perf以计算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
Run Code Online (Sandbox Code Playgroud)

该循环在Skylake(i7-6700k)上每次迭代运行约1.25个循环,最大化每个时钟4个uop的前端吞吐量.5个融合域uops(uops_issued.any),6个未融合域uops(uops_executed.thread).因此,在movhps没有任何SSE/AVX问题的情况下,微型融合肯定会发生.

改变它以将其vmovhps xmm0, xmm0, [rdi+rbx*8+8]减慢到每次迭代1.50个循环,现在6个融合域,但仍然是相同的6个未融合域uops.

如果上半部分ymm0movhps xmm0, [mem]运行时变脏,则没有额外的uop .我通过评论来测试vmovq.但是vmovq改为movq 确实会导致额外的uop:movq变成一个微融合的加载+合并,它取代了低64位(并且仍然将xmm0的高64位归零,所以它并不完全movlps).


另请注意,pinsrq xmm0, [mem], 1即使没有VEX ,也无法进行微型保险丝.但是对于VEX,你应该更喜欢vmovhps代码大小的原因.

您的编译器可能想要"优化" movhps整数数据的内在函数vpinsrq,但是,我没有检查.


Pet*_*des 5

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

然后将其展开两个,并使用mod1 << ans和将ans和ans2并行运行mod3 << ans2。注释您的代码,并使用英文文本解释大图中发生的事情!这只是Eratosthenes常规筛网的位设置循环的非常复杂的实现。(因此,如果问题首先是这样说,那就太好了。)

与多个开始/跨步并行展开是一个非常好的优化,因此通常在L1d中仍处于高温状态的高速缓存行中设置多个位。 一次针对较少的不同因素进行缓存阻塞具有类似的优势。在移至下一个之前,针对多个因素(跨步)重复遍历相同的8kiB或16kiB内存块。为2个不同的步幅分别展开4个偏移量可能是创建更多ILP的好方法。

但是,并行运行的步幅越多,第一次触摸新的缓存行时就越慢。(提供高速缓存/ TLB预取空间可避免最初的停顿)。因此,缓存阻止并不能消除多步执行的所有优势。


迈步<256的可能特殊情况

单个256位向量加载/ VPOR /存储可以设置多个位。诀窍是创建一个矢量常量或一组矢量常量,并将位放在正确的位置。不过,重复模式有点LCM(256, bit_stride)长。例如,stride = 3将以3个向量长的模式重复。除非有更聪明的方法,否则这很快就变得无法用于奇数/素数的大步走了:(

64位标量很有趣,因为可以使用按位旋转来创建一系列模式,但是在SnB系列CPU上进行可变计数旋转的成本为2微秒。

您可能会做更多的事情。也许未对齐的负载可能会有所帮助。

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


通过更窄的负载/存储减少冲突

当您只设置一个位时,就不必加载/修改/存储64位块。您的RMW操作范围越窄,位索引就可以越接近而不会发生冲突。

(但是您在同一位置没有长循环承载的dep链;您将继续前进,直到OoO exec停顿,等待长链末端的重新加载。因此,如果冲突不是正确性问题,那就是不太可能在这里产生很大的性能差异。与位图直方图或附近容易发生长串重复命中的东西不同。

32位元素将是显而易见的选择。x86可以有效地将dword加载/存储到SIMD寄存器以及从中存储dword和标量。(标量字节操作也很有效,但是来自SIMD寄存器的字节存储总是需要使用多个pextrb。)

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

使用64位SIMD元素的主要优点是,如果您要计算的是指针而不是索引。如果您可以限制sieveX为32位,则仍然可以执行此操作。例如mmap(..., MAP_32BIT|MAP_ANONYMOUS, ...)在Linux上分配。 这是假设您不需要2 ^ 32位(512MiB)的筛子空间,因此位索引始终适合32位元素。 如果不是这种情况,那么到目前为止,您仍然可以使用32位元素向量,然后将当前循环用于较高部分。

如果您使用32位SIMD元素而不将sieveX限制为32位点指针,则必须放弃使用SIMD指针计算,而只提取位索引,或者仍然将SIMD拆分为idx/ bit并同时提取。

(对于32位元素,基于存储/重新加载的SIMD->标量策略看起来更具吸引力,但是在C语言中,这主要取决于编译器。)

如果您手动收集到32位元素,你不能使用movhps。您必须将pinsrd/ pextrd用于高3个元素,而那些绝不会微熔丝/总是需要SnB系列上的port5 uop。(与movhps纯粹的商店不同)。但这意味着vpinsrd仍具有索引寻址模式的2微码。您仍然可以vmovhps将元素2用于(然后用覆盖向量的顶部双字vpinsrd);未对齐的载荷很便宜,可以覆盖下一个元素。但是您不能做movhps商店,那真的很不错。


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

显然,您有时会将此元素与mod1mod3being一起使用0,从而导致完全无用的工作浪费,[mem] |= 0为这些进步做着努力。

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

这似乎很笨拙。一个元素到达终点的更好策略可能是用标量完成该向量中的其余三个,一次一次,然后运行剩余的单个__m256i向量。如果所有步幅都在附近,则可能会获得良好的缓存位置。


标量便宜,或者可能仍然是SIMD,但仅提取位索引

对于标量“或”情况,使用SIMD将位索引分为qword索引和位掩码,然后分别提取它们,这会花费大量的运算量:如此之多,即使您有1个时钟存储吞吐量,也不会出现瓶颈我的分散/聚集答案中的所有优化。(有时,高速缓存未命中可能会减慢此速度,但是较少的前端uops意味着较大的乱序窗口可发现并行性并保持运行中的更多内存操作。)

如果我们可以使编译器制作好的标量代码来拆分位索引,则可以考虑使用纯标量。或者至少仅提取位索引并跳过SIMD移位/掩码内容。

太糟糕了,标量存储目标bts不是很快。 bts [rdi], rax会在位字符串中设置该位,即使该位不在所选择的dword之外[rdi]。(这种疯狂的CISC行为就是为什么它不快的原因,就像Skylake上的10微秒。)

但是,纯标量可能并不理想。我在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
}
Run Code Online (Sandbox Code Playgroud)

请注意,我是如何替换外部的if()以对cur1,cur2进行排序的。

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

如果的索引寻址模式or [mem], reg不花费额外的成本,那就太好了。

问题是您需要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
Run Code Online (Sandbox Code Playgroud)

因此,给定GP寄存器中的位索引,设置内存中的位为4微秒。请注意,load和store都使用mov,因此索引寻址模式对Haswell及更高版本没有影响。

但是,我认为使用shlx / shr /可以使编译器获得的最佳性能是5 or [mem], reg。(在索引寻址模式下,该or值为3而不是2。)

我认为,如果您愿意使用手写的asm,则可以使用此标量和完全放弃SIMD来更快。冲突绝不是正确的问题。

也许您甚至可以让编译器发出可比的结果,但是每个展开的RMW甚至一个额外的uop都是一件大事。