Best way to handle a list of arrays (Asking for recommendations)

lets say i have a list of coordinates in vectors eg

v1 = [1,2,3]
v2 = [4,5,6]
v3 = [7,8,9]

(list is much longer)

I can imagine two ways to organise these:

  1. Put them in a large matrix as done for example by eigvecs,

    V = [1 4 7 …;
         2 5 8 …;
         3 6 9 …]
    
    • Advantages: Easy to grab a vector of one coordinate i by v[i,:], probably a better use of memory since everything is on one place in one piece
    • Disadvantage: More difficult to create a Vector of a function of the vectors, for example the norm of all vectors
  2. Put them in a Vector of Vectors as done for example by Solvers of DifferentialEquations,

    V = [[1,2,3],
         [4,5,6],
         [7,8,9]]
    
    • Advantage: Easy to do a norm of all vectors by norm.(V)
    • Disadvantages: More difficult to grab a vector of one coordinate, maybe not continuous in memory and thus slower (is that true?)

So is there a general recommendation/best pattern?

I would highly recommend using a Vector of SVectors (statically-sized immutable vectors from StaticArrays.jl). This gives you basically the best of both worlds:

julia> using StaticArrays

julia> v1 = SVector(1,2,3); v2 = SVector(4,5,6); v3 = SVector(7,8,9)
3-element SArray{Tuple{3},Int64,1,3}:
 7
 8
 9

julia> V = [v1, v2, v3]
3-element Array{SArray{Tuple{3},Int64,1,3},1}:
 [1, 2, 3]
 [4, 5, 6]
 [7, 8, 9]

Accessing each element is easy:

julia> V[2]
3-element SArray{Tuple{3},Int64,1,3}:
 4
 5
 6

Applying a function to each vector is easy:

julia> using LinearAlgebra

julia> norm.(V)
3-element Array{Float64,1}:
  3.7416573867739413                                                                                                                                                                                                                                                                                                                  
  8.774964387392123                                                                                                                                                                                                                                                                                                                   
 13.92838827718412 

But the vectors are actually stored contiguously in memory (Julia guarantees this for bitstypes), so lookup is extremely efficient. You can even reinterpret the vector-of-svectors as a matrix with no copying:

julia> reshape(reinterpret(Int, V), (3, :))
3×3 reshape(reinterpret(Int64, ::Array{SArray{Tuple{3},Int64,1,3},1}), 3, 3) with eltype Int64:
 1  4  7
 2  5  8
 3  6  9

And, even better, constructing an SVector is essentially free:

julia> @btime SVector(1,2,3)
  1.240 ns (0 allocations: 0 bytes)
7 Likes

Oh, and one more advantage: it’s extremely clear to users of your code (or to future you) how to access each element (V[2] for the second vector, for example). If you use a matrix for storage, you will forever need to remember whether you need to do V[2, :] or V[:, 2]

2 Likes

Thank you
The tip of using StaticArrays is nice, it at least solves the problem of memory usage

EDIT:
I only noticed it now the reinterpret thing also gives you also the possibility to extract the vectors of components easily. that’s nice!