机械工程控制理论MATLAB上机实验指导1.docx
《机械工程控制理论MATLAB上机实验指导1.docx》由会员分享，可在线阅读，更多相关《机械工程控制理论MATLAB上机实验指导1.docx（22页珍藏版）》请在冰点文库上搜索。
机械工程控制理论MATLAB上机实验指导1
《机械工程控制理论》
MATLAB实践指导书
《机械工程控制理论》MATLAB实践指导书
1MATLAB概述
MATLAB是MATrixLABoratory的缩写，早期主要用于现代控制中复杂的矩阵、向量的各种运算。
由于MATLAB提供了强大的矩阵处理和绘图功能，很多专家因此在自己擅长的领域用它编写了许多专门的MATLAB工具包（toolbox），如控制系统工具包（controlsystemstoolbox）；系统辨识工具包（systemidentificationtoolbox）；信号处理工具包（signalprocessingtoolbox）；鲁棒控制工具包（robustcontroltoolbox）；最优化工具包（optimizationtoolbox）等等。
由于MATLAB功能的不断扩展，所以现在的MATLAB已不仅仅局限与现代控制系统分析和综合应用，它已是一种包罗众多学科的功能强大的“技术计算语言（TheLanguageofTechnicalComputing）”。
MathWorks公司于1992年推出了具有划时代意义的MATLAB4.0版本，并推出了交互式模型输入与仿真系统SIMULINK，它使得控制系统的仿真与CAD应用更加方便、快捷，用户可以方便地在计算机上建模和仿真实验。
1997年MathWorks推出的MATLAB5.0版允许了更多的数据结构，1999年初推出的MATLAB5.3版在很多方面又进一步改进了MATLAB语言的功能。
2000年底推出的MATLAB6.0。
最新版本是MATLAB7.0。
MATLAB以矩阵作为基本编程单元，它提供了各种矩阵的运算与操作，并有较强的绘图功能。
MATLAB集科学计算、图像处理、声音处理于一身，是一个高度的集成系统，有良好的用户界面，并有良好的帮助功能。
MATLAB不仅流行于控制界，在机械工程、生物工程、语音处理、图像处理、信号分析、计算机技术等各行各业中都有极广泛的应用。
2如何获得MATLAB帮助
在MATLAB主窗口中键入help，即可获得第一层帮助：
help％加重型内容为用户键入的内容，其它为执行后显示的内容。
HELPtopics:
matlab\general------------Generalpurposecommands.
matlab\ops------------------Operatorsandspecialcharacters.
matlab\lang-----------------Programminglanguageconstructs.
matlab\elmat---------------Elementarymatricesandmatrixmanipulation.
matlab\elfun----------------Elementarymathfunctions.
matlab\specfun-----------Specializedmathfunctions.
matlab\matfun-------------Matrixfunctions-numericallinearalgebra.
simulink\simulink-------Simulink
simulink\blocks----------Simulinkblocklibrary.
simulink\simdemos---Simulink3demonstrationsandsamples.
simulink\dee--------------DifferentialEquationEditor
MATLABR11\work------（Notableofcontentsfile）
toolbox\local--------------Preferences.
如果用户对MATLAB的语言结构lang感兴趣，想进一步了解，则键入：
helplang
Programminglanguageconstructs.
Controlflow.
if-----------Conditionallyexecutestatements.
else------IFstatementcondition.
elseif----IFstatementcondition.
end------TerminatescopeofFOR,WHILE,SWITCH,TRYandIFstatements.
for--------Repeatstatementsaspecificnumberoftimes.
while---Repeatstatementsanindefinitenumberoftimes.
如果想进一步了解for语句，则键入：
helpfor
FORRepeatstatementsaspecificnumberoftimes.
ThegeneralformofaFORstatementis:
FORvariable=expr,statement,...,statementEND
Thecolumnsoftheexpressionarestoredoneatatimeinthevariableandthenthefollowingstatements,uptotheEND,areexecuted.……Someexamples（assumeNhasalreadybeenassignedavalue）.
FORI=1:
N,
FORJ=1:
N,
A（I,J）=1/（I+J-1）;
END
END
同样，如果想了解MATLAB中有关矩阵的操作运算函数，可以键入：
helpmatfun
Matrixfunctions-numericallinearalgebra.
Matrixanalysis.
norm---------Matrixorvectornorm.
normest----Estimatethematrix2-norm.
rank----------Matrixrank.
det------------Determinant.
trace---------Sumofdiagonalelements.
null------------Nullspace.
orth-----------Orthogonalization.
rref------------Reducedrowechelonform.
subspace--Anglebetweentwosubspaces.
Eigenvaluesandsingularvalues.
eig------------Eigenvaluesandeigenvectors.
svd-----------Singularvaluedecomposition.
gsvd---------Generalizedingularvaluedecomposition.
eigs----------Afeweigenvalues.
svds---------Afewsingularvalues.
poly----------Characteristicpolynomial.
polyeig-----Polynomialeigenvalueproblem.
condeig----Conditionnumberwithrespecttoeigenvalues.
hess---------Hessenbergform.
qz------------QZfactorizationforgeneralizedeigenvalues.
schur-------Schurdecomposition.
Matrixfunctions.
expm-------Matrixexponential.
logm--------Matrixlogarithm.
sqrtm-------Matrixsquareroot.
funm--------Evaluategeneralmatrixfunction.
上面所列的都是有关矩阵的操作函数。
如eig（A）可求出A的特征根及其特征向量，具体执行方法为：
A=[01;-6-5]％输入A矩阵
A=
0
1
-6
-5
E=eig（A）%求出方阵A的特征根E
E=
-2
-3
[V,D]=eig（A）%求出方阵A的特征向量V及其A的对角型D
V=
0.4472
-0.3162
-0.8944
0.9487
D=
-2
0
0
-3
3MATLAB基本功能
我们下面给出一些MATLAB的常用的功能，不过这只是MATLAB及其众多TOOLBOX中的极少极少部分。
用户可以参阅有关MATLAB的手册，或直接在MATLAB系统中用HELP命令查阅其它功能。
3．1MATLAB的主要线性代数运算
如表A－1所示为常用的矩阵和线性代数运算函数，用户可以用helpmatfun获得更多内容。
表A－1常用线性代数函数
B=A’
矩阵转置
C=A+B
矩阵相加
C=A*B
矩阵相乘
C=A^k
矩阵幂
C=A.*B
矩阵点乘，即两维数相同的矩阵各对应元素相乘
expm（A）
指数矩阵，也就是
inv（A）
矩阵的逆矩阵
det（A）
矩阵的行列式的值
rank（A）
计算矩阵的秩
eig（A）
矩阵的特征值
[X,D]=eig（A）
矩阵的特征向量X和以特征值为元素的对角阵D
p=poly（A）
矩阵的特征多项式
r=roots（p）
特征多项式方程的根
conv（p1,p2）
两多项式相乘
3．2常用的控制系统处理函数
（1）TF2SS将传递函数转换到状态空间表达式
[A,B,C,D]=TF2SS（NUM,DEN）
将系统：
转换成：
其中:
NUM=[bm,bm-1,…,b1,b0]，DEN=[1,an-1,an-2,…,a1,a0]
（2）ZP2SS将零极点型传递函数转换到状态空间表达式
[A,B,C,D]=ZP2SS（Z,P,K）
除了
以外，其它与TF2SS相同。
（3）SS2TF将状态空间表达式转换到传递函数
[NUM,DEN]=SS2TF（A,B,C,D,iu）
即求第iu个输入信号对输出y（t）的传递函数，即：
SS2TF的调用返回值为G（s）的分子多项式的系数矩阵NUM和分母多项式的系数向量DEN。
（4）SS2ZP将状态空间表达式到零极点形式的传递函数的转换
[Z,P,K]=SS2ZP（A,B,C,D,iu）
（5）TF2ZP一般传递函数转换到零极点型传递函数
[Z,P,K]=TF2ZP（NUM,DEN）
（6）ZP2TF零极点型传递函数转换到一般传递函数
[NUM,DEN]=ZP2TF（Z,P,K）
（7）SS2SS状态空间表达式的线性变换
[A1,B1,C1,D1]=ss2ss（A,B,C,D,T）
其中T为变换矩阵。
注意变换方程为：
X1=TX，而不是常见的X=TX1。
所以要与用户习惯的变换方程一致，我们必须用T的逆代入上式，即：
[A1,B1,C1,D1]=ss2ss（A,B,C,D,inv（T））
（8）CANON求状态空间表达式的对角标准型
[As,Bs,Cs,Ds,Ts]=canon（A,B,C,D,'mod'）
其中Ts为变换矩阵，注意变换方程为：
Xs=TsX。
（9）CTRB计算系统的可控判别矩阵M
M=ctrb（A,B）
（10）OBSV计算系统的可观判别矩阵N
N=obsv（A,C）
（11）IMPULSE求系统的单位脉冲响应
[y,x]=impulse（A,B,C,D,in,t）
[y,x]=impulse（num,den,t）
（12）STEP求系统的单位阶跃响应
[y,x]=step（A,B,C,D,in,t）
[y,x]=step（num,den,t）
（13）LSIM求系统对任意输入函数u（t）的响应
[y,x]=lsim（A,B,C,D,u,t）
[y,x]=lsim（num,den,u,t）
（14）C2D连续系统状态方程转换为离散状态方程,T为采样周期
[G,H]=c2d（A,B,T）
相关的函数还有D2C,D2D。
（15）LYAP求解如下形式的李雅普诺夫方程：
所以求解用户习惯的
李雅普诺夫方程，我们必须用A的转置A’代入，即：
P=Lyap（A’，Q）。
（16）PLACE极点配置
F=PLACE（A，B，P）
（17）PLOT绘图函数
plot（X,Y,’str’）可以用不同颜色、不同符号绘制曲线，其中‘str’可以是下列三组选项的任意组合。
y-----yellow
.------point
------solid
m----magenta
o------circle
:
-----dotted
c-----cyan
x------x-mark
-.----dashdot
r------red
+----plus
------dashed
g----green
*----star
b----blue
s---square
w---white
d---diamond
k----black
v---triangle（down）
^-----triangle（up）
<----triangle（left）
>----triangle（right）
p----pentagram
h----hexagram
4例题分析
[例1]给定某控制系统的状态空间描述为：
试求对角规范型和变换矩阵W，并根据其对角规范型绘制系统在初值
，
时的响应曲线。
解：
closeall;
%状态方程模型
A=[01-1;-6-116;-6-115];B=[0;0;1];C=[100];
%求取对角规范型
[W,lamda]=eig（A）;L=inv（W）*A*W;b=inv（W）*B;c=C*W;
%显示结果
disp（'TheDiagonalCanonicalFromofSystemis:
'）;
L
b
c
disp（'TransformationMatrixis:
'）;
W
%仿真数据初始化
t=0:
0.01:
3;x0=[5;10;15];xx0=inv（W）*x0;n=length（t）;
x=zeros（3,n）;xx=zeros（3,n）;
%求解状态变量
fori=1:
n
xx（:
i）=expm（L*t（i））*xx0;
x（:
i）=W*xx（:
i）;
end
%绘图
plot（t,x）
title（'SystemTimeResponsewithu=0'）;
xlabel（'Time-sec'）;ylabel（'Response-value'）;
text（0.8,3.7,'x
（1）'）;text（0.8,1.5,'x（3）'）;text（0.8,-0.4,'x
（2）'）;
>>TheDiagonalCanonicalFromofSystemis:
L=-1.0000-0.0000-0.0000
0.0000-2.0000-0.0000
-0.00000.0000-3.0000
b=-2.8284
-13.7477
10.8628
c=0.7071-0.2182-0.0921
TransformationMatrixis:
W=0.7071-0.2182-0.0921
0.0000-0.4364-0.5523
0.7071-0.8729-0.8285
分析：
通过对角规范型求借系统状态方程的优点是，不必计算矩阵指数，只需计算单个的标量的指数即可。
在本例具体的求解的过程中，需要注意的是求得对角规范型的状态xx的时间响应之后，还要用变换矩阵W将其转化为原状态变量x，然后画图。
从系统的响应曲线来看，各状态变量最终都趋于零。
这一方面是因为输入量u=0，另一方面是因为初值不为零，但矩阵A的特征值均为零，系统的零输入响应呈衰减趋势。
但状态变量x（3）下降幅度和速度最大，x
（1）下降幅度和速度最小。
这是因为状态变量x（3）对应的特征值是-3，而x
（1）对应的特征值是-1，其衰减速率不同。
[例2]给定线性定常系统的状态空间描述为：
试判断系统的能控性。
如果系统状态完全能控，试求其能控标准型和变换矩阵。
解：
%系统状态方程模型
A=[2001;0413;0041;0002];
B=[1;0;1;2];
C=[1100];
%系统阶次
n=length（A）;
%求解系统能控性矩阵
Q=ctrb（A,B）;
%系统能控性矩阵的秩
m=rank（Q）;
%判断系统是否状态完全能控，并求解能控标准型
ifm==n
Ac1=inv（Q）*A*Q;
Bc1=inv（Q）*B;
Cc1=C*Q;
disp（'SystemisControllable.'）;
disp（'SystemFirstControllableCanonnicalFormis:
'）;
Ac1
Bc1
Cc1
disp（'TheTransformationMatrixis:
'）;
Q
else
disp（'SystemStateVariablecannotbetotallyControlled.'）;
disp（'TherankofSystemControllableMatrixis:
'）;
m
end
运行结果：
SystemisControllable.
SystemFirstControllableCanonnicalFormis:
Ac1=0-0.0000-0.0000-64.0000
1.0000-0.0000-0.000096.0000
01.0000-0.0000-52.0000
0.00000.00001.000012.0000
Bc1=1
0
0
0
Cc1=11158268
TheTransformationMatrixis:
Q=141232
0746236
1628120
24816
[例3]对于上例的线性定常系统，试判断系统的能观性，如果系统状态完全能观的，试求其能观标准型和变换矩阵。
解：
%系统状态方程模型
A=[2001;0413;0041;0002];
B=[1;0;1;2];
C=[1100];
%系统阶次
n=length（A）;
P=obsv（A,C）;
%系统能观性矩阵的秩
m=rank（P）;
ifm==n
Ao1=P*A*inv（P）;
Bo1=P*B;
Co1=C*inv（P）;
disp（'SystemisObservable.'）;
disp（'SystemFirstObservableCanonnicalFormis:
'）;
Ao1
Bo1
Co1
disp（'TheTransformationMatrixis:
'）;
inv（P）
else
disp（'SystemStateVariablecannotbetotallyObserved.'）;
disp（'TherankofSystemObservableMatrixis:
'）;
m
end
运行结果：
ystemisObservable.
SystemFirstObservableCanonnicalFormis:
Ao1=0100
0010
0001
-6496-5212
Bo1=1
11
58
268
Co1=1000
TheTransformationMatrixis:
ans=-14.000016.0000-5.37500.5625
15.0000-16.00005.3750-0.5625
01.0000-0.75000.1250
-8.00008.0000-2.50000.2500
[例4]某控制系统的状态方程描述如下。
试判断其稳定性并绘制其时间响应来验证上述判断。
解：
z=-1
p=-1.4737+2.2638i
-1.4737-2.2638i
-0.0263+0.7399i
-0.0263-0.7399i
k=1
Systemisstable
分析：
系统的极点均在复平面的左半平面，因而系统是稳定的。
从系统的时间响应曲线上看也是如此，系统的输出振荡衰减，但可以看出，系统衰减响应曲线的振荡频率很高，衰减速率很慢。
这是因为系统有两个极点的实部为-0.0263，非常靠近虚轴，在实际过程中，这种系统属于性能非常差的一类。
稳定裕量很小，稍有干扰系统的极点就可能越过虚轴跑到复平面的右半面去。
从而造成系统的不稳定。
对于这类系统，一般都要加校正手段，使其极点远里离虚轴。
[例5]含积分环节的伺服系统设计，设对象为
设计控制器
，使闭环系统具有极点
。
解：
%Poleplacement
closeall;
a=[010;001;0-2-3];
b=[0;0;1];
c=[100];
d=[0];
disp（'Therankofcontrollabilitymatrix'）
rc=rank（ctrb（a,b））
p=[-2+2*sqrt（3）*j-2-2*sqrt（3）*j-10];
k=place（a,b,p）
a1=a-b*k;b1=b*k
（1）;c1=c;d1=d;
figure
（1）
step（a1,b1,c1,d1）
title（'StepResponeseofDesignedServoSystem'）
升级会员 