使用整数的 Numpy 点积非常慢答案

【问题标题】：Numpy dot product very slow using ints使用整数的 Numpy 点积非常慢
【发布时间】：2012-08-08 01:05:16
【问题描述】：

很抱歉有这么多问题。我在 Intel core 2 Duo 上运行 Mac OSX 10.6。我正在为我的研究运行一些基准测试，但我遇到了另一件让我困惑的事情。

如果我跑步

python -mtimeit -s 'import numpy as np; a = np.random.randn(1e3,1e3)' 'np.dot(a,a)'

我得到以下输出：10 loops, best of 3: 142 msec per loop

但是，如果我运行

python -mtimeit -s 'import numpy as np; a = np.random.randint(10,size=1e6).reshape(1e3,1e3)' 'np.dot(a,a)'

我得到以下输出：10 loops, best of 3: 7.57 sec per loop

然后我跑了

python -mtimeit -s 'import numpy as np; a = np.random.randn(1e3,1e3)' 'a*a'然后

python -mtimeit -s 'import numpy as np; a = np.random.randint(10,size=1e6).reshape(1e3,1e3)' 'a*a'

每个循环都以大约 7.6 毫秒的速度运行，所以它不是乘法。添加也有相似的速度，所以这些都不应该影响点积，对吧？那么为什么使用整数计算点积比使用浮点数慢 50 倍以上呢？

【问题讨论】：

在 Linux 上对我来说也是如此——float64 大约需要 3 秒，int32 大约需要 10 秒（这是一台较旧的机器）。不是 50 倍，但仍然很奇怪。

标签： python floating-point numpy int dot-product

【解决方案1】：

非常有趣，我很好奇它是如何实现的，所以我做了：

>>> import inspect
>>> import numpy as np
>>> inspect.getmodule(np.dot)
<module 'numpy.core._dotblas' from '/Library/Python/2.6/site-packages/numpy-1.6.1-py2.6-macosx-10.6-universal.egg/numpy/core/_dotblas.so'>
>>>

所以看起来它使用了BLAS 库。

所以：

>>> help(np.core._dotblas)

我从中找到了这个：

当 Numpy 使用 ATLAS 等加速 BLAS 构建时，这些函数被替换以使用更快的实现。越快实现仅影响 float32、float64、complex64 和 complex128 数组。此外，BLAS API 仅包含矩阵矩阵，矩阵向量和向量向量乘积。具有较大阵列的产品维度使用内置函数并且没有加速。

所以看起来 ATLAS 微调了某些功能，但它只适用于某些数据类型，非常有趣。

是的，看来我会更频繁地使用花车...

【讨论】：

@Luke 感谢 +1 为你做回溯，这也是另一种有用的方法。
很高兴知道。当我上班时，我将使用相同的方法来查看 MKL 是否也适用。谢谢您的帮助。卢克更喜欢这个答案，所以你会接受。
我喜欢答案不仅回答了问题，而且还展示了如何找到了答案，并在此过程中介绍了一些技巧。 +1。在这里需要更多这样的答案。

【解决方案2】：

使用 int vs float 数据类型会导致执行不同的代码路径：

float 的堆栈跟踪如下所示：

(gdb) backtr
#0  0x007865a0 in dgemm_ () from /usr/lib/libblas.so.3gf
#1  0x007559d5 in cblas_dgemm () from /usr/lib/libblas.so.3gf
#2  0x00744108 in dotblas_matrixproduct (__NPY_UNUSED_TAGGEDdummy=0x0, args=(<numpy.ndarray at remote 0x85d9090>, <numpy.ndarray at remote 0x85d9090>), 
kwargs=0x0) at numpy/core/blasdot/_dotblas.c:798
#3  0x08088ba1 in PyEval_EvalFrameEx ()
...

..而 int 的堆栈跟踪如下所示：

(gdb) backtr
#0  LONG_dot (ip1=0xb700a280 "\t", is1=4, ip2=0xb737dc64 "\a", is2=4000, op=0xb6496fc4 "", n=1000, __NPY_UNUSED_TAGGEDignore=0x85fa960)
at numpy/core/src/multiarray/arraytypes.c.src:3076
#1  0x00659d9d in PyArray_MatrixProduct2 (op1=<numpy.ndarray at remote 0x85dd628>, op2=<numpy.ndarray at remote 0x85dd628>, out=0x0)
at numpy/core/src/multiarray/multiarraymodule.c:847
#2  0x00742b93 in dotblas_matrixproduct (__NPY_UNUSED_TAGGEDdummy=0x0, args=(<numpy.ndarray at remote 0x85dd628>, <numpy.ndarray at remote 0x85dd628>), 
kwargs=0x0) at numpy/core/blasdot/_dotblas.c:254
#3  0x08088ba1 in PyEval_EvalFrameEx ()
...

这两个调用都导致了 dotblas_matrixproduct，但似乎 float 调用保留在 BLAS 库中（可能访问一些优化良好的代码），而 int 调用被踢回 numpy 的 PyArray_MatrixProduct2。

所以这要么是一个错误，要么是 BLAS 不支持 matrixproduct 中的整数类型（这似乎不太可能）。

这是一个简单且廉价的解决方法：

af = a.astype(float)
np.dot(af, af).astype(int)

【讨论】：

值得注意的是，如果您的数据具有非常大的值，此解决方法可能会导致错误，并且可能需要复制完整的矩阵，因此如果矩阵非常大，成本会很高。
谢谢，卢克。该解决方法确实复制了矩阵，结果非常麻烦（出于内存问题），但就时间而言，1e4x1e4 矩阵的速度要快数千倍！任何更大并且使用整数测试乘法太慢。 @Dougal 这仅适用于使用 64 位浮点数大于 2^52 的数字，对吗？数字不会大于这个数字，如果可能的话，我想利用这个加速。
@Nino 是的，就在那儿。 BLAS 库不适用于整数类型真的太糟糕了，而且 numpy 的内置 dot 慢得多。如果内存问题是个问题，你可以考虑写一个小的 ctypes 接口，在 Eigen 或类似的地方做乘法，应该更快。
BLAS 中的基本矩阵-矩阵乘法有 4 个版本：sgemm.f、dgemm.f、cgemm.f、fzgemm.f 分别用于单、双、复、双复。线性代数，包括矩阵求逆和线性方程组的解，几乎没有用于纯整数计算。 netlib.org/blas