Friday, September 14, 2012

Memory Leak in Matlab (multi-dimensional array passing into a function)

First of all, Memory Leak in Matlab maybe misleading.

There is a bug (well I think it is a bug) in Matlab, where memory garbage collection is not done quickly enough, which can result in accumulation of non-usable memory space.

This accumulation of non-usable memory can lead to crashes of the code you are running in Matlab. I have confirm this issue at least in Matlab 2009a. I don't know if this is still an issue for higher Matlab versions. It is likely that it is still an issue and style of coding in Matlab to avoid this issue may be preferred.


Let's look at a function "y = myfunc( x )" which is called multiple times (many many many times) in the whole algorithm. myfunc() take in a matrix and does some processing on it (without changing the value of the input matrix itself).

e.g.

-----problematic function-----------------
function [ y ] =  myfunc( x )

% take x and process something to get results
w = rand( 100, 1 );
y = w'*x*w;

end % end of function
------------------------------------------


Lets assume the input of the function belongs to a higher dimensional array. For example if input size of the function is [ 100, 100 ] dimensional array, the actual array size is [ 100, 100, 4 ] dimensional. We would input a portion of the multi-dimensional array into the function.

An example is shown below.

----- problematic algorithm --------------

full_x =  rand( 100, 100, 4 );

for ii=1:100000
     myfunc( full_x(:,:,1) );  % copy of full_x(:,:,1) is performed "under the hood"

     myfunc( full_x(:,:,2) );  % the copy of full_x(:,:,1) is not released and being accumulated in memory
     myfunc( full_x(:,:,3) );
     myfunc( full_x(:,:,4) );

end
------------------------------------------

With the above implementation, when passing the input argument into the function Matlab makes a temporary copy of x(:,:,k) and passes the temporary copy into the function (where it is processed). Now the problem is not only the Matlab makes a partial copy of the variable, the garbage collection is not done properly. The temporary copy is still somewhere in Matlab memory taking up space for each call of myfunc(). Because of this Matlab will crash with a "out of memory" error after some time.


A method of coding to go around this issue to pass the ENTIRE variable into the function, along with the index information in which you need to process the variable. An example implementation is shown below.


-----more efficient function-----------------
function [ y ] =  myfunc( x, idx )

% take x and only use the first two dimensions along indice "idx"
% and process something to get results
w = rand( 100, 1 );
y = w'*x(:,:,.idx)*w;

end % end of function
------------------------------------------


----- algorithm fixed ---------------------

full_x = rand( 100, 100, 4 );

for ii=1:100000
     myfunc( full_x, 1 ); % no copy of full_x

     myfunc( full_x, 2 );
     myfunc( full_x, 3 );
     myfunc( full_x, 4 );

end
------------------------------------------


As long as the input variable is not changed in value, Matlab does not copy the entire input variable when calling the function "myfunc()". No copy is performed and thus no memory leakage(?). The only draw-back is that the myfunc() now requires information of the entire variable dimensions, and also require additional information of the index of the variable in which the processing needs to be done.


In many cases this is not going to be an issue, but if handling very large data sets in thousands or millions, this is really an critical issue for stability of the entire code.


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