How to elegantly express the same operation for permutations of x, y and z?

General Usage
1

This is a situation I face very often: I have the same operation that needs to be run once for each dimension of a 3D construct.

This is complicated because the order of indices differs for each variant, and I don’t know how to express this different ordering elegantly, causing me to duplicate or triplicate lots of code.

Here’s one example: Let’s say the task is to flip two dimensions of a point, the dimension staying intact is given as an integer. This would be a solution like the ones I usually write, which are annoying me because they are so redundant and need control flow:

function flip_parts(point, dim::Int)
    if dim == 1
        (point[1], point[3], point[2])
    elseif dim == 2
        (point[3], point[2], point[1])
    elseif dim == 3
        (point[2], point[1], point[3])
    else
        error()
    end
end

I don’t know how to generically express the different ordering of elements to build the return tuples given only the dim input argument.

Who can point me to an elegant way of solving problems like this? Ideally in the most generic and performant way possible?

2

I have not thought through the whole problem, but one way of finding which the two dimensions that should be flipped are, given the fixed dimension d is using mod1:

julia> for d = 1:3
           @show d, mod1(d-1,3), mod1(d+1,3)
       end
(d, mod1(d - 1, 3), mod1(d + 1, 3)) = (1, 3, 2)
(d, mod1(d - 1, 3), mod1(d + 1, 3)) = (2, 1, 3)
(d, mod1(d - 1, 3), mod1(d + 1, 3)) = (3, 2, 1)

EDIT: It appears this snippet gives you the ordering you want:

julia> for d = 1:3
           @show circshift([d, mod1(d-1, 3), mod1(d+1, 3)], d-1)
       end
circshift([d, mod1(d - 1, 3), mod1(d + 1, 3)], d - 1) = [1, 3, 2]
circshift([d, mod1(d - 1, 3), mod1(d + 1, 3)], d - 1) = [3, 2, 1]
circshift([d, mod1(d - 1, 3), mod1(d + 1, 3)], d - 1) = [2, 1, 3]
2 Likes
3

Looks like the orderings are just the cyclic permutations of (3,2,1). So maybe circshift is indeed the way to go:

julia> circshift.(Ref(reverse(1:3)), 1:3)
3-element Array{Array{Int64,1},1}:
 [1, 3, 2]
 [2, 1, 3]
 [3, 2, 1]

You can do the operation in place with circshift! to limit memory allocation. I don’t know if this can be better than the branching though.

4

Maybe to get away from this particular problem, because it’s not about just this one, let me rephrase the question:

What is the best way to express an algorithm, in which the order of elements in something like a tuple, svector, etc. needs to be determined at runtime? Can this be done without allocating?

5

I’ve found TupleTools.jl to be useful for type stable tuple manipulations like sorting, maybe it can help here too.

1 Like
6

I’m not sure that there is a one-size-fits-all solution without seeing more examples, but this appears to be allocation-free, and rather quick:

function test_permutation(d)
    a,b,c = ((1,3,2), (2,1,3), (3,2,1))[d]
end

using BenchmarkTools
@eval @benchmark test_permutation($(rand(1:3)))

BenchmarkTools.Trial:
  memory estimate:  0 bytes
  allocs estimate:  0
  --------------
  minimum time:     0.040 ns (0.00% GC)
  median time:      0.043 ns (0.00% GC)
  mean time:        0.043 ns (0.00% GC)
  maximum time:     0.220 ns (0.00% GC)
  --------------
  samples:          10000
  evals/sample:     1000
1 Like
7

You need a Ref there. But still this solution takes “no time at all” and doesn’t allocate as you said.

julia> @btime test_permutation($(Ref(2))[])
  1.970 ns (0 allocations: 0 bytes)
8

map and ntuple are very helpful when dealing with Tuples, so (complementing the reply above) things like

perm = ((1,3,2), (2,1,3), (3,2,1))[dim]
map(i -> point[i], perm)

should be fast and non-allocating.

Alternatively, if you have a function f (or a matrix mat) such as f(i, dim) (resp. mat[i, dim]) returns the new index, you could also use ntuple:

ntuple(3) do i
    idx = f(i, dim)
    return point[idx]
end
1 Like
9

Ah the map and ntuple path could be what I was looking for!

10

I would just define tiny functions for the permutation, eg

dim2perm(t) = (t[3], t[2], t[1])

and pass these around as arguments to the relevant methods.