|
7017
|
1 ## Copyright (C) 1996, 2000, 2002, 2004, 2005, 2007 |
|
|
2 ## Auburn University. All rights reserved. |
|
3431
|
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. |
|
3431
|
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. |
|
3431
|
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/>. |
|
3431
|
19 |
|
|
20 ## -*- texinfo -*- |
|
|
21 ## @deftypefn {Function File} {} analdemo () |
|
|
22 ## Octave Controls toolbox demo: State Space analysis demo |
|
|
23 ## @end deftypefn |
|
|
24 |
|
|
25 ## Author: David Clem |
|
|
26 ## Created: August 15, 1994 |
|
|
27 ## Updated by John Ingram December 1996 |
|
|
28 |
|
|
29 function analdemo () |
|
|
30 |
|
|
31 while (1) |
|
|
32 clc |
|
|
33 k=0; |
|
|
34 while(k > 8 || k < 1) |
|
|
35 k = menu("Octave State Space Analysis Demo", ... |
|
5016
|
36 "System gramians (gram, dgram)", ... |
|
3431
|
37 "System zeros (tzero)", ... |
|
|
38 "Continuous => Discrete and Discrete => Continuous conversions (c2d,d2c)", ... |
|
|
39 "Algebraic Riccati Equation (are, dare)", ... |
|
|
40 "Balanced realizations (balreal, dbalreal)", ... |
|
|
41 "Open loop truncations via Hankel singular values (balreal, dbalreal)", ... |
|
|
42 "SISO pole placement", ... |
|
|
43 "Return to main demo menu"); |
|
|
44 endwhile |
|
|
45 if (k == 1) |
|
|
46 clc |
|
|
47 help dgram |
|
|
48 prompt |
|
|
49 |
|
|
50 clc |
|
5016
|
51 disp("System Gramians: (see Moore, IEEE T-AC, 1981) \n"); |
|
3431
|
52 disp("Example #1, consider the discrete time state space system:\n"); |
|
|
53 a=[1, 5, -8.4; 1.2, -3, 5; 1, 7, 9] |
|
|
54 b=[1, 5; 2, 6; -4.4, 5] |
|
|
55 c=[1, -1.5, 2; 6, -9.8, 1] |
|
|
56 d=0 |
|
|
57 prompt |
|
5016
|
58 disp("\nThe discrete controllability gramian is computed as follows:"); |
|
|
59 cmd = "gramian = dgram(a, b);"; |
|
3431
|
60 run_cmd; |
|
|
61 disp("Results:\n"); |
|
5016
|
62 gramian = dgram(a,b) |
|
3431
|
63 disp("Variable Description:\n"); |
|
5016
|
64 disp("gramian => discrete controllability gramian"); |
|
3431
|
65 disp("a, b => a and b matrices of discrete time system\n"); |
|
5016
|
66 disp("A dual approach may be used to compute the observability gramian."); |
|
3431
|
67 prompt |
|
|
68 clc |
|
|
69 |
|
|
70 help gram |
|
|
71 prompt |
|
|
72 clc |
|
|
73 |
|
|
74 disp("Example #2, consider the continuous state space system:\n"); |
|
4139
|
75 a=[1, 3, -10.2; 3.7, -2, 9; 1, 3, 7] |
|
3431
|
76 b=[1, 12; 6, 2; -3.8, 7] |
|
|
77 c=[1, -1.1, 7; 3, -9.8, 2] |
|
|
78 d=0 |
|
|
79 prompt |
|
5016
|
80 disp("\nThe continuous controllability gramian is computed as follows:"); |
|
|
81 cmd = "gramian = gram(a, b);"; |
|
3431
|
82 run_cmd; |
|
|
83 disp("Results:\n"); |
|
5016
|
84 gramian = gram(a,b) |
|
3431
|
85 disp("Variable Description:\n"); |
|
5016
|
86 disp("gramian => continuous controllability gramian"); |
|
3431
|
87 disp("a, b => a and b matrices of continuous time system\n"); |
|
5016
|
88 disp("A dual approach may be used to compute the observability gramian."); |
|
3431
|
89 prompt |
|
|
90 clc |
|
|
91 |
|
|
92 |
|
|
93 elseif (k == 2) |
|
|
94 clc |
|
|
95 help tzero |
|
|
96 prompt |
|
|
97 |
|
|
98 disp("System zeros (tzero) example\n"); |
|
|
99 disp("Example #1, consider the state space system:\n"); |
|
|
100 a=[0, 1, 0; -10, -2, 0; -10, 0, -8] |
|
|
101 b=[0; 1; 9] |
|
|
102 c=[-10, 0, -4] |
|
|
103 d=1 |
|
|
104 prompt |
|
|
105 disp("\nTo compute the zeros of this system, enter the following command:\n"); |
|
|
106 cmd = "zer = tzero(a,b,c,d);"; |
|
|
107 run_cmd; |
|
|
108 disp("Results:\n"); |
|
|
109 zer = tzero(a,b,c,d) |
|
|
110 disp("Variable Description:\n"); |
|
|
111 disp("zer => zeros of state space system"); |
|
|
112 disp("a, b, c, d => state space system used as input argument"); |
|
|
113 prompt |
|
|
114 clc |
|
|
115 |
|
|
116 disp("Example #2, consider the state space system from example 1 again:"); |
|
4771
|
117 cmd = "sys = ss(a,b,c,d);"; |
|
3431
|
118 disp(cmd); |
|
|
119 eval(cmd); |
|
|
120 sysout(sys); |
|
|
121 disp("\nThe zeros of this system can also be calculated directly from the"); |
|
|
122 disp("system variable:"); |
|
|
123 cmd = "zer = tzero(sys);"; |
|
|
124 run_cmd; |
|
|
125 disp("Results:\n") |
|
|
126 zer = tzero(sys) |
|
|
127 disp("Variable Description:\n"); |
|
|
128 disp("zer => zeros of state space system"); |
|
|
129 disp("sys => state space system used as input argument"); |
|
|
130 prompt |
|
|
131 clc |
|
|
132 |
|
|
133 elseif (k == 3) |
|
|
134 clc |
|
|
135 help c2d |
|
|
136 prompt |
|
|
137 |
|
|
138 clc |
|
|
139 disp("Continuous => Discrete and Discrete => Continuous conversions (c2d,d2c)"); |
|
|
140 disp("\nExample #1, consider the following continuous state space system"); |
|
4771
|
141 cmd = "sys_cont = ss([-11, 6; -15, 8], [1; 2], [2, -1], 0);"; |
|
3431
|
142 eval(cmd); |
|
|
143 disp(cmd); |
|
|
144 disp("Examine the poles and zeros of the continuous system:"); |
|
|
145 sysout(sys_cont,"all"); |
|
|
146 disp("\nTo convert this to a discrete system, a sampling time is needed:"); |
|
|
147 cmd = "Tsam = 0.5;"; |
|
|
148 run_cmd; |
|
|
149 disp("\nNow convert to a discrete system with the command:"); |
|
|
150 cmd = "sys_disc = c2d(sys_cont,Tsam);"; |
|
|
151 run_cmd; |
|
|
152 disp("Examine the poles and zeros of the discrete system:"); |
|
|
153 sysout(sys_disc,"all"); |
|
|
154 prompt |
|
|
155 clc |
|
|
156 |
|
|
157 disp("\nNow we will convert the discrete system back to a continuous system"); |
|
|
158 disp("using the d2c command:"); |
|
|
159 help d2c |
|
|
160 prompt |
|
|
161 cmd = "new_sys_cont = d2c(sys_disc);"; |
|
|
162 run_cmd; |
|
|
163 disp("\nExamine the poles and zeros of the discrete system:"); |
|
|
164 sysout(new_sys_cont,"all"); |
|
|
165 prompt |
|
|
166 |
|
|
167 elseif (k == 4) |
|
|
168 clc |
|
|
169 help are |
|
|
170 prompt |
|
|
171 clc |
|
|
172 |
|
|
173 disp("Algebraic Riccati Equation (are, dare)"); |
|
|
174 |
|
|
175 disp("\nExample #1, consider the continuous state space system:\n"); |
|
|
176 a=[1, 3, -10.2; 3.7, -2, 9; 1, 3, 7] |
|
|
177 b=[1, 12; 6, 2; -3.8, 7] |
|
|
178 c=[1, -1.1, 7; 3, -9.8, 2] |
|
|
179 d=0 |
|
|
180 prompt |
|
|
181 disp("\nThe solution to the continuous algebraic riccati equation"); |
|
|
182 disp("is computed as follows:"); |
|
|
183 cmd = "x_cont = are(a, b, c);"; |
|
|
184 run_cmd; |
|
|
185 disp("Results:\n") |
|
|
186 x_cont = are(a,b,c) |
|
|
187 disp("Variable Description:\n") |
|
|
188 disp("x_cont => solution to the continuous algebraic riccati equation"); |
|
|
189 disp("a, b, c => a, b, and c matrices of continuous time system\n"); |
|
|
190 prompt |
|
|
191 |
|
|
192 clc |
|
|
193 help dare |
|
|
194 prompt |
|
|
195 clc |
|
|
196 |
|
|
197 disp("Example #2, consider the discrete time state space system:\n"); |
|
|
198 a=[1, 5, -8.4; 1.2, -3, 5; 1, 7, 9] |
|
|
199 b=[1, 5; 2, 6; -4.4, 5] |
|
|
200 c=[1, -1.5, 2; 6, -9.8, 1] |
|
|
201 d=0 |
|
|
202 r=eye(columns(b)) |
|
|
203 prompt |
|
|
204 disp("\nThe solution to the continuous algebraic riccati equation"); |
|
|
205 disp("is computed as follows:"); |
|
|
206 cmd = "x_disc = dare(a, b, c, r);"; |
|
|
207 run_cmd; |
|
|
208 disp("Results:\n") |
|
|
209 x_disc = dare(a,b,c,r) |
|
|
210 disp("Variable Description:\n"); |
|
|
211 disp("x_disc => solution to the discrete algebraic riccati equation"); |
|
|
212 disp("a, b, c => a, b and c matrices of discrete time system\n"); |
|
|
213 prompt |
|
|
214 clc |
|
|
215 |
|
|
216 elseif (k == 5) |
|
|
217 disp("--- Balanced realization: not yet implemented") |
|
|
218 elseif (k == 6) |
|
|
219 disp("--- Open loop balanced truncation: not yet implemented") |
|
|
220 elseif (k == 7) |
|
|
221 disp("SISO pole placement example:") |
|
4771
|
222 cmd = "sys=tf(1, [1, -2, 1]);"; |
|
3431
|
223 run_cmd |
|
|
224 disp("System in zero-pole form is:") |
|
|
225 cmd = "sysout(sys,\"zp\");"; |
|
|
226 run_cmd |
|
|
227 disp("and in state space form:") |
|
|
228 cmd = "sysout(sys,\"ss\");"; |
|
|
229 run_cmd |
|
|
230 disp("Desired poles at -1, -1"); |
|
|
231 cmd = "K=place(sys, [-1, -1])"; |
|
|
232 run_cmd |
|
|
233 disp("Check results:") |
|
|
234 cmd = "[A,B] = sys2ss(sys);"; |
|
|
235 run_cmd |
|
|
236 cmd = "poles=eig(A-B*K)"; |
|
|
237 run_cmd |
|
|
238 prompt |
|
|
239 elseif (k == 8) |
|
|
240 return |
|
|
241 endif |
|
|
242 endwhile |
|
|
243 endfunction |
|
|
244 |
