# Matrix-by-(slice of)vector multiplication with limited allocation

**URL:** https://discourse.julialang.org/t/matrix-by-slice-of-vector-multiplication-with-limited-allocation/46464
**Category:** New to Julia
**Tags:** question
**Created:** [September 11, 2020, 3:24pm UTC](https://discourse.julialang.org/t/matrix-by-slice-of-vector-multiplication-with-limited-allocation/46464 "2020-09-11T15:24:37Z")
**Posts on this page:** 7
**Page:** 1

<div class="post-metadata">

### Author: ![yutaka](https://avatars.discourse-cdn.com/v4/letter/y/dc4da7/32.png) [@yutaka](https://discourse.julialang.org/u/yutaka)
#### Post date: [September 11, 2020, 3:24pm UTC](https://discourse.julialang.org/t/matrix-by-slice-of-vector-multiplication-with-limited-allocation/46464/1 "2020-09-11T15:24:37Z")

</div>

I am trying to multiply a small matrix by a slice of a vector (subvector), i.e., `z[3:7] .= A * x[3:7]`, and it allocates memory for a copy of the slice. A view is a lightweight option, but it still needs some memory and takes time when I do it many times (10 million times or more). Is there any way to perform such an operation with limited memory usage? (like GEMM or SYMM for slices?)

Examples:

```julia
julia> z = zeros(7);

julia> x = ones(7);

julia> A=rand(5,5);

# naive case
julia> @time z[3:7] .= A * x[3:7];
  0.000012 seconds (6 allocations: 384 bytes)

# view for x
julia> @time z[3:7] .= A * (@view x[3:7]);
  0.000012 seconds (6 allocations: 304 bytes)

# view for z
julia> @time vz = @view z[3:7];
  0.000002 seconds (1 allocation: 48 bytes)

julia> @time vz .= A * (@view x[3:7]);
  0.000006 seconds (2 allocations: 176 bytes)

```

---

<div class="post-metadata">

### Author: ![mcabbott](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mcabbott/32/6603_2.png) [@mcabbott](https://discourse.julialang.org/u/mcabbott)
#### Post date: [September 11, 2020, 3:32pm UTC](https://discourse.julialang.org/t/matrix-by-slice-of-vector-multiplication-with-limited-allocation/46464/2 "2020-09-11T15:32:53Z")

</div>

Some ways:

```julia
julia> @btime $z[3:7] .= $A * $x[3:7];
  131.082 ns (2 allocations: 256 bytes)

julia> using LinearAlgebra, Einsum

julia> @btime @views mul!($z[3:7], $A, $x[3:7]);
  69.175 ns (0 allocations: 0 bytes)

julia> emul!(C,A,B) = @einsum C[i] = A[i,k] * B[k]

julia> @btime @views emul!($z[3:7], $A, $x[3:7]);
  26.772 ns (0 allocations: 0 bytes)

```

---

<div class="post-metadata">

### Author: ![yutaka](https://avatars.discourse-cdn.com/v4/letter/y/dc4da7/32.png) [@yutaka](https://discourse.julialang.org/u/yutaka)
#### Post date: [September 11, 2020, 4:14pm UTC](https://discourse.julialang.org/t/matrix-by-slice-of-vector-multiplication-with-limited-allocation/46464/3 "2020-09-11T16:14:22Z")

</div>

Thank you for your quick reply! I tried it but it seems like it doesn’t allocate with `@btime` but does allocate in actual situation. I am afraid I am missing the point.

```julia
using LinearAlgebra, Einsum, TimerOutputs
function alloc_test1()
   reset_timer!()
   for i=1:10000000
      @timeit "time" @views mul!(z[3:7], A, x[3:7]);
   end
end
alloc_test1(); print_timer();

emul!(C,A,B) = @einsum C[i] = A[i,k] * B[k]
function alloc_test2()
   reset_timer!()
   for i=1:10000000
      @timeit "time" @views emul!(z[3:7], A, x[3:7])
   end
end
alloc_test2(); print_timer();

```

```julia
# alloc_test1()
 ──────────────────────────────────────────────────────────────────
                           Time Allocations
                   ────────────────────── ───────────────────────
 Tot / % measured: 2.83s / 81.2% 1.34GiB / 100%

 Section ncalls time %tot avg alloc %tot avg
 ──────────────────────────────────────────────────────────────────
 time 10.0M 2.30s 100% 230ns 1.34GiB 100% 144B
 ──────────────────────────────────────────────────────────────────

# alloc_test2()
 ──────────────────────────────────────────────────────────────────
                           Time Allocations
                   ────────────────────── ───────────────────────
 Tot / % measured: 2.21s / 76.4% 917MiB / 100%

 Section ncalls time %tot avg alloc %tot avg
 ──────────────────────────────────────────────────────────────────
 time 10.0M 1.69s 100% 169ns 917MiB 100% 96.2B
 ──────────────────────────────────────────────────────────────────

```

---

<div class="post-metadata">

### Author: ![mcabbott](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mcabbott/32/6603_2.png) [@mcabbott](https://discourse.julialang.org/u/mcabbott)
#### Post date: [September 11, 2020, 4:27pm UTC](https://discourse.julialang.org/t/matrix-by-slice-of-vector-multiplication-with-limited-allocation/46464/4 "2020-09-11T16:27:03Z")

</div>

I guess this means `@timeit` is unsuitable for measuring nanoseconds. `@btime` has its quirks but it works well, and is already running things many times for you. If you want to check:

```julia
julia> function alloc_test_n(z,A,x,n)
          for i=1:n
             @views mul!(z[3:7], A, x[3:7]);
          end
       end;

julia> alloc_test_n(z,A,x,10^4);

julia> @time alloc_test_n(z,A,x,10^5)
  0.006389 seconds (1 allocation: 16 bytes)

julia> @time alloc_test_n(z,A,x,10^6)
  0.063919 seconds (1 allocation: 16 bytes)

julia> @time alloc_test_n(z,A,x,10^7)
  0.637590 seconds (1 allocation: 16 bytes)

julia> @time alloc_test_n(z,A,x,10^8)
  6.398111 seconds (1 allocation: 16 bytes)

```

---

<div class="post-metadata">

### Author: ![stevengj](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/stevengj/32/71_2.png) [@stevengj](https://discourse.julialang.org/u/stevengj)
#### Post date: [September 11, 2020, 4:35pm UTC](https://discourse.julialang.org/t/matrix-by-slice-of-vector-multiplication-with-limited-allocation/46464/5 "2020-09-11T16:35:28Z")

</div>

> [@yutaka](#):
>
> I am trying to multiply a small matrix by a slice of a vector (subvector), i.e., `z[3:7] .= A * x[3:7]` , and it allocates memory for a copy of the slice. A view is a lightweight option, but it still needs some memory and takes time when I do it many times (10 million times or more). Is there any way to perform such an operation with limited memory usage? (like GEMM or SYMM for slices?)

If your matrices and vectors are this small (\< 10 \times 10) and you are doing millions of operations, then you are may be better off with StaticArrays.jl so that it can inline and unroll the operations (and will also allow you to eliminate allocations, assuming the sizes are fixed). BLAS `gemm` operations (which can also avoid allocations if you drop down to a low enough level) are going to be suboptimal here because they incur a lot of overhead for such small matrices.

e.g.

```julia
julia> using StaticArrays, BenchmarkTools

julia> A = SMatrix{5,5}(rand(5,5));

julia> z = zeros(7); x = ones(7);

julia> @btime $z[3:7] .= $A * SVector(ntuple(i -> $x[2+i], Val{5}()));
  6.219 ns (0 allocations: 0 bytes)

```

(The business with `ntuple` is so I can write out the equivalent of `SVector(x[3],x[4],x[5],x[6],x[7])`, i.e. `SVector(x[3:7]...)` but with the loop unrolled statically. However, if you are working with `StaticArrays` then you may want to re-think your other data structures to use `SVector` as well.)

---

<div class="post-metadata">

### Author: ![kristoffer.carlsson](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/kristoffer.carlsson/32/22_2.png) [@kristoffer.carlsson](https://discourse.julialang.org/u/kristoffer.carlsson)
#### Post date: [September 11, 2020, 5:28pm UTC](https://discourse.julialang.org/t/matrix-by-slice-of-vector-multiplication-with-limited-allocation/46464/6 "2020-09-11T17:28:42Z")

</div>

> [@mcabbott](#):
>
> I guess this means `@timeit` is unsuitable for measuring nanoseconds. ``

Yes, see [GitHub - KristofferC/TimerOutputs.jl: Formatted output of timed sections in Julia](https://github.com/KristofferC/TimerOutputs.jl#overhead).

> There is a small overhead in timing a section (0.25 μs) which means that this package is not suitable for measuring sections that finish very quickly. For proper benchmarking you want to use a more suitable tool like [_BenchmarkTools_](https://github.com/JuliaCI/BenchmarkTools.jl).

---

<div class="post-metadata">

### Author: ![yutaka](https://avatars.discourse-cdn.com/v4/letter/y/dc4da7/32.png) [@yutaka](https://discourse.julialang.org/u/yutaka)
#### Post date: [September 11, 2020, 7:21pm UTC](https://discourse.julialang.org/t/matrix-by-slice-of-vector-multiplication-with-limited-allocation/46464/7 "2020-09-11T19:21:17Z")

</div>

I appreciate all your replies.

@mcabbott Your example is helpful to understand Julia’s behavior! Defining functions often removes the extra allocations.

@stevengj Thank you for your suggestion and explanation. I am still learning a way of optimization in Julia. StaticArray.jl is a definite option for me.

@kristoffer.carlsson Good to know! Thank you for the information.
