Scipy稀疏矩阵 - 密集向量乘法性能 - 块与大矩阵

ali*_*i_p 7 python performance numpy scipy sparse-matrix

我有一些scipy稀疏矩阵(目前采用CSR格式),我需要与密集的numpy 1D向量相乘.该向量称为G:

print G.shape, G.dtype
(2097152,) complex64
Run Code Online (Sandbox Code Playgroud)

每个稀疏矩阵都有形状(16384,2097152)并且非常稀疏.密度约为4.0e-6.我有一个名为100的这些稀疏矩阵的列表spmats.

我可以很容易地将每个矩阵相乘G:

res = [spmat.dot(G) for spmat in spmats]
Run Code Online (Sandbox Code Playgroud)

这导致(16384,)如预期的形状的密集矢量列表.

我的应用程序是相当性能的,因此我尝试了替代方法,即首先将所有稀疏矩阵连接成一个大的sparce矩阵,然后只使用一次这样的调用dot():

import scipy.sparse as sp
SPMAT = sp.vstack(spmats, format='csr')
RES = SPMAT.dot(G)
Run Code Online (Sandbox Code Playgroud)

这导致一个长矢量RES具有形状(1638400,),并且是res上述所有结果矢量的连接形式,如预期的那样.我检查过结果是一样的.

也许我完全错了,但我预计第二种情况应该比第一种情况要快,因为numpy调用,内存分配,python对象的创建,python循环等等都少得多.我不关心时间需要连接稀疏矩阵,只需要计算结果的时间.%timeit然而,根据:

%timeit res = [spmat.dot(G) for spmat in spmats]
10 loops, best of 3: 91.5 ms per loop
%timeit RES = SPMAT.dot(G)
1 loops, best of 3: 389 ms per loop
Run Code Online (Sandbox Code Playgroud)

我已经检查过,在任何一个操作中我都没有内存不足,似乎没有任何可疑的事情发生.我疯了,还是真的很奇怪?这是否意味着所有稀疏矩阵矢量产品应该以块的形式完成,一次几行,以使它们更快?据我所知,具有密集向量的稀疏矩阵乘法时间在非零元素的数量上应该是线性的,在上述两种情况下不变.有什么可以产生这样的差异?

我正在使用EPD7.3在单核linux机器上运行4GB内存

编辑:

这是一个为我重现问题的小例子:

import scipy.sparse as sp
import numpy as n

G = n.random.rand(128**3) + 1.0j*n.random.rand(128**3)

spmats = [sp.rand (128**2, 128**3, density = 4e-6, format = 'csr', dtype=float64) for i in range(100)]
SPMAT = sp.vstack(spmats, format='csr')

%timeit res = [spmat.dot(G) for spmat in spmats]
%timeit RES = SPMAT.dot(G)
Run Code Online (Sandbox Code Playgroud)

我明白了:

1 loops, best of 3: 704 ms per loop
1 loops, best of 3: 1.34 s per loop
Run Code Online (Sandbox Code Playgroud)

在这种情况下,性能差异并不像我自己的具有某种结构的稀疏矩阵(可能是因为缓存)那么大,但连接矩阵的情况更糟.

我尝试了scipy 10.1和12.0.

ali*_*i_p 4

我还没有找到问题中提到的奇怪行为的原因,但是我找到了一种显着加快计算速度的方法,这可能对其他人有用。

因为在我的特定情况下,我正在计算 float32 稀疏矩阵和复数 64 稠密向量的乘积,所以我可以分别将实部和虚部相乘。这使我的速度提高了 4 倍。

这需要 2.35 秒SPMAT.shape == (16384000, 2097152)

RES = SPMAT.dot(G)
Run Code Online (Sandbox Code Playgroud)

虽然这只需要 541ms:

RES = n.zeros((SPMAT.shape[0],),dtype=complex64)
RES.real = SPMAT.dot(G.real); RES.imag = SPMAT.dot(G.imag)
Run Code Online (Sandbox Code Playgroud)

结果是相同的。我认为n.zeros预分配可能没有必要,但我不知道还能怎么做。