matlab - 在Matlab中求解稀疏三角线性系统

Question

我在Matlab中实现了一些大型稀疏矩阵的LU分解算法来求解线性系统。当我得到 L,U 矩阵时，我使用后向替换和前向替换算法来求解三角线性系统：

    %x = U\y;

for i = n : -1 : 1
    x(i,:) = (y(i,:)-U(i,:)*x)/U(i,i);
end

但我发现这段代码是瓶颈。虽然我可以使用 A\b 来获得解决方案，但是我想知道如何在 Matlab 中实现一个高效的算法来解决这个问题，例如，我可以编写矩阵乘积来模拟以下动作而无需 for 循环？

（我有一些参考书和论文，但所有代码都不是在 Matlab 中，仅用于 C++ 或 C 代码）

score 3 · Accepted Answer

首先：正确性先于速度；您发布的循环产生的结果与不同，U\y因此您可能需要先检查一下 :)

AFAIK，反斜杠对输入矩阵进行一些检查，并相应地调用最快的算法。当这些检查表明A是下三角形时，它将完全按照您所做的（但可能更有效）。

无论如何，为了加速你的代码：你应该 pre-allocate x，否则 Matlab 被迫在每次迭代时增长向量。另外，调用循环变量ii -- i是虚数单位，每次迭代的名称解析需要一些时间。所以，总结一下：

x = zeros(size(y));
for ii = n : -1 : 1
    x(ii,:) = (y(ii,:)-U(ii,:)*x)/U(ii,ii);
end

请注意，没有“矢量化”解决方案，因为下一个结果取决于前一个结果。

1 回答 1