On Sat, Nov 28, 2009 at 11:04 AM, David Bateman
<address@hidden> wrote:
Jaroslav Hajek wrote:
I see, thanks. I randomized the guess again. My problem with the previous
code was that a priori forming A'*A is quite costly, especially if you want
just low precision and expect just a dozen iterations:
A = rand (1e3);
tic; [e,c] = normest (A); toc
tic; [e,c] = normest (A, 1e-1); toc
now:
Elapsed time is 0.0229969 seconds.
Elapsed time is 0.0183249 seconds.
previously:
Elapsed time is 0.220653 seconds.
Elapsed time is 0.210241 seconds.
Yes, I noticed that it starts to pay off in case the number of iterations gets comparable to the order of the matrix.
But then, normest times comparably to plain "norm". Besides, it's trivial to do `sqrt(normest(A'*A))', getting the speed-up back.
I suppose the timings used to be different, because prior to 3.2.x, A'*x was much slower, doing a physical matrix conjugate transpose.
In that case a similar issue probably exists in the code bicgstab where is two preconditioning matrices are given then
precon = @(x) M2 \ (M1 \ x);
if M1 and M2 are sparse or
M = M1 \ M2;
precon @(x) M \ x;
otherwise. Though bicgstab usually has many more iterations than a simple power method as in normest, so perhaps the existing code is the right ting to do..
It seems OK. Usually preconditioners are matrices with fast left division, so they're likely to be either sparse (permuted) triangular or diagonal matrices. Dense triangular preconditioners are unlikely. So, if M1 and M2 aren't both sparse, probably one is diagonal (or both are), in which case M1*M2 is efficient and does not adversely impact sparsity.