1 ## Copyright (C) 2009, 2010, 2011 Lukas F. Reichlin
3 ## This file is part of LTI Syncope.
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.
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.
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/>.
19 ## @deftypefn {Function File} {@var{sys} =} minreal (@var{sys})
20 ## @deftypefnx {Function File} {@var{sys} =} minreal (@var{sys}, @var{tol})
21 ## Minimal realization or zero-pole cancellation of LTI models.
24 ## Author: Lukas Reichlin <lukas.reichlin@gmail.com>
25 ## Created: October 2009
28 function sys = minreal (sys, tol = "def")
30 if (nargin > 2) # nargin == 0 not possible because minreal is inside @lti
34 if (! is_real_scalar (tol) && tol != "def")
35 error ("minreal: second argument must be a real scalar");
38 sys = __minreal__ (sys, tol);
43 ## ss: minreal (SLICOT TB01PD)
46 %! A = ss (-2, 3, 4, 5);
48 %! C = minreal (B, 1e-15);
70 %! sys = ss (A, B, C, D, "scaled", true);
71 %! sysmin = minreal (sys, 0.0);
72 %! [Ar, Br, Cr, Dr] = ssdata (sysmin);
73 %! M = [Ar, Br; Cr, Dr];
75 %! Ae = [ 1.0000 -1.4142 1.4142
76 %! -2.8284 0.0000 1.0000
77 %! 2.8284 1.0000 0.0000 ];
83 %! Ce = [ 0.0000 0.0000 -1.4142
84 %! 0.0000 0.7071 0.7071 ];
88 %! Me = [Ae, Be; Ce, De];
90 %!assert (M, Me, 1e-4);
93 ## dss: minreal (SLICOT TG01JD)
94 ## FIXME: Test fails with larger ldwork in sltg01jd.cc
95 %!shared Ar, Br, Cr, Dr, Er, Ae, Be, Ce, De, Ee
96 %! A = [ -2 -3 0 0 0 0 0 0 0
98 %! 0 0 -2 -3 0 0 0 0 0
104 %! 0 0 0 0 0 0 0 0 1 ];
106 %! E = [ 1 0 0 0 0 0 0 0 0
114 %! 0 0 0 0 0 0 0 1 0 ];
126 %! C = [ 1 0 1 -3 0 1 0 2 0
127 %! 0 1 1 3 0 1 0 0 1 ];
131 %! sys = dss (A, B, C, D, E, "scaled", true);
132 %! sysmin = minreal (sys, 0.0);
133 %! [Ar, Br, Cr, Dr, Er] = dssdata (sysmin);
135 %! Ae = [ 1.0000 -0.0393 -0.0980 -0.1066 0.0781 -0.2330 0.0777
136 %! 0.0000 1.0312 0.2717 0.2609 -0.1533 0.6758 -0.3553
137 %! 0.0000 0.0000 1.3887 0.6699 -0.4281 1.6389 -0.7615
138 %! 0.0000 0.0000 0.0000 -1.2147 0.2423 -0.9792 0.4788
139 %! 0.0000 0.0000 0.0000 0.0000 -1.0545 0.5035 -0.2788
140 %! 0.0000 0.0000 0.0000 0.0000 0.0000 1.6355 -0.4323
141 %! 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 1.0000 ];
143 %! Ee = [ 0.4100 0.2590 0.5080 -0.3109 0.0705 0.1429 -0.1477
144 %! -0.7629 -0.3464 0.0992 -0.3007 0.0619 0.2483 -0.0152
145 %! 0.1120 -0.2124 -0.4184 -0.1288 0.0569 -0.4213 -0.6182
146 %! 0.0000 0.1122 -0.0039 0.2771 -0.0758 0.0975 0.3923
147 %! 0.0000 0.0000 0.3708 -0.4290 0.1006 0.1402 -0.2699
148 %! 0.0000 0.0000 0.0000 0.0000 0.9458 -0.2211 0.2378
149 %! 0.0000 0.0000 0.0000 0.5711 0.2648 0.5948 -0.5000 ];
151 %! Be = [ -0.5597 0.2363
159 %! Ce = [ 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 4.0000
160 %! 0.0000 0.0000 0.0000 0.0000 0.0000 3.1524 -1.7500 ];
162 %! De = zeros (2, 2);
164 %!assert (Ar, Ae, 1e-4);
165 %!assert (Br, Be, 1e-4);
166 %!assert (Cr, Ce, 1e-4);
167 %!assert (Dr, De, 1e-4);
168 %!assert (Er, Ee, 1e-4);
174 %! G1 = (s+1)*s*5/(s+1)/(s^2+s+1);
175 %! G2 = tf ([1, 1, 1], [2, 2, 2]);
176 %! G1min = minreal (G1);
177 %! G2min = minreal (G2);
178 %! a = G1min.num{1, 1};
179 %! b = G1min.den{1, 1};
180 %! c = G2min.num{1, 1};
181 %! d = G2min.den{1, 1};
182 %!assert (a, [5, 0], 1e-4);
183 %!assert (b, [1, 1, 1], 1e-4);
184 %!assert (c, 0.5, 1e-4);
185 %!assert (d, 1, 1e-4);