optimizing software in c++




用于测试Collat z猜想的C++代码比手写程序集更快-为什么? (8)

我为 Project Euler Q14 编写了这两个解决方案,在汇编和C ++中。 它们是用于测试 Collat​​z猜想的 相同蛮力方法。 装配解决方案与组装

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微体系结构)。


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

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

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


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

test rax, 1
jpe even

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

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


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

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

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

    1. [甚至] [奇数]
    2. [奇偶]
    3. [甚至] [偶数]

    飞跃过去,这些2种元素的方法来计算 (N >> 1) + N + 1((N << 1) + N + 1) >> 1 并且 N >> 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 。 右移第一个结果给了我们我们想要的东西。

    QED : (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

见C和改进/ bugfixed版本的ASM的彼得·柯德斯 在Godbolt 。 (编者注:很抱歉把我的东西放在你的答案中,但是我的答案达到了Godbolt链接+文字的30k char限制!)


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

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

  1. 确保您的系统处于正常/空闲状态。 停止所有正在运行的进程(应用程序)或者密集使用CPU(或通过网络轮询)。
  2. 您的数据大小必须更大。
  3. 您的测试必须运行超过5-10秒。
  4. 不要只依赖一个样本。 进行N次测试。 收集结果并计算结果的均值或中位数。

简单的答案:

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

  • ADD 1可能比INC快

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

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

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


为了获得更好的性能:一个简单的变化是观察到在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真的很大。


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

请参阅 Agner Fog的优化装配指南 ,了解如何编写高效的asm。 他还有针对特定CPU的具体细节的指令表和微指南指南。 有关更多perf链接,另请参阅 x86 标记wiki。

另请参阅关于使用手写asm击败编译器的这个更一般的问题: 内联汇编语言是否比本机C ++代码慢? 。 TL:DR:是的,如果你做错了(就像这个问题)。

通常你很好让编译器做它的事情,特别是如果你 试图编写可以有效编译的C ++ 。 还看到 汇编比编译语言快吗? 。 其中一个答案链接到 这些简洁的幻灯片, 显示各种C编译器如何使用很酷的技巧优化一些非常简单的功能。

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

在Intel Haswell上, div r64 为36 div r64 延迟为 32-96个 周期 ,吞吐量为每21-74个周期一个。 (加上2 uop来设置RBX和零RDX,但是乱序执行可以提前运行)。 像DIV这样的高uop计数指令是微编码的,这也可能导致前端瓶颈。 在这种情况下,延迟是最相关的因素,因为它是循环携带的依赖链的一部分。

shr rax, 1 执行相同的无符号除法:它是1 uop,具有1c延迟 ,并且每个时钟周期可以运行2。

相比之下,32位除法速度更快,但与换档相比仍然很糟糕。 idiv r32 是9 uops,22-29c延迟,Haswell每8-11c吞吐量一个。

正如你从gcc的 -O0 asm输出( Godbolt编译器浏览器 )中看到的 那样 ,它只使用移位指令 。 clang -O0 确实像你想象的那样天真地编译,甚至使用64位IDIV两次。 (当进行优化时,编译器会在源进行除法时使用IDIV的两个输出,并使用相同的操作数使用模数,如果它们完全使用IDIV的话)

海湾合作委员会没有完全天真的模式; 它总是通过GIMPLE进行转换,这意味着无法禁用某些“优化” 。 这包括识别 div_by_13 除法和使用移位(2的幂)或 固定点乘法逆 (2的非幂)以避免IDIV(参见上面的godbolt链接中的div_by_13)。

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

帮助编译器

(本案例摘要:使用 uint64_t n

首先,看一下优化的编译器输出是很有意思的。 ( -O3 )。 -O0 速度基本上没有意义。

看看你的asm输出(在Godbolt上,或者看看 如何从GCC / clang组件输出中删除“噪音”? )。 当编译器首先没有制作最佳代码时: 以引导编译器制作更好代码的方式编写C / C ++源代码通常是最好的方法 。 你必须知道asm,并知道什么是有效的,但你间接地应用这些知识。 编译器也是一个很好的想法来源:有时clang会做一些很酷的事情,你可以手动执行gcc做同样的事情:看看 这个答案 以及我在@ Veedrac的代码中对非展开循环所做的事情。)

这种方法是可移植的,并且在未来20年内,一些未来的编译器可以将其编译为在未来硬件上有效的任何内容(x86或不使用),可能使用新的ISA扩展或自动向量化。 从15年前手写的x86-64 asm通常不会为Skylake进行最佳调整。 例如,比较和分支宏观融合当时不存在。 对于一个微体系结构,手工制作的asm现在最理想的可能不是其他当前和未来CPU的最佳选择。 关于@ johnfound的回答的评论 讨论了AMD Bulldozer和Intel Haswell之间的主要差异,这对这段代码产生了很大的影响。 但理论上, g++ -O3 -march=bdver3g++ -O3 -march=skylake 将做正确的事情。 (或者 -march=native 。)或者 -mtune=... 只是调整,而不使用其他CPU可能不支持的指令。

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

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

在这种情况下 :你的 n 有一个带符号的类型,gcc使用SAR / SHR / ADD序列给出正确的舍入。 (对于负输入,IDIV和算术移位“舍入”不同,请参阅 SAR insn set ref手册输入 )。 (IDK如果gcc尝试并且未能证明 n 不能为负,或者是什么。签名溢出是未定义的行为,所以应该能够。)

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

顺便说一下, gcc的 优化 asm输出看起来相当不错(使用 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)。

总计:每次迭代5个周期,延迟瓶颈 。 乱序执行与此并行处理其他所有事情(理论上:我没有使用perf计数器进行测试,看它是否真的以5c / iter运行)。

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

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

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

打败编译器

海湾合作委员会在这里做得很好。 它可以通过使用 inc edx 而不是 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的另一个聪明诀窍的答案:通过分支SHR的标志结果以及将其用于CMOV来删除CMP:仅当n为1(或0)时才为零。 (有趣的事实: 如果你读取标志结果,那么Nehalem或更早的SHR会计数!= 1会导致失速 。这就是他们如何使它成为单uop。但是,1转换特殊编码很好。)

避免MOV对Haswell的延迟没有帮助( Can x86的MOV真的可以“免费”吗?为什么我根本不能重现这个? )。 它对英特尔pre-IvB和AMD Bulldozer系列等CPU有 很大 帮助,其中MOV不是零延迟。 编译器浪费的MOV指令会影响关键路径。 BD的复杂LEA和CMOV都具有较低的延迟(分别为2c和1c),因此它是延迟的一小部分。 此外,吞吐量瓶颈成为一个问题,因为它只有两个整数ALU管道。 请参阅@ johnfound的答案 ,他有AMD CPU的计时结果。

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

LEA的延迟取决于 Intel SnB系列CPU 上的寻址模式 。 3c为3个组件( [base+idx+const] ,需要两个单独的添加),但只有1c有2个或更少的组件(一个添加)。 有些CPU(如Core2)在一个周期内甚至可以执行3分量LEA,但SnB系列则不行。 更糟糕的是, 英特尔SnB系列标准化延迟,因此没有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上收支平衡,在Core2上更糟糕。 它确实需要额外的uop,这通常不值得节省1c延迟,但延迟是这里的主要瓶颈而Haswell有足够宽的管道来处理额外的uop吞吐量。

gcc,icc和clang(在godbolt上)都没有使用SHR的CF输出,总是使用AND或TEST 。 愚蠢的编译器。 :P它们是复杂机械的重要组成部分,但聪明的人通常可以在小规模问题上击败它们。 (当然,考虑它的时间要长几千到几百万倍!编译器不会使用详尽的算法来搜索每种可能的做事方式,因为在优化大量内联代码时需要花费太长时间,这就是他们做得最好。他们也没有在目标微体系结构中对管道进行建模,至少与 IACA 或其他静态分析工具没有相同的细节;他们只是使用了一些启发式方法。)

简单的循环展开无济于事 ; 这个循环瓶颈是循环传输依赖链的延迟,而不是循环开销/吞吐量。 这意味着它可以很好地处理超线程(或任何其他类型的SMT),因为CPU有很多时间来交错来自两个线程的指令。 这意味着在 main 并行化循环,但这很好,因为每个线程只能检查 n 值的范围并因此产生一对整数。

在单个线程内手工交错也是可行的 。 也许并行计算一对数字的序列,因为每个数字只需要几个寄存器,并且它们都可以更新相同的 max / maxi 。 这创建了更多的 指令级并行性

诀窍是决定是否要等到所有 n 值都达到 1 才能获得另一对起始 n 值,或者是否突破并为只有达到结束条件的一个获得新的起点,而不触及寄存器其他顺序。 可能最好让每个链都处理有用的数据,否则你必须有条件地增加它的计数器。

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

我认为检测 1 “粘性”的最佳策略是屏蔽你添加的所有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 次迭代。 这可能比SIMD矢量化更好,因为没有SSE或AVX指令可以做到这一点。 但它仍然兼容在不同的整数寄存器中并行执行多个标量 n s。

所以循环可能如下所示:

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的答案通过延迟tzcnt / shift进行多次迭代来克服所有这些,这非常有效(见下文)。

我们可以安全地交替使用 BSF TZCNT ,因为 n 在此时永远不会为零。 TZCNT的机器代码在不支持BMI1的CPU上解码为BSF。 (忽略无意义的前缀,因此REP BSF作为BSF运行)。

TZCNT在支持它的AMD CPU上比BSF表现要好得多,所以使用 REP BSF 也是个好主意,即使输入为零而不是输出也不关心设置ZF。 当你使用 __builtin_ctzll 甚至使用 -mno-bmi 时,一些编译器会这样做。

它们在Intel CPU上执行相同的操作,因此只需保存字节,如果这一切都很重要。 英特尔(pre-Skylake)上的TZCNT仍然假定依赖于所谓的只写输出操作数,就像BSF一样,以支持输入= 0的BSF未修改的BSF的未记录行为。 因此,除非仅针对Skylake进行优化,否则您需要解决这个问题,因此额外的REP字节无法获得任何好处。 (英特尔经常超出x86 ISA手册所要求的范围,以避免破坏广泛使用的代码,这些代码依赖于它不应该的东西,或者是追溯不允许的。例如, Windows 9x假定没有推测预取TLB条目 ,这是安全的代码编写时, 在英特尔更新TLB管理规则之前 。)

无论如何,Haswell的LZCNT / TZCNT与POPCNT具有相同的假设:请参阅 此问答 。 这就是为什么在@ Veedrac的代码的gcc的asm输出中,当你不使用dst = src时,你会看到它在将要用作TZCNT的目标的寄存器上 用xor- zero 打破dep链 。 由于TZCNT / LZCNT / POPCNT永远不会使其目标未定义或未修改,因此对Intel CPU上的输出的这种错误依赖纯粹是性能错误/限制。 据推测,值得一些晶体管/功率使它们像其他uops一样运行到同一个执行单元。 唯一的软件可见上行是与另一个微体系结构限制的交互: 它们可以 在Haswell上 使用索引寻址模式 内存操作数进行微融合 ,但在Skylake上,英特尔取消了对LZCNT / TZCNT的错误依赖,他们“取消层压”索引寻址模式,而POPCNT仍然可以微融合任何地址模式。

从其他答案改进想法/代码:

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

我整理/改进了C(它指导编译器产生更好的asm),并在Godbolt上测试+工作更快的asm(在C下方的注释中):请参阅 @ hidefromkgb的答案中 的链接。 (这个答案达到了大型Godbolt网址的30k char限制,但是 短链接可能会腐烂 ,而且对于goo.gl来说太长了。)

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

@ Veedrac的代码

我从正确的移位中得到了一个非常小的加速,我们 知道 需要做的事情,并检查继续循环。 从CoresDuo(Merom)的7.5s for limit = 1e8下降到7.275s,展开因子为16。

关于Godbolt的 代码+评论。 请勿将此版本与clang一起使用; 它使用延迟循环做一些愚蠢的事情。 使用tmp计数器 k ,然后将其添加到 count 以后更改clang的作用,但这会 轻微 伤害gcc。

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


来自评论:

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

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

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

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

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

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

我们可以生产2 ^ 64吗?如果起始编号为 0x5555555555555555 , 则为 奇数,则下一个编号为3n + 1,即 0xFFFFFFFFFFFFFFFF + 1 = 0 。理论上处于未定义的算法状态,但johnfound的优化答案将通过退出ZF = 1来恢复。的 cmp rax,1 彼得科尔德的 将在无限循环结束 (QED变体1,“小气鬼”通过未定义 0 数量)。

一些更复杂的数字怎么样,哪个会创造没有循环 0 ? 坦率地说,我不确定,我的数学理论太朦胧,无法得到任何认真的想法,如何以严肃的方式处理它。 但直觉上我会说每个数字系列会收敛到1:0 <数字,因为3n + 1公式会慢慢将原始数字(或中间)的每个非2素数因子转换为2的幂,迟早。 所以我们不需要担心原始系列的无限循环,只有溢出才会妨碍我们。

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

有三个值溢出到 02271708585 直接进入 0 ,另外两个朝向前进 85 )。

但是没有创造循环溢出种子的价值。

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

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

不幸的是,我不能让这个系列的任何深刻的把握,为什么它收敛,什么是截断他们的影响 ķ 位,但与 cmp number,1 中止条件的确可以把算法与特定输入值结束的无限循环 0 后截断。

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







x86