Mercurial > hg > octave-avbm
changeset 6953:4567a35e0777
[project @ 2007-10-04 02:53:11 by jwe]
author | jwe |
---|---|
date | Thu, 04 Oct 2007 02:54:00 +0000 |
parents | 09a89fb42c09 |
children | 9dabcb305dda |
files | scripts/ChangeLog scripts/linear-algebra/__norm__.m scripts/linear-algebra/norm.m src/ChangeLog src/data.cc |
diffstat | 5 files changed, 231 insertions(+), 211 deletions(-) [+] |
line wrap: on
line diff
--- a/scripts/ChangeLog +++ b/scripts/ChangeLog @@ -1,3 +1,9 @@ +2007-10-03 John W. Eaton <jwe@octave.org> + + * linear-algebra/Makefile.in (SOURCES): Rename norm.m to __norm__.m. + * linear-algebra/__norm__.m: Rename from norm.m. Eliminate + special for __vnorm__. + 2007-10-03 Quentin Spencer <qspencer@ieee.org> * linear-algebra/norm.m: Special case vector 1-norm and 2-norm.
new file mode 100644 --- /dev/null +++ b/scripts/linear-algebra/__norm__.m @@ -0,0 +1,93 @@ +## Copyright (C) 1996, 1997 John W. Eaton +## +## This file is part of Octave. +## +## Octave is free software; you can redistribute it and/or modify it +## under the terms of the GNU General Public License as published by +## the Free Software Foundation; either version 2, or (at your option) +## any later version. +## +## Octave is distributed in the hope that it will be useful, but +## WITHOUT ANY WARRANTY; without even the implied warranty of +## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +## General Public License for more details. +## +## You should have received a copy of the GNU General Public License +## along with Octave; see the file COPYING. If not, write to the Free +## Software Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA +## 02110-1301, USA. + +## Undocumented internal function. + +## Author: jwe + +function retval = __norm__ (x, p) + + if (nargin < 1 || nargin > 2) + print_usage (); + endif + + if (isempty (x)) + retval = []; + return; + endif + + if (ndims (x) > 2) + error ("norm: only valid on 2-D objects") + endif + + if (nargin == 1) + p = 2; + endif + + ## Do we have a vector or matrix as the first argument? + + if (ndims(x) == 2 && (rows (x) == 1 || columns (x) == 1)) + if (ischar (p)) + if (strcmp (p, "fro")) + retval = sqrt (sum (abs (x) .^ 2)); + elseif (strcmp (p, "inf")) + retval = max (abs (x)); + else + error ("norm: unrecognized norm"); + endif + else + if (p == Inf) + retval = max (abs (x)); + elseif (p == -Inf) + retval = min (abs (x)); + elseif (p == 1) + retval = sum (abs (x)); + elseif (p == 2) + if (iscomplex (x)) + x = abs (x); + endif + retval = sqrt (sumsq (x)); + else + retval = sum (abs (x) .^ p) ^ (1/p); + endif + endif + else + if (ischar (p)) + if (strcmp (p, "fro")) + retval = sqrt (sum (sum (abs (x) .^ 2))); + elseif (strcmp (p, "inf")) + retval = max (sum (abs (x'))); + else + error ("norm: unrecognized vector norm"); + endif + else + if (p == 1) + retval = max (sum (abs (x))); + elseif (p == 2) + s = svd (x); + retval = s (1); + elseif (p == Inf) + retval = max (sum (abs (x'))); + else + error ("norm: unrecognized matrix norm"); + endif + endif + endif + +endfunction
deleted file mode 100644 --- a/scripts/linear-algebra/norm.m +++ /dev/null @@ -1,152 +0,0 @@ -## Copyright (C) 1996, 1997 John W. Eaton -## -## This file is part of Octave. -## -## Octave is free software; you can redistribute it and/or modify it -## under the terms of the GNU General Public License as published by -## the Free Software Foundation; either version 2, or (at your option) -## any later version. -## -## Octave is distributed in the hope that it will be useful, but -## WITHOUT ANY WARRANTY; without even the implied warranty of -## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU -## General Public License for more details. -## -## You should have received a copy of the GNU General Public License -## along with Octave; see the file COPYING. If not, write to the Free -## Software Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA -## 02110-1301, USA. - -## -*- texinfo -*- -## @deftypefn {Function File} {} norm (@var{a}, @var{p}) -## Compute the p-norm of the matrix @var{a}. If the second argument is -## missing, @code{p = 2} is assumed. -## -## If @var{a} is a matrix: -## -## @table @asis -## @item @var{p} = @code{1} -## 1-norm, the largest column sum of the absolute values of @var{a}. -## -## @item @var{p} = @code{2} -## Largest singular value of @var{a}. -## -## @item @var{p} = @code{Inf} -## @cindex infinity norm -## Infinity norm, the largest row sum of the absolute values of @var{a}. -## -## @item @var{p} = @code{"fro"} -## @cindex Frobenius norm -## Frobenius norm of @var{a}, @code{sqrt (sum (diag (@var{a}' * @var{a})))}. -## @end table -## -## If @var{a} is a vector or a scalar: -## -## @table @asis -## @item @var{p} = @code{Inf} -## @code{max (abs (@var{a}))}. -## -## @item @var{p} = @code{-Inf} -## @code{min (abs (@var{a}))}. -## -## @item other -## p-norm of @var{a}, @code{(sum (abs (@var{a}) .^ @var{p})) ^ (1/@var{p})}. -## @end table -## @seealso{cond, svd} -## @end deftypefn - -## Author: jwe - -function retval = norm (x, p) - - if (nargin < 1 || nargin > 2) - print_usage (); - endif - - if (isempty (x)) - retval = []; - return; - endif - - if (ndims (x) > 2) - error ("norm: only valid on 2-D objects") - endif - - if (nargin == 1) - p = 2; - endif - - ## Do we have a vector or matrix as the first argument? - - if (ndims(x) == 2 && (rows (x) == 1 || columns (x) == 1)) - if (isinteger (x) || issparse (x)) - if (ischar (p)) - if (strcmp (p, "fro")) - retval = sqrt (sum (abs (x) .^ 2)); - elseif (strcmp (p, "inf")) - retval = max (abs (x)); - else - error ("norm: unrecognized norm"); - endif - else - if (p == Inf) - retval = max (abs (x)); - elseif (p == -Inf) - retval = min (abs (x)); - elseif (p == 1) - retval = sum (abs (x)); - elseif (p == 2) - if (iscomplex (x)) - x = abs (x); - endif - retval = sqrt (sum (x .* x)); - else - retval = sum (abs (x) .^ p) ^ (1/p); - endif - endif - else - retval = __vnorm__ (x, p); - endif - else - if (ischar (p)) - if (strcmp (p, "fro")) - retval = sqrt (sum (sum (abs (x) .^ 2))); - elseif (strcmp (p, "inf")) - retval = max (sum (abs (x'))); - else - error ("norm: unrecognized vector norm"); - endif - else - if (p == 1) - retval = max (sum (abs (x))); - elseif (p == 2) - s = svd (x); - retval = s (1); - elseif (p == Inf) - retval = max (sum (abs (x'))); - else - error ("norm: unrecognized matrix norm"); - endif - endif - endif - -endfunction - -%!shared x -%! x = [1, -3, 4, 5, -7]; -%!assert(norm(x,1), 20); -%!assert(norm(x,2), 10); -%!assert(norm(x,3), 8.24257059961711, -4*eps); -%!assert(norm(x,Inf), 7); -%!assert(norm(x,-Inf), 1); -%!assert(norm(x,"inf"), 7); -%!assert(norm(x,"fro"), 10); -%!assert(norm(x), 10); -%!assert(norm([1e200, 1]), 1e200); -%!assert(norm([3+4i, 3-4i, sqrt(31)]), 9, -4*eps); -%!shared m -%! m = magic (4); -%!assert(norm(m,1), 34); -%!assert(norm(m,2), 34); -%!assert(norm(m,Inf), 34); -%!assert(norm(m,"inf"), 34);
--- a/src/ChangeLog +++ b/src/ChangeLog @@ -1,3 +1,8 @@ +2007-10-03 John W. Eaton <jwe@octave.org> + + * data.cc (Fnorm): New function. + (F__vnorm__): Delete. + 2007-10-03 Michael Goffioul <michael.goffioul@gmail.com> * DLD-FUNCTIONS/typecast.cc: Include <algorithm>.
--- a/src/data.cc +++ b/src/data.cc @@ -34,18 +34,19 @@ #include "str-vec.h" #include "quit.h" +#include "Cell.h" #include "defun.h" #include "error.h" #include "gripes.h" +#include "oct-map.h" +#include "oct-obj.h" #include "ov.h" #include "ov-complex.h" #include "ov-cx-mat.h" -#include "variables.h" -#include "oct-obj.h" +#include "parse.h" +#include "pt-mat.h" #include "utils.h" -#include "Cell.h" -#include "oct-map.h" -#include "pt-mat.h" +#include "variables.h" #define ANY_ALL(FCN) \ \ @@ -2614,74 +2615,141 @@ return retval; } +/* +%!shared x +%! x = [1, -3, 4, 5, -7]; +%!assert(norm(x,1), 20); +%!assert(norm(x,2), 10); +%!assert(norm(x,3), 8.24257059961711, -4*eps); +%!assert(norm(x,Inf), 7); +%!assert(norm(x,-Inf), 1); +%!assert(norm(x,"inf"), 7); +%!assert(norm(x,"fro"), 10); +%!assert(norm(x), 10); +%!assert(norm([1e200, 1]), 1e200); +%!assert(norm([3+4i, 3-4i, sqrt(31)]), 9, -4*eps); +%!shared m +%! m = magic (4); +%!assert(norm(m,1), 34); +%!assert(norm(m,2), 34); +%!assert(norm(m,Inf), 34); +%!assert(norm(m,"inf"), 34); +*/ + // Compute various norms of the vector X. -DEFUN (__vnorm__, args, , +DEFUN (norm, args, , "-*- texinfo -*-\n\ -@deftypefn {Built-in Function} {} __vnorm__ (@var{x}, @var{p})\n\ -Undocumented internal function.\n\ +@deftypefn {Function File} {} norm (@var{a}, @var{p})\n\ +Compute the p-norm of the matrix @var{a}. If the second argument is\n\ +missing, @code{p = 2} is assumed.\n\ +\n\ +If @var{a} is a matrix:\n\ +\n\ +@table @asis\n\ +@item @var{p} = @code{1}\n\ +1-norm, the largest column sum of the absolute values of @var{a}.\n\ +\n\ +@item @var{p} = @code{2}\n\ +Largest singular value of @var{a}.\n\ +\n\ +@item @var{p} = @code{Inf}\n\ +@cindex infinity norm\n\ +Infinity norm, the largest row sum of the absolute values of @var{a}.\n\ +\n\ +@item @var{p} = @code{\"fro\"}\n\ +@cindex Frobenius norm\n\ +Frobenius norm of @var{a}, @code{sqrt (sum (diag (@var{a}' * @var{a})))}.\n\ +@end table\n\ +\n\ +If @var{a} is a vector or a scalar:\n\ +\n\ +@table @asis\n\ +@item @var{p} = @code{Inf}\n\ +@code{max (abs (@var{a}))}.\n\ +\n\ +@item @var{p} = @code{-Inf}\n\ +@code{min (abs (@var{a}))}.\n\ +\n\ +@item other\n\ +p-norm of @var{a}, @code{(sum (abs (@var{a}) .^ @var{p})) ^ (1/@var{p})}.\n\ +@end table\n\ +@seealso{cond, svd}\n\ @end deftypefn") { - octave_value retval; + // Currently only handles vector norms for full double/complex + // vectors internally. Other cases are handled by __norm__.m. + + octave_value_list retval; int nargin = args.length (); if (nargin == 1 || nargin == 2) { - double p_val; - - octave_value p_arg; - - if (nargin == 1) - p_arg = 2; - else - p_arg = args(1); - - if (p_arg.is_string ()) + octave_value x_arg = args(0); + + if (x_arg.is_empty ()) + retval(0) = 0.0; + else if (x_arg.ndims () == 2) { - std::string p = args(1).string_value (); - - if (p == "inf") - p_val = octave_Inf; - else if (p == "fro") - p_val = -1; + if ((x_arg.rows () == 1 || x_arg.columns () == 1) + && ! (x_arg.is_sparse_type () || x_arg.is_integer_type ())) + { + double p_val; + + octave_value p_arg; + + if (nargin == 1) + p_arg = 2; + else + p_arg = args(1); + + if (p_arg.is_string ()) + { + std::string p = args(1).string_value (); + + if (p == "inf") + p_val = octave_Inf; + else if (p == "fro") + p_val = -1; + else + error ("norm: unrecognized norm `%s'", p.c_str ()); + } + else + { + p_val = p_arg.double_value (); + + if (error_state) + error ("norm: unrecognized norm value"); + } + + if (! error_state) + { + if (x_arg.is_real_type ()) + { + MArray<double> x (x_arg.array_value ()); + + if (! error_state) + retval(0) = x.norm (p_val); + else + error ("norm: expecting real vector"); + } + else + { + MArray<Complex> x (x_arg.complex_array_value ()); + + if (! error_state) + retval(0) = x.norm (p_val); + else + error ("norm: expecting complex vector"); + } + } + } else - { - error ("norm: unrecognized norm `%s'", p.c_str ()); - return retval; - } + retval = feval ("__norm__", args); } else - { - p_val = args(1).double_value (); - - if (error_state) - { - error ("norm: unrecognized norm value"); - return retval; - } - } - - octave_value x_arg = args(0); - - if (x_arg.is_real_type ()) - { - ColumnVector x (x_arg.vector_value ()); - - if (! error_state) - retval = x.norm (p_val); - else - error ("norm: expecting real vector"); - } - else - { - ComplexColumnVector x (x_arg.complex_vector_value ()); - - if (! error_state) - retval = x.norm (p_val); - else - error ("norm: expecting complex vector"); - } + error ("norm: only valid for 2-D objects"); } else print_usage ();