控制工程基础经典控制部分MATLAB分析.docx
《控制工程基础经典控制部分MATLAB分析.docx》由会员分享，可在线阅读，更多相关《控制工程基础经典控制部分MATLAB分析.docx（95页珍藏版）》请在冰点文库上搜索。
控制工程基础经典控制部分MATLAB分析
控制工程基础（经典控制部分）》
的MATLAB分析
第一章MATLAB的基本使用
§1-1MATLAB语言简介
MATLAB是一种高级矩阵语言，它由MathWorks公司于1984年正式推出，它的基本处理对象是矩阵，即使是一个标量纯数，MATLAB也认为它是只有一个元素的矩阵。
随着MATLAB的发展，特别是它所包含的大量工具箱（应用程序集）的集结，使MATLAB已经成为带有独特数据结构、输入输出、流程控制语句和函数、并且面向对象的高级语言。
MATLAB语言被称为一种“演算纸式的科学计算语言”，它在数值计算、符号运算、数据处理、自动控制、信号处理、神经网络、优化计算、模糊逻辑、系统辨识、小波分析、图象处理、统计分析、甚至于金融财会等广大领域有着十分广泛的用途。
MATLAB语言在工程计算与分析方面具有无可比拟的优异性能。
它集计算、数据可视化和程序设计于一体，并能将问题和解决方案以使用者所熟悉的数学符号或图形表示出来。
MATLAB语言和C语言的关系与C语言和汇编语言的关系类似。
例如当我们需要求一个矩阵的特征值时,在MATLAB下只需由几个字符组成的一条指令即可得出结果，而不必去考虑用什么算法以及如何实现这些算法等低级问题，也不必深入了解相应算法的具体内容。
就象在C语言下不必象汇编语言中去探究乘法
是怎样实现的，而只需要采用乘积的结果就可以了。
MATLAB语言还有一个巨大的优点是其高度的可靠性。
例如对于一个病态矩阵的处理，MATLAB不会得出错误的结果，而用C或其它高级语言编写出来的程序可能会得出错误的结果。
这是因为MATLAB函数集及其工具箱都是由一些在该领域卓有研究成果，造诣很深的权威学者经过反复比较所得出来的最优方法，而且经过多年的实践检验被证明是正确可靠的。
§1-2MATLAB的工作窗口
下面以MATLAB6.1为例介绍。
从实用的角度MATLAB的工作窗口包括命令窗口、M文件编辑器窗口、图形编辑窗口、数学函数庫、应用程序接口及在线窗口。
下面首先介绍MATLAB的命令窗口及M文件编辑器。
、命令窗口
启动MATLAB之后，屏幕上自动出现命令窗口MATLAB,它是MATLAB提供给用户的操作界面，用户可以在命令窗口内提示符“>>”之后（有的MATLAB版
本命令窗口没有提示符）键入MATLAB命令，回车即获得该命令的答案。
命令窗口内有File、Edit、View、Web、Window、Help等菜单条。
、M文件编辑窗口
M文件是MATLAB语言所特有的文件。
用户可以在M文件编辑窗口内，编写一
段程序，调试，运行并存盘，所保存的用户程序即是用户自己的M文件。
MATLAB
工具箱中大量的应用程序也是以M文件的形式出现的，这些M文件可以打开来
阅读，甚至修改，但应注意，不可改动工具箱中的M文件！
1.进入M文件窗口有两种方法
1）命令窗口一File—New—M-File;
2）命令窗口一点击“File”字样下面的口图标。
M文件编辑窗口的标记是“Untitled”（无标题的）。
当用户编写的程序要存盘时，
Untitled作为默认文件名提供给用户，自然，用户可以，也应当自己命名。
若用户不自己命名，贝UMATLAB会对Untitled进行编号。
2.M文件的执行：
返回命令窗口，在当前目录（CurrentDirectory）内选择所要运行的M文件的目录,在命令窗口提示符“>>”后，直接键入文件名（不加后缀）即可运行。
注意：
（1）机器默认路径为一级子目录MATLAB6p1\work;
（2）MATLAB6.1以前的版本，运行M文件的方法稍有不同，它必须在File菜单下，打开“RunScript••”子菜单，键入需要运行的文件路径及名称再回车，在这种情况下，work作为根目录对待，不出现在M文件的路径之中。
本讲义的参考程序都是在M文件窗口下编制的。
三、在线帮助窗口
在命令窗口中键入Help（空格）函数名，可以立即获得该函数的使用方法。
§1-3MATLAB最基本的矩阵操作作为命令窗口及M文件编辑器的应用实例，介绍几个最基本的矩阵运算命令。
一、矩阵的输入
在方括号内依次按行键入矩阵元素，在一行内的各元素之间用空格或逗号分开，每行之间用分号分开。
例如，在命令窗内输入
A=[223;454;789]（注意：
方括号，分号为矩阵行标记）
B=[1,3,5;6,-4,2;3,5,1]
（逗号与空格功能相同
）
A=2
2
3
B=135
4
5
4
6-4
2
7
8
9
3
51
同理：
输入
A1=
2
4
6]得到行矢量，
输入
A2=
2；
4；
6]得到列矢量，
于是，当输入
C=[A；
A1]
有
C=1
2
3
4
5
6
7
8
9
2
4
6
A1作为矩阵C的最后一行，C和A相比，增加了一行。
二、矩阵的转置
矩阵A的转置用A'表示，显然，A1与A2互为转置，即A1/会得到以2,4,6为元素的列矢量。
思考一下输入C1=：
AA2:
C2=[AA1/]
有什么结果？
而输入[A；A1＇]有无意义？
三、矩阵的四则运算
1．矩阵的加减法：
当两个矩阵维数相同时可以直接进行“＋”或“－”运算。
如D1=A+B,D2=A—B
2．矩阵的乘法：
当矩阵A，B维数相容时
C3=A*B:
普通意义下的矩阵相乘
C4=A.*B:
矩阵A与B的对应元素相乘
显然，A*B工B*A（一般情况），而A.*B=B.*A。
A.*B称为数列型乘法，它要求参加运算的矩阵或数列具有相同的行列数，这是MATLAB语言中的一种特殊运算，它在今后求取函数值等运算时是很重要的。
实际上，前面所述的矩阵加、减法就是一种数列型运算。
3．矩阵的除法
D4=A\B:
表示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文件编辑窗口，将上述命令依次键入，得到fanli001如下：
参考程序fanli001:
矩阵的四则运算
A=[223;454;789]%三阶矩阵输入
B=[1,3,5;6,-4,2;3,5,1]%三阶矩阵输入
A1=[246]%行向量A2=[2;4;6]%列向量C=[A;A1]%矩阵A增加一行C1=[AA2]%矩阵A增加一列C2=[AA1']%矩阵A增加一列D1=[A+B]%矩阵相加D2=A-B%矩阵相减C3=A*B%矩阵与矩阵相乘C4=A.*B%矩阵的对应元素相乘D3=A\B%A的逆左乘BD4=B/A%A的逆右乘BD6=A.\B%B的各元素被A的对应元素除D7=A./B%A的各元素被B的对应元素除D8=inv（A）%A的逆矩阵语句后面的%为语句说明符。
MATLAB中矩阵运算的其它主要命令可通过在线帮助获得。
§1-4MATLAB的符号运算操作
一、进入符号运算功能在命令窗口键入symsxyzt
此后，即可以使用x，y，z，t等作自变量定义函数。
symsxyztreal规定所定义的变量为实型。
二、代数方程求解
使用命令solve可以求解代数方程，如求下例方程组
厂2
x+xy+y-3=0
<
2
/—4x—2y十3=0的解，命令为
[x,y]=solve（‘xA2+x*y-+y0,xA2-4*x-2*y+3=0'）
程序见范例程序fanli002。
参考程序fanli002:
代数方程求解
symsxy%进入符号运算功能
f1=xA2+x*y+y-3%函数fl
f2=xA2-4*x-2*y+3%函数f2
[x,y]=solve（'xA2+x*y+y-3=0','xA2-4*x-2*y+3=0'）
%求解方程组
f1a=simplify（subs（f1））%用求解出的x,y检验方程1
f2a=simplify（subs（f2））%用求解出的x,y检验方程2
x=double（x）%将符号变量转换成浮点数
y=double（y）%将符号变量转换成浮点数
f1=subs（f1）%用浮点数x,y检验方程1
f2=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
等命令都支持符号运算。
_sint—cost]
A二
设[costsint
计算其特征值eig（A），逆矩阵inv（A），程序见fanli003。
参考程序fanli003:
符号矩阵的特征值
symstreal%定义为实型变量
A=[sin（t）-cos（t）;cos（t）sin（t）]%定义矩阵A
B1=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对自变量x的二阶导数；
int（f）:
求函数f的不定积分
例1.1，设
f-
fanli004。
5+4cosx
试计算其一阶，二阶导数，积分运算，并作出函数图象，见范例
参考程序fanli004:
函数的微分与积分
symsxf=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二重积分的曲线e=f-g%二阶导数的二重积分与原函数的差e=simple（e）
figure,ezplot（e）
程序中
ezplot（f）:
作函数f=f（X）的图形，x的取值范围默认值为-2nfanli004的函数曲线
2
4/（5+4cos（x））sin（x）
x
fanli004的一阶导函数曲线
fanli004的二阶导函数曲线
x
fanliOO4的二阶导函数的二重积分曲线
ezplot（f）是一个很有用的作图命令，它的其它应用形式，请查在线帮助。
细心的读者会发现，一个函数求二阶导数后再对二阶导数进行二重积分，其结果与原函数相差一个常数。
相当于纵坐标发生平移。
第二章系统的时域特性
§2-1传递函数
、传递函数的两种形式
传递函数通常表达成s的有理分式形式及零极点增益形式。
im-iS
设传递函数
G（s）二
i£
n
―an_jS
j=0
G（s）二
Gi（s）二
numden
心【（sZi）
i2
n
i【（sPj）
j
2s1
s32s22s1
1.有理分式形式
0补齐。
分别将分子、分母中s多项式的系数按降幕排列成行矢量，缺项的系数用
上述函数可表示为
num仁[21]%（注意：
方括号，同一行的各元素间留空格或逗号）。
den1=[1221]
syss1=tf（num1,den1）
运行后，返回传递函数G1（s）的形式。
这种形式不能直接进行符号运算!
2．零极点增益形式
[Z,P,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=[21]%传递函数的分子系数向量den1=[1221]%传递函数的分母系数向量sys1=tf（num1,den1）%传递函数的有理分式模型[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）
%状态空间模型转换成零极点增益模型
程序中，命令tf2ss,zp2ss及ss2tf,ss2zp是状态空间模型与有理分式及零、极点、增益模型之间的相互转换。
、传递函数框图的处理
用框图可以方便地表示传递函数的并联，串联及反馈。
为简洁，仅以有理分式模
型为例。
1.并联
sysp=parallel（sys1,sys2）
[num,den]=parallel（num1,den1,num2,den2）
2.串联
syss=series（sys1,sys2）
[nums,dens]=series（num1,den1,num2,den2）
3.反馈
G1（S）G2（S）
1+G1（s）G2（s）G3（s）
sysc=feedback（syss,sys3,±1）%默认值（-1）[nume,dene]=feedback（nums,dens,num3,den3）
Gl（S）G2（S）
1+Gl（S）G2（S）
sysd=feedback（syss,1）
[numd,dend]=feedback（nums,dens,1,1）%（单位反馈）
上面给出了同一指令的两种形式，相当于两套平行指令。
对于零极点增益形式，书写稍复杂一些，可先用zpk转换成系统形式，或用zp2tf
转折换成有理分式形式后再进行框图化简操作。
三、简单函数的拉普拉斯变换
在MATLAB的符号功能中，可以对简单函数进行拉普拉斯正、逆变换。
拉氏正变换：
laplace（f（t））
拉氏逆变换：
ilaplace（L（s））
其中f（t）为原函数，L（s）为象函数。
命令格式参见fanli007。
参考程序fanli007:
拉普拉斯变换
symsstwabc
f1=sqrt（（b-a）A2+wA2）/w*exp（-a*t）*sin（w*t+atan（w/（b-a）））%原函数fl
L1=laplace（f1）
%f1的拉氏变换（象函数）
L1=simple（L1）
%化简
f2=ilaplace（L1）
%L1的拉氏逆变换
f2=simple（f2）%化简
在MATLAB中使用laplace及ilaplace命令时，要注意象、原函数的符号，特别是对初相不等于零的振荡系统，运行结果常常同手册上的结果相差一个符号，这要注意函数表达式成立的条件。
保险的办法是再使用拉氏变换的初值定理确定象、原函数的符号。
2-2系统时域特性曲线
在MATLAB中，当传递函数已知时，可以方便地求出系统的单位脉冲响应、单位阶跃响应等曲线。
一、系统的单位阶跃响应step
step有以下几种格式
step（sys）：
直接作出sys的单位阶跃响应曲线。
其中sys=tf（num,den）或sys=zpk（z,p,k），MATLAB自动决定响应时间。
step（sys,t）
设定响应时间的单位阶跃响应。
t可以设定为最大响应时间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后面跟的两个参数横坐标在前，纵坐标在后。
参考程序见fanli008：
参考程序fanli008:
系统的单位阶跃响应
num1=[42]
den1=[2814114]
sys1=tf（num1,den1）%系统G1（s）num2=[21]
den2=[14673]
sys2=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）在自选时间（20秒）内的单位阶跃响曲线
figure,plot（t1,y1）%系统G1（s）的单位阶跃响应曲线
figure,plot（t2,y2）%系统G2（s）的单位阶跃响应曲线
fanli008:
step（sys1,sys2,t）单位阶跃曲线
二、系统的单位脉冲响应
impulse命令格式与单位阶跃响应step的命令格式完全相同，只需将语句中的
step用impulse代替即可。
针对同样的系统，其单位脉冲响应的参考程序见
fanli009。
参考程序fanli009:
系统的单位脉冲响应
num仁[42]
den仁[2814114]
sys仁tf（num1,den1）
num2=[21]
den2=[14673]
sys2=tf（num2,den2）
[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）的单位脉冲响应曲线
holdon,step（sys2）%系统G2（s）的单位阶跃和单位脉冲响应曲线
Time（sec）
fanli009:
impulse（sys1,sys2）单位脉冲响应曲线
程序的最后一句holdon是当前图形保护模式。
当要将新图形作在当前图形上时，必须使用holdon。
而figure的含意是另开一个新的图形窗口，如果不用figure或holdon，则新的图形会占用原图形窗口，始终只保留一个最新的图形窗口。
三、一阶系统及二阶系统的时域特性
一阶系统及二阶系统是最基本也是最重要的系统，高阶系统总可以
视为由若干个一阶和（或）二阶系统组合构成。
1.一阶系统（设增益为1）
1
G（s厂
Ts+1
影响系统特性的参数是其时间常数T,T越大，系统惯性越大，响应越慢。
参考程序fanli010给出了T=0.4,1.2,2.0,2.8,3.6,4,4六条单位阶跃响应曲线。
参考程序fanli010:
—阶系统的单位阶跃响应曲线
num=1;i=1;
fordel=0.1:
0.2:
1.1%一阶系统时间常数递增间隔
den=[4*del1];%一阶系统分母向量
step（tf（num,den））%一阶系统单位阶跃响应曲线
holdon,%不同时间常数的一阶系统单位阶跃响应曲线簇
i=i+1;
end
同理，可以作出对应的单位脉冲响应曲线，参考程序fanli011。
参考程序fanli011:
—阶系统的单位脉冲响应曲线
num=1;i=1;
StepResponse
edur-rpmA
20
25
fanli010—阶系统时间常数对单位阶跃响应的影响
fanli011—阶系统时间常数对单位脉冲响应的影响
fordel=0.1:
0.2:
1.1
den=[4*del1];
impulse（tf（num,den）,10）%一阶系统单位阶脉冲应曲线
holdon,%不同时间常数的一阶系统单位脉冲响应曲线簇
i=i+1;
end
注意MATLAB
中for语句的结构。
读者可以改变不同的增益,
看看图形有何变化。
2.二阶系统（设0设二阶系统为
G（s）
二阶系统的特征参数为固有频率''n及阻尼比E。
当'n增大，
升级会员 