0

我开始研究 Fortran 库 (BLAS/LAPACK) 的 SciPy 接口,可以在这里看到:Calling BLAS / LAPACK directly using the SciPy interface and Cython并想出了一个解决方案,但不得不求助于使用numpy.zeros它实际上杀死了任何速度直接调用 Fortran 代码的好处。问题是 Fortran 代码需要一个 0 值的输出矩阵(它在内存中就地处理矩阵)以匹配 Numpy 版本(np.outer)。

所以我有点困惑,因为 Python 中的 1000x1000 零矩阵只需要 8us(使用 %timeit 或 0.008ms)那么为什么添加 Cython 代码会杀死运行时,注意到我也在 memoryview 上创建它?(基本上 1000 x 1000 矩阵乘法从 3ms 到 8ms 左右)。然后在 SO 上搜索后,我在其他地方找到了一个使用memset最快的数组更改机制的解决方案。因此,我将引用帖子中的整个代码重写为更新的memoryview格式,并得到了类似这样的内容(Cythoncyblas.pyx文件):

import cython
import numpy as np
cimport numpy as np
from libc.string cimport memset #faster than np.zeros

REAL = np.float64
ctypedef np.float64_t REAL_t
ctypedef np.uint64_t  INT_t

cdef int ONE = 1
cdef REAL_t ONEF = <REAL_t>1.0

@cython.boundscheck(False)
@cython.wraparound(False)
@cython.cdivision(True)
cpdef outer_prod(double[::1] _x, double[::1] _y, double[:, ::1] _output):

    cdef int M = _y.shape[0]
    cdef int N = _x.shape[0]
    memset(&_output[0,0], 0, M*N)

    with nogil:
        dger(&M, &N, &ONEF, &_y[0], &ONE, &_x[0], &ONE, &_output[0,0], &M)

测试脚本

import numpy as np;
from cyblas import outer_prod;
a=np.random.randint(0,100, 1000);
b=np.random.randint(0,100, 1000);
a=a.astype(np.float64)
b=b.astype(np.float64)
cy_outer=np.zeros((a.shape[0],b.shape[0]));
np_outer=np.zeros((a.shape[0],b.shape[0]));

%timeit outer_prod(a,b, cy_outer)
%timeit np.outer(a,b, np_outer)

因此,这解决了我将输出矩阵值重置为第 125 行的问题,该问题仍然存在(已解决见下文)。我认为将memset长度参数设置为 M*N 会清除内存中的 1000*1000 ,但显然不是。

有人知道如何使用 将整个输出矩阵重置为 0memset吗?非常感激。

4

1 回答 1

3

[更新 - 修复是:它需要#bytes 而不仅仅是数组大小M*N,即M*N*variable_bytes作为长度输入。前面的结果可以看出这里的逻辑:第 125 行是它停止的位置 *8 字节 = 1000 行,因此回答了这个问题。仍然指标不是很好:100 个循环,最好的 3:每个循环 5.41 毫秒(cython) 100 个循环,最好的 3:每个循环 3.95 毫秒(numpy)但仍然解决了。对上述代码的更改是添加:cdef variable_bytes = np.dtype(REAL).itemsize #added to get bytes for memset, after REAL is defined, in this case 8 bytes 然后当您调用 memset 时:memset(&_output[0,0], 0, M*N*variable_bytes) # gives the same output as np.zeros function我现在可以看到的唯一加快速度的地方是prange在大型矩阵上使用 OpenMP 语句,但还有待观察。

于 2017-07-02T00:54:41.407 回答