[Octave-bug-tracker] [bug #29768] logm fix

[Jaroslav Hajek]

Follow-up Comment #5, bug #29768 (project octave):

I checked in the following patch:
http://hg.savannah.gnu.org/hgweb/octave/rev/f9860b622680

Summary of changes:

I wrote an in-place column-major (i.e. cache-aligned) version of the triple
loop from Higham's paper for computing the upper triangular sqrtm. The
in-place version reduces memory usage when applied to temporaries (such as the
result of a real->complex conversion).

For diagonal matrices, sqrtm will simply reduce to sqrt.
For upper triangular matrices, the triple loop will be used directly. Lower
triangular matrices are handled via double transpose. If the matrix has a
positive diagonal, the sqrtm will be done in real arithmetic.

Here goes a benchmark:

n = 500;

disp ("diagonal matrix");
x = diag (rand (n, 1)); x = x*x;
tic; sqrtm (x); toc

disp ("upper triangular matrix");
x = triu (rand (n) + eye (n)); x = x * x;
tic; sqrtm (x); toc

disp ("lower triangular matrix");
x = tril (rand (n) + eye (n)); x = x * x;
tic; sqrtm (x); toc

disp ("full matrix");
x = rand (n) + eye (n); x = x * x;
tic; sqrtm (x); toc

previously, I got:

diagonal matrix
Elapsed time is 0.562037 seconds.
upper triangular matrix
Elapsed time is 0.554102 seconds.
lower triangular matrix
Elapsed time is 0.577609 seconds.
full matrix
Elapsed time is 1.78512 seconds.

with the new patch, I get:

diagonal matrix
Elapsed time is 0.000648975 seconds.
upper triangular matrix
Elapsed time is 0.0239501 seconds.
lower triangular matrix
Elapsed time is 0.026382 seconds.
full matrix
Elapsed time is 1.54298 seconds.

we still don't have the error estimating code (Matlab has it).

    _______________________________________________________

Reply to this item at:

  <http://savannah.gnu.org/bugs/?29768>

_______________________________________________
  Message sent via/by Savannah
  http://savannah.gnu.org/