# Avoiding allocations of small but non-trivial arrays (work array alternative?)

**URL:** <https://discourse.julialang.org/t/avoiding-allocations-of-small-but-non-trivial-arrays-work-array-alternative/90084>\
**Category:** Performance\
**Tags:** question\
**Created:** [November 11, 2022, 12:17am UTC](https://discourse.julialang.org/t/avoiding-allocations-of-small-but-non-trivial-arrays-work-array-alternative/90084 "2022-11-11T00:17:02Z")\
**Posts on this page:** 19\
**Page:** 2

<div class="post-metadata">

**Author:** ![jmair](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/jmair/32/35117_2.png) [@jmair](https://discourse.julialang.org/u/jmair)\
**Post date:** [November 12, 2022, 6:58am UTC](https://discourse.julialang.org/t/avoiding-allocations-of-small-but-non-trivial-arrays-work-array-alternative/90084/21 "2022-11-12T06:58:12Z")

</div>

Strange that it doesn’t work, I am only trying on my phone (Arm), but it seems to run fine for me. Maybe try:

```julia
for k in eachindex(xnodes)
      x = Int(k==j)
      barw_j *= (xnodes[j] - xnodes[k]) + x
 end

```

Note that when I tried it, I also used `eachindex` instead of `1:n` but I’m not sure if it made a difference.

Or maybe

```julia
for k in eachindex(xnodes)
      x = ifelse(k==j, 1.0,0.0)
      barw_j *= (xnodes[j] - xnodes[k]) + x
 end

```

---

<div class="post-metadata">

**Author:** ![DNF](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dnf/32/10191_2.png) [@DNF](https://discourse.julialang.org/u/DNF)\
**Post date:** [November 12, 2022, 7:11am UTC](https://discourse.julialang.org/t/avoiding-allocations-of-small-but-non-trivial-arrays-work-array-alternative/90084/22 "2022-11-12T07:11:26Z")

</div>

> [@martin.d.maas](#):
>
> Yes, the `MArrays` are working fine for me now, in combination with Value-types.

As I mentioned upthread, you don’t need the `Val` type. Passing the length of an array as a separate argument does not look good, imo. If you know `n` statically, why are you not encoding it in the type of `xnodes` instead, as a static array?

---

<div class="post-metadata">

**Author:** ![martin.d.maas](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/martin.d.maas/32/50964_2.png) [@martin.d.maas](https://discourse.julialang.org/u/martin.d.maas)\
**Post date:** [November 12, 2022, 7:19am UTC](https://discourse.julialang.org/t/avoiding-allocations-of-small-but-non-trivial-arrays-work-array-alternative/90084/23 "2022-11-12T07:19:46Z")

</div>

> [@DNF](#):
>
> As I mentioned upthread, you don’t need the `Val` type. Passing the length of an array as a separate argument does not look good, imo. If you know `n` statically, why are you not encoding it in the type of `xnodes` instead, as a static array?

I agree that passing the length of an array as a separate argument doesn’t look good. It’d be nice to find a way to get rid of the extra function argument, I just don’t know how to do that in my real-life scenario, where `xnodes` is actually a `@view` to a small subset of a large, standard Array. I’m not sure if I could re-cast that as a StaticArray…

---

<div class="post-metadata">

**Author:** ![DNF](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dnf/32/10191_2.png) [@DNF](https://discourse.julialang.org/u/DNF)\
**Post date:** [November 12, 2022, 7:53am UTC](https://discourse.julialang.org/t/avoiding-allocations-of-small-but-non-trivial-arrays-work-array-alternative/90084/24 "2022-11-12T07:53:34Z")

</div>

Yeah, so this isn’t all that obvious, but I think the key is indexing into the large array with an `SVector` as index. Maybe there’s a more elegant solution, but this works for me:

```julia
function staticslice(x, start, ::Val{N}) where {N}
    ind = SVector{N}(start:start+N-1)
    return x[ind]
end

julia> x = rand(10_000);

julia> @btime staticslice($x, 200, Val(4))
  2.100 ns (0 allocations: 0 bytes)
4-element SVector{4, Float64} with indices SOneTo(4):
 0.521546493646345
 0.33167560079710856
 0.08923575401846218
 0.8255242560839948

```

This way you also get possible benefits of `xnodes` being an SVector.

---

<div class="post-metadata">

**Author:** ![lrnv](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/lrnv/32/19373_2.png) [@lrnv](https://discourse.julialang.org/u/lrnv)\
**Post date:** [November 12, 2022, 11:00am UTC](https://discourse.julialang.org/t/avoiding-allocations-of-small-but-non-trivial-arrays-work-array-alternative/90084/25 "2022-11-12T11:00:01Z")

</div>

Moving to Svector’s, switching a few things inside your function to reduce divisions, merging the two loops on j, moving the stopping condition upper in the function, and a few other modifications, gave me a 4x to 5x speedup :

```julia
using LoopVectorization
using StaticArrays
using BenchmarkTools

function Interp1d_no_array( xnodes, fvals, t )
    n = length(xnodes)
    numt, denomt = 0, 0

    @inbounds for j ∈ eachindex(xnodes)
        barw_j = 1.0
        for k ∈ 1:n
            if k != j
            barw_j *= (xnodes[j] - xnodes[k])
            end
        end
        barw_j = 1.0 / barw_j
        tdiff = t - xnodes[j]
        numt = numt + barw_j / tdiff * fvals[j]
        denomt = denomt + barw_j / tdiff
        if ( abs(tdiff) < 1e-14 )
            numt = fvals[j]
            denomt = 1.0
            break
        end
    end
    return numt / denomt
end
function run_Interp1d_no_array()
    xnodes = LinRange(0,1,4)
    fvals = xnodes.^4
    t = 0.3
    
    for i ∈ 1:100000
        u = Interp1d_no_array( xnodes, fvals, t )
    end
end

function Interp1d_MVector( xnodes, fvals, t, ::Val{n} ) where n 
    
    barw = zeros(MVector{n,Float64})
    for j ∈ 1:n
        barw[j] = 1.0
        for k ∈ 1:n
            if k != j
                barw[j] = barw[j] * (xnodes[j] - xnodes[k])
            end
        end
        barw[j] = 1.0/barw[j]
    end

    numt, denomt = 0, 0
    @inbounds for j ∈ eachindex(xnodes)
        tdiff = t - xnodes[j]
        numt = numt + barw[j] / tdiff * fvals[j]
        denomt = denomt + barw[j] / tdiff
        if ( abs(tdiff) < 1e-14 )
            numt = fvals[j]
            denomt = 1.0
            break
        end
    end
    return numt / denomt
end
function run_Interp1d_MVector()
    xnodes = LinRange(0,1,4)
    fvals = xnodes.^4
    t = 0.3
    
    n = Val(length(xnodes))
    for i ∈ 1:100000
        u = Interp1d_MVector( xnodes, fvals, t, n )
    end
end

function Interp1d_Hardcoded( xn, fvals, t )
    
    barw = SA[1.0/((xn[1]-xn[2])*(xn[1]-xn[3])*(xn[1]-xn[4])), 1.0/((xn[2]-xn[1])*(xn[2]-xn[3])*(xn[2]-xn[4])), 1.0/((xn[3]-xn[1])*(xn[3]-xn[2])*(xn[3]-xn[4])), 1.0/((xn[4]-xn[1])*(xn[4]-xn[2])*(xn[4]-xn[3])) ]

    numt, denomt = 0, 0
    @inbounds for j ∈ eachindex(xn)
        tdiff = t - xn[j]
        numt = numt + barw[j] / tdiff * fvals[j]
        denomt = denomt + barw[j] / tdiff
        if ( abs(tdiff) < 1e-14 )
            numt = fvals[j]
            denomt = 1.0
            break
        end
    end
    return numt / denomt
end
function run_Interp1d_Hardcoded()
    xnodes = LinRange(0,1,4)
    fvals = xnodes.^4
    t = 0.3
    
    for i ∈ 1:100000
        u = Interp1d_Hardcoded( xnodes, fvals, t )
    end
end

function lrnv(xnodes,fvals, t) 

    numt, denomt = 0.0, 0.0
    n = length(xnodes)

    @inbounds for j ∈ eachindex(xnodes)
        tdiff = t - xnodes[j]
        if (abs(tdiff) < 1e-14)
            return fvals[j]
        end
    end

    @inbounds for j ∈ eachindex(xnodes)
        
        xⱼ = xnodes[j]
        bⱼ = (t-xⱼ)
        for k ∈ 1:n
            if k != j
                bⱼ *= (xⱼ - xnodes[k])
            end
        end
        inv_bⱼ = 1 / bⱼ
        numt += inv_bⱼ * fvals[j]
        denomt += inv_bⱼ
    end
    return numt / denomt
end
function run_lrnv()
    xnodes = SVector{4}(LinRange(0,1,4))
    fvals = xnodes.^4
    t = 0.3
    u = 0.0
    for i ∈ 1:100000
        u += lrnv( xnodes,fvals, t)
    end
    return u
end

@btime run_Interp1d_no_array() # 4.984 ms on my machine, 1 allocation
@btime run_Interp1d_MVector() # 2.846 ms on my machine, 1 allocation
@btime run_Interp1d_Hardcoded() # 3.709 ms on my machine, 1 allocation
@btime run_lrnv() # 0.728 ms on my machine, 0 allocation

```

I appologies for the @turbo non-working proposal, i was on my phone and did not ran the code. I tried my best but the `k != j` is a deal bearker for @turbo, sadly.

There is probably something smart to do to enable @simd or @turbo to work, but i did not find it.

---

<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:** [November 12, 2022, 11:31am UTC](https://discourse.julialang.org/t/avoiding-allocations-of-small-but-non-trivial-arrays-work-array-alternative/90084/26 "2022-11-12T11:31:13Z")

</div>

In all cases, it is not clear if you know the size of `xnodes` at compile time. We are all assuming it is known to be `4`. If we remove that restriction, by passing that size as a parameter, the allocations come back, even with the MVector version.

One solution is to introduce a function barrier on the loop:

```julia
using StaticArrays

function Interp1d_MVector( xnodes, fvals, t, ::Val{n} ) where n 
    
    barw = zeros(MVector{n,Float64})
    for j ∈ 1:n
        barw[j] = 1.0
        for k ∈ 1:n
            if k != j
                barw[j] = barw[j] * (xnodes[j] - xnodes[k])
            end
        end
        barw[j] = 1.0/barw[j]
    end

    numt, denomt = 0, 0
    @inbounds for j ∈ eachindex(xnodes)
        tdiff = t - xnodes[j]
        numt = numt + barw[j] / tdiff * fvals[j]
        denomt = denomt + barw[j] / tdiff
        if ( abs(tdiff) < 1e-14 )
            numt = fvals[j]
            denomt = 1.0
            break
        end
    end
    return numt / denomt
end

function the_loop(xnodes, fvals, t, n)
    for i ∈ 1:100000
        u = Interp1d_MVector(xnodes, fvals, t, n)
    end
end

function run_Interp1d_MVector(N)
    xnodes = LinRange(0,1,N)
    fvals = xnodes.^4
    t = 0.3
    n = Val(length(xnodes))
    the_loop(xnodes, fvals, t, n)
end

```

which gives:

```julia
julia> @btime run_Interp1d_MVector(4)
  3.586 ms (2 allocations: 144 bytes)

```

We can make that faster by removing some small type instability (you had initialized `numt` and `denomt` as integers), and also by making `xnodes` a static array:

```julia
using StaticArrays

function Interp1d_MVector2(xnodes::SVector{N,T}, fvals::SVector{N,T}, t) where {N,T}
    barw = zeros(MVector{N,T})
    for j ∈ 1:N
        barw[j] = one(T)
        for k ∈ 1:N
            if k != j
                barw[j] = barw[j] * (xnodes[j] - xnodes[k])
            end
        end
        barw[j] = inv(barw[j])
    end

    numt, denomt = zero(T), zero(T)
    @inbounds for j ∈ eachindex(xnodes)
        tdiff = t - xnodes[j]
        numt = numt + barw[j] / tdiff * fvals[j]
        denomt = denomt + barw[j] / tdiff
        if ( abs(tdiff) < 1e-14 )
            numt = fvals[j]
            denomt = one(T)
            break
        end
    end
    return numt / denomt
end

function the_loop(xnodes, fvals, t)
    for i ∈ 1:100000
        u = Interp1d_MVector2(xnodes, fvals, t)
    end
end

function run_Interp1d_MVector2(N)
    xnodes = SVector{N,Float64}(LinRange(0,1,4))
    fvals = xnodes.^4
    t = 0.3
    the_loop(xnodes, fvals, t)
end

```

which gives:

```julia
julia> @btime run_Interp1d_MVector2(4)
1.036 ms (16 allocations: 800 bytes)                                                                                                         

```

Probably with some other tweaks proposed above, it can be even faster. But you won´t get that complete code to be completely non-allocating if one assumes that the length of `xnodes` cannot be known at compile time. What is important is to isolate that from the hot loop.

edit: Seems you can use `@turbo` elimnating the branch with `ifelse`:

```julia
    @turbo for j ∈ 1:N
        barw[j] = one(T)
        for k ∈ 1:N
            barw[j] = ifelse(k == j , barw[j] , barw[j] * (xnodes[j] - xnodes[k]) )
        end
        barw[j] = inv(barw[j])
    end

```

With which I get:

```julia
julia> @btime run_Interp1d_MVector2(4)
584.132 μs (16 allocations: 800 bytes)

```

(I’m not checking the results)

The whole run becomes non-allocating if `N` is hard coded, but the time does not change much.

---

<div class="post-metadata">

**Author:** ![Elrod](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/elrod/32/22461_2.png) [@Elrod](https://discourse.julialang.org/u/Elrod)\
**Post date:** [November 12, 2022, 11:45am UTC](https://discourse.julialang.org/t/avoiding-allocations-of-small-but-non-trivial-arrays-work-array-alternative/90084/27 "2022-11-12T11:45:23Z")

</div>

> [@lmiq](#):
>
> ```julia
> @turbo for j ∈ 1:N
> barw[j] = one(T)
> for k ∈ 1:N
> barw[j] = ifelse(k == j , barw[j] , barw[j] * (xnodes[j] - xnodes[k]) )
> end
> barw[j] = inv(barw[j])
> end
> 
> ```

Try replacing `1:N` with `eachindex(barw)`.  
That should give LV some more information.

---

<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:** [November 12, 2022, 11:46am UTC](https://discourse.julialang.org/t/avoiding-allocations-of-small-but-non-trivial-arrays-work-array-alternative/90084/28 "2022-11-12T11:46:40Z")

</div>

> [@Elrod](#):
>
> Try replacing `1:N` with `eachindex(barw)`.

```julia
julia> @btime run_Interp1d_MVector2(4)
  515.025 μs (16 allocations: 800 bytes)

```

Thus, this is the fastest so far for me:

> **Code**
>
> ```julia
> using StaticArrays
> using LoopVectorization
> 
> function Interp1d_MVector2(xnodes::SVector{N,T}, fvals::SVector{N,T}, t) where {N,T}
> barw = zeros(MVector{N,T})
> @turbo for j ∈ eachindex(barw)
> barw[j] = one(T)
> for k ∈ 1:N
> barw[j] = ifelse(k == j , barw[j] , barw[j] * (xnodes[j] - xnodes[k]) )
> end
> barw[j] = inv(barw[j])
> end
> 
> numt, denomt = zero(T), zero(T)
> for j ∈ eachindex(xnodes)
> tdiff = t - xnodes[j]
> numt = numt + barw[j] / tdiff * fvals[j]
> denomt = denomt + barw[j] / tdiff
> if ( abs(tdiff) < 1e-14 )
> numt = fvals[j]
> denomt = one(T)
> break
> end
> end
> return numt / denomt
> end
> 
> function the_loop(xnodes, fvals, t)
> for i ∈ 1:100000
> u = Interp1d_MVector2(xnodes, fvals, t)
> end
> end
> 
> function run_Interp1d_MVector2(N)
> xnodes = SVector{N,Float64}(LinRange(0,1,4))
> fvals = xnodes.^4
> t = 0.3
> the_loop(xnodes, fvals, t)
> end
> 
> ```

ps: replacing `1:N` in the inner loop does not work, any particular reason for that?

---

<div class="post-metadata">

**Author:** ![Elrod](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/elrod/32/22461_2.png) [@Elrod](https://discourse.julialang.org/u/Elrod)\
**Post date:** [November 12, 2022, 11:49am UTC](https://discourse.julialang.org/t/avoiding-allocations-of-small-but-non-trivial-arrays-work-array-alternative/90084/29 "2022-11-12T11:49:52Z")

</div>

> [@lmiq](#):
>
> ```julia
> @turbo for j ∈ eachindex(barw)
> barw[j] = one(T)
> for k ∈ 1:N
> 
> ```

You can replace the second, too.

---

<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:** [November 12, 2022, 11:50am UTC](https://discourse.julialang.org/t/avoiding-allocations-of-small-but-non-trivial-arrays-work-array-alternative/90084/30 "2022-11-12T11:50:31Z")

</div>

Then I get:

```julia
ERROR: "Reduction not found."
Stacktrace:
  [1] reduce_to_onevecunroll

```

---

<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:** [November 12, 2022, 12:08pm UTC](https://discourse.julialang.org/t/avoiding-allocations-of-small-but-non-trivial-arrays-work-array-alternative/90084/31 "2022-11-12T12:08:59Z")

</div>

Curiously, if I actually return something from the loop (computing the sum of `u` in it), the function becomes zero-allocating:

```julia
using StaticArrays
using LoopVectorization

function Interp1d_MVector2(xnodes::SVector{N,T}, fvals::SVector{N,T}, t) where {N,T}
    barw = zeros(MVector{N,T})
    @turbo for j ∈ eachindex(barw)
        barw[j] = one(T)
        for k ∈ 1:N
            barw[j] = ifelse(k == j , barw[j] , barw[j] * (xnodes[j] - xnodes[k]) )
        end
        barw[j] = inv(barw[j])
    end

    numt, denomt = zero(T), zero(T)
    for j ∈ eachindex(xnodes)
        tdiff = t - xnodes[j]
        numt = numt + barw[j] / tdiff * fvals[j]
        denomt = denomt + barw[j] / tdiff
        if ( abs(tdiff) < 1e-14 )
            numt = fvals[j]
            denomt = one(T)
            break
        end
    end
    return numt / denomt
end

function the_loop(xnodes, fvals, t)
    u = zero(eltype(xnodes))
    for i ∈ 1:100000
        u += Interp1d_MVector2(xnodes, fvals, t)
    end
    return u
end

function run_Interp1d_MVector2(N)
    xnodes = SVector{N,Float64}(LinRange(0,1,N))
    fvals = xnodes.^4
    t = 0.3
    the_loop(xnodes, fvals, t)
end

```

```julia
julia> @btime run_Interp1d_MVector2(4)
  538.863 μs (0 allocations: 0 bytes)
7584.269662929191

```

This is also probably more realistic, otherwise the function does nothing overall.

---

<div class="post-metadata">

**Author:** ![lrnv](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/lrnv/32/19373_2.png) [@lrnv](https://discourse.julialang.org/u/lrnv)\
**Post date:** [November 12, 2022, 12:38pm UTC](https://discourse.julialang.org/t/avoiding-allocations-of-small-but-non-trivial-arrays-work-array-alternative/90084/32 "2022-11-12T12:38:02Z")

</div>

Thanks for the `ifelse`, this is perfect. Moving the computation out of the `ifelse` statement ensures it is done in every case and allows `@turbo` to perform a little better:

```julia
using LoopVectorization
using StaticArrays
using BenchmarkTools

function v2(x,fvals, t) 
    @inbounds for j ∈ eachindex(x)
        tdiff = t - x[j]
        if abs(tdiff) < 1e-14
            return fvals[j]
        end
    end

    nₜ = 0.0
    dₜ = 0.0
    @turbo for j ∈ eachindex(x)
        xⱼ = x[j]
        bⱼ = t - xⱼ
        for k ∈ eachindex(x)
            rₖ = xⱼ - x[k]
            bⱼ = bⱼ * ifelse(k==j, 1.0, rₖ)
        end
        ibⱼ = inv(bⱼ)
        nₜ += ibⱼ * fvals[j]
        dₜ += ibⱼ
    end
    return nₜ / dₜ
end
function run_v2()
    xnodes = SVector{4,Float64}(LinRange(0,1,4))
    fvals = xnodes.^4
    t = 0.3
    u = 0.0
    for i ∈ 1:100000
        u += v2(xnodes,fvals, t)
    end
    return u
end

```

Giving :

```julia
@btime run_v2() # 0.498 ms, 0 allocation

```

---

<div class="post-metadata">

**Author:** ![Elrod](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/elrod/32/22461_2.png) [@Elrod](https://discourse.julialang.org/u/Elrod)\
**Post date:** [November 12, 2022, 1:13pm UTC](https://discourse.julialang.org/t/avoiding-allocations-of-small-but-non-trivial-arrays-work-array-alternative/90084/33 "2022-11-12T13:13:13Z")

</div>

=( Feel free to file an issue, but I hope to have LoopModels working as a replacement by JuliaCon 2023, so I’m unlikely to spend any time on LoopVectorization.

---

<div class="post-metadata">

**Author:** ![martin.d.maas](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/martin.d.maas/32/50964_2.png) [@martin.d.maas](https://discourse.julialang.org/u/martin.d.maas)\
**Post date:** [November 14, 2022, 6:36am UTC](https://discourse.julialang.org/t/avoiding-allocations-of-small-but-non-trivial-arrays-work-array-alternative/90084/36 "2022-11-14T06:36:48Z")

</div>

> [@lrnv](#):
>
> Thanks for the `ifelse`, this is perfect. Moving the computation out of the `ifelse` statement ensures it is done in every case and allows `@turbo` to perform a little better:

Oh, wow, this is fast!

Casting a `@view` as a StaticArray ended up being completely straightforward (given that the size is known at compile-time)… I just do this:

```julia
xnodes = SVector{N}(@view(A[i:i+N-1]))

```

where if I don’t, and just pass the @view, the function is also non-allocating (because that was solved by unrwaping the vector `barw`, but the same holds for the `MVector` version), but there is a penalty of about 2X in performance.

```julia
function run_v2_view()
    A = LinRange(0,1,1_000_000)
    F = A.^4
    t = 0.3
    order = 6
    for i ∈ 1:10000
        xnodes = @view(A[i:i+order-1])
        fvals = @view(F[i:i+order-1])
        u = v2(xnodes,fvals, t)
    end
end

function run_v2_static_view()
    A = LinRange(0,1,1_000_000)
    F = A.^4
    t = 0.3
    order = 6
    for i ∈ 1:10000
        xnodes = SVector{order}(@view(A[i:i+order-1]))
        fvals = SVector{order}(@view(F[i:i+order-1]))
        u = v2(xnodes,fvals, t)
    end
end

@btime run_v2_view()
11.757 ms (2 allocations: 7.63 MiB)

@btime run_v2_static_view()
6.128 ms (2 allocations: 7.63 MiB)

```

(the allocations are outside of the hot loop, so no problem with that)

---

<div class="post-metadata">

**Author:** ![martin.d.maas](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/martin.d.maas/32/50964_2.png) [@martin.d.maas](https://discourse.julialang.org/u/martin.d.maas)\
**Post date:** [November 14, 2022, 6:53am UTC](https://discourse.julialang.org/t/avoiding-allocations-of-small-but-non-trivial-arrays-work-array-alternative/90084/37 "2022-11-14T06:53:33Z")

</div>

> [@martin.d.maas](#):
>
> `SVector{order}(@view(A[i:i+order-1]))`

Hmmm, I don’t understand, I tested this in my real code, and these kinds of statements do seem to be triggering allocations… but on another level.

However, the new code with SIMD working and the original allocation problem solved seems to be the fastest by a good margin!

---

<div class="post-metadata">

**Author:** ![martin.d.maas](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/martin.d.maas/32/50964_2.png) [@martin.d.maas](https://discourse.julialang.org/u/martin.d.maas)\
**Post date:** [November 14, 2022, 6:56am UTC](https://discourse.julialang.org/t/avoiding-allocations-of-small-but-non-trivial-arrays-work-array-alternative/90084/38 "2022-11-14T06:56:28Z")

</div>

Ahh, this changes when using Threads:

```julia
function run_v2_view()
    A = LinRange(0,1,1_000_000)
    F = A.^4
    t = 0.3
    order = 6
    Threads.@threads for i ∈ 1:10000
        xnodes = @view(A[i:i+order-1])
        fvals = @view(F[i:i+order-1])
        u = v2(xnodes,fvals, t)
    end
end

function run_v2_static_view()
    A = LinRange(0,1,1_000_000)
    F = A.^4
    t = 0.3
    order = 6
    Threads.@threads for i ∈ 1:10000
        xnodes = SVector{order}(@view(A[i:i+order-1]))
        fvals = SVector{order}(@view(F[i:i+order-1]))
        u = v2(xnodes,fvals, t)
    end
end

@btime run_v2_view()
11.486 ms (51 allocations: 7.63 MiB)

@btime run_v2_static_view()
16.378 ms (240051 allocations: 19.84 MiB)

```

---

<div class="post-metadata">

**Author:** ![Mason](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mason/32/2423_2.png) [@Mason](https://discourse.julialang.org/u/Mason)\
**Post date:** [November 15, 2022, 12:09am UTC](https://discourse.julialang.org/t/avoiding-allocations-of-small-but-non-trivial-arrays-work-array-alternative/90084/39 "2022-11-15T00:09:54Z")

</div>

Looks like the multithreading is causing a type instability since `order` isn’t being constant-propagated across the task boundary. I’d rewrite it like this:

```julia
function run_v2_static_view(::Val{order}) where {order}
    A = LinRange(0,1,1_000_000)
    F = A.^4
    t = 0.3
    Threads.@threads for i ∈ 1:10000
        xnodes = SVector{order}(@view(A[i:i+order-1]))
        fvals = SVector{order}(@view(F[i:i+order-1]))
        u = v2(xnodes,fvals, t)
    end
end

```

```julia
julia> @btime run_v2_static_view(Val(6))
  3.371 ms (39 allocations: 7.63 MiB)

```

---

<div class="post-metadata">

**Author:** ![DNF](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/dnf/32/10191_2.png) [@DNF](https://discourse.julialang.org/u/DNF)\
**Post date:** [November 15, 2022, 6:18am UTC](https://discourse.julialang.org/t/avoiding-allocations-of-small-but-non-trivial-arrays-work-array-alternative/90084/40 "2022-11-15T06:18:47Z")

</div>

> [@Mason](#):
>
> ```julia
> function run_v2_static_view(::Val{order}) where {order}
> A = LinRange(0,1,1_000_000)
> F = A.^4
> t = 0.3
> Threads.@threads for i ∈ 1:10000
> xnodes = SVector{order}(@view(A[i:i+order-1]))
> fvals = SVector{order}(@view(F[i:i+order-1]))
> u = v2(xnodes,fvals, t)
> end
> end
> 
> ```

What happens if you do

```julia
function run_v2_static_index()
    A = LinRange(0,1,1_000_000)
    F = A.^4
    t = 0.3
    inds = SVector{6}(0:5)
    Threads.@threads for i ∈ 1:10000
        xnodes = A[inds.+i]
        fvals = F[inds.+i]
        u = v2(xnodes,fvals, t)
    end
end

```

?

---

<div class="post-metadata">

**Author:** ![martin.d.maas](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/martin.d.maas/32/50964_2.png) [@martin.d.maas](https://discourse.julialang.org/u/martin.d.maas)\
**Post date:** [November 17, 2022, 12:52am UTC](https://discourse.julialang.org/t/avoiding-allocations-of-small-but-non-trivial-arrays-work-array-alternative/90084/41 "2022-11-17T00:52:01Z")

</div>

That works too, however it seems to be only slightly slower than the previous alternative

[Previous page](https://discourse.julialang.org/t/avoiding-allocations-of-small-but-non-trivial-arrays-work-array-alternative/90084.md?page=1)
