控制系统建模与仿真-基于MATLABSimulink的分析与实现（第2版）课件 第8讲MATLAB在状态空间分析中的应用

第8讲MATLAB在状态空间分析中的应用主讲内容一、

极点配置的概念二、系统的可控性与可观测性三、极点配置与状态反馈四、最优二次型设计一、

极点配置的概念1、什么是极点配置

对线性定常系统，系统的稳定性和各种性能指标，主要是由闭环系统的极点位置所决定的。极点配置方法就是把系统的闭环极点配置到希望的极点位置上，使得闭环系统的极点位于s平面上所期望的值，从而获得良好的性能指标。线性系统的动态性能，如系统稳定性，时域分析中的超调量、过渡时间等指标，主要取决于系统的极点位置。因此，作为系统性能指标的一种形式，往往是给定一组期望极点，或者根据时域指标转换成一组等价的期望极点。极点配置问题，就是通过选择反馈增益矩阵，将闭环系统的极点恰好配置在根平面上所期望的位置，以获得所希望的动态性能。对于SISO系统来说，极点配置问题描述如下：设控制系统的状态方程模型为：2、状态空间模型式中：u—系统的控制输入向量；X—系统的状态向量；Y—系统的输出向量A—控制对象状态矩阵B—控制矩阵；C—输出矩阵；当引入状态反馈后，系统的控制信号为，这里是系统外部的输入；是行向量，称为反馈矩阵；此时闭环系统状态方程模型为：

加入状态反馈增益阵K=[k0,k1…..kn-1],可求得对的闭环状态方程表达式

状态反馈系数矩阵二、系统的可控性与可观测性1、可控性系统可控性指的是控制作用对被控系统的状态和输出进行控制的可能性。2、可观测性可观测性反映由能直接测量的输入输出的量测值来确定反映系统内部动态特性的状态的可能性。状态状态状态可控性反映输入u(t)对状态x(t)的控制能力。状态可观测性反映系统外部可直接或间接测量的输出y(t)来确定或反映系统状态的能力。三、极点配置与状态反馈1、极点配置必要条件控制系统中对动态性能的要求可通过设计一组希望的闭环极点（主导极点）来实现。如果把一个n阶系统中的n个状态变量作为系统的反馈信号，在系统可控的条件下，就能实现对系统极点的任意配置。

由于线性定常系统的特征多项式为实系数多项式，因此考虑到问题的可解性，对期望的极点的选择应注意下列问题:（1）对于

n阶系统，必须给出n个期望的极点；（2）期望极点必须是实数或成对出现的共轭复数；（3）期望极点必须体现对闭环系统性能品质指标要求完全可控的充分必要条件：状态阵A满秩，即：对矩阵A的所有的特征值λi=[i=1,2,3……n]

有

其中:n为矩阵

A的维数。

也称为系统的可控性判别阵。

2、用MATLAB进行极点配置的步骤（1）检查系统的可控性，当时系统可控；在MATLAB中可用B=rank(ctrb(A,b))实现；（2）由希望的闭环极点得到希望的闭环特征方程，确定系统矩阵A的特征多项式系数，在MATLAB中，可用poly函数实现；可控性判别Matlab函数如下：functionm=controllble(A,B)ctrl=rank(ctrb(A,B));n=length(A);ifn==ctrldisp(‘系统可控‘)else

disp(‘系统不可控‘)end（3）确定变换矩阵T,T=B*w,B为能控性判别矩阵，在MATLAB中可用ctrb(a,b)实现，w在MATLAB中可用hankel()函数来实现。（4）用引入反馈系数后的特征方程与希望的特征方程对比,确定期望特征多项式系数（5）求增益矩阵K的系数：（6）确定输入变换器L状态反馈极点配置：通过状态反馈矩阵K的选取，使闭环系统的极点，即A-BK的特征值恰好处于所希望的一组给定闭环极点的位置上。

此时，状态反馈律u=-Kx+v/L作用下的闭环系统如图所示：

系统状态反馈后，需要输入变换器1/L进行配准，配准后的控制阵B1=L.*B，L取值为极点配置（状态阵A1=A-BK）后传递函数s=0的值，以消除极点配置后的稳态误差，即：L的取值是加入状态反馈系数矩阵后的传递函数转换成多项式传递函数，再令s为零，即：程序命令：[num1,den1]=ss2tf(A-B*K,B,C,D);L=polyval(den1,0)/polyval(num1,0)MATLAB可以用polyval()函数取函数某零点的值得到，相应的实现程序如下：（1）acker函数：功能：SISO极点配置格式：K=acker(a,b,p)%只适用于SISO系统说明：acker函数利用Ackermann公式计算反馈增益矩阵K，使采用全反馈的单输入系统具有指定的闭环极点，即（2）place函数功能：MIMO极点配置格式：K=place(A,B,p)[K,prec,message]=place(A,B,p)

3、acker函数和place函数

MATLAB工具箱提供了place和acker函数进行极点配置设计，可直接采用这两个函数进行设计。说明：place函数利用Ackermann公式计算反馈增益矩阵K，使采用全反馈的多输入多输出系统具有指定的闭环极点，即[K,prec,message]=place(A,B,p)同时返回系统闭环极点与希望极点p的接近程度prec。prec中的每个量的值为匹配的位数。如果系统闭环极点的实际极点偏离希望极点10%以上，则message将给出警告信息。观测器极点为重根时，提示出错信息，place()命令不适于进行极点为重根的极点配置，对于有重根的情况可直接用acker()命令计算。设系统的传递函数为：

试设计状态反馈控制器，使闭环系统的极点为

P=[-2，-10，-10]。例1：程序代码：（法一）

num=[10];den=conv([10],conv([11],[12]));[A,B,C,D]=tf2ss(num,den)Nctr=rank(ctrb(A,B));n=length(A);ifn==Nctr

disp('该系统是可控的');p=[-2-10-10];K=acker(A,B,p)%使用place函数出错，因为p有重根else

message('系统是不可控的，无法进行极点配置')endeig(A-B*K)结果该系统是可控的K=19138200ans=-10.0000-10.0000-2.0000num=[10];den=conv([10],conv([11],[12]));[A,B,C,D]=tf2ss(num,den)Nctr=rank(ctrb(A,B));n=length(A);ifn==Nctr

disp('该系统是可控的');fq=poly([-2,-10,-10])%计算配置极点的特征多项式thta=polyvalm(fq,A)k=[001]*inv(ctrb(A,B))*thatelse

message('系统是不可控的，无法进行极点配置')endeig(A-B*K)程序代码：（法二）已知系统的传递函数为：希望极点为，试设计状态反馈矩阵K。程序代码：num=[1];den=[118720];[A,B,C,D]=tf2ss(num,den);Nctr=rank(ctrb(A,B));n=length(A);ifn==Nctrdisp('该系统是可控的');p=[-2+2j-2-2j-10];K=place(A,B,p)elsemessage('系统是不可控的，无法进行极点配置')end例2：

结果：

该系统是可控的K=-4-2480ans=-10.0000+0.0000i-2.0000+2.0000i-2.0000-2.0000ieig(A-B*K)%计算特征值4、系统可观测系统能设计状态观测器的前提条件是系统必须可观，在设计状态观测器前必须判断系统的可观性，判断系统可观性程序如下:functionm=observable(a,c)obsvb=rank(obsv(a,c));n=size(a,1);ifn==obsvbdisp(['系统是可观测的！'])elsedisp(['系统是不可观测的！'])end若线性系统为：给定观测器极点为，计算观测器反馈矩阵G例3：A=[01;-2-3];B=[0;1];C=[20];D=0Nobsv=rank(obsv(A,C));n=length(A);ifn==Nobsvdisp('该系统是可观测的');p=[-10-10];G=acker(A',C',p)elsemessage('系统是不可观测的，无法进行极点配置')end程序代码（法一）结果G=8.500023.5000A=[01;-2-3];B=[0;1];C=[20];D=0Nobsv=rank(obsv(A,C));n=length(A);ifn==Nobsvdisp('该系统是可观测的');fq=poly([-10-10])%计算配置极点的特征多项式pthi=polyvalm(fq,A)G=[01]*inv(obsv(A‘,C))*pthielse？？？message('系统是不可观测的，无法进行极点配置')end程序代码（法二）由观测器的状态方程可以看到，观测器的系数矩阵的(A-GC)计算是状态反馈矩阵(A-BK)的转置计算，即观测器方程的对偶系统的状态反馈，因此MATLAB函数place()和acker()同样可以用于状态观测器反馈矩阵G的计算。结论：C=[100]y=cx

已知系统状态方程，若期望特征值为p=[-2+2j,-2-2j,-10]，判断系统是否可控？若完全可控，求状态增益矩阵K。判断系统是否可观测？若完全能观测，求观测矩阵G。例4：程序代码：A=[010;001;-1-5-6];B=[0;0;1];C=[100];Nctr=rank(ctrb(A,B));Nobsv=rank(obsv(A,C));n=length(A);ifn==Nctrdisp('该系统是可控的');p=[-2+2j-2-2j-10];K=acker(A,B,p)ifn==Nobsvdisp('该系统是可观测的');G=acker(A',C',p)elsedisp('该系统是不可观测的');endelsedisp('该系统是不可控的');disp('该系统也是不可观测的')end结果：该系统是可控的K=[79438]该系统是可观测的G=[8-569]设受控系统的传递函数为:用状态反馈将闭环系统极点配置为设计实现上述反馈的全维观测器，设其极点为(-10，-10)。例5：程序代码：num=[1];den=conv([10],[16]);[A,B,C,D]=tf2ss(num,den);Nctr=rank(ctrb(A,B));Nobsv=rank(obsv(A,C));n=length(A);ifNctr==n&Nobsv==ndisp('该系统是可控的，也是可观测的');p1=[-4+6i-4-6i];p2=[-10-10];K=acker(A,B,p1)G=acker(A',C',p2)elsemessage(‘系统不可控，无法进行极点配置和计算观测阵')end结果:该系统是可控的，也是可观测的K=252G=1614已知单倒立摆系统摆长为，重量为m的摆用铰链安装在质量为M的小车上，小车受执行电机操纵，在水平方向施加控制力u，相对参考系产生位移x，重力加速度为g，当摆角为θ时，保持保持在垂直位置，列出系统的状态空间表达式为：设M=1Kg,m=0.1Kg,L=1m,g=9.8m/s2,用状态反馈配置系统极点，运用全维状态观测器实现状态反馈。例6：程序代码：M=1;m=0.1;l=1;g=9.8A=[0100;00-(m*g/M)0;0001;00(M+m)*g/(M*l)0];B=[0;1/M;0;-(1/M*l)];C=[1000];Nctr=rank(ctrb(A,B));Nobsv=rank(obsv(A,C));n=length(A);ifNctr==n&Nobsv==ndisp('该系统是可控的，也是可观测的');p=[-1-2-1+j-1-j]p2=[-2-3-2+j-2-j]K=acker(A,B,p)G=acker(A',C',p2)'elsemessage('系统是不可控的，无法运用全维状态观测器')endK=-0.4082-1.0204-21.1882-6.0204G=9.000041.7800-149.0000-490.1922结果给定一个n阶线性控制对象，其状态方程是

寻求最优控制u(t)，使系统能够满足某种最优的要求，最简单地方法是引入线性二次型(LinearQuadratic，简称LQ)最优控制。四、最优二次型设计1、连续系统的LQ调节器设计格式：[K,P,E]=lqr(A,B,Q,R)[K,P,E]=lqr(A,B,Q,R,N)说明：A为系统的状态矩阵；B为系统的输出矩,阵；Q为给定的半正定实对称常数矩阵；R为给定的正定实对称常数矩阵；K为最优反馈增益矩阵；P为对应Riccati方程的唯一正定解。（若矩阵A-BK是稳定矩阵，则总有正定解P存在）；E为矩阵A-BK的闭环特征值。设计线性二次型最优控制器关键是选择加权矩阵Q。一般说来，Q选择的越大，系统达到稳态所需的时间越短，当然，还要实际的系统允许。基于最小值原理的线性二次型最优控制，通过求解代数Riccati方程，得到的状态反馈控制律K，可以使系统的各状态获得渐近稳定特性。它的不足之处在于，加权矩阵Q、R的值与系统响应性能之间的关系是定性的，往往不能一次得到满意的结果，需要多次调整它们的值得到满意的系统响应性能。已知系统状态方程为:R=1，求最优二次型解。例7：程序代码：A=[010;001;-35-27-9];B=[0;0;1];Q=eye(3);R=1;[k,p,e]=lqr(A,B,Q,R)

结果：k= 0.01430.11070.0676p= 4.26252.49570.01432.49572.81500.11070.01430.11070.0676e= -5.0958+0.0000i-1.9859+1.7110i-1.9859-1.7110i得出最优反馈增益矩阵k，闭环特征值e和Riccati方程的正定矩阵解。采用输入反馈，系统的性能指标为：R=1设计LQ最优控制器，计算最优状态反馈矩阵K=[k1,k2,k3]，并对闭环系统进行单位阶跃仿真。例8：A=[0,1,0;1,0,1;-5,-14,-19];B=[0,0,1]';C=[1,0,0];D=0;Q=diag([1,1,1]);R=1;G0=ss(A,B,C,D);p=eig(A)figure(1);step(G0)K=lqr(A,B,Q,R)k1=K(1);Ac=A-B*K;p1=eig(Ac)Bc=B*k1;G=ss(Ac,Bc,C,D);figure(2);step(G);程序代码：控制前后的单位阶跃响应曲线如图结论：从控制前系统的特征值看出，系统有个右极点，是不稳定的。经最优输出反馈后系统稳定，闭环系统阶跃响应曲线没有超调量,稳态时间也比较小。

某倒立单摆系统如图所示，其状态方程为：例9：其中：

小车的质量为M=0.5kg,摆长倒立单摆的质量m=0.2Kg

小车的摩擦系数b=0.1N/m

倒立单摆的惯量为I=0.006Kgm2

小车的坐标为x，倒立单摆的垂直角度、。使用最优的规则设计线性二次型最优控制器，要求：（1）、输出量x和的过渡过程时间小于2s。（2）、输出量的上升时间小于0.5s（3）、输出量的超调显小于20oclear;clc;M=0.5;m=0.2;b=0.1;i=0.006;g=9.8;l=0.3;p=i*(M+m)+M*m*l^2;A=[0100;0-(i+m*l^2)*b/p(m^2*g*l^2)/p0;0001;0-(m*l*b)/pm*g*l*(M+m)/p0];B=[0;(i+m*l^2)/p;0;m*l/p];C=[1000;0010];D=[0;0];p=eig(A)G0=ss(A,B,C,D);figure(1);step(G0)Q=[1000;0000;0010;0001];R=1;K=lqr(A,B,Q,R)Ac=A-B*K;p1=eig(Ac)G1=ss(Ac,B,C,D)；figure(2);step(G1)程序代码：p=05.5651-0.1428-5.6041结果系统系统的特征值有一个右半平面的极点，因此不稳定，必须加入校正装置。添加最优控制器后，系统的特征值：p1=-8.2782+0.0000i-3.8096+0.0000i-0.8428+0.8380i-0.8428-0.8380i系统稳定。结果：从系统的阶跃响应曲线中可以看到，系统超调量基本满足要求，但一方面系统的稳态值与期望相差很远(小车坐标的响应曲线稳态值为负值)，另一方面过渡过程时间和上升时间都很大，必须重新校正。方法就是加大加权矩阵Q的值。结论：MATLAB的工具箱中提供了Kalman()函数来求解系统的Kalman滤波器。其调用格式为：对于一个给定系统sys，噪声协方差Q，R，N函数返回一个Kalman滤波器的状态空间模型kest，其滤波器反馈增益为L，状态估计误差的协方差为P。

2、Kalman滤波器MATLAB实现

Kalman滤波器就是最优观测器，能够抑制或滤掉噪声对系统的干扰和影响。利用Kalman滤波器对系统进行最优控制是非常有效的。

已知系统的状态方程为：用Q=0.001、R=0.1，设计Kalman滤波器的增益矩阵与估计误差的协方差。例10：clc;A=[-1,0,1;1,0,0;-3,7,-2];B=[6,1,1]';C=[0,0,1];D=0;G=ss(A,B,C,D);Q=0.001;R=0.1;[kest,L,P]=kalman(S,Q,R)结果：L=1.01501.20561.8469P=0.06800.07220.