# Sparse Array Compression and Decompression upgrade to v0.7+

Summary: I can “compress” and “decompress” a sparse array using deprecated functions. I want to upgrade the code to v0.7+. (I want to use this compression/decompression to numerically evaluate the Jacobian of a function of a vector of size ~106.)

1. I have a function `f(x)`, and a boolean sparse matrix `AM` that represents the Jacobian’s sparsity pattern. `AM` is a sort of an adjacency matrix (thus the name `AM`) of dependencies between different indices of `x`. `AM` has non-zero sub indices `i` and `j`, and non-zero linear indices `ij`:

``````i, j, _ = findnz(AM)
ij = find(AM)
``````
2. I also have a home-made vector of colors `col_idx` which groups columns of the matrix `AM` by color. (The vector `col_idx` has integer values.) Two columns (say, `j1` and `j2`) that share the same color (i.e., `col_idx[j1] == col_idx[j2]`) do not have a common non-zero. (I think the exact term is structural orthogonality.) This allows me to construct the corresponding indices of my compressed objects:

``````icol = i
jcol = col_idx[j]
ijcol = sub2ind(icol, jcol)
``````
3. Once these 6 linear and sub indices (`i`, `j`, `ij`, `icol`, `jcol`, and `ijcol`) are created, I can easily

• compress using

``````compress(Jac, icol, jcol, ij) = sparse(icol, jcol, Jac[ij])
``````
• and decompress using

``````decompress(smallJac, i, j, ijcol) = sparse(i, j, smallJac[ijcol])
``````

Now, with Julia v0.7+, the functions `findnz`, `find`, and `sub2ind` are deprecated, so I am trying to upgrade this code to v0.7+. Deprecation of those function suggest using the `CartesianIndex` and `LinearIndex` types, but those are not straightforward to use with sparse matrices. I have been scratching my head for almost a week on related topics now, so any help is much appreciated!

The `findnz` function is not deprecated but part of the `SparseArrays` stdlib, so I think this becomes:

``````AM = sparse([1,2,3],[1,2,3],ones(3))
i, j, _ = findnz(AM)
LI = LinearIndices(AM)
ij = LI[i,j]
icol = i
jcol = col_idx[j]
ijcol = LI[icol,jcol]
``````
1 Like

is not the correct answer because `LI[i,j]` is a 3x3 array while I want a 3x1 array like `ind2sub(i,j)`.

``````ij = LI[CartesianIndex.(i, j)]
``````

works. This `ij` can then also be used for getting at the values of a sparse matrix `S` via `S[ij]`.

Thanks for your reply though as it helped be find the solution! I’ll report back if this works all the way through More details: I think the issue stemmed from the suggestion in `ind2sub(i,j)`'s deprecation warning. `ind2sub(i,j)` may be wrongly interpreted because the suggestion:

``````julia> ij = sub2ind(i,j)
┌ Warning: `sub2ind(A::AbstractArray, I...)` is deprecated, use `(LinearIndices(A))[I...]` instead.
│   caller = top-level scope at none:0
└ @ Core none:0
``````

is what you suggested (i.e., `ij = LI[i,j]`) and does not correspond.
Instead, it should have been the one that comes when using the size of the matrix in the arguments:

``````julia> ij = sub2ind(size(AM), i, j)
┌ Warning: `sub2ind(inds::Union{DimsInteger, Indices}, I1::AbstractVector{T}, I::AbstractVector{T}...) where T <: Integer` is deprecated, use `(LinearIndices(inds))[CartesianIndex.(I1, I...)]` instead.
│   caller = top-level scope at none:0
└ @ Core none:0
``````
2 Likes

So after a bit of work, I finally found a decent solution, posting if anyone is interested… Although the solution above with `ij = LI[CartesianIndex.(i, j)]` worked, I thought I should note that once the compressed sparse matrix `smallJac` is built, accessing it via `smallJac[ijcol]` can be quite slow for large sizes. This is because CSC sparse matrices are stored by columns, so that the order of indices are perturbed by compression and decompression. It is thus faster (better?) to figure out the permutation of indices that occurs during compression and use it, via

``````_, _, perm = sparse(icol, jcol, 1:length(icol))
iperm, jperm = i[perm], j[perm]
decompress(smallJac, iperm, jperm) = sparse(iperm, jperm, smallJac.nzval)
``````

This is faster because accessing the non-zero values of `smallJac` via its `.nzval` field is fastest, right? (`Jac.nzval` can also be used for the compression part, which might be faster and avoids linear indices altogether.)

I don’t know if it’s applicable for your use case, but maybe you could also write an `AbstractArray` implementation that wraps your compressed matrix, eliminating the need for compression/decompression altogether and doing the index translation on the fly?

That would be nice but I would need some help! Could you give me some hints on how to do that?

Basically you just need to define a new type that holds your compressed array, make it inherit from `AbstractArray` and then implement the 4 functions described here:

https://docs.julialang.org/en/stable/manual/interfaces/#man-interface-array-1

1 Like