Immutable field slicing

[c42f]

When using arrays of small immutable types, I commonly find myself wanting to slice along the "field name dimension". Say I have an array of boring old Vec3 type:

using FixedSizeArrays

immutable Vec3{T} <: FixedVectorNoTuple{3, T}
    x::T
    y::T
    z::T
end

v = [Vec3(i,j,0) for i = 1:4, j = 1:4]

Now I want all the x coordinates... darn, am I going to need a comprehension or something equally verbose and inflexible? (Possibly not, I'm a newcomer to FixedSizeArrays as of today.)

What I'd really like is to "slice along the x field dimension", thus: v[:x, :, :]. Here's an interesting proof of concept which seems to work (please excuse the hacky reshape junk, it's just a proof of concept):

function Base.getindex{T<:FixedVectorNoTuple}(a::Array{T}, fieldname::Symbol, inds...)
    ty = eltype(a)
    fieldindex = find(fieldnames(ty) .== fieldname)
    squeeze(reinterpret(fieldtype(ty, fieldname), a, (nfields(ty), size(a)...))[fieldindex, inds...], 1)
end

function Base.setindex!{T<:FixedVectorNoTuple}(a::Array{T}, value, fieldname::Symbol, inds...)
    ty = eltype(a)
    fieldindex = find(fieldnames(ty) .== fieldname)
    reinterpret(fieldtype(ty, fieldname), a, (nfields(ty), size(a)...))[fieldindex, inds...] = value
end

Now I can do things like

julia> v = [Vec3(i,j,0) for i = 1:4, j = 1:4]
4x4 Array{Vec3{Int64},2}:
 Vec3(1,1,0)  Vec3(1,2,0)  Vec3(1,3,0)  Vec3(1,4,0)
 Vec3(2,1,0)  Vec3(2,2,0)  Vec3(2,3,0)  Vec3(2,4,0)
 Vec3(3,1,0)  Vec3(3,2,0)  Vec3(3,3,0)  Vec3(3,4,0)
 Vec3(4,1,0)  Vec3(4,2,0)  Vec3(4,3,0)  Vec3(4,4,0)

julia> v[:x, :, :]
4x4 Array{Int64,2}:
 1  1  1  1
 2  2  2  2
 3  3  3  3
 4  4  4  4

julia> v[:z, 2:3, 2:3] = 100
100

julia> v
4x4 Array{Vec3{Int64},2}:
 Vec3(1,1,0)  Vec3(1,2,0)    Vec3(1,3,0)    Vec3(1,4,0)
 Vec3(2,1,0)  Vec3(2,2,100)  Vec3(2,3,100)  Vec3(2,4,0)
 Vec3(3,1,0)  Vec3(3,2,100)  Vec3(3,3,100)  Vec3(3,4,0)
 Vec3(4,1,0)  Vec3(4,2,0)    Vec3(4,3,0)    Vec3(4,4,0)

Good idea or terrible idea? It sure seems super handy.

[c42f]

@PaulBellette behold, field slicing crazyness. It's not quite the syntax v[2:3, 2:3].z which seems more intuitive at first glance, but it's a lot more implementable, and probably more self-consistent.

[SimonDanisch]

Ah yeah, I definitely want to have this as well. Especially, when we decide to use FixedSizeArrays for ColorTypes again.

Cool, this is faster than it looks :)
Didn't do extensive benchmarks, but without investing too much time, only explicitely doing [vec.x for i...] was faster.

[c42f]

Hah, please don't take the implementation above too seriously :-) Now that I know it makes sense to someone other than myself, I'm a bit keener to spend time making it actually work. Some thoughts/questions:

• Is the syntax above the right thing? I do like the slicing syntax for mutating via setindex(), since it emphasizes the array (as the mutable container for immutable types) is being mutated, rather than the individual immutable values themselves. I guess this is arguably a bit nasty, but it's so incredibly handy that I want some convenient syntax.
• What can we do for slicing along the indexed dimensions of a FixedVector? The above syntax only works nicely for a FixedVectorNoTuple since typeof(:1) is Int rather than Symbol. It would seem sensible to use Val{1} in which case we'd be able to make the implementation more efficient by dispatching on the index type, but syntax like v[Val{1}, :, :] is less compelling. (IIRC there was some talk of making {...} return a type tuple - that might work acceptably v[{1}, :, :]. Ah, here was the issue what should we use { } for? JuliaLang/julia#8470)
• Could FixedVector slicing be made harmonious with the usual rules for "fancy indexing" - indexing with arrays or BitVectors etc? If the type tuple syntax is accepted we might be able to have v[{1,2}, :, :] which would be pretty neat.

[c42f]

Taking the Val approach for now, the following versions work for FixedVector:

function Base.getindex{V<:FixedVector,I}(a::Array{V}, ::Type{Val{I}}, inds...)
    squeeze(reinterpret(eltype(V), a, (length(V), size(a)...))[I, inds...], 1)
end

function Base.setindex!{V<:FixedVector,I}(a::Array{V}, rhs, ::Type{Val{I}}, inds...)
    reinterpret(eltype(V), a, (length(V), size(a)...))[I, inds...] = rhs
end

Example usage:

julia> v = [Vec(i,j,0) for i = 1:4, j = 1:4]
4x4 Array{FixedSizeArrays.Vec{3,Int64},2}:
 Vec(1,1,0)  Vec(1,2,0)  Vec(1,3,0)  Vec(1,4,0)
 Vec(2,1,0)  Vec(2,2,0)  Vec(2,3,0)  Vec(2,4,0)
 Vec(3,1,0)  Vec(3,2,0)  Vec(3,3,0)  Vec(3,4,0)
 Vec(4,1,0)  Vec(4,2,0)  Vec(4,3,0)  Vec(4,4,0)

julia> v[Val{1}, :, :] = 1
1

julia> v
4x4 Array{FixedSizeArrays.Vec{3,Int64},2}:
 Vec(1,1,0)  Vec(1,2,0)  Vec(1,3,0)  Vec(1,4,0)
 Vec(1,1,0)  Vec(1,2,0)  Vec(1,3,0)  Vec(1,4,0)
 Vec(1,1,0)  Vec(1,2,0)  Vec(1,3,0)  Vec(1,4,0)
 Vec(1,1,0)  Vec(1,2,0)  Vec(1,3,0)  Vec(1,4,0)

[SimonDanisch]

Pretty cool. Yeah Val was one option, using types would be another.

using FixedSizeArrays
abstract Unit{T}
immutable X{T} <: Unit{T}
val::T
end
immutable Y{T} <: Unit{T}
val::T
end
immutable Z{T} <: Unit{T}
val::T
end
type2index{T <: X}(::Type{T}) = 1
type2index{T <: Y}(::Type{T}) = 2
type2index{T <: Z}(::Type{T}) = 3
function Base.getindex{V<:FixedVector,I<:Unit}(a::Array{V}, ::Type{I}, inds...)
    squeeze(reinterpret(eltype(V), a, (length(V), size(a)...))[type2index(I), inds...], 1)
end
function Base.setindex!{V<:FixedVector,I<:Unit}(a::Array{V}, rhs::I, inds...)
    reinterpret(eltype(V), a, (length(V), size(a)...))[type2index(I), inds...] = rhs.val
end
v = [Vec(i,j,0) for i = 1:4, j = 1:4]
v[:, :] = X(1)
4x4 Array{FixedSizeArrays.Vec{3,Int64},2}:
 Vec(1,1,0)  Vec(1,2,0)  Vec(1,3,0)  Vec(1,4,0)
 Vec(1,1,0)  Vec(1,2,0)  Vec(1,3,0)  Vec(1,4,0)
 Vec(1,1,0)  Vec(1,2,0)  Vec(1,3,0)  Vec(1,4,0)
 Vec(1,1,0)  Vec(1,2,0)  Vec(1,3,0)  Vec(1,4,0)

v[X, :, :]
4x4 Array{Int64,2}:
 1  1  1  1
 1  1  1  1
 1  1  1  1
 1  1  1  1

The last case should return Array{X{Int}}, I was just too lazy to handle the cases with different eltypes. We could also allow v[X{Float32}, :,:] which would return an Array{X{Float32}}.
This would be nice for visualizations, since I have more information (red channel -> red image, array with Z values --> heightfield).

[c42f]

I had a bit more of a think about the right way to approach this today. With respect to FixedVector, all the above options involving special types (including my suggestions of using Val and Tuple) to distinguish between indexing the fixed vs variable dimensions seem pretty hacky to me.

The general question I'd really like a good answer to, is how to consistently slice an Array{FixedArray{T, N, SIZE},M}. This requires N+M indices/ranges, with the first N indices being fixed size and the last M being variable. Perhaps a careful implementation of getindex with the usual index semantics and types would sort this out. Might need to be @generated, I haven't had time to prototype something yet.

Or maybe the idea of slicing arrays of immutable types to modify them is just inconsistent (they're meant to be immutable after all, they're not value types even though I keep thinking of them that way...) Is there some need for some kind of hybrid ArrayWithLeadingFixedDimensions container (obviously with a better name)?

[c42f]

I just came across JuliaLang/julia#11902, which is obviously relevant to the design problems here. In particular, Keno's suggestions around a new setindex address the ickyness of changing immutables in place.

[SimonDanisch]

Yup, that's one of the most wanted features for FixedSizeArrays! ;)
This will also make most of the @generated functions unnecessary and it will be a lot easier to inherit from AbstractArray!

[c42f]

Right, makes sense. I tried implementing a getindex for arrays of Mat which takes normal integers as indices, but immediately ran into conflicts with the default getindex implementation, as you might expect.

Unfortunately I think there may be a need for some syntax to separate the fixed dimensions from the usual array dimensions. Using a macro I managed to abuse a semicolon for this:

a = [Mat(((i,i), (j,j))) for i=1:2, j=1:2]

@fslice a[:,:; 1,1] = 10

After some thought I quite like the semicolon for this, other than its visual similarity to the colon. If slicing of general immutables became a standard feature, maybe the parser would need some syntactical way of figuring out when to lower the ref syntax to setindex! vs the immutable setindexversion proposed in JuliaLang/julia#11902

Implementation:

# Place holder type for dispatch to function which slices along dimensions of
# fixed size elements of an Array.
immutable SliceElements; end

# Turn A[1,2; 3,4] into A[SliceElements(), 1,2,3,4]
function wrapinds(expr)
    if expr.head == :(=)
        return Expr(expr.head, wrapinds(expr.args[1]), expr.args[2:end]...)
    end
    @assert expr.head == :typed_vcat
    @assert expr.args[2].head == :parameters # what?
    extrainds = expr.args[3:end]
    normalinds = expr.args[2].args
    Expr(:ref, expr.args[1], SliceElements(), extrainds..., normalinds...)
end

"""
Slice along dimensions of fixed sized arrays inside an array container, using
the syntax

    @fslice a[i,j; k,l]

where a is an array of `F<:FixedArray`.  i and j slice the dimensions of F; k
and l slice along the usual dimensions of the container a.
"""
macro fslice(expr)
    wrapinds(expr)
end


function Base.getindex{F<:FixedArray}(a::Array{F}, ::SliceElements, inds...)
    SIZE = size(F)
    T = eltype(F)
    reinterpret(T, a, (SIZE..., size(a)...))[inds...]
end

function Base.setindex!{F<:FixedArray}(a::Array{F}, val, ::SliceElements, inds...)
    SIZE = size(F)
    T = eltype(F)
    reinterpret(T, a, (SIZE..., size(a)...))[inds...] = val
end

[c42f]

Darn, fslice above doesn't work for fixed vectors due to the fact that I'm abusing the parsing of typed_vcat. Actually new syntax for distinguishing the fixed dimensions would only be needed(?) in the context of special setindex lowering in the 11902 julep, it's a bit useless here.

[c42f]

Ok, I think I'm converging on something useful - I'll prepare a proper PR for it, hopefully tonight.

[c42f]

In #62 I tried an @fslice macro as syntax for extending slicing across fixed dimensions.

It turns out that it's possible without a macro, but getting the dispatch working correctly seems to require some nasty hacks. Basically I'd like to get julia's dispatch system to call my special purpose implementation of getindex{F<:FixedArray,N}(a::Array{F,N}, inds...) only when ndims(F) + N == length(inds). I'd love to be wrong, but I don't think this constraint is possible at the language level. It can be managed with the following thunk but it's rather disgusting and depends on Base implementation details. I also ran into some unfortunate method ambiguity warnings which must be silenced by creating several additional redundant methods.

@generated function getindex_impl{F<:FixedArray,N}(a::Array{F,N}, inds...)
    if length(inds) == ndims(F) + N
        quote
            destructure(a)[inds...]
        end
    else
        quote
            # Hack: Need to call Base implementation detail to avoid infinite
            # recursion of getindex(::Array{F<:FixedArray}) calling itself.
            Base._getindex(linearindexing(a), a, inds...)
        end
    end
end

# Hack: several definitions to avoid ambiguity warnings vs underlying Array implementation
Base.getindex{F<:FixedArray}(a::Array{F}, ind::Real) = getindex_impl(a, ind)
Base.getindex{F<:FixedArray}(a::Array{F}, ind::UnitRange{Int64}) = getindex_impl(a, ind)
Base.getindex{F<:FixedArray}(a::Array{F}, ind::Colon) = getindex_impl(a, ind)
Base.getindex{F<:FixedArray,T<:Real}(a::Array{F}, ind::Range{T}) = getindex_impl(a, ind)

Base.getindex{F<:FixedArray}(a::Array{F}, inds...) = getindex_impl(a, inds...)

[c42f]

Hmm, I'm still not fully satisfied by any of the above (particularly the reshape() stuff which is far from a generic solution). I think I'll try an alternative approach: implement destructure() as a function which wraps the provided array inside a DestructuredArray type, and see where that leads.

[c42f]

I guess #62 addresses this well enough, we'll see :)