|
7017
|
1 ## Copyright (C) 1996, 2000, 2004, 2005, 2007 |
|
|
2 ## Auburn University. All rights reserved. |
|
3430
|
3 ## |
|
|
4 ## This file is part of Octave. |
|
|
5 ## |
|
|
6 ## Octave is free software; you can redistribute it and/or modify it |
|
7016
|
7 ## under the terms of the GNU General Public License as published by |
|
|
8 ## the Free Software Foundation; either version 3 of the License, or (at |
|
|
9 ## your option) any later version. |
|
3430
|
10 ## |
|
7016
|
11 ## Octave is distributed in the hope that it will be useful, but |
|
|
12 ## WITHOUT ANY WARRANTY; without even the implied warranty of |
|
|
13 ## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU |
|
|
14 ## General Public License for more details. |
|
3430
|
15 ## |
|
|
16 ## You should have received a copy of the GNU General Public License |
|
7016
|
17 ## along with Octave; see the file COPYING. If not, see |
|
|
18 ## <http://www.gnu.org/licenses/>. |
|
3430
|
19 |
|
|
20 ## -*- texinfo -*- |
|
5016
|
21 ## @deftypefn {Function File} {@var{K} =} hinf_ctr (@var{dgs}, @var{f}, @var{h}, @var{z}, @var{g}) |
|
|
22 ## Called by @code{hinfsyn} to compute the |
|
|
23 ## @iftex |
|
|
24 ## @tex |
|
|
25 ## $ { \cal H }_\infty $ |
|
|
26 ## @end tex |
|
|
27 ## @end iftex |
|
|
28 ## @ifinfo |
|
|
29 ## H-infinity |
|
|
30 ## @end ifinfo |
|
|
31 ## optimal controller. |
|
3430
|
32 ## |
|
|
33 ## @strong{Inputs} |
|
|
34 ## @table @var |
|
|
35 ## @item dgs |
|
|
36 ## data structure returned by @code{is_dgkf} |
|
3502
|
37 ## @item f |
|
|
38 ## @itemx h |
|
3430
|
39 ## feedback and filter gain (not partitioned) |
|
|
40 ## @item g |
|
|
41 ## final gamma value |
|
|
42 ## @end table |
|
|
43 ## @strong{Outputs} |
|
5016
|
44 ## @table @var |
|
|
45 ## @item K |
|
3502
|
46 ## controller (system data structure) |
|
5016
|
47 ## @end table |
|
3430
|
48 ## |
|
|
49 ## Do not attempt to use this at home; no argument checking performed. |
|
|
50 ## @end deftypefn |
|
|
51 |
|
|
52 ## Author: A. S. Hodel <a.s.hodel@eng.auburn.edu> |
|
|
53 ## Created: August 1995 |
|
|
54 ## Revised by Kai P. Mueller April 1998 to solve the general H_infinity |
|
|
55 ## problem using unitary transformations Q (on w and z) |
|
|
56 ## and non-singular transformations R (on u and y). |
|
|
57 |
|
|
58 function K = hinf_ctr (dgs, F, H, Z, g) |
|
|
59 |
|
7125
|
60 if (nargin != 5) |
|
|
61 print_usage (); |
|
|
62 endif |
|
|
63 |
|
3430
|
64 nw = dgs.nw; |
|
|
65 nu = dgs.nu; |
|
|
66 nz = dgs.nz; |
|
|
67 ny = dgs.ny; |
|
|
68 d22nz = dgs.Dyu_nz; |
|
|
69 |
|
|
70 B1 = dgs.Bw; |
|
|
71 B2 = dgs.Bu; |
|
|
72 C1 = dgs.Cz; |
|
|
73 C2 = dgs.Cy; |
|
|
74 C = [C1; C2]; |
|
|
75 D11 = dgs.Dzw; |
|
|
76 D12 = dgs.Dzu; |
|
|
77 D21 = dgs.Dyw; |
|
|
78 D22 = dgs.Dyu; |
|
|
79 A = dgs.A; |
|
|
80 Ru = dgs.Ru; |
|
|
81 Ry = dgs.Ry; |
|
|
82 |
|
|
83 nout = nz + ny; |
|
|
84 nin = nw + nu; |
|
7125
|
85 nstates = size (A, 1); |
|
3430
|
86 |
|
|
87 F11 = F(1:(nw-ny),:); |
|
|
88 F12 = F((nw-ny+1):nw,:); |
|
|
89 F2 = F((nw+1):nin,:); |
|
|
90 H11 = H(:,1:(nz-nu)); |
|
|
91 H12 = H(:,(nz-nu+1):nz); |
|
|
92 H2 = H(:,(nz+1):nout); |
|
|
93 |
|
|
94 ## D11 partitions |
|
|
95 D1111 = D11(1:(nz-nu),1:(nw-ny)); |
|
|
96 D1112 = D11(1:(nz-nu),(nw-ny+1):nw); |
|
|
97 D1121 = D11((nz-nu+1):nz,1:(nw-ny)); |
|
|
98 D1122 = D11((nz-nu+1):nz,(nw-ny+1):nw); |
|
|
99 |
|
|
100 ## D11ik may be the empty matrix, don't calculate with empty matrices |
|
7125
|
101 [nd1111, md1111] = size (D1111); |
|
|
102 md1112 = length (D1112); |
|
|
103 md1121 = length (D1121); |
|
3430
|
104 |
|
7125
|
105 if (nd1111 == 0) || md1112 == 0) |
|
3430
|
106 d11hat = -D1122; |
|
|
107 else |
|
7125
|
108 xx = inv (g*g*eye(nz-nu) - D1111*D1111'); |
|
3430
|
109 d11hat = -D1121*D1111'*xx*D1112 - D1122; |
|
|
110 endif |
|
|
111 if (md1112 == 0) |
|
7125
|
112 d21hat = eye (ny); |
|
3430
|
113 elseif (nd1111 == 0) |
|
7125
|
114 d21hat = chol (eye(ny) - D1112'*D1112/g/g); |
|
3430
|
115 else |
|
7125
|
116 xx = inv (g*g*eye(nz-nu) - D1111*D1111'); |
|
|
117 xx = eye (ny) - D1112'*xx*D1112; |
|
|
118 d21hat = chol (xx); |
|
3430
|
119 endif |
|
|
120 if (md1121 == 0) |
|
7125
|
121 d12hat = eye (nu); |
|
3430
|
122 elseif (md1111 == 0) |
|
7125
|
123 d12hat = chol (eye(nu) - D1121*D1121'/g/g)'; |
|
3430
|
124 else |
|
7125
|
125 xx = inv (g*g*eye(nw-ny) - D1111'*D1111); |
|
|
126 xx = eye (nu)-D1121*xx*D1121'; |
|
|
127 d12hat = chol (xx)'; |
|
3430
|
128 endif |
|
|
129 |
|
|
130 b2hat = (B2+H12)*d12hat; |
|
|
131 c2hat = -d21hat*(C2+F12)*Z; |
|
|
132 b1hat = -H2 + (b2hat/d12hat)*d11hat; |
|
|
133 c1hat = F2*Z + (d11hat/d21hat)*c2hat; |
|
|
134 ahat = A + H*C + (b2hat/d12hat)*c1hat; |
|
|
135 |
|
|
136 ## rescale controller by Ru and Ry |
|
|
137 b1hat = b1hat/Ry; |
|
|
138 c1hat = Ru\c1hat; |
|
|
139 bhat = [b1hat, b2hat]; |
|
|
140 chat = [c1hat; c2hat]; |
|
|
141 dhat = [Ru\d11hat/Ry, Ru\d12hat; d21hat/Ry, 0*d11hat']; |
|
|
142 |
|
|
143 ## non-zero D22 is a special case |
|
|
144 if (d22nz) |
|
7125
|
145 if (rank (eye(nu) + d11hat*D22) < nu) |
|
7131
|
146 error (" *** cannot compute controller for D22 non-zero."); |
|
3430
|
147 endif |
|
|
148 |
|
|
149 d22new = [D22, zeros(ny,ny); zeros(nu,nu), 0*D22']; |
|
7125
|
150 xx = inv (eye(nu+ny) + d22new*dhat); |
|
|
151 mhat = inv (eye(nu+ny) + dhat*d22new); |
|
3430
|
152 ahat = ahat - bhat*((eye(nu+ny)-xx)/dhat)*chat; |
|
|
153 bhat = bhat*xx; |
|
|
154 chat = mhat*chat; |
|
|
155 dhat = dhat*xx; |
|
|
156 |
|
|
157 endif |
|
|
158 |
|
7125
|
159 K = ss (ahat, bhat(:,1:ny), chat(1:nu,:), dhat(1:nu,1:ny)); |
|
3430
|
160 |
|
|
161 endfunction |
