No BSD License  

Highlights from
Fast sparse matrix vector product

3.8 | 5 ratings Rate this file 17 Downloads (last 30 days) File Size: 3.15 KB File ID: #10556

Fast sparse matrix vector product


Darren Engwirda (view profile)


28 Mar 2006 (Updated )

sparse matrix vector product

| Watch this File

File Information

This is a dedicated mex-file to efficiently compute the Sparse-Matrix-Vector-Product (smvp).

x = smvp(A,b);

This is generally 4 times faster than x = A*b for large A (On my machine).

Type "help smvp" and "smvp_demo" for more information.

NOTE: You will need either windows or the MATLAB compiler to use this function. A pre-compiled version is included, but you can re-compile via:

mex -O smvp.c

You may also find my efficient LU substitution function useful:

P.S This again can be useful in the numerical solution of PDE's.

NOTE: The zip file should now work (03/04/2006) but I have removed the demo

MATLAB release MATLAB 6.5 (R13)
Tags for This File   Please login to tag files.
Please login to add a comment or rating.
Comments and Ratings (7)
21 Jul 2012 Christopher Kanan

I'm not sure why, but it is much slower on R2012a. I ran the example and the inbuilt required 0.01 s and this version required 10 seconds.

I doubt MATLAB could improve that much, so I'm guessing some sort of optimization used in the mex code is no longer functional.

13 Mar 2009 Stephen Becker

Stephen Becker (view profile)

There's actually a very important point here. Because of how Matlab stores a sparse matrix (see Tim Davis' book, or read the help on, say, mxGetIc), and because a major bottleneck for these computations is loading data into the CPU's cache, it is MUCH faster to store a matrix by row, instead of column, when doing matrix-vector multiples.

In Matlab, this is easy to remedy: simply take the transpose.

So, this cool trick can save you time. If you want to do this repeatedly:
>> y = A*b
instead, do this:
>> At = A.'; % a one-time cost
>> y = At.'b

But, this trick only works on newer versions of Matlab. On, say, R2006b, a call like
>> y = At.'b
will actually calculate the transpose of At, so this is very slow!

On older versions of Matlab, it's likely that this SMVP code will really help, but it might not be as helpful on R2008, for example.

On linux, non-multithreaded R2006a, I found that SMVP took about 60% the time of Matlab's own multiply. So, good work!

21 Feb 2007 Tim Davis

Regarding the comment on checking for NULL return of mxCreateDoubleMatrix: it's not necessary. If the allocation fails, MATLAB will automatically terminate the mexFunction and clean up any allocated memory. If it was malloc, then yes, you would need to check for NULL. This is a fine example of how to use mxCreateDoubleMatrix.

The only design issue that should be corrected is the "nlhs!=1" test (a minor point). If you simply do

>> smvp(A,b)

with no output arguments, then nlhs is set to zero. You can still write to plhs[0], and the result will go in the "ans" variable.

Comment only
20 Feb 2007 Martin S.

Works but offers almost very minor 8% speedup as of Matlab 2006.
I multiplied 700,000 x 700,000 matrix with 2,000,000 nonzero entries times a 700,000 vector and got only 8% faster than Matlab 2006. Not worth using unless there is more speedup in other settings.

2) The routine fails to check that the matrix allociation succeeded (may run out of memory). Need to add a NULL check to

plhs[0] = mxCreateDoubleMatrix(mxGetM(prhs[0]), 1, mxREAL);

3) Change if (nlhs!=1) to if (nlhs > 1)

Summary: works but very limited speedup.

20 Dec 2006 Shelley Lee

It would be excellent if it can support
complex arithmetics.

30 Apr 2006 Robert Beardmore

I had no problems using this file and it is around 5 times quicker than matlab's A*b for the problems I tried.

01 Apr 2006 Michal Kvasnicka

The zip file is corrupted. The matrix.mat file can not be extracted.

Comment only

Contact us