Sparse Cholesky of Gram Matrix (SuiteSparse?)

[stevengj]

@guille

Yes, sure, it is possible to do it with two factorizations; one which initially solves the underdetermined system (e.g., via a Moore-Penrose Pseudoinverse or other least-norm type solution) and one which solves the resulting least squares problem, but just computing the two factorizations would override any performance gains from simply factorizing the Gram matrix and reusing this factorization

As I noted in my comment, you don’t need two factorizations, you only need one. The QR factorization of A (which you use to solve the least-squares problem) automatically gives you the LQ factorization of A^T (which you can use to solve the underdetermined problem).

Unfortunately, qr(A)' \ c is not implemented in Julia, and I was too lazy to dig through the QR factorization object F = qr(A) to figure out how to use its transpose for the underdetermined problem. (There are some permutation vectors inside F that you need to decipher.) In principle this shouldn’t be too much work, however, which is why I filed missing qr(A)' \ b for sparse · Issue #115 · JuliaSparse/SparseArrays.jl · GitHub.