8

I'm fairly new to MATLAB. Normal matrix multiplication of a M x K matrix by an K x N matrix -- C = A * B -- has c_ij = sum(a_ik * b_kj, k = 1:K). What if I want this to be instead c_ij = sum(op(a_ik, b_kj), k = 1:K) for some simple binary operation op? Is there any nice way to vectorize this in MATLAB (or maybe even a built-in function)?

EDIT: This is currently the best I can do.

% A is M x K, B is K x N
% op is min
C = zeros(M, N);
for i = 1:M:
    C(i, :) = sum(bsxfun(@min, A(i, :)', B));
end
Shai
  • 111,146
  • 38
  • 238
  • 371
Nick
  • 2,821
  • 5
  • 30
  • 35

4 Answers4

2

Listed in this post is a vectorized approach that persists with bsxfun by using permute to create singleton dimensions as needed by bsxfun to let the singleton-expansion do its work and thus essentially replacing the loop in the original post. Please be reminded that bsxfun is a memory hungry implementation, so expect speedup with it only until it is stretched too far. Here's the final solution code -

op = @min;   %// Edit this with your own function/ operation
C = sum(bsxfun(op, permute(A,[1 3 2]),permute(B,[3 2 1])),3)

NB - The above solution was inspired by Removing four nested loops in Matlab.

Community
  • 1
  • 1
Divakar
  • 218,885
  • 19
  • 262
  • 358
  • The documentation for bsxfun does not indicate it is memory hungry. How did you determine that fact? I was unaware of that issue with bsxfun. – John Mar 25 '15 at 02:35
  • 2
    @John Well [`bsxfun`](http://in.mathworks.com/help/matlab/ref/bsxfun.html) by its very definition expands on singleton dimensions internally and then does elementwise operations. That expansion essentially needs such a large contiguous area on memory. This search for contiguous memory means more runtime. This "memory hungry" nature could also be confirmed by watching system monitor as the bsxfun based codes are run. But the upside of using `bsxfun` is that once the memory is allocated, the operations are done in a super fast way. Hope that made some sense! – Divakar Mar 25 '15 at 04:09
1

if the operator can operate element-by-element (like .*):

if(size(A,2)~=size(B,1))
    error(blah, blah, blah...);
end

C = zeros(size(A,1),size(B,2));
for i = 1:size(A,1)
    for j = 1:size(B,2)
        C(i,j) = sum(binaryOp(A(i,:)',B(:,j)));
    end
end
Cavaz
  • 2,996
  • 24
  • 38
  • I actually have a way to do it with one for-loop (using bsxfun). I was wondering if there was a purely vectorized way to do it. – Nick Nov 09 '11 at 02:10
  • I added my current algorithm to the original post. – Nick Nov 09 '11 at 07:37
0

Depending on your specific needs, you may be able to use bsxfun in 3D to trick the binary operator. See this answer for more infos: https://stackoverflow.com/a/23808285/1121352 Another alternative would be to use cellfun with a custom function: http://matlabgeeks.com/tips-tutorials/computation-using-cellfun/

Community
  • 1
  • 1
gaborous
  • 15,832
  • 10
  • 83
  • 102
  • Thanks for getting my attention to this problem! – Divakar Jun 14 '14 at 17:14
  • @Divakar Thank's to you for your great answers! If you are interested by another similar challenge, I have posted a question about Generalized Matrix Multiplication, the goal being to go further in the optimization: http://stackoverflow.com/questions/24245225/matlab-generalized-matrix-multiplication – gaborous Jun 16 '14 at 13:45
0

You can always write the loops yourself:

A = rand(2,3);
B = rand(3,4);

op = @times;            %# use your own function here
C = zeros(size(A,1),size(B,2));
for i=1:size(A,1)
    for j=1:size(B,2)
        for k=1:size(A,2)
            C(i,j) = C(i,j) + op(A(i,k),B(k,j));
        end
    end
end

isequal(C,A*B)
Amro
  • 123,847
  • 25
  • 243
  • 454