控制工程基础经典控制部分MATLAB分析.docx

1、控制工程基础经典控制部分MATLAB分析控制工程基础（经典控制部分） 的 MATLAB 分析第一章 MATLAB 的基本使用 1-1 MATLAB 语言简介MATLAB 是一种高级矩阵语言，它由 Math Works 公司于 1984 年正式推出， 它的基本处理对象是矩阵， 即使是一个标量纯数， MATLAB 也认为它是只有一个 元素的矩阵。随着 MA TLAB 的发展，特别是它所包含的大量工具箱（应用程序 集）的集结，使 MATLAB 已经成为带有独特数据结构、输入输出、流程控制语 句和函数、并且面向对象的高级语言。MATLAB 语言被称为一种“演算纸式的科学计算语言” ，它在数值计算、符

2、号运算、数据处理、自动控制、信号处理、神经网络、优化计算、模糊逻辑、系 统辨识、小波分析、图象处理、统计分析、甚至于金融财会等广大领域有着十分 广泛的用途。MATLAB 语言在工程计算与分析方面具有无可比拟的优异性能。它集计算、 数据可视化和程序设计于一体，并能将问题和解决方案以使用者所熟悉的数学符 号或图形表示出来。MATLAB 语言和 C 语言的关系与 C 语言和汇编语言的关系类似。例如当我 们需要求一个矩阵的特征值时 ,在 MATLAB 下只需由几个字符组成的一条指令即 可得出结果，而不必去考虑用什么算法以及如何实现这些算法等低级问题，也不 必深入了解相应算法的具体内容。就象在 C 语言

3、下不必象汇编语言中去探究乘法是怎样实现的，而只需要采用乘积的结果就可以了。MATLAB 语言还有一个巨大的优点是其高度的可靠性。 例如对于一个病态矩 阵的处理， MATLAB 不会得出错误的结果， 而用 C 或其它高级语言编写出来的程 序可能会得出错误的结果。这是因为 MATLAB 函数集及其工具箱都是由一些在 该领域卓有研究成果，造诣很深的权威学者经过反复比较所得出来的最优方法， 而且经过多年的实践检验被证明是正确可靠的。 1-2 MATLAB 的工作窗口下面以 MATLAB6.1 为例介绍。从实用的角度 MATLAB 的工作窗口包括命令窗口、 M 文件编辑器窗口、图形编 辑窗口、数学函数庫

4、、应用程序接口及在线窗口。下面首先介绍 MATLAB 的命 令窗口及 M 文件编辑器。、命令窗口启动MATLAB之后，屏幕上自动出现命令窗口 MATLAB ,它是MATLAB提供给 用户的操作界面，用户可以在命令窗口内提示符“ ”之后（有的 MATLAB版本命令窗口没有提示符）键入 MATLAB命令，回车即获得该命令的答案。命令窗口内有 File、Edit、View、Web、Window、Help 等菜单条。、M文件编辑窗口M文件是MATLAB语言所特有的文件。用户可以在 M文件编辑窗口内，编写一段程序，调试，运行并存盘，所保存的用户程序即是用户自己的 M文件。MATLAB工具箱中大量的应用程

5、序也是以 M文件的形式出现的，这些 M文件可以打开来阅读，甚至修改，但应注意，不可改动工具箱中的 M文件！1 .进入M文件窗口有两种方法1） 命令窗口 一 File New M-File ;2） 命令窗口 一点击“ File”字样下面的 口图标。M文件编辑窗口的标记是“ Un titled ”（无标题的）。当用户编写的程序要存盘时，Un titled作为默认文件名提供给用户，自然，用户可以，也应当自己命名。若用户 不自己命名，贝U MATLAB会对Un titled进行编号。2.M文件的执行：返回命令窗口，在当前目录（Current Directory ）内选择所要运行的 M文件的目录, 在命

6、令窗口提示符“ ”后，直接键入文件名（不加后缀）即可运行。注意：（1）机器默认路径为一级子目录 MATLAB6 p1work ;（2） MATLAB 6.1以前的版本，运行 M文件的方法稍有不同，它必 须在File菜单下，打开“ Run Script ”子菜单，键入需要运行的文件路径及名称 再回车，在这种情况下， work作为根目录对待，不出现在 M文件的路径之中。本讲义的参考程序都是在 M文件窗口下编制的。三、在线帮助窗口在命令窗口中键入 Help（ 空格 ） 函数名，可以立即获得该函数的使用方法。 1-3 MATLAB 最基本的矩阵操作 作为命令窗口及 M 文件编辑器的应用实例，介绍几个最

7、基本的矩阵运算命令。一、矩阵的输入在方括号内依次按行键入矩阵元素，在一行内的各元素之间用空格或逗号分开， 每行之间用分号分开。例如，在命令窗内输入A=2 2 3;4 5 4;7 8 9 （ 注意：方括号，分号为矩阵行标记 ）B=1,3,5;6,-4,2;3,5,1（ 逗号与空格功能相同)A = 223B = 1 3 54546 -4278935 1同理：输入A1 =246 得到行矢量，输入A2 =2；4；6 得到列矢量，于是，当输入C=A；A1 有C= 123456789246A1 作为矩阵 C 的最后一行， C 和 A 相比，增加了一行。二、矩阵的转置矩阵A的转置用A 表示，显然，A1与A2

8、互为转置，即 A1 /会得到以2,4 , 6 为元素的列矢量。思考一下输入 C1 =： A A2 :C2 = A A1 /有什么结果？而输入 A ；A1 有无意义？三、矩阵的四则运算1矩阵的加减法：当两个矩阵维数相同时可以直接进行“”或“”运算。如 D1 = A + B , D2 = A B2 矩阵的乘法：当矩阵 A， B 维数相容时C3 = A * B :普通意义下的矩阵相乘C4 = A . * B :矩阵A与B的对应元素相乘显然，A * B工B * A (一般情况)，而A . * B = B . * A。A . * B称为数列型乘法， 它要求参加运算的矩阵或数列具有相同的行列数，这是 MA

9、 TLAB 语言中的一种 特殊运算，它在今后求取函数值等运算时是很重要的。实际上，前面所述的矩阵 加、减法就是一种数列型运算。3 矩阵的除法D4 = AB :表示A-1 * B或inv(A)*B，即A的逆矩阵左乘矩阵 B。D5 = B/A :表示B * A-1或B * inv(A)，即A的逆矩阵右乘 B。D6 = A .B :表示B的每一个元素被 A的对应元素除。D7 = A ./B :表示A的每一个元素被 B的对应元素除。显然， A .B与A ./B的各对应元素互为倒数。读者可以思考一下， D6 . * D7 等于什么？D8 = inv(A) : A的逆矩阵。打开 M 文件编辑窗口，将上述命

10、令依次键入，得到 fanli001 如下：参考程序 fanli001: 矩阵的四则运算A=2 2 3;4 5 4;7 8 9 % 三阶矩阵输入B=1,3,5;6,-4,2;3,5,1 % 三阶矩阵输入A1=2 4 6 % 行向量 A2=2;4;6 % 列向量 C=A;A1 % 矩阵 A 增加一行 C1=A A2 % 矩阵 A 增加一列 C2=A A1 % 矩阵 A 增加一列 D1=A+B % 矩阵相加 D2=A-B % 矩阵相减 C3=A*B % 矩阵与矩阵相乘 C4=A.*B % 矩阵的对应元素相乘 D3=AB %A 的逆左乘 B D4=B/A %A 的逆右乘 B D6=A.B %B 的各元

11、素被 A 的对应元素除 D7=A./B %A 的各元素被 B 的对应元素除 D8=inv(A) %A 的逆矩阵 语句后面的 %为语句说明符。MATLAB 中矩阵运算的其它主要命令可通过在线帮助获得。 1-4 MATLAB 的符号运算操作一、进入符号运算功能 在命令窗口键入 syms x y z t此后，即可以使用 x，y，z，t 等作自变量定义函数。 syms x y z t real 规定所定义的变量为实型。二、 代数方程求解使用命令solve可以求解代数方程，如求下例方程组厂 2x +xy + y -3 = 02/ 4x2y十3 = 0的解，命令为x,y=solve( xA2+x*y-+y

12、0,xA2-4*x- 2*y+3=0 )程序见范例程序 fanli002。参考程序fanli002:代数方程求解syms x y %进入符号运算功能f1=xA2+x*y+y-3 % 函数 flf2=xA2-4*x-2*y+3 % 函数 f2x,y=solve(xA2+x*y+y-3=0,xA2-4*x-2*y+3=0)%求解方程组f1a=simplify(subs(f1) % 用求解出的x,y检验方程1f2a=simplify(subs(f2) % 用求解出的x,y检验方程2x=double(x) %将符号变量转换成浮点数y=double(y) %将符号变量转换成浮点数f1=subs(f1)

13、%用浮点数x,y检验方程1f2=subs(f2) %用浮点数x,y检验方程2由solve求出的根是根的符号表达形式，是准确解。命令simplify(f)表示化简，subs(f1)表示将求解出x，y代回fl中；double是将符号变量转换成浮点数，是准确解 x，y的近似值。这从程序运行 f1a=0 , f2a=0 , fl工0, f2工0可以确认。三、 符号矩阵运算符号矩阵可以和数值矩阵一样进行运算， 例如：求矩阵特征值eig，求矩阵的逆inv等命令都支持符号运算。_si nt costA 二设 cost sint计算其特征值 eig(A)，逆矩阵inv(A)，程序见fanli003 。参考程序

14、fanli003:符号矩阵的特征值syms t real %定义为实型变量A=sin(t) -cos(t);cos(t) sin(t) % 定义矩阵 AB1=eig(A) % 求矩阵A的特征值B仁simple(BI) %化简A的特征值表达式B2=inv(A) % 求矩阵A的逆矩阵B2=simple(B2) %化简A的逆矩阵表达式C1=A*B2 %检验A的逆矩阵C1=simple(C1) %C1 为单位矩阵 注意函数的输入方法，自变量用圆括号括起来。四、微积分运算设函数f二f(X)，则MATLAB中微积分运算命令为fiff(f):求函数f对自变量x的一阶导数；diff(f，2):求函数f对自变量

15、x的二阶导数；int(f):求函数f的不定积分例1.1，设f -fanli004。5+4 cos x试计算其一阶，二阶导数，积分运算，并作出函数图象，见范例参考程序fanli004:函数的微分与积分syms x f=1/(5+4*cos(x)ezplot(f) %函数f的曲线f1=diff(f) %函数f的一阶导数 figure,ezplot(f1) % 函数f一阶导数的曲线 f2=diff(f,2) % 函数f的二阶导数 figure,ezplot(f2) % 函数f二阶导数的曲线 g=int(int(f2) %函数f的二阶导数f2的二重积分 figure,ezplot(g) % 函数f2二

16、重积分的曲线 e=f-g %二阶导数的二重积分与原函数的差 e=simple(e)figure,ezplot(e)程序中ezplot(f):作函数f = f(X)的图形，x的取值范围默认值为-2 n x 2 nfanli004的函数曲线24/(5+4 cos(x) sin(x)xfanli004的一阶导函数曲线fanli004的二阶导函数曲线xfanliOO4的二阶导函数的二重积分曲线ezplot(f)是一个很有用的作图命令，它的其它应用形式，请查在线帮助。细心的读者会发现，一个函数求二阶导数后再对二阶导数进行二重 积分，其结果与原函数相差一个常数。相当于纵坐标发生平移。第二章 系统的时域特性

17、 2-1 传递函数、传递函数的两种形式传递函数通常表达成 s的有理分式形式及零极点增益形式。i m-iS设传递函数G(s)二i nan _j Sj=0G(s)二Gi(s)二num den心【(s Zi)i 2ni【(s Pj)j2s 1s3 2s2 2s 11.有理分式形式0补齐。分别将分子、分母中s多项式的系数按降幕排列成行矢量， 缺项的系数用上述函数可表示为num仁2 1 % (注意：方括号，同一行的各元素间留空格或逗号) 。den 1=1 2 2 1syss1=tf( nu m1,de n1)运行后，返回传递函数 G1 (s)的形式。这种形式不能直接进行符号运算!2零极点增益形式Z,P,

18、K = tf2zp(num1,den1)sys2 = zpk(Z,P ,K)返回零、极点、增益表达式，其 Z， P 分别将零点和极点表示成列向量，若无零点或极点用 (空矩阵)代替。运行得到 G1(s) 的零点 Z = -0.5极点 P= -1，-0.5 j0.866增益 K = 2指令zp2tf(Z , P , K)将零极点增益变换成有理分式形式，见程序 fanli005 。参考程序 fanli005: 传递函数的有理分式及零极点增益模型 num1=2 1 % 传递函数的分子系数向量 den1=1 2 2 1 % 传递函数的分母系数向量 sys1=tf(num1,den1) % 传递函数的有理

19、分式模型 Z,P,K=tf2zp(num1,den1)% 有理分式模型转换成零极点增益模型 num2,den2=zp2tf(Z,P ,K)% 零极点增益模型转换成有理分式模型 sys2=zpk(Z,P ,K) % 传递函数的零极点增益模型 A1,B1,C1,D1=tf2ss(num1,den1)% 有理分式模型转换成状态空间模型 A2,B2,C2,D2=zp2ss(Z,P ,K)% 零极点及增益模型转换成状态空间模型 num1,den1=ss2tf(A1,B1,C1,D1)% 状态空间模型转换成有理分式模型 Z,P,K=ss2zp(A2,B2,C2,D2)% 状态空间模型转换成零极点增益模型程

20、序中，命令tf2ss , zp2ss及ss2tf , ss2zp是状态空间模型与有理分式及零、极 点、增益模型之间的相互转换。、传递函数框图的处理用框图可以方便地表示传递函数的并联，串联及反馈。为简洁，仅以有理分式模型为例。1.并联sysp = parallel(sys1,sys2)nu m,de n=parallel( nu m1,de n1, nu m2,de n2)2.串联syss = series(sys1,sys2)nu ms,de ns = series( nu m1,de n1,nu m2,de n2)3.反馈G1(S)G2(S)1+G1(s)G2(s) G3(s)sysc=fe

21、edback(syss,sys3, 1) %默认值(-1) nu me, dene = feedback (nu ms, dens, nu m3, den3)Gl(S)G2(S)1+Gl(S)G2(S)sysd = feedback(syss, 1)numd, dend = feedback(nums, dens, 1, 1) %(单位反馈)上面给出了同一指令的两种形式，相当于两套平行指令。对于零极点增益形式， 书写稍复杂一些，可先用 zpk转换成系统形式，或用 zp2tf转折换成有理分式形式后再进行框图化简操作。三、简单函数的拉普拉斯变换在MATLAB的符号功能中，可以对简单函数进行拉普拉斯

22、正、逆变换。 拉氏正变换：laplace(f(t)拉氏逆变换：ilaplace(L(s)其中f (t)为原函数，L(s)为象函数。命令格式参见 fanli007。参考程序fanli007 :拉普拉斯变换syms s t w a b cf1=sqrt(b-a)A2+wA2)/w*exp(-a*t)*sin(w*t+atan(w/(b-a)% 原函数 flL1=laplace(f1)%f1的拉氏变换(象函数)L1=simple(L1)%化简f2=il aplace(L1)%L1的拉氏逆变换f2=simple(f2) % 化简在MATLAB中使用laplace及ilaplace命令时，要注意象、原函

23、数的符号，特别 是对初相不等于零的振荡系统，运行结果常常同手册上的结果相差一个符号，这 要注意函数表达式成立的条件。保险的办法是再使用拉氏变换的初值定理确定象、 原函数的符号。2-2 系统时域特性曲线在 MATLAB 中，当传递函数已知时，可以方便地求出系统的单位脉冲响应、单 位阶跃响应等曲线。一、系统的单位阶跃响应 stepstep 有以下几种格式step(sys) ：直接作出 sys 的单位阶跃响应曲线。 其中 sys = tf(num, den) 或 sys = zpk(z, p, k) ， MATLAB 自动决定响应时间。step(sys, t)设定响应时间的单位阶跃响应。 t 可以设

24、定为最大响应时间 t =t 终值 (秒 ) ，也可以设置为一个向量t = 0 : t: t 终值 注意冒号的使用。它产生一个从 0 到 t 终值的行矢量，元素之间的间隔为 t。step(sys1, sys2, ,sysn)在同一幅图上画出几个系统的单位阶跃响应。y, t = step(sys) ；命令输出对应时刻 t 的各个单位阶跃响应值，不画图。语句后的分号控制数据的 屏幕显示。如果要查看机器计算了多少个数据，可以使用命令size(y)得出的结果也表明数据作为列矢量的行数。要将计算出的 y, t 作成曲线，使用一 般的作图命令plot(t, y)plot 后面跟的两个参数横坐标在前，纵坐标在

25、后。参考程序见 fanli008 ：参考程序 fanli008: 系统的单位阶跃响应num1=4 2den1=2 8 14 11 4sys1=tf(num1,den1) % 系统 G1(s) num2=2 1den 2=1 4 6 7 3sys2=tf(num2,den2) % 系统 G2(s)y1,t1=step(sys1); % 系统G1(s)的单位阶跃响应数据y2,t2=step(sys2); % 系统G2(s)的单位阶跃响应数据 step(sys1,sys2) % 系统G1(s)、G2(s)的单位阶跃响曲线 figure,step(sys1,sys2,20)%系统G1(s)、G2(s)

26、在自选时间(20秒)内的单位阶跃响曲线figure,plot(t1,y1) % 系统G1(s)的单位阶跃响应曲线figure,plot(t2,y2) % 系统G2(s)的单位阶跃响应曲线fan li008:step(sys1,sys2,t)单位阶跃曲线二、系统的单位脉冲响应impulse命令格式与单位阶跃响应 step的命令格式完全相同，只需将语句中的step用impulse 代替即可。针对同样的系统，其单位脉冲响应的参考程序见fanli009 。参考程序fanli009:系统的单位脉冲响应num仁4 2den仁2 8 14 11 4sys仁tf( nu m1,de n1)num2=2 1de

27、n 2=1 4 6 7 3sys2=tf( num2,de n2)y1,t1=impulse(sys1); % 系统G1(s)的单位脉冲响应数据y2,t2=impulse(sys2); % 系统G2(s)的单位脉冲响应数据 impulse(sys1,sys2) % 系统G1(s)、G2(s)的单位脉冲响应曲线 figure,impulse(sys1,sys2,20)%系统G1(s)、G2(s)在自选时间(20秒)内的单位脉冲响应曲线figure,plot(t1,y1) % 系统G1(s)的单位脉冲响应曲线figure,plot(t2,y2) % 系统G2(s)的单位脉冲响应曲线hold on,

28、step(sys2) % 系统G2(s)的单位阶跃和单位脉冲响应曲线Time (sec)fanli009:impulse(sys1,sys2)单位脉冲响应曲线程序的最后一句hold on是当前图形保护模式。当要将新图形作在当前图形 上时，必须使用 hold on。而figure的含意是另开一个新的图形窗口，如果不用 figure或hold on，则新的图形会占用原图形窗口，始终只保留一个最新的图形窗 口。三、一阶系统及二阶系统的时域特性一阶系统及二阶系统是最基本也是最重要的系统，高阶系统总可以视为由若干个一阶和(或)二阶系统组合构成。1. 一阶系统(设增益为1)1G(s厂Ts+1影响系统特性的

29、参数是其时间常数 T, T越大，系统惯性越大，响应越慢。参考程序fanli010给出了 T= 0.4, 1.2, 2.0, 2.8, 3.6, 4,4六条单位阶跃响应曲线。参考程序fanli010: 阶系统的单位阶跃响应曲线num=1;i=1;for del=0.1:0.2:1.1 % 一阶系统时间常数递增间隔den=4*del 1; % 一阶系统分母向量step(tf( nu m,de n) % 一阶系统单位阶跃响应曲线hold on, %不同时间常数的一阶系统单位阶跃响应曲线簇i=i+1;end同理，可以作出对应的单位脉冲响应曲线，参考程序 fanl i011。参考程序fanli011:

30、阶系统的单位脉冲响应曲线num=1;i=1;Step Resp onseedur-rpmA2025fan li010 阶系统时间常数对单位阶跃响应的影响fanl i011 阶系统时间常数对单位脉冲响应的影响for del=0.1:0.2:1.1den= 4*del 1;impulse(tf( num,de n),10) % 一阶系统单位阶脉冲应曲线hold on, %不同时间常数的一阶系统单位脉冲响应曲线簇i=i+1;end注意MATLAB中for语句的结构。读者可以改变不同的增益,看看图形有何变化。2.二阶系统(设0 E 1)设二阶系统为G(s)二阶系统的特征参数为固有频率 n及阻尼比E。当 n增大，