[eigen] Solving positive semidefinite systems

[ Thread Index | Date Index | More lists.tuxfamily.org/eigen Archives ]


Hi,

First of all, I'd like to thank you all for the wonderful library!

In trying to convert some particular old code to use Eigen, I've run into a case where Eigen's Cholesky solvers don't have a away to ignore A(k,k) if it's 0 (or sufficiently close to 0). Instead, they either silently give me bad results or give up altogether (depending on if it's the sparse or dense variants).

Golub and Van Loan's Matrix Computations, chapter 4.2.8 (The semdidefinite case) gives a very easy modification to the regular Cholesky decomposition that allows it to be used for positive semidefinite matrices. All it does is skip the update iteration if A(k,k) is found to be 0. In practice, the implementation that I'm used to checks if A(k,k) is less than some epsilon, and if it is, sets the current diagonal entry and the rest of the row (or column) to be exactly 0 and then skip the rest of the loop. Then during back substitution, it checks whether the diagonal value is 0, and if so, then it sets the corresponding component of the solution to 0 instead of dividing by it. This trick is similarly used in the LU matrix solvers as well with the exception that we just set the 0 diagonal entries to 0 during factorization and simply skip the divide by 0's during back substitution.

Can such a modification be considered for Eigen's dense/sparse Cholesky/LU solvers? I have managed to work around the problem by using other Eigen solvers like FullPivLU (for dense matrices) or SparseQR (for sparse matrices). However, they're much slower.

Thanks,
-Edward



Mail converted by MHonArc 2.6.19+ http://listengine.tuxfamily.org/