《现代控制理论》MATLAB实践指导书.docx

《《现代控制理论》MATLAB实践指导书.docx》由会员分享，可在线阅读，更多相关《《现代控制理论》MATLAB实践指导书.docx（15页珍藏版）》请在课桌文档上搜索。

1、现代控制理卷富“2名实践播导零现代控制理论MATLAB实践指导书1MATLAB概述MATLAB是MATriXLABOratory的缩写，早期主要用于现代控制中复杂的矩阵、向量的各种运算。由于MATLAB提供了强大的矩阵处理和绘图功能，很多专家因此在自己擅长的领域用它编写了许多专门的MATLAB工具包(toolbox),如控制系统工具包(COntrOlsystemstoolbox)；系统辨识工具包(SyStemidentifiCationtoOlbox)；信号处理工具包(Signalprocessingtoolbox)；鲁棒控制工具包(robustcontroltoolbox)；最优化工具包(O

2、PtimiZatiOntoolbox)等等由于MATLAB功能的不断扩展，所以现在的MATLAB已不仅仅局限与现代控制系统分析和综合应用，它己是一种包罗众多学科的功能强大的“技术计算语言(TheLanguageOfTechnicalComputing)”。MathWorks公司于1992年推出了具有划时代意义的MATLAB4.0版本，并推出了交互式模型输入与仿真系统Simulink,它使得控制系统的仿真与CAD应用更加方便、快捷，用户可以方便地在计算机上建模和仿真实验。1997年MathWOrkS推出的MATLAB5.0版允许了更多的数据结构，1999年初推出的MATLAB5.3版在很多方面又

3、进一步改进了MATLAB语言的功能。2000年底推出的MATLAB6.0o最新版本是MATLAB7.0。MATLAB以矩阵作为基本编程单元，它提供了各种矩阵的运算与操作，并有较强的绘图功能。MATLAB集科学计算、图像处理、声音处理于一身，是一个高度的集成系统，有良好的用户界面，并有良好的帮助功能。MATLAB不仅流行于控制界，在机械工程、生物工程、语音处理、图像处理、信号分析、计算机技术等各行各业中都有极广泛的应用。2如何获得MATLAB帮助在MATLAB主窗口中键入help,即可获得第一层帮助：help%加重型内容为用户键入的内容，其它为执行后显示的内容。HELPtopics:matlab

4、generalGeneralpurposecommands.matlabopsOperatorsandspecialcharacters.matlablangProgramminglanguageconstructs.matlabelmatElementarymatricesandmatrixmanipulation.matlabelfunElementarymathfunctions.matlabspecfunSpecializedmathfunctions.matlabmatfunMatrixfunctions-numericallinearalgebra.simulinksimulink

5、SimulinksimulinkblocksSimulinkblocklibrary.simulinksimdemos-Simulink3demonstrationsandsamples,simulinkdeeDifferentialEquationEditorMATLABR11work(Notableofcontentsfile)toolboxlocalPreferences.如果用户对MATLAB的语言结构Iang感兴趣，想进一步了解，则键入:helpIangProgramminglanguageconstructs.Controlflow.ifConditionallyexecutest

6、atements.elseIFstatementcondition.elseif-IFstatementcondition.endTerminatescopeofFOR,WHILE,SWITCH,TRYandIFstatements,forRepeatstatementsaspecificnumberoftimes.while-Repeatstatementsanindefinitenumberoftimes.如果想进一步了解for语句，则键入：helpforFORRepeatstatementsaspecificnumberoftimes.ThegeneralformofaFORstatem

7、entis:FORvariable=expr,statement,.,statementENDThecolumnsoftheexpressionarestoredoneatatimeinthevariableandthenthefollowingstatements,uptotheEND,areexecutedSomeexamples(assumeNhasalreadybeenassignedavalue).FORI=1:N,FORJ=1:N,A(U)=1(I+J-1);ENDEND同样，如果想了解MATLAB中有关矩阵的操作运算函数，可以键入：helpmatfunMatrixfunction

8、s-numericallinearalgebra.Matrixanalysis.normMatrixorvectornorm.normestEstimatethematrix2-norm.rankMatrixrank.detDeterminant.traceSumofdiagonalelements.nullNullspace.orthOrthogonalization.rrefReducedrowechelonform.subspaceAnglebetweentwosubspaces.Eigenvaluesandsingularvalues.eigEigenvaluesandeigenvec

9、tors.svdSingularvaluedecomposition.gsvdGeneralizedingularvaluedecomposition.eigsAfeweigenvalues.svdsAfewsingularvalues.polyCharacteristicpolynomial.polyeigPolynomialeigenvalueproblem.condeigConditionnumberwithrespecttoeigenvalues.hessHessenbergform.qzQZfactorizationforgeneralizedeigenvalues.schurSch

10、urdecomposition.Matrixfunctions.expmMatrixexponential.IogmMatrixlogarithm.sqrtmMatrixsquareroot.funmEvaluategeneralmatrixfunction.上面所列的都是有关矩阵的操作函数。如eig（八）可求出A的特征根及其特征向量，具体执行方法为：A=01;-6-5%输入A矩阵A=O1-6-5E=eig（八）%求出方阵A的特征根EE=-2-3V,D=eig（八）%求出方阵A的特征向量V及其A的对角型DV=0.4472-0.3162-0.89440.9487D=-200-33MATLAB基本

11、功能我们下面给出一些MATLAB的常用的功能,不过这只是MATLAB及其众多TOOLBOX中的极少极少部分。用户可以参阅有关MATLAB的手册,或直接在MATLAB系统中用HELP命令查阅其它功能。3. 1MATLAB的主要线性代数运算如表A-I所示为常用的矩阵和线性代数运算函数，用户可以用helpmatfun获得更多内容。表A1常用线性代数函数B=A,矩阵转置C=A+B矩阵相加C=A*B矩阵相乘C=Ak矩阵哥C=A.*B矩阵点乘，即两维数相同的矩阵各对应元素相乘expm（八）指数矩阵，也就是Vinv（八）矩阵的逆矩阵det（八）矩阵的行列式的值rank（八）计算矩阵的秩eig（八）矩阵的特征

12、值X,D=eig（八）矩阵的特征向量X和以特征值为元素的对角阵Dp=poly（八）矩阵的特征多项式r=roots(p)特征多项式方程的根conv(pl,p2)两多项式相乘3.2常用的控制系统处理函数(1) TF2SS将传递函数转换到状态空间表达式A,B,C,D=TF2SS(NUM,DEN)将系统：=NUMG)=U/+b7+4$+%DEN(三)sn+a“_ST+qs+4转换成：X=AX+BUY=CX+DU其中：NUM=bm,bm-1,bl,b，DEN=fl,an-l,an-2,.,a1,aan-1an-20一。.0ao0-10A=:1:B=0.00000100C=KT%D=b(2) ZP2SS将

13、零极点型传递函数转换到状态空间表达式A,B,C,D=ZP2SS(Z)HK)除了G(三)=K(S+Z)(s+z?)(s+z”)以外,其它与tf2SS相同。(s+p)(s+%)(s+P)(3) SS2TF将状态空间表达式转换到传递函数NUM,DEN=SS2TF(ABC,D,iu)即求第iu个输入信号对输出y(t)的传递函数，即：gcv)=2W =m u，S) DENG)= C(Sl-A) B +D =NUM (S)s 1 + + qs + &SS2TF的调用返回值为G(三)的分子多项式的系数矩阵NUM和分母多项式的系数向量DENo(4) SS2ZP将状态空间表达式到零极点形式的传递函数的转换Z,P

14、,K=SS2ZP(A,B,C,D,iu)(5) )TF2ZP一般传递函数转换到零极点型传递函数Z,P,K=TF2ZP(NUM,DEN)(6) ZP2TF零极点型传递函数转换到一般传递函数NUM,DENJ=ZP2TF(Z,RK)(7) SS2SS状态空间表达式的线性变换A1,BLCI,D1=ss2ss(A,B,C,D,T)其中T为变换矩阵。注意变换方程为：Xl=TX,而不是常见的X=TX1。所以要与用户习惯的变换方程一致，我们必须用T的逆代入上式，即：Al,BLCl,Dl=ss2ss(A,B,C,Djnv(T)(8) CANON求状态空间表达式的对角标准型As,Bs,Cs,Ds,Ts=canon

15、(A,B,C,Dmod)其中TS为变换矩阵，注意变换方程为：Xs=TsXo(9) CTRB计算系统的可控判别矩阵MM=ctrb(A,B)(10) OBSV计算系统的可观判别矩阵NN=obsv(A,C)(IDIMPULSE求系统的单位脉冲响应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)(14) C2

16、D连续系统状态方程转换为离散状态方程,T为采样周期QH=c2d(A,B,T)相关的函数还有D2CQ2D。(15) LYAP求解如下形式的李雅普诺夫方程：AP+PA=-Q所以求解用户习惯的AP+PA=-Q李雅普诺夫方程，我们必须用A的转置A，代入，即：P=Lyap(A,Q)O(16) PLACE极点配置F=PLACE(A,B,P)(17) PLOT绘图函数PlOt(X,YStr)可以用不同颜色、不同符号绘制曲线，其中st可以是下列三组选项的任意组合。yyellow.pointsolidm-magentaocircle:dottedccyanXx-mark-.-dashdotrred+-plusd

17、ashedg-green*-starb-blues-squarew-whited-diamondkblackv-triangle(down)triangle(up)- triangle(left)- -triangle(right)ppentagramh-hexagram4例题分析例11给定某控制系统的状态空间描述为：- o1-11oX=-6-116x+Ow- 6-115J1y=l00卜试求对角规范型和变换矩阵W,并根据其对角规范型绘制系统在初值X=51015,W=O时的响应曲线。解：closeall;%状态方程模型A=014;-6-116;-6-H5;B=0;0;l;C=l00;%求取对角规

18、范型W,lamda=eig（八）;L=inv(W)*A*W;b=inv(W)*B;c=C*W;%显示结果disp(TheDiagonalCanonicalFromofSystemis:);Lbcdisp(TransformationMatrixis:);W%仿真数据初始化t=0:0.01:3;x0=5J0J5:xx0=inv(W)*x0;n=length(t);x=zeros(3,n);xx=zeros(3,n);%求解状态变量fori=l:nxx(:,i)=expm(L*t(i)*xxO;x(:,i)=W*xx(:,i);end%绘图PIOt(I,)title(SystemTimeRespo

19、nsewithu=0);XIabe1(Time-Sec);ylabel(Response-value);IeXI(0.8,3.7,x(1);IeXl(0.8,15x(3);text(0.8,-0.4,x(2);TheDiagonalCanonicalFromofSystemis:L=-LOOOO-0.0000-0.00000.0000-2.0000-0.0000-0.00000.0000-3.0000b=-2.8284-13.747710.8628c=0.7071-0.2182-0.0921TransformationMatrixis:W=0.7071-0.2182-0.09210.0000-

20、0.4364-0.55230.7071-0.8729-0.8285System Time Response with u=016111114 r3n-g,asuod00:分析：通过对角规范型求借系统状态方程的优点是，不必计算矩阵指数，只需计算单个的标量的指数即可。在本例具体的求解的过程中，需要注意的是求得对角规范型的状态XX的时间响应之后，还要用变换矩阵W将其转化为原状态变量X,然后画图。从系统的响应曲线来看，各状态变量最终都趋于零。这一方面是因为输入量U=0,另一方面是因为初值不为零，但矩阵A的特征值均为零，系统的零输入响应呈衰减趋势。但状态变量X(3)下降幅度和速度最大，X(1)下降幅度和

21、速度最小。这是因为状态变量X(3)对应的特征值是-3,而X(1)对应的特征值是-1,其衰减速率不同。例2给定线性定常系统的状态空间描述为：2001104130X=X+U004110002J_2_=l10Ok试判断系统的能控性。如果系统状态完全能控，试求其能控标准型和变换矩阵。解：%系统状态方程模型A=2001;0413;0041;0002;B=l;0；l；2;C=1100:%系统阶次n=length（八）;%求解系统能控性矩阵Q=ctrb(A,B);%系统能控性矩阵的秩m=rank(Q);%判断系统是否状态完全能控，并求解能控标准型ifm=nAc1=inv(Q)*A*Q;Bcl=inv(Q)*

22、B;Ccl=C*Q;disp(SystemisControllable.);disp(SystemFirstControllableCanonnicalFonnis:);AclBclCcldisp(TheTransfonnationMatrixis:);Qelsedisp(SystemStateVariablecannotbetotallyControlled.);disp(TherankofSystemControllableMatrixis:);mend运行结果：SystemisControllable.SystemFirstControllableCanonnicalFormis:Acl

23、=O-0.0000-0.0000-64.00001.0000-0.0000-0.000096.000001.0000-0.0000-52.00000.00000.00001.000012.0000Bcl=1000Ccl=11158268TheTransformationMatrixis:Q=1412320746236162812024816例3对于上例的线性定常系统，试判断系统的能观性，如果系统状态完全能观的，试求其能观标准型和变换矩阵。解：系统状态方程模型A=200l;04I3;0041;0002;B=l;0;l;2;C=l100;%系统阶次n=length（八）;P=obsv(A,C);%

24、系统能观性矩阵的秩m=rank(P);ifm=nAo1=P*A*inv(P);Bol=P*B;Col=C*inv(P);disp(SystemisObservable.);disp(SystemFirstObservableCanonnicalFormis:*);AolBolColdisp(TheTransformationMatrixis:)；inv(P)elseenddisp(,System State Variable cannot be totally Observed.); disp(The rank of System Observable Matrix is:):运行结果:yst

25、em is Observable.System First Observable Canonnical Fonn is:Aol =0-649612Il58268ColTheTransformationMatrixans = -14.000016.0000-5.37500.562515.0000 -16.00005.3750-0.562501.0000-0.75000.1250-8.00008.0000-2.50000.2500例4某控制系统的状态方程描述如下。试判断其稳定性并绘制其时间响应来验证上述判断。-3-8-2-4100001000010Lly=解：z=-lp=-1.4737+2.263

26、8i-1.4737-2.2638i-0.0263+0.7399i-0.0263-0.7399ik=lSystemisstableSystemTimeRespone0.611；：分析：系统的极点均在复平面的左半平面，因而系统是稳定的。从系统的时间响应曲线上看也是如此，系统的输出振荡衰减，但可以看出，系统衰减响应曲线的振荡频率很高，衰减速率很慢。这是因为系统有两个极点的实部为-0.0263,非常靠近虚轴，在实际过程中，这种系统属于性能非常差的一类。稳定裕量很小，稍有干扰系统的极点就可能越过虚轴跑到复平面的右半面去。从而造成系统的不稳定。对于这类系统，一般都要加校正手段，使其极点远里离虚轴。例5含积

27、分环节的伺服系统设计，设对象为G(三)=!s(s+l)(s+2)设计控制器“=-Kx+kr,使闭环系统具有极点-223-10。解：%Poleplacementcloseall;a=0I0;00l;0-2-3;b=0;0;l;c=l00J;d=0;disp(Therankofcontrollabilitymatrix1)rc=rank(ctrb(a,b)p=-2+2*sqrt(3)*j-2-2*sqrt(3)*j-10;k=place(a,b)a1=a-b*k;b1=b*k(l);cl=c;dI=d;figurestep(al,bl,cl,dl)title(,StepResponeseofDes

28、ignedServoSystem)figurey,x,t=ste(al,bl,cl,dl);pot(t,x)title(,StepResponeseCurvesforxl,x2,x3)gridonTherankofcontrollabilitymatrixrc=3k=160.000054.000011.0000StepResponeseofDesignedServoSystem1.411.2OL-1-O0.611.522.63TimA伍m例6开环系统的状态方程为-010比=OOl-6-11-6y=10Ok0+ 0设计全维状态观测器，使观测器的闭环极点为-2/2百,-5o解：状态观测器的结构见教

29、材图55,为求出状态观测器的反馈增益阵K,先为原系统构造一对偶系统p=a,p+c,vq=p然后采用极点配置方法对对偶系统进行闭环极点培植，得到反馈增益阵k,从而可由对偶原理得到原系统的状态观测器增益阵勺=/。程序如下:%Designofafull-orderstateobservera=010;00-6;b=0;0;l;c=l00;%Checkobservabilitydisp(Therankofobservabilitymatrix)rt)=rank(obsv(a,c)%Constructthedaulsystem(al,bl,cl)al=a;bl=c,;cl=b;%sloveKusingp

30、olesplacementp=-2+2*sqrt(3)*j-2-2*sqrt(3)*j-5;k=acker(aLbl,p);ke=k运行结果：Therankofobservabilitymatrixr=3ke=3.00007.0000-1.00005上机习题5-1假设系统由两个Gl（三）和Gz（s）串联连接而成，已知（S+1FlL5+35+5G（5）JLG2（三）-7752+3s+4.V4+4.v3+3.S2+2s+1若想求出总系统的状态方程模型，请在MATLAB下比较下面两种方法有何不同结果：（1） 将两个传递函数模型进行串联连接，然后求出整个系统的状态方程模型。（2） 求出两个模型的状态方

31、程表示，然后求出整个系统的状态方程模型。5-2某位置控制系统方块图如下图所示，求其状态空间表达式。53试将5-2题所得动态方程转换成Jordan标准型,设k=2,k 1=-0.05, k2=0.16, k3=.24,并试用MATLAB语言求解。5-4设k=2,kl=-0.05,k2=O.I6,k3=24,用MATLAB求题5-2所示系统的单位阶跃响应。55判断卜.列系统的能控性和能观性，并将系统按能控性和能观性进行分解(1)11-1-2-21y = l -1 Ok(2)-3-1001-5-2-25-6某亚音速战斗机状态空间模型为：A=-O.OO8I9,-25.70839,0,-32.17095

32、；-0.00019,-1.27626,1.0000,0；0.00069,1.02176,-2.40523,0；0,0,1.0000,0B=L-6.80939;-0.14968；-14.06111；0C=l,0,0,0；0,1,0,0；0,0,1,0；0,0,0,1试判断系统的稳定性、能控性、能观性，可借用MATLAB编程。5-7试绘制下列系统模型的单位脉冲响应和状态变量的脉冲响应。-2119-20019-2120X +140-40-402y = l O 2X58对给定的对象模型2OOO12OOOO-1O-1O111请设计出一个状态反馈向量K,使得闭环系统的极点配置到(-2,-2,-1,-1)位置,另外，如果我们想将系统的所有极点均配置到-2,这样的配置是否可行，请解释原因。-O101ri59被控对象X=x+U,y=l0L20.6OjL1J设计调节器使闭环极点为-1.82.4,而且状态不可测，因此设计状态观测器使其闭环极点为-8,-8。510给定线性系统-l-2-22X=0-11X+0，y=1IOk10-IjH(1) 确定它的全维状态观测器，并使观测器的极点为刀=42=-3,43=-4。(2) 在前面所求得的状态观测器的基础上，设计一个观测器一一状态反馈控制系统，使期望的闭环极点m=22=T,l3=-2,并画出系统的单位阶跃响应曲线。