octave_packages/control-2.3.52/@lti/mpower.m

   1 ## Copyright (C) 2009   Lukas F. Reichlin
   2 ##
   3 ## This file is part of LTI Syncope.
   4 ##
   5 ## LTI Syncope is free software: you can redistribute it and/or modify
   6 ## it under the terms of the GNU General Public License as published by
   7 ## the Free Software Foundation, either version 3 of the License, or
   8 ## (at your option) any later version.
   9 ##
  10 ## LTI Syncope is distributed in the hope that it will be useful,
  11 ## but WITHOUT ANY WARRANTY; without even the implied warranty of
  12 ## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  13 ## GNU General Public License for more details.
  14 ##
  15 ## You should have received a copy of the GNU General Public License
  16 ## along with LTI Syncope.  If not, see <http://www.gnu.org/licenses/>.
  17
  18 ## -*- texinfo -*-
  19 ## Matrix power of LTI objects.  The exponent must be an integer.
  20 ## Used by Octave for "sys^int".
  21
  22 ## Author: Lukas Reichlin <lukas.reichlin@gmail.com>
  23 ## Created: October 2009
  24 ## Version: 0.1
  25
  26 function retsys = mpower (sys, e)
  27
  28   if (nargin != 2)       # prevent sys = mpower (a, b, c, ...)
  29     error ("lti: mpower: this is a binary operator");
  30   endif
  31
  32   if (! is_real_scalar (e) || e != round (e))
  33     error ("lti: mpower: exponent must be an integer");
  34   endif
  35
  36   [p, m] = size (sys);
  37
  38   if (p != m)
  39     error ("lti: mpower: system must be square");
  40   endif
  41
  42   ex = round (abs (e));  # make sure ex is a positive integer
  43
  44   switch (sign (e))
  45     case -1              # lti^-ex
  46       sys = inv (sys);
  47       retsys = sys;
  48     case 0               # lti^0
  49       retsys = eye (p);
  50       return;
  51     case 1               # lti^ex
  52       retsys = sys;
  53   endswitch
  54
  55   for k = 2 : ex
  56     retsys = retsys * sys;
  57   endfor
  58
  59 endfunction