Row vs. column major for arrays of coordinate vectors

General Usage Performance
1

There is one advantage that I can think of with column major storage for geometry (coordinate storage). 1D, 2D and 3D coordinates can be represented without hassle. A row-major format will require variable stride (1 in 1D, 2 in 2D and 3 in 3D) to access the coordinate value. But it is always fixed (stride of 1) in the case of column major formats.

2

No, this is the same for both formats — it simply depends on whether the coordinate vectors are stored in the rows or columns.

In any case, in Julia I would typically recommend using StaticArrays.jl for 1d/2d/3d coordinate vectors, in which case you store an array of coordinates as a 1d array of SVectors instead of a 2d array.

3

Is this not row vs. column major?

x_{i,j} = \{x_{i}\}_j = x(i,j) (i=dim, j = #of points)

May not be optimal in many situations (the AoS vs. SoA debate)

4

Why not just tuples?

1 Like
5

StaticArrays are wrappers around Tuples that have extra array-like behavior defined for them, eg. vector addition, matrix multiplication etc. For coordinates, to take an example, this makes it easy to transform between reference frames.

2 Likes
6

You could, but you would need to add your operations, e.g. addition, as tuples do not have it defined. It could work, but you might/likely will leave out operations you do now know about (ok if you don’t need them, until you do).

julia> supertype(SVector)
StaticVector (alias for StaticArray{Tuple{S}, T, 1} where {S, T})

julia> supertype(StaticArray)
AbstractArray

vs.

julia> supertype(NTuple)  # also for Tuple
Any

It’s useful for all Arrays and Vectors to derive from AbstractArray (rather than only Any). I trust what the operations you need are as fast as possible with StaticArrays.jl, you wouldn’t do better, the only slight downside to defining more operations in the package than you think you need is startup-overhead:

julia> @time using StaticArrays
  0.957173 seconds (2.18 M allocations: 156.695 MiB, 1.31% compilation time)  # in 1.8, down from 2.05 sec in 1.6.0

It seems it could do better with startup overhead, and there are ways for that.

7

No. You can choose to put your 3d coordinate vectors in A[i,:] (e.g. A is N\times3) or A[:,i] (3\times N). Row vs. column major is about how the resulting arrays are stored in memory. So, if you want the components to be consecutive in row-major, you use A[i,:], whereas if you want them to be consecutive in column-major, you use A[:,i].

Typical applications using 1d/2d/3d coordinate vectors access multiple components together frequently (e.g. for performing vector operations), in which case an array of SVector (AoS) is a good choice.

8

With the latest release of StructArrays, you can have have a StructArray{SVector} which behaves like Array{SVector} (v[i] gives you the i-th SVector) but is stored as a SoA. If you need both to access the “component vectors” as well as the individual SVectors in an efficient way, it can probably be a good compromise.

Just realized the above sounds a bit abstract. More concretely you can do things like:

julia> using StructArrays, StaticArrays, BenchmarkTools

julia> s = StructArray(SVector(1, 2) for _ in 1:10);

julia> @btime StructArrays.components($s)
  1.300 ns (0 allocations: 0 bytes)
([1, 1, 1, 1, 1, 1, 1, 1, 1, 1], [2, 2, 2, 2, 2, 2, 2, 2, 2, 2])

julia> @btime $s[3]
  1.600 ns (0 allocations: 0 bytes)
2-element SVector{2, Int64} with indices SOneTo(2):
 1
 2
5 Likes