What is the most efficient way to fill a sparse matrix? I know that sparse matrixes are CSC, so I expected it to be fast to fill them column by column like
using SparseArrays M = 100 N = 1000 sparr = spzeros(Float64, M, N) for n = 1:N # do some math here idx = <<boolean index array of nonzero elements in the nth column>> vals = <<values of nonzero elements in the nth column>> sparr[idx, n] = vals end However, I find that this scales very poorly with N. Is there a better way to fill the array? Or perhaps, I should not bother with filling the array and instead initialize the matrix differently?
11 Answer
You can do sparse(I, J, V, M, N) directly:
julia> using SparseArrays julia> M = 100; julia> N = 1000; julia> nz = 2000; # number of nonzeros julia> I = rand(1:M, nz); # dummy I indices julia> J = rand(1:N, nz); # dummy J indices julia> V = randn(nz); # dummy matrix values julia> sparse(I, J, V, M, N) 100×1000 SparseMatrixCSC{Float64, Int64} with 1982 stored entries: ⣻⣿⣿⣿⣿⡿⣾⣿⣿⣿⣿⣿⣿⣷⣾⣽⣿⢿⢿⣿⣿⣿⢿⣿⣾⣿⣽⣿⣿⣾⣿⣿⣿⣿⣿⣿⣿⣿⣻⣿ ⣼⣿⣿⡿⣿⣿⡽⣿⣿⣿⣿⣿⣿⣿⣿⣿⣿⣿⣿⣻⣿⡿⣿⣿⣿⡿⣿⡿⣯⢿⣿⠾⣿⣿⡿⢿⣿⣻⡿⣾ which should scale decently with size. For more expert use, you could directly construct the SparseMatrixCSC object.
EDIT:
Note that if you have to stick with the pseudo code you gave with a for-loop for column indices and values, you could simply concatenate them and create I, J, and V:
I = Int[] J = Int[] V = Float64[] for n = 1:N # do some math here idx = <<boolean index array of nonzero elements in the nth column>> vals = <<values of nonzero elements in the nth column>> I = [I; idx] J = [J; fill(n, length(I))] V = [V; vals] end but that'd be slower I think.