首页 文章

用于测试Collatz猜想的C代码比手写程序集更快 - 为什么?

提问于
浏览
754

我在Project Euler Q14中为汇编和C编写了这两个解决方案 . 它们是用于测试Collatz conjecture的相同蛮力方法 . 组装解决方案与组装

nasm -felf64 p14.asm && gcc p14.o -o p14

C编译了

g++ p14.cpp -o p14

大会, p14.asm

section .data
    fmt db "%d", 10, 0

global main
extern printf

section .text

main:
    mov rcx, 1000000
    xor rdi, rdi        ; max i
    xor rsi, rsi        ; i

l1:
    dec rcx
    xor r10, r10        ; count
    mov rax, rcx

l2:
    test rax, 1
    jpe even

    mov rbx, 3
    mul rbx
    inc rax
    jmp c1

even:
    mov rbx, 2
    xor rdx, rdx
    div rbx

c1:
    inc r10
    cmp rax, 1
    jne l2

    cmp rdi, r10
    cmovl rdi, r10
    cmovl rsi, rcx

    cmp rcx, 2
    jne l1

    mov rdi, fmt
    xor rax, rax
    call printf
    ret

C,p14.cpp

#include <iostream>

using namespace std;

int sequence(long n) {
    int count = 1;
    while (n != 1) {
        if (n % 2 == 0)
            n /= 2;
        else
            n = n*3 + 1;

        ++count;
    }

    return count;
}

int main() {
    int max = 0, maxi;
    for (int i = 999999; i > 0; --i) {
        int s = sequence(i);
        if (s > max) {
            max = s;
            maxi = i;
        }
    }

    cout << maxi << endl;
}

我知道编译器优化以提高速度和一切,但我没有看到很多方法来进一步优化我的程序集解决方案(以编程方式而不是数学方式) .

C代码每个项都有模数,每个偶数项都有一个除法,其中汇编每个偶数项只有一个除法 .

但是组装平均比C解决方案长1秒 . 为什么是这样?我主要是好奇地问 .

执行时间

我的系统:1.4 GHz Linux 1.4 GHz Intel Celeron 2955U(Haswell微体系结构) .

11 回答

  • 20

    即使不考虑装配,最明显的原因是 /= 2 可能被优化为 >>=1 并且许多处理器具有非常快速的移位操作 . 但即使处理器没有移位操作,整数除法也比浮点除法快 .

    Edit: 您的milage可能会因上述"integer division is faster than floating point division"声明而异 . 下面的评论表明,现代处理器优先考虑优化fp除以整数除法 . 因此,如果有人正在寻找这个线程问题所要求的加速的最可能原因,那么编译器优化 /=2 作为 >>=1 将是最好看的第一名 .


    在一个不相关的注释中,如果 n 是奇数,则表达式 n*3+1 将始终为偶数 . 所以没有必要检查 . 您可以将该分支更改为

    {
       n = (n*3+1) >> 1;
       count += 2;
    }
    

    那么整个陈述就是这样

    if (n & 1)
    {
        n = (n*3 + 1) >> 1;
        count += 2;
    }
    else
    {
        n >>= 1;
        ++count;
    }
    
  • 4

    作为一般性答案,并非专门针对此任务:在许多情况下,您可以通过在高级别上进行改进来显着加快任何程序 . 就像计算一次而不是多次数据一样,完全避免不必要的工作,以最好的方式使用缓存,等等 . 这些东西在高级语言中更容易实现 .

    编写汇编程序代码,可以改进优化编译器的功能,但这很难 . 一旦完成,您的代码就很难修改,因此添加算法改进要困难得多 . 有时处理器具有您无法从高级语言使用的功能,内联汇编在这些情况下通常很有用,并且仍然允许您使用高级语言 .

    在Euler问题中,大多数时候你通过构建一些东西来成功,找到它缓慢的原因,更好地构建一些东西,找到它为什么慢,等等 . 使用汇编程序非常非常困难 . 以一半可能的速度运行的更好的算法通常会在全速运行时击败更差的算法,并且在汇编程序中获得全速并非易事 .

  • 5

    简单的答案:

    • 做MOV RBX,3和MUL RBX很贵;只需添加RBX,RBX两次

    • ADD 1可能比INC更快

    • MOV 2和DIV非常昂贵;向右转

    • 64位代码通常明显慢于32位代码,并且对齐问题更复杂;对于像这样的小程序,你必须打包它们,这样你就可以进行并行计算,以便有机会比32位代码更快

    如果为C程序生成程序集列表,则可以看到它与程序集的不同之处 .

  • -2

    来自评论:

    但是,这段代码永远不会停止(因为整数溢出)!?! Yves Daoust

    对于许多数字,它不会溢出 .

    如果它会溢出 - 对于那些不幸的初始种子之一,溢出的数字很可能会收敛到1而没有另一个溢出 .

    这仍然是一个有趣的问题,是否有一些溢出循环的种子数?

    任何简单的最终收敛序列都以两个值的幂开始(显而易见?) .

    2 ^ 64将溢出为零,这是根据算法未定义的无限循环(仅以1结尾),但由于 shr rax 产生ZF = 1,因此最佳的答案解决方案将完成 .

    我们可以 生产环境 2 ^ 64吗?如果起始编号为 0x5555555555555555 ,则为奇数,下一个编号为3n 1,即 0xFFFFFFFFFFFFFFFF + 1 = 0 . 理论上处于未定义的算法状态,但johnfound的优化答案将通过退出ZF = 1来恢复 . Peter Cordes的 cmp rax,1 will end in infinite loop (QED变体1,"cheapo"到undefined 0 数字) .

    如何更复杂的数字,这将创建没有 0 的循环?坦率地说,我需要担心原始系列的无限循环,只有溢出才能阻碍我们 .

    所以我只是将几个数字放入表中并查看了8位截断数字 .

    有三个值溢出到 02271708585 直接转到 0 ,其他两个转向 85 ) .

    但是没有创造循环溢出种子的 Value .

    有趣的是,我做了一个检查,这是第一个遭受8位截断的数字,已经 27 受到影响!它在正确的非截断系列中达到值 9232 (第12步中第一个截断值为 322 ),并且非截断方式中任何2-255输入数字达到的最大值为 13120 (对于 255 本身),收敛到 1 的最大步数约为 128 (-2,不确定"1"是否计数等等) .

    有趣的是(对我来说)数字 9232 对于许多其他源数字来说是最大的,它有什么特别之处呢? :-O 9232 = 0x2410 ......嗯......不知道 .

    遗憾的是我无法深入掌握这个系列,它为什么会聚并且将它们截断为k位有什么含义,但是在 cmp number,1 终止条件下,当然可以将算法置于无限循环中,特定输入值结束为截断后 0 .

    但是对于8位情况溢出的值 27 有点警报,这看起来如果你计算达到值 1 的步数,你将得到整个k位整数集中的大多数数字的错误结果 . 对于8位整数,256个中的146个数字通过截断影响了系列(其中一些可能仍然意外地达到了正确的步数,我懒得检查) .

  • 3

    对于Collatz问题,您可以通过缓存"tails"来显着提升性能 . 这是时间/记忆的权衡 . 请参阅:memoization(https://en.wikipedia.org/wiki/Memoization) . 您还可以研究动态编程解决方案以进行其他时间/内存权衡 .

    示例python实现:

    import sys
    
    inner_loop = 0
    
    def collatz_sequence(N, cache):
        global inner_loop
    
        l = [ ]
        stop = False
        n = N
    
        tails = [ ]
    
        while not stop:
            inner_loop += 1
            tmp = n
            l.append(n)
            if n <= 1:
                stop = True  
            elif n in cache:
                stop = True
            elif n % 2:
                n = 3*n + 1
            else:
                n = n // 2
            tails.append((tmp, len(l)))
    
        for key, offset in tails:
            if not key in cache:
                cache[key] = l[offset:]
    
        return l
    
    def gen_sequence(l, cache):
        for elem in l:
            yield elem
            if elem in cache:
                yield from gen_sequence(cache[elem], cache)
                raise StopIteration
    
    if __name__ == "__main__":
        le_cache = {}
    
        for n in range(1, 4711, 5):
            l = collatz_sequence(n, le_cache)
            print("{}: {}".format(n, len(list(gen_sequence(l, le_cache)))))
    
        print("inner_loop = {}".format(inner_loop))
    
  • 1783

    声称C编译器可以生成比合格的汇编语言程序员更优的代码是一个非常糟糕的错误 . 特别是在这种情况下 . 人总是可以使编码器能够更好地编写代码,并且这种特殊情况很好地说明了这种说法 .

    您所看到的时序差异是因为问题中的汇编代码在内部循环中远非最佳 .

    (以下代码为32位,但可以很容易地转换为64位)

    例如,序列函数可以优化为仅5条指令:

    .seq:
            inc     esi                 ; counter
            lea     edx, [3*eax+1]      ; edx = 3*n+1
            shr     eax, 1              ; eax = n/2
            cmovc   eax, edx            ; if CF eax = edx
            jnz     .seq                ; jmp if n<>1
    

    整个代码看起来像:

    include "%lib%/freshlib.inc"
    @BinaryType console, compact
    options.DebugMode = 1
    include "%lib%/freshlib.asm"
    
    start:
            InitializeAll
            mov ecx, 999999
            xor edi, edi        ; max
            xor ebx, ebx        ; max i
    
        .main_loop:
    
            xor     esi, esi
            mov     eax, ecx
    
        .seq:
            inc     esi                 ; counter
            lea     edx, [3*eax+1]      ; edx = 3*n+1
            shr     eax, 1              ; eax = n/2
            cmovc   eax, edx            ; if CF eax = edx
            jnz     .seq                ; jmp if n<>1
    
            cmp     edi, esi
            cmovb   edi, esi
            cmovb   ebx, ecx
    
            dec     ecx
            jnz     .main_loop
    
            OutputValue "Max sequence: ", edi, 10, -1
            OutputValue "Max index: ", ebx, 10, -1
    
            FinalizeAll
            stdcall TerminateAll, 0
    

    为了编译此代码,需要FreshLib .

    在我的测试中,(1 GHz AMD A4-1200处理器),上述代码比问题中的C代码快了大约四倍(使用 -O0 编译时:430 ms与1900 ms),速度提高了两倍多(当使用 -O3 编译C代码时,430 ms与830 ms) .

    两个程序的输出相同:i = 837799时最大序列= 525 .

  • 3

    你没有发布编译器生成的代码,所以这里有一些猜测,但即使没有看到它,也可以这样说:

    test rax, 1
    jpe even
    

    ...有50%的机会错误预测分支机构,这将会变得昂贵 .

    编译器几乎肯定会进行两种计算(由于div / mod的延迟相当长,所以成本可忽略不计,所以乘法加法是"free")并跟随CMOV . 当然,这被错误预测的可能性为零 .

  • 4

    如果您认为64位DIV指令是除以2的好方法,那么难怪编译器的asm输出会超过您的手写代码,即使使用 -O0 (编译速度快,没有额外的优化,并且存储/重新加载到内存之后) /在每个C语句之前,所以调试器可以修改变量) .

    请参阅Agner Fog's Optimizing Assembly guide以了解如何编写高效的asm . 他还提供了指令表和微指南,以获取特定CPU的具体细节 . 有关更多perf链接,另请参阅x86标记wiki .

    另请参阅这个关于使用手写asm击败编译器的更一般的问题:Is inline assembly language slower than native C++ code? . TL:DR:是的,如果你做错了(比如这个问题) .

    通常你很好让编译器做它的事情,特别是如果你 try to write C++ that can compile efficiently . 另见is assembly faster than compiled languages? . 其中一个答案链接到these neat slides,显示各种C编译器如何使用很酷的技巧优化一些非常简单的功能 .


    even:
        mov rbx, 2
        xor rdx, rdx
        div rbx
    

    在Intel Haswell上, div r64 是36 uop,具有 latency of 32-96 cycles ,每21-74个周期的吞吐量为1 . (加上2个uop来设置RBX和零RDX,但是乱序执行可以提前运行) . High-uop-count instructions like DIV are microcoded, which can also cause front-end bottlenecks.在这种情况下,延迟是最相关的因素,因为它是循环携带的依赖链的一部分 .

    shr rax, 1 does the same unsigned division: It's 1 uop, with 1c latency ,每个时钟周期可以运行2个 .

    相比之下,32位除法速度更快,但对于换档仍然很可怕 . idiv r32 是9 uops,22-29c延迟,每8-11c一个Haswell的吞吐量 .


    As you can see from looking at gcc's -O0 asm output (Godbolt compiler explorer), it only uses shifts instructions . clang -O0 确实像你想象的那样天真地编译,即使使用64位IDIV两次也是如此 . (优化时,编译器在源进行除法时使用IDIV的两个输出,并使用相同的操作数使用模数,如果它们使用IDIV的话)

    海湾合作委员会没有完全天真的模式; it always transforms through GIMPLE, which means some "optimizations" can't be disabled . 这包括识别逐次除法和使用移位(2的幂)或a fixed-point multiplicative inverse(2的非幂)来避免IDIV(参见上面的godbolt链接中的 div_by_13 ) .

    gcc -Os (针对大小进行优化)确实将IDIV用于非2次幂除法,不幸的是,即使在乘法逆码仅略大但速度更快的情况下也是如此 .


    帮助编译器

    (本案例摘要:使用 uint64_t n

    首先,看一下优化的编译器输出是很有意思的 . ( -O3 ) . -O0 speed is basically meaningless.

    看看你的asm输出(在Godbolt上,或者参见How to remove "noise" from GCC/clang assembly output?) . 当编译器首先没有制作最佳代码时: Writing your C/C++ source in a way that guides the compiler into making better code is usually the best approach . 你必须知道asm,并知道什么是有效的,但你间接地应用这些知识 . 编译器也是一个很好的想法来源:有时clang会做一些很酷的事情,你可以手动执行gcc做同样的事情:请参阅this answer以及我在@ Veedrac的代码中对未展开的循环做了什么 . )

    这种方法是可移植的,并且在未来20年内,一些未来的编译器可以将其编译为在未来硬件上有效的任何内容(x86或不是x86),可能使用新的ISA扩展或自动矢量化 . 从15年前手写的x86-64 asm通常不会为Skylake进行最佳调整 . 例如当时并不存在比较和分支宏观融合 . What's optimal now for hand-crafted asm for one microarchitecture might not be optimal for other current and future CPUs. Comments on @johnfound's answer讨论AMD Bulldozer和Intel Haswell之间的主要差异,这对这段代码有很大的影响 . 但理论上, g++ -O3 -march=bdver3g++ -O3 -march=skylake 会做正确的事情 . (或 -march=native . )或 -mtune=... 只是调整,而不使用其他CPU可能不支持的指令 .

    我的感觉是指导编译器asm对你关心的当前CPU有好处,对于未来的编译器来说应该不是问题 . 他们希望在寻找转换代码的方法方面比现有的编译器更好,并且可以找到适用于未来CPU的方法 . 无论如何,未来的x86对于当前x86上的任何好处都可能不会很糟糕,并且未来的编译器会在实现类似C源代码的数据移动时避免任何特定于asm的陷阱,如果它没有看到更好的东西 .

    手写asm是优化器的黑盒子,因此当内联使输入成为编译时常量时,常量传播不起作用 . 其他优化也会受到影响 . 在使用asm之前阅读https://gcc.gnu.org/wiki/DontUseInlineAsm . (并避免使用MSVC样式的内联asm:输入/输出必须通过内存which adds overhead . )

    In this case :你的 n 有一个带符号的类型,gcc使用SAR / SHR / ADD序列给出正确的舍入 . (对于负输入,IDIV和算术移位"round"不同,请参阅SAR insn set ref manual entry) . (IDK如果gcc尝试并且未能证明 n 不能为负,或者是什么 . 签名溢出是未定义的行为,所以应该可以 . )

    你应该使用 uint64_t n ,所以它只能SHR . 因此它可以移植到 long 只有32位的系统(例如x86-64 Windows) .


    顺便说一句, gcc's optimized asm output looks pretty good (using unsigned long n) :它内联到 main() 的内部循环这样做:

    # from gcc5.4 -O3  plus my comments
    
     # edx= count=1
     # rax= uint64_t n
    
    .L9:                   # do{
        lea    rcx, [rax+1+rax*2]   # rcx = 3*n + 1
        mov    rdi, rax
        shr    rdi         # rdi = n>>1;
        test   al, 1       # set flags based on n%2 (aka n&1)
        mov    rax, rcx
        cmove  rax, rdi    # n= (n%2) ? 3*n+1 : n/2;
        add    edx, 1      # ++count;
        cmp    rax, 1
        jne   .L9          #}while(n!=1)
    
      cmp/branch to update max and maxi, and then do the next n
    

    内循环是无分支的,循环携带依赖链的关键路径是:

    • 3组分LEA(3个循环)

    • cmov(Haswell上2个循环,Broadwell上或之后1c) .

    Total: 5 cycle per iteration, latency bottleneck . 乱序执行与此并行处理其他所有事情(理论上:我没有使用perf计数器进行测试,看它是否真的以5c / iter运行) .

    cmov (由TEST生成)的FLAGS输入比RAX输入(来自LEA-> MOV)生成得更快,因此它不在关键路径上 .

    类似地,产生CMOV的RDI输入的MOV-> SHR不在关键路径上,因为它也比LEA快 . IvyBridge上的MOV以及之后的延迟为零(在寄存器重命名时间处理) . (它仍然需要一个uop和一个插槽,因此它不是免费的,只是零延迟) . LEA dep链中的额外MOV是其他CPU瓶颈的一部分 .

    cmp / jne也不是关键路径的一部分:它不是循环携带的,因为控制依赖性是通过分支预测推测执行来处理的,这与关键路径上的数据依赖性不同 .


    击败编译器

    海湾合作委员会在这里做得很好 . 它可以通过使用inc edx instead of add edx, 1来保存一个代码字节,因为没有人关心P4及其对部分标志修改指令的错误依赖性 .

    它还可以保存所有MOV指令,并且TEST:SHR设置CF =移出的位,因此我们可以使用 cmovc 而不是 test / cmovz .

    ### Hand-optimized version of what gcc does
    .L9:                       #do{
        lea     rcx, [rax+1+rax*2] # rcx = 3*n + 1
        shr     rax, 1         # n>>=1;    CF = n&1 = n%2
        cmovc   rax, rcx       # n= (n&1) ? 3*n+1 : n/2;
        inc     edx            # ++count;
        cmp     rax, 1
        jne     .L9            #}while(n!=1)
    

    请参阅@johnfound 's answer for another clever trick: remove the CMP by branching on SHR' s标志结果以及将其用于CMOV:仅当n为1(或0)时才为零 . (有趣的事实:SHR with count != 1 on Nehalem or earlier causes a stall if you read the flag results . 这就是他们如何让它成为单一的但是,按1转换特殊编码很好 . )

    避免MOV对Haswell(Can x86's MOV really be "free"? Why can't I reproduce this at all?)的延迟没有帮助 . 它确实对像英特尔pre-IvB和AMD Bulldozer系列这样的CPU有很大帮助,其中MOV不是零延迟 . 编译器's wasted MOV instructions do affect the critical path. BD'的复合LEA和CMOV都具有较低的延迟(分别为2c和1c),因此它是延迟的一小部分 . 此外,吞吐量瓶颈成为一个问题,因为它只有两个整数ALU管道 . See @johnfound's answer,他有AMD CPU的计时结果 .

    即使在Haswell上,这个版本可能会有所帮助,避免一些偶然的延迟,其中非关键uop从关键路径上的一个执行端口窃取执行端口,将执行延迟1个周期 . (这称为资源冲突) . 它还保存了一个寄存器,这可能有助于在交错循环中并行执行多个 n 值(见下文) .

    LEA's latency depends on the addressing mode ,在Intel SnB系列CPU上 . 3c为3个组件( [base+idx+const] ,需要两个单独的添加),但只有1c有2个或更少的组件(一个添加) . 有些CPU(如Core2)在一个周期内甚至可以执行3分量LEA,但SnB系列则不行 . 更糟糕的是,Intel SnB-family standardizes latencies so there are no 2c uops,否则3分量LEA将只有2c像Bulldozer . (3组件LEA在AMD上也较慢,只是没有那么多) .

    所以 lea rcx, [rax + rax*2] / inc rcx 在像Haswell这样的Intel SnB系列CPU上只有2c延迟,比 lea rcx, [rax + rax*2 + 1] 快 . 在BD上收支 balancer ,在Core2上更糟糕 . 它确实需要花费额外的uop,这通常不值得节省1c延迟,但延迟是这里的主要瓶颈,Haswell有足够宽的管道来处理额外的uop吞吐量 .

    Neither gcc, icc, nor clang (on godbolt) used SHR's CF output, always using an AND or TEST . 愚蠢的编译器 . :P他们使用详尽的算法来搜索每种可能的做事方式,因为在优化大量内联代码时需要花费太长时间,这是他们最擅长的 . 他们也没有在目标微体系结构中对管道进行建模,至少与IACA或其他静态分析工具的细节不同;他们只是使用一些启发式方法 . )


    Simple loop unrolling won't help ;这个循环瓶颈是循环携带的依赖链的延迟,而不是循环开销/吞吐量 . 这意味着它可以很好地处理超线程(或任何其他类型的SMT),因为CPU有很多时间来交错来自两个线程的指令 . 这意味着在 main 中并行循环,但这很好,因为每个线程只能检查一系列 n 值并产生一对整数 .

    Interleaving by hand within a single thread might be viable, too . 也许并行计算一对数字的序列,因为每个数字只需要几个寄存器,并且它们都可以更新相同的 max / maxi . 这会创建更多instruction-level parallelism .

    诀窍是决定是否等到所有 n 值都达到 1 之后才获得另一对起始 n 值,或者是否突破并获得一个达到结束条件的新起点,而不触及寄存器其他顺序 . 可能它必须有条件地增加其计数器 .


    您甚至可以使用SSE打包比较来有条件地增加 n 尚未达到 1 的向量元素的计数器 . 然后为了隐藏SIMD条件增量实现的更长延迟,您需要保留更多 n 值的向量 . 也许只值得256b矢量(4x uint64_t ) .

    我认为检测 1 "sticky"的最佳策略是屏蔽你添加的所有1的向量以递增计数器 . 因此,在元素中看到 1 之后,增量向量将为零,而= 0为无操作 .

    手工矢量化未经测试的想法

    # starting with YMM0 = [ n_d, n_c, n_b, n_a ]  (64-bit elements)
    # ymm4 = _mm256_set1_epi64x(1):  increment vector
    # ymm5 = all-zeros:  count vector
    
    .inner_loop:
        vpaddq    ymm1, ymm0, xmm0
        vpaddq    ymm1, ymm1, xmm0
        vpaddq    ymm1, ymm1, set1_epi64(1)     # ymm1= 3*n + 1.  Maybe could do this more efficiently?
    
        vprllq    ymm3, ymm0, 63                # shift bit 1 to the sign bit
    
        vpsrlq    ymm0, ymm0, 1                 # n /= 2
    
        # There may be a better way to do this blend, avoiding the bypass delay for an FP blend between integer insns, not sure.  Probably worth it
        vpblendvpd ymm0, ymm0, ymm1, ymm3       # variable blend controlled by the sign bit of each 64-bit element.  I might have the source operands backwards, I always have to look this up.
    
        # ymm0 = updated n  in each element.
    
        vpcmpeqq ymm1, ymm0, set1_epi64(1)
        vpandn   ymm4, ymm1, ymm4         # zero out elements of ymm4 where the compare was true
    
        vpaddq   ymm5, ymm5, ymm4         # count++ in elements where n has never been == 1
    
        vptest   ymm4, ymm4
        jnz  .inner_loop
        # Fall through when all the n values have reached 1 at some point, and our increment vector is all-zero
    
        vextracti128 ymm0, ymm5, 1
        vpmaxq .... crap this doesn't exist
        # Actually just delay doing a horizontal max until the very very end.  But you need some way to record max and maxi.
    

    您可以而且应该使用内在函数来实现它,而不是手写的asm .


    算法/实施改进:

    除了使用更高效的asm实现相同的逻辑之外,还要寻找简化逻辑或避免冗余工作的方法 . 例如memoize以检测序列的常见结尾 . 或者甚至更好,一次查看8个尾随位(gnasher的答案)

    @EOF指出 tzcnt (或 bsf )可用于在一个步骤中进行多次 n/=2 迭代 . 但是,'s probably better than SIMD vectorizing, because no SSE or AVX instruction can do that. It'仍然兼容在不同的整数寄存器中并行执行多个标量 n .

    所以循环可能如下所示:

    goto loop_entry;  // C++ structured like the asm, for illustration only
    do {
       n = n*3 + 1;
      loop_entry:
       shift = _tzcnt_u64(n);
       n >>= shift;
       count += shift;
    } while(n != 1);
    

    这可能会显着减少迭代次数,但在没有BMI2的Intel SnB系列CPU上,可变计数移位速度很慢 . 3次uops,2c延迟 . (它们对FLAGS有输入依赖性,因为count = 0表示标志未经修改 . 它们将此处理为数据依赖性,并且因为uop只能有2个输入(无论如何都是HSW / BDW),因此需要多次uop) . 这是人们抱怨x86的疯狂CISC设计所指的那种 . 它使得x86 CPU比现在从头开始设计ISA的速度慢,即使是以大致相似的方式 . (即这是“x86税”的一部分)这会花费速度/功率 . )SHRX / SHLX / SARX(BMI2)是一个巨大的胜利(1 uop / 1c延迟) .

    它还在关键路径上放置了tzcnt(Haswell及更高版本的3c),因此它显着延长了循环携带的依赖链的总延迟 . 但它确实消除了对CMOV的任何需要,或者准备了一个保存 n>>1 的寄存器 . @Veedrac's answer overcomes all this by deferring the tzcnt/shift for multiple iterations, which is highly effective (see below).

    我们可以安全地交替使用BSFTZCNT,因为 n 在此时永远不会为零 . TZCNT 's machine-code decodes as BSF on CPUs that don' t支持BMI1 . (忽略无意义的前缀,因此REP BSF作为BSF运行) .

    TZCNT在支持它的AMD CPU上比BSF表现要好得多,所以使用 REP BSF 是个好主意,即使你不关心如果输入为零而不是输出设置ZF . 当您使用 __builtin_ctzll 时,某些编译器会执行此操作,即使使用 -mno-bmi 也是如此 .

    它们在Intel CPU上执行相同的操作,因此如果从额外的REP字节中无法获得's all that matters. TZCNT on Intel (pre-Skylake) still has a false-dependency on the supposedly write-only output operand, just like BSF, to support the undocumented behaviour that BSF with input = 0 leaves its destination unmodified. So you need to work around that unless optimizing only for Skylake, so there',则只需保存该字节 . (英特尔经常超出x86 ISA手册所要求的范围,以避免破坏广泛使用的代码,这些代码依赖于它不应该的东西,或者追溯不允许的东西 . 例如Windows 9x's assumes no speculative prefetching of TLB entries,这在编写代码时是安全的,before Intel updated the TLB management rules . )

    无论如何,Haswell的LZCNT / TZCNT与POPCNT具有相同的假设:见this Q&A . 这就是为什么在gcc 's asm output for @Veedrac'代码中,你看到breaking the dep chain with xor-zeroing在寄存器's about to use as TZCNT'的目的地,当它没有值得一些晶体管/电源使它们像其他uops一样运行到同一个执行单元 . 唯一的软件可见上行是与另一个微体系结构限制的交互:在Haswell上,但在Skylake上,英特尔删除了对LZCNT / TZCNT的错误依赖性,他们"un-laminate"索引寻址模式,而POPCNT仍然可以微融合任何地址模式 .


    来自其他答案的想法/代码的改进:

    @hidefromkgb's answer 有一个很好的观察,你可以保证在3n之后能够做一次右移 . 你可以更有效地计算它,而不仅仅是在步骤之间省略检查 . 然而,该答案中的asm实现被破坏了(它取决于OF,在SHRD之后未定义,计数> 1),并且缓慢: ROR rdi,2SHRD rdi,rdi,2 快,并且在关键路径上使用两个CMOV指令比额外的TEST可以并行运行 .

    我整理/改进了C(指导编译器生成更好的asm),并测试了更快的工作速度asm(在C中的评论中)在Godbolt上:请参阅@hidefromkgb's answer中的链接 . (这个答案达到了大型Godbolt网址的30k字符限制,但是shortlinks can rot并且对于goo.gl而言太长了 . )

    还改进了输出打印以转换为字符串并使一个 write() 而不是一次写一个字符 . 这样可以最大限度地减少对整个程序的计时影响 perf stat ./collatz (记录性能计数器),并且我对一些非关键的asm进行了去混淆 .


    @Veedrac's code

    我根据我们知道的需要做了一次非常小的加速,并检查继续循环 . 从CoresDuo(Merom)的7.5s for limit = 1e8下降到7.275s,展开因子为16 .

    代码注释on Godbolt . 不要与clang一起使用这个版本;它使用延迟循环做一些愚蠢的事情 . 使用tmp计数器 k 然后将其添加到 count 以后更改clang的作用,但这会轻微伤害gcc .

    请参阅注释中的讨论:Veedrac的代码在具有BMI1的CPU上非常出色(即不是赛扬/奔腾)

  • 4

    在从源代码生成机器代码期间,C程序被转换为汇编程序 . 说装配比C慢,这几乎是错误的 . 而且,生成的二进制代码因编译器而异 . 因此,智能C编译器可以生成比哑组装器代码更优化和更高效的二进制代码 .

    但是我相信你的分析方法有一定的缺陷 . 以下是分析的一般准则:

    • 确保您的系统处于正常/空闲状态 . 停止所有正在运行的进程(应用程序)或者密集使用CPU(或通过网络轮询) .

    • 您的数据大小必须更大 .

    • 您的测试必须运行超过5-10秒 .

    • 不要只依赖一个样本 . 进行N次测试 . 收集结果并计算结果的平均值或中位数 .

  • 92

    为了获得更好的性能:一个简单的变化是观察到在n = 3n 1之后,n将是偶数,所以你可以立即除以2 . n不会是1,所以你不需要测试它 . 所以你可以保存一些if语句并写:

    while (n % 2 == 0) n /= 2;
    if (n > 1) for (;;) {
        n = (3*n + 1) / 2;
        if (n % 2 == 0) {
            do n /= 2; while (n % 2 == 0);
            if (n == 1) break;
        }
    }
    

    这是一个很大的胜利:如果你看一下n的最低8位,所有步骤直到你除以2 8次完全取决于这8位 . 例如,如果最后8位是0x01,那就是二进制,你的数字是???? 0000 0001然后接下来的步骤是:

    3n+1 -> ???? 0000 0100
    / 2  -> ???? ?000 0010
    / 2  -> ???? ??00 0001
    3n+1 -> ???? ??00 0100
    / 2  -> ???? ???0 0010
    / 2  -> ???? ???? 0001
    3n+1 -> ???? ???? 0100
    / 2  -> ???? ???? ?010
    / 2  -> ???? ???? ??01
    3n+1 -> ???? ???? ??00
    / 2  -> ???? ???? ???0
    / 2  -> ???? ???? ????
    

    所以这些步骤都可以预测,256k 1被替换为81k 1.所有类似的东西都会发生组合 . 所以你可以使用一个大的switch语句创建一个循环:

    k = n / 256;
    m = n % 256;
    
    switch (m) {
        case 0: n = 1 * k + 0; break;
        case 1: n = 81 * k + 1; break; 
        case 2: n = 81 * k + 1; break; 
        ...
        case 155: n = 729 * k + 425; break;
        ...
    }
    

    运行循环直到n≤128,因为在那个点上n可以变为1而少于8个除以2,并且一次做8个或更多个步骤将使你错过第一次达到1的点 . 然后继续“正常”循环 - 或准备一个表格,告诉您需要多少步骤才能达到1 .

    PS . 我强烈怀疑Peter Cordes的建议会让它变得更快 . 除了一个之外,根本没有条件分支,除非循环实际结束,否则将正确预测一个分支 . 所以代码就像是

    static const unsigned int multipliers [256] = { ... }
    static const unsigned int adders [256] = { ... }
    
    while (n > 128) {
        size_t lastBits = n % 256;
        n = (n >> 8) * multipliers [lastBits] + adders [lastBits];
    }
    

    在实践中,您将测量一次处理n的最后9,10,11,12位是否会更快 . 对于每个位,表中的条目数将加倍,并且当表不再适合L1高速缓存时,我的速度减慢 .

    PPS . 如果你需要操作次数:在每次迭代中我们完成八个除以2和一个可变数量的(3n 1)操作,因此计算操作的一个明显方法是另一个数组 . 但我们实际上可以计算步数(基于循环的迭代次数) .

    我们可以稍微重新定义问题:如果奇数将n替换为(3n 1)/ 2,如果是偶数则用n / 2替换n . 然后每次迭代都会完成8个步骤,但你可以考虑作弊:-)所以假设有r个操作n < - 3n 1和s个操作n < - n / 2 . 结果将非常精确地为n'= n * 3 ^ r / 2 ^ s,因为n < - 3n 1意味着n < - 3n *(1 1 / 3n) . 取对数,我们发现r =(s log2(n'/ n))/ log2(3) .

    如果我们执行循环直到n≤1,000,000并且具有预先计算的表,从任何起始点n≤1,000,000需要多少次迭代然后如上所述计算r,四舍五入到最接近的整数,将给出正确的结果,除非s真的很大 .

  • 17

    在一个相当无关的说明:更多性能黑客!

    [第一个«猜想»最终被@ShreevatsaR揭穿;删除]

    • 遍历序列时,我们只能在当前元素 N 的2邻域中获得3个可能的情况(如前所示):

    • [偶] [奇数]

    • [单数] [偶数]

    • [偶数] [偶数]

    跳过这两个元素意味着分别计算 (N >> 1) + N + 1((N << 1) + N + 1) >> 1N >> 2 .

    让我们证明对于两种情况(1)和(2)都可以使用第一个公式 (N >> 1) + N + 1 .

    案例(1)显而易见 . 情况(2)暗示 (N & 1) == 1 ,所以如果我们假设(不失一般性)N是2位长并且其位从最重要到最不重要是 ba ,那么 a = 1 ,并且以下成立:

    (N << 1) + N + 1:     (N >> 1) + N + 1:
    
            b10                    b1
             b1                     b
           +  1                   + 1
           ----                   ---
           bBb0                   bBb
    

    其中 B = !b . 右移第一个结果给了我们我们想要的东西 .

    Q.E.D . : (N & 1) == 1 ⇒ (N >> 1) + N + 1 == ((N << 1) + N + 1) >> 1 .

    事实证明,我们可以使用单个三元运算一次遍历序列2个元素 . 另外2倍的时间缩短 .

    生成的算法如下所示:

    uint64_t sequence(uint64_t size, uint64_t *path) {
        uint64_t n, i, c, maxi = 0, maxc = 0;
    
        for (n = i = (size - 1) | 1; i > 2; n = i -= 2) {
            c = 2;
            while ((n = ((n & 3)? (n >> 1) + n + 1 : (n >> 2))) > 2)
                c += 2;
            if (n == 2)
                c++;
            if (c > maxc) {
                maxi = i;
                maxc = c;
            }
        }
        *path = maxc;
        return maxi;
    }
    
    int main() {
        uint64_t maxi, maxc;
    
        maxi = sequence(1000000, &maxc);
        printf("%llu, %llu\n", maxi, maxc);
        return 0;
    }
    

    这里我们比较 n > 2 因为如果序列的总长度是奇数,则进程可以在2而不是1处停止 .

    [编辑:]

    让我们把它翻译成大会!

    MOV RCX, 1000000;
    
    
    
    DEC RCX;
    AND RCX, -2;
    XOR RAX, RAX;
    MOV RBX, RAX;
    
    @main:
      XOR RSI, RSI;
      LEA RDI, [RCX + 1];
    
      @loop:
        ADD RSI, 2;
        LEA RDX, [RDI + RDI*2 + 2];
        SHR RDX, 1;
        SHRD RDI, RDI, 2;    ror rdi,2   would do the same thing
        CMOVL RDI, RDX;      Note that SHRD leaves OF = undefined with count>1, and this doesn't work on all CPUs.
        CMOVS RDI, RDX;
        CMP RDI, 2;
      JA @loop;
    
      LEA RDX, [RSI + 1];
      CMOVE RSI, RDX;
    
      CMP RAX, RSI;
      CMOVB RAX, RSI;
      CMOVB RBX, RCX;
    
      SUB RCX, 2;
    JA @main;
    
    
    
    MOV RDI, RCX;
    ADD RCX, 10;
    PUSH RDI;
    PUSH RCX;
    
    @itoa:
      XOR RDX, RDX;
      DIV RCX;
      ADD RDX, '0';
      PUSH RDX;
      TEST RAX, RAX;
    JNE @itoa;
    
      PUSH RCX;
      LEA RAX, [RBX + 1];
      TEST RBX, RBX;
      MOV RBX, RDI;
    JNE @itoa;
    
    POP RCX;
    INC RDI;
    MOV RDX, RDI;
    
    @outp:
      MOV RSI, RSP;
      MOV RAX, RDI;
      SYSCALL;
      POP RAX;
      TEST RAX, RAX;
    JNE @outp;
    
    LEA RAX, [RDI + 59];
    DEC RDI;
    SYSCALL;
    

    使用以下命令编译:

    nasm -f elf64 file.asm
    ld -o file file.o
    

    通过Peter Cordes on Godbolt查看C和asm的改进/ bugfixed版本 . (编者注:很抱歉把我的东西放在你的答案中,但我的答案达到Godbolt链接文字的30k char限制!)

相关问题