3

I have the following large, very inefficient loop.

P is a [2000 x 200 x 5] matrix
D is a [2000 x 200 x 5] matrix
S is a [200 x 1005] matrix
PS is a [2000 x 1000 x 5] matrix

I want to compute the following loop:

for k=1:2000
   for n=1:200
      for t=1:5
          P(k,n,t) = sum(S(n,t+1:t+1000) .* PS(k,1:1000,t));
      end
   end
end

Obviously this is very inefficient. I tried parfor, but I would rather a vectorized solution. I tried couple of things with bsxfun, but also never managed to get it working.

Thank you.

Divakar
  • 218,885
  • 19
  • 262
  • 358
phdstudent
  • 1,060
  • 20
  • 41
  • Could you a small sample data and the expected result,,please? It would help to verify any answers. – kkuilla Dec 07 '15 at 09:01
  • The data comes from a long code. I think for the purposes at hand any random data should do the trick: P = rand( 2000,200,5); D = rand( 2000,200,5); S = rand( 200,1005); PS = rand( 2000,1000,5); Thank you! – phdstudent Dec 07 '15 at 09:03
  • @volcompt please add data (or random data) to your question (i.e. not in the comments) and perhapd change the scale of the problem for the example or else replace the `1000`, `200` and `5` with variables we can alter to work on a smaller problem. Also, you should post your `bsxfun` attempts and tell us where they go wrong. Lastly, your current code completely ignores the first column of `S`, is that correct? – Dan Dec 07 '15 at 09:15

1 Answers1

6

Here's an almost (almost because we still have a loop, but with only 5 iterations) vectorized approach using powerful matrix-multiplication -

out = zeros(2000,200,5);
for t=1:size(P,3) %// size(P,3) = 5
    out(:,:,t) = PS(:,:,t)*S(:,t+1:t+1000).';
end

Runtime tests and verify output -

%// Inputs
D = rand(2000,200,5);
S = rand(200,1005);
PS = rand(2000,1000,5);

disp('--------------------- No Matrix-mult-fun')
tic
P = zeros(2000,200,5);
for k=1:2000
   for n=1:200
      for t=1:5
          P(k,n,t) = sum(S(n,t+1:t+1000) .* PS(k,1:1000,t));
      end
   end
end
toc

disp('--------------------- Fun fun Matrix-mult-fun')
tic
out = zeros(2000,200,5);
for t=1:size(P,3) %// size(P,3) = 5
    out(:,:,t) = PS(:,:,t)*S(:,t+1:t+1000).';
end
toc

error_val = max(abs(P(:)-out(:)))

Output -

--------------------- No Matrix-mult-fun
Elapsed time is 70.223008 seconds.
--------------------- Fun fun Matrix-mult-fun
Elapsed time is 0.624308 seconds.
error_val =
     1.08e-12
Community
  • 1
  • 1
Divakar
  • 218,885
  • 19
  • 262
  • 358
  • 1
    Of, course... Whom else....? :-) The summation is done by the matrix multiplication, I presume? – kkuilla Dec 07 '15 at 09:24
  • That was impressive. Thank you Divakar! – phdstudent Dec 07 '15 at 09:35
  • @Divakar - on a sidenote that is unrelated to your code: you can add `` before `code` sections where you don't need syntax highlighting to happen (such as in the _Output_ box in your answer). – Dev-iL Dec 07 '15 at 12:29
  • @Dev-iL Hmm, would see how to use it in future posts, thanks for chipping in with that info! – Divakar Dec 07 '15 at 12:50