Fan*_*ang 48 c++ performance cmath pow
我一直在做一些项目Euler练习,以提高我对C++的知识.
我写了以下函数:
int a = 0,b = 0,c = 0;
for (a = 1; a <= SUMTOTAL; a++)
{
for (b = a+1; b <= SUMTOTAL-a; b++)
{
c = SUMTOTAL-(a+b);
if (c == sqrt(pow(a,2)+pow(b,2)) && b < c)
{
std::cout << "a: " << a << " b: " << b << " c: "<< c << std::endl;
std::cout << a * b * c << std::endl;
}
}
}
Run Code Online (Sandbox Code Playgroud)
这计算在17毫秒.
但是,如果我改变了这条线
if (c == sqrt(pow(a,2)+pow(b,2)) && b < c)
Run Code Online (Sandbox Code Playgroud)
至
if (c == sqrt((a*a)+(b*b)) && b < c)
Run Code Online (Sandbox Code Playgroud)
计算在2毫秒内完成.是否有一些明显的实现细节pow(int, int)
我缺少哪个使得第一个表达式计算速度慢得多?
Rin*_*g Ø 68
pow()
使用真实的浮点数并在引擎盖下使用公式
pow(x,y) = e^(y log(x))
Run Code Online (Sandbox Code Playgroud)
计算x^y
.将int
被转换为double
之前调用pow
.(log
是自然对数,基于e)
x^2
pow()
因此使用比慢x*x
.
根据相关评论进行编辑
pow
整数指数甚至可能会产生不正确的结果(PaulMcKenzie)pow
还有函数调用(而x*x
不是)(jtbandes)Pet*_*des 37
你已经选择了一种最慢的方法来检查
c*c == a*a + b*b // assuming c is non-negative
Run Code Online (Sandbox Code Playgroud)
这会编译为三个整数乘法(其中一个可以从循环中提升).即使没有pow()
,你仍然会转换double
并采用平方根,这对吞吐量来说很糟糕.(还有延迟,但分支预测+现代CPU上的推测性执行意味着延迟不是这里的一个因素).
Intel Haswell的SQRTSD指令的吞吐量为每8-14个周期一个(源:Agner Fog的指令表),所以即使你的sqrt()
版本保持FP sqrt执行单元饱和,它仍然比我发出的gcc慢4倍(下面).
当b < c
条件的一部分变为false时,您还可以优化循环条件以摆脱循环,因此编译器只需执行该检查的一个版本.
void foo_optimized()
{
for (int a = 1; a <= SUMTOTAL; a++) {
for (int b = a+1; b < SUMTOTAL-a-b; b++) {
// int c = SUMTOTAL-(a+b); // gcc won't always transform signed-integer math, so this prevents hoisting (SUMTOTAL-a) :(
int c = (SUMTOTAL-a) - b;
// if (b >= c) break; // just changed the loop condition instead
// the compiler can hoist a*a out of the loop for us
if (/* b < c && */ c*c == a*a + b*b) {
// Just print a newline. std::endl also flushes, which bloats the asm
std::cout << "a: " << a << " b: " << b << " c: "<< c << '\n';
std::cout << a * b * c << '\n';
}
}
}
}
Run Code Online (Sandbox Code Playgroud)
这-O3 -mtune=haswell
用内循环编译(使用gcc6.2 )代码.查看Godbolt编译器资源管理器的完整代码.
# a*a is hoisted out of the loop. It's in r15d
.L6:
add ebp, 1 # b++
sub ebx, 1 # c--
add r12d, r14d # ivtmp.36, ivtmp.43 # not sure what this is or why it's in the loop, would have to look again at the asm outside
cmp ebp, ebx # b, _39
jg .L13 ## This is the loop-exit branch, not-taken until the end
## .L13 is the rest of the outer loop.
## It sets up for the next entry to this inner loop.
.L8:
mov eax, ebp # multiply a copy of the counters
mov edx, ebx
imul eax, ebp # b*b
imul edx, ebx # c*c
add eax, r15d # a*a + b*b
cmp edx, eax # tmp137, tmp139
jne .L6
## Fall-through into the cout print code when we find a match
## extremely rare, so should predict near-perfectly
Run Code Online (Sandbox Code Playgroud)
在Intel Haswell上,所有这些指令均为1 uop.(并且cmp/jcc将宏 - 融合成比较和分支的uops.)因此,这是10个融合域uops,它可以每2.5个周期在一次迭代中发出.
Haswell imul r32, r32
以每个时钟一次迭代的吞吐量运行,因此内部循环内的两个乘法不会使端口1饱和,每2.5c两次乘以.这留下了空间来吸收来自ADD和SUB窃取端口1的不可避免的资源冲突.
我们甚至没有接近任何其他执行端口瓶颈,因此前端瓶颈是唯一的问题,这应该在Intel Haswell及更高版本上每2.5个周期运行一次.
循环展开可以帮助减少每次检查的uop数量.例如,用于lea ecx, [rbx+1]
计算下一次迭代的b + 1,因此我们可以imul ebx, ebx
不使用MOV使其具有非破坏性.
强度降低也是可能的:鉴于b*b
我们可以尝试在(b-1) * (b-1)
没有IMUL的情况下进行计算. (b-1) * (b-1) = b*b - 2*b + 1
,也许我们可以做一个lea ecx, [rbx*2 - 1]
然后从中减去b*b
.(没有寻址模式可以减去而不是添加.嗯,也许我们可以保留-b
在寄存器中,并向上计数到零,因此我们可以使用ECX中的lea ecx, [rcx + rbx*2 - 1]
更新b*b
,-b
在EBX中给出).
除非你真的遇到了IMUL吞吐量的瓶颈,否则这可能最终会带来更多的麻烦,而不是一场胜利.看看编译器在C++源代码中实现这种强度降低的效果可能会很有趣.
您也可以使用SSE或AVX对其进行矢量化b
,并行检查4个或8个连续值.由于命中率非常罕见,因此您只需检查8中的任何一个是否有命中,然后在极少数情况下找出匹配中的哪一个.
有关更多优化内容,另请参阅x86标记wiki.