# HG changeset patch # User Marco Caliari # Date 1509960429 -3600 # Node ID 3d96400df7132e98a63417017c1fc706515e9565 # Parent 3d26a174e143a9ceceb62d57ffac7a8ffa94d14f Add function vecnorm (bug #52342). * scripts/linear-algebra/vecnorm.m: New function. * script/linear-algebra/module.mk: Add vecnorm function to build system. * NEWS: Announce new vecnorm function. * linalg.txi: Add vecnorm to manual. diff -r 3d26a174e143 -r 3d96400df713 NEWS --- a/NEWS Wed Nov 08 18:46:24 2017 -0500 +++ b/NEWS Mon Nov 06 10:27:09 2017 +0100 @@ -89,6 +89,7 @@ repelem rticks thetaticks + vecnorm xticklabels xticks yticklabels diff -r 3d26a174e143 -r 3d96400df713 doc/interpreter/linalg.txi --- a/doc/interpreter/linalg.txi Wed Nov 08 18:46:24 2017 -0500 +++ b/doc/interpreter/linalg.txi Mon Nov 06 10:27:09 2017 +0100 @@ -137,6 +137,8 @@ @DOCSTRING(rref) +@DOCSTRING(vecnorm) + @node Matrix Factorizations @section Matrix Factorizations @cindex matrix factorizations diff -r 3d26a174e143 -r 3d96400df713 libinterp/corefcn/data.cc --- a/libinterp/corefcn/data.cc Wed Nov 08 18:46:24 2017 -0500 +++ b/libinterp/corefcn/data.cc Mon Nov 06 10:27:09 2017 +0100 @@ -5581,7 +5581,7 @@ compute its norm. The result is returned as a column vector. Similarly, if @var{opt} is @qcode{"columns"} or @qcode{"cols"} then compute the norms of each column and return a row vector. -@seealso{normest, normest1, cond, svd} +@seealso{normest, normest1, vecnorm, cond, svd} @end deftypefn */) { int nargin = args.length (); diff -r 3d26a174e143 -r 3d96400df713 scripts/linear-algebra/module.mk --- a/scripts/linear-algebra/module.mk Wed Nov 08 18:46:24 2017 -0500 +++ b/scripts/linear-algebra/module.mk Mon Nov 06 10:27:09 2017 +0100 @@ -30,7 +30,8 @@ %reldir%/rref.m \ %reldir%/subspace.m \ %reldir%/trace.m \ - %reldir%/vech.m + %reldir%/vech.m \ + %reldir%/vecnorm.m %canon_reldir%dir = $(fcnfiledir)/linear-algebra diff -r 3d26a174e143 -r 3d96400df713 scripts/linear-algebra/vecnorm.m --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/scripts/linear-algebra/vecnorm.m Mon Nov 06 10:27:09 2017 +0100 @@ -0,0 +1,116 @@ +## Copyright (C) 2017 Marco Caliari +## +## 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 3 of the License, 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, see +## . + +## -*- texinfo -*- +## @deftypefn {} {@var{n} =} vecnorm (@var{A}) +## @deftypefnx {} {@var{n} =} vecnorm (@var{A}, @var{p}) +## @deftypefnx {} {@var{n} =} vecnorm (@var{A}, @var{p}, @var{dim}) +## Return the p-norm of the elements of @var{A} along dimension @var{dim}. +## +## The p-norm of a vector is defined as +## +## @tex +## $$ {\Vert A \Vert}_p = \left[ \sum_{i=1}^N {| A_i |}^p \right] ^ {1/p} $$ +## @end tex +## @ifnottex +## +## @example +## @var{p-norm} (@var{A}, @var{p}) = sum (abs (@var{A}) .^ @var{p})) ^ (1/@var{p}) +## @end example +## +## @end ifnottex +## If @var{p} is omitted it defaults to 2 (Euclidean norm). @var{p} can be +## @code{Inf} (absolute value of largest element). +## +## If @var{dim} is omitted the first non-singleton dimension is used. +## +## @seealso{norm} +## @end deftypefn + +function n = vecnorm (A, p = 2, dim) + + if (nargin < 1 || nargin > 3) + print_usage (); + endif + + if (! isscalar (p) || ! isreal (p) || (p <= 0)) + error ("vecnorm: P must be positive real scalar or Inf"); + endif + + sz = size (A); + if (nargin <= 2) + ## Find the first non-singleton dimension. + (dim = find (sz > 1, 1)) || (dim = 1); + elseif (! (isscalar (dim) && dim == fix (dim) && dim > 0)) + error ("vecnorm: DIM must be an integer and a valid dimension"); + endif + + ## Calculate norm using the value of p to accelerate special cases + switch (p) + case {1} + n = sum (abs (A), dim); + + case {2} + n = sqrt (sumsq (A, dim)); + + case {Inf} + n = max (abs (A), [], dim); + + otherwise + if (rem (p,2) == 0) + n = (sum ((real (A).^2 + imag (A).^2) .^ (p/2), dim)) .^ (1 / p); + else + n = (sum (abs (A) .^ p, dim)) .^ (1 / p); + endif + + endswitch + +endfunction + + +%!test +%! A = [0 1 2; 3 4 5]; +%! c = vecnorm (A); +%! r = vecnorm (A, 2, 2); +%! i = vecnorm (A, Inf); +%! assert (c, [3.0000, 4.1231, 5.3852], 1e-4); +%! assert (r, [2.2361; 7.0711], 1e-4); +%! assert (i, [3, 4, 5]); +%!test +%! A = [1, 2]; +%! assert (vecnorm (A), 2.2361, 1e-4); +%!test +%! A(:, :, 1) = [1, 2]; +%! A(:, :, 2) = [3, 4]; +%! A(:, :, 3) = [5, 6]; +%! ret(:, :, 1) = 2.2361; +%! ret(:, :, 2) = 5; +%! ret(:, :, 3) = 7.8102; +%! assert (vecnorm (A), ret, 1e-4); + +## Test input validation +%!error vecnorm () +%!error vecnorm (1,2,3,4) +%!error

vecnorm (1, [1 2]) +%!error

vecnorm (1, i) +%!error

vecnorm (1, -1) +%!error

vecnorm (1, 0) +%!error vecnorm (1, 2, [1 2]) +%!error vecnorm (1, 2, [1 2]) +%!error vecnorm (1, 2, 0) +%!error vecnorm (1, 2, -1)