python - 求解大量小型线性系统

Question

我需要用 Python 求解大量 3x3 对称正定系统。到目前为止，我做到了

res = numpy.zeros(n)
for k, obj in enumerate(data_array):
    # construct A, rhs, idx from obj
    res[idx] += numpy.linalg.solve(A, rhs)

这会产生正确的结果，但是如果 n 很大，它也会很慢。 (嗯...是的。)也许 3x3 不是问题大小，调用 solve() 很有意义。

有什么提示吗？

Answer 1

在 NumPy 1.8 及更高版本中，numpy.linalg.solve actually broadcasts .对于 numpy.linalg.solve(a, b)，如果 b.ndim == a.ndim - 1，它将执行广播矩阵向量求解；否则，它将进行广播矩阵-矩阵求解。 (没有记录此决策标准；我必须查看来源。)

如果您可以有效地构造一个 A 和 rhs 的堆栈，您可以调用一次 solve 并避免 Python 循环。