Cython比纯Python略快或略慢

4
我正在使用几种技术(NumPy、Weave和Cython)来执行Python性能基准测试。代码在数学上的基本操作是C = AB,其中A、B和C是N x N矩阵(注意:这是矩阵乘积而不是逐元素相乘)。
我编写了5个不同的实现:
1.纯Python(循环遍历2D Python列表) 2.NumPy(2D NumPy数组的点积) 3.Weave inline(C++循环遍历2D数组) 4.Cython(循环遍历2D Python列表+静态类型) 5.Cython-Numpy(循环遍历2D NumPy数组+静态类型)
我的期望是实现2到5比实现1快得多。然而,我的结果表明相反。这些是相对于纯Python实现的归一化加速结果: python_list: 1.00 numpy_array: 330.09 weave_inline: 30.72 cython_list: 2.80 cython_array: 0.14
我对NumPy的性能非常满意,但对Weave的性能不太热衷,而Cython的性能让我哭泣。我的整个代码分为两个文件。一切都是自动化的,您只需运行第一个文件即可查看所有结果。请问有人可以帮助我指出我可以做些什么以获得更好的结果吗?
matmul.py:
import time

import numpy as np
from scipy import weave
from scipy.weave import converters

import pyximport
pyximport.install()
import cython_matmul as cml


def python_list_matmul(A, B):
    C = np.zeros(A.shape, dtype=float).tolist()
    A = A.tolist()
    B = B.tolist()
    for k in xrange(len(A)):
        for i in xrange(len(A)):
            for j in xrange(len(A)):
                C[i][k] += A[i][j] * B[j][k]
    return C


def numpy_array_matmul(A, B):
    return np.dot(A, B)


def weave_inline_matmul(A, B):
    code = """
       int i, j, k;
       for (k = 0; k < N; ++k)
       {
           for (i = 0; i < N; ++i)
           {
               for (j = 0; j < N; ++j)
               {
                   C(i, k) += A(i, j) * B(j, k);
               }
           }
       }
       """

    C = np.zeros(A.shape, dtype=float)
    weave.inline(code, ['A', 'B', 'C', 'N'], type_converters=converters.blitz, compiler='gcc')
    return C


N = 100
A = np.random.rand(N, N)
B = np.random.rand(N, N)

function = []
function.append([python_list_matmul, 'python_list'])
function.append([numpy_array_matmul, 'numpy_array'])
function.append([weave_inline_matmul, 'weave_inline'])
function.append([cml.cython_list_matmul, 'cython_list'])
function.append([cml.cython_array_matmul, 'cython_array'])

t = []
for i in xrange(len(function)):
    t1 = time.time()
    C = function[i][0](A, B)
    t2 = time.time()
    t.append(t2 - t1)
    print function[i][1] + ' \t: ' + '{:10.6f}'.format(t[0] / t[-1])

cython_matmul.pyx:

import numpy as np
cimport numpy as np

import cython
cimport cython

DTYPE = np.float
ctypedef np.float_t DTYPE_t


@cython.boundscheck(False)
@cython.wraparound(False)
@cython.nonecheck(False)
cpdef cython_list_matmul(A, B):

    cdef int i, j, k
    cdef int N = len(A)

    A = A.tolist()
    B = B.tolist()
    C = np.zeros([N, N]).tolist()

    for k in xrange(N):
        for i in xrange(N):
            for j in xrange(N):
                C[i][k] += A[i][j] * B[j][k]
    return C


@cython.boundscheck(False)
@cython.wraparound(False)
@cython.nonecheck(False)
cpdef cython_array_matmul(np.ndarray[DTYPE_t, ndim=2] A, np.ndarray[DTYPE_t, ndim=2] B):

    cdef int i, j, k, N = A.shape[0]
    cdef np.ndarray[DTYPE_t, ndim=2] C = np.zeros([N, N], dtype=DTYPE)

    for k in xrange(N):
        for i in xrange(N):
            for j in xrange(N):
                C[i][k] += A[i][j] * B[j][k]
    return C

据我所知,您只运行了一次代码。这会对像weave这样的东西产生偏差,因为在第一次计算期间,weave需要实际编译内联代码——随后的调用可能会绕过这个非常昂贵的步骤,因为我认为内联代码将被缓存。 - mgilson
是的,编写的代码确实只执行一次,但是当我更改某些内容时,我会确保多次运行代码。 - Aeronaelius
为什么你在Cython-Numpy中使用显式循环?如果我没记错的话,包括np.import_array()可以让你使用Numpy的C API,并调用np.PyArray_MatrixProduct2(A, B, C)而不是执行循环。 (或者这是否使您的cython函数变得无用?)使用np.matrix而不是np.ndarray也可能会影响性能,但我不确定有多大影响。 - JAB
到目前为止,我只使用过基本的Python和基本的Numpy。我对Python的C API和Numpy没有任何熟悉度。这就是我的原因。但是感谢您指出这一点。我可能会开始深入挖掘,并通过尽可能地挤取所有内容来使自己更加熟练。 - Aeronaelius
1个回答

11

Python列表和高性能数学不兼容,忘记cython_list_matmul。

你的cython_array_matmul唯一的问题是索引使用不正确。应该是

C[i,k] += A[i,j] * B[j,k]

这是Python中numpy数组的索引方式,也是Cython优化的语法。通过这种改变,您应该可以获得不错的性能表现。
Cython的注释功能在发现此类优化问题方面非常有帮助。您可能会注意到A[i][j]产生了大量的Python API调用,而A[i,j]则没有任何调用。
此外,如果您手动初始化所有条目,则np.emptynp.zeros更合适。

Nikita - 我本以为自己是一个NumPy专家,但我不知道np.empty。谢谢! - Rick
谢谢Nikita,这个方法很有效。使用cython_array_matmul,我加速了约70倍。有没有办法让它更快?也许可以使用指针算术?@Rick,不管你是多么专业,你总能学到新东西,这太棒了 :)。 - Aeronaelius
强制使用连续缓冲区([mode='c'])可能会稍微提高一点性能,但此时您应该切换到现代的缓存无关算法,采用多线程并重写为SIMD汇编代码......或者只需使用优化的BLAS实现(例如Intel MKL)。 - Nikita Nemkin
1
OpenBLAS 是免费的,并且在速度方面与 MKL 相当。 - ali_m

网页内容由stack overflow 提供, 点击上面的
可以查看英文原文,
原文链接