home *** CD-ROM | disk | FTP | other *** search
/ Liren Large Software Subsidy 10 / 10.iso / l / l455 / 10.ddi / CONTROL.DI$ / JETDEMO.M < prev    next >
Encoding:
Text File  |  1993-03-11  |  8.1 KB  |  213 lines
  1. echo on
  2.  
  3. % This file demonstrates MATLAB's ability for classical control system
  4. % design by going through the design of a YAW DAMPER for a Jet Transport
  5. % aircraft.
  6.  
  7. echo off
  8. %    Copyright (c) 1986-93 by the MathWorks, Inc.
  9. echo on
  10.  
  11. pause % Press any key to continue ...
  12.  
  13.  
  14. % Define Jet Transport model MACH=0.8 H=40,000ft
  15.  
  16. A=[-.0558 -.9968 .0802 .0415
  17.     .598 -.115 -.0318 0
  18.    -3.05 .388 -.4650 0
  19.    0 0.0805 1 0];
  20. B=[.0729  .0001
  21.    -4.75   1.23
  22.     1.53   10.63
  23.     0      0];
  24. C=[0 1 0 0
  25.    0 0 0 1];
  26. D=[0 0
  27.    0 0];
  28. states='beta yaw roll phi';
  29. inputs='rudder aileron';
  30. outputs='yaw-rate bank-angle';
  31. printsys(A,B,C,D,inputs,outputs,states)
  32.  
  33. % These are the state space matrices for a Jet Transport during cruise flight.
  34. % The model has two inputs and two outputs.  The units are radians for beta
  35. % (sideslip angle) and phi (bank angle) and radians/sec for yaw (yaw rate) and
  36. % roll (roll rate).  The rudder and aileron deflections are in degrees.
  37.  
  38. pause % Press any key to continue ...
  39.  
  40. % This model has one set of eigenvalues that are lightly damped.  They 
  41. % correspond to what is called the Dutch Roll Mode.  We need to design
  42. % a compensator that increases the damping of these poles.  
  43.  
  44. disp('Open Loop Eigenvalues'), damp(A); 
  45. subplot(111)
  46. pzmap(A,[],[],[]); pause % Press any key after plot ...
  47.  
  48. % Our design criteria is to provide damping ratio, zeta > 0.35, with natural
  49. % frequency Wn < 1.0 rad/sec.  We want to design the compensator using
  50. % classical methods. 
  51.  
  52. % Let's do some open loop analysis to determine possible control strategies.
  53.  
  54. % Time response (we could use STEP or IMPULSE here)
  55. impulse(A,B,C,D); pause % Press any key after plot ...
  56.  
  57. % The time responses show that the system is indeed lightly damped.  But the
  58. % time frame is much to long.  Let's look at the response over a smaller
  59. % time frame.  Define the time vector T before invoking IMPULSE.
  60.  
  61. T = 0:0.2:20;    % Define time vector from 0 to 20 secs in steps of 0.2
  62.  
  63. pause % Press any key to continue ...
  64. % Plot responses as separate graphs
  65. %subplot(221), impulse(A,B,C(1,:),D(1,:),1,T); title('Input 1 Output 1')
  66. %subplot(222), impulse(A,B,C(2,:),D(2,:),1,T); title('Input 1 Output 2')
  67. %subplot(223), impulse(A,B,C(1,:),D(1,:),2,T); title('Input 2 Output 1')
  68. %subplot(224), impulse(A,B,C(2,:),D(2,:),2,T); title('Input 2 Output 2')
  69. %pause % Press any key after plot ...
  70. echo off
  71. subplot(221), impulse(A,B,C(1,:),D(1,:),1,T); title('Input 1 Output 1')
  72. subplot(222), impulse(A,B,C(2,:),D(2,:),1,T); title('Input 1 Output 2')
  73. subplot(223), impulse(A,B,C(1,:),D(1,:),2,T); title('Input 2 Output 1')
  74. subplot(224), impulse(A,B,C(2,:),D(2,:),2,T); title('Input 2 Output 2')
  75. pause
  76. echo on
  77.  
  78. % Look at the plot from aileron (input 2) to bank_angle (output 2).  The
  79. % aircraft is oscillating around a non-zero bank angle.  Thus the aircraft
  80. % turns in response to an aileron impulse.  This behavior will be important
  81. % later.
  82.  
  83. pause % Press any key to continue ...
  84.  
  85. % Typically yaw dampers are designed using yaw-rate as the sensed output
  86. % and rudder as the input.  Let's look at that frequency response.
  87. bode(A,B,C(1,:),D(1,:),1);pause % Press any key after plot ...
  88.  
  89. % From this frequency responses we see that the rudder has much effect around
  90. % the lightly damped Dutch roll mode (at 1 rad/sec). 
  91.  
  92. % To make the design easier, extract of the subsystem from rudder to yaw_rate.
  93. [a,b,c,d] = ssselect(A,B,C,D,1,1); % Extract system with input 1 and output 1
  94.  
  95. % Let's do some designs.  The simplest compensator is a gain.  We can determine
  96. % values for this gain using the root locus
  97. clg
  98. rlocus(a,b,c,d); pause % Press any key after plot ...
  99.  
  100. % Oops, looks like we need positive feedback.
  101. rlocus(a,b,-c,-d); sgrid, pause % Press any key after plot ...
  102.  
  103. % That looks better.  So using just simple feedback we can achieve a 
  104. % damping ratio of zeta=0.45.
  105.  
  106. % Now its your turn, Use RLOCFIND to select the point on the root locus
  107. % with maximum damping.
  108.  
  109. pause % Press any key and then select point on the plot ...
  110. [k,poles] = rlocfind(a,b,-c,-d);
  111.  
  112. disp(['You chose gain: ',num2str(k)]), damp(esort(poles));
  113.  
  114. pause % Press any key to continue ...
  115.  
  116. % Let's form the closed loop system so that we can analyze the design.
  117. [ac,bc,cc,dc] = feedback(a,b,c,d,[],[],[],-k);
  118.  
  119. % These eigenvalues should match the ones you chose.
  120. disp('Closed loop eigenvalues'), damp(ac);
  121.  
  122. % Time response using our time vector T
  123. impulse(ac,bc,cc,dc,1,T);  pause % Press any key after plot ...
  124.  
  125. % So the response looks pretty good.  Let's close the loop on the original
  126. % model and see how the response from the aileron looks.  Feedback using
  127. % input 1 and output 1 of plant.
  128.  
  129. % Feedback with selection vectors assumes positive feedback
  130. [Ac,Bc,Cc,Dc] = feedback(A,B,C,D,[],[],[],k,[1],[1]);
  131. disp('Closed loop eigenvalues'), damp(Ac);
  132.  
  133. % Time response
  134. %subplot(221), impulse(Ac,Bc,Cc(1,:),Dc(1,:),1,T); title('Input 1 Output 1')
  135. %subplot(222), impulse(Ac,Bc,Cc(2,:),Dc(2,:),1,T); title('Input 1 Output 2')
  136. %subplot(223), impulse(Ac,Bc,Cc(1,:),Dc(1,:),2,T); title('Input 2 Output 1')
  137. %subplot(224), impulse(Ac,Bc,Cc(2,:),Dc(2,:),2,T); title('Input 2 Output 2')
  138. %pause % Press any key after plot ...
  139. echo off
  140. subplot(221), impulse(Ac,Bc,Cc(1,:),Dc(1,:),1,T); title('Input 1 Output 1')
  141. subplot(222), impulse(Ac,Bc,Cc(2,:),Dc(2,:),1,T); title('Input 1 Output 2')
  142. subplot(223), impulse(Ac,Bc,Cc(1,:),Dc(1,:),2,T); title('Input 2 Output 1')
  143. subplot(224), impulse(Ac,Bc,Cc(2,:),Dc(2,:),2,T); title('Input 2 Output 2')
  144. pause
  145. echo on
  146.  
  147. % Look at the plot from aileron (input 2) to bank_angle (output 2).  When
  148. % we move the aileron the system no longer continues to bank like a normal
  149. % aircraft.  We have over-stabilized the spiral mode.  The spiral mode is
  150. % typically a very slow mode and allows the aircraft to bank and turn without
  151. % constant aileron input.  Pilots are used to this behavior and will not like
  152. % our design if it doesn't allow them to fly normally.  Our design has moved
  153. % the spiral mode so that it has a faster frequency.
  154.  
  155. pause % Press any key to continue ...
  156.  
  157. % What we need to do is make sure the spiral mode doesn't move farther into
  158. % the left half plane when we close the loop.  One way flight control designers
  159. % have fixed this problem is to use a washout filter, i.e.
  160. %              Ks
  161. %    H(s) = --------
  162. %            (s + a)
  163.  
  164. % Choosing a = 0.333 for a time constant of 3 seconds, form the washout
  165. [aw,bw,cw,dw] = zp2ss([0],[-.333],1);
  166.  
  167. % Connect the washout in series with our design model
  168. [a,b,c,d] = series(a,b,c,d,aw,bw,cw,dw);
  169.  
  170. % Do another root locus
  171. clg
  172. rlocus(a,b,-c,-d); sgrid, pause % Press any key after plot
  173.  
  174. % Now the maximum damping is zeta = 0.25  Using RLOCFIND chose the gain for
  175. % maximum damping:
  176.  
  177. pause % Press any key and then select point on plot ...
  178.  
  179. [k,poles] = rlocfind(a,b,-c,-d);
  180. disp(['You chose gain: ',num2str(k)]), damp(esort(poles));
  181.  
  182. % Look at the closed loop response
  183. [ac,bc,cc,dc] = feedback(a,b,c,d,[],[],[],-k);
  184. impulse(ac,bc,cc,dc,1,T);  pause % Press any key to continue ...
  185.  
  186. % Now form the controller (washout + gain) 
  187. [aw,bw,cw,dw] = series(aw,bw,cw,dw,[],[],[],k);
  188.  
  189. % Close the loop with the original model
  190. [Ac,Bc,Cc,Dc] = feedback(A,B,C,D,aw,bw,cw,dw,[1],[1]);
  191.  
  192. % Final closed-loop response
  193. %subplot(221), impulse(Ac,Bc,Cc(1,:),Dc(1,:),1,T); title('Input 1 Output 1')
  194. %subplot(222), impulse(Ac,Bc,Cc(2,:),Dc(2,:),1,T); title('Input 1 Output 2')
  195. %subplot(223), impulse(Ac,Bc,Cc(1,:),Dc(1,:),2,T); title('Input 2 Output 1')
  196. %subplot(224), impulse(Ac,Bc,Cc(2,:),Dc(2,:),2,T); title('Input 2 Output 2')
  197. %pause % Press any key after plot
  198. echo off
  199. subplot(221), impulse(Ac,Bc,Cc(1,:),Dc(1,:),1,T); title('Input 1 Output 1')
  200. subplot(222), impulse(Ac,Bc,Cc(2,:),Dc(2,:),1,T); title('Input 1 Output 2')
  201. subplot(223), impulse(Ac,Bc,Cc(1,:),Dc(1,:),2,T); title('Input 2 Output 1')
  202. subplot(224), impulse(Ac,Bc,Cc(2,:),Dc(2,:),2,T); title('Input 2 Output 2')
  203. pause % Press any key to continue ...
  204. echo on
  205.  
  206. % Although we didn't quite meet the criteria, our design increased the damping
  207. % of the system substantially and the design does allow the pilots to fly the
  208. % aircraft normally.
  209.  
  210. echo off
  211.  
  212.  
  213.