Missing matrix-by-sparse-matrix multiplication specializations?

General Usage
1

It looks like the sparse matrix * array specializations are missing a simple-to-add case:

A = sparse (n-by-n, say)
X = wide (k-by-n)
Y = wide (k-by-n)

Y = X*A or mul!(Y,X,A)

seems to fall back to the generic case.

I did a quick search and couldn’t find any issues addressing this… happy to contribute (or work on) a fix if it isn’t know.

Possibly related… https://github.com/JuliaLang/julia/issues/32172

Code to verify the performance difference is below.

julia> using LinearAlgebra, SparseArrays; n = 1000; A = sprand(n,n,5/n); X = rand(2,n); Y = copy(X); @time mul!(Y,X,A); x = rand(n,2); y = copy(x); @time mul!(y,A,x);
  0.647105 seconds (1.27 M allocations: 58.626 MiB, 13.17% gc time)
  0.024182 seconds (82.77 k allocations: 4.418 MiB)

julia> using LinearAlgebra, SparseArrays; n = 1000; A = sprand(n,n,5/n); X = rand(2,n); Y = copy(X); @time mul!(Y,X,A); x = rand(n,2); y = copy(x); @time mul!(y,A,x);
  0.010937 seconds (10 allocations: 496 bytes)
  0.000035 seconds (4 allocations: 160 bytes)
2

This seems to be fixed on master:

julia> versioninfo()
Julia Version 1.4.0-DEV.86
Commit 0849742f73 (2019-09-04 18:28 UTC)
Platform Info:
  OS: Linux (x86_64-linux-gnu)
  CPU: Intel(R) Core(TM) i9-9820X CPU @ 3.30GHz
  WORD_SIZE: 64
  LIBM: libopenlibm
  LLVM: libLLVM-6.0.1 (ORCJIT, skylake)
Environment:
  JULIA_PKG_DEVDIR = /home/patrick/julia-dev

julia> using LinearAlgebra, SparseArrays; n = 1000; A = sprand(n,n,5/n); X = rand(2,n); Y = copy(X); @time mul!(Y,X,A); x = rand(n,2); y = copy(x); @time mul!(y,A,x);
  0.023139 seconds (75.52 k allocations: 3.984 MiB)
  0.021833 seconds (63.96 k allocations: 3.303 MiB)

julia> using LinearAlgebra, SparseArrays; n = 1000; A = sprand(n,n,5/n); X = rand(2,n); Y = copy(X); @time mul!(Y,X,A); x = rand(n,2); y = copy(x); @time mul!(y,A,x);
  0.000073 seconds (4 allocations: 160 bytes)
  0.000069 seconds (4 allocations: 160 bytes)

I didn’t look too closely but it might have been fixed as a result of (https://github.com/JuliaLang/julia/pull/32195), which was written to fix https://github.com/JuliaLang/julia/issues/32172.

I didn’t check to see if this is going to be in 1.3 or any future 1.2.x releases (if there will be any of those?). I would guess it will be in 1.3 but it is probably worth checking and posting an issue if it isn’t already slated for inclusion in 1.3.

1 Like
3

Thanks – for some reason I would have thought that would have made it into 1.2, but I guess I was mistaken. I just checked on 1.3rc1 and it is fixed there!