cp7_2.m

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

%%%%%%%%%%% Comprehensive Problem 7.2 %%%%%%%%%%%
%   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                 %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Electric Power Generation System

clear all, close all
disp('Comprehensive Problem CP7.2')

dpower2    % read data

disp('*****>'), pause
disp('discretized model with a sampling period of 0.01 s')
Ts = 0.01
Gz = c2d(Gs,Ts) 
disp('Note that most of the zero entries of the A & B matrices')
disp('  of the continuous time model have become nonzero in the')
disp('  discrete time models. The C & D matrcies are unchanged.')
disp('*****>'), pause

disp('System model with Vterminal as the only output')
G1 = Gz(1,1)
disp('*****>'), pause

disp('Applying PI controller to Vterminal control')
disp('PI controller is of form')
disp('             (1+K_I T_s) z - 1  ')
disp('K_V(z) = K_P(-----------------) ')
disp('                    z-1         ')

disp('Varying K_I and K_P')
disp('Integral gain K_I')
K_I = [1:1:10]*0.05
len_KI = length(K_I);
disp('*****>'), pause
disp('Proportional gain') 
K_P = [10:5:40]
len_KP = length(K_P);

disp('Sweeping of parameters ...')
km_mat = zeros(len_KP,len_KI);
pm_mat = km_mat;
Mo_mat = km_mat;
tr_mat = km_mat;
ts_mat = km_mat;
t = [0:0.01:30]';
for ii = 1:len_KI
 for jj = 1:len_KP
  K_V = K_P(jj)*tf([(1+K_I(ii)*Ts) -1],[1 -1],Ts);
  [km,pm,gw,pw] = margin(G1*K_V);
  G1cl = feedback(G1*K_V,1);
  y = step(G1cl,t);
  y = 0.05*y;
  [Mo,tp,tr,ts,ess] = kstats(t,y,0.05);
  if isempty(ts) == 1
      ts = max(t);
  end
  if isempty(Mo) == 1
      Mo = 0;
  end
  if isempty(tr) == 1
      tr = 0;
  end
  km_mat(jj,ii) = km;
  pm_mat(jj,ii) = pm;
  Mo_mat(jj,ii) = Mo;
  tr_mat(jj,ii) = tr;
  ts_mat(jj,ii) = ts;
 end
end
disp('*****>'), pause

figure, mesh(K_I,K_P,km_mat)
xlabel('K_I'), ylabel('K_P')
zlabel('Gain margin (dB)') 
title('Gain margin k_m vs K_I and K_P')
disp('Plotting gain margin as a function of gains.')
disp('*****>'), pause


figure, mesh(K_I,K_P,pm_mat)
xlabel('K_I'), ylabel('K_P')
zlabel('Phase margin (deg)') 
title('Phase margin p_m vs K_I and K_P')
disp('Plotting phase margin as a function of gains.')
disp('*****>'), pause

figure, mesh(K_I,K_P,Mo_mat)
xlabel('K_I'), ylabel('K_P')
zlabel('Overshoot (%)') 
title('Overshoot M_o vs K_I and K_P')
disp('Plotting overshoot as a function of gains.')
disp('*****>'), pause

figure, mesh(K_I,K_P,tr_mat)
xlabel('K_I'), ylabel('K_P')
zlabel('Rise time (s)') 
title('Rise time t_r vs K_I and K_P')
disp('Plotting rise time as a function of gains.')
disp('*****>'), pause

figure, mesh(K_I,K_P,ts_mat)
xlabel('K_I'), ylabel('K_P')
zlabel('Settling time (s)') 
title('Settling time t_s vs K_I and K_P')
disp('Plotting settling time as a function of gains.')
disp('*****>'), pause

% selection of gain

figure
[cMo,hMo] = contour(K_I,K_P,Mo_mat,'--'); % plot level contour lines
clabel(cMo,hMo);                          % label contour lines
grid
xlabel('K_I')
ylabel('K_P')
hold on
[ctr,htr] = contour(K_I,K_P,tr_mat);   % plot level contour lines
clabel(ctr,htr);                       % label contour lines
title('t_r in s (solid), M_o in percent (dashed)')
hold off
disp('Contour plot of settling time and overshoot vs K_I and K_P')
disp('*****>'), pause

disp('From the contour plot, to achieve t_r < 1 s and M_o < 10%,')
disp('we need to select K_P = 20 and K_I = 0.35')
disp(' ')

disp('end of Comprehensive Problem 7.2')
%%%%%%%%%%