实验二利用MATLAB进行时域分析报告

1、实验二 利用 MATLAB 进行时域分析本实验容包含以下三个部分：基于MATLAB的线性系统稳定性分析、基于MATLAB勺线性 系统动态性能分析、和 MATALB进行控制系统时域分析的一些其它实例。一、 基于MATLAB的线性系统稳定性分析线性系统稳定的充要条件是系统的特征根均位于S平面的左半部分。系统的零极点模型可以直接被用来判断系统的稳定性。另外，MATLAB语言中提供了有关多项式的操作函数，也可以用于系统的分析和计算。(1) 直接求特征多项式的根设p为特征多项式的系数向量，则MATLAB数roots()可以直接求出方程 p=0在复数围的解 v, 该函数的调用格式为：v=roots(p)例

2、 3.1 已知系统的特征多项式为：x5 3x3 2x2 x 1特征方程的解可由下面的MATLAB命令得出。>> p=1,0,3,2,1,1;v=roots(p)结果显示：v =0.3202 + 1.7042i0.3202 - 1.7042i-0.72090.0402 + 0.6780i0.0402 - 0.6780i利用多项式求根函数 roots(), 可以很方便的求出系统的零点和极点，然后根据零极点 分析系统稳定性和其它性能。(2) 由根创建多项式如果已知多项式的因式分解式或特征根，可由MATLABI数poly()直接得出特征多项式系数向量，其调用格式为：p=poly(v)如上例

3、中：v=0.3202+1.7042i;0.3202-1.7042i; -0.7209;0.0402+0.6780i; 0.0402-0.6780i;>> p=poly(v)结果显示p =1.0000 0.00013.0000 2.00010.9998 0.9999由此可见，函数roots()与函数poly()是互为逆运算的。(3) 多项式求值在MATLAB过函数polyval()可以求得多项式在给定点的值，该函数的调用格式为:polyval(p,v)对于上例中的p值，求取多项式在x点的值，可输入如下命令：>> p=1,0,3,2,1,1;x=1polyval(p,x)结

4、果显示x =1ans =8(4) 部分分式展开 考虑下列传递函数：M (s)numb0snb1sn1bnN(s)dena0sna1sn1an式中a。0，但是ai和bj中某些量可能为零。MATLAB函数可将展开成部分分式，直接求出展开式中的留数、极点和余项。该N(s)函数的调用格式为：r, p, kresidue(num, den)则器的部分分式展开由下式给出:k(s) s P(n)M (s)r(1)r(2)N(s) s p(1) s p(2)式中p(1)P1, p(2)P2,，p( n)Pn,为极点,r(1)r1 , r(2)r?,r(n)rn为各极点的留数，k(s)为余项。例3.2设传递函数

5、为:32G(s)2s 5s 3s 6s3 6s211s 6该传递函数的部分分式展开由以下命令获得:>> nu m=2,5,3,6;den=1,6,11,6;r,p,k=residue( nu m,de n) 命令窗口中显示如下结果r =-6.0000-4.00003.0000P =-3.0000-2.0000-1.00002中留数为列向量r，极点为列向量 p,余项为行向量 k。由此可得出部分分式展开式：G(s) 6432s 3 s 2 s 1该函数也可以逆向调用，把部分分式展开转变回多项式M(s)之比的形式，命令格式为:N(s)nu m,de n=residue(r,p,k)对上例

6、有：>> nu m,de n=residue(r,p,k)结果显示num =2.0000 5.0000 3.0000 6.0000den =1.0000 6.0000 11.0000 6.0000应当指出，如果p(j)=p(j+1)=p(j+m-1),则极点p(j)是一个m重极点。在这种情况下，部分分式展开式将包括下列诸项：r(j)s p( j)r(j 1)s p( j)r(j m 1) s p(j)例3.3设传递函数为:G(s)s2 2s 3(s 1)3S2 2s 3s3 3s2 3s 1则部分分式展开由以下命令获得:>> v=-1,-1,-1num=0,1,2,3;

7、den=poly(v);r,p,k=residue( nu m,de n)结果显示v =-1 -1 -11.00000.00002.0000P =-1.0000-1.0000-1.0000k =其中由poly()命令将分母化为标准降幕排列多项式系数向量den, k=为空矩阵。由上可得展开式为：G(s)0(s 1)22(s 1)3(5) 由传递函数求零点和极点。在MATLABS制系统工具箱中，给出了由传递函数对象G求出系统零点和极点的函数，其调用格式分别为：Z=tzero(G)P=G .P 1注意：式19中要求的G必须是零极点模型对象，且出现了矩阵的点运算 “.”和大括号訂表示的矩阵元素，详细容

8、参阅后面章节。例3.4已知传递函数为：2、6.8s61.2s 95.2G(s) 3s4 7.5s3 22s219.5s输入如下命令：num=6.8,61.2,95.2;den=1,7.5,22,19.5,0;G=tf( nu m,de n);G1=zpk(G);Z=tzero(G)P=G1.P1结果显示Z =-7-20-3.0000 + 2.0000i-3.0000 - 2.0000i-1.5000(6) 零极点分布图在MATLAB中，可利用pzmap()函数绘制连续系统的零、极点图，从而分析系统的稳 定性，该函数调用格式为：pzmap( nu m,de n)例3.5给定传递函数：3s4 2s

9、3 5s24s 6G(s) 5432s 3s 4s 2s 7s 2利用下列命令可自动打开一个图形窗口，显示该系统的零、极点分布图，如图3-1所示。>> nu m=3,2,5,4,6;den=1,3,4,2,7,2;pzmap (nu m,de n)title(1 Pole-Zero Map 1) % 图形标题。5 0 5o o- cfxa wanraa-1-1.5-2-1.5Pole-Zero Map-1-0.5Real Axis0.5系统动态特性分析。控制系统软件包提供了控制系统工程需要的基本的时域与频域分析工具函数。ImpulseStepLsimBode连续时间系统分析函数脉冲

10、响应阶跃响应任意输入的仿真 波特图Nyquist奈奎斯特图Lyap雅普诺夫方程Gram可控性与可观性(1)时域响应解析算法部分分式展开法用拉氏变换法求系统的单位阶跃响应，可直接得出输出c(t)随时间t变化的规律，对于高阶系统，输出的拉氏变换象函数为：1num 1C(s) G(s)sden s对函数c(s)进行部分分式展开，我们可以用num,den,0来表示c(s)的分子和分母。例3.6给定系统的传递函数：s3 7s224s24432s 10s35s50s 24用以下命令对g©进行部分分式展开。s>> num,den,0num=1,7,24,24den=1,10,35,50

11、,24r,p,k=residue( nu m,de n, 0) 输出结果为num =172424den =110355024-1.00002.0000-1.0000-1.00001.0000P =-4.0000-3.0000-2.0000-1.00000k =输出函数c(s)为：C (s) =num/ (den*s )num =172424 ans =1103550240拉氏变换得：L=laplace(sA3+7*s+24*s+24)/(sA4+10*sA3+35*sA2+50*s+24),s,t)(2 )单位阶跃响应的求法：如果已知控制系统工具箱中给出了一个函数step()来直接求取线性系统

12、的阶跃响应,传递函数为：G(s)numden则该函数可有以下几种调用格式:(b)step (nu m,de n) step (nu m,de n,t)step(G)(c)step(G,t)(d)该函数将绘制出系统在单位阶跃输入条件下的动态响应图，同时给出稳态值。对于式(b)和(d)，t为图像显示的时间长度，是用户指定的时间向量。式(a)和(c)的显示时间由系统根据输出曲线的形状自行设定。如果需要将输出结果返回到MATLABT作空间中，则采用以下调用格式：c=step(G)此时，屏上不会显示响应曲线，必须利用plot()命令去查看响应曲线。plot可以根据两个或多个给定的向量绘制二维图形，详细介

13、绍可以查阅后面的章节。例3.7已知传递函数为：G(s)25s2 4s 25利用以下MATLAB命令可得阶跃响应曲线如图3-2所示。unit-siqu恐讯曲篦or g【s：户筑 锢沪普)图3-2 MATLAB绘制的响应曲线>> nu m=0,0,25;den=1,4,25;step( nu m,de n)grid %绘制网格线title( 1Unit-Step Response of G(s)=25/(sA2+4s+25)1) % 图像标题我们还可以用下面的语句来得出阶跃响应曲线>> G=tf(0,0,25,1,4,25);t=0:0.1:5;%从0到5每隔0.1取一个值

14、。c=step(G,t);%动态响应的幅值赋给变量cplot(t,c)%绘二维图形，横坐标取t,纵坐标取c。Css=dcgai n(G) %求取稳态值。系统显示的图形类似于上一个例子，在命令窗口中显示了如下结果Css=1(3 )求阶跃响应的性能指标MATLAB提供了强大的绘图计算功能，可以用多种方法求取系统的动态响应指标。我们首先介绍一种最简单的方法一一游动鼠标法。对于例16,在程序运行完毕后，用鼠标左键点击时域响应图线任意一点， 系统会自动跳出一个小方框， 小方框显示了这一点的横坐标 (时 间)和纵坐标(幅值)。按住鼠标左键在曲线上移动，可以找到曲线幅值最大的一点一一即 曲线最大峰值, 此时

15、小方框中显示的时间就是此二阶系统的峰值时间, 根据观察到的稳态值 和峰值可以计算出系统的超调量。 系统的上升时间和稳态响应时间可以依此类推。 这种方法 简单易用,但同时应注意它不适用于用 plot() 命令画出的图形。另一种比较常用的方法就是用编程方式求取时域响应的各项性能指标。与上一段介绍的游动鼠标法相比,编程方法稍微复杂,但通过下面的学习,读者可以掌握一定的编程技巧, 能够将控制原理知识和编程方法相结合, 自己编写一些程序, 获取一些较为复杂的性能指标。通过前面的学习,我们已经可以用阶跃响应函数 step( ) 获得系统输出量,若将输出量 返回到变量 y 中,可以调用如下格式y,t=ste

16、p(G)该函数还同时返回了自动生成的时间变量t，对返回的这一对变量 y和t的值进行计算，可以得到时域性能指标。 峰值时间 (timetopeak) 可由以下命令获得：Y,k=max(y);timetopeak=t(k)应用取最大值函数 max()求出y的峰值及相应的时间，并存于变量Y和k中。然后在变量 t 中取出峰值时间,并将它赋给变量 timetopeak 。 最大 ( 百分比 ) 超调量 (percentovershoot) 可由以下命令得到：C=dcgain(G);Y,k=max(y);percentovershoot=100*(Y-C)/Cdcga in()函数用于求取系统的终值，将终

17、值赋给变量C，然后依据超调量的定义，由 Y 和 C 计算出百分比超调量。 上升时间(risetime) 可利用MATLAB控制语句编制M文件来获得。首先简单介绍 一下循环语句 while 的使用。while 循环语句的一般格式为：while< 循环判断语句 >循环体end其中 ,循环判断语句为某种形式的逻辑判断表达式。当表达式的逻辑值为真时, 就执行循环体的语句； 当表达式的逻辑值为假时, 就退出当 前的循环体。 如果循环判断语句为矩阵时, 当且仅当所有的矩阵元素非零时, 逻辑表达式的 值为真。为避免循环语句陷入死循环,在语句必须有可以自动修改循环控制变量的命令。要求出上升时间,可

18、以用 while 语句编写以下程序得到：C=dcgain(G);n=1;while y(n)<Cn=n+1;endrisetime=t(n)在阶跃输入条件下，y的值由零逐渐增大，当以上循环满足y=C时，退出循环，此时对应的时刻,即为上升时间。对于输出无超调的系统响应，上升时间定义为输出从稳态值的10%上升到90%所需时间,则计算程序如下：C=dcgai n( G);n=1;while y( n)<0.1*Cn=n+1;endm=1;while y( n)<0.9*Cm=m+1;endrisetime=t(m)-t( n) 调节时间(setllingtime)可由while语句

19、编程得到：C=dcgai n( G);i=le ngth(t);while(y(i)>0.98*C)&(y(i)<1.02*C)i=i-1;endsetlli ngtime=t(i)用向量长度函数length() 可求得t序列的长度，将其设定为变量i的上限值。例3.8已知二阶系统传递函数为：G(s)(s 13i)(s 13i)利用下面的stepanalysis.m程序可得到阶跃响应如图3-3及性能指标数据。G(S)=3/(SA2+2S+10)>> nu m=0,0,3;den=1,2,10;step( nu m,de n)gridtitle(1Un it-Ste

20、p Respo nse of G(s)=3/(sA2+2s+10)1)或者G=tf(0,0,3,1,2,10);t=0:0.1:5; c=step(G,t);plot(t,c)>> Css=dcga in(G)Css =0.3000G=tf(0,0,3,1,2,10);>> y,t=step(G)00.00440.01680.03610.06100.09020.12230.15620.19070.22480.25740.28790.31550.33980.36040.37710.38980.39870.40370.40530.40360.39920.39230.3835

21、0.37330.36200.35000.33790.32600.31460.30390.29430.28580.27860.27280.26830.26530.26360.26310.26370.26530.26770.27080.27440.27840.28260.28680.29100.29500.29870.30210.30510.30760.30960.31120.31220.31280.31300.31270.31220.31130.31020.30890.30750.30610.30460.30310.30170.30040.29920.29820.29730.29660.2961

22、0.29570.29550.29550.29550.29570.29600.29640.29690.29740.29790.29840.29890.29940.29990.30030.30060.30100.30120.30140.30150.30160.30160.30160.30150.30140.30120.30110.30090.30070.30060.30040.30020.30000.29990.299800.05520.11040.16560.22090.27610.33130.38650.44170.49690.55210.60740.66260.71780.77300.82820.88340.93860.99391.04911.10431.15951.21471.26991.32521.38041.43561.49081.54601.60121.65641.71171.76691.82211.87731.93251.98772.04292.09822.15342.20862.26382.31902.37422