我有一个用常规编写的函数numpy ndarray
,另一个用typed memoryview
. 但是,我无法让该memoryview
版本比普通版本运行得更快(不像许多博客,例如memoryview benchmarks)。
任何提高 memoryview 代码速度与 numpy 替代方案的指针/建议将不胜感激!...或者...如果有人能指出任何明显的原因,为什么 memoryview 版本并不比常规 numpy 版本快多少
在下面的代码中有两个函数,它们都接受两个向量bi
并xi
返回一个矩阵。第一个功能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% 的速度提升(不像博客中的巨大速度提升)。