octave_packages/m/linear-algebra/isdefinite.m

   1 ## Copyright (C) 2003-2012 Gabriele Pannocchia
   2 ##
   3 ## This file is part of Octave.
   4 ##
   5 ## Octave is free software; you can redistribute it and/or modify it
   6 ## under the terms of the GNU General Public License as published by
   7 ## the Free Software Foundation; either version 3 of the License, or (at
   8 ## your option) any later version.
   9 ##
  10 ## Octave is distributed in the hope that it will be useful, but
  11 ## WITHOUT ANY WARRANTY; without even the implied warranty of
  12 ## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  13 ## General Public License for more details.
  14 ##
  15 ## You should have received a copy of the GNU General Public License
  16 ## along with Octave; see the file COPYING.  If not, see
  17 ## <http://www.gnu.org/licenses/>.
  18
  19 ## -*- texinfo -*-
  20 ## @deftypefn  {Function File} {} isdefinite (@var{x})
  21 ## @deftypefnx {Function File} {} isdefinite (@var{x}, @var{tol})
  22 ## Return 1 if @var{x} is symmetric positive definite within the
  23 ## tolerance specified by @var{tol} or 0 if @var{x} is symmetric
  24 ## positive semidefinite.  Otherwise, return -1.  If @var{tol}
  25 ## is omitted, use a tolerance of
  26 ## @code{100 * eps * norm (@var{x}, "fro")}
  27 ## @seealso{issymmetric, ishermitian}
  28 ## @end deftypefn
  29
  30 ## Author: Gabriele Pannocchia <g.pannocchia@ing.unipi.it>
  31 ## Created: November 2003
  32 ## Adapted-By: jwe
  33
  34 function retval = isdefinite (x, tol)
  35
  36   if (nargin < 1 || nargin > 2)
  37     print_usage ();
  38   endif
  39
  40   if (! isfloat (x))
  41     x = double (x);
  42   endif
  43
  44   if (nargin == 1)
  45     tol = 100 * eps (class (x)) * norm (x, "fro");
  46   endif
  47
  48   if (! ishermitian (x, tol))
  49     error ("isdefinite: X must be a Hermitian matrix");
  50   endif
  51
  52   e = tol * eye (rows (x));
  53   [r, p] = chol (x - e);
  54   if (p == 0)
  55     retval = 1;
  56   else
  57     [r, p] = chol (x + e);
  58     if (p == 0)
  59       retval = 0;
  60     else
  61       retval = -1;
  62     endif
  63   endif
  64
  65 endfunction
  66
  67
  68 %!test
  69 %! A = [-1 0; 0 -1];
  70 %! assert (isdefinite (A), -1)
  71
  72 %!test
  73 %! A = [1 0; 0 1];
  74 %! assert (isdefinite (A), 1)
  75
  76 %!test
  77 %! A = [2 -1 0; -1 2 -1; 0 -1 2];
  78 %! assert (isdefinite (A), 1)
  79
  80 %!test
  81 %! A = [1 0; 0 0];
  82 %! assert (isdefinite (A), 0)
  83
  84 %!error isdefinite ()
  85 %!error isdefinite (1,2,3)
  86 %!error <X must be a Hermitian matrix> isdefinite ([1 2; 3 4])
  87