如您所知,黄线表示发生了与python的某些交互,即使用了python功能而不是原始的c功能,您可以查看生成的代码,看看会发生什么,以及是否可以/应该解决/避免. .
并非每次与python的交互都意味着(可衡量的)速度下降.
让我们看一下这个简化的函数:
%%cython
cimport numpy as np
def use_slices(np.ndarray[np.double_t] a):
a[0:len(a)]=0.0
当我们查看生成的代码时,我们看到了(我只保留了重要部分):
__pyx_t_1 = PyObject_Length(((PyObject *)__pyx_v_a));
__pyx_t_2 = PyInt_FromSsize_t(__pyx_t_1);
__pyx_t_3 = PySlice_New(__pyx_int_0, __pyx_t_2, Py_None);
PyObject_SetItem(((PyObject *)__pyx_v_a)
因此,基本上我们得到了一个新的切片(它是一个numpy数组),然后使用numpy的功能(PyObject_SetItem)将所有元素设置为0.0,这是引擎盖下的C代码.
让我们看一下手写的for循环的版本:
cimport numpy as np
def use_for(np.ndarray[np.double_t] a):
cdef int i
for i in range(len(a)):
a[i]=0.0
它仍然使用PyObject_Length(由于长度)和边界检查,但是否则它是C代码.当我们比较时间时:
>>> import numpy as np
>>> a=np.ones((500,))
>>> %timeit use_slices(a)
100000 loops, best of 3: 1.85 ?s per loop
>>> %timeit use_for(a)
1000000 loops, best of 3: 1.42 ?s per loop
>>> b=np.ones((250000,))
>>> %timeit use_slices(b)
10000 loops, best of 3: 104 ?s per loop
>>> %timeit use_for(b)
1000 loops, best of 3: 364 ?s per loop
您可以看到为小尺寸创建切片的额外开销,但是for-version中的额外检查意味着从长远来看,它会有更多开销.
让我们禁用这些检查:
%%cython
cimport cython
cimport numpy as np
@cython.boundscheck(False)
@cython.wraparound(False)
def use_for_no_checks(np.ndarray[np.double_t] a):
cdef int i
for i in range(len(a)):
a[i]=0.0
在产生的html中,我们可以看到a [i]变得像它一样简单:
__pyx_t_3 = __pyx_v_i;
*__Pyx_BufPtrStrided1d(__pyx_t_5numpy_double_t *, __pyx_pybuffernd_a.rcbuffer->pybuffer.buf, __pyx_t_3, __pyx_pybuffernd_a.diminfo[0].strides) = 0.0;
}
__Pyx_BufPtrStrided1d(type,buf,i0,s0)是为(type)((char *)buf i0 * s0)定义的.
现在:
>>> %timeit use_for_no_checks(a)
1000000 loops, best of 3: 1.17 ?s per loop
>>> %timeit use_for_no_checks(b)
1000 loops, best of 3: 246 ?s per loop
我们可以通过在for循环中释放gil来进一步改善它:
%%cython
cimport cython
cimport numpy as np
@cython.boundscheck(False)
@cython.wraparound(False)
def use_for_no_checks_no_gil(np.ndarray[np.double_t] a):
cdef int i
cdef int n=len(a)
with nogil:
for i in range(n):
a[i]=0.0
现在:
>>> %timeit use_for_no_checks_no_gil(a)
1000000 loops, best of 3: 1.07 ?s per loop
>>> %timeit use_for_no_checks_no_gil(b)
10000 loops, best of 3: 166 ?s per loop
因此它速度更快,但是对于较大的数组,您仍然无法击败numpy.
我认为有两点可借鉴:
> Cython不会将切片转换为通过for循环进行的访问,因此必须使用Python功能.
>开销很小,但是它只调用numpy功能,大部分工作是通过numpy代码完成的,因此无法通过Cython加快速度.
最后尝试使用memset函数:
%%cython
from libc.string cimport memset
cimport numpy as np
def use_memset(np.ndarray[np.double_t] a):
memset(&a[0], 0, len(a)*sizeof(np.double_t))
我们得到:
>>> %timeit use_memset(a)
1000000 loops, best of 3: 821 ns per loop
>>> %timeit use_memset(b)
10000 loops, best of 3: 102 ?s per loop
对于大型数组,它也与numpy代码一样快.
正如DavidW所建议的,可以尝试使用内存视图:
%%cython
cimport numpy as np
def use_slices_memview(double[::1] a):
a[0:len(a)]=0.0
导致小数组的代码速度稍快,但大数组的代码相似(与numpy-slices相比):
>>> %timeit use_slices_memview(a)
1000000 loops, best of 3: 1.52 ?s per loop
>>> %timeit use_slices_memview(b)
10000 loops, best of 3: 105 ?s per loop
这意味着内存视图切片比numpy切片具有更少的开销.这是产生的代码:
__pyx_t_1 = __Pyx_MemoryView_Len(__pyx_v_a);
__pyx_t_2.data = __pyx_v_a.data;
__pyx_t_2.memview = __pyx_v_a.memview;
__PYX_INC_MEMVIEW(&__pyx_t_2, 0);
__pyx_t_3 = -1;
if (unlikely(__pyx_memoryview_slice_memviewslice(
&__pyx_t_2,
__pyx_v_a.shape[0], __pyx_v_a.strides[0], __pyx_v_a.suboffsets[0],
0,
0,
&__pyx_t_3,
0,
__pyx_t_1,
0,
1,
1,
0,
1) < 0))
{
__PYX_ERR(0, 27, __pyx_L1_error)
}
{
double __pyx_temp_scalar = 0.0;
{
Py_ssize_t __pyx_temp_extent = __pyx_t_2.shape[0];
Py_ssize_t __pyx_temp_idx;
double *__pyx_temp_pointer = (double *) __pyx_t_2.data;
for (__pyx_temp_idx = 0; __pyx_temp_idx < __pyx_temp_extent; __pyx_temp_idx++) {
*((double *) __pyx_temp_pointer) = __pyx_temp_scalar;
__pyx_temp_pointer += 1;
}
}
}
__PYX_XDEC_MEMVIEW(&__pyx_t_2, 1);
__pyx_t_2.memview = NULL;
__pyx_t_2.data = NULL;
我认为最重要的部分是:此代码不会创建其他临时对象-它为切片重用了现有的内存视图.
如果使用内存视图,则编译器(至少对于我的机器而言)会生成稍快的代码.不确定是否值得调查.乍一看,每个迭代步骤的区别是:
# created code for memview-slices:
*((double *) __pyx_temp_pointer) = __pyx_temp_scalar;
__pyx_temp_pointer += 1;
#created code for memview-for-loop:
__pyx_v_i = __pyx_t_3;
__pyx_t_4 = __pyx_v_i;
*((double *) ( /* dim=0 */ ((char *) (((double *) data) + __pyx_t_4)) )) = 0.0;
我希望不同的编译器以不同的方式处理此代码.但是很明显,第一个版本更容易优化.
正如Behzad Jamali所指出的,double [:] a和double [:: 1] a之间是有区别的.使用切片的第二个版本在我的计算机上快约20%.区别在于,在编译期间,对于double [:: 1]版本知道,内存访问将是连续的,并且可以将其用于优化.在具有double [:]的版本中,直到运行时我们才对步幅一无所知.