I would appreciate advise how to best implement the following:
I have an N^d array A, where the dimension d depends on the
application (d \in \{1, 2, 3\}). Somewhere else I have a list L which is a
d x M array or integers corresponding to points/elements in this array,
e.g. if
L[:, 3] == [5, 2, 8]
then A would be N^3 and
M[L[1,3], L[2,3], L[3,3]]
would read
A[5, 2, 8]
The difficulty is that I would like to read from A in a *fast*
dimension-independent way. My only idea was to define a type and overload
getindex, setindex! ; see below. If I tested this correctly, then I only
loose a factor of two in terms of both speed and memory.
I am slightly worried that it is a non-standard access to a standard array.
+ I'd love to get that factor 2 back. Hence, I'd like to know whether
there are faster / more elegant alternatives, or just other ways of doing
this?
module damod
export darray
immutable type darray
data
dim
end
function darray(A)
return darray(A, length(size(A)))
end
export getindex
function getindex(A::darray, ii)
if A.dim == 2
return A.data[ii[1], ii[2]]
elseif A.dim == 3
return A.data[ii[1], ii[2], ii[3]]
end
end
end
using damod
N = 100
a = darray(rand(N, N, N))
@time for n = 1:N, m=1:N, k=1:N; dummy = a[[n,m,k]] end
@time for n = 1:N, m=1:N, k=1:N; dummy = a.data[n,m,k] end
@time for n = 1:N, m=1:N, k=1:N; dummy = a[[n,m,k]] end
@time for n = 1:N, m=1:N, k=1:N; dummy = a.data[n,m,k] end
@time for n = 1:N, m=1:N, k=1:N; dummy = a[[n,m,k]] end
@time for n = 1:N, m=1:N, k=1:N; dummy = a.data[n,m,k] end
elapsed time: 0.470031905 seconds (128633244 bytes allocated, 11.90% gc time)
elapsed time: 0.209358484 seconds (48565624 bytes allocated, 8.98% gc time)
elapsed time: 0.42617227 seconds (128641396 bytes allocated, 16.12% gc time)
elapsed time: 0.21529278 seconds (48565624 bytes allocated, 10.65% gc time)
elapsed time: 0.457094657 seconds (128565624 bytes allocated, 14.37% gc time)
elapsed time: 0.201969234 seconds (48565624 bytes allocated, 10.24% gc time)
