# Filling a preallocated array without allocations

Hello, does anybody know if it is possible to fill a preallocated array without allocations?

In other words, is there a way of performing this function `loop` without allocations

``````mutable struct KF
F::Array{Float64, 3}
function KF()
new(rand(2, 2, 100))
end
end

function loop(kfaux::KF)
@inbounds for t in 40:100
kfaux.F[:, :, t] = kfaux.F[:, :, t-1]
end
return
end

kfaux = KF()
using BenchmarkTools
@benchmark loop(\$kfaux)
``````

I got

``````julia> @benchmark loop(\$kfaux)
BenchmarkTools.Trial:
memory estimate:  6.67 KiB
allocs estimate:  61
--------------
minimum time:     3.361 ÎĽs (0.00% GC)
median time:      3.467 ÎĽs (0.00% GC)
mean time:        4.745 ÎĽs (15.47% GC)
maximum time:     860.028 ÎĽs (98.64% GC)
--------------
samples:          10000
evals/sample:     8
``````

Sureâ€“one way to do it is just to write out the loop:

``````julia> function loop(kfaux::KF)
for t in 40:100
for j in axes(kfaux.F, 2)
for i in axes(kfaux.F, 1)
kfaux.F[i, j, j] = kfaux.F[i, j, t - 1]
end
end
end
end
loop (generic function with 1 method)

julia> @btime loop(\$kfaux)
496.222 ns (0 allocations: 0 bytes)
``````
4 Likes

In addition, if the dimensions of your example is representative of what you are really trying to do, maybe it makes sense to swap `Array{T, 3}` to `Vector{SMatrix}` using StaticArrays

``````using StaticArrays

mutable struct KFS
F::Vector{SMatrix{2,2,Float64,4}}
function KFS()
new([@SMatrix rand(2,2) for _ in 1:100])
end
end

function loop(kfaux::KFS)
F = kfaux.F
@inbounds for t in 40:100
F[t] = F[t-1]
end
return
end
``````

Benchmarks:

``````julia> using BenchmarkTools

julia> @btime loop(\$kfaux) # Original version
2.503 ÎĽs (61 allocations: 6.67 KiB)

julia> @btime loopâ€˛(\$kfaux) # Robin's version
459.447 ns (0 allocations: 0 bytes)

julia> @btime loop(\$kfauxs)
18.602 ns (0 allocations: 0 bytes)
``````
5 Likes

Can you explain why this is allocation less? Would like to learn

Kind regards

Sure!

When you write `x[:, :, t]`, you are making a copy of that slice of the array, so this line:

``````x[:, :, t] = x[:, :, t - 1]
``````

makes a copy of `x[:, :, t - 1]` and then assigns that to `x[:, :, t]`.

You can use a view instead:

``````x[:, :, t] .= @view x[:, :, t - 1]
``````

but it does not make a huge difference in your particular case because the sizes of the first two dimensions of your array are so small.

A loop doesnâ€™t require creating a copy or a view at all, so itâ€™s a nice solution. I would also second @fredrikekreâ€™s suggestion: if your data structure is actually a collection of N 2x2 matrices, then itâ€™s much nicer to work with it in that way, and using StaticArrays makes that very efficient.

2 Likes

Thanks, so loops are actually pretty smart, nice to know. Used to that vectorizing would be faster, so that is why I sometimes have a hard time grasping these kind of things, but I really liked your explanation.

Kind regards

â€śVectorizingâ€ť in that sense is only necessary in languages where loops are slow (matlab, python, etc.). After all, a vectorized operation in numpy or matlab is actually implemented via a loop in a fast language (like C). Julia is a fast language, so you can choose a vectorized operation or a loop based on whichever one makes your particular task easier to write, read, or understand.

6 Likes

Thanks! It helped a lot!