Mercurial > octave-nkf

comparison scripts/control/hinf/h2syn.m @ 7131:a184bc985c37

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
[project @ 2007-11-08 15:55:02 by jwe]
author jwe
date Thu, 08 Nov 2007 15:55:02 +0000
parents a1dbe9d80eee
children eb7bdde776f2
comparison
equal deleted inserted replaced
7130:5eeb46c784d7 7131:a184bc985c37
77 77
78 ## Updated for System structure December 1996 by John Ingram 78 ## Updated for System structure December 1996 by John Ingram
79 79
80 function [K, gain, Kc, Kf, Pc, Pf] = h2syn (Asys, nu, ny, tol) 80 function [K, gain, Kc, Kf, Pc, Pf] = h2syn (Asys, nu, ny, tol)
81 81
82 if ((nargin < 3) | (nargin > 4)) 82 if (nargin < 3 || nargin > 4)
83 print_usage (); 83 print_usage ();
84 elseif(nargin == 3 ) 84 elseif (nargin == 3)
85 [chkdgkf,dgs] = is_dgkf(Asys,nu,ny); 85 [chkdgkf, dgs] = is_dgkf (Asys, nu, ny);
86 elseif(nargin == 4) 86 elseif (nargin == 4)
87 [chkdgkf,dgs] = is_dgkf(Asys,nu,ny,tol); 87 [chkdgkf, dgs] = is_dgkf (Asys, nu, ny, tol);
88 endif 88 endif
89 89
90 if (!chkdgkf ) 90 if (! chkdgkf)
91 disp("h2syn: system does not meet required assumptions") 91 error ("h2syn: system does not meet required assumptions")
92 help is_dgkf
93 error("h2syn: exit");
94 endif 92 endif
95 93
96 ## extract dgs information 94 ## extract dgs information
97 nw = dgs.nw; nu = dgs.nu; 95 nw = dgs.nw;
98 A = dgs.A; Bw = dgs.Bw; Bu = dgs.Bu; 96 nu = dgs.nu;
99 Cz = dgs.Cz; Dzw = dgs.Dzw; Dzu = dgs.Dzu; nz = dgs.nz; 97 nz = dgs.nz;
100 Cy = dgs.Cy; Dyw = dgs.Dyw; Dyu = dgs.Dyu; ny = dgs.ny; 98 ny = dgs.ny;
99
100 A = dgs.A;
101
102 Bw = dgs.Bw;
103 Bu = dgs.Bu;
104
105 Cz = dgs.Cz;
106 Cy = dgs.Cy;
107
108 Dzw = dgs.Dzw;
109 Dzu = dgs.Dzu;
110
111 Dyw = dgs.Dyw;
112 Dyu = dgs.Dyu;
113
101 d22nz = dgs.Dyu_nz; 114 d22nz = dgs.Dyu_nz;
115
102 dflg = dgs.dflg; 116 dflg = dgs.dflg;
103 117
104 if(norm(Dzw,Inf) > norm([Dzw, Dzu ; Dyw, Dyu],Inf)*1e-12) 118 if (norm (Dzw, Inf) > norm ([Dzw, Dzu; Dyw, Dyu], Inf)*1e-12)
105 warning("h2syn: Dzw nonzero; feedforward not implemented") 119 warning ("h2syn: Dzw nonzero; feedforward not implemented")
106 Dzw 120 Dzw
107 D = [Dzw, Dzu ; Dyw, Dyu] 121 D = [Dzw, Dzu ; Dyw, Dyu]
108 endif 122 endif
109 123
110 ## recover i/o transformations 124 ## recover i/o transformations
111 Ru = dgs.Ru; Ry = dgs.Ry; 125 Ru = dgs.Ru;
112 [ncstates, ndstates, nout, nin] = sysdimensions(Asys); 126 Ry = dgs.Ry;
113 Atsam = sysgettsam(Asys); 127
114 [Ast, Ain, Aout] = sysgetsignals(Asys); 128 [ncstates, ndstates, nout, nin] = sysdimensions (Asys);
115 129 Atsam = sysgettsam (Asys);
116 if(dgs.dflg == 0) 130 [Ast, Ain, Aout] = sysgetsignals (Asys);
117 Pc = are(A,Bu*Bu',Cz'*Cz); # solve control, filtering ARE's 131
118 Pf = are(A',Cy'*Cy,Bw*Bw'); 132 if (dgs.dflg == 0)
133 Pc = are (A, Bu*Bu', Cz'*Cz); # solve control, filtering ARE's
134 Pf = are(A', Cy'*Cy, Bw*Bw');
119 F2 = -Bu'*Pc; # calculate feedback gains 135 F2 = -Bu'*Pc; # calculate feedback gains
120 L2 = -Pf*Cy'; 136 L2 = -Pf*Cy';
121 137
122 AF2 = A + Bu*F2; 138 AF2 = A + Bu*F2;
123 AL2 = A + L2*Cy; 139 AL2 = A + L2*Cy;
124 CzF2 = Cz + (Dzu/Ru)*F2; 140 CzF2 = Cz + (Dzu/Ru)*F2;
125 BwL2 = Bw+L2*(Ry\Dyw); 141 BwL2 = Bw+L2*(Ry\Dyw);
126 142
127 else 143 else
128 ## discrete time solution 144 ## discrete time solution
129 error("h2syn: discrete-time case not yet implemented") 145 error ("h2syn: discrete-time case not yet implemented")
130 Pc = dare(A,Bu*Bu',Cz'*Cz); 146 Pc = dare (A, Bu*Bu', Cz'*Cz);
131 Pf = dare(A',Cy'*Cy,Bw*Bw'); 147 Pf = dare (A', Cy'*Cy, Bw*Bw');
132 endif 148 endif
133 149
134 nn = ncstates + ndstates; 150 nn = ncstates + ndstates;
135 In = eye(nn); 151 In = eye (nn);
136 KA = A + Bu*F2 + L2*Cy; 152 KA = A + Bu*F2 + L2*Cy;
137 Kc1 = ss(AF2,Bw,CzF2,zeros(nz,nw)); 153 Kc1 = ss (AF2, Bw, CzF2, zeros (nz, nw));
138 Kf1 = ss(AL2,BwL2,F2,zeros(nu,nw)); 154 Kf1 = ss (AL2, BwL2, F2, zeros (nu, nw));
139 155
140 g1 = h2norm(Kc1); 156 g1 = h2norm (Kc1);
141 g2 = h2norm(Kf1); 157 g2 = h2norm (Kf1);
142 158
143 ## compute optimal closed loop gain 159 ## compute optimal closed loop gain
144 gain = sqrt ( g1*g1 + g2*g2 ); 160 gain = sqrt (g1*g1 + g2*g2);
145 161
146 if(nargout) 162 if (nargout)
147 Kst = strappend(Ast,"_K"); 163 Kst = strappend (Ast, "_K");
148 Kin = strappend(Aout((nout-ny+1):(nout)),"_K"); 164 Kin = strappend (Aout((nout-ny+1):(nout)), "_K");
149 Kout = strappend(Ain((nin-nu+1):(nin)),"_K"); 165 Kout = strappend (Ain((nin-nu+1):(nin)), "_K");
150 166
151 ## compute systems for return 167 ## compute systems for return
152 K = ss(KA,-L2/Ru,Ry\F2,zeros(nu,ny),Atsam,ncstates,ndstates,Kst,Kin,Kout); 168 K = ss (KA, -L2/Ru, Ry\F2, zeros(nu,ny), Atsam, ncstates,
169 ndstates, Kst, Kin, Kout);
153 endif 170 endif
154 171
155 if (nargout > 2) 172 if (nargout > 2)
156 ## system full information control state names 173 ## system full information control state names
157 stname2 = strappend(Ast,"_FI"); 174 stname2 = strappend (Ast, "_FI");
158 175
159 ## system full information control input names 176 ## system full information control input names
160 inname2 = strappend(Ast,"_FI_in"); 177 inname2 = strappend (Ast, "_FI_in");
161 178
162 ## system full information control output names 179 ## system full information control output names
163 outname2 = strappend(Aout(1:(nout-ny)),"_FI_out"); 180 outname2 = strappend (Aout(1:(nout-ny)), "_FI_out");
164 181
165 nz = rows (Cz); 182 nz = rows (Cz);
166 nw = columns (Bw); 183 nw = columns (Bw);
167 184
168 Kc = ss(AF2, In, CzF2, zeros(nz,nn), Atsam, ... 185 Kc = ss (AF2, In, CzF2, zeros(nz,nn), Atsam,
169 ncstates, ndstates, stname2, inname2, outname2); 186 ncstates, ndstates, stname2, inname2, outname2);
170 endif 187 endif
171 188
172 if (nargout >3) 189 if (nargout >3)
173 ## fix system state estimator state names 190 ## fix system state estimator state names
174 stname3 = strappend(Ast,"_Kf"); 191 stname3 = strappend (Ast, "_Kf");
175 192
176 ## fix system state estimator input names 193 ## fix system state estimator input names
177 inname3 = strappend(Ast,"_Kf_noise"); 194 inname3 = strappend (Ast, "_Kf_noise");
178 195
179 ## fix system state estimator output names 196 ## fix system state estimator output names
180 outname3 = strappend(Ast,"_est"); 197 outname3 = strappend (Ast, "_est");
181 198
182 Kf = ss(AL2, BwL2, In, zeros(nn,nw),Atsam, ... 199 Kf = ss (AL2, BwL2, In, zeros(nn,nw),Atsam,
183 ncstates, ndstates, stname3, inname3,outname3); 200 ncstates, ndstates, stname3, inname3, outname3);
184 endif 201 endif
185 202
186 endfunction 203 endfunction