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

Question

我有许多scipy稀疏矩阵（目前为CSR格式），需要与密集的numpy 1D向量相乘。该向量称为G：

print G.shape, G.dtype
(2097152,) complex64

每个稀疏矩阵都具有形状（163842097152），并且非常稀疏。密度约为4.0e-6。我有一个包含100个稀疏矩阵的列表，称为spmats。

我可以轻松地将每个矩阵与G相乘，如下所示：

res = [spmat.dot(G) for spmat in spmats]

这将产生一个形状密集向量列表（16384，）。

我的应用程序对性能相当关键，所以我尝试了另一种方法，即首先将所有稀疏矩阵连接到一个大的稀疏矩阵中，然后像这样只使用dot（）的一个调用：

import scipy.sparse as sp
SPMAT = sp.vstack(spmats, format='csr')
RES = SPMAT.dot(G)

这将产生一个长向量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

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

我运行在一台单核linux机器上，使用EPD7，内存为4GB。3.

编辑：

下面是一个小示例，它为我再现了这个问题：

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)

我得到：

1 loops, best of 3: 704 ms per loop
1 loops, best of 3: 1.34 s per loop

这种情况下的性能差异没有我自己的稀疏矩阵有一定的结构（可能是因为缓存）那么大，但串联矩阵更糟糕。

我已经尝试了与这两个10.1和12.0。

Answer 1

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

因为在我的特殊情况下，我计算float32稀疏矩阵和complex64密集向量的乘积，我可以分别乘以实部和虚部。这为我提供了4倍的加速。

这需要2.35s（使用<代码>SPMAT）。形状==（163840002097152）：

RES = SPMAT.dot(G)

虽然这只需要541ms：

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

结果是一样的。我想可能不需要预先分配n.zero，但我不知道如何才能做到这一点。