Re: [eigen] [patch] LDLt decomposition with rank-deficient matrices

[Benoit Jacob <jacob.benoit.1@xxxxxxxxx>]

2010/2/24 Gael Guennebaud <gael.guennebaud@xxxxxxxxx>:
> Hi,
>
> note that the LDLT decomp does pivoting, so it is normal that you don't get
> the same diagonal matrix. If you reconstruct the matrix from the decomp you
> will see that in this case the decomp is correct.
>
> In the devel branch I've just added a reconstructedMatrix() to LDLT for that
> precise purpose, so you can do:
>
> cout << (A - A_factored.reconstructedMatrix()) << endl;

Good idea, but for homogeneity, we should then do that also for all
other decompositions. Do you agree?

Benoit

> As a side note,
>
> A.diagonal().setOnes();
>
> works as expected ;) (simpler than a home made for loop)
>
> gael.
>
>
> On Wed, Feb 24, 2010 at 3:29 AM, Ben Goodrich <bgokgm@xxxxxxxxxxxxxx> wrote:
>>
>> Hi,
>>
>> I have found another robustness problem with ldlt() in the current
>> development branch. Some combination of rank-deficiency, constant
>> diagonal, pivoting, and double precision can make it impossible to
>> roundtrip a matrix through the LDLt decomposition. I am attaching
>> another test case that yields the output below.
>>
>> Thank you,
>> Ben
>>
>> A = LDL' has a unit diagonal (apart from numerical noise)
>> 1
>> 1
>> 1
>> 1
>> 1
>> This is a numerical disaster (compare the second cells).
>> Here is the true D
>> 1
>> 0.915092
>> 0
>> 0
>> 0
>> Here is the calculated D
>> 1
>> 0.0785762
>> 2.01553e-16
>> 7.97973e-17
>> 1.39591e-15
>> This is good (because the diagonal of A is forced to have 1.0 in all its
>> cells).
>> Here is the true D
>> 1
>> 0.915092
>> 0
>> 0
>> 0
>> Here is the calculated D
>> 1
>> 0.915092
>> -6.06069e-16
>> -8.0296e-16
>> -8.4893e-16
>
>