I am trying to permute/shuffle columns of a multidimensional array, an example is here:

using BenchmarkTools
const arra = zeros(Int, 4, 3, 2)
arra[:,:, 1] = 1:12
arra[:,:, 2] = 13:24
shuffle_col = [3, 1, 1]
@benchmark arra[:, shuffle_col, :] #solution 1 - works but allocates
@benchmark @view arra[:, shuffle_col, :] #solution 2 - much faster but arra is now a subarray

Both solutions work but the first one allocates, and the second one forces me to continue to work with a subarray (which I would like to avoid). I tried to put everything in a function in a loop, but my solution takes longer than the two above for bigger arrays (and needs a copy to work to begin with):

function permute_col_array!(array::AbstractArray, permutation::Vector{Int}) #, dimension::Int
temp = copy(array)
for (i, permute) in enumerate(permutation)
array[:,i,:] = temp[:,permute,:]
end
end
@benchmark permute_col_array!(arra, shuffle_col)

Does anyone here know if I can permute multidimensional arrays without allocation? I would gain significant improvements in functions if I could do that. There is a base function permute! in Julia but it only works on vector elements, not on whole columns.

I’ve been stuck trying to solve this problem too! I have a use case for this kind of mechanism that I haven’t been able to crack: Complex Numbers Wrappers

AFAIK permuting arrays in place in general is not a trivial problem. If allocation (and not copying) is the costly part per se, perhaps you should just preallocate for output.

Preallocation is unfortunately not possible as I need to shuffle the array in a loop and the object of interest is the array itself ( a filtering problem).

Could you elaborate on this ? In place permutations are really difficult to optimize. In 3D you have to compute some permutation cycles… In general out of place (with pre-allocation + ptr swap) permutations are much easier to implement and perform better. If it is simple permutedims you can tile.

#Input
observations = Matrix(T,P)
#Algorithm
function myproblem(observations, k, ....)
#initialize trajectory
trajectories = Matrix(T,M,N)
for i in 1:T
#get column order
weights = somefunction(observations[1:i,:], trajectories[:,(i-k):i,:] )
#Reshuffle trajectories
trajectories = trajectories[:,weights,:]
end
return trajectories
end

I can do preallocation for trajectories, but so far I have not found a way for trajectories[:,weights,:] to not allocate a new array at each iteration unfortunately. I also cannot simply make a output container as I only know the output after looping through all columns in this case.

Also, looping column by column performs better here than row by row, as I fill the trajectories via columns in this case.

Thank you a lot! I tried to get it to run in Julia 1.2, but failed so far. I will think about how I could extend it to more dimension without copying at least p.