1、第 1 章 绪论一、例题例 1线性系统的建模仿真:开环控制系统;闭环控制系统。解 开环控制系统运行后可得下图:0 1 2 3 4 5 6 7 8 9 10-5-4-3-2-1012345闭环控制系统运行后得下图:0 1 2 3 4 5 6 7 8 9 10-5-4-3-2-1012345例 2非线性系统的建模仿真:开环控制系统;闭环控制系统。 、解 开环控制系统运行后得下图:0 1 2 3 4 5 6 7 8 9 10-5-4-3-2-1012345闭环控制系统运行后得下图:0 1 2 3 4 5 6 7 8 9 10-505二、仿真下图为在 Simulink 工具里面的搭建的仿真模块,实现控
2、制的稳定性。图 1.1 控制系统结构模型图对模型中的数据进行合理的设计,运行图形如下:图 1.2 控制系统结构波形图分析:由图示结果看出较为稳定,超调量小,调节时间也很短。在 t=0.2s时基本达到稳定。第 2 章 自动控制系统的数学模型一、例题例 12 两个子系统为13()4Gs24()3sG将两个系统按并联方式连接,可输入:num1=3;den1=1,4;num2=2,4;den2=1,2,3;num,den=parallel(num1,den1,num2,den2)则得num =0 5 18 25den =1 6 11 12因此 2123518()()6sGs例 13 两个子系统为251
3、()3sG5(2)10sH将两个系统按反馈方式连接,可输入numg=2 5 1;deng=1 2 3;numh=5 10;denh=1 10;num,den=feedback(numg,deng,numh,denh)则得num =2 25 51 10den =11 57 78 40因此闭环系统的传递函数为 32()510()1784cnumssGde二、仿真系统 1 为: , 系统 2 为求按串联、并联、正反馈、负反馈连接时的系统状态方程及系统 1 按单位负反馈连接时的状态方程。编写程序如下:clca1=0 1;-1 -2;b1=0;1;c1=1 3;d1=1;a2=0 1;-1 -3;b2=
4、0;1;c2=1 4;d2=0;a,b,c,d=series(a1,b1,c1,d1,a2,b2,c2,d2) %串联连接a,b,c,d=parallel(a1,b1,c1,d1,a2,b2,c2,d2) %并联连接a,b,c,d=feedback(a1,b1,c1,d1,a2,b2,c2,d2,+1) %正反馈连接a,b,c,d=feedback(a1,b1,c1,d1,a2,b2,c2,d2) %负反馈连接a,b,c,d=cloop(a1,b1,c1,d1) %单位负反馈连接运行结果如下(仅列出串联连接的结果):串联连接a = b =0 1 0 0 0-1 -3 1 3 10 0 0 1
5、00 0 -1 -2 1c = 1 4 0 0d =0分析:可见在 MATLAB 中,可以用程序建立各种数学模型,而且可以进行各类数学模型见的转换,因此利用 MATLAB 建立数学模型应用较为广泛。第 3 章 时域分析法要判断系统的稳定性,只需要确定系统闭环极点在 s 平面上的分布。利用MATLAB 命令可以快速求出闭环系统零极点并绘制其零极点图,也可以方便绘出系统的时间响应曲线。因此,利用 MATLAB 可以可以方便、快捷地对控制系统进行时域分析。例 18 已知连续系统的传递函数为432546()7ssG要求:求出该系统的零、极点及增益。绘出起零、极点图,判断系统稳定性。解 可执行如下程序:
6、%This program creates a transfer function and then finds/displays its poles、zerosand % gainnum=3,2,5,4,6;den=1,3,4,2,7,2;z,p,k=tf2zp(num,den);pzmap(num,den);title(Poles and zeros map)程序执行结果如下:z=0.4019+101965i p= -1.7680+1.2673i0.4019-101965i 1.7680-1.2673i-0.7352+0.8455i 0.4176+1.1130i-0.7352-0.8455
7、i 0.4176 -1.1130i-0.2991K=3同时屏幕上显示系统的零极点分布图(如图所示)-2 -1.5 -1 -0.5 0 0.5-1.5-1-0.500.511.5 0.140.280.420.560.70.820.910.9750.140.280.420.560.70.820.910.9750.250.50.7511.251.51.752Poles and zeros mapReal AxisImaginaryAxis分析:无论是从所求的系统零、极点,还是绘制出的零、极点图,都可以看到系统中有两个极点位于 s 的右半平面,因此,该系统不稳定。例 19 已知典型二阶系统的传递函数为
8、22()nGss式中 =6,绘图系统在 =0.1,0.2,1.0,2.0 时的单位阶跃响应。n解 可执行如下程序:%This program plots a curve of step responsewn=6;kosi=0.1,0.2,1.0,2.0;figure(1)hold onfor kos=kosinum=wn.2;den=1,2*kos*wn,wn.2;step(num,den)endtitle(Step Response)hold off0 1 2 3 4 5 6 7 8 9 1000.20.40.60.811.21.41.61.8 step responseTime (sec)
9、Amplitude从图中可以看出,在过阻尼和临界阻尼曲线中,临界阻尼的响应具有最短的上升时间,响应速度最快;在欠阻尼的的响应曲线中,阻尼系数越小,超调量越大,上升时间越短。例 20 已知三阶系统的传递函数为3210()().4.0.4sGs绘制系统的单位阶跃响应和单位脉冲响应曲线。解 可执行如下程序:% This program plots a curve of step response and step impulse for three order % systemclfnum=100 200;den=1 1.4 100.44 100.04;h=tf(num,den);y,t,x=ste
10、p(h)y1,t1,x1=impulse(h)subplot(211),plot(t ,y)title(Step Response)xlabel(time),ylabel(amplitude)subplot(212),plot(t1,y1)title(impulse response)xlabel(time),ylabel(amplitude)0 5 10 15 20 250123 Step Responsetimeamplitude0 5 10 15 20 25 30-1001020 impulse responsetimeamplitude第 4 章 根轨迹法一、例题例 9 已知单位反馈系
11、统的开环传递函数为2(3)()56)sGsK试在系统的闭环的根轨迹图上选择一点,求出该点的增益 K 及其系统的闭环极点位置,并判定在该点系统的闭环稳定性。解 调用 rlocfind( )函数,Matlab 程序为:num=1 3;den=conv(conv(conv(1 0,1 5),1 6),1 2 2);sys=tf(num,den);rlocus(sys)k,poles=rlocfind(sys)title(根轨迹分析)xlabel(实轴)ylabel(虚轴)-15 -10 -5 0 5 10-10-8-6-4-20246810 冲冲冲冲冲冲冲冲冲执行程序后用光标在根轨迹图上选一点,可得
12、到相应的该点的系统的增益和其闭环极点:k =82.3756poles =-5.6903 + 1.2000i-5.6903 - 1.2000i-2.2680 0.3243 + 1.7654i0.3243 - 1.7654i二,仿真例 1. 某开环系统传递函数如 ,要求绘制系统的闭环22)34()sksGo根轨迹,分析其稳定性,并绘制出当 k=55 和 k=56 时系统的闭环冲激响应。解:可执行如下程序:clcnumo=1 2;den=1 4 3;deno=conv(den,den);figure(1)k=0:0.1:150;rlocus(numo,deno,k)title(root locus)
13、p,z=pzmap(numo,deno);k,p1=rlocfind(numo,deno); %求出系统临界稳定增益kfigure(2) %验证系统的稳定性subplot(211)k=55;num2=k*1 2;den=1 4 3;den2=conv(den,den);numc,denc=cloop(num2,den2,-1);impulse(numc,denc)title(impulse response k=55);subplot(212)k=56;num3=k*1 2;den=1 4 3;den3=conv(den,den);numcc,dencc=cloop(num3,den3,-1)
14、;impulse(numcc,dencc)title(impulse response k=56);程序执行结果如下图所示:运行后系统的闭环根轨迹如图 4.1 所示:root locusReal AxisImaginaryAxis-8 -7 -6 -5 -4 -3 -2 -1 0 1-5-4-3-2-1012345图 4.1.闭环系统根轨迹图执行程序后,用光标在根轨迹图上选一点,可得相应的该点的系统增益。如:selected_point =0.1789 - 3.4627i k =72.2648selected_point = -0.1836 + 2.7795i k =39.5736(注:运行了
15、两次,选取了两个不同的值。K=72.2648 时系统不稳定;k=39.5736 时,系统稳定)同时,通过绘制 k=55 和 k=56 系统的闭环冲激响应曲线,验证其稳定性。所得图如下:0 200 400 600 800 1000 1200-4-2024 impulse response k=55Time (sec)Amplitude0 50 100 150 200 250 300 350 400-40-2002040 impulse response k=56Time (sec)Amplitude图 4.2.系统的闭环冲激响应曲线分析:当 k=55 时,系统根轨迹处于 s 左半平面,即其所有闭
16、环极点的实部均为负值,所以在该点处,闭环系统是稳定的,如图 4.2 上所示。当 k=56 时,系统根轨迹处于 s 右半平面,其闭环极点的实部有正值,所以在该点处,闭环系统是不稳定的,如图 4.2 下所示。第 5 章 频域分析法例 11 有一个二阶系统,其自然频率 =1,阻尼因子 =0.2,要绘n制出系统的幅频和相频曲线,可输入:a,b,c,d=ord(1,0.2);bode(a,b,c,d);title(Bode Plot)执行后得到下图:-40-30-20-10010Magnitude (dB)10-1 100 101-180-135-90-450Phase (deg)Bode PlotFr
17、equency (rad/sec)例 12 典型二阶系统:2()nkGss绘制出 取不同值时的 Bode 图。解 取 =6, 取0.1:1.0 时二阶系统的 Bode 图可直接采用 Bode 得到。nMATLAB 程序为%Example 5.1%wn=6;kosi=0.1:1.0;w=logspace(-1,1,100);figure(1)num=wn.2;for kos=kosiden=1 2*kos*wn wn.2;mag,pha,w1=bode(num,den,w);subplot(2,1,1);hold onsemilogx(w1,mag);subplot(2,1,2);hold on
18、semilogx(w1,pha);endsubplot(2,1,1);grid ontitle(Bode Plot);xlabel(Frequency(rad/sec);ylabel(Gain dB);subplot(2,1,2);grid onxlabel(Frequency(rad/sec);ylabel(Phase deg);hold off执行后得下图:0 1 2 3 4 5 6 7 8 9 100246 Bode PlotFrequency(rad/sec)GaindB0 1 2 3 4 5 6 7 8 9 10-200-150-100-500Frequency(rad/sec)Ph
19、ase deg例 13 有系统:210(4)().5ksGs绘制出系统的 Bode 图。解 MATLAB 程序为:%Example 5.2%k=100;z=-4;p=0 -0.5 -50 -50;num,den=zp2tf(z,p,k);bode(num,den);title(Bode Plot)执行后得下图:-150-100-50050Magnitude(dB)10-2 10-1 100 101 102 103-270-225-180-135-90Phase (deg)Bode PlotFrequency (rad/sec)例 14 有二阶系统:251()3ksG现要得到系统的 Nyquis
20、t 曲线,可输入:num=2 5 1;den=1 2 3;nyquist(num,den);title(Nyquist Plot)执行后可得到下图:-1 -0.5 0 0.5 1 1.5 2 2.5 3-2-1.5-1-0.500.511.52 Nyquist PlotReal AxisImaginaryAxis由于曲线没有包围-1+j0 点且 p=0,所以 G(s)单位负反馈构成的闭环系统稳定。例 15 开环系统: 50()(2)kGss解 MATLAB 程序为:%Example 5.3%k=50;z= ;p=-5 2;num,den=zp2tf(z,p,k);figure(1)nyquis
21、t(num,den)title(Nyquist Plot);figure(2)num1,den1=cloop(num,den);impulse(num1,den1);title(Impulse Response)执行后得下图:-5 -4.5 -4 -3.5 -3 -2.5 -2 -1.5 -1 -0.5 0-1.5-1-0.500.511.5 Nyquist PlotReal AxisImaginary Axis0 0.5 1 1.5 2 2.5 3 3.5 4-3-2-10123456 Impulse ResponseTime (sec)Amplitude从图中可以看出,系统 Nyquist
22、 曲线按逆时针方向包围(-1,j0)点 1 圈,而开环系统包含右半 S 平面上的 1 个极点,因此闭环系统稳定。例 16 开环系统:50()1()2kGss绘制系统的 Nyquist 曲线,判断闭环系统稳定性,绘制出闭环系统的单位冲激响应。解 MATLAB 程序为:%Example 5.4%k=50;z= ;p=-1 -5 2;num,den=zp2tf(z,p,k);figure(1)nyquist(num,den)title(Nyquist Plot);figure(2)num1,den1=cloop(num,den);impulse(num1,den1);title(Impulse Re
23、sponse)执行后得下图:-5 -4.5 -4 -3.5 -3 -2.5 -2 -1.5 -1 -0.5 0-2-1.5-1-0.500.511.52 Nyquist PlotReal AxisImaginary Axis0 0.1 0.2 0.3 0.4 0.500.511.522.533.54 Impulse ResponseTime (sec)Amplitude第 6 章 控制系统的综合与校正6.12 设单位反馈系统的开环传递函数为: )4(10)(ssG使用 Bode 设计法设计滞后超前校正装置,使校正后的系统能满足如下的性能指标: 在单位斜坡信号作用下,系统的速度误差系数 。10s
24、Kv 系统校正后的剪切频率 。1/5.sradWc 系统校正后相角稳定裕度 。40 校正后系统时域性能指标: 。sttp6,2%,3MATLAB 命令:k0=30;n1=1;d1=conv(conv(1 0,0.1 1),0.2 1);mag,phase,w=bode(k0*n1,d1);figure(1);margin(mag,phase,w);hold onfigure(2);s1=tf(k0*n1,d1);sys=feedback(s1,1);step(sys)执行后得到:(1)未校正系统的 BODE 图:-150-100-50050Magnitude (dB)10-1 100 101
25、102 103-270-225-180-135-90Phase (deg)Bode DiagramGm = -6.02 dB (at 7.07 rad/sec) , Pm = -17.2 deg (at 9.77 rad/sec)Frequency (rad/sec)(2)未校正系统的阶跃响应曲线:0 5 10 15-6-4-202468x 106 Step ResponseTime (sec)Amplitude%求滞后校正器的传递函数:MATLAB 命令:wc=1.5;k0=40;n1=1;d1=conv(conv(1 0,1 1),1 4);beta=9.5;T=1/(0.1*wc);be
26、tat=beta*T;Gc1=tf(T 1,betat 1)执行后所得结果:Transfer function:6.667 s + 163.33 s + 1%求超前校正器的传递函数:MATLAB 命令:n1=conv(0 40,6.667 1);d1=conv(conv(conv(1 0,1 1),1 4),63.33 1);sope=tf(n1,d1);wc=1.5;num=sope.num1;den=sope.den1;na=polyval(num,j*wc);da=polyval(den,j*wc);g=na/da;g1=abs(g);h=20*log10(g1);a=10(h/10);
27、wm=wc;T=1/(wm*(a)(1/2);alphat=a*T;Gc=tf(T 1,alphat 1)执行后所得结果:Transfer function:1.82 s + 10.2442 s + 1%校验:MATLAB 命令:n1=40;d1=conv(conv(1 0,1 1),1 4);s1=tf(n1,d1);s2=tf(6.667 1,63.33 1);s3=tf(1.82 1,0.2442 1);sope=s1*s2*s3;mag,phase,w=bode(sope);margin(mag,phase,w)执行后所得结果(校正后的系统 Bode 图)为:-100-50050100
28、Magnitude (dB)10-3 10-2 10-1 100 101 102-270-225-180-135-90Phase (deg)Bode DiagramGm = 14 dB (at 4.34 rad/sec) , Pm = 57.8 deg (at 1.5 rad/sec)Frequency (rad/sec)%校验后性能指标及阶跃响应:MATLAB 命令:global y t;k0=30;n1=40;d1=conv(conv(1 0,1 1),1,4);s1=tf(n1,d1);s2=tf(6.667 1,63.33 1);s3=tf(1.82 1,0.2442 1);sope=
29、s1*s2*s3;sys=feedback(sope,1);step(sys)y,t=step(sys);运行后,得到校正后的单位阶跃响应曲线如图示:0 5 10 15 20 25 30 3500.20.40.60.811.21.4 Step ResponseTime (sec)Amplitude第 7 章 离散控制系统一、例题例 21 某二阶系统:23.415()608zG要求其阶跃响应,可输入:num=2 -3.4 1.5;den=1 -1.6 0.8;dstep(num,den)title(Discrete Step Response)执行后可得下图:0 10 20 30 40 50 6
30、0-0.500.511.52 Discrete Step ResponseTime (sec)Amplitude例 23 有系统:23.415()608zG可输入:num=2 -3.4 1.5;den=1 -1.6 0.8;axis(square)zgrid(new)rlocus(num,den);title(Root Locus)执行后可得下图:-1 -0.8 -0.6 -0.4 -0.2 0 0.2 0.4 0.6 0.8 1-1-0.8-0.6-0.4-0.200.20.40.60.810.1/T0.2/T0.3/T0.4/T0.5/T0.6/T0.7/T0.8/T0.9/T/T0.1/
31、T0.2/T0.3/T0.4/T0.5/T0.6/T0.7/T0.8/T0.9/T/T0.10.20.30.40.50.60.70.80.9Root LocusReal AxisImaginaryAxis例 24 已知离散系统: ,绘制出系统的20.69()17583HzzNyquist 曲线,判别闭环系统的稳定性,并绘制出闭环系统的单位冲激响应。解 MATLAB 程序如下:num=0.692;den=1,-1.758,0.375;z,p,k=tf2zp(num,den);pfigure(1)dnyquist(num,den,0.1)title(离散 Nyquist 曲线图);xlabel(实
32、数轴 );ylabel(虚数轴);figure(2)num1,den1=cloop(num,den);dimpulse(num1,den1);title(离散冲激响应);xlabel(时间);ylabel( 振幅);运行程序可得下图:-2 -1.5 -1 -0.5 0 0.5-0.4-0.3-0.2-0.100.10.20.30.4 冲冲Nyquist冲冲冲冲冲冲冲冲冲0 50 100 150 200 250 300 350 400 450-2-1.5-1-0.500.511.522.5x 106 冲冲冲冲冲冲冲冲 (sec)冲冲p =1.50960.2484由仿真图可知该离散系统是发散的。二
33、、仿真例 1.已知一个离散系统如图所示,其中采样周期 TS=1s.,对象模型,零阶保持器 ,试求开环增益的稳定范围。)1()skGp sesT1)(G0图 7.1 系统模型图解:执行如下程序:num=1;den=1 1 0;sys=tf(num,den); %连续系统传递函数c2d (sys,1) %离散系统传递函数运行结果如下:Transfer function:0.3679 z + 0.2642-z2 - 1.368 z + 0.3679Sampling time: 1继续编写程序,求取该离散系统的根轨迹图:num=0.3679 0.2642;den=1 -1.368 0.3679;G=t
34、f(num,den,-1);rlocus(G)k,poles=rlocfind(G)运行结果如下:用鼠标单击根轨迹与单位圆的交点,可以得到交点的极点坐标以及交点处的开环增益 K 值,如图 7.2 所示。运行程序如下:selected_point =1.0118 - 0.0047i selected_point = -0.7227 + 0.0047ik =0.0127 k =778.5564poles = poles =0.9873 -284.33820.3761 -0.7247selected_point = 0.2536 - 0.9643ik = 2.3513poles =0.2515 +
35、0.9622i0.2515 - 0.9622i-2.5 -2 -1.5 -1 -0.5 0 0.5 1 1.5-1.5-1-0.500.511.5 Root LocusReal AxisImaginary Axis图 7.2 离散系统的根轨迹图分析:在程序运行中选取了三个点,极点、零点以及根轨迹与单位园的交点,可得在整个根轨迹图中,k=0778.5564 。在离散系统根轨迹图上,虚线表示的是单位元,由理论分析可知,系统闭环传递函数的所有极点位于 z 平面的单位圆内时,该离散系统是稳定的。从根轨迹的分布图上可以看出,当0k2.3513 时,该离散系统是稳定的。令 k=1,编写程序,还可以求得该离
36、散系统的单位冲激响应曲线,进一步验证该离散系统是否稳定。程序如下:num=0.3679 0.2642;den=1 -1.368 0.3679;num1,den1=cloop(num,den);dimpulse(num1,den1)运行程序,得到离散闭环系统的单位冲激响应如图 7.3,可见该离散系统在k=1 时是稳定的。0 5 10 15 20 25-0.3-0.2-0.100.10.20.30.40.50.60.7 Impulse ResponseTime (sec)Amplitude第 8 章 控制系统状态空间分析与综合例 19 已知线性定常系统如图所示:试求系统的的状态方程,选择正定的实对
37、称矩阵 Q 后计算李雅普诺夫方程的解并利用李雅普诺夫函数确定系统的稳定性。解 选择正定的实对称矩阵10Q编写如下程序:n1=5;d1=1 1;s1=tf(n1,d1);n2=1;d2=1 2;s2=tf(n2,d2);n3=1;d3=1 0;s3=tf(n3,d3);s123=s1*s2*s3;sys=feedback(s123,1);A B C D=tf2ss(sys.num1,sys,den1);A B C D=tf2ss(sys.num1,sys.den1);q=1 0 0;0 1 0;0 0 1;if det(A)=0p=lyap(A,q)det1=det(p(1,1)det2=det
38、(p(1:2,1:2)detp=det(p)end程序运行结果为:p =23.0000 -0.5000 -13.5000-0.5000 13.5000 -0.5000-13.5000 -0.5000 8.2000det1 =23.0000det2 =310.2500detp =71.1750因为10Q故有: 是负定的。22()(13)TVxx从运行结果可以看出:对各阶主子行列式(det1,det2,det3)进行校验说明 p阵确实是正定阵,因此系统在坐标原点的平衡状态时稳定的,而且是大范围渐进稳定的。第 9 章 非线性控制系统例 19 试用描述函数法分析下图非线性系统的稳定性:解 程序如下:k
39、=input(k=)num=0 0 0 k;den=1 6 5 0;w=0.1:0.1:100;re,im,w=nyquist(num,den,w);v=-4 4,-5 5;axis(v);plot(re,im);title(Curves of - 1/N(X) and G(jw);xlabel(Re);ylabel(I);grid on;hold on;x=0.1:0.1:50;z=-(pi/4)*x;plot(real(z),imag(z)仿真结果为:-0.4 -0.2 0 0.2 0.4 0.6-0.25-0.2-0.15-0.1-0.0500.050.1 Curves of - 1/N
40、(X) and G(jw)ReI图 9.1 K=5 时-0.15 -0.1 -0.05 0 0.05 0.1-0.015-0.01-0.00500.0050.010.015 Curves of - 1/N(X) and G(jw)ReI图 9.2 K=2 时仿真结果分析:由图 9.1 和图 9.2 可以看出,当曲线跟曲线相交时,非线性系统是不稳定的;当曲线和曲线没有相交点时,非线性系统是稳定的。MATLAB 在电力系统中的应用实例:短路故障分析(1)单相接地短路在 0-0.03s 时线路工作在稳定状态,三相电流、电压对称。在 0.03s 时发生 A 相接地短路,A 相电压基本为 0,B、C 相
41、电压也相对减小;故障相 A 相电流迅速上升为短路电流,B、C 相电流也相对增大;三相电压、电流不再对称,说明单相接地短路为不对称短路。在 0.08s 时,故障切除,三相电压电流经暂态后达到新的稳定状态,并且重新恢复三相对称运行的工作状态(2)两相接地短路在 0-0.03s 时线路工作在稳定状态,三相电流、电压对称。在 0.03s 时发生 B、C 两相接地短路,B、C 两相电压基本为 0,A 相电压略有增大;故障相B、C 两相电流迅速上升为短路电流,C 相电流相对增大;三相电压、电流不再对称,说明两相短路为不对称短路。在 0.08s 时,故障切除,三相电压电流经暂态后达到新的稳定状态,并且重新恢复三相对称运行的工作状态。(3)三相短路