Pagewise multiplication along a dimension without using a for loop?

Hello, is there a way code pagewise multiplication along a single dimension of a matrix without using a for loop?

In the example below, I have a single 4D matrix where I use pagewise multiplication along its 4th dimension in a for loop to achieve a 3D matrix.

Is there a way to remove the for loop from this?

X = randi([1 10],2,2,10,5);
Y = pagemtimes(X(:,:,:,1),X(:,:,:,2));
for i = 3:length(X(1,1,1,:))
    Y = pagemtimes(Y,X(:,:,:,i));
end

Answer 1

[Matt J]

No, but looking at your code, I would guess that the speed is being encumbered much more by the operation length(X(1,1,1,:)) than by anything else. You should instead do,

for i = 3:size(X,4)

Y = pagemtimes(Y,X(:,:,:,i));

end

You would need an incredibly large size(X,4) relative to the other dimensions for the for-loop itself to be a significant bottleneck. If that's the case though, you might be able to benefit from a parfor loop, as demonstrated below,

X = randi([1 10],2,2,10,1e6);
Y = pagemtimes(X(:,:,:,1),X(:,:,:,2));

tic;
parfor i = 3:size(X,4)
    Y = pagemtimes(Y,X(:,:,:,i));
end
toc

Elapsed time is 0.784269 seconds.

tic;
for i = 3:size(X,4)
    Y = pagemtimes(Y,X(:,:,:,i));
end
toc

Elapsed time is 2.254356 seconds.

[Nathan Zechar]

Thank you for your answer. I've played around with parfor in the past, and it has generally been slower for what I've applied it to. And it is in the case as well.

The last dimension you've changed from a length of 5, to a length of 1e6. In which case, the regular for loop is significantly faster.

[Matt J]

The last dimension you've changed from a length of 5, to a length of 1e6

Yes, but if size(X,4)=5 is what you have in practical examples, it is not worth worrying about for-loop overhead. As I said, it has to be much larger for the loop to be the problem.

[Nathan Zechar]

X could range from 5 to maybe as high as 100 for my application,

The output gives an analytical solution used to minimize data to find an unknown. So process has to be called several times to generate central difference derivatives within a loop, and that process repeats until convergence critera has been reached. It adds up. But aside from that, I was just curious if something like this could be done within matlab.

[Matt J]

The Optimization Toolbox solvers let yoy parallelize central difference computations in their iterative loops. That might be somethign to consider.

Answer 2

[Matt J]

If you can build X in cell array form, instead of 4D form, you can save some time as well:

X = randi([1 10],2,2,10,5);

Xc = num2cell( X , 1:3  );

timeit(@()version1(X))*1e6

ans = 18.1973

timeit(@()version2(Xc))*1e6

ans = 6.8640

function version1(X)
 Y = pagemtimes(X(:,:,:,1),X(:,:,:,2));
 for i = 3:size(X,4)
    Y = pagemtimes(Y,X(:,:,:,i));
 end
end

function version2(Xc)
 Y = pagemtimes(Xc{1},Xc{2});
 for i = 3:numel(Xc)
    Y = pagemtimes(Y,Xc{i});
 end
end