python中稀疏矩阵的矩阵幂

JRu*_*Run 6 python numpy linear-algebra scipy sparse-matrix

我试图找到一种方法来对稀疏矩阵 M 进行矩阵幂: M^k = M*...*M k 次,其中 * 是矩阵乘法(numpy.dot),而不是逐元素乘法

我知道如何对普通矩阵执行此操作:

import numpy as np
import scipy as sp
N=100
k=3
M=(sp.sparse.spdiags(np.ones(N), 0, N, N)-sp.sparse.spdiags(np.ones(N), 2, N, N)).toarray()
np.matrix_power(M,k)
Run Code Online (Sandbox Code Playgroud)

对于稀疏 M 我该如何做:

M=(sp.sparse.spdiags(np.ones(N), 0, N, N)-sp.sparse.spdiags(np.ones(N), 2, N, N))
Run Code Online (Sandbox Code Playgroud)

当然,我可以通过递归乘法来做到这一点,但我想知道 scipy 中是否有类似用于稀疏矩阵的矩阵幂之类的功能。非常感谢任何帮助。提前致谢。

hpa*_*ulj 6

**已实施csr_matrix. 有一个__pow__方法。

在处理了一些特殊情况后,这__pow__会:

            tmp = self.__pow__(other//2)
            if (other % 2):
                return self * tmp * tmp
            else:
                return tmp * tmp
Run Code Online (Sandbox Code Playgroud)

对于稀疏矩阵,*是矩阵乘积(dot对于 ndarray)。所以它正在做递归乘法。


如前所述mathnp.matrix也将**( __pow__) 实现为矩阵幂。事实上它最终会调用np.linalg.matrix_power.

np.linalg.matrix_power(M, n)是用Python编写的,所以你可以很容易地看到它做了什么。

因为n<=3只是重复dot.

对于较大的n,它会进行二元分解以减少 s 的总数dot。我认为这意味着n=4

result = np.dot(M,M)
result = np.dot(result,result)
Run Code Online (Sandbox Code Playgroud)

稀疏版本并不那么通用。它只能处理正整数幂。

您不能指望numpy在备用矩阵上运行的函数。真正起作用的是将操作传递给数组自己的方法。例如np.sum(A)打电话A.sum()