changeset 7265:7da4a5262e2e

[project @ 2007-12-06 19:16:47 by jwe]
author jwe
date Thu, 06 Dec 2007 19:16:48 +0000
parents fac10884ddd4
children b42f8f3527a5
files liboctave/CMatrix.cc liboctave/ChangeLog liboctave/dMatrix.cc scripts/ChangeLog scripts/general/Makefile.in scripts/general/ishermitian.m scripts/general/issymmetric.m
diffstat 7 files changed, 122 insertions(+), 12 deletions(-) [+]
line wrap: on
line diff
--- a/liboctave/CMatrix.cc	Wed Dec 05 22:09:53 2007 +0000
+++ b/liboctave/CMatrix.cc	Thu Dec 06 19:16:48 2007 +0000
@@ -2865,11 +2865,16 @@
       for (octave_idx_type i = 0; i < nc; i++)
 	{
 	  octave_idx_type k = i * nc + i;
-	  pnpp [k] = pnpp [k] + padec [j];
-	  pdpp [k] = pdpp [k] + minus_one_j * padec [j];
+	  pnpp[k] += padec[j];
+	  pdpp[k] += minus_one_j * padec[j];
 	}      
+
       npp = m * npp;
+      pnpp = npp.fortran_vec ();
+
       dpp = m * dpp;
+      pdpp = dpp.fortran_vec ();
+
       minus_one_j *= -1;
     }
 
--- a/liboctave/ChangeLog	Wed Dec 05 22:09:53 2007 +0000
+++ b/liboctave/ChangeLog	Thu Dec 06 19:16:48 2007 +0000
@@ -1,3 +1,9 @@
+2007-12-06  John W. Eaton  <jwe@octave.org>
+
+	* CMatrix.cc (ComplexMatrix::expm): Update pointers to internal
+	data for npp and dpp after assignments.
+	* dMatrix.cc (Matrix::expm): Use same method as ComplexMatrix::expm.
+
 2007-12-04  John W. Eaton  <jwe@octave.org>
 
 	* Sparse.cc (assign (Sparse<LT>&,  const Sparse<RT>&)):
--- a/liboctave/dMatrix.cc	Wed Dec 05 22:09:53 2007 +0000
+++ b/liboctave/dMatrix.cc	Thu Dec 06 19:16:48 2007 +0000
@@ -2473,15 +2473,28 @@
   // npp, dpp: pade' approx polynomial matrices.
   
   Matrix npp (nc, nc, 0.0);
+  double *pnpp = npp.fortran_vec ();
   Matrix dpp = npp;
+  double *pdpp = dpp.fortran_vec ();
   
   // Now powers a^8 ... a^1.
   
   octave_idx_type minus_one_j = -1;
   for (octave_idx_type j = 7; j >= 0; j--)
     {
-      npp = m * npp + padec[j] * m;
-      dpp = m * dpp + (minus_one_j * padec[j]) * m;
+      for (octave_idx_type i = 0; i < nc; i++)
+	{
+	  octave_idx_type k = i * nc + i;
+	  pnpp[k] += padec[j];
+	  pdpp[k] += minus_one_j * padec[j];
+	}      
+
+      npp = m * npp;
+      pnpp = npp.fortran_vec ();
+
+      dpp = m * dpp;
+      pdpp = dpp.fortran_vec ();
+
       minus_one_j *= -1;
     }
   
--- a/scripts/ChangeLog	Wed Dec 05 22:09:53 2007 +0000
+++ b/scripts/ChangeLog	Thu Dec 06 19:16:48 2007 +0000
@@ -1,3 +1,15 @@
+2007-12-06  John W. Eaton  <jwe@octave.org>
+
+	* general/issymmetric.m: Move tests here from test/test_number.m
+
+2007-12-06  Jason Riedy  <ejr@cs.berkeley.edu>
+
+	* general/issymmetric.m: To keep its argument sparse and the
+	function quick, use the infinity norm rather than the 2-norm.
+	Also measure the symmetric part rather than the Hermitian part.
+	* general/ishermitian.m: New file.  Measure the Hermitian part.
+	* general/Makefile.in: Add ishermitian.m to SOURCES.
+
 2007-12-04  John W. Eaton  <jwe@octave.org>
 
 	* plot/__go_draw_axes__.m: Omit "font \"NAME,SIZE\"" in gnuplot
--- a/scripts/general/Makefile.in	Wed Dec 05 22:09:53 2007 +0000
+++ b/scripts/general/Makefile.in	Thu Dec 06 19:16:48 2007 +0000
@@ -38,10 +38,10 @@
   celldisp.m circshift.m common_size.m cplxpair.m cumtrapz.m deal.m del2.m \
   diff.m flipdim.m fliplr.m flipud.m gradient.m ind2sub.m int2str.m interp1.m \
   interp2.m interp3.m interpn.m interpft.m is_duplicate_entry.m isa.m \
-  isdefinite.m isdir.m isequal.m isequalwithequalnans.m isscalar.m \
-  issquare.m issymmetric.m isvector.m logical.m logspace.m lookup.m mod.m \
-  nargchk.m nextpow2.m nthroot.m num2str.m perror.m pol2cart.m \
-  polyarea.m postpad.m prepad.m quadl.m randperm.m rat.m rem.m \
+  isdefinite.m isdir.m isequal.m isequalwithequalnans.m ishermitian.m \
+  isscalar.m issquare.m issymmetric.m isvector.m logical.m logspace.m \
+  lookup.m mod.m nargchk.m nextpow2.m nthroot.m num2str.m perror.m \
+  pol2cart.m polyarea.m postpad.m prepad.m quadl.m randperm.m rat.m rem.m \
   repmat.m rot90.m rotdim.m shift.m shiftdim.m sortrows.m \
   sph2cart.m strerror.m structfun.m sub2ind.m trapz.m tril.m triu.m
 
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/scripts/general/ishermitian.m	Thu Dec 06 19:16:48 2007 +0000
@@ -0,0 +1,56 @@
+## Copyright (C) 1996, 1997, 2002, 2003, 2004, 2005, 2006, 2007
+##               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 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
+## <http://www.gnu.org/licenses/>.
+
+## -*- texinfo -*-
+## @deftypefn {Function File} {} ishermitian (@var{x}, @var{tol})
+## If @var{x} is Hermitian within the tolerance specified by @var{tol},
+## then return the dimension of @var{x}.  Otherwise, return 0.  If
+## @var{tol} is omitted, use a tolerance equal to the machine precision.
+## Matrix @var{x} is considered symmetric if
+## @code{norm (@var{x} - @var{x}', inf) / norm (@var{x}, inf) < @var{tol}}.
+## @seealso{size, rows, columns, length, ismatrix, isscalar,
+## issquare, issymmetric, isvector}
+## @end deftypefn
+
+## Author: A. S. Hodel <scotte@eng.auburn.edu>
+## Created: August 1993
+## Adapted-By: jwe
+
+function retval = ishermitian (x, tol)
+
+  if (nargin == 1 || nargin == 2)
+    retval = issquare (x);
+    if (retval != 0)
+      if (nargin == 1)
+        tol = eps;
+      endif
+      norm_x = norm (x, inf);
+      if (norm_x != 0 && norm (x - x', inf) / norm_x > tol)
+        retval = 0;
+      endif
+    endif
+  else
+    print_usage ();
+  endif
+
+endfunction
+
+%!assert(ishermitian ([1, 2i; -2i, 1]) == 2);
+%!assert(!ishermitian ([1, 2i; 2i, 1]));
+%!assert(ishermitian ([1, 2.1i; -2i, 1.1], 0.2) == 2);
--- a/scripts/general/issymmetric.m	Wed Dec 05 22:09:53 2007 +0000
+++ b/scripts/general/issymmetric.m	Thu Dec 06 19:16:48 2007 +0000
@@ -22,7 +22,9 @@
 ## If @var{x} is symmetric within the tolerance specified by @var{tol},
 ## then return the dimension of @var{x}.  Otherwise, return 0.  If
 ## @var{tol} is omitted, use a tolerance equal to the machine precision.
-## @seealso{size, rows, columns, length, ismatrix, isscalar,
+## Matrix @var{x} is considered symmetric if
+## @code{norm (@var{x} - @var{x}.', inf) / norm (@var{x}, inf) < @var{tol}}.
+## @seealso{size, rows, columns, length, ishermitian, ismatrix, isscalar,
 ## issquare, isvector}
 ## @end deftypefn
 
@@ -30,7 +32,7 @@
 ## Created: August 1993
 ## Adapted-By: jwe
 
-function retval = issymmetric (x,tol)
+function retval = issymmetric (x, tol)
 
   if (nargin == 1 || nargin == 2)
     retval = issquare (x);
@@ -38,8 +40,8 @@
       if (nargin == 1)
         tol = eps;
       endif
-      norm_x = norm (x);
-      if (norm_x != 0 && norm (x - x') / norm_x > tol)
+      norm_x = norm (x, inf);
+      if (norm_x != 0 && norm (x - x.', inf) / norm_x > tol)
         retval = 0;
       endif
     endif
@@ -48,3 +50,19 @@
   endif
 
 endfunction
+
+%!assert(issymmetric (1));
+%!assert(!(issymmetric ([1, 2])));
+%!assert(!(issymmetric ([])));
+%!assert(issymmetric ([1, 2; 2, 1]) == 2);
+%!assert(!(issymmetric ("test")));
+%!assert(issymmetric ([1, 2.1; 2, 1.1], 0.2) == 2);
+%!assert(!issymmetric ([1, 2i; -2i, 1]));
+%!assert(!(issymmetric ("t")));
+%!assert(!(issymmetric (["te"; "et"])));
+%!error issymmetric ([1, 2; 2, 1], 0, 0);
+%!error issymmetric ();
+
+%!test
+%! s.a = 1;
+%! assert(!(issymmetric (s)));