Cholesky for complex *symmetric*

[phil.cummins@]

Hi,

I am working with a Galerkin FEM implementation of an elastodynamic problem in the frequency domain. For the purely elastic case, this results in a symmetric, positive definite linear system that is efficiently solved by Cholesky decomposition. In order to consider anelasticity, however, I make the elastic moduli complex, and in order to avoid "wrap-around" problems in the time domain, the frequency is also often made complex. This results in a complex symmetric, not Hermitian symmetric, linear system. My question is, what is an efficient way to solve this system? (decomposition is normally preferred, since I need to solve for many r.h.s.s). It seems to me that the Cholesky algorithm as formulated for real, symmetric matrices should still work. But all the canned Cholesky routines I'v seen for complex matrices require Hermitian symmetry. Is there any reason why Cholesky can't be used for complex symmetric, as opposed to Hermitian, matrices?

Any suggestions appreciated.

 

[hotvette]

Looks like you may need LDL^T which is 4x the cost of Cholesky.

http://www.nag.co.uk/numeric/fl/manual/pdf/F07/f07nsf.pdf

 

[phil.cummins@]

Thanks, that sounds right. I am now looking at some of the routines developed by Tim Davis et al. (e.g. CHOLMOD, http://www.cise.ufl.edu/research/sparse/cholmod/ )

 

[AlephZero]

Cholesky only works on positive definite Hermitian matrices. Real symmetric matrices are Hermitian, of course.

Freqency domain problems with any form of damping produce complex symmetric matrices which are not hermitian. Don't waste your time trying to find a Cholesky solver for them, there isn't one.

Often, if you are working in the frequency domain, you can reduce the model size by orders of magnitude by working in modal coordinates. The eignensolution for the undamped modes will be a hermitian matrix problem. If your model has more than say 1000 modes in the frequency range of interest, you are probably doing something not very clever (fewer than 100 modes would be more typical), but solving even a 1000x1000 dense matrix with no special properties is not a big deal compared with an unreduced FE model with maybe 10^6 degrees of freedom.