状态反馈的倒立摆PID控制

计算机控制理论与设计-1-基于LMI的单级倒立摆的状态反馈PID控制摘要本文以单级倒立摆为研究对象，通过物理规律得到系统的传递函数和状态方程，并结合状态反馈的概念，将约束条件转化为线性矩阵不等式求解反馈矩阵并得到PID控制器，最后进行MATLAB仿真得到仿真曲线进行对比分析指出该方法的优缺点。关键字单级倒立摆PIDMATLABLMI状态反馈前言倒立摆系统是一种典型的非线性的，不稳定的复杂系统。是控制理论教学与科研中研究诸如鲁棒问题、非线性系统的控制等问题的良好实验对象。同时，倒立摆系统作为机器人行走中平衡控制、火箭垂直姿态控制和卫星飞行中姿态控制的最简单模型在航空航天以及军工等领域有着广泛的用处。倒立摆可以根据摆杆数量的不同分为一级、二级和三级等，多级摆杆间采用自由连接。一级倒立摆的仿真与控制已广泛应用于教学科研，而二级倒立摆也已在大部分实验室中实现，至于三级倒立摆的控制问题则是国际上公认的难题。然而我国学者李洪兴教授在2002年实现了国际上首次四级倒立摆实物系统的控制，这是我国学者采用自己提出的理论完成世界性难题的重大科学成就。本文中以单级倒立摆为研究对象，根据物理定律进行建模得到数学模型，在此基础上进行PID控制，并通过MATLAB仿真对比加入PID控制器前后的响应曲线进行分析。数学模型的建立直线一级倒立摆是由小车、摆杆等部件构成。现假设有一个一级倒立摆系统，其中摆杆的长度为L，2Ll，质量为m，小车质量为M，u为作用在小车上的外力，为摆杆与垂直向上方向的夹角，x为小车的水平位移。若不考虑小车与导轨，摆杆与小车之间的摩擦并且不计各种空气阻力，则倒立摆系统的受力分析如图1所示。计算机控制理论与设计-2-Figure1规定摆杆重心的坐标为,GGxy，则sincosGGxxlyl采用隔离法，对小车有212dxuTMdt(1)其中1T为摆杆所受沿水平方向的力。对于摆杆，在水平方向上有212sindxlTmdt(2)在竖直方向上有222cosdlmgTmdt(3)其中2T为摆杆所受沿垂直方向的力。以摆杆质心为旋转中心建立动力学方程有2212sincosdTlTlJdt(4)其中2112JmL为摆杆绕其质心的转动惯量。此时，发现(2)(3)(4)均为非线性方程，为了简化计算并且便于进行控制器设计，要对上述方程进行近似线性化处理。由于控制目标是使摆杆处于垂直稳定状态，即尽可能的小，甚至接近于0，因此可以认为sin0、cos1、20、0。这里需要注意计算机控制理论与设计-3-222sinsincosddt，222coscossinddt因此方程(2)(3)(4)化为1Tmxml(5)20mgT(6)21TlTlJ(7)消去中间变量1T、2T，由(1)(5)式得到Mmxmlu(8)由(6)(7)式得到2()Jmlmlxmgl(9)由(8)(9)得到单级倒立摆方程22()()()mMmglmluMmJMmlMmJMml(10)22222()()mglJmlxuMmJMmlMmJMml(11)设小车的质量M为2Kg，摆长L为1m，摆杆质量m为0.1Kg，重力加速度g取210ms带入上述方程，可以得到角度与外力之间的传递函数Gs，20.7415.56sGsUss(12)由于系统有位于复平面右半平面的极点显然其不稳定，在阶跃函数作用下的响应的响应曲线如下图所示，同样说明系统是不稳定的，因此需要加入PID控制来改变系统的性能。Figure2计算机控制理论与设计-4-状态反馈PID设计令1x，2x，3xx，4xx建立状态方程，则由(10)(11)可得倒立摆系统的状态方程为13240100000000010000ttxxuAxbutt(13)其中12()()mMmgltMmJMml，2222()mgltMmJMml，32()mltMmJMml，242()JmltMmJMml接下来我们要用状态反馈的方法设计此倒立摆系统的PID控制器。采用文献[4]的方法，在控制输入点之前引入一个积分环节zutut，其中0ut给定。记xztu，00AbA，00001Tb因此得到增广系统的状态方程zzAzbu(14)设状态反馈zutKzt(14)1212zxtututKKKxtKutut(16)因为1rankb，则Tbb可逆，从而可以将控制输入ut用状态及状态的导数表示，即1TTutbbbxAx(17)将上式代入(16)式中，得到11112212TTTTTTutKxKbbbxAxKbbbxKKbbbAx(18)记12TTKbbb为3K，112TTKKbbbA为4K，在零初始条件下对上式积分得到计算机控制理论与设计-5-340tutKxKxd(19)可以看出这是状态反馈加状态积分反馈的形式。如果能够找到一组状态变量，使得状态反馈就是比例积分加微分，那么该状态反馈策略就是PID控制。对于n阶系统，微分部分包括了系统从1阶到1n阶的微分；若是二阶系统，该部分就是一个工程中常用的PID控制器。现在讨论反馈增益矩阵K，使得增广系统满足稳定的指标。即设计反馈矩阵K使得闭环系统稳定。此时闭环系统的状态方程为ztAbKzt(20)根据Lyapunov理论可得，存在对称正定矩阵P使得0TTTAPKbPPAPbK(21)在上式左右两边左乘1P，右乘1TP，并令1LKP则得到110TTTPALbAPbL(22)这是关于L与1P的线性矩阵不等式，又由于P正定则1P正定，因此上述LMI可以改写为111000TTTPALbAPbLP(23)通过MATLAB可以求得L与1P，由此可以解出对称正定矩阵P与反馈增益K。下面通过MATLAB进行仿真和计算。仍然采用之前仿真时使用的小车质量与摆杆质量等数据，可以得到增广系统的状态方程的系数为0100015.560000.74000100.370000.4900000A，00001b通过MATLAB计算可得反馈增益K为1126.8281.413.635.223.7，因此3022.25014.73K，4775.14281.3913.6235.19K，代入控制器u的表达式中得到340tutKxKxd，仿真曲线图如下所示。计算机控制理论与设计-6-Figure3Figure4Figure5计算机控制理论与设计-7-Figure6Figure7结论由图2可知系统在没有控制其作用下是不稳定的，通过状态反馈PID控制之后的系统状态的响应曲线分别为图3~6，可以看出，只在稳定的约束下，系统的四个状态可以在18秒左右达到稳定，证明方法的可行性。同时也发现了只对稳定性约束时的不足，即系统达到稳定所需的时间过长，为此可以在考虑稳定性的条件下加入其他指标进行约束以达到预期的效果，例如极点指标、方差指标、H指标等，甚至可以进行多指标约束，例如极点指标加方差指标等。求解时需要将这些指标化为线性矩阵不等式组，运用MATLAB求解满足需求的解。计算机控制理论与设计-8-讨论与体会在我之前的认识中，对于PID只理解到自动控制原理中通过对超调量等指标的约束画出伯德图进行求解参数的阶段。后来学习了状态空间描述系统时，就思考如何在状态方程的基础上进行PID控制器设计，因此借此计控作业的机会查找资料学习了一下。在这个作业中我还用到了一些这学期学到知识，例如线性矩阵不等式的求解和系统稳定性分析。在学习了线性矩阵不等式后也借此机会实践了一下其中的变量替换和MATLAB求解,同时也是第一次听到满意控制这个概念。因此可以说这次作业不仅仅是计算机控制的作业也从一定程度上说是对这学期的小小总结，同时在这之中也发现还有很多需要学习的地方。参考文献1.张嗣瀛，高立群．现代控制理论[M]．北京：清华大学出版社，20062.关守平，尤富强．计算机控制理论与设计[M]．北京：机械工业出版社，20123.张庆灵，张雪峰，翟丁．控制理论基础[M]．北京：高等教育出版社，20084.马建伟，李银伢．满意PID控制设计理论与方法[M]．北京：科学出版社，20075.周俊．PID控制在单级倒立摆系统中的分析与应用[J]．硅谷，2010，6：20-216.王慧萍，孔庆忠．基于MATLAB的直线一级倒立摆的PID控制研究[J]．机械工程与自动化，2015，5:179-1827.薛定宇．控制系统计算机辅助设计[M]．北京：清华大学出版社，20158.王建辉，顾树生．自动控制原理[M]．北京：清华大学出版社，2013附录:代码clearallnum=[-0.74];den=[10-15.56];G=tf(num,den);t=[0:0.1:20];[y,x,t]=step(num,den,t);figure(1);plot(t,y);grid;计算机控制理论与设计-9-xlabel('t');ylabel('y');title('单位阶跃响应曲线')%角度对输入的传递函数的阶跃响应曲线a=[0,1,0,0,0;15.56,0,0,0,-0.74;0,0,0,1,0;-0.37,0,0,0,0.49;0,0,0,0,0];b=[0;0;0;0;1];setlmis([])%开启LMI求解Q=lmivar(1,[5,1]);L=lmivar(2,[1,5]);%定义矩阵变量lmiterm([111Q],a,1,'s');lmiterm([111L],b,1,'s');lmiterm([122Q],-1,1);%描述线性矩阵不等式lmisys=getlmis;[tmin,xfeas]=feasp(lmisys);QQ=dec2mat(lmisys,xfeas,Q);LL=dec2mat(lmisys,xfeas,L);K=LL*inv(QQ)%解出反馈矩阵K%状态方程c=[0,0,0,0,0];d=[0];g=ss((a+b*K),b,c,d);t=[0:0.1:20];x0=[0,0,1,0,0];[y,t,x]=initial(g,x0,t);figure(2)plot(t,x);holdongridon%得到增广系统的状态响应曲线figure(3)plot(t,x(:,1)),grid;xlabel('t'),ylabel('角度'),title('摆角随时间变化的曲线');figure(4)计算机控制理论与设计-10-plot(t,x(:,2)),grid;xlabel('t'),ylabel('角速度'),title('角速度随时间变化的曲线');figure(5)plot(t,x(:,3)),grid;xlabel('t'),ylabel('位移'),title('小车位移随时间变化的曲线');figure(6)plot(t,x(:,4)),grid;xlabel('t'),ylabel('速度'),title('小车速度随时间变化的曲线');figure(7)plot(t,x(:,5)),grid;xlabel('t'),title('控制器变化曲线');