[jpivarski]

> I don't get what's missing? Seems like that's the list in the documentation.

The essential part of this is the fact that it is a columnar implementation. Although that's not interface-level, it's not a detail: the time complexities and degree of memory sharing between input and output are different for these algorithms. The video link above is fast-forwarded to the 3.5 minute section that describes this in detail. (I've talked about it in about a dozen talks, but that one is a pretty good presentation.)

Keep in mind, I think Julia is great! The composability and generality of the array abstract type make short work of a wide variety of non-columnar algorithms. I wish we were using it more in high energy physics. As I said above, I'm doubly interested if someone is thinking about doing columnar algorithms in Julia. (The multiple dispatch would make it a bit different: you probably wouldn't want a high-level wrapper like ak.Array, just a uniform interface on all the node types, since they don't have class methods to expose internals.)

The columnar approach has two main advantages, only one of which is to make up for Python's necessarily slow implementation. The other is that it improves opportunities for vectorization, the hardware equivalent of what we're doing for Python. Doing it in Julia would not be superfluous!

[simeonschaub]

In Julia you would just use an array of views for this, which can represent slices of arrays without any copies:

julia> struct Muon p_T::Float64; phi::Float64; eta::Float64; end

julia> a = reinterpret(Muon, rand(3*7))

7-element reinterpret(Muon, ::Vector{Float64}):

Muon(0.5512393381972832, 0.9349789151451744, 0.006690464595502932) Muon(0.5856015732294971, 0.19023473269375601, 0.40764209748521973) Muon(0.14872954753560852, 0.12281085717049867, 0.9307398048388644) Muon(0.7885776521084014, 0.1392696530731592, 0.4054805743644767) Muon(0.4841152655677211, 0.053858886714772236, 0.9556610184833677) Muon(0.5325190758093583, 0.31100637434877343, 0.4364100043728055) Muon(0.8697751162452897, 0.07683143115108726, 0.49822326551511953)

julia> jagged = [view(a, idx) for idx in [1:3, 4:4, 5:5, 6:7]]

4-element Vector{SubArray{Muon, 1, Base.ReinterpretArray{Muon, 1, Float64, Vector{Float64}, false}, Tuple{UnitRange{Int64}}, true}}:

[Muon(0.5512393381972832, 0.9349789151451744, 0.006690464595502932), Muon(0.5856015732294971, 0.19023473269375601, 0.40764209748521973), Muon(0.14872954753560852, 0.12281085717049867, 0.9307398048388644)] [Muon(0.7885776521084014, 0.1392696530731592, 0.4054805743644767)] [Muon(0.4841152655677211, 0.053858886714772236, 0.9556610184833677)] [Muon(0.5325190758093583, 0.31100637434877343, 0.4364100043728055), Muon(0.8697751162452897, 0.07683143115108726, 0.49822326551511953)]

Also note how I was able to tell Julia to just reinterpret a bunch of contiguous floating point values as objects of type `Muon`, which produced a `ReinterpretArray`. Nowhere in there did I ever copy any data from the original array produced by the `rand(3*7)` call.*

[ChrisRackauckas]

As mentioned in my post, ComponentArrays.jl will represent it columnarly and uses that for doing this efficiently with GPUs and automatic differentiation.