1 typedef Matrix
<double, 5, 3> Matrix5x3
;
2 typedef Matrix
<double, 5, 5> Matrix5x5
;
3 Matrix5x3 m
= Matrix5x3::Random();
4 cout
<< "Here is the matrix m:" << endl
<< m
<< endl
;
5 Eigen::FullPivLU
<Matrix5x3
> lu(m
);
6 cout
<< "Here is, up to permutations, its LU decomposition matrix:"
7 << endl
<< lu
.matrixLU() << endl
;
8 cout
<< "Here is the L part:" << endl
;
9 Matrix5x5 l
= Matrix5x5::Identity();
10 l
.block
<5,3>(0,0).triangularView
<StrictlyLower
>() = lu
.matrixLU();
12 cout
<< "Here is the U part:" << endl
;
13 Matrix5x3 u
= lu
.matrixLU().triangularView
<Upper
>();
15 cout
<< "Let us now reconstruct the original matrix m:" << endl
;
16 cout
<< lu
.permutationP().inverse() * l
* u
* lu
.permutationQ().inverse() << endl
;