我正在对一个科学应用程序进行数值优化。我注意到的一件事是,GCC将通过将调用pow(a,2)编译为a*a来优化它,但调用pov(a,6)并没有优化,实际上会调用库函数pow,这会大大降低性能。(相比之下,可执行icc的“英特尔C++编译器”将消除对pow(a,6)的库调用。)

我好奇的是,当我使用GCC 4.5.1和选项“-O3-lm-funroll-loops-msse4”将pow(a,6)替换为a*a*a*a*a*a时,它使用了5条多指令:

movapd  %xmm14, %xmm13
mulsd   %xmm14, %xmm13
mulsd   %xmm14, %xmm13
mulsd   %xmm14, %xmm13
mulsd   %xmm14, %xmm13
mulsd   %xmm14, %xmm13

而如果我写(a*a*a)*(a*a*a),它将产生

movapd  %xmm14, %xmm13
mulsd   %xmm14, %xmm13
mulsd   %xmm14, %xmm13
mulsd   %xmm13, %xmm13

这将乘法指令的数量减少到3。icc也有类似的行为。

为什么编译器不认识这种优化技巧?


当前回答

Lambdageek正确地指出,由于关联性不适用于浮点数,所以a*a*a*a*a*a到(a*a**a)*(a*a*a)的“优化”可能会改变值。这就是C99不允许它的原因(除非用户通过编译器标志或pragma特别允许)。一般来说,假设程序员写的东西是有原因的,编译器应该尊重这一点。如果你想要(a*a*a)*(a*a*a),请写下。

不过,写起来可能会很痛苦;当你使用pow(a,6)时,编译器为什么不能做正确的事情?因为这样做是错误的。在一个拥有良好数学库的平台上,pow(a,6)比a*a*a*a*a*a或(a*a**a)*(a*a*a)要准确得多。为了提供一些数据,我在我的Mac Pro上做了一个小实验,测量了在评估所有[1,2)之间的单精度浮点数的^6时的最差误差:

worst relative error using    powf(a, 6.f): 5.96e-08
worst relative error using (a*a*a)*(a*a*a): 2.94e-07
worst relative error using     a*a*a*a*a*a: 2.58e-07

使用pow而不是乘法树可以将误差范围减少4倍。编译器不应该(通常也不会)进行增加错误的“优化”,除非获得用户的许可(例如通过-fast math)。

注意,GCC提供__builtin_powi(x,n)作为pow()的替代方案,pow应该生成内联乘法树。如果您想以精度换取性能,但不想启用快速数学,请使用该选项。

其他回答

另一个类似的情况是:大多数编译器不会将a+b+c+d优化为(a+b)+(c+d)(这是一种优化,因为第二个表达式可以更好地进行流水线处理),并按照给定的方式对其求值(即(((a+c)+d))。这也是因为角落案例:

float a = 1e35, b = 1e-5, c = -1e35, d = 1e-5;
printf("%e %e\n", a + b + c + d, (a + b) + (c + d));

这将输出1.00000e-05 0.000000e+00

gcc实际上可以进行这种优化,即使对于浮点数也是如此。例如

double foo(double a) {
  return a*a*a*a*a*a;
}

变成

foo(double):
    mulsd   %xmm0, %xmm0
    movapd  %xmm0, %xmm1
    mulsd   %xmm0, %xmm1
    mulsd   %xmm1, %xmm0
    ret

使用-O-funcafe数学优化。但是,这种重新排序违反了IEEE-754,因此需要标记。

正如Peter Cordes在一篇评论中指出的,有符号整数可以在没有funsafe数学优化的情况下进行这种优化,因为它恰好在没有溢出的情况下有效,如果有溢出,则会出现未定义的行为。所以你得到

foo(long):
    movq    %rdi, %rax
    imulq   %rdi, %rax
    imulq   %rdi, %rax
    imulq   %rax, %rax
    ret

只需-O。对于无符号整数,这更容易,因为它们是2的模幂,因此即使在溢出的情况下也可以自由地重新排序。

Lambdageek正确地指出,由于关联性不适用于浮点数,所以a*a*a*a*a*a到(a*a**a)*(a*a*a)的“优化”可能会改变值。这就是C99不允许它的原因(除非用户通过编译器标志或pragma特别允许)。一般来说,假设程序员写的东西是有原因的,编译器应该尊重这一点。如果你想要(a*a*a)*(a*a*a),请写下。

不过,写起来可能会很痛苦;当你使用pow(a,6)时,编译器为什么不能做正确的事情?因为这样做是错误的。在一个拥有良好数学库的平台上,pow(a,6)比a*a*a*a*a*a或(a*a**a)*(a*a*a)要准确得多。为了提供一些数据,我在我的Mac Pro上做了一个小实验,测量了在评估所有[1,2)之间的单精度浮点数的^6时的最差误差:

worst relative error using    powf(a, 6.f): 5.96e-08
worst relative error using (a*a*a)*(a*a*a): 2.94e-07
worst relative error using     a*a*a*a*a*a: 2.58e-07

使用pow而不是乘法树可以将误差范围减少4倍。编译器不应该(通常也不会)进行增加错误的“优化”,除非获得用户的许可(例如通过-fast math)。

注意,GCC提供__builtin_powi(x,n)作为pow()的替代方案,pow应该生成内联乘法树。如果您想以精度换取性能,但不想启用快速数学,请使用该选项。

正如Lambdageek指出的那样,浮点乘法是不相关的,你可以得到更少的精度,但当获得更好的精度时,你可以反对优化,因为你想要一个确定性的应用程序。例如,在游戏模拟客户端/服务器中,每个客户端都必须模拟相同的世界,您希望浮点计算具有确定性。

像“pow”这样的库函数通常经过精心设计,以产生最小可能的错误(在一般情况下)。这通常是用样条逼近函数实现的(根据Pascal的评论,最常见的实现似乎是使用Remez算法)

基本上是以下操作:

pow(x,y);

具有与任何单个乘法或除法中的误差大致相同大小的固有误差。

执行以下操作时:

float a=someValue;
float b=a*a*a*a*a*a;

其固有误差大于单个乘法或除法的误差的5倍(因为您组合了5个乘法)。

编译器应该非常小心它正在进行的优化:

如果将pow(a,6)优化为a*a*a*a*a*a,可能会提高性能,但会大大降低浮点数的精度。如果将a*a*a*a*a*a优化为pow(a,6),实际上可能会降低精度,因为“a”是一个特殊的值,它允许无误差的乘法(2的幂或一些小整数)如果将pow(a,6)优化为(a*a*a)*(a*a*a)或(a*a)*。

一般来说,您知道对于任意浮点值,“pow”的精度比您最终可以编写的任何函数都要高,但在某些特殊情况下,多次乘法可能具有更好的精度和性能,这取决于开发人员选择更合适的方法,最终对代码进行注释,以便其他人不会“优化”该代码。

唯一有意义的优化(个人观点,显然是GCC中没有任何特定优化或编译器标志的选择)应该是将“pow(a,2)”替换为“a*a”。这将是编译器供应商应该做的唯一明智的事情。