Mercurial > octave-nkf
annotate scripts/linear-algebra/normest.m @ 9871:87595d714005
move normest to linear-algebra, improve it
author | Jaroslav Hajek <highegg@gmail.com> |
---|---|
date | Thu, 26 Nov 2009 07:22:02 +0100 |
parents | scripts/sparse/normest.m@1bf0ce0930be |
children | 47f36dd27203 |
rev | line source |
---|---|
8920 | 1 ## Copyright (C) 2006, 2007, 2008, 2009 David Bateman and Marco Caliari |
9871
87595d714005
move normest to linear-algebra, improve it
Jaroslav Hajek <highegg@gmail.com>
parents:
9051
diff
changeset
|
2 ## Copyright (C) 2009 VZLU Prague |
6212 | 3 ## |
7016 | 4 ## This file is part of Octave. |
6212 | 5 ## |
7016 | 6 ## Octave is free software; you can redistribute it and/or modify it |
7 ## under the terms of the GNU General Public License as published by | |
8 ## the Free Software Foundation; either version 3 of the License, or (at | |
9 ## your option) any later version. | |
10 ## | |
11 ## Octave is distributed in the hope that it will be useful, but | |
12 ## WITHOUT ANY WARRANTY; without even the implied warranty of | |
13 ## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU | |
14 ## General Public License for more details. | |
6212 | 15 ## |
16 ## You should have received a copy of the GNU General Public License | |
7016 | 17 ## along with Octave; see the file COPYING. If not, see |
18 ## <http://www.gnu.org/licenses/>. | |
6212 | 19 |
20 ## -*- texinfo -*- | |
21 ## @deftypefn {Function File} {[@var{n}, @var{c}] =} normest (@var{a}, @var{tol}) | |
6222 | 22 ## Estimate the 2-norm of the matrix @var{a} using a power series |
9051
1bf0ce0930be
Grammar check TexInfo in all .m files
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
23 ## analysis. This is typically used for large matrices, where the cost |
6222 | 24 ## of calculating the @code{norm (@var{a})} is prohibitive and an approximation |
6212 | 25 ## to the 2-norm is acceptable. |
26 ## | |
9051
1bf0ce0930be
Grammar check TexInfo in all .m files
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
27 ## @var{tol} is the tolerance to which the 2-norm is calculated. By default |
1bf0ce0930be
Grammar check TexInfo in all .m files
Rik <rdrider0-list@yahoo.com>
parents:
8920
diff
changeset
|
28 ## @var{tol} is 1e-6. @var{c} returns the number of iterations needed for |
6212 | 29 ## @code{normest} to converge. |
30 ## @end deftypefn | |
31 | |
9871
87595d714005
move normest to linear-algebra, improve it
Jaroslav Hajek <highegg@gmail.com>
parents:
9051
diff
changeset
|
32 function [e, c] = normest (A, tol = 1e-6) |
87595d714005
move normest to linear-algebra, improve it
Jaroslav Hajek <highegg@gmail.com>
parents:
9051
diff
changeset
|
33 if (! ismatrix (A) || ndims (A) != 2) |
6498 | 34 error ("normest: A must be a matrix"); |
6212 | 35 endif |
9871
87595d714005
move normest to linear-algebra, improve it
Jaroslav Hajek <highegg@gmail.com>
parents:
9051
diff
changeset
|
36 if (! isfloat (A)) |
87595d714005
move normest to linear-algebra, improve it
Jaroslav Hajek <highegg@gmail.com>
parents:
9051
diff
changeset
|
37 A = double (A); |
87595d714005
move normest to linear-algebra, improve it
Jaroslav Hajek <highegg@gmail.com>
parents:
9051
diff
changeset
|
38 endif |
87595d714005
move normest to linear-algebra, improve it
Jaroslav Hajek <highegg@gmail.com>
parents:
9051
diff
changeset
|
39 tol = max (tol, eps (class (A))); |
6212 | 40 c = 0; |
9871
87595d714005
move normest to linear-algebra, improve it
Jaroslav Hajek <highegg@gmail.com>
parents:
9051
diff
changeset
|
41 x = norm (A, "columns").'; |
87595d714005
move normest to linear-algebra, improve it
Jaroslav Hajek <highegg@gmail.com>
parents:
9051
diff
changeset
|
42 e = norm (x); |
87595d714005
move normest to linear-algebra, improve it
Jaroslav Hajek <highegg@gmail.com>
parents:
9051
diff
changeset
|
43 if (e > 0) |
6498 | 44 [m, n] = size (A); |
9871
87595d714005
move normest to linear-algebra, improve it
Jaroslav Hajek <highegg@gmail.com>
parents:
9051
diff
changeset
|
45 x /= e; |
6212 | 46 e0 = 0; |
9871
87595d714005
move normest to linear-algebra, improve it
Jaroslav Hajek <highegg@gmail.com>
parents:
9051
diff
changeset
|
47 while (abs (e - e0) > tol * e) |
87595d714005
move normest to linear-algebra, improve it
Jaroslav Hajek <highegg@gmail.com>
parents:
9051
diff
changeset
|
48 e0 = e; |
87595d714005
move normest to linear-algebra, improve it
Jaroslav Hajek <highegg@gmail.com>
parents:
9051
diff
changeset
|
49 y = A*x; |
87595d714005
move normest to linear-algebra, improve it
Jaroslav Hajek <highegg@gmail.com>
parents:
9051
diff
changeset
|
50 e = norm (y); |
87595d714005
move normest to linear-algebra, improve it
Jaroslav Hajek <highegg@gmail.com>
parents:
9051
diff
changeset
|
51 if (e == 0) |
87595d714005
move normest to linear-algebra, improve it
Jaroslav Hajek <highegg@gmail.com>
parents:
9051
diff
changeset
|
52 x = rand (n, 1); |
87595d714005
move normest to linear-algebra, improve it
Jaroslav Hajek <highegg@gmail.com>
parents:
9051
diff
changeset
|
53 else |
87595d714005
move normest to linear-algebra, improve it
Jaroslav Hajek <highegg@gmail.com>
parents:
9051
diff
changeset
|
54 x = A'*(y / e); |
87595d714005
move normest to linear-algebra, improve it
Jaroslav Hajek <highegg@gmail.com>
parents:
9051
diff
changeset
|
55 endif |
87595d714005
move normest to linear-algebra, improve it
Jaroslav Hajek <highegg@gmail.com>
parents:
9051
diff
changeset
|
56 x /= norm (x); |
87595d714005
move normest to linear-algebra, improve it
Jaroslav Hajek <highegg@gmail.com>
parents:
9051
diff
changeset
|
57 c += 1; |
87595d714005
move normest to linear-algebra, improve it
Jaroslav Hajek <highegg@gmail.com>
parents:
9051
diff
changeset
|
58 endwhile |
6212 | 59 endif |
60 endfunction |