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

   1 ## Copyright (C) 2009, 2010, 2011   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 ## @deftypefn {Function File} {@var{sys} =} sminreal (@var{sys})
  20 ## @deftypefnx {Function File} {@var{sys} =} sminreal (@var{sys}, @var{tol})
  21 ## Perform state-space model reduction based on structure.
  22 ## Remove states which have no influence on the input-output behaviour.
  23 ## The physical meaning of the states is retained.
  24 ##
  25 ## @strong{Inputs}
  26 ## @table @var
  27 ## @item sys
  28 ## State-space model.
  29 ## @item tol
  30 ## Optional tolerance for controllability and observability.
  31 ## Entries of the state-space matrices whose moduli are less or equal to @var{tol}
  32 ## are assumed to be zero.  Default value is 0.
  33 ## @end table
  34 ##
  35 ## @strong{Outputs}
  36 ## @table @var
  37 ## @item sys
  38 ## Reduced state-space model.
  39 ## @end table
  40 ##
  41 ## @seealso{minreal}
  42 ## @end deftypefn
  43
  44 ## Author: Lukas Reichlin <lukas.reichlin@gmail.com>
  45 ## Created: October 2009
  46 ## Version: 0.4
  47
  48 function sys = sminreal (sys, tol = 0)
  49
  50   if (nargin > 2)            # sminreal () not possible (inside @lti)
  51     print_usage ();
  52   endif
  53
  54   if (! isa (sys, "ss"))
  55     warning ("sminreal: system not in state-space form");
  56     sys = ss (sys);          # needed by __sys_prune__
  57   endif
  58
  59   if (! (is_real_scalar (tol) && tol >= 0))
  60     error ("sminreal: second argument is not a valid tolerance");
  61   endif
  62
  63   [a, b, c, d, e] = dssdata (sys, []);
  64
  65   a = abs (a) > tol;
  66   b = abs (b) > tol;
  67   c = abs (c) > tol;
  68
  69   if (! isempty (e))
  70     e = abs (e) > tol;
  71     a = a | e;
  72   endif
  73
  74   co_idx = __controllable_states__ (a, b);
  75   ob_idx = __controllable_states__ (a.', c.');
  76
  77   st_idx = intersect (co_idx, ob_idx);
  78
  79   sys = __sys_prune__ (sys, ":", ":", st_idx);
  80
  81 endfunction
  82
  83
  84 function c_idx = __controllable_states__ (a, b)
  85
  86   n = rows (a);              # number of states
  87   a = a & ! eye (n);         # set diagonal entries to zero
  88
  89   c_vec = any (b, 2);        # states directly controllable
  90   c_idx = find (c_vec);      # indices of directly controllable states
  91   c_idx_new = 0;             # any vector of length > 0 possible
  92
  93   while (all (length (c_idx) != [0, n]) && length(c_idx_new) != 0)
  94
  95     u_idx = find (! c_vec);  # indices of uncontrollable states
  96
  97     #{
  98     ## debug code
  99     a(u_idx, :)
 100     repmat (c_vec.', length (u_idx), 1)
 101     a(u_idx, :) & repmat (c_vec.', length (u_idx), 1)
 102     any (a(u_idx, :) & repmat (c_vec.', length (u_idx), 1), 2)
 103     find (any (a(u_idx, :) & repmat (c_vec.', length (u_idx), 1), 2))
 104     #}
 105
 106     c_idx_new = u_idx (find (any (a(u_idx, :) & repmat (c_vec.', length (u_idx), 1), 2)));
 107     c_idx = union (c_idx, c_idx_new);
 108     c_vec(c_idx_new) = 1;
 109
 110   endwhile
 111
 112 endfunction
 113
 114
 115 ## ss: sminreal
 116 %!shared B, C
 117 %!
 118 %! A = ss (-2, 3, 4, 5);
 119 %! B = A / A;
 120 %! C = sminreal (B);  # no states should be removed
 121 %!
 122 %!assert (C.a, B.a);
 123 %!assert (C.b, B.b);
 124 %!assert (C.c, B.c);
 125 %!assert (C.d, B.d);
 126
 127 %!shared A, B, D, E
 128 %!
 129 %! A = ss (-1, 1, 1, 0);
 130 %! B = ss (-2, 3, 4, 5);
 131 %! C = [A, B];
 132 %! D = sminreal (C(:, 1));
 133 %! E = sminreal (C(:, 2));
 134 %!
 135 %!assert (D.a, A.a);
 136 %!assert (D.b, A.b);
 137 %!assert (D.c, A.c);
 138 %!assert (D.d, A.d);
 139 %!assert (E.a, B.a);
 140 %!assert (E.b, B.b);
 141 %!assert (E.c, B.c);
 142 %!assert (E.d, B.d);