1、 基于 MATLAB 的控制系统频域设计 姓名： 学院： 专业： 班级： 学号： 基于 MATLAB 的控制系统频域设计 一 实验目的 1. 利用计算机作出开环系统的波特图 2. 观察记录控制系统的开环频率特性 3. 控制系统的开环频率特性分析 二 预习要点 1. 预习 Bode 图和 Nyquist 图的画法； 2. Nyquist 稳定性判据内容。 三 实验方法 1、奈奎斯特图（幅相频率特性图） 对于频率特性函 数 G(jw)，给出 w 从负无穷到正无穷的一系列数值，分别求出 Im(G(jw)和 Re(G(jw)。以 Re(G(jw) 为横坐标， Im(G(jw) 为纵坐标绘制成为极坐标频
2、率特性图。 MATLAB 提供了函数 nyquist()来绘制系统的极坐标图，其用法如下： nyquist(a,b,c,d)：绘制出系统的一组 Nyquist 曲线，每条曲线相应于连续状态空间系统 a,b,c,d的输入 /输出组合对。其中频率范围由函数自动选取，而且在响应快速变化的位置会自动采用更多取样点。 nyquist(a,b,c,d,iu)：可得到从系统第 iu 个输入到所有输出的极坐标图。 nyquist(num,den)：可绘制出以连续时间多项式传递函数表示的系统的极坐标图。 nyquist(a,b,c,d,iu,w)或 nyquist(num,den,w)：可利用指定的角频率矢量绘
3、制出系统的极坐标图。 当不带返回参数时，直接在屏幕上绘制出系统的极坐标图（图上用箭头表示 w 的变化方向，负无穷到正无穷） 。当带输出变量 re,im,w引用函数时，可得到系统频率特性函数的实部 re 和虚部 im 及角频率点 w 矢量（为正的部分）。可以用 plot(re,im)绘制出对应 w 从负无穷到零变化的部分。 2、对数频率特性图（波特图） 对数频率特性图包括了对数幅频特性图和对数相频特性图。横坐标为频率 w，采用对数分度，单位为弧度 /秒；纵坐标均匀分度，分别为幅值函数 20lgA(w)，以 dB 表示；相角，以度表示。 MATLAB 提供了函数 bode()来绘制系统的波特图，其
4、用法如下： bode(a,b,c,d,iu)：可得到从系统第 iu 个输入到所有输出的波特图。 bode(a,求取系统对数频率特性图（波特图）： bode() 求取系统奈奎斯特图（幅相曲线图或极坐标图）： nyquist() b,c,d)：自动绘制出系统的一组 Bode 图，它们是针对连续状态空间系统 a,b,c,d的每个输入的 Bode 图。其中频率范围由函数自动选取，而且在响应快速变化的位置会自动采用更多取样点。 bode(num,den)：可绘制出以连续时间多项式传递函数表示的系统的波特图。 bode(a,b,c,d,iu,w)或 bode(num,den,w)：可利用指定的角频率矢量绘
5、制出系统的波特图。 当带输出变量 mag,pha,w或 mag,pha引用函数时，可得到系统波特图相应的幅值 mag、相角 pha 及 角频率点 w 矢量或只是返回幅值与相角。相角以度为单位，幅值可转换为分贝单位： magdb=20 log10(mag) 四 实验内容 下面举例说明 用 MATLAB 对控制系统频域设计 1用 Matlab 作 Bode 图 .画出对应 Bode 图 , 并加标题 . （ 1）25425)( 2 sssGnum=25;den=1 4 25; G=tf(num,den); figure(1) margin(G); figure(2) nichols(G); axi
6、s(-207 0 -40 40);ngrid figure(3) nyquist(G); axis equal B o d e D i a g r a mG m = In f d B ( a t In f r a d / s e c ) , P m = 6 8 . 9 d e g ( a t 5 . 8 3 r a d / s e c )F r e q u e n c y ( r a d / s e c )10-1100101102- 1 8 0- 1 3 5- 9 0- 4 50S y s t e m : GF r e q u e n c y ( r a d / s e c ) : 5 .
7、9 5P h a s e ( d e g ) : - 1 1 3Phase(deg)- 6 0- 4 0- 2 0020S y s t e m : GF r e q u e n c y ( r a d / s e c ) : 5 . 9 2M a g n i t u d e ( d B ) : - 0 . 3 2 1S y s t e m : GF r e q u e n c y ( r a d / s e c ) : 1 . 4 5M a g n i t u d e ( d B ) : 0 . 5 0 3Magnitude(dB)（ 2）)92.1( )12.0(9)( 22 sss sss
8、Gnum=conv(0 1,1 0.2 1);den=conv(1 0,1 1.2 9); G=tf(num,den); figure(1) margin(G); figure(2) nichols(G); axis(-207 0 -40 40);ngrid figure(3) nyquist(G); axis equal - 4 0- 3 0- 2 0- 1 0010Magnitude(dB)10-1100101102- 9 0- 4 504590Phase(deg)B o d e D i a g r a mG m = In f , P m = 9 0 . 4 d e g ( a t 0 .
9、 1 1 r a d / s e c )F r e q u e n c y ( r a d / s e c )2用 Matlab 作 Nyquist 图 .画对应 Nyquist 图 ,并加网格标题 . 18.01)( 2 sssG num=1;den=1 0.8 1; G=tf(num,den); figure(1) margin(G); figure(2) nichols(G); axis(-207 0 -40 40);ngrid figure(3) nyquist(G); axis equal- 1 . 5 -1 - 0 . 5 0 0 . 5 1 1 . 5 2- 1 . 5-1- 0
10、 . 500 . 511 . 50 d B- 2 0 d B- 1 0 d B- 6 d B- 4 d B- 2 d B2 0 d B1 0 d B6 d B4 d B2 d BN y q u i s t D i a g r a mR e a l A x i sImaginaryAxis3典型二阶系统2222)(nnn sssG ，试绘制 取不同值时的 Bode 图。取0.1:1.0:1.0,6 n 。 当 w=6, =0.1 时 num=36;den=1 1.2 36; G=tf(num,den); figure(1) margin(G); - 6 0- 4 0- 2 0020S y s t
11、 e m : GF r e q u e n c y ( r a d / s e c ) : 5 . 7 8M a g n i t u d e ( d B ) : 1 3 . 6Magnitude(dB)100101102- 1 8 0- 1 3 5- 9 0- 4 50S y s t e m : GF r e q u e n c y ( r a d / s e c ) : 8 . 4 1P h a s e ( d e g ) : - 1 6 4Phase(deg)B o d e D i a g r a mG m = In f d B ( a t In f r a d / s e c ) , P
12、 m = 1 6 . 3 d e g ( a t 8 . 4 r a d / s e c )F r e q u e n c y ( r a d / s e c )当 w=6, =1.0 时 num=36;den=1 12 36; G=tf(num,den); figure(1) margin(G); - 5 0- 4 0- 3 0- 2 0- 1 00Magnitude(dB)10-1100101102- 1 8 0- 1 3 5- 9 0- 4 50Phase(deg)B o d e D i a g r a mG m = In f d B ( a t In f r a d / s e c )
13、 , P m = - 1 8 0 d e g ( a t 0 r a d / s e c )F r e q u e n c y ( r a d / s e c )num=50;den=conv(1 5,1 -2); G=tf(num,den); figure(1) margin(G); figure(2) nichols(G); axis(-207 0 -40 40);ngrid figure(3) nyquist(G); axis equal hold on impulse(G) 4某开环传函为：)2)(5( 50)( sssG，试绘制系统的 Nyquist 曲线，并判断闭环系统稳定性，最后
14、求出闭环系统的单位脉冲响应。 num=50;den=conv(1 5,1 -2); G=tf(num,den); figure(1) margin(G); figure(2) nichols(G); axis(-207 0 -40 40);ngrid figure(3) nyquist(G); axis equal -5 - 4 . 5 -4 - 3 . 5 -3 - 2 . 5 -2 - 1 . 5 -1 - 0 . 5 0- 1 . 5-1- 0 . 500 . 511 . 50 d B- 2 0 d B- 1 0 d B- 6 d B- 4 d B- 2 d B2 0 d B1 0 d
15、B6 d B4 d B2 d BN y q u i s t D i a g r a mR e a l A x i sImaginaryAxis有奈奎斯特图可知它有左半平面的开环极点，也可看出他包围了（ -1， j0），所以系统不稳定。 0 1 2 3 4 5 6 7 801234567x 1 07单位脉冲响应T i m e ( s e c )Amplitude5 1.0,5.0,1,2 1.0,12122 TTssTsG当 T=0.1， =2 时 num=1;den=0.01 0.4 1; G=tf(num,den); figure(1) margin(G); figure(2) nichol
16、s(G); axis(-207 0 -40 40);ngrid figure(3) nyquist(G); axis equal title(波特图 ) -1 - 0 . 8 - 0 . 6 - 0 . 4 - 0 . 2 0 0 . 2 0 . 4 0 . 6 0 . 8 1- 0 . 8- 0 . 6- 0 . 4- 0 . 200 . 20 . 40 . 60 . 80 d B- 2 0 d B- 1 0 d B- 6 d B- 4 d B- 2 d B2 0 d B1 0 d B6 d B4 d B 2 d B波特图R e a l A x i sImaginaryAxis当 T=0.1
17、， =1 时 -1 - 0 . 8 - 0 . 6 - 0 . 4 - 0 . 2 0 0 . 2 0 . 4 0 . 6 0 . 8 1- 0 . 8- 0 . 6- 0 . 4- 0 . 200 . 20 . 40 . 60 . 80 d B- 2 0 d B- 1 0 d B- 6 d B- 4 d B- 2 d B2 0 d B1 0 d B6 d B4 d B 2 d B波特图R e a l A x i sImaginaryAxis当 T=0.1， =0.5 时 - 1 . 5 -1 - 0 . 5 0 0 . 5 1 1 . 5- 1 . 5-1- 0 . 500 . 511 .
18、50 d B- 2 0 d B- 1 0 d B- 6 d B- 4 d B- 2 d B2 0 d B1 0 d B6 d B4 d B2 d B波特图R e a l A x i sImaginaryAxis当 T=0.1， =0.1 时； num=1;den=0.01 0.02 1; G=tf(num,den); figure(1) margin(G); figure(2) nichols(G); axis(-207 0 -40 40);ngrid figure(3) nyquist(G); axis equal title(波 特图 ) -6 -4 -2 0 2 4 6-6-4-2024
19、60 d B- 1 0 d B- 6 d B- 4 d B- 2 d B1 0 d B6 d B4 d B2 d B波特图R e a l A x i sImaginaryAxis6 11.0101.0 6.31 ssssG要求： (a) 作波特图 num=31.6;den=conv(1 0,0.01 1)； den=conv(den,0.1 1); G=tf(num,den); figure(1) margin(G); - 1 5 0- 1 0 0- 5 0050100Magnitude(dB)10-1100101102103104- 2 7 0- 2 2 5- 1 8 0- 1 3 5- 9
20、 0Phase(deg)B o d e D i a g r a mG m = 1 0 . 8 d B ( a t 3 1 . 6 r a d / s e c ) , P m = 2 2 . 3 d e g ( a t 1 6 . 3 r a d / s e c )F r e q u e n c y ( r a d / s e c )(b) 由稳定裕度命令计算系统的稳定裕度 gL 和 c ，并确定系统的稳定性 K=1/0.1; G0=zpk(,0 -100 -10,K); lg, c,wx,wc=margin(G0) lg =1.1000e+004 c= 89.9370 wx = 31.6228 wc = 0.0100 - 3 . 5 -3 - 2 . 5 -2 - 1 . 5 -1 - 0 . 5 0- 5 0- 4 0- 3 0- 2 0- 1 001020304050N y q u i s t D i a g r a mR e a l A x i sImaginaryAxis
