为什么用于测试Collatz猜想的C++代码比手写汇编运行得更快?

我为Project Euler Q14编写了这两个解决方案,在组装和C++中。它们实现了相同的蛮力方法来测试Collatz猜想。组装解决方案组装在一起:

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

C++汇编如下:

g++ p14.cpp -o p14

大会,p14.asm

section .datafmt db "%d", 10, 0
global mainextern printf
section .text
main:mov rcx, 1000000xor rdi, rdi        ; max ixor rsi, rsi        ; i
l1:dec rcxxor r10, r10        ; countmov rax, rcx
l2:test rax, 1jpe even
mov rbx, 3mul rbxinc raxjmp c1
even:mov rbx, 2xor rdx, rdxdiv rbx
c1:inc r10cmp rax, 1jne l2
cmp rdi, r10cmovl rdi, r10cmovl rsi, rcx
cmp rcx, 2jne l1
mov rdi, fmtxor rax, raxcall printfret

C++,p14.cpp

#include <iostream>
int sequence(long n) {int count = 1;while (n != 1) {if (n % 2 == 0)n /= 2;elsen = 3*n + 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;}}std::cout << maxi << std::endl;}

我知道编译器优化可以提高速度和一切,但我看不到有多少方法可以进一步优化我的汇编解决方案(以编程方式,而不是数学方式)。

C++代码每项使用模,每项使用除法,而汇编代码每项只使用一个除法。

但是组装的时间平均比C++解决方案长1秒。这是为什么?我主要是出于好奇问。

执行时间

我的系统:64位Linux1.4 GHz Intel Celeron 2955U(Haswell微架构)。

  • g++(未优化): avg 1272 ms.
  • g++ -O3: avg 578 ms.
  • asm (div)(原始):avg 2650 ms。
  • asm (shr): avg 679 ms.
  • @刘德华(与NASM组装):avg 501 ms。
  • @陈志立: avg 200 ms.
  • @VeedracC++-O3平均81毫秒,-O0平均305毫秒。
172646 次浏览

C++程序在从源代码生成机器代码的过程中被翻译成汇编程序。说汇编比C++慢实际上是错误的。此外,生成的二进制代码因编译器而异。因此,智能C++编译器可能产生的二进制代码比哑汇编程序代码更最优和高效。

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

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

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

请参阅Agner Fog的优化装配指南了解如何编写高效的asm。他还有指令表和微架构指南,了解特定CPU的特定详细信息。另请参阅标签wiki以获取更多perf链接。

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

通常你可以让编译器做它的事情,特别是如果你尝试编写可以有效编译的C++。另见汇编语言比编译语言快吗?。其中一个答案链接到这些整齐的幻灯片,展示了各种C编译器如何用很酷的技巧优化一些非常简单的函数。马特·戈德博尔特的CppCon2017演讲"我的编译器最近为我做了什么?解开编译器的盖子"是类似的。


even:mov rbx, 2xor rdx, rdxdiv rbx

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

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

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


从gcc的#0 asm输出可以看出(戈德博尔特),它只使用换档指令. clang-O0确实像你想象的那样天真地编译,即使使用了两次64位IDIV。(优化时,当源使用相同的操作数进行除法和取模时,编译器确实会使用IDIV的两个输出,如果他们根本使用IDIV)

GCC没有完全简单的模式;它总是通过GIMTH进行转换,这意味着某些“优化”不能被禁用。这包括识别除以常数并使用移位(2的幂)或定点乘法逆(2的非幂)来避免IDIV(参见上面的gobolt链接中的div_by_13)。

gcc -Os(优化大小)确实使用IDIV进行非2次幂除法,不幸的是,即使在乘法逆码只是稍大但快得多的情况下。


帮助编译器

(本案例摘要:使用uint64_t n

首先,查看优化的编译器输出才有趣。(-O3)。
#0速度基本上没有意义。<强>

看看你的ASM输出(在戈德博尔特上,或参见如何从GCC/clang组件输出中消除“噪音”?)。当编译器一开始就没有做出最佳代码时:以引导编译器编写更好代码的方式编写C/C++源代码通常是最好的方法。你必须知道ASM,知道什么是高效的,但你间接应用了这些知识。编译器也是一个很好的想法来源:有时clang会做一些很酷的事情,你可以手持gcc做同样的事情:参见这个答案和我在下面@Veedrac的代码中对非展开循环所做的事情。)

这种方法是可移植的,在20年后,一些未来的编译器可以将其编译为未来硬件上有效的任何东西(x86与否),也许使用新的ISA扩展或自动矢量化。15年前的手写x86-64 asm通常不会针对Skylake进行最佳调整。例如,当时还不存在比较和分支宏融合。现在对于一个微架构的手工制作的ASM来说是最佳的,但对于其他当前和未来的CPU来说可能不是最佳的。评论@johnfind的回答讨论了AMD Bulldozer和Intel Haswell之间的主要区别,这对这段代码有很大影响。但理论上,g++ -O3 -march=bdver3g++ -O3 -march=skylake会做正确的事情。(或-march=native。)或者-mtune=...只是调整,不使用其他CPU可能不支持的指令。

我的感觉是,引导编译器使用对你关心的当前CPU有益的asm对未来的编译器来说应该不是问题。他们希望在找到转换代码的方法方面比当前的编译器更好,并且能够找到一种适用于未来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不能为负,或者什么。Signed-overflow是未定义的行为,所以它应该能够。)

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


顺便说一句,gcc的优化的asm输出看起来不错(使用#0):它内联到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 + 1mov    rdi, raxshr    rdi         # rdi = n>>1;test   al, 1       # set flags based on n%2 (aka n&1)mov    rax, rcxcmove  rax, rdi    # n= (n%2) ? 3*n+1 : n/2;add    edx, 1      # ++count;cmp    rax, 1jne   .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也不是关键路径的一部分:它不是循环承载的,因为控制依赖关系是通过分支预测+推测执行处理的,这与关键路径上的数据依赖关系不同。


击败编译器

GCC在这方面做得很好。使用#0而不是#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 + 1shr     rax, 1         # n>>=1;    CF = n&1 = n%2cmovc   rax, rcx       # n= (n&1) ? 3*n+1 : n/2;inc     edx            # ++count;cmp     rax, 1jne     .L9            #}while(n!=1)

另一个聪明的技巧请参见@johnfind的回答:通过在SHR的标志结果上分支并将其用于CMOV:零来删除CMP,仅当n从一开始就是1(或0)时。(有趣的事实:如果读取标志结果,则Nehalem或更早版本上的计数为!=1的SHR会导致失速。这就是他们如何使其成为单uop的。不过,Shive-by-1特殊编码很好。)

避免MOV对Haswell(x86的MOV真的可以“免费”吗?为什么我完全不能复制这个?)的延迟没有任何帮助。它确实在像Intel pre-IvB和AMD Bulldozer系列这样的CPU上帮助了显着,其中MOV不是零延迟(以及带有更新微代码的Ice Lake)。编译器浪费的MOV指令确实影响了关键路径。BD的复杂LEA和CMOV都有较低的延迟(分别为2c和1c),所以它的延迟比例更大。此外,吞吐量瓶颈成为一个问题,因为它只有两个整数ALU管道。参见@johnfind的回答,他有AMD CPU的计时结果。

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

LEA的延迟取决于寻址模式,在英特尔SnB系列CPU上。3个组件的3c([base+idx+const],需要两个单独的添加),但只有2个或更少组件的1c(一个添加)。一些CPU(如Core2)甚至在一个周期内完成3个组件的LEA,但SnB系列没有。更糟糕的是,英特尔SnB系列标准化延迟,因此没有2c uops,否则3个组件的LEA将像推土机一样只有2c。(3个组件的LEA在AMD上也更慢,只是没有那么多)。

所以在像Haswell这样的英特尔SnB系列CPU上,lea rcx, [rax + rax*2]/inc rcx只有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向量(4xuint64_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, xmm0vpaddq    ymm1, ymm1, xmm0vpaddq    ymm1, ymm1, set1_epi64(1)     # ymm1= 3*n + 1.  Maybe could do this more efficiently?
vpsllq    ymm3, ymm0, 63                # shift bit 1 to the sign bit
vpsrlq    ymm0, ymm0, 1                 # n /= 2
# FP blend between integer insns may cost extra bypass latency, but integer blends don't have 1 bit controlling a whole qword.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, ymm4jnz  .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, 1vpmaxq .... 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实现相同的逻辑之外,还要寻找简化逻辑或避免冗余工作的方法。例如memoze以检测序列的常见结尾。或者更好的是,一次查看8个尾随位(gnasher的答案)

@EOF指出tzcnt(或bsf)可以用于在一步中执行多个n/=2迭代。这可能比SIMD向量化要好;没有SSE或AVX指令可以做到这一点。不过,它仍然兼容在不同的整数寄存器中并行执行多个标量n

所以循环可能看起来像这样:

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

这可能会做更少的迭代,但是在没有BMI2的英特尔SnB系列CPU上,变量计数移位很慢。3个uops,2c延迟。(它们对FLAGS有输入依赖关系,因为coun=0意味着标志是未修改的。它们将此作为数据依赖关系处理,并需要多个uops,因为一个uop只能有2个输入(无论如何,在HSW/BDW之前)。)。这就是人们抱怨x86疯狂的CISC设计所指的那种。这使得x86 CPU比今天从头开始设计ISA的情况下要慢,即使是以基本相似的方式。(即,这是成本速度/功耗的“x86税”的一部分。)SHRX/SHLX/SARX(BMI2)是一个巨大的胜利(1 uop/1c延迟)。

它还将tzcnt(Haswell及更高版本上的3c)放在关键路径上,因此它显着延长了循环携带的依赖链的总延迟。不过,它确实消除了对CMOV或准备寄存器保存n>>1的任何需求。@Veedrac的答案通过推迟多次迭代的tzcnt/Shift来克服这一切,这是非常有效的(见下文)。

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

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

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

无论如何,Haswell上的LZCNT/TZCNT与POPCNT具有相同的false dep:请参阅本次问答。这就是为什么在@Veedrac代码的gcc asm输出中,您会在寄存器上看到它用异零法打破dep链,当它不使用dst=src时,它将用作TZCNT的目的地。由于TZCNT/LZCNT/POPCNT永远不会让它们的目的地未定义或未修改,因此这种对英特尔CPU输出的错误依赖是性能bug/限制。大概值得一些晶体管/电源让它们像其他前往同一执行单元的uops一样运行。唯一的优势是与另一个UARC限制的交互:他们可以用索引寻址模式微融合内存操作数在Haswell上,但在Skylake上,英特尔删除了LZCNT/TZCNT的错误dep,他们“取消层压”索引寻址模式,而POPCNT仍然可以微融合任何addr模式。


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

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

我把整理/改进的C(它引导编译器产生更好的asm),并测试+更快地工作asm(在C下面的注释中)在Godbolt上:请参阅@hidefromkgb的回答中的链接。(这个答案达到了大型Godbolt URL的30kchar限制,但是短链接会腐烂,无论如何对goo.gl来说太长了。)

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


@Veedrac的代码

在Core2Duo(Merom)上,我们知道需要做的右移和检查以继续循环,我得到了一个小小的加速。从7.5秒的限制=1e8下降到7.275秒,展开因子为16。

代码+评论

声称C++编译器可以产生比合格的汇编语言程序员更多的最佳代码是一个非常糟糕的错误。尤其是在这种情况下。人类总是可以比编译器更好地编写代码,这种特殊情况很好地说明了这种说法。

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

(以下代码为32位,但可以轻松转换为64位)

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

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

整个代码看起来像:

include "%lib%/freshlib.inc"@BinaryType console, compactoptions.DebugMode = 1include "%lib%/freshlib.asm"
start:InitializeAllmov ecx, 999999xor edi, edi        ; maxxor ebx, ebx        ; max i
.main_loop:
xor     esi, esimov     eax, ecx
.seq:inc     esi                 ; counterlea     edx, [3*eax+1]      ; edx = 3*n+1shr     eax, 1              ; eax = n/2cmovc   eax, edx            ; if CF eax = edxjnz     .seq                ; jmp if n<>1
cmp     edi, esicmovb   edi, esicmovb   ebx, ecx
dec     ecxjnz     .main_loop
OutputValue "Max sequence: ", edi, 10, -1OutputValue "Max index: ", ebx, 10, -1
FinalizeAllstdcall TerminateAll, 0

为了编译此代码,需要FreshLib

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

两个程序的输出是相同的:最大序列=525 on i=837799。

来自评论:

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

对于许多数字,它将没有溢出。

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

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

任何简单的最终收敛级数都以两个值的幂开始(足够明显吗?)。

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

我们能产生2^64吗?如果开始的数字是0x5555555555555555,它是奇数,下一个数字是3n+1,即0xFFFFFFFFFFFFFFFF + 1=0。理论上处于算法的未定义状态,但约翰芬的优化答案将通过退出ZF=1来恢复。彼得·科德斯将以无限循环结束cmp rax,1(QED变体1,通过未定义的0数字的“Cheapo”)。

一些更复杂的数字如何,这将创建没有0的循环?坦率地说,我不确定我的数学理论是如此模糊,无法得到任何严肃的想法,如何严肃地处理它。但是直觉上我会说,对于0<的每个数,级数将收敛到1,因为3n+1公式迟早会将原始数(或中间数)的每个非2素数因子慢慢变成2的幂。所以我们不需要担心原始级数的无限循环,只有溢出会阻碍我们。

所以我只是把几个数字放入工作表中,并查看了8位截断的数字。

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

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

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

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

不幸的是,我无法深入了解这个系列,为什么它会收敛,以及将它们截断为k位的含义是什么,但是在cmp number,1终止条件下,当然可以将算法放入无限循环中,特定的输入值在截断后以0结尾。

但是对于8位的情况,值27溢出是一种警告,这看起来像是如果你计算达到值1的步数,你将从总的k位整数集中得到错误的结果。对于8位整数,256个数字中的146个数字被截断影响了序列(其中一些可能仍然偶然击中了正确的步数,也许,我懒得检查)。

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

  • [第一个“猜想”终于被@ShreevatsaR揭穿;删除]

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

    1. [偶数][奇数]
    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                    b1b1                     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 thingCMOVL 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.asmld -o file file.o

请参阅C和Peter Cordes的改进/错误修复版本

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

test rax, 1jpe even

…有50%的可能性会误判分枝,而且代价高昂。

编译器几乎肯定会同时进行这两种计算(由于div/mod的延迟相当长,因此multiply-add是“免费的”,因此成本可以忽略不计)并跟进CMOV。当然,这有0%的可能性被错误预测。

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

编辑:你的里程可能会因为上面的“整数除法比浮点除法快”语句而有所不同。下面的评论显示,现代处理器已经优先优化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;}

为了获得更多性能:一个简单的变化是观察到在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位是0x01,那就是二进制,你的数字是????0000 0001那么接下来的步骤是:

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

所以所有这些步骤都可以预测,256k+1被替换为81k+1。所有组合都会发生类似的情况。所以你可以用一个大开关语句做一个循环:

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可以用少于8个除以2变成1,并且一次执行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。如果你需要操作数:在每次迭代中,我们只做8除以2,以及可变数量的(3n+1)操作,所以计算操作的明显方法是另一个数组。但我们实际上可以计算步数(基于循环的迭代次数)。

我们可以稍微重新定义这个问题:如果奇数,将n替换为(3n+1)/2,如果偶数,将n替换为n/2。然后每次迭代都将准确地执行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真的很大。

作为一个通用的答案,不专门针对此任务:在许多情况下,通过在高层进行改进,可以显着加快任何程序。例如计算数据一次而不是多次,完全避免不必要的工作,以最佳方式使用缓存等。这些事情在高级语言中要容易得多。

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

在欧拉问题中,大多数时候你通过构建一些东西、找到它慢的原因、构建更好的东西、找到它慢的原因等等来获得成功。使用汇编程序非常非常困难。以一半可能速度的更好算法通常会以全速击败更差的算法,在汇编程序中获得全速并不容易。

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

python实现示例:

import sys
inner_loop = 0
def collatz_sequence(N, cache):global inner_loop
l = [ ]stop = Falsen = N
tails = [ ]
while not stop:inner_loop += 1tmp = nl.append(n)if n <= 1:stop = Trueelif n in cache:stop = Trueelif n % 2:n = 3*n + 1else:n = n // 2tails.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 elemif 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))

简单的答案:

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

  • 这里的ADD 1可能比INC快

  • MOV 2和DIV非常昂贵;只需向右移动

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

如果您为C++程序生成程序集列表,您可以看到它与您的程序集有何不同。