Reinterpret nested array as concatenated arrays

The generic method for matrix multiplication works already in that case:

julia> using StaticArrays

julia> A = [SMatrix{2, 2}(1, 2, 3, 4) for i in 1:3, j in 1:3];

julia> v = [SVector(1, 2) for i in 1:3];

julia> A * v
3-element Vector{SVector{2, Int64}}:
 [21, 30]
 [21, 30]
 [21, 30]

The other approach (reshaping into 2n x 2n matrix and 2n vector) may be more efficient, because it calls an optimized BLAS routine (I haven’t benchmarked the difference). However, you would need to create new temp arrays of that shape (the memory layout is different). It would looks something like:

julia> to_vector(nested::AbstractVector{T}) where {T <: StaticVector} = vec(reinterpret(reshape, eltype(T), nested))
to_vector (generic function with 1 method)

julia> function to_matrix(nested::AbstractMatrix{T}) where T <: StaticMatrix
           sz_in, sz_out = size(T), size(nested)
           sz = map(*, sz_in, sz_out)
           m0 = reinterpret(reshape, eltype(T), nested)
           m1 = reshape(m0, sz_in..., sz_out...)
           m2 = permutedims(m1, (1, 3, 2, 4))
           return reshape(m2, sz)
       end

where the permutedims step is allocating, but I don’t think there is a way to avoid that.

1 Like