# Inexplicable allocations when summing \`StaticArrays\`

**URL:** <https://discourse.julialang.org/t/inexplicable-allocations-when-summing-staticarrays/18012>\
**Category:** Performance\
**Created:** [November 26, 2018, 11:58am UTC](https://discourse.julialang.org/t/inexplicable-allocations-when-summing-staticarrays/18012 "2018-11-26T11:58:54Z")\
**Posts on this page:** 9\
**Page:** 1

<div class="post-metadata">

**Author:** ![zygmuntszpak](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/zygmuntszpak/32/2591_2.png) [@zygmuntszpak](https://discourse.julialang.org/u/zygmuntszpak)\
**Post date:** [November 26, 2018, 11:58am UTC](https://discourse.julialang.org/t/inexplicable-allocations-when-summing-staticarrays/18012/1 "2018-11-26T11:58:54Z")

</div>

This is a cross-post from an issue I opened in the `StaticArrays` package: [https://github.com/JuliaArrays/StaticArrays.jl/issues/537](https://github.com/JuliaArrays/StaticArrays.jl/issues/537)

I have been working on a pull-request to fix an issue I discovered with `StaticArrays`. I discovered that when one entity that happens to be a column vector (single column of matrix), and the other a row vector, and you compute the outer product of these two vectors, the returned array is not a `StaticArray` but a regular allocated array. I stumbled upon this issue whilst implementing a particular cost function. I have stripped the cost function and given it as an example below.

I came up with a proposed fix by adding the following two definitions to `matrix_multiply.jl` in the [https://github.com/JuliaArrays/StaticArrays.jl](https://github.com/JuliaArrays/StaticArrays.jl) package:

```julia
@inline *(A::StaticMatrix, B::Adjoint{<:Any, <:StaticVector}) = *(reshape(A, Size(Size(A)[1],)), B) 
@inline mul!(dest::StaticVecOrMat, A::StaticMatrix, B::Adjoint{<:Any, <:StaticVector}) = mul!(dest, reshape(A, Size(Size(A)[1],)), B) 

```

I have verified that adding these additional dispatch rules fixes the original problem (the returned array is now a `StaticArray`). However, when I am summing the result of these `StaticArray` types I am suddenly getting a lot of allocations.

You can reproduce the problem by running the following code snippet:

```julia
using StaticArrays, BenchmarkTools, LinearAlgebra

function hom(v::SVector)
    push(v,1)
end

function T(𝛉::AbstractArray, 𝒞::Tuple{AbstractArray, Vararg{AbstractArray}}, 𝒟::Tuple{AbstractArray, Vararg{AbstractArray}})
    ⊗ = kron
    l = 9
    𝐈ₗ = SMatrix{9,9}(1.0I)
    𝐈ₘ = SMatrix{1,1}(1.0I)
    𝐓 = @SMatrix zeros(9,9)
    N = length(𝒟[1])
    ℳ, ℳʹ = 𝒟
    Λ₁, Λ₂ = 𝒞
    𝚲ₙ = @MMatrix zeros(4,4)
    𝐞₁ = @SMatrix [1.0; 0.0; 0.0]
    𝐞₂ = @SMatrix [0.0; 1.0; 0.0]
    for n = 1: N
        index = SVector(1,2)
        𝚲ₙ[1:2,1:2] .= Λ₁[n][index,index]
        𝚲ₙ[3:4,3:4] .= Λ₂[n][index,index]
        𝐦 = hom(ℳ[n])
        𝐦ʹ= hom(ℳʹ[n])
        𝐔ₙ = (𝐦 ⊗ 𝐦ʹ)
        ∂ₓ𝐮ₙ = [(𝐞₁ ⊗ 𝐦ʹ) (𝐞₂ ⊗ 𝐦ʹ) (𝐦 ⊗ 𝐞₁) (𝐦 ⊗ 𝐞₂)]
        𝐁ₙ = ∂ₓ𝐮ₙ * 𝚲ₙ * ∂ₓ𝐮ₙ'
        𝚺ₙ = 𝛉' * 𝐁ₙ * 𝛉
        𝚺ₙ⁻¹ = inv(𝚺ₙ)
        𝐓₁ = @SMatrix zeros(Float64,9,9)
        for k = 1:l
            𝐞ₖ = 𝐈ₗ[:,k]
            ∂𝐞ₖ𝚺ₙ = (𝐈ₘ ⊗ 𝐞ₖ') * 𝐁ₙ * (𝐈ₘ ⊗ 𝛉) + (𝐈ₘ ⊗ 𝛉') * 𝐁ₙ * (𝐈ₘ ⊗ 𝐞ₖ)
            # Accumulating the result in 𝐓₁ allocates memory, even though
            # the two terms in the summation are both SArrays.
            𝐓₁ = 𝐓₁ + 𝐔ₙ * 𝚺ₙ⁻¹ * (∂𝐞ₖ𝚺ₙ) * 𝚺ₙ⁻¹ * 𝐔ₙ' * 𝛉 * 𝐞ₖ'
        end
        𝐓 = 𝐓 + 𝐓₁
    end
    𝐓
end

# Some sample data
N = 300
ℳ = [@SVector rand(2) for i = 1:N]
ℳʹ = [@SVector rand(2) for i = 1:N]
Λ₁ = [SMatrix{3,3}(Matrix(Diagonal([1.0,1.0,0.0]))) for i = 1:length(ℳ)]
Λ₂ = [SMatrix{3,3}(Matrix(Diagonal([1.0,1.0,0.0]))) for i = 1:length(ℳ)]
F = @SMatrix rand(3,3)
𝒞 = (Λ₁,Λ₂)
𝒟 = (ℳ, ℳʹ)

T(vec(F),𝒞,𝒟)
@btime T(vec($F),$𝒞,$𝒟) # 682.152 μs (6002 allocations: 3.85 MiB

```

I tried `@code_warntype T(vec(F),𝒞,𝒟)` which produced an output too long to list here. However, a particular line jumped out with `Any` :

```julia
1230 ─ %5524 = invoke Base.afoldl(%5506::typeof(*), %5517::SArray{Tuple{9,1},Float64,2,9}, %5475::Adjoint{Float64,SArray{Tuple{9},Float64,1,9}}, _2::SArray{Tuple{9},Float64,1,9}, %5476::Adjoint{Float64,SArray{Tuple{9},Float64,1,9}})::Any

```

I don’t know how to proceed from here and would appreciate any advice. In particular, I am not sure whether my “fix” is missing something, or whether I have stumbled upon a different “bug”.

---

<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:** [November 27, 2018, 9:40am UTC](https://discourse.julialang.org/t/inexplicable-allocations-when-summing-staticarrays/18012/2 "2018-11-27T09:40:54Z")

</div>

The main thing I’d suggest is to try to strip it down to bare essentials. While the example is lovely, there’s a lot of computation here that is presumably irrelevant to the inference problem. If you can digest it down to a couple of key lines it will be much easier to make progress.

---

<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:** [November 27, 2018, 1:52pm UTC](https://discourse.julialang.org/t/inexplicable-allocations-when-summing-staticarrays/18012/3 "2018-11-27T13:52:15Z")

</div>

Rewriting `∂ₓ𝐮ₙ * 𝚲ₙ * ∂ₓ𝐮ₙ'` as something like `∂ₓ𝐮ₙ * (𝚲ₙ * ∂ₓ𝐮ₙ')` might work around it.

---

<div class="post-metadata">

**Author:** ![zygmuntszpak](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/zygmuntszpak/32/2591_2.png) [@zygmuntszpak](https://discourse.julialang.org/u/zygmuntszpak)\
**Post date:** [November 29, 2018, 6:52am UTC](https://discourse.julialang.org/t/inexplicable-allocations-when-summing-staticarrays/18012/4 "2018-11-29T06:52:43Z")

</div>

I’ve managed to construct a minimal example which reproduces the allocation problem:

```julia
using BenchmarkTools, StaticArrays

function test_mem_ok(P,Q)
    return P[1]*Q[1] * P[2]*Q[2] * P[3]*Q[3] * P[4]*Q[4] * P[5]*Q[5] * P[6]*Q[6] * P[7]*Q[7] * P[8]*Q[8] * P[9]*Q[9]
end

function test_mem_bad(P,Q)
    return P[1]*Q[1] * P[2]*Q[2] * P[3]*Q[3] * P[4]*Q[4] * P[5]*Q[5] * P[6]*Q[6] * P[7]*Q[7] * P[8]*Q[8] * P[9]*Q[9] * P[10]*Q[10]
end

function test_mem_fix(P,Q)
    return (P[1]*Q[1] * P[2]*Q[2] * P[3]*Q[3] * P[4]*Q[4] * P[5]*Q[5]) * P[6]*Q[6] * P[7]*Q[7] * P[8]*Q[8] * P[9]*Q[9] * P[10]*Q[10]
end

function test_ok()
    A = @SVector [@SMatrix [1.0] for i = 1:10]
    C = @SVector [adjoint(SVector(1.0)) for i = 1:10]
    test_mem_ok(A,C)
end

function test_bad()
    A = @SVector [@SMatrix [1.0] for i = 1:10]
    C = @SVector [adjoint(SVector(1.0)) for i = 1:10]
    test_mem_bad(A,C)
end

function test_fix()
    A = @SVector [@SMatrix [1.0] for i = 1:10]
    C = @SVector [adjoint(SVector(1.0)) for i = 1:10]
    test_mem_fix(A,C)
end

test_ok()
@time test_ok() # 0.000002 seconds (5 allocations: 176 bytes)
@code_warntype test_ok()

test_bad()
@time test_bad() # 0.000003 seconds (26 allocations: 512 bytes)
@code_warntype test_bad()

test_fix()
@time test_fix() # 0.000001 seconds (5 allocations: 176 bytes)
@code_warntype test_fix()

```

This seems to be a case of [afoldl causes many memory allocations](https://github.com/JuliaLang/julia/issues/13724).

[There was a performance tip warning in the documentation about this issue which was subsequently removed.](https://github.com/JuliaLang/julia/issues/14821)

Evizero had some further insights that he made in [the issue that I originally opened.](https://github.com/JuliaArrays/StaticArrays.jl/issues/537)

Is this something I should report as an issue on Julialang?

---

<div class="post-metadata">

**Author:** ![Ajaychat3](https://avatars.discourse-cdn.com/v4/letter/a/ecd19e/32.png) [@Ajaychat3](https://discourse.julialang.org/u/Ajaychat3)\
**Post date:** [November 29, 2018, 7:18am UTC](https://discourse.julialang.org/t/inexplicable-allocations-when-summing-staticarrays/18012/5 "2018-11-29T07:18:18Z")

</div>

Out of curiousity, what is @time showing allocation but @btime not

```julia
using BenchmarkTools
const k = zeros(20)
function test_mem()
    c=0
    for i in 1: 10
        c += (k[i] *2)  
    end
    c
end

function test(n::Int64)
    ret = 0
    for i = 1:n
        ret += test_mem()
    end
    ret
end
@btime test(100000000)
@time test(100000000)

926.639 ms (0 allocations: 0 bytes)
 0.941509 seconds (238 allocations: 14.453 KiB)

```

---

<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:** [December 2, 2018, 6:48pm UTC](https://discourse.julialang.org/t/inexplicable-allocations-when-summing-staticarrays/18012/6 "2018-12-02T18:48:50Z")

</div>

> [@Ajaychat3](#):
>
> Out of curiousity, what is @time showing allocation but @btime not

Run it twice

```julia
julia> @time test(100000000)
  1.350088 seconds (238 allocations: 14.453 KiB)
0.0

julia> @time test(100000000)
  1.309178 seconds (5 allocations: 176 bytes)
0.0

```

The 5 allocations are from `@time` itself.

---

<div class="post-metadata">

**Author:** ![Ajaychat3](https://avatars.discourse-cdn.com/v4/letter/a/ecd19e/32.png) [@Ajaychat3](https://discourse.julialang.org/u/Ajaychat3)\
**Post date:** [December 3, 2018, 1:19am UTC](https://discourse.julialang.org/t/inexplicable-allocations-when-summing-staticarrays/18012/7 "2018-12-03T01:19:13Z")

</div>

Thanks. Any specific reason we need to run @time twice but @btime just once to get the correct resultant memory allocation? I understand the function got compiled when I first ran @btime in sequence.

Sorry if this appears to be a basic level question.

---

<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:** [December 3, 2018, 1:23am UTC](https://discourse.julialang.org/t/inexplicable-allocations-when-summing-staticarrays/18012/8 "2018-12-03T01:23:06Z")

</div>

`@btime` runs the function many time and returns the result with the shortest run time.

---

<div class="post-metadata">

**Author:** ![Ajaychat3](https://avatars.discourse-cdn.com/v4/letter/a/ecd19e/32.png) [@Ajaychat3](https://discourse.julialang.org/u/Ajaychat3)\
**Post date:** [December 3, 2018, 1:25am UTC](https://discourse.julialang.org/t/inexplicable-allocations-when-summing-staticarrays/18012/9 "2018-12-03T01:25:08Z")

</div>

Thanks for clarification.