As I understand it, you have 2 dimensions of data. You have not specified how you want to sort them (e.g. "along dimension 1"). If, in fact, these dimensions are not relevant to the sorting (i.e. you want to sort all values into one single list), then you can effectively vectorise this 2D matrix using the reshape() operator.
More importantly, your indexing is very excessive in its use of dimensions. You can easily store all the indexing information in just one dimension by encoding each [i1,i2,...,i10] into a single number. After all, there are only 5^10 (i.e. fewer than 10 million) possible index combinations - so we can easily fit this into one number. Here is an example:
close all; clear all; clc;
% N example index values in [1,M]
N = 10;
M = 5;
indexes = randi(M,N,1);
% Encode indices as a single number
encoder = M.^[0:N-1];
code = encoder*indexes
% Decode
decoded = zeros(N,1);
for i = 1:N-1
next = ceil(code/encoder(i+1)-1)*encoder(i+1);
decoded(i) = (code - next) / encoder(i);
code = next;
end
decoded(N) = code/encoder(end);
% Check the code/decode works
assert(isequal(indexes, decoded));
This will, at worst, leave you with a 3D array, which should be straightforward to handle.
