C数组与NumPy数组

Question

C数组与NumPy数组

ric*_*izy 12 c python arrays numpy cython

在性能方面(代数运算,查找,缓存等),C数组(可以作为C数组公开,或者cython.view.array[Cython数组],或上述两者的内存视图)和NumPy数组(在Cython中应该没有Python开销)

编辑:

我应该提一下,在NumPy数组中使用Cython进行静态类型化,dtypes是NumPy编译时数据类型(例如cdef np.int_t或者cdef np.float32_t),C语言中的类型是C等价物(cdef int_t和cdef float)

EDIT2:

以下是Cython Memoryview文档中的示例,以进一步说明我的问题:

from cython.view cimport array as cvarray
import numpy as np

# Memoryview on a NumPy array
narr = np.arange(27, dtype=np.dtype("i")).reshape((3, 3, 3))
cdef int [:, :, :] narr_view = narr

# Memoryview on a C array
cdef int carr[3][3][3]
cdef int [:, :, :] carr_view = carr

# Memoryview on a Cython array
cyarr = cvarray(shape=(3, 3, 3), itemsize=sizeof(int), format="i")
cdef int [:, :, :] cyarr_view = cyarr

Run Code Online (Sandbox Code Playgroud)

坚持C arrayvs a Cython arrayvs vs NumPy array？之间有什么区别吗？

Answer 1

Ian*_*anH 21

我对此的了解仍然不完善,但这可能会有所帮助.我运行了一些非正式的基准来展示每种阵列类型的优点,并对我发现的内容感兴趣.

虽然这些数组类型在很多方面都有所不同,但如果您使用大型数组进行大量计算,则应该能够从其中任何一个获得类似的性能,因为逐项访问应该在整个板上大致相同.

NumPy数组是使用Python的C API实现的Python对象.NumPy数组确实提供了C级的API,但它们不能独立于Python解释器创建.它们特别有用,因为NumPy和SciPy中提供了所有不同的数组操作例程.

Cython内存视图也是一个Python对象,但它是一个Cython扩展类型.它似乎不是设计用于纯Python,因为它不是Cython的一部分,可以直接从Python导入,但您可以从Cython函数返回Python的视图.您可以在https://github.com/cython/cython/blob/master/Cython/Utility/MemoryView.pyx上查看实现.

AC数组是C语言中的本机类型.它像索引一样被索引,但数组和指针是不同的.在http://c-faq.com/aryptr/index.html上对此有一些很好的讨论. 它们可以在堆栈上分配,并且C编译器更容易优化,但是它们将更难以访问用Cython.我知道你可以从已经由其他程序动态分配的内存中创建一个NumPy数组,但是这样看起来要困难得多.Travis Oliphant在http://blog.enthought.com/python/numpy-arrays-with-pre-allocated-memory/上发布了一个例子. 如果你在程序中使用C数组或指针进行临时存储,他们应该工作得非常好对你好它们不会像切片或任何其他类型的矢量化计算那样方便,因为你必须使用显式循环来完成所有事情,但它们应该更快地分配和释放,并且应该为速度提供良好的基线.

Cython还提供了一个数组类.它看起来像是为内部使用而设计的.复制内存视图时会创建实例.请参阅http://docs.cython.org/src/userguide/memoryviews.html#view-cython-arrays

在Cython中,您还可以分配内存并索引指针,以便像分组一样处理分配的内存.请参阅http://docs.cython.org/src/tutorial/memory_allocation.html

以下是一些基准测试,表明索引大型数组的性能有些类似.这是Cython文件.

from numpy cimport ndarray as ar, uint64_t
cimport cython
import numpy as np

@cython.boundscheck(False)
@cython.wraparound(False)
def ndarr_time(uint64_t n=1000000, uint64_t size=10000):
    cdef:
        ar[uint64_t] A = np.empty(n, dtype=np.uint64)
        uint64_t i, j
    for i in range(n):
        for j in range(size):
            A[j] = n

def carr_time(uint64_t n=1000000):
    cdef:
        ar[uint64_t] A = np.empty(n, dtype=np.uint64)
        uint64_t AC[10000]
        uint64_t a
        int i, j
    for i in range(n):
        for j in range(10000):
            AC[j] = n

@cython.boundscheck(False)
@cython.wraparound(False)
def ptr_time(uint64_t n=1000000, uint64_t size=10000):
    cdef:
        ar[uint64_t] A = np.empty(n, dtype=np.uint64)
        uint64_t* AP = &A[0]
        uint64_t a
        int i, j
    for i in range(n):
        for j in range(size):
            AP[j] = n

@cython.boundscheck(False)
@cython.wraparound(False)
def view_time(uint64_t n=1000000, uint64_t size=10000):
    cdef:
        ar[uint64_t] A = np.empty(n, dtype=np.uint64)
        uint64_t[:] AV = A
        uint64_t i, j
    for i in range(n):
        for j in range(size):
            AV[j] = n

Run Code Online (Sandbox Code Playgroud)

使用我们获得的IPython来计时

%timeit -n 10 ndarr_time()
%timeit -n 10 carr_time()
%timeit -n 10 ptr_time()
%timeit -n 10 view_time()

10 loops, best of 3: 6.33 s per loop
10 loops, best of 3: 3.12 s per loop
10 loops, best of 3: 6.26 s per loop
10 loops, best of 3: 3.74 s per loop

Run Code Online (Sandbox Code Playgroud)

这些结果让我觉得有点奇怪,考虑到效率:数组与指针,数组不太可能明显快于指针.似乎某种编译器优化使得纯C数组和类型化内存视图更快.我尝试关闭C编译器上的所有优化标志并获得时间

1 loops, best of 3: 25.1 s per loop
1 loops, best of 3: 25.5 s per loop
1 loops, best of 3: 32 s per loop
1 loops, best of 3: 28.4 s per loop

Run Code Online (Sandbox Code Playgroud)

在我看来,除了C数组和Cython内存视图似乎更容易让编译器进行优化之外,逐项访问几乎完全相同.

我在前面发现的这两篇博文中可以看到更多评论:http: //jakevdp.github.io/blog/2012/08/08/memoryview-benchmarks/ http://jakevdp.github.io /博客/ 2012/08/16/memoryview-基准-2 /

在第二篇博客文章中,他评论了如果内存视图切片被内联,它们可以提供类似于指针算术的速度.我在一些自己的测试中注意到,显式内联使用内存视图切片的函数并不总是必要的.作为一个例子,我将计算两行数组的每个组合的内积.

from numpy cimport ndarray as ar
cimport cython
from numpy import empty

# An inlined dot product
@cython.boundscheck(False)
@cython.wraparound(False)
cdef inline double dot_product(double[:] a, double[:] b, int size):
    cdef int i
    cdef double tot = 0.
    for i in range(size):
        tot += a[i] * b[i]
    return tot

# non-inlined dot-product
@cython.boundscheck(False)
@cython.wraparound(False)
cdef double dot_product_no_inline(double[:] a, double[:] b, int size):
    cdef int i
    cdef double tot = 0.
    for i in range(size):
        tot += a[i] * b[i]
    return tot

# function calling inlined dot product
@cython.boundscheck(False)
@cython.wraparound(False)
def dot_rows_slicing(ar[double,ndim=2] A):
    cdef:
        double[:,:] Aview = A
        ar[double,ndim=2] res = empty((A.shape[0], A.shape[0]))
        int i, j
    for i in range(A.shape[0]):
        for j in range(A.shape[0]):
            res[i,j] = dot_product(Aview[i], Aview[j], A.shape[1])
    return res

# function calling non-inlined version
@cython.boundscheck(False)
@cython.wraparound(False)
def dot_rows_slicing_no_inline(ar[double,ndim=2] A):
    cdef:
        double[:,:] Aview = A
        ar[double,ndim=2] res = empty((A.shape[0], A.shape[0]))
        int i, j
    for i in range(A.shape[0]):
        for j in range(A.shape[0]):
            res[i,j] = dot_product_no_inline(Aview[i], Aview[j], A.shape[1])
    return res

# inlined dot product using numpy arrays
@cython.boundscheck(False)
@cython.boundscheck(False)
cdef inline double ndarr_dot_product(ar[double] a, ar[double] b):
    cdef int i
    cdef double tot = 0.
    for i in range(a.size):
        tot += a[i] * b[i]
    return tot

# non-inlined dot product using numpy arrays
@cython.boundscheck(False)
@cython.boundscheck(False)
cdef double ndarr_dot_product_no_inline(ar[double] a, ar[double] b):
    cdef int i
    cdef double tot = 0.
    for i in range(a.size):
        tot += a[i] * b[i]
    return tot

# function calling inlined numpy array dot product
@cython.boundscheck(False)
@cython.wraparound(False)
def ndarr_dot_rows_slicing(ar[double,ndim=2] A):
    cdef:
        ar[double,ndim=2] res = empty((A.shape[0], A.shape[0]))
        int i, j
    for i in range(A.shape[0]):
        for j in range(A.shape[0]):
            res[i,j] = ndarr_dot_product(A[i], A[j])
    return res

# function calling nun-inlined version for numpy arrays
@cython.boundscheck(False)
@cython.wraparound(False)
def ndarr_dot_rows_slicing_no_inline(ar[double,ndim=2] A):
    cdef:
        ar[double,ndim=2] res = empty((A.shape[0], A.shape[0]))
        int i, j
    for i in range(A.shape[0]):
        for j in range(A.shape[0]):
            res[i,j] = ndarr_dot_product(A[i], A[j])
    return res

# Version with explicit looping and item-by-item access.
@cython.boundscheck(False)
@cython.wraparound(False)
def dot_rows_loops(ar[double,ndim=2] A):
    cdef:
        ar[double,ndim=2] res = empty((A.shape[0], A.shape[0]))
        int i, j, k
        double tot
    for i in range(A.shape[0]):
        for j in range(A.shape[0]):
            tot = 0.
            for k in range(A.shape[1]):
                tot += A[i,k] * A[j,k]
            res[i,j] = tot
    return res

Run Code Online (Sandbox Code Playgroud)

我们看到这些时间

A = rand(1000, 1000)
%timeit dot_rows_slicing(A)
%timeit dot_rows_slicing_no_inline(A)
%timeit ndarr_dot_rows_slicing(A)
%timeit ndarr_dot_rows_slicing_no_inline(A)
%timeit dot_rows_loops(A)

1 loops, best of 3: 1.02 s per loop
1 loops, best of 3: 1.02 s per loop
1 loops, best of 3: 3.65 s per loop
1 loops, best of 3: 3.66 s per loop
1 loops, best of 3: 1.04 s per loop

Run Code Online (Sandbox Code Playgroud)

结果与显式内联一样快,没有它.在这两种情况下,类型化的内存视图都与没有切片的函数版本相当.

在博客文章中,他必须编写一个特定的示例来强制编译器不要内联函数.似乎一个体面的C编译器(我正在使用MinGW)能够处理这些优化而不被告知内联某些功能.对于在Cython模块中的函数之间传递数组切片,Memoryview可以更快,即使没有显式内联也是如此.

然而,在这种特殊情况下,即使将环路推到C,也不会真正达到通过正确使用矩阵乘法可以实现的速度.BLAS仍然是做这样事情的最佳方式.

%timeit A.dot(A.T)
10 loops, best of 3: 25.7 ms per loop

Run Code Online (Sandbox Code Playgroud)

还有从NumPy数组到内存视图的自动转换,如

cimport cython

@cython.boundscheck(False)
@cython.wraparound(False)
def cysum(double[:] A):
    cdef tot = 0.
    cdef int i
    for i in range(A.size):
        tot += A[i]
    return tot

Run Code Online (Sandbox Code Playgroud)

唯一的问题是,如果您希望函数返回NumPy数组,则必须使用np.asarray将内存视图对象再次转换为NumPy数组.这是一个相对便宜的操作,因为内存视图符合http://www.python.org/dev/peps/pep-3118/

结论

类型化的内存视图似乎是NumPy数组的可行替代方案,供Cython模块内部使用.对于内存视图,数组切片会更快,但是没有像NumPy数组那样为内存视图编写的函数和方法.如果您不需要调用一堆NumPy数组方法并希望轻松进行数组切片,则可以使用内存视图代替NumPy数组.如果您需要给定数组的数组切片和 NumPy功能,则可以创建一个指向与NumPy数组相同内存的内存视图.然后,您可以使用视图在函数和数组之间传递切片以调用NumPy函数.这种方法仍然有限,但如果您使用单个阵列进行大部分处理,它将很有效.

C数组和/或动态分配的内存块对于中间计算可能很有用,但它们不容易传回Python以便在那里使用.在我看来,动态分配多维C数组也更麻烦.我所知道的最好的方法是分配一个大的内存块然后使用整数运算来索引它,就好像它是一个多维数组.如果您想在运行中轻松分配数组,这可能是一个问题.另一方面,C阵列的分配时间可能要快得多.其他数组类型的设计几乎同样快,更方便,所以除非有令人信服的理由,否则我建议使用它们.

更新:正如@Veedrac的回答中提到的,你仍然可以将Cython内存视图传递给大多数NumPy函数.当你这样做时,NumPy通常必须创建一个新的NumPy数组对象来处理内存视图,所以这会有点慢.对于大型阵列,效果可以忽略不计.np.asarray无论数组大小如何,对内存视图的调用都会相对较快.但是,为了证明这种效果,这是另一个基准:

Cython文件:

def npy_call_on_view(npy_func, double[:] A, int n):
    cdef int i
    for i in range(n):
        npy_func(A)

def npy_call_on_arr(npy_func, ar[double] A, int n):
    cdef int i
    for i in range(n):
        npy_func(A)

Run Code Online (Sandbox Code Playgroud)

在IPython中:

from numpy.random import rand
A = rand(1)
%timeit npy_call_on_view(np.amin, A, 10000)
%timeit npy_call_on_arr(np.amin, A, 10000)

Run Code Online (Sandbox Code Playgroud)

输出:

10 loops, best of 3: 282 ms per loop
10 loops, best of 3: 35.9 ms per loop

Run Code Online (Sandbox Code Playgroud)

我试图选择一个可以很好地显示这种效果的例子.除非涉及对相对较小的数组进行许多NumPy函数调用,否则这不应该改变整个时间.请记住,无论我们调用NumPy的方式如何,都会发生Python函数调用.

这仅适用于NumPy中的功能.大多数数组方法都不适用于内存视图(某些属性仍然是,如size和,shape和T).例如A.dot(A.T),NumPy阵列就会变成np.dot(A, A.T).

归档时间：	11 年，6 月前
查看次数：	5471 次
最近记录：	9 年，2 月前