为什么x ** 3比x * x * x慢？

Question

在NumPy中，x * x * x比x ** 3或甚至np.power（x，3）快一个数量级。

x = np.random.rand(1e6) 
%timeit x**3 
100 loops, best of 3: 7.07 ms per loop 

%timeit x*x*x 
10000 loops, best of 3: 163 µs per loop 

%timeit np.power(x, 3) 
100 loops, best of 3: 7.15 ms per loop 

有关为何发生此行为的任何想法？据我可以告诉所有三个产量相同的输出（使用np.allclose检查）。

Answer 1

按照this answer，这是因为幂乘的实现有一些开销，乘法不是。然而，随着指数的增加，初始乘法会变得越来越慢。经验证明：

In [3]: x = np.random.rand(1e6) 

In [15]: %timeit x**2 
100 loops, best of 3: 11.9 ms per loop 

In [16]: %timeit x*x 
100 loops, best of 3: 12.7 ms per loop 

In [17]: %timeit x**3 
10 loops, best of 3: 132 ms per loop 

In [18]: %timeit x*x*x 
10 loops, best of 3: 27.2 ms per loop 

In [19]: %timeit x**4 
10 loops, best of 3: 132 ms per loop 

In [20]: %timeit x*x*x*x 
10 loops, best of 3: 42.4 ms per loop 

In [21]: %timeit x**10 
10 loops, best of 3: 132 ms per loop 

In [22]: %timeit x*x*x*x*x*x*x*x*x*x 
10 loops, best of 3: 137 ms per loop 

In [24]: %timeit x**15 
10 loops, best of 3: 132 ms per loop 

In [25]: %timeit x*x*x*x*x*x*x*x*x*x*x*x*x*x*x 
1 loops, best of 3: 212 ms per loop 

注幂停留的时间除外x**2情况下，我怀疑或多或少的恒定，是特例，而乘法变得越来越慢。看来你可以利用这个来获得更快的整数幂...例如：

In [26]: %timeit x**16 
10 loops, best of 3: 132 ms per loop 

In [27]: %timeit x*x*x*x*x*x*x*x*x*x*x*x*x*x*x*x 
1 loops, best of 3: 225 ms per loop 

In [28]: def tosixteenth(x): 
    ....:  x2 = x*x 
    ....:  x4 = x2*x2 
    ....:  x8 = x4*x4 
    ....:  x16 = x8*x8 
    ....:  return x16 
    ....: 

In [29]: %timeit tosixteenth(x) 
10 loops, best of 3: 49.5 ms per loop 

看来，你可以通过拆分任意整数到的两个大国的总和，计算每两个功率一般应用该技术如上所述，且求和：

In [93]: %paste 
def smartintexp(x, exp): 
    result = np.ones(len(x)) 
    curexp = np.array(x) 
    while True: 
     if exp%2 == 1: 
      result *= curexp 
     exp >>= 1 
     if not exp: break 
     curexp *= curexp 
    return result 
## -- End pasted text -- 

In [94]: x 
Out[94]: 
array([ 0.0163407 , 0.57694587, 0.47336487, ..., 0.70255032, 
     0.62043303, 0.0796748 ]) 

In [99]: x**21 
Out[99]: 
array([ 3.01080670e-38, 9.63466181e-06, 1.51048544e-07, ..., 
     6.02873388e-04, 4.43193256e-05, 8.46721060e-24]) 

In [100]: smartintexp(x, 21) 
Out[100]: 
array([ 3.01080670e-38, 9.63466181e-06, 1.51048544e-07, ..., 
     6.02873388e-04, 4.43193256e-05, 8.46721060e-24]) 

In [101]: %timeit x**21 
10 loops, best of 3: 132 ms per loop 

In [102]: %timeit smartintexp(x, 21) 
10 loops, best of 3: 70.7 ms per loop 

它速度快两小甚至权力：

In [106]: %timeit x**32 
10 loops, best of 3: 131 ms per loop 

In [107]: %timeit smartintexp(x, 32) 
10 loops, best of 3: 57.4 ms per loop 

但随着指数变大变慢：

In [97]: %timeit x**63 
10 loops, best of 3: 133 ms per loop 

In [98]: %timeit smartintexp(x, 63) 
10 loops, best of 3: 110 ms per loop 

而对于大型最坏情况下并不快：

In [115]: %timeit x**511 
10 loops, best of 3: 135 ms per loop 

In [114]: %timeit smartintexp(x, 511) 
10 loops, best of 3: 192 ms per loop 

Answer 2

这是因为python中的权力是作为浮动操作执行的（对于numpy也是如此，因为它使用C）。

在C中，pow function提供3种方法：

双POW（双X，双Y）

长POWL（长双X，长双Y）

浮子函数powf（ float x，float y）

这些都是浮点运算。

Answer 3

我认为这是因为x**y必须处理通用情况，其中x和y都是浮点数。在数学上，我们可以写x**y = exp(y*log(x))。以你的例子我发现

x = np.random.rand(1e6) 
%timeit x**3 
10 loops, best of 3: 178 ms per loop 

%timeit np.exp(3*np.log(x)) 
10 loops, best of 3: 176 ms per loop 

我还没有检查实际的numpy代码，但它必须在内部做这样的事情。

Answer 4

timeit np.multiply(np.multiply(x,x),x) 

倍一样x*x*x。我的猜测是np.multiply正在使用像BLAS这样的快速Fortran线性代数包。我从另一个问题知道numpy.dot在某些情况下使用BLAS。

---

我必须回来。 np.dot(x,x)比np.sum(x*x)快3倍。所以np.multiply的速度优势与使用BLAS不一致。

---

与我numpy的（次将与机器和可用的库会发生变化）

np.power(x,3.1) 
np.exp(3.1*np.log(x)) 

花费大约在同一时间，但

np.power(x,3) 

是2X一样快。没有像x*x*x那么快，但仍然比一般电力更快。所以它正在利用整数功率。

Answer 5

作为一个说明，如果你正在计算的权力和担心速度：

x = np.random.rand(5e7) 

%timeit x*x*x 
1 loops, best of 3: 522 ms per loop 

%timeit np.einsum('i,i,i->i',x,x,x) 
1 loops, best of 3: 288 ms per loop 

为什么einsum更快仍然是mine一个悬而未决的问题。虽然其类似于einsum能够使用SSE2，而numpy的ufuncs不会到1.8。

在地方甚至更快：

def calc_power(arr): 
    for x in xrange(arr.shape[0]): 
     arr[x]=arr[x]*arr[x]*arr[x] 
numba_power = autojit(calc_power) 

%timeit numba_power(x) 
10 loops, best of 3: 51.5 ms per loop 

%timeit np.einsum('i,i,i->i',x,x,x,out=x) 
10 loops, best of 3: 111 ms per loop 

%timeit np.power(x,3,out=x) 
1 loops, best of 3: 609 ms per loop 

Answer 6

按照spec：

两个参数形式POW（X，Y）是相当于使用功率 操作者：X * * Y。

参数必须有数字类型。对于混合操作数类型，适用于二元算术运算符的 强制规则。

换句话说：由于x是浮子，该指数是从int转换为浮点数，并且执行通用浮点功率运行。在内部，这通常被改写为：

x**y = 2**(y*lg(x)) 

2**a和lg a（基地a 2对数）是现代处理器一条指令，但它仍然需要更长的时间超过两乘法。