小编dwf*_*wfm的帖子

我试图弄清楚如何有效地解决稀疏的三角形系统,在scipy稀疏中Au*x = b.

例如,我们可以构造一个稀疏的上三角矩阵,Au和一个右手边b,其中:

import scipy.sparse as sp
import scipy.sparse.linalg as sla
import numpy as np

n = 2000
A = sp.rand(n, n, density=0.4) + sp.eye(n)
Au = sp.triu(A).tocsr()
b = np.random.normal(size=(n))

我们可以使用spsolve来解决问题,但很明显三角形结构没有被利用.这可以通过对解决方案进行计时并将其与splu中的求解方法进行比较来证明.(这里使用iPython的%时间魔法)

%time x1 = sla.spsolve(Au,b)
CPU times: user 3.63 s, sys: 79.1 ms, total: 3.71 s
Wall time: 1.1 s

%time Au_lu = sla.splu(Au)
CPU times: user 3.61 s, sys: 62.2 ms, total: 3.67 s
Wall time: 1.08 s

%time x2 = Au_lu.solve(b)
CPU times: user 25 ms, …