Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

In-place Cholesky (LDLt) solve #51

Open
mfalt opened this issue Feb 20, 2018 · 2 comments
Open

In-place Cholesky (LDLt) solve #51

mfalt opened this issue Feb 20, 2018 · 2 comments

Comments

@mfalt
Copy link
Collaborator

mfalt commented Feb 20, 2018

I implemented a wapper to do in-place solves for sparse Cholesky factorizations (i.e A_ldiv_B!) in https://github.com/mfalt/CholmodSolve2.jl (which is now registered in METADATA)
This could be be used in LeastSquaresDirect here: https://github.com/kul-forbes/ProximalOperators.jl/blob/07df4b3b391d5dac71cfde61470e8919f090baec/src/functions/leastSquaresDirect.jl#L103
and probably here: https://github.com/kul-forbes/ProximalOperators.jl/blob/07df4b3b391d5dac71cfde61470e8919f090baec/src/functions/leastSquaresDirect.jl#L107
to avoid allocations completely (well O(1) allocations even on multiple calls).
The same is true for QuadraticDirect. Are there any other places in the code where we use sparse Cholesky?
EDIT: I guess IndGraphSparse

It seems like epicompose.jl, indGraphFat.jl and indGraphSkinny.jl only has dense implementations, did we have a reason for this? could also have a sparse implementation.

@lostella
Copy link
Member

You mention that in CholmodSolve2 there is no support for complex right-hand side vectors: is that a limitation of SuiteSparse?

@mfalt
Copy link
Collaborator Author

mfalt commented Mar 12, 2018

SuiteSparse solve2 should support Float64, Complex{Float32} and "zomplex" Complex{Float64} but I didn't need them so I only implemented CholmodSolve2.jl for Float64. It should be a relatively trivial to update the code to support these too, but I would need to properly test that all vectors are allocated to proper size.

It seems like julia only wraps Float64 and Complex{Float64} in general. To update the code for this, it should be sufficient to: add where Tv <: VTypes on most of the functions, potentially add a new set of buffers, and verify that the lengths are double where appropriate.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Projects
None yet
Development

No branches or pull requests

2 participants