>我在Cython中使用了很多3D内存视图,例如
cython.declare(a='double[:, :, ::1]')
a = np.empty((10, 20, 30), dtype='double')
我经常想遍历a
的所有元素。我可以使用三重循环来做到这一点,例如
for i in range(a.shape[0]):
for j in range(a.shape[1]):
for k in range(a.shape[2]):
a[i, j, k] = ...
如果我不关心索引i
、j
和k
,做一个平面循环更有效,比如
cython.declare(a_ptr='double*')
a_ptr = cython.address(a[0, 0, 0])
for i in range(size):
a_ptr[i] = ...
在这里我需要知道数组中元素(size
)的数量。这是由shape
属性中元素的乘积给出的,即size = a.shape[0]*a.shape[1]*a.shape[2]
,或者更一般地说size = np.prod(np.asarray(a).shape)
。我发现这两个写起来都很丑陋,而且(尽管很小)计算开销困扰着我。这样做的好方法是使用 memoryviews 的内置size
属性,size = a.size
。但是,由于我无法理解的原因,这会导致未优化的C代码,从Cython生成的注释html文件中可以明显看出。具体来说,size = a.shape[0]*a.shape[1]*a.shape[2]
生成的 C 代码就是
__pyx_v_size = (((__pyx_v_a.shape[0]) * (__pyx_v_a.shape[1])) * (__pyx_v_a.shape[2]));
其中从size = a.size
生成的 C 代码是
__pyx_t_10 = __pyx_memoryview_fromslice(__pyx_v_a, 3, (PyObject *(*)(char *)) __pyx_memview_get_double, (int (*)(char *, PyObject *)) __pyx_memview_set_double, 0);; if (unlikely(!__pyx_t_10)) __PYX_ERR(0, 2238, __pyx_L1_error)
__Pyx_GOTREF(__pyx_t_10);
__pyx_t_14 = __Pyx_PyObject_GetAttrStr(__pyx_t_10, __pyx_n_s_size); if (unlikely(!__pyx_t_14)) __PYX_ERR(0, 2238, __pyx_L1_error)
__Pyx_GOTREF(__pyx_t_14);
__Pyx_DECREF(__pyx_t_10); __pyx_t_10 = 0;
__pyx_t_7 = __Pyx_PyIndex_AsSsize_t(__pyx_t_14); if (unlikely((__pyx_t_7 == (Py_ssize_t)-1) && PyErr_Occurred())) __PYX_ERR(0, 2238, __pyx_L1_error)
__Pyx_DECREF(__pyx_t_14); __pyx_t_14 = 0;
__pyx_v_size = __pyx_t_7;
为了生成上面的代码,我通过编译器指令启用了所有可能的优化,这意味着a.size
生成的笨拙的 C 代码无法优化。在我看来,size
"属性"并不是真正预先计算的属性,而是实际上在查找时执行计算。此外,这种计算比简单地将乘积超过shape
属性要复杂得多。我在文档中找不到任何解释的提示。
这种行为的解释是什么,如果我真的关心这种微优化,我有比写出a.shape[0]*a.shape[1]*a.shape[2]
更好的选择吗?
通过查看生成的 C 代码,您已经可以看到size
是一个属性,而不是一个简单的 C 成员。以下是用于内存视图的原始 Cython 代码:
@cname('__pyx_memoryview')
cdef class memoryview(object):
...
cdef object _size
...
@property
def size(self):
if self._size is None:
result = 1
for length in self.view.shape[:self.view.ndim]:
result *= length
self._size = result
return self._size
很容易看出,产品只计算一次,然后缓存。显然,它对 3 维数组没有太大作用,但对于更高数量的维度,缓存可能会变得非常重要(正如我们将看到的,最多有 8 个维度,所以它不是那么明确,这个缓存是否真的值得)。
人们可以理解懒惰计算size
的决定——毕竟,size
并不总是需要/使用,也不想为此付费。显然,如果你经常使用这种size
,这种懒惰是要付出代价的——这就是cython所做的权衡。
我不会在调用a.size
的开销上赘述太久 - 与从 python 调用 cython 函数的开销相比,这算不了什么。
例如,@danny的测量仅测量此 python 调用开销,而不是不同方法的实际性能。为了说明这一点,我将第三个函数放入其中:
%%cython
...
def both():
a.size+a.shape[0]*a.shape[1]*a.shape[2]
这做了两倍的工作量,但是
>>> %timeit mv_size
22.5 ns ± 0.0864 ns per loop (mean ± std. dev. of 7 runs, 10000000 loops each)
>>> %timeit mv_product
20.7 ns ± 0.087 ns per loop (mean ± std. dev. of 7 runs, 10000000 loops each)
>>>%timeit both
21 ns ± 0.39 ns per loop (mean ± std. dev. of 7 runs, 10000000 loops each)
同样快。另一方面:
%%cython
...
def nothing():
pass
不是更快:
%timeit nothing
24.3 ns ± 0.854 ns per loop (mean ± std. dev. of 7 runs, 10000000 loops each)
简而言之:由于可读性,我会使用a.size
,假设优化它不会加快我的应用程序速度,除非分析证明不同的东西。
整个故事:变量a
属于__Pyx_memviewslice
类型,而不是人们想象的__pyx_memoryview
类型。结构__Pyx_memviewslice
具有以下定义:
struct __pyx_memoryview_obj;
typedef struct {
struct __pyx_memoryview_obj *memview;
char *data;
Py_ssize_t shape[8];
Py_ssize_t strides[8];
Py_ssize_t suboffsets[8];
} __Pyx_memviewslice;
这意味着,Cython代码可以非常有效地访问shape
,因为它是一个简单的C数组(顺便说一句。我问自己,如果维度超过8个,会发生什么?- 答案是:你不能超过 8 个维度)。
成员memview
是保存内存的地方,__pyx_memoryview_obj
是C扩展,它由我们在上面看到的密码代码生成,如下所示:
/* "View.MemoryView":328
*
* @cname('__pyx_memoryview')
* cdef class memoryview(object): # <<<<<<<<<<<<<<
*
* cdef object obj
*/
struct __pyx_memoryview_obj {
PyObject_HEAD
struct __pyx_vtabstruct_memoryview *__pyx_vtab;
PyObject *obj;
PyObject *_size;
PyObject *_array_interface;
PyThread_type_lock lock;
__pyx_atomic_int acquisition_count[2];
__pyx_atomic_int *acquisition_count_aligned_p;
Py_buffer view;
int flags;
int dtype_is_object;
__Pyx_TypeInfo *typeinfo;
};
因此,Pyx_memviewslice
并不是真正的Python对象 - 它是一种方便的包装器,它缓存重要数据,如shape
和stride
,因此可以快速,廉价地访问这些信息。
当我们打电话给a.size
时会发生什么? 首先,调用__pyx_memoryview_fromslice
,它执行一些额外的引用计数和一些进一步的内容,并从__Pyx_memviewslice
-对象返回成员memview
。
然后在此返回的 memoryview 上调用属性size
,该视图访问_size
中的缓存值,如上面的 Cython 代码所示。
看起来python程序员为shape
,strides
和suboffsets
等重要信息引入了快捷方式,但对于可能不那么重要的size
则没有 - 这就是在shape
的情况下更干净的C代码的原因。
为a.size
生成的 C 代码看起来不错。
它必须与 Python 接口,因为内存视图是 Python 扩展类型。 内存视图上的size
是一个 python 属性,并转换为ssize_t
。这就是 C 代码所做的一切。通过将size
变量键入为Py_ssize_t
而不是ssize_t
来避免转换。
因此,C 代码中没有任何内容看起来未经优化 - 它只是在 python 对象上查找属性,在这种情况下在内存视图上
查找大小。以下是两种方法的微观基准测试结果。
设置:
cimport numpy as np
import numpy as np
cimport cython
cython.declare(a='double[:, :, ::1]')
a = np.empty((10, 20, 30), dtype='double')
def mv_size():
return a.size
def mv_product():
return a.shape[0]*a.shape[1]*a.shape[2]
结果:
%timeit mv_size
10000000 loops, best of 3: 23.4 ns per loop
%timeit mv_product
10000000 loops, best of 3: 23.4 ns per loop
性能几乎相同。
产品方法是纯 C 代码,如果需要并行执行,这很重要,但除此之外,与内存视图size
相比没有性能优势。