Eigen C++求解线性系统时如何重用系数矩阵分解结果

计算科学 线性求解器 C++ 效率 本征
2021-11-26 22:29:59

我的问题需要解决密集的对称线性系统,例如:

A x = b, A y = x, A z = y+x,... 按顺序。

在 Eigen C++ 中,如果我利用以下对称性A

x=A.ldlt().solve(b);
y=A.ldlt().solve(x);
z=A.ldlt().solve(y+x);

同一个矩阵A必须Cholesky多次分解。

householderQRLU算法可能存在类似的问题。

如何为其他系统重用系数矩阵分解结果?

这似乎相当于询问是否有一个简单的上三角线性系统的反向替换求解器。

2个回答

您需要像这样声明一个 LDLT 对象:

LDLT<MatrixXd> ldlt(A);
x = ldlt.solve(b);
y = ldlt.solve(x);
...

鉴于接受的答案有错别字并且没有给出完整的代码,我决定发布我的。

#include <iostream>
#include <Eigen/Dense>
using namespace std;
using namespace Eigen;
int main()
{
    MatrixXd A;
    A.resize(3,3);
    VectorXd b;
    b.resize(3);

    A <<
      13, 5, 7 ,
      5 , 9, 3 ,
      7 , 3, 11;
    b << 3, 3, 4;

    cout << "Here is the matrix A:\n" << A << endl;
    cout << "Here is the vector b:\n" << b << endl;
    LDLT<MatrixXd> L = A.ldlt();
    VectorXd x = L.solve(b);
    cout << "The solution is:\n" << x << endl;
    cout << "reconstructed b is:\n"  << A*x << endl; 

}