首页 文章

Cython通过阵列广播加速循环

提问于
浏览 次
3

Summary:

你们太棒了......我的真实代码正常运行 . 我接受了JoshAdel的建议，即：

1）将所有ndarray更改为键入的内存视图2）手动展开所有numpy数组计算3）为索引使用静态定义的unsigned int 4）禁用boundscheck和wraparound

而且，非常感谢Veedrac的洞察力！

Original post:

我知道python做这些代码真的很慢：

import numpy as np

def func0():
    x = 0.
    for i in range(1000):
        x += 1.
    return

如果我将其更改为Cython，它可以更快：

import numpy as np
cimport numpy as np

def func0():
    cdef double x = 0.
    for i in range(1000):
        x += 1.
    return

这是结果：

# Python
10000 loops, best of 3: 58.9 µs per loop
# Cython
10000000 loops, best of 3: 66.8 ns per loop

但是，现在我有这种代码，它不是单个数字的循环，而是数组的循环 . （是的......我正在解决PDE，所以这种情况发生了） .

我知道以下示例是愚蠢的，但只是尝试了解计算类型：

纯蟒蛇：

def func1():
    array1 = np.random.rand(50000, 4)
    array2 = np.random.rand(50000)
    for i in range(1000):
        array1[:, 0] += array2
        array1[:, 1] += array2
        array1[:, 2] += array2
        array1[:, 3] += array2
    return

用Cython：

def func1():
    cdef np.ndarray[np.double_t, ndim=2] array1 = np.random.rand(50000, 4)
    cdef np.ndarray[np.double_t, ndim=1] array2 = np.random.rand(50000)
    for i in range(1000):
        array1[:, 0] += array2
        array1[:, 1] += array2
        array1[:, 2] += array2
        array1[:, 3] += array2
    return

而且几乎没有改善 . 与此同时，我知道Python由于开销很大而不善于处理这些巨大的循环 .

# Python
1 loops, best of 3: 299 ms per loop
# Cython
1 loops, best of 3: 300 ms per loop

关于如何改进这些代码的任何建议？

python arrays numpy cython

2 回答

  • 1

    在这两个其他实现中，我玩过

    • 使用cython编译器指令删除numpy通常必须进行的一些检查

    • 使用类型化的内存视图，以便我可以指定内存布局（有时它们通常比旧的缓冲区接口更快）

    • 展开循环，以便我们不要切片机：

    否则，你只是通过cython使用numpy，而numpy已经在相当快的c代码下实现了这个 .

    方法：

    import numpy as np
cimport numpy as np

cimport cython

def func1():
    cdef np.ndarray[np.double_t, ndim=2] array1 = np.random.rand(50000, 4)
    cdef np.ndarray[np.double_t, ndim=1] array2 = np.random.rand(50000)
    cdef unsigned int i
    for i in range(1000):
        array1[:, 0] += array2
        array1[:, 1] += array2
        array1[:, 2] += array2
        array1[:, 3] += array2
    return

@cython.boundscheck(False)
@cython.wraparound(False)
def func2():
    cdef np.ndarray[np.double_t, ndim=2] array1 = np.random.rand(50000, 4)
    cdef np.ndarray[np.double_t, ndim=1] array2 = np.random.rand(50000)
    cdef unsigned int i, k
    for i in range(1000):
        for k in xrange(50000):
            array1[k, 0] += array2[k]
            array1[k, 1] += array2[k]
            array1[k, 2] += array2[k]
            array1[k, 3] += array2[k]
    return


@cython.boundscheck(False)
@cython.wraparound(False)
def func3():
    cdef np.ndarray[np.double_t, ndim=2] array1 = np.random.rand(50000, 4)
    cdef np.ndarray[np.double_t, ndim=1] array2 = np.random.rand(50000)
    cdef np.double_t[::1] a2 = array2
    cdef np.double_t[:,::1] a1 = array1
    cdef unsigned int i, k

    for i in range(1000):
        for k in xrange(50000):
            a1[k, 0] += a2[k]
            a1[k, 1] += a2[k]
            a1[k, 2] += a2[k]
            a1[k, 3] += a2[k]
    return

    时间（在我的机器上 - 编译器和硬件当然可能影响时间）：

    • Pure numpy：1个循环，每个循环最好为3：419 ms

    • 您的原始cython函数，键入 i ：1个循环，最佳3：428 ms每个循环

    • func2：1循环，最佳3：336 ms每循环

    • func3：1循环，每循环最佳3：206 ms

    回复于
  • 4

    很大程度上，你错了 . Python非常适合这些类型的循环 .

    这是因为他们是重量级的 . 正如Josh Adel所展示的那样，一个有效的纯C变种只能实现2到4倍的加速 .

    最好的办法是寻找一些类型的低效率：

    • 很多切片（例如切片很多小行）

    这可以通过手动迭代并移动到 Cython 或 Numba （或另一个JITing编译器）来解决 .

    • Full-on Numpy数组计算速度不够快

    Numba 和 Parakeet 可以提供帮助，但你真的想 Theano 或 numexpr

    • 更复杂，非JITable，不可表达的缓慢的东西

    Cython 或 C 是您唯一的选择，所以转向他们

    See also the first part of this answer.

    回复于

相关问题