Bug in sparse matrix .+= scalar?


#1

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)

#2

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


#3

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?


#4

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


#5

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?


#6

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”.


#7

See e.g.


https://github.com/JuliaLang/julia/issues/7010


#8

Nice!