北京科技大学控制实验报告3

实验三 利用MATLAB进行时域分析

一、实验目的

(1) 学会使用MATLAB编程绘制控制系统的单位阶跃响应曲线；

(2) 研究二阶控制系统中，x、wn对系统动态特性和时域指标的影响； (3) 掌握准确读取动态特性指标的方法；

(4) 分析二阶系统闭环极点和闭环零点对系统动态性能的影响；

(5) 研究三阶系统单位阶跃响应及其动态性能指标与其闭环极点的关系； (6) 研究闭环极点和闭环零点对高阶系统动态性能的影响； (7) 了解高阶系统中主导极点与偶极子的作用；

(8) 了解系统阶跃响应、脉冲响应和斜坡响应输出曲线之间的联系与差别。 二、实验原理及内容 1. 求系统的特征根

若已知系统的特征多项式D (s)，利用roots ( ) 函数可以求其特征根。若已知系统的传递函数，利用eig ( ) 函数可以直接求出系统的特征根。 2、求系统的闭环根、ζ和ωn

函数damp ( ) 可以计算出系统的闭环根、ζ和wn。 3、零极点分布图

可利用pzmap()函数绘制连续系统的零、极点图，从而分析系统的稳定性，调用格式为： pzmap(num,den)

【范例3-1】给定传递函数:

3s42s35s24s6 G(s)5 432s3s4s2s7s2 利用下列命令可自动打开一个图形窗口，显示该系统的零、极点分布图，如图3- 所示。

>> num=[3,2,5,4,6]; den=[1,3,4,2,7,2];

pzmap(num,den)

title(¹Pole-Zero Map¹) % 图形标题。

图3- 1 MATLAB函数零、极点分布图

4、求系统的单位阶跃响应

step ( ) 函数可以计算连续系统单位阶跃响应（impulse( ) 函数可以计算连续系统单位脉冲响应）：

step (sys) 或step ( sys , t ) 或step (num , den)

函数在当前图形窗口中直接绘制出系统的单位阶跃响应曲线，对象sys可以由tf ( )，zpk ( ) 函数中任何一个建立的系统模型。第二种格式中t可以指定一个仿真终止时间，也可以设置为一个时间矢量（如t＝0 : dt : Tfinal，即dt是步长，Tfinal是终止时刻）。

如果需要将输出结果返回到MATLAB工作空间中，则采用以下调用格式： c=step(sys)

此时，屏上不会显示响应曲线，必须利用plot()命令查看响应曲线。plot 可以根据两个或多个给定的向量绘制二维图形。 【范例3-2】已知传递函数为：G(s)25

s24s25利用以下MATLAB命令可得阶跃响应曲线如图3- 所示。 >> num=[0,0,25]; den=[1,4,25];

step(num,den)

grid % 绘制网格线。

title(¹Unit-Step Response of G(s)=25/(s^2+4s+25) ¹) % 图像标题

图3- 2 MATLAB绘制的响应曲线

还可以用下面的语句来得出阶跃响应曲线 >> G=tf([0,0,25],[1,4,25]);

t=0:0.1:5; % 从0到5每隔0.1取一个值。 c=step(G,t); % 动态响应的幅值赋给变量c

plot(t,c) % 绘二维图形，横坐标取t，纵坐标取c。 Css=dcgain(G) % 求取稳态值。

系统显示的图形类似于上一个例子，在命令窗口中显示了如下结果 Css= 1

5、求阶跃响应的性能指标

MATLAB提供了强大的绘图计算功能，可以用多种方法求取系统的动态响应指标。首先介绍一种最简单的方法――游动鼠标法。对于例2，在程序运行完毕后，在曲线中空白区域，单击鼠标右键，在快捷菜单中选择”characteristics”，包含：Peak response (峰值)； settling time (调节时间)；Rise time(上升时间)；steady state(稳态值)；在相应位置出现相应点，用鼠标单击后，相应性能值就显示出来。用鼠标左键点击时域响应曲线任意一点，系统会自动跳出一个小方框，小方框显示了这一点的横坐标（时间）和纵坐标（幅值）。这种方法简单易用，但同时应注意它不适用于用plot()命令画出的图形。 【自我实践1】若已知单位负反馈前向通道的传递函数为：G(s)100，试作出其单位2s5s阶跃响应曲线，准确读出其动态性能指标，并记录数据。

另一种比较常用的方法就是用编程方式求取时域响应的各项性能指标。与游动鼠标法相比，编程方法稍微复杂，但可以获取一些较为复杂的性能指标。

若将阶跃响应函数step( )获得系统输出量返回到变量y中，可以调用如下格式 [y,t]=step(G)

该函数还同时返回了自动生成的时间变量t，对返回变量y和t进行计算，可以得到时域性能指标。

① 峰值时间(timetopeak)可由以下命令获得：

[Y,k]=max(y); timetopeak=t(k)

② 最大(百分比)超调量(percentovershoot)可由以下命令得到：

C=dcgain(G); [Y,k]=max(y);

percentovershoot=100*(Y-C)/C dcgain( )函数用于求取系统的终值。

③ 上升时间(risetime)可利用MATLAB中控制语句编制M文件来获得。

要求出上升时间，可以用while语句编写以下程序得到： C=dcgain(G); n=1;

while y(n)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(n)<0.9*C m=m+1; end

risetime=t(m)-t(n)

④ 调节时间(setllingtime)可由while语句编程得到：

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的上限值。 自我检测1：

>> G1=tf([0,0,100],[1,5,0]); G2=1;

G=feedback(G1,G2)

Transfer function: 100

--------------- s^2 + 5 s + 100

>> num=[0,0,100]; >> den=[1,5,100]; >> step(num,den) >> grid

>> title('Unit-Step Response') >>

上升时间：0.129s；峰值时间：1.42s；调节时间：1.41s；超调：42% 稳态值：1

【范例3-3】已知二阶系统传递函数为：G(s)3

(s13i)(s13i)利用下面的stepanalysis.m程序可得到阶跃响应如图3-及性能指标数据。 >> G=zpk([ ],[-1+3*i,-1-3*i ],3); % 计算最大峰值时间和超调量。 C=dcgain(G) [y,t]=step(G);

plot(t,y) grid

[Y,k]=max(y); timetopeak=t(k)

percentovershoot=100*(Y-C)/C % 计算上升时间。 n=1;

while y(n)risetime=t(n) % 计算调节时间。 i=length(t);

while(y(i)>0.98*C)&(y(i)<1.02*C) i=i-1; end

setllingtime=t(i)

运行后的响应图如图3-，命令窗口中显示的结果为

C = timetopeak = 0.3000 1.0491 percentovershoot = risetime = 35.0914 0.6626 setllingtime = 3.5337

图3-3 二阶系统阶跃响应

用游动鼠标法求取此二阶系统的各项性能指标与本例是一致的。

6、分析ωn不变时，改变阻尼比x，观察闭环极点的变化及其阶跃响应的变化。

【自我实践2】二阶系统，ωn=10，当x＝0，0.25，0.5，0.75，1，1.25时，求对应系统的闭环极点、自然振荡频率及阶跃响应曲线；并分析x对系统性能的影响。 参考程序：

阶跃响应曲线：

阻尼比不同时的阶跃响应曲线2ξ=01.81.6ξ=0.251.41.2ξ=0.5ξ=0.75ξ=1ξ=1.25Amplitude10.80.60.40.2000.20.40.60.81Time (sec)1.21.41.61.82

自我实践2： >> num=100;i=0;

for sigma = 0:0.25:1.25 den = [1 2*sigma*10 100]; damp(den)

sys = tf(num,den); i=i+1;

step(sys,2) hold on end

Eigenvalue Damping Freq. (rad/s) 0.00e+000 + 1.00e+001i 0.00e+000 1.00e+001 0.00e+000 - 1.00e+001i 0.00e+000 1.00e+001 Eigenvalue Damping Freq. (rad/s) -2.50e+000 + 9.68e+000i 2.50e-001 1.00e+001 -2.50e+000 - 9.68e+000i 2.50e-001 1.00e+001 Eigenvalue Damping Freq. (rad/s) -5.00e+000 + 8.66e+000i 5.00e-001 1.00e+001 -5.00e+000 - 8.66e+000i 5.00e-001 1.00e+001 Eigenvalue Damping Freq. (rad/s) -7.50e+000 + 6.61e+000i 7.50e-001 1.00e+001 -7.50e+000 - 6.61e+000i 7.50e-001 1.00e+001 Eigenvalue Damping Freq. (rad/s) -1.00e+001 1.00e+000 1.00e+001 -1.00e+001 1.00e+000 1.00e+001 Eigenvalue Damping Freq. (rad/s) -2.00e+001 1.00e+000 2.00e+001 -5.00e+000 1.00e+000 5.00e+000 >> grid hold off

title('阻尼比不同时的阶跃响应曲线')

lab1='zunibi = 0';text(0.3,1.9,lab1), lab2='zunibi = 0.25';text(0.3,1.5,lab2), >> lab3='zunibi = 0.5';text(0.3,1.2,lab3), lab4='zunibi = 0.75';text(0,3,1.05,lab4), lab5='zunibi = 1';text(0,35,0.9,lab5), lab6='zunibi = 1.25';text(0,35,0.8,lab6) >>

ωn不变时，改变阻尼比x，当ξ>1时，系统为过阻尼系统，系统的阶跃响应为非震荡过程，瞬态特性为单调变化曲线，无超调和震荡；当0<ξ<1时，系统为欠阻尼系统，系统的阶跃响应为非震荡过程，ξ越小，超调量越大，震荡次数越多，调节时间越长。当ξ=0时，系统为零阻尼系统，系统的阶跃响应为持续的等幅震荡。当ξ<0时，输出量做发散震荡。 7、保持x＝0.25不变，分析ωn变化时，闭环极点对系统单位阶跃响应的影响。

【自我实践3】二阶系统，x＝0.25，当ωn=10，30，50时，求系统的阶跃响应曲线；并分析ωn对系统性能的影响。 参考程序：

阶跃响应曲线：

wn 变化时系统的阶跃响应曲线1.5wn=10wn=30wn=501Amplitude0.5000.20.40.60.81Time (sec)1.21.41.61.82

自我实践3：

>> sgma=0.25;i=0; >> for wn=10:20:50

num=wn^2;den=[1,2*sgma*wn,wn^2]; sys = tf(num,den); i=i+1;

step(sys,2) hold on,grid end

>> hold off

>> title('wn变化时系统的阶跃响应曲线') >> lab1='wn=10';text(0.35,1.4,lab1), >> lab2='wn=30';text(0.12,1.3,lab2), >> lab3='wn=50';text(0.05,1.2,lab3) >>

结论：当ξ不变时，ωn 越大，峰值时间越短，调节时间越短，上升时间越短，超调量不变。

【综合实践】通过分别改变典型二阶系统的ξ和ωn，观察系统在脉冲、阶跃、斜坡信号作用下的响应特性，求时域指标，总结参数对系统性能影响的规律。 ξ ξ>1 设ξ=2 ωn 0.2 1 0.2 1 0.2 1 0.2 ts（s） tp（s） 78.3 15.7 46.4 9.29 4.26 0.852 σ% tr（s） 响应曲线 参数影响说明 当ξ>1时，系统为过阻尼系统，系统的阶跃响应为非震荡过程，瞬态特性为单调变化曲线，无超调和震荡；当0<ξ<1时，系统为欠阻尼系统，系统的阶跃响应为非震荡过程，ξ越0<ξ<1 6.36 1.27 脉设ξ=0.5 冲 ξ=0 -1<ξ<0

设ξ=-0.5 ξ<-1 设ξ=-2 1 0.2 1 0.2 ξ>1 1 0.2 0<ξ<1 1 8.08 3.6 16.3 1.64 14.9 40.4 18 16.3 8.23 8.22 74.4 41.2 小，超调量越大，震荡次数越多，调节时间越长。当ξ=0时，系统为零阻尼系统，系统的阶跃响应为持续的等幅震荡。当ξ<0时，输出量做发散震荡。当ξ不变时，ωn 越大，峰值时间越短，调节时间越短，上升时间越短，超调量不1 0.2 -1<ξ<0 1 0.2 ξ<-1 1 变。 阶跃 0.2 ξ=0

响应曲线：脉冲 ξ=2 ωn = 0.2

ξ=2 ωn = 1

ξ=0.5 ωn = 0.2

ξ=0.5 ωn = 1

ξ=0 ωn = 0.2

ξ=0 ωn = 1

ξ=-0.5 ωn =0.2

ξ=-0.5 ωn =1

ξ=-2 ωn =0.2

ξ=-2 ωn =1

ξ=2 ωn =0.2

ξ=2， ωn =1

ξ=0.5 ωn =0.2

ξ=0.5 ωn =1

ξ=0 ωn =0.2

ξ=0 ωn =1

ξ=-0.5 ωn =0.2

ξ=-0.5 ωn =1

ξ=-2 ωn =0.2

ξ=-2 ωn =1

8、分析系统零、极点对系统阶跃响应的影响。

【自我实践4】试作出以下系统的阶跃响应，并与原系统G(s)线进行比较，作出实验结果分析。 num=[0,0,10]; den=[1,2,10]; step(num,den)

10的阶跃响应曲

s22s10

1）系统有零点情况：z=-5，即G1(s)2s10； 2s2s10num=[0,2,10]; den=[1,2,10]; step(num,den)

结果分析：增加零点，阶跃响应上升时间变短，峰值时间变短，超调量变小，调节时间变短，系统稳态值不变。

s20.5s102) 分子与分母多项式阶数相等：n=m=2，G2(s)2；

s2s10num=[1,0.5,10]; den=[1,2,10]; step(num,den)

结果分析：分子与分母多项式阶数相等，系统从稳态震荡，峰值时间和调节时间均减小，稳态值不变。

s20.5s3) 分子多项式零次项系数为0，G3(s)2；

s2s10num=[1,0.5,0]; den=[1,2,10]; step(num,den)

结果分析：零初值的系统，零时刻对应的纵坐标为零，非零初值系统零时刻对应的纵坐标不为零。

4) 原系统的微分响应，微分系数为1/10，G4(s)s。 2s2s10num=[0,1,0]; den=[1,2,10]; step(num,den)

结果分析：若没有零点，对于单位阶跃响应的稳态值为1，而分子部分出现零值时，单位阶跃响应的稳态值为0。

【综合实践】附加零点的影响。

设原系统为：

R(s) Y(s)

W0(s) -

图

3-4 原系统结构图 附加开环零点的情况：

R(s) - P(s) Y(s) W0(s)

图3-5 附加开环零点的系统结构图

附加闭环零点的情况为：

R(s) P(s) —Y(s) W0(s) 图3-6 附加闭环零点的系统结构图

2n其中 P(S)=(TS+1)、W0(s)2。这里取ωn=1、ξ=0.5

s2ns 零点参数 闭环传递函数 单位阶跃响应曲线 说明零点的影响 T=0.2 附加开环零点 T=1 0.2 s + 1 --------- s^2 + 1.2 s + 1 s + 1 ------------- s^2 + 2 s + 1 0.2 s + 1 ----------- s^2 + s + 1 s + 1 ----------- s^2 + s + 1 T=0.2 附加闭环零点 T=1 增加系统零点会使系统稳定性增加，但会使调节时间变长。系统零点并不改变原有系统的稳定性，但会影响系统的动态参数。零点越大，上升时间和峰值时间以及调节时间越短，超调量越大，稳态值不变。 p=tf([0.2,1],[1]); W=tf([1],[1,1,0]); G=feedback(p*W,1) step(G)

Transfer function: 0.2 s + 1 --------------- s^2 + 1.2 s + 1

p=tf([1,1],[1]); W=tf([1],[1,1,0]); G=feedback(p*W,1) step(G)

Transfer function: s + 1 ------------- s^2 + 2 s + 1

clc

p=tf([0.2,1],[1]); W=tf([1],[1,1,0]); G=feedback(W,1); G1=p*G step(p*G)

Transfer function: 0.2 s + 1 ----------- s^2 + s + 1

clc

p=tf([1,1],[1]); W=tf([1],[1,1,0]); G=feedback(W,1); G1=p*G step(p*G)

Transfer function: s + 1 ----------- s^2 + s + 1 >>

【综合实践】附加极点的影响。

当附加零点中的函数变为：P(S)=1/(TS+1)，则上图3-4、3-5、3-6成为附加极点的情况。假设取ωn=1、ξ=0.5，用单位阶跃信号作为系统输入，按照下表的要求输入参数，记录仿真曲线于下表。 极点参数 闭环传递函数 1 ------------------- 0.2 s^3 + 1.2 s^2 + s + 1 1 ------------------- s^3 + 2 s^2 + s + 1 响应曲线 说明极点的影响 T=0.2 附加开环极点 T=1 系统的稳定性不变，上升时间、超调时间、峰值时间均变长，超调量减小。

1 ------------------- 0.2 s^3 + 1.2 s^2 + 1.2 s + 1 1 ------------------- s^3 + 2 s^2 + 2 s + 1 T=0.2 附加闭环极点 T=1 clc;

p=tf([1],[0.2,1]); W=tf([1],[1,1,0]); G=feedback(p*W,1) step(G)

Transfer function: 1

------------------------- 0.2 s^3 + 1.2 s^2 + s + 1 >>

clc;

p=tf([1],[1,1]); W=tf([1],[1,1,0]); G=feedback(p*W,1) step(G)

Transfer function: 1

------------------- s^3 + 2 s^2 + s + 1

clc

p=tf([1],[0.2,1]); W=tf([1],[1,1,0]); G=feedback(W,1); G1=p*G step(p*G)

Transfer function: 1

----------------------------- 0.2 s^3 + 1.2 s^2 + 1.2 s + 1 >>

clc

p=tf([1],[1,1]); W=tf([1],[1,1,0]); G=feedback(W,1); G1=p*G step(p*G)

Transfer function: 1

--------------------- s^3 + 2 s^2 + 2 s + 1

9、三阶系统的单位阶跃响应分析

研究三阶系统单位阶跃响应及其动态性能指标与其闭环极点的关系。 【自我实践5】已知三阶系统闭环传递函数为(s)5(s2)(s3)， 2(s4)(s2s2)

① 编写MATLAB程序，求取系统闭环极点及其单位阶跃响应，读取动态性能指标。

clc; Z=[-2,-3];

P=[-4,-1+j,-1-j]; K=5;

G=zpk(Z,P,K) eig(G) step(G)

Zero/pole/gain: 5 (s+2) (s+3) -------------------- (s+4) (s^2 + 2s + 2)

ans =

-4.0000 -1.0000 + 1.0000i -1.0000

- 1.0000i

② 改变系统闭环极点的位置， (s)0.625(s2)(s3)，即将原极点s = - 4改成2(s0.5)(s2s2)s= - 0.5，使闭环极点靠近虚轴，观察单位阶跃响应和动态性能指标变化。

clc; Z=[-2,-3];

P=[-0.5,-1+j,-1-j];

K=0.625; G=zpk(Z,P,K) eig(G) step(G)

Zero/pole/gain: 0.625 (s+2) (s+3) ---------------------- (s+0.5) (s^2 + 2s + 2)

ans =

-0.5000 -1.0000 + 1.0000i -1.0000 - 1.0000i

③ 改变系统闭环零点的位置，(s)10(s1)(s3)，即将原零点s = - 2改成s = - 2(s4)(s2s2)

1，观察单位阶跃响应及其动态性能指标的变化。

clc; Z=[-1,-3];

P=[-4,-1+j,-1-j]; K=10;

G=zpk(Z,P,K) eig(G) step(G)

Zero/pole/gain: 10 (s+1) (s+3) -------------------- (s+4) (s^2 + 2s + 2)

ans =

-4.0000 -1.0000 + 1.0000i -1.0000 - 1.0000i

④ 分析上面实验，给出结论。

如果系统的闭环极点离虚轴远，则相应的瞬态分量就衰减得快，系统的调节时间短。当系统闭环极点接近虚轴时，系统阻尼增加，系统超调量减小，调节时间变长。当闭环零点减小后，相当于减小了系统阻尼，系统峰值时间减小，超调时间减小，超调量增大。 10、高阶系统的单位阶跃响应分析

【自我实践6】已知控制系统的闭环传递函数

1(s)1.05(0.4762s1)

(0.125s1)(0.5s1)(s2s1)① 用MATLAB软件分析该系统的单位阶跃响应及其动态性能指标。

clc;

num=[0,0,0,0.50001,1.05];

den=[0.0625,0.6875,1.6875,1.625,1]; G=tf(num,den) zpk(G) step(G)

Transfer function:

0.5 s + 1.05

------------------------------------------------- 0.0625 s^4 + 0.6875 s^3 + 1.688 s^2 + 1.625 s + 1

Zero/pole/gain: 8.0002 (s+2.1)

------------------------- (s+8) (s+2) (s^2 + s + 1) >>

② 将该系统的阶跃响应与二阶系统2(s)主导极点的特点及作用，

1.05的单位阶跃响应比较分析闭环系统

s2s1clc;

num=[0,0,1.05]; den=[1,1,1];

G=tf(num,den) zpk(G) step(G)

Transfer function: 1.05 ----------- s^2 + s + 1

Zero/pole/gain: 1.05

------------- (s^2 + s + 1)

主导极点的特点及作用：极点的实部的模大，它所决定的瞬态分量不仅持续时间最长，其初始值也最大，充分体现了它在系统响应中的主导作用，而其他闭环极点产生的响应分量随时间的推移迅速衰减，对系统响应过程影响甚微。 ③ 比较系统3(s)1.05和 2(0.125s1)(0.5s1)(ss1)

4(s)1.05(s1)的单位阶跃响应及其动态性能指标，观察闭环零2(0.125s1)(0.5s1)(ss1)点对系统动态性能产生的影响有哪些？

clc;

num=[0,0,0,0,1.05];

den=[0.0625,0.6875,1.6875,1.625,1]; G=tf(num,den) zpk(G) step(G)

Transfer function:

1.05

------------------------------------------------- 0.0625 s^4 + 0.6875 s^3 + 1.688 s^2 + 1.625 s + 1

Zero/pole/gain: 16.8

------------------------- (s+8) (s+2) (s^2 + s + 1)

clc;

num=[0,0,0,1.05,1.05];

den=[0.0625,0.6875,1.6875,1.625,1]; G=tf(num,den) zpk(G) step(G)

Transfer function:

1.05 s + 1.05

------------------------------------------------- 0.0625 s^4 + 0.6875 s^3 + 1.688 s^2 + 1.625 s + 1

Zero/pole/gain: 16.8 (s+1)

-------------------------

(s+8) (s+2) (s^2 + s + 1)

闭环零点对系统的影响：减小峰值时间，使响应速度加快，超调量增大。零点会减小系统阻尼，并且这种作用将随闭环零点接近虚轴而加剧。， ④ 比较系统1(s)和5(s)1.05(0.4762s1)的单位阶跃响应及其动态性能指标，分

(0.5s1)(s2s1)析非主导极点对系统动态性能的影响及作用。

clc;

num=[0,0,0.50001,1.05]; den=[0.5,1.5,1.5,1]; G=tf(num,den) zpk(G) step(G)

Transfer function: 0.5 s + 1.05

-----------------------------

0.5 s^3 + 1.5 s^2 + 1.5 s + 1

Zero/pole/gain: (s+2.1)

------------------- (s+2) (s^2 + s + 1)

非主导极点对系统动态性能的影响及作用：增大峰值时间，使系统响应速度变缓，超调量减小。闭环非主导极点可以增大系统的阻尼，这种作用随闭环极点接近虚轴而加剧。 ⑤ 比较系统2(s)和5(s)的单位阶跃响应及其动态性能指标，分析偶极子对系统动态性能的影响及作用。

偶极子远离虚轴时，其对系统的动态性能的影响可以忽略；靠近虚轴时，其对系统的动态性能的影响不能忽略。