## Archive for the 'matlab' Category

### Mess up Matlab codes and outputs in LaTeX

You do not want to mess up, right? When writing a LaTeX document, you may once in a while want to include some Matlab codes and/or outputs (preferably typeset using typewriter font if you have the same taste as me) during the course of your writing. What I used to do was to copy and paste the Matlab codes into my LaTeX file, execute the codes in Matlab, then do another copy and paste to place the results in my LaTeX file, and finally decorate them in a verbatim block or something like that. Guess what, Matlab provides a command, called publish, that helps you do all these in a simpler way.

In a nutshell, the way to use publish is to first type in the texts (as your usually LaTeX editing), including Matlab codes, in a single .m file. Let’s say the file name is example.m. Then, in Matlab, you issue the command

publish('example.m', struct('format','latex','outputDir','ltx-src'));

It means that you want Matlab to process example.m and output a LaTeX file example.tex (that you can compile to get pdf) in the sub-directory ltx-src. This is it. Instead of writing example.tex, you write a file example.m.

So, how should I write example.m? It is best to give an example. See the following:

%% % <latex> % The eigenvalues of a circulant matrix can be % obtained by performing FFT on the first column % of the matrix. First, let us construct a % $5\times5$ circulant matrix \verb|C| whose first % column \verb|c| is generated with random input: % </latex> 

c = rand(5,1); % sad that Matlab does not provide a circulant() % command... C = toeplitz(c, c([1 end:-1:2]))

%% % <latex> % The eigenvalues of \verb|C| are nothing but % </latex>

lambda = fft(c)

%% % <latex> % Check it out! The output is the same as using % the \verb|eig| command: % </latex>

eig(C)

%% % Fun, isn't it? 

It is nothing but a script file that Matlab can execute, right? The trick part is that all the texts and LaTeX markups are buried in comment blocks. How the Matlab command publish makes a LaTeX output is that whenever it meets a whole block of comments starting with ‘%%‘, it strips the comment signs and decorates the whole block using the pair \begin{par} and \end{par}. On the other hand, whenever it meets a block of codes that does not start with ‘%%‘, Matlab knows that they are executable commands. Matlab uses \begin{verbatim} and \end{verbatim} to typeset these command texts, and automatically add the Matlab outputs of the commands, which are also decorated by the \begin{verbatim} and \end{verbatim} pair, in the LaTeX file. Something I am not satisfied is that Matlab does not recognize LaTeX commands such as \verb||. I have to put <latex></latex> so that Matlab can do a verbatim copy of \verb||, instead of expand the text \verb|| in some weird way, in the output LaTeX file.

It is time to try the above example yourself. Have fun.

### Call by value or call by reference?

In C, if you pass in a pointer to a function, the content of the involved variable may be changed. In Matlab, however, if you pass in an array as an argument, the content of the array will never be changed after the function call, even if the variable with the same name is modified within the function, unless the variable is also output. Ok, now comes the question: What if the function called by Matlab is written in C as a mex function? Haha, have fun with the following three tests, t1, t2, and t3.

/* File: t1.c */ #include 

void foo(double *a) { a[0] = 1; }

 void main(void) { double a[2] = {0}; foo(a); printf("first element of array a = %g\n", a[0]); } 

% File: t2.m function t2()

a = zeros(2,1); foo(a); fprintf(1, 'first element of array a = %g\n', a(1));

function foo(a)

 a(1) = 1; 

% File: t3.m function t3()

 a = zeros(2,1); foo_mex(a); fprintf(1, 'first element of array a = %g\n', a(1)); 

/* File: foo_mex.c */ #include "mex.h"

 void mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[]) { double *a = mxGetPr(prhs[0]); a[0] = 1; } 

### A sample MATLAB code for defining functions with cases

The code should be self-explaning..

function y = f(x) y1 = x.* (x<0); y2 = x.^2 .* (x>=0) .* (x<2); y3 = 4 .* (x>=2); y = y1 + y2 + y3; 

### Find all the empty cells?

Suppose I have a cell array which looks like this:

 C = [] [] [] [] [] [3x3 double] [] [] []

How can I get all the empty cells using as few commands as possible? Well, the answer is …

cellfun(@isempty, C)

So, cellfun is a function that applies a function (isempty here) to each cell of a cell array (C here). If you enjoy general Matlab functions that implicitly apply to each element of an array (say, sqrt(M) does the square roots to all the elements of M), then you must also enjoy this cellfun function that applies an arbitrary function to each element of a cell array. Check it out!

### Noninteractive Matlab

The following command runs a matlab script under Linux in a batch mode:

jchen@peach \$ matlab -nodesktop -r my_script

Remember to put the exit command at the end of the script my_script.m to enforce Matlab to quit. This noninterative way is useful when you are submitting a job to a machine that you cannot interact with. Usually this machine is a super computer..

### MATLAB sparse matrix: vec==0 vs. length(find(vec))==0

Remember that in Matlab sparse matrices are stored column-wise. Hence it is much more efficient to do column operations than to do row operations on a sparse matrix. When possible, always have two copies of the matrix A: A and AT=A’. This is important for time critical experiments.

Now let’s restrict to column operations. How do we find zero columns of a sparse matrix? There are two (or many more?) possible ways: (1) Use the A(:,i) == 0 test; (2) Use the length(find(A(:,i))) == 0 test. It turns out that the former method is way slower than the latter. Try it! Use the following codes:

m = 10000; n = 10000; density = 1e-4; n_zero_col = 100;

 % generate sparse matrix A A = sprand(m, n, density); 

 % zero out some columns idx = ceil(rand(n_zero_col,1)*n); A(:,idx) = 0; 

 % two ways to find zero columns % way 1: the ==0' construct tic; zero_cols = []; for i = 1:n if A(:,i)==0 zero_cols = [zero_cols i]; end end toc 

 % way 2: the find' construct tic; zero_cols = []; for i = 1:n if length(find(A(:,i)))==0 zero_cols = [zero_cols i]; end end toc

### Optimize Matlab Codes

When computing the inner product of two vectors, a'*b is SIX times faster than dot(a,b). Turn on profile and try the follow script:

n = 1e3; iter = 1e5; a = rand(n,1); b = rand(n,1); for i = 1:iter a'*b; dot(a,b); end

• 261,714 hits