This bugzilla service is closed. All entries have been migrated to https://gitlab.com/libeigen/eigen

Bug 1268

Summary: Add LDLT with block diagonal
Product: Eigen Reporter: w.h.greene
Component: CholeskyAssignee: Nobody <eigen.nobody>
Status: CONFIRMED ---    
Severity: Feature Request CC: allanmulin, chtz, gael.guennebaud, w.h.greene
Priority: Normal    
Version: unspecified   
Hardware: All   
OS: All   
Whiteboard:
Attachments:
Description Flags
simple test showing failure
none
Detect faillure in LDLT none

Description w.h.greene 2016-07-30 16:00:16 UTC
Created attachment 719 [details]
simple test showing failure

I have a simple symmetric matrix that I believe LDLT should solve. 
Instead it produces an incorrect solution. As the example code shows,
PartialPivLU produces a correct solution.

The incorrect solution is not reported; as can be seen info()
reports success.

Based on debugging, I can see that realAkk is flagged as zero during
the decomposition but this should not be the case. Note that
LU reports full-rank for this matrix.

In cases where a singularity is *really* encountered, I believe that
info() should report this rather than failing silently.
Comment 1 Christoph Hertzberg 2016-07-30 18:12:34 UTC
Ok, first of all, your matrix is indefinite -- you see this if you make an eigenvalue decomposition. Unfortunately, our ldlt decomposition often fails in that case.
Also realAkk==0 can happen for full-rank (but indefinite) matrices.
The easiest example is [0 1; 1 0]

I agree however, that silently returning Success in that case is not a good idea.
Perhaps, if realAkk==0 happens once, but we get abs(realAkk)>0 afterwards, declare the decomposition as failed (I have not thought through all cases ...)

Ideally, we'd support indefinite matrices as well, of course (would require allowing 2x2 blocks on the diagonal, i.e. an API change, though).
Comment 2 Gael Guennebaud 2016-08-23 17:24:10 UTC
This is rather a request for another LDLT variant.
Comment 3 Christoph Hertzberg 2016-08-23 17:29:09 UTC
I agree, but we definitely must not return Success, if the decomposition actually failed.
Comment 4 Gael Guennebaud 2016-08-23 20:46:28 UTC
Created attachment 722 [details]
Detect faillure in LDLT

right, I've been too fast sorry.

Here is an initial attempt. With this patch, LDLT returns NumericalIssue for a zero matrix, though the factorization is exact. To fix this we would need to explicitly check for a zero input in the case of a zero diagonal.
Comment 5 Gael Guennebaud 2016-08-23 21:17:41 UTC
ok, finally I've implemented a special path for a zero diagonal, and pushed the changes:

https://bitbucket.org/eigen/eigen/commits/c6f5c3c6ef20/
Summary:     Bug 1268: detect failure in LDLT and report them through info()
Comment 6 Allan Leal 2016-12-05 13:06:28 UTC
Out of curiosity, are there any plans on implementing LDLT solvers for indefinite matrices (i.e. when D has 1x1 and 2x2 blocks depending on the signs of the eigenvalues)?
Comment 7 Nobody 2019-12-04 16:04:20 UTC
-- GitLab Migration Automatic Message --

This bug has been migrated to gitlab.com's GitLab instance and has been closed from further activity.

You can subscribe and participate further through the new bug through this link to our GitLab instance: https://gitlab.com/libeigen/eigen/issues/1268.