犰狳 C++ LU 分解

Question

我正在使用 Armadillo C++ 库来求解中/大尺寸的线性系统（1000-5000 个方程）。

因为我必须解决不同的线性系统

AX=b

其中 A 始终相同而 B 发生变化，我只想对 A 进行一次 LU 因式分解，然后用不同的 b 重复使用 LU 分解。不幸的是，我不知道如何在犰狳中执行这种操作。

我所做的只是 A 矩阵的 LU 分解：

arma::mat A;
// ... fill the A matrix ...
arma::mat P,L,U;
arma::lu(L, U, P, A);

但现在我想使用矩阵 P、L 和 U 来求解具有不同 b 向量的几个线性系统。

你能帮帮我吗？

Answer 1

由于A = P.t()*L*U（由于舍入误差，相等仅是近似值），求解P.t()*L*U*x = b 中的x 需要置换B 的行并执行正向和反向替换：

x = solve(trimatu(U), solve(trimatl(L), P*b) );

由于犰狳中缺乏真正的三角求解器，以及执行行置换的快速方法，因此相对于直接调用相关的计算 LAPACK 子例程而言，此过程效率不高。

一般建议是避免在更高级别的库（如犰狳）中进行显式 LU 分解。

1. 如果同时知道所有不同的b，则将它们作为列存储在矩形矩阵B 和X = solve(A,B); 中
2. 如果不同的b 一次是已知的，那么如果不同r.h.s 的数量不同，预计算AINV = A.i(); 和x = AINV*b; 会更有效。向量足够大。看到这个answer 到一个类似的question。