|
3440
|
1 ## Copyright (C) 1996, 1998 Auburn University. All rights reserved. |
|
|
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 the |
|
|
7 ## Free Software Foundation; either version 2, or (at your option) any |
|
|
8 ## later version. |
|
|
9 ## |
|
|
10 ## Octave is distributed in the hope that it will be useful, but WITHOUT |
|
|
11 ## ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or |
|
|
12 ## FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License |
|
|
13 ## 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, write to the Free |
|
|
17 ## Software Foundation, 59 Temple Place, Suite 330, Boston, MA 02111 USA. |
|
|
18 |
|
|
19 ## -*- texinfo -*- |
|
3502
|
20 ## @deftypefn {Function File} {[@var{retval}, @var{dgkf_struct} ] =} is_dgkf (@var{asys}, @var{nu}, @var{ny}, @var{tol} ) |
|
3440
|
21 ## Determine whether a continuous time state space system meets |
|
5016
|
22 ## assumptions of @acronym{DGKF} algorithm. |
|
3440
|
23 ## Partitions system into: |
|
|
24 ## @example |
|
5016
|
25 ## [dx/dt] [A | Bw Bu ][w] |
|
|
26 ## [ z ] = [Cz | Dzw Dzu ][u] |
|
3440
|
27 ## [ y ] [Cy | Dyw Dyu ] |
|
|
28 ## @end example |
|
|
29 ## or similar discrete-time system. |
|
3502
|
30 ## If necessary, orthogonal transformations @var{qw}, @var{qz} and nonsingular |
|
|
31 ## transformations @var{ru}, @var{ry} are applied to respective vectors |
|
5016
|
32 ## @var{w}, @var{z}, @var{u}, @var{y} in order to satisfy @acronym{DGKF} assumptions. |
|
3502
|
33 ## Loop shifting is used if @var{dyu} block is nonzero. |
|
3440
|
34 ## |
|
|
35 ## @strong{Inputs} |
|
|
36 ## @table @var |
|
3502
|
37 ## @item asys |
|
3440
|
38 ## system data structure |
|
|
39 ## @item nu |
|
|
40 ## number of controlled inputs |
|
|
41 ## @item ny |
|
|
42 ## number of measured outputs |
|
|
43 ## @item tol |
|
5016
|
44 ## threshold for 0; default: 200*@code{eps}. |
|
3440
|
45 ## @end table |
|
|
46 ## @strong{Outputs} |
|
|
47 ## @table @var |
|
|
48 ## @item retval |
|
|
49 ## true(1) if system passes check, false(0) otherwise |
|
|
50 ## @item dgkf_struct |
|
5016
|
51 ## data structure of @command{is_dgkf} results. Entries: |
|
3440
|
52 ## @table @var |
|
|
53 ## @item nw |
|
|
54 ## @itemx nz |
|
|
55 ## dimensions of @var{w}, @var{z} |
|
3502
|
56 ## @item a |
|
|
57 ## system @math{A} matrix |
|
|
58 ## @item bw |
|
|
59 ## (@var{n} x @var{nw}) @var{qw}-transformed disturbance input matrix |
|
|
60 ## @item bu |
|
|
61 ## (@var{n} x @var{nu}) @var{ru}-transformed controlled input matrix; |
|
3440
|
62 ## |
|
4946
|
63 ## @math{B = [Bw Bu]} |
|
3502
|
64 ## @item cz |
|
3440
|
65 ## (@var{nz} x @var{n}) Qz-transformed error output matrix |
|
3502
|
66 ## @item cy |
|
|
67 ## (@var{ny} x @var{n}) @var{ry}-transformed measured output matrix |
|
3440
|
68 ## |
|
4946
|
69 ## @math{C = [Cz; Cy]} |
|
3502
|
70 ## @item dzu |
|
|
71 ## @item dyw |
|
|
72 ## off-diagonal blocks of transformed system @math{D} matrix that enter |
|
3440
|
73 ## @var{z}, @var{y} from @var{u}, @var{w} respectively |
|
3502
|
74 ## @item ru |
|
3440
|
75 ## controlled input transformation matrix |
|
3502
|
76 ## @item ry |
|
3440
|
77 ## observed output transformation matrix |
|
3502
|
78 ## @item dyu_nz |
|
|
79 ## nonzero if the @var{dyu} block is nonzero. |
|
|
80 ## @item dyu |
|
|
81 ## untransformed @var{dyu} block |
|
3440
|
82 ## @item dflg |
|
|
83 ## nonzero if the system is discrete-time |
|
|
84 ## @end table |
|
|
85 ## @end table |
|
|
86 ## @code{is_dgkf} exits with an error if the system is mixed |
|
5016
|
87 ## discrete/continuous. |
|
3440
|
88 ## |
|
|
89 ## @strong{References} |
|
|
90 ## @table @strong |
|
|
91 ## @item [1] |
|
5016
|
92 ## Doyle, Glover, Khargonekar, Francis, @cite{State Space Solutions to Standard} |
|
|
93 ## @iftex |
|
|
94 ## @tex |
|
|
95 ## $ { \cal H }_2 $ @cite{and} $ { \cal H }_\infty $ |
|
|
96 ## @end tex |
|
|
97 ## @end iftex |
|
|
98 ## @ifinfo |
|
|
99 ## @cite{H-2 and H-infinity} |
|
|
100 ## @end ifinfo |
|
|
101 ## @cite{Control Problems}, @acronym{IEEE} @acronym{TAC} August 1989. |
|
3440
|
102 ## @item [2] |
|
5016
|
103 ## Maciejowksi, J.M., @cite{Multivariable Feedback Design}, Addison-Wesley, 1989. |
|
3440
|
104 ## @end table |
|
|
105 ## @end deftypefn |
|
|
106 |
|
|
107 ## Author: A. S. Hodel <a.s.hodel@eng.auburn.edu> |
|
|
108 ## Updated by John Ingram July 1996 to accept structured systems |
|
|
109 |
|
|
110 ## Revised by Kai P. Mueller April 1998 to solve the general H_infinity |
|
|
111 ## problem using unitary transformations Q (on w and z) |
|
|
112 ## and non-singular transformations R (on u and y) such |
|
|
113 ## that the Dzu and Dyw matrices of the transformed plant |
|
|
114 ## |
|
|
115 ## ~ |
|
|
116 ## P (the variable Asys here) |
|
|
117 ## |
|
|
118 ## become |
|
|
119 ## |
|
|
120 ## ~ -1 T |
|
|
121 ## D = Q D R = [ 0 I ] or [ I ], |
|
|
122 ## 12 12 12 12 |
|
|
123 ## |
|
|
124 ## ~ T |
|
|
125 ## D = R D Q = [ 0 I ] or [ I ]. |
|
|
126 ## 21 21 21 21 |
|
|
127 ## |
|
|
128 ## This transformation together with the algorithm in [1] solves |
|
|
129 ## the general problem (see [2] for example). |
|
|
130 |
|
|
131 function [retval, dgkf_struct] = is_dgkf (Asys, nu, ny, tol) |
|
|
132 |
|
|
133 if (nargin < 3) | (nargin > 4) |
|
|
134 usage("[retval,dgkf_struct] = is_dgkf(Asys,nu,ny{,tol})"); |
|
4030
|
135 elseif (! isscalar(nu) | ! isscalar(ny) ) |
|
3440
|
136 error("is_dgkf: arguments 2 and 3 must be scalars") |
|
4030
|
137 elseif (! isstruct(Asys) ) |
|
3440
|
138 error("Argument 1 must be a system data structure"); |
|
|
139 endif |
|
|
140 if(nargin < 4) |
|
|
141 tol = 200*eps; |
|
|
142 elseif( !is_sample(tol) ) |
|
|
143 error("is_dgkf: tol must be a positive scalar") |
|
|
144 endif |
|
|
145 |
|
|
146 retval = 1; # assume passes test |
|
|
147 |
|
|
148 dflg = is_digital(Asys); |
|
|
149 [Anc, Anz, nin, nout ] = sysdimensions(Asys); |
|
|
150 |
|
|
151 if( Anz == 0 & Anc == 0 ) |
|
|
152 error("is_dgkf: no system states"); |
|
|
153 elseif( nu >= nin ) |
|
|
154 error("is_dgkf: insufficient number of disturbance inputs"); |
|
|
155 elseif( ny >= nout ) |
|
|
156 error("is_dgkf: insufficient number of regulated outputs"); |
|
|
157 endif |
|
|
158 |
|
|
159 nw = nin - nu; nw1 = nw + 1; |
|
|
160 nz = nout - ny; nz1 = nz + 1; |
|
|
161 |
|
|
162 [A,B,C,D] = sys2ss(Asys); |
|
|
163 ## scale input/output for numerical reasons |
|
|
164 if(norm (C, "fro") * norm (B, "fro") == 0) |
|
|
165 error("||C||*||B|| = 0; no dynamic connnection from inputs to outputs"); |
|
|
166 endif |
|
|
167 xx = sqrt(norm(B, Inf) / norm(C, Inf)); |
|
|
168 B = B / xx; C = C * xx; |
|
|
169 |
|
|
170 ## partition matrices |
|
|
171 Bw = B(:,1:nw); Bu = B(:,nw1:nin); |
|
|
172 Cz = C(1:nz,:); Dzw = D(1:nz,1:nw); Dzu = D(1:nz,nw1:nin); |
|
|
173 Cy = C(nz1:nout,:); Dyw = D(nz1:nout,1:nw); Dyu = D(nz1:nout,nw1:nin); |
|
|
174 |
|
|
175 ## Check for loopo shifting |
|
|
176 Dyu_nz = (norm(Dyu,Inf) != 0); |
|
|
177 if (Dyu_nz) |
|
|
178 warning("is_dgkf: D22 nonzero; performing loop shifting"); |
|
|
179 endif |
|
|
180 |
|
|
181 ## 12 - rank condition at w = 0 |
|
|
182 xx =[A, Bu; Cz, Dzu]; |
|
|
183 [nr, nc] = size(xx); |
|
|
184 irank = rank(xx); |
|
|
185 if (irank != nc) |
|
|
186 retval = 0; |
|
|
187 warning(sprintf("rank([A Bu; Cz Dzu]) = %d, need %d; n=%d, nz=%d, nu=%d", ... |
|
|
188 irank,nc,(Anc+Anz),nz,nu)); |
|
|
189 warning(" *** 12-rank condition violated at w = 0."); |
|
|
190 endif |
|
|
191 |
|
|
192 ## 21 - rank condition at w = 0 |
|
|
193 xx =[A, Bw; Cy, Dyw]; |
|
|
194 [nr, nc] = size(xx); |
|
|
195 irank = rank(xx); |
|
|
196 if (irank != nr) |
|
|
197 retval = 0; |
|
|
198 warning(sprintf("rank([A Bw; Cy Dyw]) = %d, need %d; n=%d, ny=%d, nw=%d", ... |
|
|
199 irank,nr,(Anc+Anz),ny,nw)); |
|
|
200 warning(" *** 21-rank condition violated at w = 0."); |
|
|
201 endif |
|
|
202 |
|
|
203 ## can Dzu be transformed to become [0 I]' or [I]? |
|
|
204 ## This ensures a normalized weight |
|
|
205 [Qz, Ru] = qr(Dzu); |
|
|
206 irank = rank(Ru); |
|
|
207 if (irank != nu) |
|
|
208 retval = 0; |
|
|
209 warning(sprintf("*** rank(Dzu(%d x %d) = %d", nz, nu, irank)); |
|
|
210 warning(" *** Dzu does not have full column rank."); |
|
|
211 endif |
|
|
212 if (nu >= nz) |
|
|
213 Qz = Qz(:,1:nu)'; |
|
|
214 else |
|
|
215 Qz = [Qz(:,(nu+1):nz), Qz(:,1:nu)]'; |
|
|
216 endif |
|
|
217 Ru = Ru(1:nu,:); |
|
|
218 |
|
|
219 ## can Dyw be transformed to become [0 I] or [I]? |
|
|
220 ## This ensures a normalized weight |
|
|
221 [Qw, Ry] = qr(Dyw'); |
|
|
222 irank = rank(Ry); |
|
|
223 if (irank != ny) |
|
|
224 retval = 0; |
|
|
225 warning(sprintf("*** rank(Dyw(%d x %d) = %d", ny, nw, irank)); |
|
|
226 warning(" *** Dyw does not have full row rank."); |
|
|
227 endif |
|
|
228 |
|
|
229 if (ny >= nw) |
|
|
230 Qw = Qw(:,1:ny); |
|
|
231 else |
|
|
232 Qw = [Qw(:,(ny+1):nw), Qw(:,1:ny)]; |
|
|
233 endif |
|
|
234 Ry = Ry(1:ny,:)'; |
|
|
235 |
|
|
236 ## transform P by Qz/Ru and Qw/Ry |
|
|
237 Bw = Bw*Qw; |
|
|
238 Bu = Bu/Ru; |
|
|
239 B = [Bw, Bu]; |
|
|
240 Cz = Qz*Cz; |
|
|
241 Cy = Ry\Cy; |
|
|
242 C = [Cz; Cy]; |
|
|
243 Dzw = Qz*Dzw*Qw; |
|
|
244 Dzu = Qz*Dzu/Ru; |
|
|
245 Dyw = Ry\Dyw*Qw; |
|
|
246 |
|
|
247 ## pack the return structure |
|
|
248 dgkf_struct.nw = nw; |
|
|
249 dgkf_struct.nu = nu; |
|
|
250 dgkf_struct.nz = nz; |
|
|
251 dgkf_struct.ny = ny; |
|
|
252 dgkf_struct.A = A; |
|
|
253 dgkf_struct.Bw = Bw; |
|
|
254 dgkf_struct.Bu = Bu; |
|
|
255 dgkf_struct.Cz = Cz; |
|
|
256 dgkf_struct.Cy = Cy; |
|
|
257 dgkf_struct.Dzw = Dzw; |
|
|
258 dgkf_struct.Dzu = Dzu; |
|
|
259 dgkf_struct.Dyw = Dyw; |
|
|
260 dgkf_struct.Dyu = Dyu; |
|
|
261 dgkf_struct.Ru = Ru; |
|
|
262 dgkf_struct.Ry = Ry; |
|
|
263 dgkf_struct.Dyu_nz = Dyu_nz; |
|
|
264 dgkf_struct.dflg = dflg; |
|
|
265 |
|
|
266 endfunction |
