Bug in sparse matrix .+= scalar?

pablosanjose · August 1, 2017, 6:29pm

I guess this is a bug, but I’d like some feedback before opening an issue in github, in case there is already an open discussion somewhere

In v0.5.1 (correct behavior)

julia> sp = sprand(10,10,.05)
10×10 sparse matrix with 5 Float64 nonzero entries:
	[8 ,  1]  =  0.384875
	[7 ,  6]  =  0.57818
	[5 ,  9]  =  0.721343
	[6 ,  9]  =  0.776902
	[4 , 10]  =  0.420691

julia> sp .+= 1.0
10×10 sparse matrix with 100 Float64 nonzero entries:
	[1 ,  1]  =  1.0
	[2 ,  1]  =  1.0
	[3 ,  1]  =  1.0
	[4 ,  1]  =  1.0
	[5 ,  1]  =  1.0
	[6 ,  1]  =  1.0
	[7 ,  1]  =  1.0
	⋮
	[3 , 10]  =  1.0
	[4 , 10]  =  1.42069
	[5 , 10]  =  1.0
	[6 , 10]  =  1.0
	[7 , 10]  =  1.0
	[8 , 10]  =  1.0
	[9 , 10]  =  1.0
	[10, 10]  =  1.0

However, I see some wacko behaviour in current master (also in v0.6)

julia> sp = sprand(10,10,.05)
10×10 SparseMatrixCSC{Float64,Int64} with 7 stored entries:
  [9 ,  2]  =  0.116704
  [8 ,  3]  =  0.852258
  [9 ,  3]  =  0.949641
  [2 ,  5]  =  0.998075
  [3 ,  6]  =  0.912836
  [2 ,  7]  =  0.304016
  [3 , 10]  =  0.896151

julia> sp .+= 1.0
10×10 SparseMatrixCSC{Float64,Int64} with 100 stored entries:
  [1 ,  1]  =  2.0
  [2 ,  1]  =  2.0
  [3 ,  1]  =  2.0
  [4 ,  1]  =  2.0
  [5 ,  1]  =  2.0
  [6 ,  1]  =  2.0
  [7 ,  1]  =  2.0
  [8 ,  1]  =  2.0
  [9 ,  1]  =  2.0
  ⋮
  [2 , 10]  =  2.0
  [3 , 10]  =  2.0
  [4 , 10]  =  2.0
  [5 , 10]  =  2.0
  [6 , 10]  =  2.0
  [7 , 10]  =  2.0
  [8 , 10]  =  2.0
  [9 , 10]  =  2.0
  [10, 10]  =  2.0

The 2.0 is twice the scalar provided in sp .+= scalar. All previous entries are overwritten.

Version information:

julia> versioninfo()
Julia Version 0.7.0-DEV.1165
Commit 1a43098cf7 (2017-07-31 03:33 UTC)
Platform Info:
  OS: macOS (x86_64-apple-darwin13.4.0)
  CPU: Intel(R) Core(TM) i7-2600 CPU @ 3.40GHz
  WORD_SIZE: 64
  BLAS: libopenblas (USE64BITINT DYNAMIC_ARCH NO_AFFINITY Sandybridge)
  LAPACK: libopenblas64_
  LIBM: libopenlibm
  LLVM: libLLVM-3.9.1 (ORCJIT, sandybridge)

kristoffer.carlsson · August 1, 2017, 7:08pm

Ref: https://github.com/JuliaLang/julia/issues/22183, https://github.com/JuliaLang/julia/issues/21693

pablosanjose · August 1, 2017, 8:12pm

Thanks Kristoffer, I suspected this might be a known issue. However, it seems the idea is to make .+= an error on sparse matrices? I find this mildly alarming.

Actually, coming from Mathematica, I had expected sparse .+= scalar to not be an issue, really. In Mathematica sparse matrices have a “default value” for non-structural entries that can be zero or non-zero. Hence sparse .+= scalar just modifies the default value, as well as the structural entries, which does not change the sparcity pattern of the matrix at all. Do you know if this kind of design has been discussed in Julia?

ChrisRackauckas · August 1, 2017, 8:17pm

Yes. And I think everyone agrees we should do it, just someone needs to do it.

pablosanjose · August 1, 2017, 8:36pm

Hi Chris, can you point us to any particular issue or thread where we can read about this discussion, and the implications of such a redesign?

ChrisRackauckas · August 1, 2017, 8:42pm

I don’t know the issue off the top of my head, but it’s one of the sparse broadcast ones. But we were talking about it in the Slack just like 2 days ago. Might need to be “formalized”.

stevengj · August 1, 2017, 10:10pm

See e.g.

github.com/JuliaLang/julia

Feature request: default (non-zero) values for sparse matrices

opened 05:56PM - 05 Mar 15 UTC

closed 08:50PM - 05 Mar 15 UTC

sbromberger

sparse

Is there any practical way to achieve this? Consider: I have a very large matrix… describing edge weights between all pairs of vertices. This is an NxN matrix and is impractical for graphs with orders of 1e7 without using sparse matrices. If I want the default edge weight of 1.0, I have to code specially around this to make a sparse matrix entry of "0" really be "1" (which is ok in this specific case because an edge weight of 0 is not defined). It would be much easier if I could specify "create this sparse matrix but any unspecified value should be X". Is there already an easy way to do this that I've missed?

github.com/JuliaSparse/SparseArrays.jl

Map for sparse matrices

opened 01:54PM - 28 May 14 UTC

closed 12:49AM - 17 Jan 22 UTC

jumutc

Please consider the code below where `julia` is converting sparse matrix into de…nse one and doing `map` on the top. Should we have `map` defined only for non-zero elements in sparse matrix instead? Otherwise the slow down is an issue and one cannot apply effectively λ-calculus in this case! ``` julia julia> test = sprand(10^3,10^3,.01) 1000x1000 sparse matrix with 10144 Float64 entries: [80 , 1] = 0.993039 [128 , 1] = 0.117601 [152 , 1] = 0.974119 [259 , 1] = 0.0362442 [289 , 1] = 0.621536 [371 , 1] = 0.653076 [631 , 1] = 0.131718 ⋮ [439 , 1000] = 0.0621062 [538 , 1000] = 0.109039 [613 , 1000] = 0.212955 [620 , 1000] = 0.147798 [640 , 1000] = 0.479203 [702 , 1000] = 0.88309 [884 , 1000] = 0.780324 [892 , 1000] = 0.0164652 julia> f = v -> v+1 (anonymous function) julia> @time map(f,test) elapsed time: 14.162436633 seconds (89388704 bytes allocated) 1000x1000 sparse matrix with 1000000 Float64 entries: [1 , 1] = 1.0 [2 , 1] = 1.0 [3 , 1] = 1.0 [4 , 1] = 1.0 [5 , 1] = 1.0 [6 , 1] = 1.0 [7 , 1] = 1.0 ⋮ [993 , 1000] = 1.0 [994 , 1000] = 1.0 [995 , 1000] = 1.0 [996 , 1000] = 1.0 [997 , 1000] = 1.0 [998 , 1000] = 1.0 [999 , 1000] = 1.0 [1000, 1000] = 1.0 julia> @time (I,J,V)=findnz(test); sparse(I,J,map(f,V)) elapsed time: 0.000144994 seconds (243784 bytes allocated) 1000x1000 sparse matrix with 10144 Float64 entries: [80 , 1] = 1.99304 [128 , 1] = 1.1176 [152 , 1] = 1.97412 [259 , 1] = 1.03624 [289 , 1] = 1.62154 [371 , 1] = 1.65308 [631 , 1] = 1.13172 ⋮ [439 , 1000] = 1.06211 [538 , 1000] = 1.10904 [613 , 1000] = 1.21296 [620 , 1000] = 1.1478 [640 , 1000] = 1.4792 [702 , 1000] = 1.88309 [884 , 1000] = 1.78032 [892 , 1000] = 1.01647 ```

github.com/JuliaLang/julia

`map` on sparse arrays does not work with non-numeric data

opened 03:44AM - 12 Dec 16 UTC

closed 05:39AM - 22 May 17 UTC

amitmurthy

regression sparse

``` julia> map(x-> x != 0.0 ? "Hello" : x, sparse(eye(4))) ERROR: MethodError:… no method matching zero(::Type{Any}) ``` Also a problem when mixing ints and floats. ``` julia> map(x-> x != 0.0 ? 1 : x, sparse(eye(4))) ERROR: MethodError: no method matching zero(::Type{Any}) ``` If sparse arrays are not meant to be used with non-numeric data, we should at least throw a better error message.

pablosanjose · August 1, 2017, 10:52pm

Nice!