VIEWS: 48 PAGES: 13 POSTED ON: 5/12/2010 Public Domain
Exercise 18 WB 2421 2004: SISO system IMPORTANT: To understand how a particular Matlab command works, type help followed by the name of the command. Example: help plot. USE M.file hinf_siso.m PART 1: Representation of systems in Matlab. One way to represent a Single Input Single Output system is through its transfer function, i.e. given its numerator and denominator polynomials. A polynomial in Matlab is represented by the vector of its coefficients, and the multiplication of polynomials is done with the command conv. Example: To represent the system s 2 10s 13 P( s ) ( s 1 )( s 2 ) use the commands: num=[1 10 13]; % numerator of P den=conv([1 1],[1 2]); % denominator of P P=nd2sys(num,den); The command nd2sys converts a SISO transfer function representation into the system matrix, which is the structure used in the Mu-Toolbox to represent a system. Now you do: Consider the model of the motion system: x2 c F M1 M2 H x2 2 b g ds c 0.02 d b 2 g b g b F s J1 J 2 s J1 J 2 ds J1 J 2 c g where J1 J 2 4.8e 6 c 0.22 d 1.e 4 . Build the system G in Matlab, determine the zero-pole-gain form, the transfer function form, and the state-space realization of G. numG=0.02*([d c]); denG=[(J1*J2) d*(J1+J2) c*(J1+J2) 0 0]; Gmu=nd2sys(numG,denG); % representation according to the Mu- % Tools G=tf(numG,denG) % transfer function form of G ZP_G=zpk(tf_G) % zero-pole-gain form of G 1 [A,B,C,D]=unpck(G) % state-space realization of G % Calculate the poles and the zeros and examine the properties of the system G plant_poles=spoles(G) % poles of the system plant_zeros=szeros(G) % zeros of the system minfo(G) % Returns structure information about G seesys(G) % Returns the state-space realization of G PART 2: Calculate and plot frequency responses of the system G. In order to calculate a frequency response of a system, first a frequency vector can be defined. This is done with the command LOGSPACE. Once the frequency vector has been defined, the frequency response can be calculated with the command FRSP. Notice that this command requires a frequency vector in rad/sec and not in Hertz! Also the commands bode and nyquist can be used directly! % Perform the following: hz=logspace(-5,5,200); % Frequency vector in Hz of 200 5 5 % points from 10 Hz to 10 Hz wz=2*pi*hz; % Frequency vector in rad/s ghat=frsp(G,wz); % varying frequency response % matrix of G figure(1); vplot('bode',ghat);grid on; % Bode plot of G subplot(211);grid on;title('Bode plot of G'); % Now suppose that we close the loop around G with a unitary negative feedback. Can you say whether the closed loop system is stable or not? ref + out G - To do this plot de Nyquist plot of the open loop system G and check the Nyquist stability requirement figure(2) subplot(211); vplot('nyq',ghat);grid on; % Nyquist plot of G title('Nyquist plot of G'); subplot(212); 2 nyquist(tf_G,{1,1000});grid on % plots the Nyquist plot between % 1 [rad/s] and 1000 [rad/s]. % zoom on to see a part of the plot near % the (–1)-point. PART 3: Calculate and plot time response of the system G In order to calculate a time response of a system, first a time vector and an input signal should be defined. The time can be defined with the command LINSPACE. The input signal is defined with the command SIGGEN. Then the time response of the system for the defined input is obtained using the command TRSP and can be plotted using VPLOT with the option 'iv,m'. Also the commands step, gensig and lsim can be used. % Plot the response of G to a step input: t=linspace(0,5,2000); % Time vector of 2000 points from % 0 to 5 sec. ustep=siggen('1',t); % Define step input. ystep=trsp(G,ustep,5); % step response of G within 5[sec] figure(3) subplot(211); vplot('iv,m',ystep);grid on; title('step response') % Plot the response of G to a sine of frequency of 5 Hz: usin=siggen('sin(2*pi*5*t)',t); %Define sine input. ysin=trsp(G,usin,5); %sine response of G subplot(212); vplot('iv,m',ysin);grid on; title('response to sinusoidal input') PART 4: Analysis of the closed loop system Consider the classical closed loop scheme of figure 1. ref + out K G - Figure 1 Now you do: Based on the analysis of the Bode and Nyquist plots of PART 2, try to determine a PD controller K for which the closed loop system is stable with a phase 3 margin of at least 30 [deg] and with a peak of the sensitivity function not bigger than 2. To do this, first determine the interconnected system in open loop form of figure 2, and then close the loop with the controller. z z2 1 w + y u G _ z1 w z2 P y u Figure 2: generalized plant structure % Import the following to build the interconnection: systemnames='Gmu'; % Defines the system blocks in the structure % of figure 2 inputvar='[w;u]'; % Defines the input variables in the structure outputvar='[w-Gmu;Gmu;w-Gmu]'; % Defines the output variables z1,z2,y in the % structure input_to_Gmu='[u]'; % Specifies the input of the plant G sysoutname='Pmu'; % Gives the name of the interconnection cleanupsysic='yes'; sysic; % Runs the interconnection routine seesys(Pmu) % Returns the structure of Pmu minfo(Pmu) % Convert to lti format [a,b,c,d]=unpck(Pmu); P=ss(a,b,c,d); P.InputName={'w';'u'}; P.OutputName={'z1';'z2';'y'}; 4 % Close the loop with a PD controller K as shown in figure 4. z1 w z2 P y u K %Figure 4: closed loop system. % Import the following: % Put here your choice of K=Kp+Kv*s %Kp=?; %Kv=?; Kp=0.004; Kv=0.00016; numk=[Kv Kp]; denk=[1.e-8 1]; K=tf(numk,denk); Pcl=lft(P,K); % Closed loop system matrix % the commando lft is the Linear Fractional % transformation of the plant P with controller K % Type "help lft" for more information. % For our specific interconnection we have % where S is sensitivity and T is % complementary sensitivity % Calculate the poles of the closed loop system: % Is the closed loop stable? closed_loop_poles=eig(Pcl) % Determine and plot the frequency responses of the following functions of the closed loop system: The complementary sensitivity function: T The sensitivity function: S T=Pcl('z2','w'); % The commando selects subsystems from a Pcl 5 S=Pcl('z1','w'); That=frsp(T,wz); % frequency response of T Shat=frsp(S,wz); % frequency response of S figure(4); subplot(311); vplot('liv,lm',That);grid on; title('Complementary sensitivity: T') subplot(312);grid on; vplot('liv,lm',Shat);grid on; title('Sensitivity: S') % Determine the peak values of the sensitivity and the complementary sensitivity functions: Tpeak=norm(T,inf) Speak=norm(S,inf) % Determine and plot the step response of the closed loop system. Is there a steady-state error? yclstep=trsp(T,ustep,5); % step response of the closed loop subplot(313); vplot('iv,m',yclstep);grid on; title('closed_loop step response') PART 5: Robust stability analysis for the PD controller. Remark: The H-infinity norm is a suitable measure to assess robustness of a control system. In Matlab it can be calculated using the command HINFNORM or NORM. Suppose that the plant is not perfectly known and that we can model it as the nominal value G with an output multiplicative uncertainty as: b g G p ( s ) 1 Wu ( s )( s ) G( s ) where Wu ( s ) gives the profile of the uncertainty as a function of the frequency. Suppose: Wu ( s ) b 10 s 0.01 g s 100 Determine whether the controlled system is robustly stable against the defined uncertainty? For robust stability the following H-infinity condition should be satisfied: 6 Wu T 1 The following block scheme gives the weighted plant without controller: z2 z1 Wu w + y u G _ Figure 5 ~ % Determine the weighted plant P of figure 5: % Import the weight Wu numwu=10*[1 0.01]; denwu=[1 100]; wu=nd2sys(numwu,denwu); % Determine the weighted complementary sensitivity function and its H -norm: norm(wu*T,inf) figure(5) subplot(211); vplot('liv,lm',That);grid on;hold on; title('Unweighted complementary sensitivity & Inverse weighting') vplot('liv,lm',wu_invhat,'r');grid on; subplot(212); vplot('liv,lm',wuThat);grid on; title('Weighed complementary sensitivity wuT') % Is the controlled system robustly % stable? 7 % Are the plots in agreement with the % H-infinity norm? PART 6: Robust performance analysis for the PD controller. Consider the uncertain system of PART5 and assume that the performance of the closed loop system is specified in terms of the weighted sensitivity function. Take as the performance weight Wp the inverse of the sensitivity function computed in PART4. Is the closed loop robustly performing for every uncertainty such that 1? For robust performance the following H-infinity condition should hold: Wu T 1. Wp S The following block scheme gives the weighed plant with Wp and Wu: z1 z2 Wp Wu w + y u G _ Figure 6 % Determine the weighted plant tP of figure 6, close the loop with the constant % controller, and Determine the H norm of the controlled sytem: wp=minv(S); wp_inv=S; wp_invhat=Shat; Wm=daug(wp,wu,1); % The daug command groups systems in diagonal form. tP=mmult(Wm,P); seesys(tP) % Close the loop 8 tN=starp(tP,K); seesys(tN); tNhat=frsp(tN,wz); hinfnorm(tNhat) % Is the closed loop system robustly performing? PART 7: S-design The performance that can be achieved by a simple constant controller is very poor. For instance the bandwidth of the closed loop system is very low. Shaping the sensitivity function S(s) by H control is a suitable design tool to improve performance. Typical specifications that are related to the sensitivity function are the bandwidth of the closed loop system, the tracking error and disturbance attenuation. Consider a performance weighting function of the form: 2 2 s 2ns n Wp ( s a ) ( s b) Select appropriate parameters n , , a and b in order to satisfy the following specifications: 1. Bandwidth between 5 and 15 rad/s. (as an example) 2. Disturbance attenuation at low frequencies. 3. Zero steady-state error of the step response. Perform an S-design that achieves the above specifications: The block scheme on which the S-design is based is given in figure 7: z1 z2 Wp w + y u G _ Figure 7 Import the following in Matlab: %select the appropriate values of the parameters: 9 wb=?; % Desired bandwidth z=?; a=?; b=?; numwp=[1 2*z*wb wb^2]; denwp=conv([1 a],[1 b]); Wp=tf(numwp,denwp); % Plot the weight function Wp: wphat=frsp(wp,w); figure(6); vplot('liv,lm',wphat);grid on; title('Weight function Wp') % Determine the weighted plant tP W=append(Wp,1,1); WP= W*P; [a,b,c,d]=ssdata(WP); WPmu=pck(a,b,c,d); seesys(WPmu) seesys(tP) % Given the weighted plant WP, compute the Hinf controller K(s) which stabilizes the plan tP and which yields a closed-loop peak smaller than gopt. gopt is the smallest achievable peak value: [gopt,Kmu] = hinflmi(MPmu,[1 1],2,1e-3); % The command hinflmi computes the % H-infinity contoller according to the % LMI-Tools. % tP : Weighted open-loop system % [1,1]: 1=Number of measured variables % 1=Number of controlled signals % 2 : Try to achieve peak value smaller % than 2 % 1e-3 : Accuracy % Determine the state-space realization of the controller K, the transfer function, % the zero-pole-gain form, and the poles of K. [Ak,Bk,Ck,Dk]=unpck(Kmu); % state-space realization K=ss(Ak,Bk,Ck,Dk); [numk,denk]=SS2TF(Ak,Bk,Ck,Dk); % numerator and denominator tf_k=tf(numk,denk) ZP_kS=zpk(tf_k) K_poles=spoles(K) % Plot the magnitude of the controller and interpret what kind of action it 10 % performs. Plot also the magnitudes of the model G, and of the open loop system GK: Khat=frsp(K,wz); GK=mmult(G,K); GKhat=frsp(GK,wz); figure(7) vplot('liv,lm',ghat); grid on;hold on; vplot('liv,lm',Khat,'r');hold on vplot('liv,lm',GKhat,'b'); title('Model G, controller K, open loop GK : S-design') % Analysis of controlled system % Close the loop of the plant P (not weighted) with the H-infinity controller K. % Calculate the closed loop poles Pcl2=lft(P,K); seesys(Pcl2) cl_poles=spoles(Pcl2) % Plot the frequency response of the following: % * Sensitivity function S % * Inverse of weight Wp % * Complementary sensitivity T % * Inverse of weight Wu % * Weighted sensitivity WpS % * Weighted complementary sensitivity % * Achieved performance gopt % Is the controlled system robustly stable against the multiplicative % uncertainty of PART 5? norm(Wu*Pcl2('z2','w'),inf) % Plot the step response of the closed loop system: figure(9);grid on step(Pcl('z2','w'),5) title('closed_loop step response: S-design') % Is the steady state error % equal to zero? % Is the steady state error % equal to zero? PART 8: S/T-design: The mixed S/T design allows to combine performance requirements, expresse d in terms of the sensitivity function, and robustness requirements, expressed in terms of the complementary sensitivity function. In this design, the objective is the minimization of: 11 F SI W GT J HK W p u where Wu reflects the size of the multiplicative uncertainty. Assume Wu as given in PART 5 and Wu as in PART 7. Perform an S/T-design that achieves the specifications in PART 7: The mixed S/T-design is based on the following scheme: z1 z2 Wp Wu w + y u G _ Figure(8) % Determine the weighted plant WmPmu with Wp and Wu as in figure 8 Wm=append(Wp,Wu,1); WmP=Wm*P; [a,b,c,d]=ssdata(WmP); WmPmu=pck(a,b,c,d); seesys(WmPmu) % Compute the Hinfinity controller K [gopt,K2mu] = hinflmi(WmPmu,[1 1],2,1e-3); [Ak2,Bk2,Ck2,Dk2]=unpck(K2mu); % state-space realization K2=ss(Ak2,Bk2,Ck2,Dk2); [numk2,denk2]=tfdata(tf(K2)); % numerator and denominator [K2_poles,K2_zeros]=pzmap(K2) % Plot the magnitude of the controller and interpret what kind of action it performs. Plot also the magnitudes of the model G, and of the open loop system GK: Khat=frsp(K2,wz); GK2=mmult(G,K2); GKhat=frsp(GK2,wz); 12 figure(10) vplot('liv,lm',ghat); grid on;hold on; vplot('liv,lm',Khat,'r');hold on vplot('liv,lm',GKhat,'b'); title('Model G, controller K, open loop Gk : S/T-design') % What is the difference between the controller of the S-design and of % the S/T-design? ( use e.g. magnitude plots and/or the zero-pole form of % the controllers) % Analysis of controlled system % Close the loop of the plant P (not weighted) with the H-infinity controller K. % Calculate the closed loop poles Pcl3=lft(P,K2); seesys(Pcl3) cl_poles=spoles(Pcl3) % Plot the frequency response of the following: % * Sensitivity function S % * Inverse of weight Wp % * Complementary sensitivity T % * Inverse of weight Wu % * The weighted sensitivity WpS % * The weighted complementary sensitivity WuT % * The achieved performance gopt % Plot the step response of the closed loop system: % Examine robust stability and robust performance % Compare the performance and robustness aspects of S-design and S/T-design % Plot the frequency response of KS, the transfer function from the reference % signal to the output of the controller. % Draw conclusions 13