dgleich
September 16, 2019, 4:03pm
#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)
```

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

dgleich
September 16, 2019, 5:26pm
#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!