SVDS-- what is the point of using it and is it ever faster than SVD??

My program finds the eigenvalues of a bunch of matrices and adds them together. But the matrices can get really big and sparse, and this adds to my computation time. I'm using the function SVD to get the eigenvalues, so I tried using the sparse version, SVDS, since I only need the biggest eigenvalues anyway. This wasn't any faster. In fact, I think it's SLOWER! If this is the case, I don't see why SVDS exists. I need to find a sparse matrix operation that will reduce my computation time. Thanks!

[Jan]

What exactly is "really big" and do you use the sparse type for the matrix or does it only contains a lot of zeros?

[Clare]

Really big means millions by millions. I use the sparse type for the matrix, but still it doesn't seem to speed things up.

[dpb]

Well, while sparse will save storage, there's always an overhead associated w/ the storage scheme of retrieving values and there's nothing that guarantees a performance improvement for any particular arrangement of sparsity while for some other it might make a large improvement.

Only an analysis of the specifics of the sparsity in a particular problem space might lead to an improved handling over the general-purpose case that Matlab implements.

[Clare]

So you're saying if my matrix has more nonzero elements, SVDS could be faster?

[Muthu Annamalai]

It seems to me one woman's food is another man's poison, or no one-size fits all :-)

[dpb]

If it were more sparse or if there was a different regularity in the existing sparsity either could make a difference, yes. All in all, it'd probably be more likely to make a difference if there were some recognizable/usable pattern in the sparsity over just numbers. What, if any, analysis SVDS does internally on that I don't know.

[Clare]

When you say usable pattern, do you mean usable to do other operations? I'm wondering if there is a workaround to use instead of SVDS.

[Jan]

@Clare: When the programmer knows, that a matrix is a block-sparse matrix with all non-zero elements in square blocks of a growing size at the main diagonal, the SVD can be calculated much faster. It is like the application of a Gauss decomposition for an upper triangle matrix. Then you do not have to search for pivot elements under the main diagonal.

An algorithm can be much faster when a known pattern can be exploited, but the programming and testing of an individually adjusted solver will take much time. So it is worth to check, if you need 2 weeks for programming only to save 4 seconds of run time.

Answer 1

[Matt J]

Really big means millions by millions. I use the sparse type for the matrix, but still it doesn't seem to speed things up.

I assume the comparison between SVD and SVDS that you're performing is on a scaled down version of this matrix. There is no way a "millions by millions" matrix can be stored in anything but sparse type and there is no way you can be running SVD on a sparse type matrix (it would return an error).

If your tests are based on small scale versions of your actual matrix, we need to ask how small. Definitely the gains in SVDS over SVD won't be seen if your examples are too small, e.g.,

>> N=100; As=speye(N);A=eye(N); tic; svd(A); toc ; tic; svds(A); toc
Elapsed time is 0.001177 seconds.
Elapsed time is 0.006077 seconds.

but there are clearly also examples where svds is better, e.g.,

>> N=5000; As=speye(N);A=eye(N); tic; svd(A); toc ; tic; svds(A); toc
Elapsed time is 47.876229 seconds.
Elapsed time is 0.324260 seconds.

[Clare]

Yes, I meant that they will eventually be millions by millions, thanks. Good to see that SVDS is clearly faster with the identity matrix.

Answer 2

[Jan]

"Millions" means something >= 2e6. Then the array has 32 TB when it is full. Do you have enough memory for this?

I tried this under Matlab 2009a/64/Win7/4 GB:

x = rand(1e4);
x(rand(size(x))<0.9) = 0;
tic; S = svd(x); toc
% 1313.510373 sec

x = sparse(x);
tic; S = svds(x); toc
% 54.953447 sec

This looks like a good argument for SVDS. Now tell us please something about the sparsity of your array and the number of eigenvalues you need. And please explain, which timings you expect for [2e6 x 2e6] matrices.

[Jan]

I haven't seen Matt's answer when I started to create mine. I've read other message while the computer was busy for 1300 seconds. The SVD and SVDS are not multi-threaded in 2009a, so the processor load was 50% only on my Dual-Core.

[Clare]

Thanks for the information and testing. I will eventually have 0.5 TB to work with. But of course I'd like to speed up things by only working with the nonzero elements of my huge matrices.

But here is my question-- if I use S = sparse(i,j,s) to construct my sparse matrix, and the sparsity is around 0.3 (meaning 30% nonzero elements, 70% zeros), will I save any memory? For each element I must store its i and j coordinates, which multiplies the storage by 3. Isn't this the same as storing the entire matrix if it's 0.3 sparse? Thanks again!

[Jan]

A good question. The sparse matrix stores the data and the indices in three different blocks, while the full matrix requires one large contiguous block of free memory. Therefore full matrices cause Out-Of-Memory errors earlier.

If the indices of sparse arrays are stored as UINT32, less memory is consumed. But I do not know the type of the indices.

[Cedric]

A very good question indeed. A sparsity of 0.3 is almost not sparse (you can find a document that describes the internals here).. and in R2012b and earlier versions at least, we cannot use UINT32 as indices.

[Matt J]

But here is my question-- if I use S = sparse(i,j,s) to construct my sparse matrix, and the sparsity is around 0.3 (meaning 30% nonzero elements, 70% zeros), will I save any memory?

@Clare - even if your sparsity is zero, it is still possible to waste memory as the following example shows,

 >> Asparse=spalloc(1,1e6,1); 
 >> Afull=full(Asparse);

 >> whos Asparse Afull

Name         Size                   Bytes  Class     Attributes

Afull        1x1000000            8000000  double              
Asparse      1x1000000            8000024  double    sparse