Hi,
I have two 3d matrices: S of size [m, n, n] and W of size [n, n, n] and I need to compute the matrix X of size [m, n, n^2] such that
X(:, :, n * (i-1) + j) = squeeze(W(:, :, i)) * squeeze(S(:, :, j))
I would like to avoid doing a double for loop on i and j, in order to keep things as much efficient as possible. I was able to avoid one of the two for loops by doing something like
X(:, :, n*(i-1)+idx) = pagemtimes(S, squeeze(W(:, :, i)));
However, profiling the code I noticed that most of the time is spent in memory management. Indeed, if rewrite the code as
aux = pagemtimes(S, squeeze(W(:, :, i)));
X(:, :, n*(i-1)+idx) = aux;
25% of the time is effectively spent in computing the product, while 75% of the time is spent in assigning X(:, :, n*(i-1)+idx) = aux. In case it is useful, in my case m = 100, n = 10, but I get a similar ratio also if I use m = 10000, n = 10. I am not sure whether this can be avoided, but as a possible solution I was trying to assemble the matrix X all together. The only idea that came to my mind was to do something like
i = repmat(1:n, [n, 1]); i = i(:)';
X = pagemtimes(S(:, :, j), W_hat(:, :, i));
which (I believe) copies the matrices S and W under the hood, and is not much more efficient for this reason. Is there a way to avoid execute this operation without duplicating memory and without the slow memory management I discussed above?
Thanks,
Ivan
