Sparse to dense matrix

For some numerical optimization I need to get the inverse of a sparse matrix. How do I compute it efficiently? Julia throws me an error when I try to use either inv() or A\B and tells me to convert my sparse matrix first to a dense one. How can I do that?

2 Likes

full or Matrix can solve your problem.

julia> a = speye(5,5)
5×5 SparseMatrixCSC{Float64,Int64} with 5 stored entries:
  [1, 1]  =  1.0
  [2, 2]  =  1.0
  [3, 3]  =  1.0
  [4, 4]  =  1.0
  [5, 5]  =  1.0

julia> b = full(a)
5×5 Array{Float64,2}:
 1.0  0.0  0.0  0.0  0.0
 0.0  1.0  0.0  0.0  0.0
 0.0  0.0  1.0  0.0  0.0
 0.0  0.0  0.0  1.0  0.0
 0.0  0.0  0.0  0.0  1.0

julia> c = Matrix(a)
5×5 Array{Float64,2}:
 1.0  0.0  0.0  0.0  0.0
 0.0  1.0  0.0  0.0  0.0
 0.0  0.0  1.0  0.0  0.0
 0.0  0.0  0.0  1.0  0.0
 0.0  0.0  0.0  0.0  1.0
3 Likes

Also take a look at
https://github.com/JuliaMath/IterativeSolvers.jl

2 Likes

Thanks a lot, I’ll have a look! :slight_smile:

Are you sure you do need the inverse? Textbook equations often use an inverse where the numerical implementation does without one. Examples include linear regression (b=(X'X)^{-1}X'y, but you use a factorization and never form the inverse), BFGS update (work with the inverse of the approximate Hessian instead), etc.

(in case you have given this question due consideration and I am preaching to the choir, please ignore)

1 Like

You want to use an iterative solver or a sparse factorization on sparse matrices. I’m surprised \ doesn’t automatically fallback to some polyalgorithm over SuiteSparse for that?

:+1: for factorizations. For efficiency and stability purposes, you should also have a look at the documentation of some factorization routines available in julia by issuing ?lufact, ?qrfact and ?schurfact, the first two of which already have methods defined for sparse matrices. For the rest, you should check an iterative solver suitable for your purpose as mentioned by Chris and Sebastian above.

Thank you guys for the many answers. Indeed, if you use A\B, A might be sparse. My problem was instead that A is actually singular and thus does not allow me to calculate the inverse. I am looking for the mistake on my side in calculating A, but if anyone of you have experience using CompEcon.jl, please let me know!

It might be better to open a new thread specific to your issue, which preferably also includes a minimal working example. In this way, you would be highlighting the problem better (in a more concrete way), and people here would be helping you more easily. Otherwise, we tend to guess the issue ourselves :slight_smile:

3 Likes

You can instead use svdfact! for the pseudoinverse of the matrix if that’s necessary. Of course, first check analytically if it’s supposed to be singular before going down this route.

1 Like

That’s not available for sparse matrices in v0.6 (and even the v0.7 I have now), but qrfact is available on both and is equally insensitive to A being singular.

1 Like

full was not defined in Julia 1.0.2

1 Like

Did you try the other thing that was suggested?

1 Like

There is no full in Julia 1.0. Use Matrix(...) instead to do the conversion.

5 Likes

Array(a) or Matrix(a) do the trick.

1 Like

Hi all, reviving this old thread to ask if there is a generic way of converting sparse arrays to dense, in particular when using sparse gpu arrays (e.g. CuSparseMatrixCSRCuMatrix). Thanks!