小论文——倒立摆与SI鲁棒极点配置

1、鲁棒极点配置概念在单输入系统控制中的应用李昊1 陈雯柏2 恒庆海2（1. 北京机械工业学院 计算机与自动化系 北京 100085；2. 北京信息工程学院 信息与通信工程系 北京 100101）摘 要：利用矩阵特征值一阶摄动理论解释了鲁棒极点配置的原始目标，指出以特征向量广义夹角加权和最大化为目标的一些方法存在认识偏差，并根据控制问题的特点重新规定了适用于评估闭环稳定性的评判数；在MATLAB环境下编制了可直接使用拉格朗日方程的程序，求出了二级倒立摆稳定位置附近局部线性化模型，以该评判数为基础在GIP-300-VPPA-L型倒立摆上实现了对二级倒立摆的稳定控制。建模程序经适应性修改后可用于对任意

2、可用拉格朗日方程描述的系统进行建模。关键字： 特征值一阶摄动 状态反馈 鲁棒极点配置 拉格朗日方程Robust pole assignment concept used on the SI plantLi Hao1,Chen wen-bai2, Heng qing-hai2(1. Department of Computer Scicnce & Automation, Beijing Institute of Machinery, Beijing 100085, China2. Department of Information & Communication Engineer

3、ing, Beijing Information Technology Institute, Beijing 100101, China)Abstract: The concept of robust pole assignment was re-explained using eigenvalue perturbation theory and some misunderstandings were pointed out. A new criterion for measuring the stable margin of state-space closed loop system wa

4、s given. Then the inverted pendulum is stabilized by the controller, which was selected by the criteria.Key words: eigenvalue; state-feedback; robust pole assignment; Lagrange equation引言应用线性系统理论的极点配置方法解决控制问题的时候，在证明一个以状态空间模型描述的线性（子）系统为能控后，常基于纯粹的性能因素指定目标极点，为了在参数发生摄动的情况下使极点配置的结果具有鲁棒性，CAVIN等人1提出了鲁棒极点配置（

5、RPA）概念，即依靠MI系统控制器相对目标极点的自由度设计具有鲁棒能力的控制器，保证闭环极点对系统参数摄动具有较低的敏感度。为此目标，CAVIN1提出了第一种解决方案，实质上以特征向量广义夹角加权和最大化为目标（以下简称MaxA）；后来KAUTSKY等人2以使闭环系统矩阵关于特征值问题的谱条件数最小为目标，以鲁棒特征结构配置（以下简称RFSP）的名目完成了第二种方法，通过使各极点对应的特征向量矩阵尽量正交化，具体的执行过程可形容为“循环掰”。根据CAVIN的方案还衍生了一系列类似的方法3、4、5。但这些算法都缺乏对不同闭环极点的区别对待，因为在所有极点具有同样摄动敏感度时，显然距离虚轴更远的比

6、较近的对系统稳定性威胁更小；MaxA类型的方法试图使各右特征向量内积之和尽量小，客观上相当于要各右特征向量之间互相正交化（90°是单角最大），在输入比状态数小不多时容易利用自由度向这一目的靠拢5，但杨亚光等人3指出，这种方法仅在右特征向量间内积绝对值之和的两倍小于n/(n-1)时才有效。问题描述由特征值一阶摄动理论6（以下简称EPT）得知，设某一n维方阵有一组相异的特征值，是一个很小的常数，为与维数相同且具有代表摄动形式特征的矩阵，在受到摄动变为后，其单根在影响下的摄动大小等于，与分别是对应的单位化了的左、右特征向量。由于实际系统的复杂性、以及在控制问题中因状态反馈下闭环矩阵具有的形

7、式，闭环摄动阵不仅与、，甚至还与有关，故经常很难对进行正交化分解和确定元素的比例大小。但即使在无法确定形式的情况下，特征值摄动的大小也仍然与成反比，并且该组值仅取决于名义闭环系统矩阵，可见，各特征值的摄动灵敏度判据之间是独立的，可能出现同一个闭环系统的两个特征值摄动灵敏度相差悬殊的情况。因此可以以使的某元素尽可能大作为目标（以下简称EPT目标）提高对应闭环极点的抗摄动能力。RFSP法试图通过使的右特征向量矩阵尽量为正交阵（即尽量大）来达到RPA目的，只使用了右特征向量系，似乎与同时需要左、右特征向量系的EPT理论不同，但实际上由于成立，即任一特征值的左特征向量与其它特征值的右特征向量正交，故任

8、一特征值的左特征向量肯定在由其它特征值的右特征向量展成空间的零空间中，反之亦然。考虑到无重根情况下特征值的左特征向量系和右特征向量系都肯定是满秩的，即的列向量之间线性无关，的行向量之间也线性无关。由于线性无关弱于正交，即一特征值的右（左）特征向量一般不在其它特征值的右（左）特征向量转置展成空间的零空间中而是与之呈一定夹角，所以一般有，如果正好在该零空间中就有成立。这样，RFSP的做法就符合了基于EPT理论产生的使的目标。而除了完全正交的情况，各特征向量广义夹角和最大与以上目标并不重合，甚至不能排除各特征向量位于一个广义平面内的情况，二阶系统时MaxA方法与RFSP法目标一致，而对三阶及更高阶系

9、统，该情况意味着所有特征值毫无抗摄动能力。算例1将以一个三阶系统为例说明这一问题，同时验证小摄动下EPT理论的正确性。主要结果控制问题更关心特征值摄动后是否越过性能界或稳定界，而不应该泛泛地要求所有特征值具有接近的摄动灵敏度，故以更强调特征值整体抗摄动能力的作为评判数就不太合适，从稳定性来说，对距虚轴较近的点显然应该有更严格的摄动限制，应该把一组目标极点对应的最小作为更有效率的评判标准，该值越大越好。建模得到的难免存在不精确现象，由于在中可以看出与由的不确定造成的闭环下摄动成正比，大的造成闭环系统矩阵中更大的参数摄动，从而有时一些看似摄动影响非常小的目标极点对应的控制结果反而不好；又由于过强的

10、输入会激发系统内部包括非线性在内的更复杂的问题，所以在用于实际问题时不能仅使用一个评判标准，而应该既保证较大，又限制对应的控制器反馈倍数。闭环目标极点确定对应某一右特征向量空间而不受其它闭环极点影响2，但其左特征向量却完全取决于其它闭环极点而不受对应极点的影响。在MI情况下，即使已经确定，由于特征向量矩阵的列向量拥有自由度，通过选择不同组合可产生一系列。算例2将利用倒立摆验证上述目标的效果，但在SI系统中闭环极点决定右特征向量，闭环目标极点组也与一一对应，虽然可以简化设计过程，但因没有MI系统可用的自由度，无法直接使用EPT方法的结论，要使用EPT目标，只能采取在一系列满足性能要求的目标极点组

11、中循环的办法，即验算哪些组的较大，以此确定闭环目标极点，同时得到对应的的集合。在用于SI系统中时，仍需进行RFSP法2所要求的前期准备工作。对SI系统的具体步骤如下，由于使用状态反馈，所以不写出输出方程。a) 对进行QR分解，有，为一个n维列向量，为一个标量。b) 从目标性能出发，划定一个满足性能要求的目标极点区域。由于最大目标极点实部小的闭环系统在同样偏离平衡点时一般需要更大的控制输入，容易导致更多问题，所以要权衡选取。c) 从RFSP理论得知，当某一闭环特征值重数大于输入维数时，会导致X奇异，即使能够得出名义控制器，但配置在这里的几重目标极点对参数摄动均毫无抵抗力，极点配置的目标更是工程问

12、题而非数学问题，因此指定SI系统目标极点互异，计算对应这些目标极点的并单位化得到，将各按的顺序排成一个n阶方阵。由于无重根，故可逆，求其逆，单位化的各行向量，得到。而后以作为该组的评判数，该指标同样越大越好，再考虑控制器反馈倍数进行取舍。采用遍历循环的办法，在目标极点配置区域内按某一策略选取很多组目标极点，在所有循环到的极点集中选出一（多）组符合评判标准的作为最终的配置极点，并同时得到多个控制器，这样做的好处是可以有更多的检验机会。由于SI系统状态静态反馈闭环极点与控制器之间一一对应，确定了后可以用任何方法求取。算例1假设三阶系统的目标极点为-1 2 -3，假设仅能将单位化前的右特征向量矩阵X

13、配置成或（实际上把X1的0.1变成0后，就是三个共面且互相成120°（60°）夹角的矢量，加0.1纯粹是为了可以求逆），对该两矩阵进行列单位化（为了EPT理论计算）。根据相似变换，不管是否对X1、X2进行列单位化，由极点对角阵和该两右特征向量矩阵总能得到和，分别记为A1、A2，两者的MaxA指标1分别为1.4901和5.3531（越小越好），注意A1的指标勉强小于3/（3-1）1.5，而A2的指标远超过该值。同时， A1的三个特征值各自的EPT衡量标准值分别为0.0995、0.0993、0.0993，而A2的为0.2075、0.2873、0.2873（越大越好），MaxA与

14、EPT的衡量结果相反。后两个特征值的EPT标准相等是因为X1、X2的列向量在三维空间对应的矢量有对称关系，当然，在对称阵I-XTX的对角线同侧也会出现两个相同的值，即该阵中会有四个一样的值，但EPT与MaxA的成因不同。另外，当X发生变化时，A也跟着变化，下面为了简化，只写出X1、X2的变化，讨论时保证对应关系。A1、A2两矩阵各元素均受到在±0.005范围内均匀分布的摄动时，1000次同样随机摄动尝试下的最大特征值结果如图1所示，其中标“”的是A1在摄动下的表现，标“o”的是A2的。其它随机尝试的结果均与其大同小异，并且所有摄动影响下的实际特征值都是纯实数。图1 同强度摄动下A1（

15、）与A2（o）特征值最大实部受到的影响可见MaxA指标较小的A1特征值受到影响的程度反而要大得多，特征值变化更符合EPT比较的结果。当X1中的0.1项减小时，其MaxA指标虽始终小于1.5，但摄动对A1特征值的影响迅速增大。当X1的0.1和X2的0.3均变成10时，MaxA指标分别减小到0.5099 和0.0398，而对应特征值1的EPT指标分别增加到0.9950和0.9901，同样摄动下特征值最大实部的情况如图2。图2 参数变成10后同强度摄动下A1与A2特征值最大实部受到的影响可见EPT指标比MaxA更接近实际情况，且趋势也一致。继续增大那三个值时情况也类似。所以，杨亚光等人的证明3只对M

16、axA本身有意义，而MaxA却与RPA没有直接联系。本质上是因为其与EPT理论不重合。算例2使用固高公司的GIP-300-VPPA-L型倒立摆，用计算机的数字运算、执行器控制连续系统，采样间隔0.005s，理论计算中的连续系统要用数字控制器来控制，从采样/保持信号与模拟信号之间的锯齿形差异和时间延迟、以及电机皮带系统的响应动态来看，也是一种加入的时变不确定因素，客观上可以检验理论的适应能力。GIP-300-VPPA-L在接成直线二级倒立摆状态下的物理模型如图3。使用的摆杆为：摆杆1为标准配置的短杆，摆杆2为标准配置的长杆。广义坐标方向即为该坐标的正向，小车质量M1.32，均质摆杆1质量m1=0

17、.04，半长l1=0.09，杆相对其质心的转动惯量J1=1/3*m1*l12=0.0001，摆杆2质量m2=0.132，半长l2=0.27，杆相对其质心的转动惯量J2=1/3*m2*l22=0.0032，质量块（实际是摆杆2的角度传感器）质量m3=0.2086，忽略所有摩擦。因为仅考虑上平衡点附近的稳定性问题，所以仅建立此处的局部线性化模型。本文利用MATLAB的符号运算功能处理第二类拉格朗日方程，编程计算装置数学模型的符号解和数值解，程序如下（状态变量顺序为，图3中的方向即为其正向。标“”的行在经过适当更改后，该程序可用于建立更多可使用拉格朗日方程描述的局部线性模型）：clearn=6; %

18、 syms M m1 m2 m3 l1 l2 L1 J1 J2 g positive; % syms ddx dx x ddsita1 dsita1 sita1 ddsita2 dsita2 sita2 F L real; % temp1=x dx sita1 dsita1 sita2 dsita2; % 状态变量temp2=dx ddx dsita1 ddsita1 dsita2 ddsita2; % 状态变量微分temp3=ddx ddsita1 ddsita2; % （角）加速度temp4=M m1 m2 m3 l1 l2 L1 J1 J2 g; % 物理参数init=sym(1.32

19、.04 .132 .208 .09 .27 2*.09 1/3*.04*.092 1/3*.132*.272 9.8); % 物理参数值L=1/2*M*dx2+1/2*m1*(dx-l1*dsita1)2+1/2*m2*(dx-L1*dsita1-l2*dsita2)2+1/2*J1*dsita12+1/2*J2*dsita22+1/2*m3*(dx-L1*dsita1)2+m1*g*l1*sita12/2+m2*g*(L1*sita12/2+l2*sita22/2)+m3*g*L1*sita12/2; % 已经局部线性化后的拉格朗日函数sysA=subs(diff(L,dx),temp1,t

20、emp2)-diff(L,x);subs(diff(L,dsita1),temp1,temp2)-diff(L,sita1);subs(diff(L,dsita2),temp1,temp2)-diff(L,sita2); % 拉格朗日方程左端项。sysB=F;0;0; % 拉格朗日方程右端项。E=jacobian(sysA,temp3);% 求加速度项的导数，仅使用temp4中参数。R=sysB-simple(sysA-E*temp3');% 系统方程为E*dx=R。sys=ER;% 系统方程为dx=sys。preA=jacobian(sys,temp1);% 求s ds sita1

21、dsita1 sita2 dsita2 的具体系数。preB=diff(sys,F);% F的系数。post,sg = subexpr(preA preB,'sg');% 求符号解preA和preB的简化形式。partA=double(subs(preA,temp4,init);% 符号解转变成数值解。partB=double(subs(preB,temp4,init);A=zeros(n,n);b=zeros(n,1);for i=1:n/2 A(2*i-1,2*i)=1; A(2*i,:)=partA(i,:); b(2*i)=partB(i);end% 最终的模型数值解完

22、成。数值解如下：，。从其形式可以看出，无论物理参数如何摄动，最终的闭环系统矩阵的摄动只发生在第二、四、六行，也就是说的第一、三、五分量不管多大，对最后的特征值摄动都没有影响，考虑到，因此这三个分量大，意味着会产生影响的第二、四、六分量小，因此用更符合实际情况的代替稍粗糙的作为评判数。设定目标极点在、所包围的矩形区域中，以1的步长进行取点，并遍历寻求所有可能的组合，满足评判数大于0.08、最大反馈倍数小于200，得到8个控制器，有6个可以实现稳定，其中一个的目标极点为，为-9.6149 -7.7542 -89.9405 -3.5447 119.5499 18.7666。两个不能保证闭环稳定的解并

23、不能证明基于EPT的算法不能成立，而只能说明抗摄动的方向与实际摄动方向相差太大，且实际摄动范围也大，这表明 EPT法也存在一定的保守性。当使用两个标准的中等长度摆杆时，按文献6中参数外推，摆杆质量m1=m2=0.083，半长l1= l2=0.17，杆相对其重心的转动惯量J1=J2=1/3*m2*l22=0.0008，其它不变，将程序中参数初始化命令init行改为“init=sym(1.32 .083 .083 .208 .17 .17 2*.17 1/3*.083*.172 1/3*.083*.172 9.8);”得到的7个控制器都可以使新的二级对象稳定，其中一个是 -4.8677 6.728

24、1 126.2604 4.9033 157.4546 21.053，目标极点为。后来的实验还发现上述能达到稳定的闭环系统对多种情况有一定适应能力，鼓励有条件的读者检验，但请注意坐标方向和变量顺序。结束语基于EPT理论证明了RFSP方法的正确性，指出了MaxA类型目标的错误，并改进得到了更有针对性的RPA评判数。该方法简单方便，物理意义明确，无需求解Riccati方程，并且文中所取的评判数门槛还有相当大的裕量。在二级倒立摆实物上实验的成功显示了新评判数的效果，而且分离了的各目标极点鲁棒性标准也为将来更细致的区别对待创造了条件。同时这种鲁棒极点配置评判数的选取也可以用于如观测器设计等场合。参考文献1 CAVIN R K,et al., Robust and Well-conditioned Eigen-structure Assignment Via Sylvesters Equations, American Control Conference, 1982, 1053-1057.2 KAUTSKY J,NICHOLS N K,VAN DOOREN P.Robu