Thread Subject:

Huge time differences on inversion

Dear all,

I did a post yesterday, but no replies came and I wanted to repeat the
post. I need the inverse of a sparse 5544X5544 matrix, where I tried some
different options in MATLAB but the results are so different which raised
my attention, especially the one on inv. I used cholesky, lu
decompositions and direct inv in MATLAB,

here is the code to test the inversion timing

clc, clear all;
% test 1
load Mss;
tic;
invMss = inv(Mss);
toc;
% test2
clear Mss, invMss;
load Mss;
tic;
L = chol(Mss, 'lower');
invChol = L'\(L\ speye(size(Mss,1)));
toc;
%
clear Mss
clear invChol
clear L;
load Mss;
tic;
[L, U] = lu(Mss);
invLu = U\(L\speye(size(Mss,1)));
toc;
%
clear Mss;
clear invLu;
clear L;
clear U;
load Mss;
tic;
[L, U, P, Q] = lu(Mss);
invLu = U\(L\speye(size(Mss,1)));
toc;

And
 
Elapsed time is 13.613957 seconds.
Elapsed time is 142.979129 seconds.
Elapsed time is 221.484886 seconds.
Elapsed time is 52.666728 seconds.

Any comments are appreciated on these results, especially the huge time
differences.

Best regards,
Umut

On 23 apr, 14:49, utab <uta...@tudelft.nl> wrote:
> Dear all,
>
> I did a post yesterday, but no replies came and I wanted to repeat the
> post. I need the inverse of a sparse 5544X5544 matrix, where I tried some
> different options in MATLAB but the results are so different which raised
> my attention, especially the one on inv. I used cholesky, lu
> decompositions and direct inv in MATLAB,

> Any comments are appreciated on these results, especially the huge time
> differences.

INV has a number of methods available for computing an inverse.
One would assume that it contains some sort of preprocessor
that analyzes the matrix, and then selects the more appropriate
method.

The other methods you tried are specific methods that haev
different properties. Specifically, the Cholesky and LU
factorizations produce dense matrices in intermediate results,
even if you start out from a sparse matrix.

All in all, you have rediscovered, for the umpteenth time, that
different algorithms for doing the same thing have different
properties and different scopes of use.

Rune

"utab" <utabak@tudelft.nl> wrote in message
news:1f4cd$4bd19764$82a112b1$31701@news1.tudelft.nl...
> Dear all,
>
> I did a post yesterday, but no replies came and I wanted to repeat the
> post. I need the inverse of a sparse 5544X5544 matrix,

DON'T DO THIS.

If you're trying to use the inverse to solve a system of equations, use
backslash instead.

Also remember that the inverse of a sparse matrix may be much more densely
populated than the original sparse matrix, even though it will be stored as
a sparse matrix.

n = 100;
S = speye(n)+sparse([ones(n, 1) zeros(n, n-1)]);
S = S+S';
Sinv = inv(S);
fprintf('S contains %d nonzeros out of %d elements for a density of %g.\n',
...
    nnz(S), numel(S), nnz(S)/numel(S))
fprintf('Sinv contains %d nonzeros out of %d elements for a density of
%g.\n', ...
    nnz(Sinv), numel(Sinv), nnz(Sinv)/numel(Sinv))

In fact, storing a densely populated matrix as a sparse matrix can require
_more_ memory than storing that same matrix as a full matrix!

SinvF = full(Sinv);
whos S Sinv SinvF

So again, do not do this unless you absolutely, positively, cannot proceed
without the inverse matrix for whatever reason -- and even then, if you
describe the problem that you're trying to solve in a reply to the
newsgroup, someone may be able to find you an INV-less alternative.

--
Steve Lord
slord@mathworks.com
comp.soft-sys.matlab (CSSM) FAQ: http://matlabwiki.mathworks.com/MATLAB_FAQ

On Fri, 23 Apr 2010 09:50:00 -0400, Steven Lord wrote:

> "utab" <utabak@tudelft.nl> wrote in message
> news:1f4cd$4bd19764$82a112b1$31701@news1.tudelft.nl...
>> Dear all,
>>
>> I did a post yesterday, but no replies came and I wanted to repeat the
>> post. I need the inverse of a sparse 5544X5544 matrix,
>
> DON'T DO THIS.
>
> If you're trying to use the inverse to solve a system of equations, use
> backslash instead.
>
> Also remember that the inverse of a sparse matrix may be much more
> densely populated than the original sparse matrix, even though it will
> be stored as a sparse matrix.
>
> n = 100;
> S = speye(n)+sparse([ones(n, 1) zeros(n, n-1)]); S = S+S';
> Sinv = inv(S);
> fprintf('S contains %d nonzeros out of %d elements for a density of
> %g.\n', ...
> nnz(S), numel(S), nnz(S)/numel(S))
> fprintf('Sinv contains %d nonzeros out of %d elements for a density of
> %g.\n', ...
> nnz(Sinv), numel(Sinv), nnz(Sinv)/numel(Sinv))
>
> In fact, storing a densely populated matrix as a sparse matrix can
> require _more_ memory than storing that same matrix as a full matrix!
>
> SinvF = full(Sinv);
> whos S Sinv SinvF
>
> So again, do not do this unless you absolutely, positively, cannot
> proceed without the inverse matrix for whatever reason -- and even then,
> if you describe the problem that you're trying to solve in a reply to
> the newsgroup, someone may be able to find you an INV-less alternative.

Hi Steve,

Thanks for the detailed information which is not new to me. I am not
trying to solve a linear system first of all. On the backslash issue, I
agree with you. I have to form a symmetrization matrix for a problem that
I am working on for a long time(not explicitly although but I need the
inverse of that matrix, after that explicitly I can get a format where I
do not need to construct the full symmetrization matrix but the inverse
of the matrix mentioned before). I know that these kinds of inverses
should be avoided as far as possible. However I need this inverse to make
my algorithm to work to symmetrize my original unsymmetric system
matrices. That is the source of the problem in brief.

Seems that I should device sth other than the inverse of the matrix or
find another way.

Best regards,
Umut

All the methods you post are not the best way to do it. You need to let chol (or lu) permute the matrix to reduce fill-in. If your run-time is much higher than x=A\b, then you are doing the wrong thing.

See the "FACTORIZE" object and the methods it contains for the best way to do this, for the sparse/dense cases, and the symmetric/unsymmetric/rectangular cases:

http://www.mathworks.com/matlabcentral/fileexchange/24119