线性控制系统的计算机辅助分析.pptx

会计学1线性控制系统线性控制系统(kn zh x tn)的计算机的计算机辅助分析辅助分析第一页，共109页。2对于离散对于离散(lsn)(lsn)时间时间系统，系统，其状态变量的解析(ji x)解为要使x(kT)有界，则要Fk有界，即F矩阵的所有特征根的模均小于1。故如果系统全部极点都位于Z平面的单位(dnwi)圆内，则系统是稳定的。若连续时间系统的全部零极点都位于S左半平面；或若离散时间系统的全部零极点都位于Z平面单位(dnwi)圆内，则系统是最小相位系统。MATLAB中，可用eig()函数直接求取系统的特征根，也可用pzmap()函数直接绘制系统的零极点。第2页/共109页第二页，共109页。3例例1 1：系统：系统(xt(xt ng)ng)传函为传函为试判断其稳定性。num=18 514 5982 36380 122664 222088 185760 40320;den=1 36 546 4536 22449 67284 118124 109584 40320;G=tf(num,den);eig(G)ans=-8.0000 -7.0000 -6.0000 -5.0000 -4.0000 -3.0000 -2.0000 -1.0000可见(kjin)，系统是稳定的。第3页/共109页第三页，共109页。4例例2 2：离散系统受控对象的传函为：离散系统受控对象的传函为控制器模型(mxng)为 试分析单位负反馈下的闭环系统稳定性。z=tf(z,0.1);G=0.00147635*(z2+3.4040929*z+0.71390672)/(z-1)*(z-0.535261429)*(z-0.951229425);Gc=1.5*(z-0.5)/(z+0.8);GG=feedback(G*Gc,1);eig(GG)abs(eig(GG)ans=-0.7991 0.7991 0.9745+0.0782i 0.9776 0.9745-0.0782i 0.9776 0.5344 0.5344 第4页/共109页第四页，共109页。5pzmap(GG)pzmap(GG)第5页/共109页第五页，共109页。64.1.2 线性系统的线性相似(xin s)变换 相似变换是对方阵进行的一种运算，它把A变成P-1AP，而可逆矩阵P称为进行这一变换的相似变换矩阵。这种变换的重要(zhngyo)意义在于简化对矩阵的各种运算，其方法是先通过相似变换，将矩阵变成与之等价的对角矩阵，再对对角矩阵进行运算，从而将比较复杂的矩阵的运算转化为比较简单的对角矩阵的运算。MATLAB控制系统工具箱提供了ss2ss()函数完成状态方程模型的相似变换：G1=ss2ss(G,T)其中，G为原始的状态方程模型，T为变换矩阵。第6页/共109页第六页，共109页。7例3：第7页/共109页第七页，共109页。8例例4 4：设系统：设系统(xt(xt ng)ng)的状态方程为的状态方程为实际应用中，若不要求将A变换为对角阵，则P也可用任意非奇异(qy)矩阵。变换矩阵P为反对角矩阵，反对角线上的元素(yun s)为1，其余元素(yun s)为0。另：首先介绍fliplr()函数，其变换矩阵行元素的左右顺序。如A=1 2 3 4 5 6 7 8 9 10 11 12fliplr(A)=4 3 2 1 8 7 6 5 12 11 10 9第8页/共109页第八页，共109页。9A=0 1 0 0;0 0 1 0;0 0 0 1;-24-50-35-10;A=0 1 0 0;0 0 1 0;0 0 0 1;-24-50-35-10;G1=ss(A,0;0;0;1,24 7 1 0,0);G1=ss(A,0;0;0;1,24 7 1 0,0);P=fliplr(eye(4);G2=ss2ss(G1,P)P=fliplr(eye(4);G2=ss2ss(G1,P)a=x1 x2 x3 x4 x1 -10 -35 -50 -24 x2 1 0 0 0 x3 0 1 0 0 x4 0 0 1 0b=u1 x1 1 x2 0 x3 0 x4 0c=x1 x2 x3 x4 y1 0 1 7 24d=u1 y1 0第9页/共109页第九页，共109页。104.1.3 4.1.3 线性系统的可控性分析线性系统的可控性分析线性系统的可控性分析线性系统的可控性分析(fnx)(fnx)现代控制理论用状态方程和输出方程描述系统，输出和输入构成系统的外部变量，而状态为系统的内部变量，系统就好比是一块集成电路芯片，内部结构可能十分复杂，物理量很多，而外部只有少数几个引脚，对电路内部物理量的控制和观测(gunc)都只能通过这为数不多的几个引脚进行。这就存在着系统内的所有状态是否都受输入控制和所有状态是否都可以从输出反映出来的问题。如果系统所有状态变量的运动都可以通过有限的控制点的输入来使其由任意的初态达到任意设定的终态，则称系统是可控的，更确切的说是状态可控的；如果系统中所有的状态都是可控的，则称该系统为完全可控的系统。否则，就称系统不可控。第10页/共109页第十页，共109页。111.1.线性系统的可控性判定线性系统的可控性判定线性系统的可控性判定线性系统的可控性判定(pndng)(pndng)可通过构造可控性判定矩阵可通过构造可控性判定矩阵若若TcTc为满秩矩阵，则系统为完全为满秩矩阵，则系统为完全(wnqun)(wnqun)可控的。如果该矩可控的。如果该矩阵不是满秩矩阵，则它的秩为系统的可控状态的个数。阵不是满秩矩阵，则它的秩为系统的可控状态的个数。可控性判定矩阵由可控性判定矩阵由Tc=ctrb(A,B)Tc=ctrb(A,B)函数构造。函数构造。rank()rank()函数可求出矩阵的秩。函数可求出矩阵的秩。例例5 5：试判断(pndun)系统的可控性。第11页/共109页第十一页，共109页。12若系统的若系统的GramGram矩阵是非奇异矩阵，则该系统是矩阵是非奇异矩阵，则该系统是完全可控的。完全可控的。Gram Gram矩阵为矩阵为GramGram矩阵是以下矩阵是以下LyapunovLyapunov方程的解方程的解求解该求解该LyapunovLyapunov方程可用，方程可用，lyap(A,B*B)lyap(A,B*B)若调用函数不能求出方程的解，则该系统若调用函数不能求出方程的解，则该系统不完全可控。不完全可控。控制系统的可控控制系统的可控GramGram矩阵还可以矩阵还可以(ky(ky)由由 Gc=gram(G,c)Gc=gram(G,c)直接求出。直接求出。第12页/共109页第十二页，共109页。13例例6 6：已知采样周期已知采样周期(zhuq)(zhuq)为为0.1s0.1s，求系统可控，求系统可控GramGram矩阵。矩阵。num=0.1324-0.5743 0.3879-0.0889;num=0.1324-0.5743 0.3879-0.0889;den=1-3.233 3.9869-2.2209 0.4723;den=1-3.233 3.9869-2.2209 0.4723;G=tf(num,den,Ts,0.1);G=tf(num,den,Ts,0.1);Lc=gram(ss(G),c);Lc=gram(ss(G),c);第13页/共109页第十三页，共109页。142.2.可控性阶梯可控性阶梯可控性阶梯可控性阶梯(jit)(jit)分解分解分解分解对于对于(duy)(duy)不完全可控的系统，可对其进行可控性阶梯不完全可控的系统，可对其进行可控性阶梯分解，即构造变换矩阵分解，即构造变换矩阵T T，将状态方程变换为，将状态方程变换为这样即将系统的可控子空间，与不可控子空间分离出来。这样即将系统的可控子空间，与不可控子空间分离出来。Ac,Bc,Cc,Tc=ctrbf(A,B,C)Ac,Bc,Cc,Tc=ctrbf(A,B,C)函数，可将系统变换为可控函数，可将系统变换为可控性阶梯模型，其中，性阶梯模型，其中，TcTc为相似变换矩阵。为相似变换矩阵。例例7 7：请将例：请将例4 4中的系统进行可控性阶梯分解。中的系统进行可控性阶梯分解。A=;B=;C=;A=;B=;C=;Ac,Bc,Cc,Tc=ctrbf(A,B,C)Ac,Bc,Cc,Tc=ctrbf(A,B,C)第14页/共109页第十四页，共109页。153.3.可控标准型及其可控标准型及其MATLABMATLAB实现实现 若系统完全可控，则可利用若系统完全可控，则可利用(lyng)(lyng)矩阵矩阵TcTc将其变将其变换为第一可控规范型换为第一可控规范型 ，其系数阵之间，其系数阵之间满足关系满足关系其中(qzhng)，第15页/共109页第十五页，共109页。16例例8 8：已知系统：已知系统(xt(xt ng)ng)的系数阵为的系数阵为试判断其可控性。若完全可控，则求其第一(dy)可控规范型。A=2 0 0;0 4 1;0 0 4;B=1;0;1;C=1 1 0;Tc=ctrb(A,B);rank(Tc)Ac=inv(Tc)*A*Tc;Bc=inv(Tc)*B;Cc=C*Tc;Ac=0 0 32 1 0 -32 0 1 10Bc=1 0 0Cc=1 3 12结果(ji gu)第16页/共109页第十六页，共109页。174.1.4 4.1.4 线性系统的可观线性系统的可观线性系统的可观线性系统的可观(kgun)(kgun)性分析性分析性分析性分析 如果系统所有的状态变量任意形式的运动均可由有限测量点的输出(shch)完全确定出来，则称系统可观测；如果系统中所有的状态都是可观测的，则称该系统为完全可观测的系统。反之，则称系统不可观测。1.线性系统的可观测性判定可通过构造可观测可通过构造可观测(gunc)(gunc)性判定矩阵性判定矩阵若若ToTo为满秩矩阵，则系统为完全可观测为满秩矩阵，则系统为完全可观测(gunc)(gunc)的。的。可观测可观测(gunc)(gunc)性判定矩阵由性判定矩阵由To=obsv(A,C)To=obsv(A,C)函数构造。函数构造。rank()rank()函数可求出矩阵的秩。函数可求出矩阵的秩。第17页/共109页第十七页，共109页。18A=4 4 4;-11-12-12;13 14 13;A=4 4 4;-11-12-12;13 14 13;BB=1;-1;0;=1;-1;0;C=1 1 1;C=1 1 1;To=obsv(A,C)To=obsv(A,C)rank(To)rank(To)例9：已知系统(xtng)的系数阵为试判断(pndun)它的可观性。To=1 1 1 6 6 5 23 22 17rank(To)=3故系统(xtng)可观。结果第18页/共109页第十八页，共109页。192.2.可观标准型及其可观标准型及其MATLABMATLAB实现实现 若系统完全可观，则可利用矩阵若系统完全可观，则可利用矩阵P P将其变换将其变换(binhun)(binhun)为第一可观规范型，其系数阵之间满足关为第一可观规范型，其系数阵之间满足关系系其中(qzhng)，第19页/共109页第十九页，共109页。20例：上例中，求其第一例：上例中，求其第一(dy)(dy)可控规范型。可控规范型。P=inv(P=inv(ToTo););Ao=inv(P)*A*P;Ao=inv(P)*A*P;Bo=inv(P)*B;Bo=inv(P)*B;Co=C*P;Co=C*P;Ao=-0.0000 1.0000 0.0000 0.0000 -0.0000 1.0000 4.0000 -8.0000 5.0000Bo=0.0000 0.0000 1.0000Co=1.0000 0.0000 -0.0000结果(ji gu)第20页/共109页第二十页，共109页。213.3.系统的可观系统的可观(kgun)(kgun)性分解性分解 对于状态不完全可观对于状态不完全可观(kgun)(kgun)的系统，同样可对其的系统，同样可对其进行可观进行可观(kgun)(kgun)性分解。在性分解。在MATLABMATLAB中可调用中可调用obsvfobsvf()()函数。函数。AO,BO,CO,T,K=obsvf(A,B,C)AO,BO,CO,T,K=obsvf(A,B,C)其中，其中，T T为相似变换阵，为相似变换阵，K K为可观为可观(kgun)(kgun)子阵的阶次子阵的阶次向量。对系统进行可观向量。对系统进行可观(kgun)(kgun)性分解后得到相应性分解后得到相应可观可观(kgun)(kgun)子系统子系统其中，其中，(Ao,Co)(Ao,Co)为可观为可观(kgun)(kgun)子对。子对。第21页/共109页第二十一页，共109页。22例例1010：试确定：试确定(qudng)(qudng)系统的可观性并进行可观性分解。系统的可观性并进行可观性分解。A=-2 1;1-2;B=1;0;C=1-1;n=length(A);Q=C;C*A;r=rank(Q);if r=n disp(system is observable.)Qelse disp(system is unobservable.)disp(rank)rend system is unobservable.rankr=1结果(ji gu)第22页/共109页第二十二页，共109页。23AO,BO,CO,T,K=obsvf(A,B,C)AO,BO,CO,T,K=obsvf(A,B,C)disp(observable submatrix)disp(observable submatrix)AO=AO(2,2)AO=AO(2,2)CO=CO(1,2)CO=CO(1,2)AO=-1.0000 0 -0.0000 -3.0000BO=0.7071 -0.7071CO=0 -1.4142T=0.7071 0.7071 -0.7071 0.7071K=1 0observable submatrixAO=CO=-3.0000 -1.4142结果(ji gu)：第23页/共109页第二十三页，共109页。244.2 控制系统(kn zh x tn)的时域分析 系统的过渡过程性能（如上升时间、调节时间、超调量及稳态误差）常用典型输入作用下的时间响应来描述(mio sh)。响应是指零初始值条件下某种典型的输入函数作用下对象的响应，其响应是时间t的函数，故称为时域响应。控制系统常用的输入函数为单位阶跃函数和脉冲激励函数（即冲激函数）。1.单位(dnwi)阶跃函数1(t)X(t)1t第24页/共109页第二十四页，共109页。252.单位冲激函数 t0X(t)t 斜坡(xip)函数B=1时称为单位斜坡(xip)函数。其拉氏变换(binhun)后的像函数为：第25页/共109页第二十五页，共109页。264.2.1 4.2.1 线性定常系统线性定常系统线性定常系统线性定常系统(xt(xt ng)ng)状态方程的解状态方程的解状态方程的解状态方程的解x(0)=x0 x(0)=x0是系统的初始状态。是系统的初始状态。问题：对给定的控制输入和初始状态，如何确定任意时刻的系统状问题：对给定的控制输入和初始状态，如何确定任意时刻的系统状态和输出态和输出(shch)(shch)；状态的变化行为。；状态的变化行为。首先考虑首先考虑两边(lingbin)同时左乘e-At,得根据矩阵微积分知识，上式进一步有：(4.1)(4.2)(4.3)(4.4)第26页/共109页第二十六页，共109页。27两边同时在t0,t区间(q jin)积分，得两边同时(tngsh)左乘eAt,并整理得即：当初始(ch sh)时刻t0=0时，上式可变为由此可知，非齐次状态方程的解由两部分组成，第一部分是在初始状态X(t0)作用下的自由运动，第二部分为在系统输入U(t)的作用下的强制运动。(4.5)(4.6)(4.8)(4.9)(4.10)第27页/共109页第二十七页，共109页。28当U(t)为几种典型的控制输入时，则有如下(rxi)形式。1.冲激函数输入，即即：(4.11)第28页/共109页第二十八页，共109页。292.阶跃信号(xnho)输入，即3.斜坡(xip)信号输入，即U(t)=kt，可以求得：(4.12)(4.13)第29页/共109页第二十九页，共109页。30例11：求以下(yxi)系统在单位阶跃函数作用下的状态响应由或得则(4.14)(4.15)第30页/共109页第三十页，共109页。31grid;t=0:0.1:10;x1=0.5-exp(-t)+0.5*exp(-2*t);x2=exp(-t)-exp(-2*t);plot(t,x1,x,t,x2,*);xlabel(时间轴)ylabel(x代表(dibio)x1，-*代表(dibio)x2)系统(xtng)状态轨迹图第31页/共109页第三十一页，共109页。32另：另：另：另：当输入当输入当输入当输入(shr)(shr)信号为信号为信号为信号为可将系统(xtng)原始的状态方程转化为其中(qzhng)，(4.16)(4.17)(4.18)第32页/共109页第三十二页，共109页。33则，系统(xtng)状态方程的解析解为以下为MATLAB编写的ss_augment()，用来求取(qi q)系统的增广状态方程模型。其中，cc=c0,c1,ck，dd=d1,d2,d3,d4(4.19)(4.20)(4.21)第33页/共109页第三十三页，共109页。34function Ga,Xa=ss_augment(G,cc,dd,X)G=ss(G);Aa=G.a;Ca=G.c;Xa=X;Ba=G.b;D=G.d;if(length(dd)0&sum(abs(dd)1e-5),if(abs(dd(4)1e-5),Aa=Aa dd(2)*Ba,dd(3)*Ba;.zeros(2,length(Aa),dd(1),-dd(4);dd(4),dd(1);Ca=Ca dd(2)*D dd(3)*D;Xa=Xa;1;0;Ba=Ba;0;0;else Aa=Aa dd(2)*B;zeros(1,length(Aa)dd(1);Ca=Ca dd(2)*D;Xa=Xa;1;Ba=B;0;endend第34页/共109页第三十四页，共109页。35if(length(cc)0&sum(abs(cc)1e-5),if(length(cc)0&sum(abs(cc)1e-5),M=length(cc);M=length(cc);Aa=Aa Ba zeros(length(Aa),M-1);Aa=Aa Ba zeros(length(Aa),M-1);zeros(M-1,length(Aa)+1)eye(M-1);zeros(M-1,length(Aa)+1)eye(M-1);zeros(1,length(Aa)+M);zeros(1,length(Aa)+M);Ca=Ca D zeros(1,M-1);Xa=Xa;cc(1);ii=1;Ca=Ca D zeros(1,M-1);Xa=Xa;cc(1);ii=1;for i=2:M for i=2:M ii=ii*i;%ii=ii*i;%阶乘阶乘(ji chn(ji chn)运算运算 Xa(length(Aa)+i)=cc(i)*ii;Xa(length(Aa)+i)=cc(i)*ii;end;end;endendGa=ss(Aa,zeros(size(Ca),Ca,D);Ga=ss(Aa,zeros(size(Ca),Ca,D);第35页/共109页第三十五页，共109页。36例12:系统(xtng)的状态方程模型为cc=2;dd=-4 0 2 3;A=-5 2 0 0;0-4 0 0;-3 2-4-1;-3 2 0-4;x0=1;2;0;1;B=1;2;3;4;C=1 1 1 1;D=0;G=ss(A,B,C,D);Ga,xx0=ss_augment(G,cc,dd,x0)其中状态变量初值xT(0)=1,2,0,1。假设系统的输入信号(xnho)u(t)=2+2e-4tsin(3t)，则由ss_augment()得系统的增广状态方程则由下面语句可直接(zhji)获得生成信号的解析解syms t;y=Ga.c*expm(Ga.a*t)*xx0;latex(y)第36页/共109页第三十六页，共109页。374.2.2 基于部分(b fen)分式展开方法求解 控制系统常用到并联系统，这时就要对系统函数进行(jnxng)分解，使其表现为一些基本控制单元的和的形式。如若系统(xtng)输出的Laplace 式子的分母多项式的根pi都是不重复的，则若其中的j个根为m重根pj，则该部分的Laplace 反变换为(4.22)(4.23)第37页/共109页第三十七页，共109页。38部分部分(b fen)(b fen)分式展开：分式展开：num=2,0,9,1;num=2,0,9,1;den=1,1,4,4;r,p,k=residue(num,den)den=1,1,4,4;r,p,k=residue(num,den)p=0.0000+2.0000i 0.0000-2.0000i -1.0000k=2r=0.0000-0.2500i 0.0000+0.2500i -2.0000结果(ji gu)表达式：第38页/共109页第三十八页，共109页。39 MATLAB环境中，函数r,p,k=residue(num,den)对两个(lin)多项式进行部分展开，以及把传函分解为微分单元的形式。向量num和den是按s的降幂排列的多项式系数。部分分式展开后，k为展开余项。例13:考虑系统的传函为系统(xtng)的输入信号为单位阶跃信号，则系统(xtng)的输出信号为将输出信号用部分分式(fnsh)展开,则有第39页/共109页第三十九页，共109页。40 num=1 7 24 24;den=1 10 35 50 24 0;r,p,k=residue(num,den);r,pans=-1.0000 -4.0000 2.0000 -3.0000 -1.0000 -2.0000 -1.0000 -1.0000 1.0000 0故,由上面结果(ji gu)写出系统输入表达式为此例为系统传函只含有实数极点的情况(qngkung)，下面的例子中系统输出还含有复数极点。第40页/共109页第四十页，共109页。41例14:考虑(kol)系统的传函为 num=1 3;den=1 2 11 18 18 0;r,p,k=residue(num,den);r,pans=0.0020+0.0255i 0.0000+3.0000i 0.0020-0.0255i 0.0000-3.0000i -0.0853+0.0088i -1.0000+1.0000i -0.0853-0.0088i -1.0000-1.0000i 0.1667 0 故第41页/共109页第四十一页，共109页。42若想获得输出(shch)的物理表达形式，则可借助(jizh)函数pfrac()function R,P,K=pfrac(num,den)R,P,K=residue(num,den);for i=1:length(R),if abs(imag(P(i)eps a=real(R(i);b=imag(R(i);R(i)=-2*sqrt(a2+b2);R(i+1)=-atan2(a,b);elseif abs(imag(P(i)eps R(i)=real(R(i);endend其中(qzhng)，第42页/共109页第四十二页，共109页。43 离散系统的解析(ji x)解离散系统的输出信号为则若考虑采样周期T的因素，则系统输入信号(xnho)的解析解为(4.24)(4.25)(4.26)第43页/共109页第四十三页，共109页。44D=conv(1-1/3,conv(1-1/4,conv(1 1/5,1-1);%分母N=0 0 conv(1-1/2,1 0);%分子前面(qin mian)要补足零N=N(end:-1:1);D=D(end:-1:1);%将z变换式逆序排列R,P,K=residue(N,D);R,P,-R./P例15：离散系统的传函为假设系统的输入为阶跃信号，其z变换为z/(z-1)，则其输出(shch)可有下列语句计算ans=12.1528 -5.0000 2.4306 35.5556 4.0000 -8.8889 -16.8750 3.0000 5.6250 -0.8333 1.0000 0.8333可写出第44页/共109页第四十四页，共109页。45下例为系统的输出(shch)含有重根的情况。例16：离散系统的传函为则其阶跃响应(xingyng)的解析解为D=conv(1-1/2,conv(1-1/2,conv(1-1/2,.conv(1,-1/3,1-1);N=0 0 0 5-2 0;R,P,K=residue(N(end:-1:1),D(end:-1:1);R,Pans=324.0000 3.0000-240.0000 2.0000 -96.0000 2.0000 192.0000 2.0000 -36.0000 1.0000于是(ysh)有第45页/共109页第四十五页，共109页。46 时间延迟系统的解析解 对于带有时间延迟的系统G(s)e-Ls和H(z)z-k，若将其直接进行部分(b fen)分式展开比较困难。可先不考虑时间延迟，得出系统的解析解y(t)和y(n)，然后分别用t-L或n-k代替得出的解析解中的t或n，即可得到时间延迟系统的解析解。例17：带有时间延迟的系统将上例的解的n替换(t hun)为n-5，即可第46页/共109页第四十六页，共109页。474.2.3 二阶系统的阶跃响应(xingyng)及阶跃响应(xingyng)指标 假设系统的开环模型为 ，由单位(dnwi)负反馈构造出的闭环系统模型为于是(ysh)，二阶系统的阶跃响应y(t)为其中，根据的取值不同，系统输入也不同。第47页/共109页第四十七页，共109页。48s2s1 当=0，y(t)=1-cos(nt)为无阻尼振荡。当011，称为(chn wi)过阻尼振荡。1S1,2S1S21第49页/共109页第四十九页，共109页。50wn=1;yy=;t=0:.1:12;zet=0:0.1:1,2,3,5;for z=zet if z=0,y=1-cos(wn*t);elseif(z0&z1,dd=sqrt(z2-1);lam1=-z-dd;lam2=-z+dd;y=1-0.5*wn*(exp(lam1*t)/lam1-exp(lam2*t)/lam2)/dd;end yy=yy;y;endplot(t,yy)当n=1rad/sec，选择不同(b tn)的，由下列命令可以得到系统在不同(b tn)阻尼比下的阶跃响应曲线。第50页/共109页第五十页，共109页。51不同(b tn)阻尼比下的阶跃响应第51页/共109页第五十一页，共109页。52 稳态值稳态值 线性系统典型的阶跃响应(xingyng)曲线，如下图所示误差带：误差带：5，2（1）延迟时间）延迟时间（2）上升时间）上升时间（3）峰值时间）峰值时间（4）调节时间调节时间调节时间调节时间其中(qzhng)，系统的稳态值可由dcgain()函数获得。（5）超调量超调量超调量超调量第52页/共109页第五十二页，共109页。534.3 4.3 线性系统的数字线性系统的数字线性系统的数字线性系统的数字(shz)(shz)仿真分析仿真分析仿真分析仿真分析4.3.1 线性系统的时域响应(xingyng)1、求取(qi q)系统单位阶跃响应 step()格式1：Y,T=step(G)格式2：Y,T=step(G,t)格式3：Y,T=step(G,iu)格式4：Y,T=step(G,iu,t)说明：G为tf(),zpk(),ss()中任一种模型。如果用户在调用step()函数时不返回任何向量，则将自动地绘出阶跃响应输出曲线。对于带有返回参数的将不绘制曲线,其中Y是输出向量,T是时间向量。t为用户设定的时间向量,一般可以由t=0:step:end等步长地产生出来。对于MIMO系统,iu表示第iu个输入到所有输出的冲激响应曲线.第53页/共109页第五十三页，共109页。54num=7 7;den=conv(conv(1 0,1 3),1 4 5);g=tf(num,den);gg=feedback(g,1,-1);y,t,x=step(gg);plot(t,y)例17 已知系统(xtng)框图其中(qzhng)，输入(shr)以下MATLAB命令，可得系统的单位阶跃响应曲线如下图所示。第54页/共109页第五十四页，共109页。55图4.1 阶跃响应(xingyng)指标显示第55页/共109页第五十五页，共109页。56例19：设连续(linx)系统的数学模型为选择采用周期T=0.1s，试比较各种离散方法的时域响应(xingyng)的差异。G=tf(1,1 3 3 1,iodelay,.5);G1=c2d(G,0.1,zoh);G2=c2d(G,0.1,foh);G3=c2d(G,0.1,tustin);step(G,-r,G1,-g,G2,:b,G3,-.y)第56页/共109页第五十六页，共109页。57第57页/共109页第五十七页，共109页。582、impulse 求连续系统的单位脉冲响应。格式(g shi)1：Y,T=impulse(G)格式(g shi)2：Y,T=impulse(G,t)格式(g shi)3：Y,T=impulse(G,iu)格式(g shi)4：Y,T=impulse(G,iu,t)说 明：impulse()函数与step()函数调用格式(g shi)完全一致。例20:系统(xtng)则其脉冲响应曲线(qxin)由以下语句获得G=zpk(-1;-2;-3,-1+i;-1-i;-3.5;-4;-5,8,iodelay,2);impulse(G,8)第58页/共109页第五十八页，共109页。59第59页/共109页第五十九页，共109页。604.3.2 4.3.2 任意输入任意输入任意输入任意输入(shr)(shr)下系统的响应下系统的响应下系统的响应下系统的响应lsim()函数求任意输入信号时系统的响应(xingyng)格式1：Y,T=lsim(sys1,u,t)格式2：Y,T=lsim(sys2,u,t,x0)说明：u为输入信号.t为等间隔时间向量.sys1为tf()或zpk()模型。sys2为ss()模型,x0为初始条件。例21：系统(xtng)其两路输入分别为第60页/共109页第六十页，共109页。61MATLAB语句(yj)为：g11=tf(0.1134,1.78,4.48,1);g12=tf(0.924,2.07,1);g21=tf(0.3378,0.361,1.09,1);g22=tf(-0.318,2.93,1);G=g11,g12;g21,g22;G.ioDelay=0.72,0;0.3,1.29;t=0:.1:15;u=1-exp(-t).*sin(3*t+1),sin(t).*cos(t+2)lsim(G,u,t);第61页/共109页第六十一页，共109页。62在两路输入下，系统的时域响应(xingyng)曲线为：第62页/共109页第六十二页，共109页。63系统(xtng)的开环传函为综合(zngh)应用（1）求单位负反馈时，其闭环传递函数，并绘制(huzh)输出量阶跃响应曲线和脉冲响应曲线。（2）选择函数的状态变量将其化为状态方程模型，并绘制(huzh)状态变量的阶跃响应曲线和脉冲响应曲线。close all;clear;numo=0 0 0 0 200;deno=1 20 140 400 384;%系统开环传函numc=numo;n=length(deno);denc=zeros(1,n);denc=numo+deno;%系统闭环传函disp(System Closed Loop Transfer Function is:)numcdenc第63页/共109页第六十三页，共109页。64t=0:.05:3;y=step(numc,denc,t);%系统(xtng)输出量的阶跃响应yy=impulse(numc,denc,t);%系统(xtng)输出量的脉冲响应plot(t,y);title(System Step Response);xlabel(Time-sec);ylabel(Response-value);grid;pauseplot(t,yy);title(system impulse respone);xlabel(time-sec);ylabel(respone-value);grid;pauseA,B,C,D=tf2ss(numc,denc);disp(system state-space model is:);第64页/共109页第六十四页，共109页。65A，B，C，Dys,x=step(A,B,C,D,1,t);%系统(xtng)状态变量的阶跃响应yys,xx=impulse(A,B,C,D,1,t);%系统(xtng)状态变量的阶跃响应plot(t,x(:,1),t,x(:,2),+,t,x(:,3),.-,t,x(:,4),.)title(system state-variables step response);xlabel(time-sec);ylabel(response-value);grid;pauseplot(t,xx(:,1),t,xx(:,2),+,t,xx(:,3),.-,t,xx(:,4),.)title(system state-variables impulse response);xlabel(time-sec);ylabel(response-value);grid;第65页/共109页第六十五页，共109页。66Notice:仿真时间t的选择：对于典型二阶系统根据其响应时间的估算公式 可以确定。对于高阶系统往往其响应时间很难估计，一般采用试探的方法，把t选大一些，看看响应曲线的结果，最后再确定其合适的仿真时间。一般来说，先不指定仿真时间，由MATLAB自己确定，然后根据结果，最后确定合适的仿真时间。在指定仿真时间时，步长的不同会影响到输出曲线的光滑(gung hu)程度，一般不易取太大。第66页/共109页第六十六页，共109页。674.4 控制系统(kn zh x tn)的根轨迹分析方法 所谓根轨迹(Root Locus)是指，当系统的某个参数(cnsh)(如开环增益K)由零连续变化到无穷大时,闭环特征根在复平面上运动的轨迹。一般来说，在无零极点对消时，闭环系统特征根就是闭环传递函数的极点。根轨迹法(Root Locus Technique)是分析和设计线性定常控制系统的图解方法，使用十分简便。利用它可以对系统进行各种性能分析。一、根轨迹(guj)分析方法的概念第67页/共109页第六十七页，共109页。68图 4-1 控制系统(kn zh x tn)框图 将图4-1所示系统(xtng)的开环传递函数转化为 其中，k=2K，式(4.27)便是绘制根轨迹(guj)所用的开环传递函数的标准形式。(4.27)下面结合图4-1所示的二阶系统的例子,介绍有关根轨迹的基本概念。第68页/共109页第六十八页，共109页。69 由式(4.27)可得两开环极点分别(fnbi)为 p1=0,p2=-2，无开环零点。将这两个开环极点绘于图4-2上，并用“”表示。即 由式(4.27)可得闭环系统(xtng)的特征方程为 得闭环系统的特征(tzhng)根(闭环极点)为 从上式看出，闭环系统极点s1,s2与标准化参数k存在一定的关系，参数k变化会引起