犰狳 C++ LU 分解

Posted 2023-02-17

【中文标题】犰狳 C++ LU 分解

【英文标题】：Armadillo C++ LU decomposition

【发布时间】：2013-07-30 07:23:56

【问题描述】：

我正在使用 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 向量的几个线性系统。

你能帮帮我吗？

【问题讨论】：

也许您可以告诉我们您已经尝试过什么？也许您可以编辑您的问题以包含您尝试的SSCCE？

感谢 Joachim，我添加了我尝试过的内容...

【参考方案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 分解。

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