Re: [eigen] Complex symmetric Cholesky factorization

[John Daily <John.Daily@xxxxxxxxxxxxxxxxx>]

Thanks Gael. 

On 16 August 2015 at 07:00, Gael Guennebaud <gael.guennebaud@xxxxxxxxx> wrote:

On Thu, Aug 13, 2015 at 9:28 PM, John Daily <John.Daily@xxxxxxxxxxxxxxxxx> wrote:

Greetings, certain symmetric complex (non-Hermitian) matrices admit an LL^T/LDL^T factorization[1] (rather than an LL^H LDL^H). I am working with sparse matrices and the code changes to support this are minimal and are all encompassed in Simplicial_Cholesky_impl.h, in the factorize method, as well as the updated solve calling transpose rather than selfadjointView. I have two questions:

Yes we are aware of that and this has been a planed feature for a long time.

 

2) If so what would be an acceptable way to do this, that would not cause surprise to users? It would seem that a separate template parameter or even function call would make this explicit, such as instead of .factorize() calling .factorizeComplexSymmetric() or something. 

We already have a SelfAdjoint and Symmetric enums to distinguish between the two cases, so perhaps an idea would be to allow to pass Symmetric to the second template parameter, like:

[Simplicial]L[D]LT<MatrixType, {{Lower,Upper,Lower|Upper}[|Symmetric] , Symmetric}>

I don't quite follow your notation here for the Symmetric enum, did you mean one of those to be SelfAdjoint as written? 

 

For consistency I'm including the dense factorizations here.

Sure, we should do both while we are at it and have them consistent. 

 

Actually, this solution is already the one we adopted for PardisoLLT<> and PardisoLDLT<>.

I'll take a look at that. I'll try and make a branch achieving this and will let you know when I have something that can be reviewed. 

 

In a second time, for consistency, we should also add a symmetricView<>() method to [Sparse]MatrixBase<> returning a SelfAdjointView<...,{Upper,Lower}[|Symmetric]>. Likewise, we could also think about adding templated variants of adjoint() / selfadjointView() to ease writing generic code with true compile-time branches...

I don't have enough familiarity with the sparse internals to know how difficult this part will be. I'll start with the above. 

Thanks,

 

cheers,

gael