Re: chol2inv for Octave?

[David Bateman]

By solve I mean the function in the Matrix and ComplexMatrix classes 
that is behind the "\" and "/" operators. In most cases there is never 
any need to calculate the inverse (I can think of a few counter examples 
though LDPC codes come to mind), as the "\" and "/" can treat multiple 
RHS by taking a matrix argument for the RHS. So ask yourself if you 
really need inv... As a
example with the octave-fore chol function installed

octave:1> a = randn(1024,1024) + 100* eye(1024,1024); a= a+a';
octave:2> tic; x1= a \ ones(1024,1); toc
ans = 0.94385
octave:3> tic; r = chol(a); x2 = r \ (r' \ ones(1024,1)); toc
ans = 0.63005
octave:4> r = chol(a); tic; x2 = r \ (r' \ ones(1024,1)); toc
ans = 0.21868
octave:5> ainv = inv(a); tic; x1 = ainv * ones(1024,1); toc
ans = 0.015854
octave:5> norm(x2 - x1)
ans =  2.9087e-16

This shows that chol is faster than "\" directly, and that even if you 
wanted to use multiple RHS but not use a matrix argument, the back 
subsitution part is relatively fast, but not as fast as "ainv * b". This 
speed difference is mainly the
fault of the home-grown back substitution code used by the octave-forge 
"chol" stuff.

You could do something like "R \ (R' \ I)" to get the inverse from the 
cholesky factoruzation, except that is a very inefficient means of 
forming the inverse of a matrix.

As for implementing different versions of inv, I hadn't intended to do 
that since as previously noted its more important that "\" and "/" are 
efficient. Although the factorization in "\" is the same as for inv, 
what's done with it after that isn't, so sharing the code between "\" 
and inv would require some reorganization.....

D.



-------------------------------------------------------------
Octave is freely available under the terms of the GNU GPL.

Octave's home on the web:  http://www.octave.org
How to fund new projects:  http://www.octave.org/funding.html
Subscription information:  http://www.octave.org/archive.html
-------------------------------------------------------------