cython memoryview 不比 ndarray 快

Posted

技术标签:

【中文标题】cython memoryview 不比 ndarray 快【英文标题】:cython memoryview not faster than ndarray 【发布时间】:2015-10-05 16:33:09 【问题描述】:

我有一个用常规numpy ndarray 编写的函数,另一个用typed memoryview 编写的函数。但是,我无法让memoryview 版本比普通版本运行得更快(不像memoryview benchmarks 等许多博客)。

任何提高 memoryview 代码速度与 numpy 替代方案相比的指针/建议将不胜感激! ...或者...如果有人能指出为什么 memoryview 版本没有比常规 numpy 版本快多少的明显原因

在下面的代码中有两个函数,它们都接受两个向量bixi并返回一个矩阵。第一个函数shrink_correl 是常规的numpy 版本,第二个函数shrink_correl2 是memoryview 替代方案(让文件为sh_cor.pyx)。

# cython: boundscheck=False
# cython: wraparound=False
# cython: cdivision=True

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

# -- ***this is the Regular Cython version*** -
cpdef ar[double, ndim=2, mode='c'] shrink_correl(ar[double, ndim=1, mode='c'] bi, ar[double, ndim=1, mode='c'] xi):
    cdef:
        int n_ = xi.shape[0]
        int n__ = int(n_*(n_-1)/2)
        ar[double, ndim=2, mode='c'] f = np.zeros([n__, n_+1])
        int x__ = 0
        ar[double, ndim=2, mode='c'] f1 = np.zeros([n_, n_+1])
        ar[double, ndim=2, mode='c'] f2 = np.zeros([n__, n_+1])
        ar[double, ndim=1, mode='c'] g = np.zeros(n_+1)
        ar[double, ndim=1, mode='c'] s = np.zeros(n__)
        ar[double, ndim=2, mode='c'] cori_ = np.zeros([n_, n_])
        Py_ssize_t j, k

    with nogil:
        for j in range(0, n_-1):
            for k in range(j+1, n_):
                x__ += 1
                f[x__-1, j] = bi[k]*xi[k]*1000
                f[x__-1, k] = bi[j]*xi[j]*1000
    f1 = np.dot(np.transpose(f), f)      
    with nogil:
        for j in range(0, n_):
            f1[n_, j] = xi[j]*1000
            f1[j, n_] = f1[n_, j]
    f2 = np.dot(f, np.linalg.inv(f1))
    with nogil:
        for j in range(0, n_):
            g[j] = -bi[j]*xi[j]*1000

    s = np.dot(f2, g)

    with nogil:
        for j in range(0, n_):
            cori_[j, j] = 1.0
    x__ = 0

    with nogil:
        for j in range(0, n_-1):
            for k in range(j+1, n_):
                x__ += 1
                cori_[j, k] = s[x__-1]
                cori_[k, j] = cori_[j, k]
    return cori_

# -- ***this is the MemoryView Cython version*** -    
cpdef ar[double, ndim=2, mode='c'] shrink_correl2(double[:] bi, double[:] xi):
    cdef:
        int n_ = xi.shape[0]
        int n__ = int(n_*(n_-1)/2)
        double[:, ::1] f = np.zeros([n__, n_+1])
        int x__ = 0
        double[:, ::1] f1 = np.zeros([n_, n_+1])
        double[:, ::1] f2 = np.zeros([n__, n_+1])
        double[:] g = np.zeros(n_+1)
        double[:] s = np.zeros(n__)
        double[:, ::1] cori_ = np.zeros([n_, n_])
        ar[double, ndim=2, mode='c'] cori__ = np.zeros([n_, n_])
        Py_ssize_t j, k
    with nogil:
        for j in range(0, n_-1):
            for k in range(j+1, n_):
                x__ += 1
                f[x__-1, j] = bi[k]*xi[k]*1000
                f[x__-1, k] = bi[j]*xi[j]*1000
    f1 = np.dot(np.transpose(f), f)      
    with nogil:
        for j in range(0, n_):
            f1[n_, j] = xi[j]*1000
            f1[j, n_] = f1[n_, j]
    f2 = np.dot(f, np.linalg.inv(f1))
    with nogil:
        for j in range(0, n_):
            g[j] = -bi[j]*xi[j]*1000

    s = np.dot(f2, g)

    with nogil:
        for j in range(0, n_):
            cori_[j, j] = 1.0
    x__ = 0

    with nogil:
        for j in range(0, n_-1):
            for k in range(j+1, n_):
                x__ += 1
                cori_[j, k] = s[x__-1]
                cori_[k, j] = cori_[j, k]
    cori__[:, :] = cori_
    return cori__

这是使用以下setup.py代码编译的

from distutils.core import setup
from distutils.extension import Extension
from Cython.Distutils import build_ext
import numpy as np
import os

ext_modules = [Extension('sh_cor', ['sh_cor.pyx'], include_dirs=[np.get_include(),
                                                                 os.path.join(np.get_include(), 'numpy')],
                         define_macros=[('NPY_NO_DEPRECATED_API', None)],
                         extra_compile_args=['-O3', '-march=native', '-ffast-math', '-flto'],
                         libraries=['m']
                         )]

setup(
    name="Sh Cor",
    cmdclass='build_ext': build_ext,
    ext_modules=ext_modules
)

用于测试速度的代码是

import numpy as np
import sh_cor  # this the library created by the setup.py file
import time

b = np.random.random(400)
b = b/np.sum(b)

x = np.random.random(400)-0.5

n = 10 

t0 = time.time()
for i in range(n):
    v1 = sh_cor.shrink_correl(b, x)
t1 = time.time()
print((t1-t0)/n)

t0 = time.time()
for i in range(n):
    v2 = sh_cor.shrink_correl2(b, x)
t1 = time.time()
print((t1-t0)/n)

我的电脑上的输出是:

0.7070999860763549   # regular numpy
0.6726999998092651   # memoryview

使用 memoryview(在上面的代码中)只给我 5% 的速度提升(不像博客中的巨大速度提升)。

【问题讨论】:

有什么不同?只是函数参数的规范?当我看到这么多没有解释的代码时,我的眼睛都呆住了。 是的,差不多。大部分功能只是循环 我对@9​​87654322@ 的印象是,当您专注于低级别操作时,您将从memoryviews 中获得最大收益。但我对你的代码研究得还不够多,无法知道你的情况。 另外 - 您链接到的博客所讨论的优势是在您创建切片时,例如double [:] x = some_2d_array[:,0],显然与ndarrays 相比效率低下。由于您在这里不这样做,因此您可能不会获得该好处。直接的元素访问可能非常相似(查看此处的 C 代码进行检查)。 我在两者之间看到的一个区别是,在您的第一个函数输入 ndarrays 以及 sgc-contiguous 而在第二个函数中它们不是。使用double [::1] 而不是double[:]。除此之外,您的大多数操作都需要相同的 numpy C-API,因此可以肯定您不会看到任何加速。如果你让 cython 吐出 html 文件 (cython -a),你基本上会看到类似的 C 代码。 【参考方案1】:

@uday 给我大约一周的时间,因为我的电脑少了,但这里是加快速度让你开始的地方:1)而不是使用 np.transpose 附加 gil 创建一个与你想要转置的内容相同的内存视图在任何循环之前(即您将变量 f 声明为不需要 gil 的内存视图,只需在该 f_t 上创建一个视图,即 cdef double[:, ::1] f_T = np.transpose(f) 或只是 =f.T

2) 这一步有点棘手,因为您需要 np.dot 的 C/C++ 样式包装器版本(所以在这种情况下,请确保对 dgemm 函数的调用在其上方是 with nogil: 并缩进下一行释放具有 4 个空格缩进 SO 的 gil 的函数需要):https://gist.github.com/pv/5437087。该示例看起来可行(尽管您必须将包含 f2pyptr.h 文件保存并将其放在正在构建项目的位置;我还怀疑您应该添加 cimport numpy as np);如果不是,它需要模组,你可以看到我在另一篇文章中所做的那样:Calling BLAS / LAPACK directly using the SciPy interface and Cython (pointer issue?)/- also how to add MKL 然后你需要在顶部添加from cython.parallel cimport prange 并将所有循环从range 更改为prange 并成为确保您的所有prange 部分都是nogil 并且所有变量都是cdef 在操作之前声明的。此外,您必须将-openmp 添加到编译器参数中的 setup.py 以及指向其包含库的链接。如果您需要澄清,请提出更多问题。这并不像应有的那么容易,但是通过一些指导变得非常简单。基本上,一旦您的 setup.py 被修改为包含它将继续工作的所有内容。

3) 尽管可能最容易修复 - 摆脱该列表。如果您需要文本和数据,请将其设为 numpy 数组或 pandas 数据框。每当我使用数据列表时,速度的下降都是令人难以置信的。

【讨论】:

以上是关于cython memoryview 不比 ndarray 快的主要内容,如果未能解决你的问题,请参考以下文章

《Cython系列》7. Cythonnumpy以及类型化memoryview

Cython:从类型化的memoryview数组中读取值的最快方法是什么?

HTML Insertar Flash deformaestándar

以类型化 memoryview 作为成员的结构定义

笔记bytes和bytearray还有memoryview

笔记bytes和bytearray还有memoryview