如何在使用C扩展numpy时考虑列连续数组

Question

我有一个C函数来规范化日志空间中的数组行(这可以防止数字下溢).

我的C函数的原型如下:

void normalize_logspace_matrix(size_t nrow, size_t ncol, double* mat);

你可以看到它需要一个指向数组的指针并在适当的位置修改它.C代码当然假设数据被保存为C连续数组,即行连续.

我使用Cython(导入并cdef extern from省略)将函数包装如下:

def normalize_logspace(np.ndarray[np.double_t, ndim=2] mat):
    cdef Py_ssize_t n, d
    n = mat.shape[0]
    d = mat.shape[1]
    normalize_logspace_matrix(n, d, <double*> mat.data)
    return mat

大多数情况下,numpy-arrays是行连续的,函数运行正常.但是,如果先前已转换了numpy-array,则不会复制数据,而只返回数据的新视图.在这种情况下,我的函数失败,因为数组不再是行连续的.

我可以通过将数组定义为具有Fortran连续顺序来解决这个问题,这样在转置后它将是C连续的:

A = np.array([some_func(d) for d in range(D)], order='F').T
A = normalize_logspace(A)

显然,这非常容易出错,用户必须注意数组的顺序是正确的,这是用户在Python中不需要关心的.

使用行和列连续数组的最佳方法是什么？我假设在Cython中进行某种数组顺序检查是可行的方法.当然,我更喜欢不需要将数据复制到新数组的解决方案,但我几乎认为这是必要的.

Answer 1

如果要在不复制的情况下支持C和Fortran顺序的数组,则C函数需要足够灵活以支持两种顺序.这可以通过将NumPy数组的步幅传递给C函数来实现:将原型更改为

void normalize_logspace_matrix(size_t nrow, size_t ncol, 
                               size_t rowstride, size_t colstride,
                               double* mat);

和Cython打电话给

def normalize_logspace(np.ndarray[np.double_t, ndim=2] mat):
    cdef Py_ssize_t n, d, rowstride, colstride
    n = mat.shape[0]
    d = mat.shape[1]
    rowstride = mat.strides[0] // mat.itemsize
    colstride = mat.strides[1] // mat.itemsize
    normalize_logspace_matrix(n, d, rowstride, colstride, <double*> mat.data)
    return mat

然后,mat[row*ncol + col]通过mat[row*rowstride + col*colstride] 替换C代码中的每个出现.