从Cython代码生成SIMD指令

Question

我需要大致了解一下在高性能数字代码中使用Cython可以获得的性能。我感兴趣的一件事是找出优化的C编译器是否可以向量化Cython生成的代码。因此，我决定编写以下小示例：

import numpy as np
cimport numpy as np
cimport cython

@cython.boundscheck(False)
@cython.wraparound(False)
cpdef int f(np.ndarray[int, ndim = 1] f):
    cdef int array_length =  f.shape[0]
    cdef int sum = 0
    cdef int k
    for k in range(array_length):
        sum += f[k]
    return sum

我知道有Numpy函数可以完成这项工作，但是我想编写一个简单的代码来了解Cython的功能。事实证明，生成的代码是：

from distutils.core import setup
from Cython.Build import cythonize

setup(ext_modules = cythonize("sum.pyx"))

并致电：

python setup.py build_ext --inplace

生成一个看起来像这样的C代码循环：

for (__pyx_t_2 = 0; __pyx_t_2 < __pyx_t_1; __pyx_t_2 += 1) {
  __pyx_v_sum = __pyx_v_sum + (*(int *)((char *) 
    __pyx_pybuffernd_f.rcbuffer->pybuffer.buf +
    __pyx_t_2 * __pyx_pybuffernd_f.diminfo[0].strides)));
}

此代码的主要问题是，编译器在编译时不知道__pyx_pybuffernd_f.diminfo[0].strides数组的元素在内存中是否紧密在一起。没有这些信息，编译器将无法有效地向量化。

Cython有办法做到这一点吗？

Answer 1

您的代码中有两个问题（使用选项-a使其可见）：

1. numpy数组的索引效率不高
2. 你忘int了cdef sum=0

考虑到这一点，我们得到：

cpdef int f(np.ndarray[np.int_t] f):  ##HERE
    assert f.dtype == np.int
    cdef int array_length =  f.shape[0]
    cdef int sum = 0                  ##HERE
    cdef int k
    for k in range(array_length):
        sum += f[k]
    return sum

对于循环，以下代码：

int __pyx_t_5;
int __pyx_t_6;
Py_ssize_t __pyx_t_7;
....
__pyx_t_5 = __pyx_v_array_length;
for (__pyx_t_6 = 0; __pyx_t_6 < __pyx_t_5; __pyx_t_6+=1) {
   __pyx_v_k = __pyx_t_6;
   __pyx_t_7 = __pyx_v_k;
   __pyx_v_sum = (__pyx_v_sum + (*__Pyx_BufPtrStrided1d(__pyx_t_5numpy_int_t *, __pyx_pybuffernd_f.rcbuffer->pybuffer.buf, __pyx_t_7, __pyx_pybuffernd_f.diminfo[0].strides)));

}

这样做还不错，但是对于优化器而言，却不如人类编写的普通代码那么容易。正如您已经指出的那样，__pyx_pybuffernd_f.diminfo[0].strides在编译时未知，这会阻止矢量化。

但是，使用类型化的内存视图时，您会得到更好的结果，即：

cpdef int mf(int[::1] f):
    cdef int array_length =  len(f)
...

这导致了不透明的C代码-至少是我的编译器可以更好地优化：

 __pyx_t_2 = __pyx_v_array_length;
  for (__pyx_t_3 = 0; __pyx_t_3 < __pyx_t_2; __pyx_t_3+=1) {
    __pyx_v_k = __pyx_t_3;
    __pyx_t_4 = __pyx_v_k;
    __pyx_v_sum = (__pyx_v_sum + (*((int *) ( /* dim=0 */ ((char *) (((int *) __pyx_v_f.data) + __pyx_t_4)) ))));
  }

这里最关键的是，我们要让cython清楚，内存是连续的，即int[::1]与int[:]numpy数组所看到的相比，stride!=1必须考虑可能的情况。

在这种情况下，用Cython生成的C代码导致同一汇编作为代码我会写。正如crisb指出的那样，相加-march=native会导致向量化，但是在这种情况下，两个函数的汇编器将再次略有不同。

但是，以我的经验，编译器通常会遇到一些问题，无法优化由cython创建的循环，并且/或者容易遗漏一个细节，从而阻止生成真正好的C代码。因此，我处理工作循环的策略是用纯C语言编写它们，并使用cython来包装/访问它们-
通常会更快一些，因为也可以使用专用的编译器标志来捕获此代码而不会影响整个Cython-模块。