Numpy vs Cython速度

Question

我有一个分析代码,使用numpy进行一些繁重的数值运算.只是为了好奇,尝试用cython编译它几乎没有变化,然后我用循环为numpy部分重写它.

令我惊讶的是,基于循环的代码要快得多(8x).我不能发布完整的代码,但我把一个非常简单的无关计算放在一起,显示出类似的行为(虽然时间差异不是很大):

版本1(没有cython)

import numpy as np

def _process(array):

    rows = array.shape[0]
    cols = array.shape[1]

    out = np.zeros((rows, cols))

    for row in range(0, rows):
        out[row, :] = np.sum(array - array[row, :], axis=0)

    return out

def main():
    data = np.load('data.npy')
    out = _process(data)
    np.save('vianumpy.npy', out)

版本2(使用cython构建模块)

import cython
cimport cython

import numpy as np
cimport numpy as np

DTYPE = np.float64
ctypedef np.float64_t DTYPE_t

@cython.boundscheck(False)
@cython.wraparound(False)
@cython.nonecheck(False)
cdef _process(np.ndarray[DTYPE_t, ndim=2] array):

    cdef unsigned int rows = array.shape[0]
    cdef unsigned int cols = array.shape[1]
    cdef unsigned int row
    cdef np.ndarray[DTYPE_t, ndim=2] out = np.zeros((rows, cols))

    for row in range(0, rows):
        out[row, :] = np.sum(array - array[row, :], axis=0)

    return out

def main():
    cdef np.ndarray[DTYPE_t, ndim=2] data
    cdef np.ndarray[DTYPE_t, ndim=2] out
    data = np.load('data.npy')
    out = _process(data)
    np.save('viacynpy.npy', out)

版本3(使用cython构建模块)

import cython
cimport cython

import numpy as np
cimport numpy as np

DTYPE = np.float64
ctypedef np.float64_t DTYPE_t

@cython.boundscheck(False)
@cython.wraparound(False)
@cython.nonecheck(False)
cdef _process(np.ndarray[DTYPE_t, ndim=2] array):

    cdef unsigned int rows = array.shape[0]
    cdef unsigned int cols = array.shape[1]
    cdef unsigned int row
    cdef np.ndarray[DTYPE_t, ndim=2] out = np.zeros((rows, cols))

    for row in range(0, rows):
        for col in range(0, cols):
            for row2 in range(0, rows):
                out[row, col] += array[row2, col] - array[row, col]

    return out

def main():
    cdef np.ndarray[DTYPE_t, ndim=2] data
    cdef np.ndarray[DTYPE_t, ndim=2] out
    data = np.load('data.npy')
    out = _process(data)
    np.save('vialoop.npy', out)

使用保存在data.npy中的10000x10矩阵,时间为:

$ python -m timeit -c "from version1 import main;main()"
10 loops, best of 3: 4.56 sec per loop

$ python -m timeit -c "from version2 import main;main()"
10 loops, best of 3: 4.57 sec per loop

$ python -m timeit -c "from version3 import main;main()"
10 loops, best of 3: 2.96 sec per loop

这是预期的还是我缺少的优化？版本1和版本2给出相同结果的事实是以某种方式预期的,但为什么版本3更快？

Ps.-这不是我需要做的计算,只是一个显示同样事物的简单例子.

Answer 1

稍作修改,版本3的速度提高了两倍:

@cython.boundscheck(False)
@cython.wraparound(False)
@cython.nonecheck(False)
def process2(np.ndarray[DTYPE_t, ndim=2] array):

    cdef unsigned int rows = array.shape[0]
    cdef unsigned int cols = array.shape[1]
    cdef unsigned int row, col, row2
    cdef np.ndarray[DTYPE_t, ndim=2] out = np.empty((rows, cols))

    for row in range(rows):
        for row2 in range(rows):
            for col in range(cols):
                out[row, col] += array[row2, col] - array[row, col]

    return out

计算中的瓶颈是内存访问.您的输入数组是C有序的,这意味着沿最后一个轴移动会使内存中的跳转最小.因此,您的内环应该沿轴1而不是轴0.进行此更改会将运行时间减少一半.

如果需要在小输入数组上使用此函数,则可以使用np.empty而不是使用来减少开销np.ones.为了减少PyArray_EMPTY从numpy C API 进一步使用的开销.

如果在非常大的输入数组(2**31)上使用此函数,则用于索引(和range函数)的整数将溢出.为安全起见:

cdef Py_ssize_t rows = array.shape[0]
cdef Py_ssize_t cols = array.shape[1]
cdef Py_ssize_t row, col, row2

代替

cdef unsigned int rows = array.shape[0]
cdef unsigned int cols = array.shape[1]
cdef unsigned int row, col, row2

定时:

In [2]: a = np.random.rand(10000, 10)
In [3]: timeit process(a)
1 loops, best of 3: 3.53 s per loop
In [4]: timeit process2(a)
1 loops, best of 3: 1.84 s per loop

process你的版本3 在哪里

Answer 2

正如其他答案中所提到的,版本2与版本1基本相同,因为cython无法深入到阵列访问运算符以优化它.这有两个原因

• 首先,与优化的C代码相比,每次调用numpy函数都会产生一定的开销.但是,如果每个操作处理大型数组,则此开销将变得不那么重要

• 其次,有中间数组的创建.如果你考虑更复杂的操作,例如,这会更清楚out[row, :] = A[row, :] + B[row, :]*C[row, :].在这种情况下,B*C必须在内存中创建整个数组,然后添加到A.这意味着CPU缓存正在被打乱,因为数据正在从内存中读取和写入,而不是保存在CPU中并立即使用.重要的是,如果您正在处理大型数组,这个问题会变得更糟.

特别是因为您声明您的实际代码比您的示例更复杂,并且它显示出更高的加速,我怀疑第二个原因可能是您的主要因素.

顺便说一句,如果你的计算足够简单,你可以通过使用numexpr来克服这种影响,虽然当然cython在很多情况下都很有用,所以它可能是更好的方法.

Answer 3

我建议使用-a标志让cython生成html文件,该文件显示正在转换为纯c与调用python API的内容:

http://docs.cython.org/src/quickstart/cythonize.html

版本2提供了与版本1几乎相同的结果,因为所有繁重的工作都是由Python API(通过numpy)完成的,并且cython没有为您做任何事情.实际上在我的机器上,numpy是针对MKL构建的,所以当我使用gcc编译cython生成的c代码时,版本3实际上比其他两个慢一点.

当您进行numpy无法以"矢量化"方式执行的数组操作时,或者当您执行内存密集型操作时,Cython会让您无法创建大型临时数组.对于我自己的一些代码,我使用cython vs numpy获得了115倍的加速:

https://github.com/synapticarbors/pylangevin-integrator

部分原因是在c代码级别调用randomkit目录而不是调用它numpy.random,但大部分是cython将计算密集型for循环转换为纯c而不调用python.