# Matrix-vector product faster than matrix addition?

**URL:** <https://discourse.julialang.org/t/matrix-vector-product-faster-than-matrix-addition/87954>\
**Category:** Performance\
**Tags:** question, linearalgebra, benchmark\
**Created:** [September 28, 2022, 9:01pm UTC](https://discourse.julialang.org/t/matrix-vector-product-faster-than-matrix-addition/87954 "2022-09-28T21:01:38Z")\
**Posts on this page:** 6\
**Page:** 1

<div class="post-metadata">

**Author:** ![goerz](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/goerz/32/3269_2.png) [@goerz](https://discourse.julialang.org/u/goerz)\
**Post date:** [September 28, 2022, 9:01pm UTC](https://discourse.julialang.org/t/matrix-vector-product-faster-than-matrix-addition/87954/1 "2022-09-28T21:01:38Z")

</div>

I’m running some benchmarks, and to my surprize I’m finding adding two matrices _much_ slower than a matrix-vector multiplication. The [relevant benchmark code is this (see link for context)](https://github.com/goerz-research/2022-09_benchmark_hamiltonian/blob/master/Benchmarks.jl#L161):

```julia
function benchmark_mv_vs_mpm()

    N = 1000
    H = random_hermitian_matrix(N)
    H0 = random_hermitian_matrix(N)
    Ψ = random_state_vector(N)
    ϕ = random_state_vector(N)
    val = 1.15

    println("*** matrix-vector product ϕ = H Ψ")
    b_mv1 = @benchmark mul!($ϕ, $H, $Ψ)
    display(b_mv1)

    println("*** matrix-vector product ϕ += v H Ψ")
    b_mv2 = @benchmark mul!($ϕ, $H, $Ψ, $val, true)
    display(b_mv2)

    println("*** matrix-vector product ϕ += H Ψ")
    b_mv3 = @benchmark mul!($ϕ, $H, $Ψ, true, true)
    display(b_mv3)

    println("*** matrix-matrix addition H += v H0")
    b_mpm1 = @benchmark axpy!($val, $H0, $H)
    display(b_mpm1)

    println("*** matrix-matrix addition H += H0")
    b_mpm2 = @benchmark axpy!(true, $H0, $H)
    display(b_mpm2)

    println("*** matrix-copy H = H0")
    b_mpm3 = @benchmark copyto!($H, $H0)
    display(b_mpm3)

end

```

with the initialization routines in [a separate file `testutils.jl`](https://github.com/goerz-research/2022-09_benchmark_hamiltonian/blob/02b7727d01e7d30bc62c4715e8970fec751dca36/testutils.jl) (not that it matters: everything is just standard dense matrices/vectors)

The resulting runtimes are this:

 ![Screen Shot 2022-09-28 at 16.40.38](https://global.discourse-cdn.com/julialang/original/3X/7/3/73fdd52c1fd36ea93daf9da555e946dc9fd7e195.jpeg)

I guess I’ll take the super-fast matrix-vector multiplication, but this just seems a little odd to me. Fundamentally, both operations should scale as N^2. I don’t think I would have expected that even just copying a matrix would be so much slower than doing a matrix-vector product.

Does anybody have any insights into this behavior?

---

<div class="post-metadata">

**Author:** ![Oscar\_Smith](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/oscar_smith/32/25343_2.png) [@Oscar\_Smith](https://discourse.julialang.org/u/Oscar_Smith)\
**Post date:** [September 28, 2022, 9:11pm UTC](https://discourse.julialang.org/t/matrix-vector-product-faster-than-matrix-addition/87954/2 "2022-09-28T21:11:16Z")

</div>

you should expect matrix addition to be roughly 3x slower because you will be bottlenecked by memory bandwidth and adding 2 matrices into a 3rd requires looking at 3x as much memory.

---

<div class="post-metadata">

**Author:** ![goerz](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/goerz/32/3269_2.png) [@goerz](https://discourse.julialang.org/u/goerz)\
**Post date:** [September 28, 2022, 9:13pm UTC](https://discourse.julialang.org/t/matrix-vector-product-faster-than-matrix-addition/87954/3 "2022-09-28T21:13:39Z")

</div>

Interesting! So this is true in general? That is in other languages, say, Fortran? (I could test this out, but writing benchmark code for Fortran would probably take me a day ;-))

---

<div class="post-metadata">

**Author:** ![Oscar\_Smith](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/oscar_smith/32/25343_2.png) [@Oscar\_Smith](https://discourse.julialang.org/u/Oscar_Smith)\
**Post date:** [September 28, 2022, 9:22pm UTC](https://discourse.julialang.org/t/matrix-vector-product-faster-than-matrix-addition/87954/4 "2022-09-28T21:22:12Z")

</div>

Yeah. Fundamentally, BLAS1 and BLAS2 for largish matrices are linear in the amount of data, so their runtime will be directly correlated with memory visited.

---

<div class="post-metadata">

**Author:** ![lmiq](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/lmiq/32/18314_2.png) [@lmiq](https://discourse.julialang.org/u/lmiq)\
**Post date:** [September 28, 2022, 10:58pm UTC](https://discourse.julialang.org/t/matrix-vector-product-faster-than-matrix-addition/87954/5 "2022-09-28T22:58:13Z")

</div>

> [@goerz](#):
>
> That is in other languages, say, Fortran?

And these are likely Fortran routines: [?axpy](https://www.intel.com/content/www/us/en/develop/documentation/onemkl-developer-reference-fortran/top/blas-and-sparse-blas-routines/blas-routines/blas-level-1-routines-and-functions/axpy.html)

---

<div class="post-metadata">

**Author:** ![goerz](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/goerz/32/3269_2.png) [@goerz](https://discourse.julialang.org/u/goerz)\
**Post date:** [September 28, 2022, 11:04pm UTC](https://discourse.julialang.org/t/matrix-vector-product-faster-than-matrix-addition/87954/6 "2022-09-28T23:04:51Z")

</div>

Oh, yeah, definitely! I’m using those in my Fortran code. I’m also using the matrix-vector products instead of matrix-matrix sums in Fortran, but in my Julia prototype, I decided to sum all my operators first, because I didn’t think it would make _that_ much of a difference and it was easier to get started.

Nobody has written `BenchmarkTools` for Fortran, which makes benchmarking a whole lot more time consuming, and I never thoroughly investigated the two alternatives 😉

In any case, this is good news, because my Julia code is going to become quite significantly faster! 🙂
