# How to re-use the coefficient matrix decomposition result when solving linear systems by Eigen C++

My problem needs to solve dense symmetric linear systems something like:

A x = b, A y = x, A z = y+x,...in sequence.

In Eigen C++, if I take advantage of the symmetry of A by using:

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


The same matrix A has to be Cholesky factorized many times.

Similar issue may exists for householderQR and LU algorithms.

How can I reuse the coefficient matrix decomposition result for other systems?

It seems this is equivalent to ask whether there is a simple back-substitution solver for upper triangular linear systems.

You need to declare a LDLT object like this:

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

• Thank you. So if I have A changed in a loop, then just re-declare it? – LCFactorization Dec 24 '13 at 15:51
• There is a compute method like this: ldlt.compute(B); – ggael Dec 25 '13 at 10:00

Given that the accepted answer has typos and does not give full code, I decided to post mine.

#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;

}