Hi,

I am trying to perform a LU decomposition of a square matrix using Eigen::FullPivLU so the original matrix m = L*U.

However, when I do that, the matrices L and U calculated, when multiplied, do not result in the original matrix m.

The code is below and matrix m (just an example) is m = [88, -57, -63, 0;
-57, 39, 0, -36;
-63, 0, 45, -46;
0, -36, -46, 1]

Code: Select all


Eigen::FullPivLU<MatrixXf> lu(tam,tam);
lu = lu.compute(m);
std::cout << "Here is, up to permutations, its LU decomposition matrix:"
<< std::endl << lu.matrixLU() << std::endl;

MatrixXf u = lu.matrixLU().triangularView<Eigen::Upper>();
std::cout << "Here is the U part:" << std::endl;
std::cout << u << std::endl;

MatrixXf l;
l.setIdentity();
l.triangularView<Eigen::StrictlyLower>() = lu.matrixLU();
std::cout << "Here is the L part:" << std::endl;
std::cout << l << std::endl;



When I multiply L*U, it does not result in m.

Anyone knows why does it happen?

As explained there: https://eigen.tuxfamily.org/dox/classEi ... PivLU.html, it compute M = P * L * U where P is a row-interchange permutation matrix.