Modifying subentries in a Vector{SVector}, @set and reinterpret

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?

2 Likes

Accessors and Base.setindex doesn’t skip instantiation methods, and even if it’s as simple as repackaging fields into new, the compiler occasionally may still decide against tweaking bytes. After all, even if the source code looks very similar to mutation before macro expansion, that doesn’t give an iron-clad guarantee to the compiler that you really don’t need the previous value after the tweak. An unsafer option would probably be needed to do the same work as real mutation.

Technically you’re not because immutables lack the sharing property of mutables. In effect, you unsafely instantiated an SVector and unsafely mutated the outer Vector by skipping the associated methods. For these particular types and an isbits element, you can probably get away with that. FWIW, setindex!(v::MArray, val, i::Int) makes a pointer to its backend Tuple field if isbits and uses unsafe_store! to tweak an offset element instead of safely instantiating the new Tuple instance.

1 Like