Previous Up Next

6  Miscellaneous Functions

6.1  Function mexConjGrad

Implementation of a conjugate gradient for solving a linear system Ax=b when A is positive definite. In some cases, it is faster than the Matlab function pcg, especially when the library uses the Intel Math Kernel Library.


% Usage:   x =mexConjGrad(A,b,x0,tol,itermax)
%
% Name: mexConjGrad
%
% Description: Conjugate gradient algorithm, sometimes faster than the 
%    equivalent Matlab function pcg. In order to solve Ax=b;
%
% Inputs: A:  double square n x n matrix. HAS TO BE POSITIVE DEFINITE
%         b:  double vector of length n.
%         x0: double vector of length n. (optional) initial guess.
%         tol: (optional) tolerance.
%         itermax: (optional) maximum number of iterations.
%
% Output: x: double vector of length n.
%
% Author: Julien Mairal, 2009

The following piece of code contains usage examples:

A=randn(5000,500);
A=A'*A;
b=ones(500,1);
x0=b;
tol=1e-4;
itermax=0.5*length(b);

tic
for ii = 1:20
x1 = mexConjGrad(A,b,x0,tol,itermax);
end
t=toc;
fprintf('
mex-file time: %fs\n',t);

tic
for
 ii = 1:20
x2 = pcg(A,b);
end
t=toc;
fprintf('Matlab time: %fs\n',t);
sum((x1(:)-x2(:)).^2)

6.2  Function mexBayer

Apply a Bayer pattern to an input image


% Usage:   x =mexConjGrad(A,b,x0,tol,itermax)
%
% Name: mexConjGrad
%
% Description: Conjugate gradient algorithm, sometimes faster than the 
%    equivalent Matlab function pcg. In order to solve Ax=b;
%
% Inputs: A:  double square n x n matrix. HAS TO BE POSITIVE DEFINITE
%         b:  double vector of length n.
%         x0: double vector of length n. (optional) initial guess.
%         tol: (optional) tolerance.
%         itermax: (optional) maximum number of iterations.
%
% Output: x: double vector of length n.
%
% Author: Julien Mairal, 2009

The following piece of code contains usage examples:

A=randn(5000,500);
A=A'*A;
b=ones(500,1);
x0=b;
tol=1e-4;
itermax=0.5*length(b);

tic
for ii = 1:20
x1 = mexConjGrad(A,b,x0,tol,itermax);
end
t=toc;
fprintf('
mex-file time: %fs\n',t);

tic
for
 ii = 1:20
x2 = pcg(A,b);
end
t=toc;
fprintf('Matlab time: %fs\n',t);
sum((x1(:)-x2(:)).^2)

6.3  Function mexCalcAAt

For an input sparse matrix A, this function returns AAT. For some reasons, when A has a lot more columns than rows, this function can be much faster than the equivalent matlab command X*A'.


% Usage:   AAt =mexCalcAAt(A);
%
% Name: mexCalcAAt
%
% Description: Compute efficiently AAt = A*A', when A is sparse 
%   and has a lot more columns than rows. In some cases, it is
%   up to 20 times faster than the equivalent Matlab expression
%   AAt=A*A';
%
% Inputs: A:  double sparse m x n matrix   
%
% Output: AAt: double m x m matrix 
%
% Author: Julien Mairal, 2009

The following piece of code contains usage examples:

A=sprand(200,200000,0.05);

tic
AAt=mexCalcAAt(A);
t=toc;
fprintf('mex-file time: %fs\n',t);

tic
AAt2=A*A';
t=toc;
fprintf('
matlab time: %fs\n',t);

sum((AAt(:)-AAt2(:)).^2)

6.4  Function mexCalcXAt

For an input sparse matrix A and a matrix X, this function returns XAT. For some reasons, when A has a lot more columns than rows, this function can be much faster than the equivalent matlab command X*A'.


% Usage:   XAt =mexCalcXAt(X,A);
%
% Name: mexCalcXAt
%
% Description: Compute efficiently XAt = X*A', when A is sparse and has a 
%   lot more columns than rows. In some cases, it is up to 20 times 
%   faster than the equivalent Matlab expression XAt=X*A';
%
% Inputs: X:  double m x n matrix
%         A:  double sparse p x n matrix   
%
% Output: XAt: double m x p matrix 
%
% Author: Julien Mairal, 2009

The following piece of code contains usage examples:

X=randn(64,200000);
A=sprand(200,200000,0.05);

tic
XAt=mexCalcXAt(X,A);
t=toc;
fprintf('mex-file time: %fs\n',t);

tic
XAt2=X*A';
t=toc;
fprintf('
mex-file time: %fs\n',t);

sum((XAt(:)-XAt2(:)).^2)

6.5  Function mexCalcXY

For two input matrices X and Y, this function returns XY.


% Usage:   Z =mexCalcXY(X,Y);
%
% Name: mexCalcXY
%
% Description: Compute Z=XY using the BLAS library used by SPAMS.
%
% Inputs: X:  double m x n matrix
%         Y:  double n x p matrix   
%
% Output: Z: double m x p matrix 
%
% Author: Julien Mairal, 2009

The following piece of code contains usage examples:

X=randn(64,200);
Y=randn(200,20000);

tic
XY=mexCalcXY(X,Y);
t=toc;
fprintf('mex-file time: %fs\n',t);

tic
XY2=X*Y;
t=toc;
fprintf('mex-file time: %fs\n',t);

sum((XY(:)-XY2(:)).^2)

6.6  Function mexCalcXYt

For two input matrices X and Y, this function returns XYT.


% Usage:   Z =mexCalcXYt(X,Y);
%
% Name: mexCalcXYt
%
% Description: Compute Z=XY' using the BLAS library used by SPAMS.
%
% Inputs: X:  double m x n matrix
%         Y:  double p x n matrix   
%
% Output: Z: double m x p matrix 
%
% Author: Julien Mairal, 2009

The following piece of code contains usage examples:

X=randn(64,200);
Y=randn(200,20000)';

tic
XYt=mexCalcXYt(X,Y);
t=toc;
fprintf('
mex-file time: %fs\n',t);


tic
XYt2=X*Y';
t=toc;
fprintf('
matlab-file time: %fs\n',t);

sum((XYt(:)-XYt2(:)).^2)

6.7  Function mexCalcXtY

For two input matrices X and Y, this function returns XTY.


% Usage:   Z =mexCalcXtY(X,Y);
%
% Name: mexCalcXtY
%
% Description: Compute Z=X'Y using the BLAS library used by SPAMS.
%
% Inputs: X:  double n x m matrix
%         Y:  double n x p matrix   
%
% Output: Z: double m x p matrix 
%
% Author: Julien Mairal, 2009

The following piece of code contains usage examples:

X=randn(64,200)';
Y=randn(200,20000);

tic
XtY=mexCalcXtY(X,Y);
t=toc;
fprintf('
mex-file time: %fs\n',t);

tic
XtY2=X'*Y;
t=toc;
fprintf('
matlab-file time: %fs\n',t);

sum((XtY(:)-XtY2(:)).^2)

6.8  Function mexInvSym

For an input symmetric matrices A in ℝn × n, this function returns A−1.


% Usage:   B =mexInvSym(A);
%
% Name: mexInvSym
%
% Description: returns the inverse of a symmetric matrix A
%
% Inputs: A:  double n x n matrix   
%
% Output: B: double n x n matrix 
%
% Author: Julien Mairal, 2009

The following piece of code contains usage examples:

A=rand(1000,1000);
A=A'*A;

tic
B=mexInvSym(A);
t=toc;
fprintf('
mex-file time: %fs\n',t);

tic
B2=inv(A);
t=toc;
fprintf('matlab-file time: %fs\n',t);

sum((B(:)-B2(:)).^2)

6.9  Function mexNormalize


% Usage:   Y =mexNormalize(X);
%
% Name: mexNormalize
%
% Description: rescale the columns of X so that they have
%        unit l2-norm.
%
% Inputs: X:  double m x n matrix   
%
% Output: Y: double m x n matrix 
%
% Author: Julien Mairal, 2010

The following piece of code contains usage examples:

X=rand(100,1000);

tic
Y=mexNormalize(X);
t=toc;

6.10  Function mexSort


% Usage:   Y =mexSort(X);
%
% Name: mexSort
%
% Description: sort the elements of X using quicksort
%
% Inputs: X:  double vector of size n
%
% Output: Y: double  vector of size n
%
% Author: Julien Mairal, 2010

The following piece of code contains usage examples:

X=rand(1,300000);

tic
Y=mexSort(X);
t=toc;
toc

tic

Y2=sort(X,'ascend');
t=toc;
toc

sum
((Y2(:)-Y(:)).^2)

6.11  Function mexDisplayPatches

Print to the screen a matrix containing as columns image patches.

6.12  Function mexCountPathsDAG

This function counts the number of paths in a DAG using dynamic programming.


% Usage:   num=mexCountPathsDAG(G);
%
% Name: mexCountPathsDAG
%
% Description: mexCountPathsDAG counts the number of paths 
%       in a DAG.
%
%       for a graph G with |V| nodes and |E| arcs,
%       G is a double sparse adjacency matrix of size |V|x|V|,
%       with |E| non-zero entries.
%       (see example in test_CountPathsDAG.m
%
%
% Inputs: G:  double sparse |V| x |V| matrix (full graph)
%
% Output: num: number of paths
%
% Author: Julien Mairal, 2012

The following piece of code contains usage examples:

% graph corresponding to figure 1 in arXiv:1204.4539v1
fprintf('test mexCountPathsDAG\n');
% this graph is a DAG
G=[0 0 0 0 0 0 0 0 0 0 0 0 0;
   1 0 0 0 0 0 0 0 0 0 0 0 0;
   1 1 0 0 0 0 0 0 0 0 0 0 0;
   1 1 0 0 0 0 0 0 0 0 0 0 0;
   0 0 0 1 0 0 0 0 0 0 0 0 0;
   0 1 1 0 1 0 0 0 0 0 0 0 0;
   0 1 0 0 1 0 0 0 0 0 0 0 0;
   0 0 0 0 0 1 1 0 0 0 0 0 0;
   1 0 0 1 0 0 0 0 0 0 0 0 0;
   1 0 0 0 0 0 0 0 1 0 0 0 0;
   0 0 0 0 0 0 0 0 1 1 0 0 0;
   0 0 0 0 0 0 0 0 1 0 1 0 0;
   0 0 0 0 1 0 0 0 1 0 0 1 0];
G=sparse(G);
num=mexCountPathsDAG(G);
fprintf('Num of paths: %d\n',num);

6.13  Function mexRemoveCyclesGraph

One heuristic to remove cycles from a graph.


% Usage:   DAG=mexRemoveCycleGraph(G);
%
% Name: mexRemoveCycleGraph
%
% Description: mexRemoveCycleGraph heuristically removes
%       arcs along cycles in the graph G to obtain a DAG.
%       the arcs of G can have weights. The heuristic will
%       remove in priority arcs with the smallest weights.
%
%       for a graph G with |V| nodes and |E| arcs,
%       G is a double sparse adjacency matrix of size |V|x|V|,
%       with |E| non-zero entries. The non-zero entries correspond
%       to the weights of the arcs.
%
%       DAG is a also double sparse adjacency matrix of size |V|x|V|,
%       but the graph is acyclic.
%
%       Note that another heuristic to obtain a DAG from a general 
%       graph consists of ordering the vertices.
%
% Inputs: G:  double sparse |V| x |V| matrix
%
% Output: DAG:  double sparse |V| x |V| matrix
%
% Author: Julien Mairal, 2012

The following piece of code contains usage examples:

fprintf('test mexRemoveCyclesGraph\n');
% this graph is not a DAG
G=[0   0   0   0   0   0   0   0   0   0   0   0   0;
   1   0   0   0.1 0   0   0   0.1 0   0   0   0   0;
   1   1   0   0   0   0.1 0   0   0   0   0   0   0;
   1   1   0   0   0   0   0   0   0   0   0   0.1 0;
   0   0   0   1   0   0   0   0   0   0   0   0   0;
   0   1   1   0   1   0   0   0   0   0   0   0   0;
   0   1   0   0   1   0   0   0   0   0   0   0   0;
   0   0   0   0   0   1   1   0   0   0   0   0   0;
   1   0   0   1   0   0   0   0   0   0   0   0   0;
   1   0   0   0   0   0   0   0   1   0   0   0   0;
   0   0   0   0   0   0   0   0   1   1   0   0   0;
   0   0   0   0   0   0   0   0   1   0   1   0   0;
   0   0   0   0   1   0   0   0   1   0   0   1   0];
G=sparse(G);
DAG=mexRemoveCyclesGraph(G);
format compact;
fprintf('Original graph:\n');
full(G)
fprintf('New graph:\n');
full(DAG)

6.14  Function mexCountConnexComponents

Count the number of connected components of a subgraph from a graph.


% Usage:   num=mexCountConnexComponents(G,N);
%
% Name: mexCountConnexComponents
%
% Description: mexCountConnexComponents counts the number of connected
%       components of the subgraph of G corresponding to set of nodes
%       in N. In other words, the subgraph of G by removing from G all
%       the nodes which are not in N.
%
%       for a graph G with |V| nodes and |E| arcs,
%       G is a double sparse adjacency matrix of size |V|x|V|,
%       with |E| non-zero entries.
%       (see example in test_CountConnexComponents.m)
%       N is a dense vector of size |V|. if  N[i] is non-zero,
%       it means that the node i is selected.
%
%
% Inputs: G:  double sparse |V| x |V| matrix (full graph)
%         N:  double full |V| vector.
%
% Output: num: number of connected components.
%
% Author: Julien Mairal, 2012

The following piece of code contains usage examples:

% graph corresponding to figure 1 in arXiv:1204.4539v1
fprintf('test mexCountConnexComponents\n');
% this graph is a DAG
G=[0 0 0 0 0 0 0 0 0 0 0 0 0;
   1 0 0 0 0 0 0 0 0 0 0 0 0;
   1 1 0 0 0 0 0 0 0 0 0 0 0;
   1 1 0 0 0 0 0 0 0 0 0 0 0;
   0 0 0 1 0 0 0 0 0 0 0 0 0;
   0 1 1 0 1 0 0 0 0 0 0 0 0;
   0 1 0 0 1 0 0 0 0 0 0 0 0;
   0 0 0 0 0 1 1 0 0 0 0 0 0;
   1 0 0 1 0 0 0 0 0 0 0 0 0;
   1 0 0 0 0 0 0 0 1 0 0 0 0;
   0 0 0 0 0 0 0 0 1 1 0 0 0;
   0 0 0 0 0 0 0 0 1 0 1 0 0;
   0 0 0 0 1 0 0 0 1 0 0 1 0];
G=sparse(G);
nodes=[0 1 1 0 0 1 0 0 1 0 1 1 0];
num=mexCountConnexComponents(G,nodes);
fprintf('Num of connected components: %d\n',num);

% this graph is a not a DAG anymore. This function works
% with general graphs.

G=[0 0 0 0 0 0 0 0 0 0 0 0 0;
   1 0 0 0 0 0 0 0 0 0 0 0 0;
   1 1 0 1 0 0 0 0 0 0 0 0 0;
   1 1 0 0 0 0 0 0 0 0 0 0 0;
   0 0 0 1 0 0 0 0 0 0 0 0 0;
   0 1 1 0 1 0 0 0 0 0 0 0 0;
   0 1 0 0 1 0 0 0 0 0 0 0 0;
   0 0 0 0 0 1 1 0 0 0 0 0 0;
   1 0 0 1 0 0 0 0 0 0 0 0 0;
   1 0 0 0 0 0 0 0 1 0 0 0 0;
   0 0 0 0 0 0 0 0 1 1 0 0 0;
   0 0 0 0 0 0 0 0 1 0 1 0 0;
   0 0 0 0 1 0 0 0 1 0 0 1 0];
nodes=[0 1 1 0 0 1 0 0 1 0 1 1 0];
G=sparse(G);
num=mexCountConnexComponents(G,nodes);
fprintf('Num of connected components: %d\n',num);

nodes=[0 1 1 1 0 1 0 0 1 0 1 1 0];
num=mexCountConnexComponents(G,nodes);
fprintf('Num of connected components: %d\n',num);

Previous Up Next