numpy ufunc/算术性能 - 整数不使用 SSE?

Posted

技术标签:

【中文标题】numpy ufunc/算术性能 - 整数不使用 SSE?【英文标题】:numpy ufunc/arithmetic performance - integer not using SSE? 【发布时间】:2015-02-18 19:39:42 【问题描述】:

考虑以下 iPython 性能测试,我们在其中创建一对 10,000 长的 32 位向量并将它们相加。先用整数运算,再用浮点运算:

from numpy.random import randint
from numpy import int32, float32

a, b = randint(255,size=10000).astype(int32), randint(255,size=10000).astype(int32)
%timeit a+b  # int32 addition, gives 20.6µs per loop

a, b = randint(255,size=10000).astype(float32), randint(255,size=10000).astype(float32)
%timeit a+b  # float32 addition, gives 3.91µs per loop

为什么浮点版本快 5 倍?

如果您使用 float64 进行相同的测试,则所需时间是 float32 的两倍,如果我们充分利用硬件,这是您所期望的。然而,对于int8int64,整数情况的时间似乎是恒定的。这一点,再加上 5 倍的减速,让我怀疑它完全无法使用 SSE。

对于int32,当a+ba & 0xffa >> 2 替换时,我观察到类似的20µs 值,这表明问题不仅限于加法。

我正在使用numpy 1.9.1,但不幸的是我不记得我是在本地编译它还是下载了二进制文件。但无论哪种方式,这种性能观察对我来说都是非常令人震惊的。我的版本怎么可能在整数运算方面如此无望?

编辑: 我还在一台类似但独立的 PC 上进行了测试,运行 numpy 1.8,我很确定它直接来自 PythonXY 二进制文件。我得到了同样的结果。

问题:其他人是否看到类似的结果,如果没有,我该怎么做才能像他们一样?

更新:我创建了a new issue on numpy's github repo。

【问题讨论】:

您可以使用np.show_config()import numpy.distutils.system_info as sysinfo; sysinfo.show_all() 收集有关如何编译numpy 的信息。见this SO question。 @unutbu - 我刚刚运行了sysinfo.show_all() - 请参阅output on pastebin - 它主要是“不可用”的列表。 您的观察是正确的:在撰写本文时,只有布尔运算和浮点运算使用 SIMD。这是通过显式内在函数调用完成的,因此它应该很少依赖于编译设置。与任何其他开源项目一样,我们始终欢迎贡献... @Jaime - 你在说这个文件吗:numpy/numpy/core/src/umath/simd.inc.src。我不声称对它在做什么有任何想法,但令我惊讶的是,一旦程序员处理了所有广播/步进逻辑,编译器就没有机会自动矢量化循环。是否存在我可以表现出兴趣的现有问题/里程碑。如果没有,您是否会考虑创建一个...我看不到自己在短期内破解 numpy 的核心。 【参考方案1】:

如果编译器支持,尚未发布的 numpy 1.10 也将向量化整数运算。 它是在此更改中添加的: https://github.com/numpy/numpy/pull/5144

例如您使用 gcc 4.8 编译的当前 git head 的测试用例导致 int 和 float 的速度相同,并且生成的代码看起来不错:

  0.04 │27b:   movdqu (%rdx,%rax,1),%xmm0
 25.33 │       add    $0x1,%r10
       │       movdqu (%r8,%rax,1),%xmm1
       │       paddd  %xmm1,%xmm0
 23.17 │       movups %xmm0,(%rcx,%rax,1)
 34.72 │       add    $0x10,%rax
 16.05 │       cmp    %r10,%rsi
       │     ↑ ja     27b

如果 cpu 支持(例如 intel haswell),可以使用 AVX2 存档额外的加速,尽管目前需要通过使用 OPT="-O3 -mavx2" 编译来完成,但在 numpy 中还没有运行时检测。

【讨论】:

好消息!各种观看/广播要求都是这样吗?还是仅限于我上面给出的简单示例? 应该涵盖所有整数 ufunc,如果遗漏了什么,请提交错误。请注意,只有连续的操作(步幅等于元素大小)才能被矢量化。此外,默认的 64 位整数通常不会从矢量化中受益。 原则上你当然可以在其他一些情况下进行矢量化,例如当 MMX 寄存器大小是 elementsize*stride 的整数倍或因子时。我想这是一个相当常见的用例,有一个 nx2 数组或一个 nd 数组,每个轴的长度为 2 的幂。 (我猜你可以更进一步,找到这两种尺寸的最小公倍数等......但我想在某些时候你只会导致代码膨胀并获得递减收益。)【参考方案2】:

在现代 CPU 上,有很多因素会影响性能。数据是整数还是浮点数只是其中之一。

诸如数据是在缓存中还是必须从 RAM 中获取(或者更糟的是从交换中获取)等因素会产生很大的影响。

用来编译numpy的编译器也会有很大的影响;使用SIMD 指令如SSE 有多好?这些可以显着加快数组操作。

我的系统(Intel Core2 Quad Q9300)的结果;

In [1]: from numpy.random import randint

In [2]: from numpy import int32, float32, float64

In [3]: a, b = randint(255,size=10000).astype(int32), randint(255,size=10000).astype(int32)

In [4]: %timeit a+b
100000 loops, best of 3: 12.9 µs per loop

In [5]: a, b = randint(255,size=10000).astype(float32), randint(255,size=10000).astype(float32)

In [6]: %timeit a+b
100000 loops, best of 3: 8.25 µs per loop

In [7]: a, b = randint(255,size=10000).astype(float64), randint(255,size=10000).astype(float64)

In [8]: %timeit a+b
100000 loops, best of 3: 13.9 µs per loop

所以在这台机器上,int32float32 之间没有五分之一的差异。在float32float64 之间也没有两倍。

从处理器利用率来看,timeit 循环仅使用四个可用内核之一。 这似乎证实了这些简单的操作不使用 BLAS 例程,因为这个 numpy 是使用并行的 openBLAS 构建的。

numpy 的编译方式也会产生重大影响。 根据对this question 的回答,我可以看到使用objdump 我的numpy 使用SSE2 指令和xmm 寄存器。

In [9]: from numpy import show_config

In [10]: show_config()
atlas_threads_info:
    library_dirs = ['/usr/local/lib']
    language = f77
    include_dirs = ['/usr/local/include']
    define_macros = [('ATLAS_INFO', '"\\"None\\""')]
    libraries = ['alapack', 'ptf77blas', 'ptcblas', 'atlas']
openblas_lapack_info:
  NOT AVAILABLE
blas_opt_info:
    library_dirs = ['/usr/local/lib']
    language = f77
    libraries = ['openblasp', 'openblasp']
mkl_info:
  NOT AVAILABLE
lapack_mkl_info:
  NOT AVAILABLE
lapack_opt_info:
    library_dirs = ['/usr/local/lib']
    language = f77
    include_dirs = ['/usr/local/include']
    define_macros = [('ATLAS_INFO', '"\\"None\\""')]
    libraries = ['alapack', 'ptf77blas', 'ptcblas', 'atlas']
openblas_info:
    library_dirs = ['/usr/local/lib']
    language = f77
    libraries = ['openblasp', 'openblasp']
blas_mkl_info:
  NOT AVAILABLE

如果您想查看您使用的 BLAS 的效果,请使用使用不同 BLAS 库编译的 numpy 运行以下程序。

from __future__ import print_function
import numpy
import sys
import timeit

try:
    import numpy.core._dotblas
    print('FAST BLAS')
except ImportError:
    print('slow blas')

print("version:", numpy.__version__)
print("maxint:", sys.maxsize)
print()

setup = "import numpy; x = numpy.random.random((1000,1000))"
count = 5

t = timeit.Timer("numpy.dot(x, x.T)", setup=setup)
print("dot:", t.timeit(count)/count, "sec")

在我的机器上我得到了;

FAST BLAS
version: 1.9.1
maxint: 9223372036854775807

dot: 0.06626860399264842 sec

根据该测试的结果,我从 ATLAS 切换到 OpenBLAS,因为它在我的机器上明显更快。

【讨论】:

鉴于两组数据具有相同的缓存占用空间:10,000 个 32 位元素连续存储,因此关于缓存的声明可能与此处无关。此外,%timeit% 应该在运行主循环之前预热缓存,因此它们都从平等的基础开始。 此外,多核利用率的缺乏并不令人惊讶,因为可能有一个很小的甜蜜点,数据适合缓存,但不会小到导致多核过度杀伤......如果数据不适合缓存,那么多核没有任何好处。

以上是关于numpy ufunc/算术性能 - 整数不使用 SSE?的主要内容,如果未能解决你的问题,请参考以下文章

Python Numpy TypeError:输入类型不支持ufunc'isfinite'

ValueError:numpy.ufunc 大小已更改,可能表示二进制不兼容。预期来自 C 标头的 216,从 PyObject 获得 192

我可以使用`xarray.apply_ufunc`并行化`numpy.bincount`吗?

为啥numpy在尝试预测回归时会引发异常错误:“ufunc'add'没有包含具有签名匹配类型的循环”?

numpy vectorize np.prod 无法构造超过 32 个操作数的 ufunc

numpy