现代控制理论实验指导书

实用文档标准文案现代控制理论实验指导书实验一：线性系统状态空间分析1、模型转换图1、模型转换示意图及所用命令传递函数一般形式：MATLAB表示为：G=tf(num,den)，其中num,den分别是上式中分子，分母系数矩阵。零极点形式：MATLAB表示为：G=zpk(Z,P,K)，其中Z，P，K分别表示上式中的零点矩阵，极点矩阵和增益。传递函数向状态空间转换：[A,B,C,D]=TF2SS(NUM,DEN)；状态空间转换向传递函数：[NUM,DEN]=SS2TF(A,B,C,D,iu)---iu表示对系统的第iu个输入量求传递函数；对单输入iu为1;验证教材P438页的例9-6。求P512的9-6题的状态空间描述。>>A=[01;0-2];>>B=[10;01];>>C=[10;01];>>D=[00;00];>>[NUM,DEN]=ss2tf(A,B,C,D,1)NUM=012000DEN=120>>[NUM,DEN]=ss2tf(A,B,C,D,2)NUM=001010DEN=120给出的结果是正确的，是没有约分过的形式P5129-6>>[A,B,C,D]=tf2ss([168],[143])A=-4-310B=10C=25D=12、状态方程求解单位阶跃输入作用下的状态响应：G=ss(A,B,C,D);[y,t,x]=step(G);plot(t,x).零输入响应[y,t,x]=initial(G,x0)其中，x0为状态初值。验证P435的例9-4，P437的例9-5。9-4A=[01;-2-3];B=[0;0];C=[00];D=[0];G=ss(A,B,C,D);[y,t,x]=initial(G,[1;2]);plot(t,x)（设初始状态为[1;2]）零输入响应9-5零输入响应A=[01;-2-3];B=[0;1];C=[00];D=[0];G=ss(A,B,C,D);[y,t,x]=initial(G,[1;2]);plot(t,x)零状态响应，阶跃信号激励下>>A=[01;-2-3];B=[0;1];C=[00];D=[0];>>G=ss(A,B,C,D);[y,t,x]=step(G);plot(t,x)总响应>>A=[01;-2-3];B=[0;1];C=[00];D=[0];G=ss(A,B,C,D);[y1,t1,x1]=step(G);[y2,t2,x2]=initial(G,[1;2]);>>x=x1+x2;>>plot(t1,x)3、系统可控性和可观测性可控性判断：首先求可控性矩阵：co=ctrb(A，B)。然后求rank(co)并比较与A的行数n的大小，若小于n则不可控，等于为可控。也可以求co的绝对值，不等于0，系统可控，否则不可控。验证P456例9-14。A=[01;-1-2];B=[1;-1];>>C=[10];>>co=ctrb(A,B)co=1-1-11>>rank(co)ans=1系统不可控可观测性判断：首先求可观测性矩阵ob=obsv(A，C),或者ob=ctrb(A'，C');然后求rank(ob)并比较与A的行数大小，若小于，为不可观测，等于则为可观测。验证P458例9-15。1>>A=[-20;0-1];B=[3;1];C=[10];>>ob=obsv(A,C)ob=10-20>>rank(ob)ans=1不可观测2、>>A=[1-1;11];B=[2-1;10];C=[10;-11];>>ob=obsv(A,C)ob=10-111-102>>rank(ob)ans=2可观测4、线性变换一个系统可以选用不同的状态变量，所以状态方程是不唯一的。但是这些方程之间是可以相互转换的。[At，Bt，Ct，Dt]=ss2ss(A，B，C，D，T)变换矩阵T不同，可得到不同的状态方程表示形式，如可控型，可观测型，Jordan标准型表示。matlab变换与控制书上讲的变换略有差别。这里是，其中x是原来的变量，z是现在的变量。书上则是。因此线性变换时，首先要对给定的变换矩阵进行逆变换，然后将其代入上面指令的T中。求对角阵（或约当阵）：MATLAB提供指令：[At，Bt，Ct，Dt，T]=canon(A，B，C，D，'modal')它可将系统完全对角化，不会出现经典控制中的约当块。>>A=[-4-3;10];B=[1;0];C=[25];D=[0];>>[At,Bt,Ct,Dt,T]=canon(A,B,C,D,'modal')At=-300-1Bt=-3.6056-2.2361Ct=-0.1387-0.6708Dt=0T=-3.6056-3.6056-2.2361-6.7082>>inv(T)ans=-0.41600.22360.1387-0.2236求可观测标准型：[At，Bt，Ct，Dt，T]=canon(A，B，C，D，'companion')>>[At,Bt,Ct,Dt,T]=canon(A,B,C,D,'companion')At=0-31-4Bt=10Ct=2-3Dt=0T=1401>>inv(T)ans=1-401求可控标准型：首先需要求可观测标准型，然后根据对偶关系求[At'，Ct'，Bt'，Dt']>>At'ans=01-3-4Ct'ans=2-3Bt'ans=10验证P512的9-6习题。5、线性定常系统的结构分解当系统是不可控的，可以进行可控性规范分解。使用[a1,b1,c1,t,k]=ctrbf(A,B,C)命令。验证P473例题9-19。当系统是不可观测的，可以进行可观测性规范分解。使用[a2,b2,c2,t,k]=obsvf(A,B,C)命令。验证P475例题9-20。>>A=[12-1;010;1-43];B=[0;0;1];C=[1-11];>>co=ctrb(A,B)co=0-1-4000138>>rank(co)ans=2系统不可控[a1,b1,c1,t,k]=ctrbf(A,B,C)a1=100-211-4-13b1=001c1=-1-11t=010-100001k=110>>ob=obsv(A,C)ob=1-112-324-74>>rank(ob)ans=2系统不可观测>>[a2,b2,c2,t,k]=obsvf(A,B,C)a2=2.00002.30944.0825-0.00000.6667-0.9428-0.00000.47142.3333b2=0.7071-0.40820.5774c2=-0.0000-0.00001.7321t=-0.7071-0.00000.7071-0.4082-0.8165-0.40820.5774-0.57740.5774k=1106、极点配置算法调用命令格式为K=place(A,B,P)。A，B为系统系数矩阵，P为配置极点，K为反馈增益矩阵。验证P484例9-21。>>A=[000;1-60;01-12];B=[1;0;0];P=[-2-1+j-1-j];>>K=place(A,B,P)K=1.0e+003*-0.01400.1860-1.2200用下列编码对状态反馈前后的输出响应进行比较（附带文件control.m）。t=0:0.01:5;U=0.025*ones(size(t));%幅值为0.025输入阶跃信号[Y1,X1]=lsim(A,B,C,D,U,t);[Y2,X2]=lsim(A-B*K,B,C,D,U,t);figure(1)plot(t,Y1);grid;title('反馈前');figure(2)plot(t,Y2);title('反馈后');grid;>>A=[000;1-60;01-12];B=[1;0;0];P=[-2-1+j-1-j];C=[100];D=[0];K=1.0e+003*[-0.01400.1860-1.2200];t=0:0.01:5;U=0.025*ones(size(t));%幅值为0.025输入阶跃信号[Y1,X1]=lsim(A,B,C,D,U,t);[Y2,X2]=lsim(A-B*K,B,C,D,U,t);figure(1)plot(t,Y1);grid;title('反馈前');figure(2)plot(t,Y2);title('反馈后');grid;7、线性定常系统稳定判据函数lyap(A,Q)求如下式的李氏方程：AP+PAT=-Q注意与教材的区别，应将给定A矩阵转置后再代入lyap函数。验证P495的例9-25。>>A=[02;1-1];>>Q=[10;01];>>lyap(A,Q)ans=-0.7500-0.2500-0.25000.2500实验二：单级倒立摆的LQR状态调节器设计有一个倒立摆小车系统如图一所示。它由质量为M的小车，长为2L的倒立摆构成，倒立摆的质量为m，铰链在小车上，小车在控制函数u的作用下，沿滑轨在x方向运动，使倒立摆在垂直平面内稳定。2L2LxMmu图一小车倒立摆系统示意图为了简单起见，设倒立摆为均匀细杆，执行机构和轴无摩擦，此时系统的动力学非线性微分方程为：其中，M=1kg，m=0.1kg，L=1m，g=9.81m/s2，f=50N/s。现设计LQR控制器，使系统倒立摆在初始条件[x(0),dx(0),(0),d(0)]T=[0.05,0,0.08,0]T下稳定。注：倒立摆偏角可以通过同轴旋转电位计送出正比的电信号。思路：先建立状态空间模型A、B、C和D矩阵，当摆杆角度很小时，可令，四个状态为[x(0),dx(0),(0),d(0)]T。采用LQR命令求最优K矩阵，定义Q和R阵，用对角阵，人工定义对角线上的加权系数，由此决定对每个状态分量的影响。如对角度侧重，则对应系数加大。>>A=[0100;0-4.972500;0001;03.729400]>>B=[0;0.9756;0;-0.7317];>>C=[1000;0010];>>D=[0;0];>>Q=[100000;0000;00100;0000];>>R=1;>>[K,P]=lqr(A,B,Q,R)K=9.8312369.6931-0.5786490.5742P=1.0e+012*0.00000.00010.00000.00010.00010.73370.00020.97830.00000.0002