崔 颖，王 曦

（1.中国航发贵州红林航空动力控制科技有限公司，贵阳550009；2.北京航空航天大学能源与动力工程学院，北京100191；3.先进航空发动机协同创新中心，北京100191）

0 引言

现代航空发动机的控制系统采用先进的全权限数字电子和多变量控制技术。为了获得优越的性能，现代航空发动机工作范围非常接近部件的机械、气动和热负荷的极限最大状态。因此，需要采用先进的多变量控制技术以保证发动机在全工作范围内的高精准稳态控制[1-2]，具有先进的伺服跟踪和抗干扰的鲁棒控制器成为现代航空发动机控制的核心技术。

在工程设计中，许多闭环反馈控制设计问题可以转化为静态输出反馈进行求解[3-4]。近年来，相关研究取得了一系列成果。文献[5]证明了对于严格真的线性系统都存在1 个稳定的静态输出反馈控制律的Riccati 方程解；文献[6]在离散域研究了参数不确定性的静态输出反馈设计问题；文献[7]提出了1 种改进的静态输出反馈线性矩阵不等式的迭代求解算法；文献[8] 通过最小化2 个决策矩阵变量的方法将非凸优化问题转化为凸优化问题进行求解，从而放宽了静态输出反馈控制的约束条件；文献[9]提出了1 种通过求解线性矩阵不等式LMI（Linear Matrix Inequality）和矩阵等式组成的方程组，获得了严格真线性状态空间模型输出反馈控制律。

LMI 是MATLAB 软件中解决鲁棒控制问题的有效工具[10-12]，以椭球法和内点法作为求解凸优化数值求解问题的基础，近年来在控制领域的应用已十分广泛[13-15]；采用双目标最优化修正的自适应控制则是近年来针对系统输入中存在不确定性问题的1 种鲁棒控制的有效方法[16]。

在稳态控制器设计中，闭环极点在复平面的不同位置决定了系统的稳定性及动态性能[17]。本文以高精准稳态控制为目标，在静态输出反馈渐近稳定条件的基础上，考虑控制系统的稳定裕度和鲁棒稳定性能，推导了极点配置圆的条件，提出了1 种多变量控制系统闭环极点配置圆的LMI 设计方法，并对涡扇发动机的双回路控制和单回路控制进行了仿真验证。

1 输出反馈PI 控制的转化

设包含执行机构动态特性的被控对象为

式中：x(t)∈Rn为系统状态向量；u(t)∈Rp为系统控制向量；y(t)∈Rm为系统输出向量。

考虑被控对象∑1的闭环输出反馈PI 控制律为

式中：Kp和Ki分别为PI 控制器的比例增益和积分增益。

以被控对象∑1的状态向量和输出向量构造另一个被控对象∑2的新的状态向量为

式中：T 为转置符号。

则存在以下关系

式中：I∈Rm×m为单位对角矩阵。

将式（5）和式（6）代入式（2），则

并定义以下矩阵

则另一个被控对象∑2可用状态向量z(t)∈Rn+m、控制向量u(t)∈Rp、输出向量q(t)∈Rm+m构造如下

由式（7）和式（9）可知，存在∑2的闭环反馈输出控制律

则

2 极点配置圆条件

∑2闭环系统等价于如下自治系统

式中：Kg为被控对象∑2的输出反馈增益矩阵。

选Lyapunov 二次函数V(z) = zTQz ＞0，则

式中：Q 为正定矩阵。

为使V˙(z)负定，令Q-1= P，则

在上述渐近稳定条件的基础上，考虑控制系统的稳定裕度和鲁棒稳定性能，设存在1 个正定矩阵P＞0和鲁棒稳定性能特征因子α ＞0, β ＞0，构造如下矩阵不等式

上述不等式（15）的第2 项P(Ag+ BgKgCg+ αI)T+(Ag+ BgKgCg+ αI)P ＜0反映了∑2闭环系统的稳定裕度，而第1 项β(Ag+ BgKgCg)P(Ag+ BgKgCg)T＜0 反映了∑2闭环系统的鲁棒稳定性能。

设ξ 为(Ag+ BgKgCg)T的特征值λ 的特征向量，对上式左乘ξ*右乘ξ，得

式中：λ*为λ 的共轭转置，ξ*为ξ 的共轭转置。即

式（17）表示闭环系统的极点矩阵(Ag+ BgKgCg)T的特征值λ 均落在如图1 所示的以(-c,0)为圆心、r 为半径的圆域内。

图1 闭环极点圆在复平面上的分布区域

对上述矩阵不等式进行变换，得

由Schur 补引理，得线性矩阵不等式

设Y = KgCgP，则

求式（21）LMI 的解，其解P 若为正定阵，可构造如下矩阵

当M 非奇异时有惟一解为

当M 奇异时，存在广义逆的逼近解为

3 多变量控制极点配置算法

综上所述，同时考虑到涡扇发动机多变量控制系统中不同物理参数量纲变化较大，应对其进行归一化处理，形成多变量控制极点配置算法：

（1）对在稳态设计点获得的涡扇发动机状态空间模型进行归一化处理，获得涡扇发动机状态空间归一化线性模型；

（2）将执行机构动态模型增广到归一化线性模型中，构建涡扇发动机增广状态空间模型∑1；

（3）将被控对象∑1转化为∑2；

（4）给定极点配置圆的几何参数c、r；

（5）由式（26）计算稳定裕度和鲁棒稳定性能特征因子α、β；

（6）用LMI 工具箱求解式（21），求得P、Y；

（7）由式（22）构造M 矩阵，由M 的奇异性，通过式（24）、（25）求得Kg；

（8）对式（11）分解，可得Kp、Ki。

4 仿真验证

为了验证上述方法的伺服跟踪和抗干扰性能，分双轴涡扇发动机的多变量和单变量控制进行仿真和分析。

4.1 双回路控制

双轴涡扇发动机状态空间模型为

式中：状态向量为x = [NLNH]T，NL为低压转子转速，NH为高压转子转速；输入向量为u = [ WfA8]T，Wf为主燃油流量，A8为尾喷口喉道面积；输出向量为y =[NHπT]T，πT为涡轮落压比。

设

则可得归一化线性模型

采用双回路控制系统结构，控制目标

设调节主燃油流量回路和尾喷口喉道面积回路的执行机构传递函数是时间常数为0.1 s 的1 阶惯性环节。

将发动机归一化线性模型与执行机构模型进行增广，增广状态空间模型的系数矩阵为

考虑到闭环极点圆的圆心位置在复平面上，若靠近虚轴，系统的动态响应会变慢；若离虚轴太远，系统的动态响应太快，与较慢的执行机构动态产生不匹配的问题。折中考虑后圆心选为（-5，0）；同时，考虑到极点不能落到复平面的右半平面内，以及极点位置不能太靠近虚轴，通过半径对极点配置圆进行约束，半径选为r = 4，可得α = 0.9，β = 0.2，按本文所述方法求解LMI，计算结果为

将上述解按归一化增广对象求得的PI 控制器进行反变换，得

对上述控制系统进行阶跃和斜波响应的双回路闭环仿真验证。仿真时间为第5～45 s，其中在第5 ～30 s 考察阶跃跟踪响应情况，在第30～45 s 考察斜波跟踪响应情况。高压转子转速参考指令如图2所示，涡轮落压比参考指令如图3 中所示。

图2 NH 阶跃指令、斜波指令和NH 响应曲线

图3 πT 阶跃指令、斜波指令和πT 响应曲线

从图2 中可见，高压转子转速参考指令在第5～10 s 保持NH= 10000 r/min 不变，在第10 s 加入ΔNH= 2800 r/min 的阶跃信号，在第10～20 s 保持NH=12800 r/min 不变，在第20 s 加入ΔNH=-2800 r/min的阶跃信号，在第20～30 s 保持NH=10000 r/min 不变，在第30～40 s 加入斜波指令信号，斜率为280 r/min/s，在第40～45 s 保持NH= 12800 r/min 不变。

从图3 中可见，涡轮落压比参考指令在第5～15 s保持πT= 8 不变，在第15 s 加入ΔπT= 0.5 的阶跃信号，在第15～25 s 保持πT=8.5 不变，在第25 s 加入ΔπT= -0.5 的阶跃信号，在第25～35 s 保持πT=8 不变，在第35～40 s 加入斜波指令信号，斜率为0.1/s，在第40～45 s 保持πT=9 不变。

高压转子转速NH和涡轮落压比πT的伺服跟踪响应曲线如图2、3 虚线所示。在第5～30 s 的仿真过程中可见，在双回路控制中，2 个回路在各自不同的阶跃输入指令下，在NH和πT动态调节过程中存在相互耦合干扰。

当加入NH阶跃指令信号时，NH响应能够伺服跟踪第10、20 s 的阶跃指令，对第15、25 s 由另一回路πT阶跃响应耦合的干扰具有抑制效果，且进入稳态后，能够无静差伺服跟踪参考指令。

当加入πT阶跃指令信号时，πT响应能够伺服跟踪在第15、25 s 的阶跃指令，但在第10 、20 s 由另一回路NH阶跃响应耦合的干扰，会由于高压转子转速NH的突变对涡轮落压比πT带来超调量约为1.3%的干扰，这种影响的动态调节时间不大于2 s，随后进入稳态后，能够无静差伺服跟踪参考指令。

对于第30～40 s 的NH斜波指令和对于第35～40 s 的πT斜波指令，2 个回路之间的干扰作用很小，NH斜波跟踪指令误差不大于1.5%，πT斜波跟踪指令误差不大于0.4%。

4.2 单回路控制

固定喷口面积的双轴涡扇发动机归一化状态空间模型为

执行机构是时间常数为0.2 s 的1 阶惯性环节。为了考察执行机构动态对控制系统的影响，分以下2种情况进行对比分析。

4.2.1 情况1：不考虑执行机构动态的设计

极点圆配置设计在以（-5，0）为圆心，r= 4 为半径的圆内，按上述方法对归一化固定喷口面积的双轴涡扇发动机求解控制器得Knp=4.8349，Kni=17.1438。其闭环系统的极点、零点如图4 所示。

图4 情况1 的闭环系统的极点、零点

从图中可见，这3 个闭环极点其中1 个极点为-1.4+0i，与1 个闭环零点距离很近，因此动态性能主要由1 对共轭主导极点-3.2±1.1i 决定，其阻尼比约为0.95。

其Nyquist 曲线如图5 所示，Bode 图曲线如图6所示。从图中可见，系统有无穷大的幅值裕度和86毅的相角裕度。

首先，进行不带执行机构模型的闭环阶跃响应仿真，参考指令如图7 中的实线所示。在稳态点分别在第2、6、10、14 s 加入4 个小阶跃和在第18 s 加入1个大阶跃参考指令信号，获得的低压转子转速响应如图7 中的虚线所示。无动态超调，调节时间约为1 s左右，稳态误差为0。

图5 情况1 的Nyquist 曲线

图6 情况1 的Bode 图曲线

图7 情况1 的不带执行机构模型的闭环NL 阶跃响应曲线

其次，将执行机构模型嵌入闭环中进行仿真，在不同阶跃幅值下的NL响应动态性能变差，在第18 s加入较大的阶跃指令幅值下的响应超调量约为11%，如图8 所示。

为了分析动态性能变差的原因，对执行机构模型增广得到的开环传递函数为

其闭环极点、零点分布如图9 所示。

图9 增广执行机构后的闭环系统极点、零点

从图中可见，这4 个闭环极点中的1 个与1 个闭环零点距离很近，动态性能主要由1 对共轭主导极点-2.3±3.2i 决定，其阻尼比约为0.55。

其Nyquist 曲线如图10 所示。从图中可见，Nyquist 频谱左向弯曲靠近（-1，0i）点。

其Bode 曲线如图11 所示。从图中可见，虽然具有无穷大的幅值裕度，但是相角裕度只剩下55毅。

可见，在不考虑执行机构动态时，闭环系统的阻尼比将减小0.4，相角裕度减少了30毅左右，这是导致系统动态性能变差的主要原因。

图11 增广执行机构模型后的Bode 图曲线

如果进一步考察传感器噪声对控制系统性能的影响，加入频率为10 Hz、幅值为±0.5 的均值为0 的高斯白噪声测量信号，带噪声的转速测量信号如图12 所示。

图12 情况1 的带有噪声NL 的测量信号

仿真结果如图13 所示。从图中可见，在不同阶跃幅值下的阶跃响应动态性能进一步变差，不仅超调增大，还出现不同程度的转速摆动现象。

图13 情况1 下闭环NL 阶跃响应曲线

4.2.2 情况2：考虑执行机构动态的设计

作为对比设计，将发动机归一化线性模型与执行机构模型进行增广，增广状态空间模型矩阵为

极点圆配置同上，按上述方法求解，得Knp=2.2048，Kni=7.0648。其闭环极点、零点分布如图14所示。其Nyquist 曲线如图15 所示。

Bode 图曲线如图16 所示。从图中可见，系统具有无穷大的幅值裕度和近78毅的相角裕度。

图14 情况2 的闭环系统的极点、零点

图15 情况2 的Nyquist 曲线

图16 情况2 的Bode 图曲线

首先，在仿真中未加入传感器测量噪声，在不同阶跃幅值下的低压转子转速响应曲线如图17 所示。NL转速无超调，调节时间为1.5 s。

图17 情况2 的闭环NL 阶跃响应曲线（传感器不带噪声）

其次，在仿真中加入频率为10 Hz、幅值为 的均值为0 的高斯白噪声测量信号，转速传感器信号如图18 所示。

带传感器噪声的不同转速阶跃幅值下的响应曲线如图19 所示。NL转速无超调，调节时间为1.5 s。

图18 情况2 的带有噪声的NL 测量信号

图19 情况2 的闭环NL 阶跃响应曲线（传感器带噪声）

进一步考虑执行机构建模的不确定性，设执行机构的实际时间常数Ta=0.35 s，加入的传感器噪声信号同前，NL测量信号如图20 所示。

仿真结果如图21 所示。从图中可见，即使在执行机构未建模动态存在0.15 s 的情况下，控制系统的动态性能未明显变差，仍具有鲁棒性能，小阶跃超调量在0.5%之内，调节时间小于2 s，大阶跃超调量在3%之内，调节时间小于3.5 s。

图20 情况2 的含噪声的反馈NL 信号

图21 情况2 的闭环NL 阶跃响应曲线（执行机构Ta=0.35 s）

5 结论

本文提出了1 种多变量控制系统闭环极点配置圆的LMI 设计方法，在双转子涡扇发动机上进行了仿真验证，得到结论如下：

（1）对双转子涡扇发动机双回路控制的仿真表明：控制系统对于高压转子转速回路和涡轮落压比回路具有伺服跟踪性能和抗回路耦合干扰性能；高压转子转速NH的阶跃突变对涡轮落压比πT超调量约1.3%的干扰，动态调节时间不大于2 s，进入稳态后，能够无静差伺服跟踪参考指令。

（2）对双转子涡扇发动机单回路控制的仿真表明：不考虑执行机构动态直接进行控制器设计，其相角裕度将减少30毅左右，导致系统的动态性能和稳定性变差。