欺骗SciPy进行稀疏对角系统的卷积

Question

我正在尝试把一些代码转换成Python，但我发现SciPy库在处理对角线系统时，稀疏对角操作有点问题。

举个例子，下面的代码可以写成逐像素的卷积，在我用C++实现的时候速度非常快。加上重叠后，差不多就是内存访问的时间。我本以为Python会知道这是对角线系统。

但是当我在Python中运行时，我的终端占用了太多系统资源，差点把我的系统搞崩了。比如，Matlab版本运行得就快得多。

import numpy as np
from scipy import sparse
print(np.version.version)
color=imread('lena.png')
gray=mean(color,2)
[h,w]=gray.shape
img=gray.flatten()
ne=h*w;
L=np.ones(ne);
M=.5*np.ones(ne);
R=np.ones(ne);
Diags=[L,M,R]
mtx = sparse.spdiags(Diags, [-1,0,1], ne, ne);
#blured=np.dot(mtx,img) Dies

我知道我可以把它重写成一个遍历像素的'for'循环，但有没有办法在保持稀疏数据结构的同时，提升性能呢？

Answer 1

使用 mtx.dot 来代替 np.dot，因为

blured = mtx.dot(img)

或者直接使用

blured = mtx * img # where mtx is sparse and img is `dense` or `sparse` 

即使 np.dot 的两个参数中有一个是 sparse（稀疏矩阵），它仍然会把这两个参数都当作 dense（密集矩阵）来处理。这样做可能会导致 MemoryError（内存错误）。