Hi,
Most of the time it makes conceptually more sense to store some data as a Vector{SVector{Float64, 3}}
v
as opposed to a Matrix{Float64}
M
of size
(3, n)
. However, this makes updating (sub)entries more difficult.
I. What would be the ‘canonical’ equivalent to M[3, 5] = 0
for v
? The best I can come up with is x = v[5]; v[5] = @set x[3] = 0
using Accessors.jl, or v[5] = setindex(v[5], 3, 0)
, but both options don’t seem particularly elegant. Unless the compiler really steps up, this will also be slower (see below).
II. In some sense there is essentially no difference between v
and M
as the data is stored in the same way. This also means we can easily convert between them: M = reinterpret(reshape, eltype(eltype(v)), v)
and v = reinterpret(SVector{3, Float64}, M)
. Now we could use this fact to perform the equivalent to M[3, 5] = 0
on v
via such a reinterpret
:
v = rand(SVector{3, Float64}, 1000)
println(v[5][3]) # e.g. 0.9433828830786585
reinterpret(reshape, eltype(eltype(v)), v)[3, 5] = 0
println(v[5][3]) # 0
But this feels wrong, as we are mutating an immutable SVector
. On the other hand, it seems more performant than the @set
or setindex
approach:
using Accessors, BenchmarkTools, StaticArrays
function fill_set!(v)
for i = eachindex(v)
for j = 1:length(eltype(v)) # (or get this from type-annotating v::AbstractVector{SVector{N, T}})
v[i] = setindex(v[i], j, j) # v[i][j] = j
end
end
end
function fill_reinterpret!(v)
M = reinterpret(reshape, eltype(eltype(v)), v)
for i = eachindex(v)
for j = 1:length(eltype(v))
M[j, i] = j
end
end
end
v = zeros(SVector{3, Float64}, 1000)
@btime fill_set!($v) # 812.791 ns (0 allocations: 0 bytes)
@btime fill_reinterpret!($v) # 687.333 ns (0 allocations: 0 bytes)
v = zeros(SVector{4, Float64}, 1000)
@btime fill_set!($v) # 1.420 μs (0 allocations: 0 bytes)
@btime fill_reinterpret!($v) # 320.482 ns (0 allocations: 0 bytes)
# Faster for 4 than 3 because of SIMD?
For completeness: updating all of v[i] in one go, fill! and (FixedSize)Matrix
using FixedSizeArrays
function fill_all!(v::AbstractVector{SVector{N, T}}) where {N, T}
for i = eachindex(v)
v[i] = SVector{N, T}(1:N)
end
end
function fill_matrix!(M)
for i = axes(M, 2)
for j = axes(M, 1)
M[j, i] = j
end
end
end
function fill_matrix_val!(M, ::Val{N}) where N
for i = axes(M, 2)
for j = 1:N
M[j, i] = j
end
end
end
v = zeros(SVector{3, Float64}, 1000)
M = zeros(3, 1000)
fM = FixedSizeArray(M)
@btime fill_all!($v); # 479.487 ns (0 allocations: 0 bytes)
@btime fill!($v, SVector{3, Float64}(1:3)); # 480.000 ns (0 allocations: 0 bytes)
@btime fill_matrix!($M) # 1.920 μs (0 allocations: 0 bytes)
@btime fill_matrix!($fM) # 1.810 μs (0 allocations: 0 bytes)
@btime fill_matrix_val!($M, Val(3)) # 495.876 ns (0 allocations: 0 bytes)
@btime fill_matrix_val!($fM, Val(3)) # 539.683 ns (0 allocations: 0 bytes)
v = zeros(SVector{4, Float64}, 1000)
M = zeros(4, 1000)
fM = FixedSizeArray(M)
@btime fill_all!($v); # 325.820 ns (0 allocations: 0 bytes)
@btime fill!($v, SVector{4, Float64}(1:4)); # 323.684 ns (0 allocations: 0 bytes)
@btime fill_matrix!($M) # 2.322 μs (0 allocations: 0 bytes)
@btime fill_matrix!($fM) # 2.233 μs (0 allocations: 0 bytes)
@btime fill_matrix_val!($M, Val(4)) # 391.089 ns (0 allocations: 0 bytes)
@btime fill_matrix_val!($fM, Val(4)) # 372.683 ns (0 allocations: 0 bytes)
Will using reinterpret
to mutate immutables in this way lead to certain problems?