Wednesday, August 8, 2012

squeeze() is slow!


Test code to see how whether to use squeeze or use direct index addressing is faster

=========

clear all;


N = 50000;
tic;
for kk=1:N
x = rand(2,2,50);
y1 = x(1,1,:);
y2 = x(1,2,:);
y3 = x(2,1,:);
y4 = x(2,2,:);

z1 = y1(1,1,:).*conj(y1(1,1,:)) + y3(1,1,:).*conj(y3(1,1,:));
z2 = y2(1,1,:).*conj(y2(1,1,:)) + y4(1,1,:).*conj(y4(1,1,:));
end
toc;

tic;
for kk=1:N
x = rand(2,2,50);
y1 = squeeze(x(1,1,:));
y2 = squeeze(x(1,2,:));
y3 = squeeze(x(2,1,:));
y4 = squeeze(x(2,2,:));

z1 = y1.*conj(y1) + y3.*conj(y3);
z2 = y2.*conj(y2) + y4.*conj(y4);
end
toc;

=======

Execution Results :

Elapsed time is 2.188846 seconds.
Elapsed time is 11.824104 seconds.

Conclusion : do not use the squeeze() function in Matlab if you want speed!


Tuesday, August 7, 2012

Vector Multiplication of Array of Matrices

Testing multiplication of vector of matrices in Matlab.

> conclusion : loop explicit unrolling to a singleton vector multiplication is the fastest!


function [t,v] = matrixMultTest()
    n = 2; m = 2; p = 1e5;
    A = rand(n,m,p);
    B = rand(m,n,p);

    %# time functions
    t = zeros(5,1);
    t(1) = timeit( @() func1(A,B,n,m,p) );
    t(2) = timeit( @() func2(A,B,n,m,p) );
    t(3) = timeit( @() func3(A,B,n,m,p) );
    t(4) = timeit( @() func4(A,B,n,m,p) );
    t(5) = timeit( @() func5(A,B,n,m,p) );

    %# check the results
    v = cell(5,1);
    v{1} = func1(A,B,n,m,p);
    v{2} = func2(A,B,n,m,p);
    v{3} = func3(A,B,n,m,p);
    v{4} = func4(A,B,n,m,p);
    v{5} = func5(A,B,n,m,p);
    assert( isequal(v{:}) )
end
%# simple FOR-loop
function C = func1(A,B,n,m,p)
    C = zeros(n,n,p);
    for k=1:p
        C(:,:,k) = A(:,:,k) * B(:,:,k);
    end
end
%# ARRAYFUN
function C = func2(A,B,n,m,p)
    C = arrayfun(@(k) A(:,:,k)*B(:,:,k), 1:p, 'UniformOutput',false);
    C = cat(3, C{:});
end
%# NUM2CELL/FOR-loop/CELL2MAT
function C = func3(A,B,n,m,p)
    Ac = num2cell(A, [1 2]);
    Bc = num2cell(B, [1 2]);
    C = cell(1,1,p);
    for k=1:p
        C{k} = Ac{k} * Bc{k};
    end;
    C = cell2mat(C);
end
%# NUM2CELL/CELLFUN/CELL2MAT
function C = func4(A,B,n,m,p)
    Ac = num2cell(A, [1 2]);
    Bc = num2cell(B, [1 2]);
    C = cellfun(@mtimes, Ac, Bc, 'UniformOutput', false);
    C = cell2mat(C);
end
%# Loop Unrolling
function C = func5(A,B,n,m,p)
    C = zeros(n,n,p);
    C(1,1,:) = A(1,1,:).*B(1,1,:) + A(1,2,:).*B(2,1,:);
    C(1,2,:) = A(1,1,:).*B(1,2,:) + A(1,2,:).*B(2,2,:);
    C(2,1,:) = A(2,1,:).*B(1,1,:) + A(2,2,:).*B(2,1,:);
    C(2,2,:) = A(2,1,:).*B(1,2,:) + A(2,2,:).*B(2,2,:);
end
The results:
>> [t,v] = matrixMultTest();
>> t
t =
      0.63633      # FOR-loop
      1.5902       # ARRAYFUN
      1.1257       # NUM2CELL/FOR-loop/CELL2MAT
      1.0759       # NUM2CELL/CELLFUN/CELL2MAT
      0.05712      # Loop Unrolling



From stackoverflow http://stackoverflow.com/questions/6580656/matlab-how-to-vector-multiply-two-arrays-of-matrices