MATLAB實現控制系統的時域分析
目的
用 MATLAB 對控制系統進行時域分析，包括典型響應、判斷系統穩 定性和分析系統的動態特性。
理論
二、典型響應及其效能分析
1、單位階躍響應
單位階躍響應呼叫格式為：
(1) step(num,den)
(2) step(num,den,t)
(3) step(G)
(4) step(G,t)
該函式將繪製出系統在單位階躍輸入條件下的動態響應圖，同時給出穩態值。其中 t 為影象顯示的時間長度，是使用者指定的時間向量。
如果需要將輸出結果返回到 MATLAB 工作空間中，則採用以下呼叫格式：
c=step(G)
2、求階躍響應的效能指標
MATLAB 提供了強大的繪圖計算功能，可以用多種方法求取系統的動態響應指標。
(1) 最簡單的方法――遊動滑鼠法
在程式執行完畢後，用滑鼠左鍵點選時域響應圖線任意一點，
系統會自動跳出一個小方框， 小方框顯示了這一點的橫座標 （時間）和縱座標 （幅值）。按住滑鼠左鍵在曲線上移動，可以找到曲線幅值最大的一點, 即曲線最大峰值，此時小方框中顯示的時間就是此二階系統的峰值時間，根據觀察到的穩態值和峰值可以計算出系統的超調量。
這種方法簡單易用，但同時應注意它不適用於用 plot()命令畫出的圖形。
(2) 用程式設計方式求取時域響應的各項效能指標
由上面內容可知, 用階躍響應函式 step( )可以獲得系統輸出量，若將輸出量
返回到變數 y 中，可以呼叫如下格式命令：
[y,t]=step(G)
該函式還同時返回了自動生成的時間變數 t，對返回的這一對變數 y 和 t 的
值進行計算，可以得到時域效能指標。
① 峰值時間
[Y,k]=max(y);timetopeak=t(k)
應用取最大值函式 max()求出 y 的峰值及相應的時間，並存於變數 Y 和 k 中。
然後在變數 t 中取出峰值時間，並將它賦給變數 timetopeak。
② 超調量
C=dcgain(G);[Y,k]=max(y);overshoot=100*(Y-C)/C
dcgain( )函式用於求取系統的終值，將終值賦給變數 C，然後依據超調量的
定義，由 Y 和 C 計算出超調量。
③ 上升時間
可以用 while 語句編寫以下程式得到：
C=dcgain(G); n=1;
while y(n)<C
n=n+1;
end
risetime=t(n)在階躍輸入條件下， y 的值由零逐漸增大，當以上迴圈滿足 y=C 時，退出
迴圈，此時對應的時刻，即為上升時間。
對於輸出無超調的系統響應，上升時間定義為輸出從穩態值的 10%上升到
90%所需時間，則計算程式如下：
C=dcgain(G);n=1;
while y(n)<0.1*C
n=n+1;
end
m=1;
while y(m)<0.9*C
m=m+1;
end
risetime=t(m)-t(n)
④ 調節時間(取誤差帶為 0.02)
C=dcgain(G);i=length(t);
while(y(i)>0.98*C)&(y(i)<1.02*C)
i=i-1;
end
setllingtime=t(i)
用向量長度函式 length( )可求得 t 序列的長度，將其設定為變數 i 的上限值。
3、脈衝響應
脈衝響應呼叫格式為：
(1) impulse (num,den)
(2) impulse (num,den,t) 時間向量 t 的範圍可以由人工給定（例如 t=0:0.1:10）
(3) [y,x]=impulse(num,den) 返回變數 y 為輸出向量， x 為狀態向量
(4) [y,x,t]=impulse(num,den,t) 向量 t 表示脈衝響應進行計算的時間
4、高階系統的降階處理
利用閉環主導極點的概念, 可將高階系統進行降階處理。
5、分析零點對系統單位階躍響應的影響
三、系統穩定性分析
MATLAB 中有以下三種方法對系統進行穩定性分析：
1、利用 pzmap 繪製連續系統的零極點圖；
2、利用 tf2zp 求出系統零極點；
3、利用 roots 求分母多項式的根來確定系統的極點
實操
2-1
clc,clear
num=[1,7,24,24];
den=[1,10,35,50,24,0];
[r,p,k]=residue(num,den) 
2-2
clc,clear
num=[1];
den=[0.2,1];
G=tf(num,den);
t=0:0.1:4;
subplot(2,2,1)
step(G,t),
grid on
xlabel('t/s'),ylabel('c(t)'),title('單位階躍響應曲線')
u1=t;
subplot(2,2,2)
lsim(G,u1,t),
grid on
xlabel('t/s'),ylabel('c(t)'),title('單位斜坡響應曲線')
u2=(1/2)*(t.^2);
subplot(2,2,3)
lsim(G,u2,t),
grid on
xlabel('t/s'),ylabel('c(t)'),title('單位加速度響應曲線')
2-3
(1) 最簡單的方法――遊動滑鼠法
num=[10];den=[1 2 10];G=tf(num,den),grid
title('Unit-Step Response of G(s)=10/(s^2+2s+10)')
(2) 用程式設計方式求取時域響應的各項效能指標
clc,clear
num=[10];den=[1 2 10];G=tf(num,den)
Wn=sqrt(10) % 求取無阻尼振盪頻率
E=1/Wn % 求取阻尼比
[r,p,k]=residue(num,den);% 求取閉環極點
p
C=dcgain(G) % 求取穩態值
[y,t]=step(G);
plot(t,y),grid
title('Unit-Step Response of G(s)=10/(s^2+2s+10)')
[Y,k]=max(y) % 計算峰值、峰值時間
timetopeak=t(k)
overshoot=100*(Y-C)/C % 計算超調量
n=1;
while y(n)<C
n=n+1;
end
risetime=t(n) % 計算上升時間
i=length(t);
while(y(i)>0.98*C)&(y(i)<1.02*C)
i=i-1;
end
setllingtime=t(i)% 計算調節時間
2-4
1.
clc,clear
syms Wn;
num1=[1];
den1 =[1 2*0 1];
den2 =[1 2*0.1 1];
den3 =[1 2*0.2 1];
den4 =[1 2*0.3 1];
den5 =[1 2*0.4 1];
den6 =[1 2*0.5 1];
den7 =[1 2*0.6 1];
den8 =[1 2*0.7 1];
den9 =[1 2*0.8 1];
den10=[1 2*1.0 1];
den11=[1 2*2.0 1];
t=[0:0.1:10];
[y1,x,t]=step(num1,den1,t);
[y2,x,t]=step(num1,den2,t);
[y3,x,t]=step(num1,den3,t);
[y4,x,t]=step(num1,den4,t);
[y5,x,t]=step(num1,den5,t);
[y6,x,t]=step(num1,den6,t);
[y7,x,t]=step(num1,den7,t);
[y8,x,t]=step(num1,den8,t);
[y9,x,t]=step(num1,den9,t);
[y10,x,t]=step(num1,den10,t);
[y11,x,t]=step(num1,den11,t);
plot(t,y1,'-',t,y2,'--',t,y3,'--',t,y4,'--',...
t,y5,'--',t,y6,':',t,y7,'--',t,y8,'+',...
t,y9,'--',t,y10,'--',t,y11,'-.'),grid
xlabel('time[sec]'),ylabel('step response')
由圖可知，當 =1rad/s， 分別為 0, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 1.0, 2.0 時的單位階躍響應曲線中，=0.4到0.8時（‘＋’線附近）效能最好。
=0時是無阻尼二階系統，單位階躍響應是一條平均值為1的餘弦形式的等幅振盪曲線，
0<<1的範圍內是欠阻尼二階系統，隨著的增大，單位階躍響應的為衰減的餘弦振盪曲線，振盪幅度逐漸減小；
=1時為臨界阻尼二階系統，單位階躍響應沒有了振盪，是一條穩態值為1的無超調單調上升的曲線；
>1時為過阻尼二階系統，隨著的增大，單位階躍響應為上升的非振盪曲線，斜率隨著的增大，逐漸減小。
在過阻尼和臨界阻尼的響應曲線中，臨界阻尼響應具有最短的上升時間，響應速度最快。在欠阻尼響應曲線中，阻尼比越小，超調量越大，上升時間越短，在=0.4-0.8時超調量適度，調節時間較短，故通常可在此範圍內的選取。
2.
clc,clear
num1=[2.^2];
num2=[4.^2];
num3=[6.^2];
num4=[8.^2];
num5=[10.^2];
num6=[12.^2];
den1 =[1 2*0.6*2 2.^2];
den2 =[1 2*0.6*4 4.^2];
den3 =[1 2*0.6*6 6.^2];
den4 =[1 2*0.6*8 8.^2];
den5 =[1 2*0.6*10 10.^2];
den6 =[1 2*0.6*12 12.^2];
t=[0:0.1:5];
[y1,x,t]=step(num1,den1,t);
[y2,x,t]=step(num2,den2,t);
[y3,x,t]=step(num3,den3,t);
[y4,x,t]=step(num4,den4,t);
[y5,x,t]=step(num5,den5,t);
[y6,x,t]=step(num6,den6,t);
plot(t,y1,'-',t,y2,'--',t,y3,'--',...
t,y4,'--', t,y5,'--',t,y6,':'),grid
xlabel('time[sec]'),ylabel('step response')
由圖可知，當 =0.6, n 分別為 2， 4， 6， 8， 10， 12rad/s 時的單位階躍響應曲線中，隨著n的增大，上升時間越短，峰值時間越快到達，調節時間越短，效能越好。
2-5
clc,clear
num1=[16.41 16.41*1.95];
den1=[1 12 38 52 32];
num2=[1];
den2=[0.5 1 1];
t=[0:0.1:10];
[y1,x,t]=step(num1,den1,t);
[y2,x,t]=step(num2,den2,t);
plot(t,y1,':',t,y2,'--'),grid
xlabel('time[sec]'),ylabel('step response')
2-6
clc,clear
num=[1 2 3 2 1];den=[1 2 6 7 3 1];
G=tf(num,den);
pzmap(G);%利用 pzmap 繪製連續系統的零極點圖；
pole=roots(den)%利用 roots 求分母多項式的根來確定系統的極點
[z,p,k]=tf2zp(num,den)%利用 tf2zp 求出系統零極點；pole =
-0.2939 + 2.1097i
-0.2939 - 2.1097i
-1.0000 + 0.0000i
-0.2061 + 0.4218i
-0.2061 - 0.4218i
z =
-0.5000 + 0.8660i
-0.5000 - 0.8660i
-0.5000 + 0.8660i
-0.5000 - 0.8660i
p =
-0.2939 + 2.1097i
-0.2939 - 2.1097i
-1.0000 + 0.0000i
-0.2061 + 0.4218i
-0.2061 - 0.4218i
k =
1
零極點圖如上圖所示，由圖和計算結果可知，該系統的所有極點均不具
有正實部，故系統穩定。