sparse arrays with algebraic, non-numerical data

[abraunst]

using SparseArrays,StaticArrays

julia> T=SMatrix{2,2,Float64};

julia> s=sparse(1:10,1:10,fill(zero(T), 10));

julia> nnz(s)
10

julia> nnz(sparse(Matrix(s)))
100

This is just an example, there are other methods that don't work quite as they should (e.g. setindex!) The problem seems to be that the code for sparse matrices is riddled with != 0 checks instead of !iszero checks. This seems related to JuliaLang/julia#19561, but in that case the problem was for types without zero(T) for which sparse does not make real sense anyway. Maybe it would be reasonable to enforce the value type <: Number for SparseMatrixCSC? Are there people using it for non-number types? In that case one could fix all the != 0 instances and see if things work...

[ViralBShah]

Happy to accept a PR and see this work.

[abraunst]

Tried, with mixed results ;-) Doing so will forbid creating sparse vectors or matrices of types T with no defined zero(T) for instance in SparseVector{Tv,Ti}(s::AbstractVector{Tv}) as it needs to determine if the values are non-zero. That would be fine I think, except that in line 1246 of stdlib/SparseArrays/test/sparse.jl there is this test

   A = sparse(["a", "b"])
   @test_throws MethodError findmin(A, dims=1)

That now fail in the creation of the vector. I don't know what to think, is SparseVector{String,Ti} something that should be supported? All remaining tests pass (both for vector and matrices).
All in all I think that the best would be to be able to specify the non-stored value for sparse arrays.

[ViralBShah]

I would be ok with forbidding creation of sparse matrices where zero(T) is not defined. Would be good to hear other thoughts.

[stevengj]

See also JuliaLang/julia#19561, JuliaLang/julia#11408.

[abraunst]

See also JuliaLang/julia#19561, JuliaLang/julia#11408.

Thanks @stevengj. I was unaware of JuliaLang/julia#11408. If there is no consensus regarding the drop of support for types T without zero(T), then maybe it is better close this (and the PR)?

[andreasnoack]

I'm in favor of restricting SparseMatrixCSC to elements with a zero but regardless of that discussion, I think JuliaLang/julia#30580 improves that the situation for StaticArrays, right?

[abraunst]

I'm in favor of restricting SparseMatrixCSC to elements with a zero but regardless of that discussion, I think JuliaLang/julia#30580 improves that the situation for StaticArrays, right?

Yes, but I admit that I found the example after I ran into the "bug" in the code. The current code uses != 0 in several places, which is of course wrong for non-numerical types that define zero; but at least it doesn't give an error for non-numerical types that do not define zero (although !=0 is always true there, so for instance the sparse-dense-sparse round-trip gives a matrix full of stored zeros). Using iszero would mean to really drop support for those types (e.g. String) even for simple storage (which I don't know if people really use, but at least it appears in some tests).

[yurivish]

D4M sometimes uses sparse arrays to store non-numeric data.

I found this issue after being tripped up by the fact that broadcast non-numeric data didn't work, and was a little surprised that you couldn't broadcast an Any[1, 2] with a sparse matrix either, though I see now why that's the case. Regardless of the conclusions here, it would be nice if that error message were improved to be a bit more specific. Thanks to everyone for taking the time to thoroughly map out the territory — I learned a lot reading through some of the issues just now.

I wonder whether Base should maybe define zero(::String) == "", which would at least make that specific case work.

[abraunst]

I wonder whether Base should maybe define zero(::String) == "", which would at least make that specific case work.

I get your point, but at least from the mathematical point of view it seems an awkward stretch IMHO: "" is the identity element of the multiplicative (non-abelian) semigroup (*,String), whereas zero(T) should be the identity of the additive group (+,T). I think all this awkwardness could be avoided if one could specify the "non-stored" element of sparse arrays. On numerical data, linear operations could be fast-tracked if the non-stored element iszero so as to retain the current performance (and a sensible default could mimic current or desired behavior). I'm not saying all this is easy to implement though of course, but it seems like a fun project.

[StefanKarpinski]

I wonder whether Base should maybe define zero(::String) == "", which would at least make that specific case work.

This should absolutely not be defined. The empty string is the multiplicative unit of the string semigroup which is why we have one(String) == "". There is no additive identity because we do not use + for concatenation. If we did, then "foo" + "bar" should produce the regular expression r"foo|bar". String and regular expressions form a semiring under the operations of * for concatenation and + for alternation. So the correct value for zero(String) if any, should be a regular expression that matches nothing, so we could have zero(::Type{String}) = r"(?!)".

[yurivish]

Right, yeah, I realized later that I was confused. Returning regular expressions seems mathematically correct but not very nice for performance (or type-stability).

[lr4d]

I ran into the following error today when trying to use spare arrays with non-numeric data types. I could create a SparseVector{Symbol} but couldn't use it as intended

julia> v = SparseArrays.sparsevec([1;22], [:a, :b])
22-element SparseVector{Symbol,Int64} with 2 stored entries:
  [1 ]  =  :a
  [22]  =  :b

julia> findfirst(isequal(:b), v)
ERROR: MethodError: no method matching zero(::Type{Symbol})

If these types do not really support working with non-numeric values, then I would also suggest to forbid their creation

[KlausC]

As far as I know, the SparseVector and SparseMatrixCSC types have minimal requirements, which are support for iszero and zero. Both are essential for the intended functionalty of the structures - and make not much sense for Symbol or String.
Neverthless it is not desirable restrict the element types to (for example) Number s, because there are other senseful types, which have both iszero and zero.

[abraunst]

I wonder if we could not simply add two new parameters SparseMatrixCSC{Tv,Ti,Z,IZ}, and make the constructors make them default to zero and iszero respectively. This would make using SparseMatrixCSC as storage for non-numerical types possible, on pair with Matrix.

[KlausC]

I would prefer not to add functions as type parameters.
Alternatively introduce new functions like sparse_unstored and is_sparse_unstored, which default to zero and iszero for Number, and which can be redefined for other types.