# Slow arithmetic on views of sparse matrices

**URL:** <https://discourse.julialang.org/t/slow-arithmetic-on-views-of-sparse-matrices/3644>\
**Category:** General Usage\
**Created:** [May 11, 2017, 8:56am UTC](https://discourse.julialang.org/t/slow-arithmetic-on-views-of-sparse-matrices/3644 "2017-05-11T08:56:19Z")\
**Posts on this page:** 13\
**Page:** 1

<div class="post-metadata">

**Author:** ![mkborregaard](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mkborregaard/32/556_2.png) [@mkborregaard](https://discourse.julialang.org/u/mkborregaard)\
**Post date:** [May 11, 2017, 8:56am UTC](https://discourse.julialang.org/t/slow-arithmetic-on-views-of-sparse-matrices/3644/1 "2017-05-11T08:56:20Z")

</div>

I was profiling some code that was surprisingly slow, and it turned out that the source of the issue was a sum over a view into a large sparse matrix. Here is an MWE that recreates it:

```julia
d = sprand(Bool,10000,10000, 0.01)
e = view(d, rand(1:10000,5000), rand(1:10000,9000))

using BenchmarkTools
@benchmark sum($d, 1)
BenchmarkTools.Trial:
  memory estimate: 78.20 KiB
  allocs estimate: 2
  --------------
  minimum time: 421.354 μs (0.00% GC)
  median time: 450.828 μs (0.00% GC)
  mean time: 463.575 μs (0.00% GC)
  maximum time: 1.226 ms (0.00% GC)
  --------------
  samples: 10000
  evals/sample: 1

@benchmark sum($e, 1)
BenchmarkTools.Trial:
  memory estimate: 585.06 KiB
  allocs estimate: 51
  --------------
  minimum time: 3.286 s (0.00% GC)
  median time: 3.312 s (0.00% GC)
  mean time: 3.312 s (0.00% GC)
  maximum time: 3.339 s (0.00% GC)
  --------------
  samples: 2
  evals/sample: 1

```

That’s a slowdown by a factor of 7000.

Is this 1) expected behaviour and I’m best off not using views into sparse matrices?; 2) am I doing something wrong?; or 3) is it a bug?

EDIT: `f = view(d, :, :)` has the same issue.

Thanks!

---

<div class="post-metadata">

**Author:** ![mauro3](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mauro3/32/292_2.png) [@mauro3](https://discourse.julialang.org/u/mauro3)\
**Post date:** [May 11, 2017, 9:44am UTC](https://discourse.julialang.org/t/slow-arithmetic-on-views-of-sparse-matrices/3644/2 "2017-05-11T09:44:27Z")

</div>

I would think this is expected: The sum of a sparse CSC matrix is trivial: just sum the internal vector storing the values of non-zeros. Your sum of the view on the other hand, uses the generic `sum` function: `sum(A::AbstractArray, region) at reducedim.jl:320`. This will look up each index to sum it, which for a CSC sparse matrix is two pointer indirections, most of which will lead to “zero”. Thus 7000x slowdown seems about right. I suspect to make a specialized `sum` function specialized to random access sparse matrix will be hard and not much more efficient (although I might be wrong).

---

<div class="post-metadata">

**Author:** ![tkelman](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/tkelman/32/692_2.png) [@tkelman](https://discourse.julialang.org/u/tkelman)\
**Post date:** [May 11, 2017, 9:48am UTC](https://discourse.julialang.org/t/slow-arithmetic-on-views-of-sparse-matrices/3644/3 "2017-05-11T09:48:48Z")

</div>

Depending on what your real use case is, if you’re trying to sum elements only within some contiguous subrange of rows and columns, it wouldn’t be hard to write a specialized routine for SparseMatrixCSC storage. Views of sparse matrices are not well optimized since they fall back to generic operations which aren’t aware of the sparsity.

---

<div class="post-metadata">

**Author:** ![mkborregaard](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mkborregaard/32/556_2.png) [@mkborregaard](https://discourse.julialang.org/u/mkborregaard)\
**Post date:** [May 11, 2017, 9:55am UTC](https://discourse.julialang.org/t/slow-arithmetic-on-views-of-sparse-matrices/3644/4 "2017-05-11T09:55:44Z")

</div>

So is the recommended approach to dispense with the view and just slice the matrix instead?

---

<div class="post-metadata">

**Author:** ![tkelman](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/tkelman/32/692_2.png) [@tkelman](https://discourse.julialang.org/u/tkelman)\
**Post date:** [May 11, 2017, 10:27am UTC](https://discourse.julialang.org/t/slow-arithmetic-on-views-of-sparse-matrices/3644/5 "2017-05-11T10:27:21Z")

</div>

Slicing would make a copy, which you don’t necessarily need. Depends how sensitive you are to the performance here vs how much special-purpose code you want to write for a particular operation.

---

<div class="post-metadata">

**Author:** ![mkborregaard](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mkborregaard/32/556_2.png) [@mkborregaard](https://discourse.julialang.org/u/mkborregaard)\
**Post date:** [May 11, 2017, 10:30am UTC](https://discourse.julialang.org/t/slow-arithmetic-on-views-of-sparse-matrices/3644/6 "2017-05-11T10:30:44Z")

</div>

Yes, I realize it’ll make a copy, but in the end I think it’ll be a lot faster. Though, even though they are sparse 1e6 filled elements in the matrix is realistic for my use. I am sensitive to the performance though, the issue is that the elements are never contiguous in my actual use case. I think writing the specialized code is beyond me.

---

<div class="post-metadata">

**Author:** ![mkborregaard](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mkborregaard/32/556_2.png) [@mkborregaard](https://discourse.julialang.org/u/mkborregaard)\
**Post date:** [May 11, 2017, 11:13am UTC](https://discourse.julialang.org/t/slow-arithmetic-on-views-of-sparse-matrices/3644/7 "2017-05-11T11:13:29Z")

</div>

Do any of you know where the specialized sum code for SparseMatrixCSC is located? I can’t seem to find it in Base.

---

<div class="post-metadata">

**Author:** ![tkelman](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/tkelman/32/692_2.png) [@tkelman](https://discourse.julialang.org/u/tkelman)\
**Post date:** [May 11, 2017, 11:19am UTC](https://discourse.julialang.org/t/slow-arithmetic-on-views-of-sparse-matrices/3644/8 "2017-05-11T11:19:02Z")

</div>

Probably calls into this? [https://github.com/JuliaLang/julia/blob/4dbfe4b31deacc7cb586ff0093525724382c42c9/base/sparse/sparsematrix.jl#L1594](https://github.com/JuliaLang/julia/blob/4dbfe4b31deacc7cb586ff0093525724382c42c9/base/sparse/sparsematrix.jl#L1594)

---

<div class="post-metadata">

**Author:** ![mauro3](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mauro3/32/292_2.png) [@mauro3](https://discourse.julialang.org/u/mauro3)\
**Post date:** [May 11, 2017, 11:28am UTC](https://discourse.julialang.org/t/slow-arithmetic-on-views-of-sparse-matrices/3644/9 "2017-05-11T11:28:55Z")

</div>

`@edit sum(sq)` to look at the code, but essentially it boils down to `sum(s.nzval)`. Anyway, have a look at the doc of `nzrange` for an example of an iteration over the non-zeros of a sparse matrix.

---

<div class="post-metadata">

**Author:** ![Per](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/per/32/10387_2.png) [@Per](https://discourse.julialang.org/u/Per)\
**Post date:** [May 11, 2017, 12:23pm UTC](https://discourse.julialang.org/t/slow-arithmetic-on-views-of-sparse-matrices/3644/10 "2017-05-11T12:23:15Z")

</div>

If you use a matrix product instead of a view to compute the sum, it’ll be fast.

```julia
d = sprand(Bool,10000,10000, 0.01)
i = rand(1:10000,5000)
j = rand(1:10000,9000)
e = view(d, i, j)

s = sum(e,1)
t = (sparse(ones(Int,length(i)), i, ones(Int,length(i)),1,10000)*d)[j]
@assert s[:] == t

```

---

<div class="post-metadata">

**Author:** ![mkborregaard](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mkborregaard/32/556_2.png) [@mkborregaard](https://discourse.julialang.org/u/mkborregaard)\
**Post date:** [May 11, 2017, 12:31pm UTC](https://discourse.julialang.org/t/slow-arithmetic-on-views-of-sparse-matrices/3644/11 "2017-05-11T12:31:21Z")

</div>

That’s plain amazing 🙂

---

<div class="post-metadata">

**Author:** ![mkborregaard](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mkborregaard/32/556_2.png) [@mkborregaard](https://discourse.julialang.org/u/mkborregaard)\
**Post date:** [May 11, 2017, 12:35pm UTC](https://discourse.julialang.org/t/slow-arithmetic-on-views-of-sparse-matrices/3644/12 "2017-05-11T12:35:52Z")

</div>

I’ve opened a feature request for the general functionality: [https://github.com/JuliaLang/julia/issues/21796](https://github.com/JuliaLang/julia/issues/21796)

---

<div class="post-metadata">

**Author:** ![tim.holy](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/tim.holy/32/52_2.png) [@tim.holy](https://discourse.julialang.org/u/tim.holy)\
**Post date:** [May 11, 2017, 1:35pm UTC](https://discourse.julialang.org/t/slow-arithmetic-on-views-of-sparse-matrices/3644/13 "2017-05-11T13:35:59Z")

</div>

The real problem is one of iteration—you want to visit just the relevant entries in an arbitrary `AbstractArray`. [https://github.com/timholy/ArrayIteration.jl](https://github.com/timholy/ArrayIteration.jl) might be a solution some day, but it needs some compiler improvements to be blazingly fast. (See [https://github.com/JuliaLang/julia/issues/21796#issuecomment-300789568](https://github.com/JuliaLang/julia/issues/21796#issuecomment-300789568) for more information.)

It still might be worth playing with—in particular, it might be worth “finishing” (I haven’t had the time lately) and registering as a package, even if it’s not as fast as it could be someday. If anyone wants to join the fun, I’d be happy to help coach folks through anything that’s not clear.