Optimal way to do in-place modification of sparse matrix?

General Usage
1

I’m looking for the most performant way to modify a large SparseMatrixCSC in-place. The matrix is initially built using the sp = sparse(I,J,K) method. Later, I need to modify only the non-zero entries in sp. I’m tempted to do

sp.nzval .= newK

where newK is the updated K vector. This is really fast, but I worry this approach is fragile (as e.g. there are no guarantees that zeros in the original K are dropped when building the original sp, or that the order of sp.nzval is the same as in K). Can this be made safe? If not, is there some recommended approach with comparable performance?

1 Like
2

There is dropzeros! to remove all the stored zeros.

3

Hi Kristoffer!
I know about dropzeros!, but my problem is the opposite. I want to keep some structural zeroes when building sp, as they may become non-zero when updating it. Moreover, to use sp.nzval .= newK safely, I need to be sure about the order of sp.nzval, so I know which entry in the matrix I’m modifying for each element of newK. I have the feeling, however, that this is just not the proper Julian way to do what I want, hence the question.

I see roughly three possible approaches to modify structural values of sp after the fact

1- In-place broadcast! directly over sp.nzval (risky!)

julia> @btime sp.nzval .= newK;
  5.358 μs (0 allocations: 0 bytes)

2- In-place setindex! over structural values

julia> updatesp!(sp, I, J, newK) = for (k,(i,j)) in enumerate(zip(I, J)); sp[i,j] = newK[k]; end
updatesp! (generic function with 1 method)

julia> @btime updatesp!(sp, I, J, newK)
  154.394 μs (0 allocations: 0 bytes)

3- Allocating a whole new sparse matrix

julia> @btime sparse(I, J, newK);
  334.870 μs (19 allocations: 547.78 KiB)

I like the speed of 1 :slight_smile:. The question is whether it is possible to do something that is “safe” with comparable performance.

EDIT: corrected typo in code.

4

In the docs for sparse (Julia 0.6) it says: “Numerical zeros in (I, J, V) are retained as structural nonzeros; to drop numerical zeros, use dropzeros!.”. So your sp.nzval .= newK is fine.

5

Ah, great. So there is a guarantee that structural zeros are not dropped. On the other hand I only now realised that the order of sp.nzval is also guaranteed to sweep columns of sp in order, so the sp.nzval .= newK is indeed safe. Thanks a lot!

1 Like