cp10_1.m

离散控制系统设计的MATLAB 代码

%%%%%%%%%%% Comprehensive Problem 10.1 %%%%%%%%%%
%   Discrete-Time Control Problems using        %
%       MATLAB and the Control System Toolbox   %
%   by J.H. Chow, D.K. Frederick, & N.W. Chbat  %
%         Brooks/Cole Publishing Company        %
%                September 2002                 %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Ball and Beam System
% System has 4 states, 1 input, and 2 outputs
% 4 states : ball's position and velocity, and
%            wheel's angular position and velocity
% 1 input  : motor voltage (u)
% 2 output : ball's position (xi), and 
%            wheel's angular position (theta)

disp('CP10.1:Ball and Beam')

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Part A : Load model
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
clear; close all

% Load continuous and discrete-time systems
if ~exist('bbeam.mat'), 
    dbbeam
elseif exist('bbeam.mat') == 2,
    load bbeam 
end

Ts = 0.02;  % (50Hz) in dbbeam already
disp('Loaded dbbeam')
disp('*******>'); pause; disp(' ')

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Part B : full state feedback control design
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
disp('Design of full-state control using pole-placement')
disp(' ')

% determining the desired pole locations
disp('The desired servo-motor transient decay should be about 0.2 s')
disp('that is, the transient should decay to 5% after 10 time periods')
disp('of 0.02 s (Ts). So the 10th root of 0.05 should be a good ')
disp('starting point.')
p_motor = 0.05^0.1
disp('*******>'); pause; disp(' ')

disp('The desired ball position step response should have a rise time')
disp('about 2.5 s and an overshoot shoot of less than 15%.')
disp('This will translate into a continuous-time oscillatory mode')
disp('with a period of 10 (4x2.5) s and a damping ratio of 0.5. ')
disp('Then the continuous complex poles can be mapped into')
disp('discrete-time poles using z = exp(s*Ts).')
w = 2*pi/10
s = -cot(acos(0.5))*w + sqrt(-1)*w
p_ball = exp(s*Ts)
disp('*******>'); pause; disp(' ')

disp('desired poles are')
p_des = [p_motor; p_motor-0.01; p_ball; conj(p_ball)];
xy2p(p_des);
disp('*******>'); pause; disp(' ')

disp('Full-state-feedback gain') 
K = place(Gz.a,Gz.b,p_des)

disp('Closed-loop system model')
Gz_fs = ss(Gz.a - Gz.b*K, Gz.b, Gz.c, Gz.d, Ts);
disp('*******>'); pause; disp(' ')

disp('initial condition response with zero input')
disp('initial condition is [0.1 0 0 0], ball is 0.1 m off center')
x0 = [0.1 0 0 0]'
t = [0:1:500]'*Ts;
y = lsim(Gz_fs,0*t,t,x0);

figure
plot(t,y(:,1),'o',t,y(:,2),'^'); grid  
drawnow
xlabel('Time (s)')
ylabel('Amplitude')
legend('ball position (m)', 'wheel angle (rad)')
title('Full state feedback system response')
disp('*****>');pause

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Part D : Observer design 
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

% assign observer eigenvalues at twice the decay rate of 
% the full-state feedback poles
pObs = [p_motor^2; (p_motor-0.01)^2; p_ball^2; conj(p_ball^2)]

disp('observer gain')
L = place(Gz.a',Gz.c',pObs);
L = L'

disp('*****>');pause
A_f = Gz.a-Gz.b*K-L*Gz.c               % controller system matrix
disp('Controller poles')
xy2p(eig(A_f));                        % controller poles
disp('Controller zeros')
xy2p(tzero(A_f,-L,-K,[0 0]));          % controller zeros
disp('*****>');pause
KObs = ss(A_f,-L,-K,0,Ts);

GK = Gz*KObs;
GObsCL = feedback(GK,[1 0;0 1]);       % closed-loop system
disp('Closed-loop poles')
xy2p(pole(GObsCL));                    % closed-loop system poles
disp('*****>');pause
disp('Initial condition for system and observer')
xObs0 = [0.1 0 0 0 0 0 0 0]'
t = [0:1:500]'*Ts;
y = lsim(GObsCL,t*[0 0],t,xObs0);
figure, plot(t,y(:,1),'o',t,y(:,2),'^'), grid
xlabel('Time (s)')
ylabel('Amplitude')
legend('ball position (m)', 'wheel angle (rad)')
title('Observer feedback system response')
% Done
disp('End of Comprehensive Problem CP10.1')
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%