2. 北京航空航天大学 能源与动力工程学院,北京 100191;
3. 先进航空发动机协同创新中心,北京 100191
2. School of Energy and Power Engineering, Beihang University, Beijing 100191, China;
3. Collaborative Innovation Center of Advanced Aero-Engine, Beijing 100191, China
高空台高空模拟腔是进行发动机全飞行包线功能、性能模拟的核心试验装置,能够进行大、中等推力发动机试验的试验设备,对于新建高空台设备需要进行可行性方案的论证,以避免试验过程中出现的不可预测风险。
据公开资料报道的美国阿诺德工程发展中心(AEDC)的J2高空试验台,通过采用掺混调温的方法实现了高空模拟腔的温度、压力同时控制[1~5],但没有公开相关控制器的设计方法;德国斯图加特高空台采用自适应控制实现了高空模拟腔的温度、压力同时控制[6],加拿大渥太华运用前馈控制器和解耦控制器实现了高空模拟腔温度、压力的同时控制[7, 8],但与美国阿诺德高空试验台相比,温度回路的控制效果较差。
国内20世纪80年代建立的高空台是继美、俄、英、法4国之后拥有此类大型试验设施的国家,其测量精度与俄罗斯、英国高空台相当[9],目前在建高空台采用了特种流量调节阀,由于试验条件的限制无法获得精准的输入输出特性,由此带来了系统建模的不确定性,这对于方案论证中需要建立系统模型,并基于模型展开控制器的设计及仿真验证带来设计上的难度。
高空模拟腔的设计要求为在飞行环境条件的干扰下能够伺服跟踪飞行环境的参考指令。传统的控制系统设计均采用单变量的PID结构进行调节,然而,这种方案是无法保证上述飞行环境模拟腔的设计要求,飞行环境模拟腔的控制参数为模拟腔温度T和压力p,这二个参数在进行发动机高空性能试验时是同时变化的。
μ综合控制由于其理论是针对不确定性被控系统的鲁棒控制设计的保守性太强而导致控制性能下降的问题而发展起来的,近年来受到工程界的广泛重视[10, 11],它的核心思想是采用了D-K迭代算法对变尺度D、K矩阵进行μ综合结构奇异值的最优寻优[12~16],在整个频域工作范围获得较为均匀的H∞抗干扰性能,以此克服保守性。
针对高空模拟腔温度、压力同时控制的论证方案,本文采用系统建模的方法,考虑了系统未建模动态不确定性,建立了由高空模拟腔、特种流量阀调节阀执行机构组成的增广被控对象不确定性模型,提出的一种积分型μ综合控制器的设计方案,通过对温度、压力回路采用不同的性能加权函数,并考虑了对控制器输出幅值采用分频加权的设计思想,采用D-K迭代算法设计了μ综合控制器,并以涡扇发动机仿真试验对算法进行了仿真验证。
2 系统不确定性建模设国内高空台高空模拟腔控制系统结构如图 1所示。
图中r为飞行环境模拟指令,e为偏差信号,u为闭环控制器输出信号,v为执行机构输出信号,y为高空模拟腔输出信号。
2.1 高空模拟腔模型高空模拟腔结构简图如图 2所示,该腔有两路进气和一路排气,其中第一路是温度为Tin1的高温热气流,第二路是温度为Tin2的低温冷气流,分别由1号、2号调节阀来控制,其温度和压力非线性微分方程为
$ \begin{array}{l} \frac{{{\rm{d}}T}}{{{\rm{d}}t}} = \frac{{RT}}{{VP\left( {{c_p} - R} \right)}}\left[ { - \left( {{h_{{\rm{out}}}} - RT} \right)\left( {{W_{{\rm{in1}}}} + {W_{{\rm{in2}}}} - {W_{{\rm{out}}}}} \right) + } \right.\\ \;\;\;\;\;\;\;\;\left. {\left( {{h_{{\rm{in1}}}}{W_{{\rm{in1}}}} + {h_{{\rm{in2}}}}{W_{{\rm{in2}}}} - {h_{{\rm{out}}}}{W_{{\rm{out}}}}} \right)} \right] \end{array} $ | (1) |
$ \begin{array}{l} \frac{{{\rm{d}}p}}{{{\rm{d}}t}} = \frac{R}{V}\left[ {\left( {T - \frac{{{h_{{\rm{out}}}} - RT}}{{{c_p} - R}}} \right)\left( {{W_{{\rm{in1}}}} + {W_{{\rm{in2}}}} - {W_{{\rm{out}}}}} \right) + } \right.\\ \left. {\frac{R}{{V\left( {{c_p} - R} \right)}}\left( {{h_{{\rm{in1}}}}{W_{{\rm{in1}}}} + {h_{{\rm{in2}}}}{W_{{\rm{in2}}}} - {h_{{\rm{out}}}}{W_{{\rm{out}}}}} \right)} \right] \end{array} $ | (2) |
其中,参数变量物理意义参见表 1。
定义状态空间模型的状态向量为
当发动机工作在某一节流状态时,利用高空模拟腔非线性模型获得稳态工作点数据如表 1所示。
在上述稳态工作点对模型进行线性化,可得线性模型为
$ \begin{array}{l} \dot x = \mathit{\boldsymbol{A}}x + {\mathit{\boldsymbol{B}}_1}u + {\mathit{\boldsymbol{B}}_2}d\\ y = \mathit{\boldsymbol{C}}x \end{array} $ | (3) |
$ \mathit{\boldsymbol{A}} = \left[ {\begin{array}{*{20}{c}} {\frac{{\partial \dot T}}{{\partial T}}}&{\frac{{\partial \dot T}}{{\partial p}}}\\ {\frac{{\partial \dot p}}{{\partial T}}}&{\frac{{\partial \dot p}}{{\partial p}}} \end{array}} \right],{\mathit{\boldsymbol{B}}_1} = \left[ {\begin{array}{*{20}{c}} {\frac{{\partial \dot T}}{{\partial {W_{{\rm{in1}}}}}}}&{\frac{{\partial \dot T}}{{\partial {W_{{\rm{in2}}}}}}}\\ {\frac{{\partial \dot p}}{{\partial {W_{{\rm{in1}}}}}}}&{\frac{{\partial \dot p}}{{\partial {W_{{\rm{in2}}}}}}} \end{array}} \right] $ |
$ {\mathit{\boldsymbol{B}}_2} = \left[ {\begin{array}{*{20}{c}} {\frac{{\partial \dot T}}{{\partial {W_{{\rm{out}}}}}}}&{\frac{{\partial \dot T}}{{\partial {T_{{\rm{in1}}}}}}}&{\frac{{\partial \dot T}}{{\partial {T_{{\rm{in2}}}}}}}\\ {\frac{{\partial \dot p}}{{\partial {W_{{\rm{out}}}}}}}&{\frac{{\partial \dot p}}{{\partial {T_{{\rm{in1}}}}}}}&{\frac{{\partial \dot p}}{{\partial {T_{{\rm{in2}}}}}}} \end{array}} \right],\mathit{\boldsymbol{C}} = \left[ {\begin{array}{*{20}{c}} 1&0\\ 0&1 \end{array}} \right] $ |
按表 1建立的稳态工作点线性模型如下
$ \mathit{\boldsymbol{A}} = \left[ {\begin{array}{*{20}{c}} { - 0.8403}&0\\ { - 295.4755}&0 \end{array}} \right],{\mathit{\boldsymbol{B}}_1} = \left[ {\begin{array}{*{20}{c}} {0.5478}&{0.0898}\\ {428.9472}&{267.8978} \end{array}} \right], $ |
$ {\mathit{\boldsymbol{B}}_2} = \left[ {\begin{array}{*{20}{c}} \begin{array}{l} - 0.2693\\ - 331.0314 \end{array}&\begin{array}{l} 0.3315\\ 116.5715 \end{array}&\begin{array}{l} 0.5102\\ 179.4130 \end{array} \end{array}} \right],\mathit{\boldsymbol{C}} = \left[ {\begin{array}{*{20}{c}} 1&0\\ 0&1 \end{array}} \right] $ |
调节阀执行机构模型由两部分组成:从控制指令输出到液压缸活塞位移的传递函数和轮盘式特种阀流量模型的传递函数。
(1) 液压缸模型
带执行机构的调节阀在图 2中用1号和2号调节阀表示,结构上为轮盘式特种阀,如图 3所示。该特种阀的直径为2m,在整体上分为三个部分:前壳体、中壳体和后壳体,由大功率液压缸来驱动。
根据实际液压缸的特性数据采用辨识方法获得液压缸由给定信号到输出信号的传递函数为
$ {g_{\rm{o}}} = \frac{{ - 1.333s + 2.665}}{{s + 2.664}} $ | (4) |
(2) 特种阀模型
轮盘式特种阀的流量计算公式为
$ W = \varphi mA\sqrt {2\rho {p_1}} $ | (5) |
式中W为阀的空气流量,m为特种阀流通面积与管道截面积之比,A为管道截面积,ρ为空气的密度,p1为阀前压力,p2为阀后压力,
$ \varphi = f\left( {m,Pr} \right),Pr = \frac{{{p_2}}}{{{p_1}}} $ | (6) |
为了运用线性系统理论进行控制器设计,需要对上述特种阀的非线性模型进行线性化,根据该特种阀在对应的表 1稳态工作点的数据,可得1号和2号调节阀液压缸位移到特种阀流量的增益传递函数分别为
$ {K_{V1}} = 1822.6,{K_{V2}} = 2311.5 $ | (7) |
另外,根据计算,在全飞行包线工作范围内KV1和KV2的变化在1300~2700,相对于稳态设计点,其乘性不确定性可按相对误差计算,KV1存在约48%的乘性不确定性Δ1,KV2存在约43.7%的乘性不确定性Δ2。
根据上述液压缸传递函数和调节阀传递函数可得调节阀执行机构模型为
$ {\mathit{\boldsymbol{G}}_{{\rm{act}}}} = \left[ {\begin{array}{*{20}{c}} {{g_0} \times {K_{V1}}\left( {1 + {\mathit{\Delta }_1}} \right)}&0\\ 0&{{g_0} \times {K_{V2}}\left( {1 + {\mathit{\Delta }_2}} \right)} \end{array}} \right] $ | (8) |
将高空模拟腔模型和调节阀执行机构模型进行增广构成增广被控对象,其频域特性由图 4所示。
由图 4(a)可知,温度控制回路开环传函截止频率为147Hz,由图 4(b)可知,压力控制回路开环传函截止频率为71.4kHz,二者的频响特性相差数百倍。
2.4 发动机稳态简化模型相对于研究对象的高空模拟腔,被试验的发动机进口空气流量对其直接构成了干扰,在设计过程中应加以考虑。
为模拟发动机在整个飞行包线内工作的情况,高空试验台装配了飞行高度H,马赫数Ma,油门杆角度PLA的操作装置,试验过程中通过改变H,Ma生成高空模拟腔的指令,同时操作员通过改变PLA模拟发动机在高空的工作情况,显然,发动机的实际进口空气流量可用H,Ma,PLA这三个参数描述。
由发动机控制计划可以获得低压转子转速N1相对于PLA的函数关系,可得
$ {N_1} = f\left( {PLA} \right) $ | (9) |
由发动机工作原理,换算转速为
$ {N_{1,{\rm{cor}}}} = f\left( {H,Ma,{N_1}} \right) $ | (10) |
根据发动机稳态工作原理,当喷口固定不变条件下,在压气机特性图上存在确定的稳态共同工作线,由该工作线与等换算转速线的交点即可获得换算流量如下
$ {W_{{\rm{a}},{\rm{cor}}}} = f\left( {{N_{1,{\rm{cor}}}}} \right) $ | (11) |
再根据发动机相似工作原理可获得高空飞行条件下的发动机实际进口空气流量
$ {W_{\rm{a}}} = f\left( {H,Ma,{W_{{\rm{a}},{\rm{cor}}}}} \right) $ | (12) |
考虑高空模拟腔双变量的伺服设计要求,针对图 1控制系统结构提出以下如图 5所示的积分型μ综合控制的设计方案。
图 5中G为高空模拟腔模型,Gact为调节阀执行机构模型,K为μ综合控制器,Wu为对控制量的加权函数,Wp为伺服性能加权函数,Int为对角型积分模块,不同物理参数{1-x}表示了变量的维数。
3.2 性能加权函数设计在前述的增广被控对象双回路开环频响特性中可以看到,截止频率相差很大,为此,对温度、压力回路采用不同加权函数的方法是解决温度/压力同时控制的有效途径,以实现后双回路的同步作用效果。
考虑图 5的积分型μ综合闭环系统,对温度伺服性能的加权含积分前加权和积分后加权,同理对压力伺服性能的加权也含积分前加权和积分后加权,首先,对于二路积分前的加权函数分别设计如下:
(1) 对温度回路偏差的加权函数设计为
$ {W_{p11}} = {10^4}\frac{{\frac{s}{3} + 1}}{{\frac{s}{{0.0002}} + 1}} $ | (13) |
(2) 对压力回路偏差的加权函数设计为
$ {W_{p22}} = {10^4}\frac{{\frac{s}{{0.04}} + 1}}{{\frac{s}{{0.0000025}} + 1}} $ | (14) |
(3) 对温度回路偏差积分的加权函数设计为
$ {W_{p33}} = 2 \times {10^4}\frac{{\frac{s}{3} + 1}}{{\frac{s}{{0.0002}} + 1}} $ | (15) |
(4) 对压力回路偏差积分的加权函数设计为
$ {W_{p44}} = {10^5}\frac{{\frac{s}{{0.04}} + 1}}{{\frac{s}{{0.00000025}} + 1}} $ | (16) |
则,性能加权函数为
$ {W_p} = \left[ {\begin{array}{*{20}{c}} {{W_{p11}}}&0&0&0\\ 0&{{W_{p22}}}&0&0\\ 0&0&{{W_{p33}}}&0\\ 0&0&0&{{W_{p44}}} \end{array}} \right] $ | (17) |
Wp的幅频响应如图 6所示。
采用低频放开限制、中频逐步限制、高频最大限制的原则设计控制量加权函数
$ {W_{u11}} = {W_{u22}} = 0.8\frac{{\frac{s}{{0.025}} + 1}}{{\frac{s}{{400}} + 1}} $ | (18) |
则,构造的双回路控制量加权函数为
$ {W_u} = \left[ {\begin{array}{*{20}{c}} {{W_{u11}}}&0\\ 0&{{W_{u22}}} \end{array}} \right] $ | (19) |
Wu幅频响应如图 7所示。
经过上述处理,至此,对图 6中进行μ综合设计前的工作已准备完毕,按下述算法构建系统的内联模块
systemnames= 'Gact G Wp Wu Int ';
inputvar= '[ref{2}; dist{3}; control{2}]';
outputvar= '[Wp; Wu; ref; -G; Int]';
input_to_Gact = '[control]';
input_to_G = '[Gact; dist]';
input_to_Wp = '[ref-G; Int]';
input_to_Wu = '[control]';
input_to_Int= '[ref(1)-G(1);ref(2)-G(2)]';
sys_ic= sysic;
并给定测量输出的数目nmeas=6和控制器输出ncont=2,采用D-K迭代算法(dksyn函数)求解,其迭代过程如表 3所示。
上述求解获得一个阶次为14阶的μ综合控制器。
温度和压力回路的开环伯德图如图 8所示,由图可知,温度回路在1.73rad/s频率处的幅值裕度为8.88dB,在0.716rad/s频率处的相角裕度为36.6°;压力回路在5.32rad/s的幅值裕度为26.1dB,在0.574 rad/s的相角裕度为71°。
温度回路和压力回路的灵敏度函数和补灵敏度函数的幅频特性图如图 9所示,由图 9(a)可知,在低频段灵敏度函数的幅值远小于0dB,系统具有抗干扰性能和伺服跟踪性能,在高频段补灵敏度函数每十倍频程幅值衰减-40dB,具有高频抗噪声的抑制能力。由图 9(b)所示,压力回路频域响应特性与温度回路类似,同样具有抗干扰性能和伺服跟踪性能以及对高频噪声的抑制能力。
为了验证算法的有效性,以亚声速大涵道比民用涡扇发动机工作状态点来设定仿真试验条件,如图 10所示。并假定在仿真试验过程中,高空模拟腔两路进气阀前的温度和压力保持不变。
对图 10发动机试验条件的描述:
A~B:在0时刻,H为5km,Ma为0.4,PLA处于慢车状态,保持10s。
B~C:从10s开始,其后20s内,H不变,Ma从0.4增到0.8,PLA由慢车变化到巡航状态。
C~D:从30s到80s,保持H为5km、Ma为0.8,PLA为45°这一试验状态不变。
D~E:从80s开始,Ma保持0.8不变,H在25s内由5km上升到8km,PLA在1s内快速地由巡航状态变化到最大状态,之后保持到仿真结束。
同时,发动机进口流量的变化曲线如图 11所示。
在上述条件下运行,高空模拟腔压力的仿真结果如图 12所示。由图 12可知,从10s~30s,H保持不变,Ma逐渐增大,发动机由慢车过渡到巡航状态过程中,发动机进口流量逐渐增大,发动机进口压力逐渐增大,高空模拟腔的压力始终能跟踪发动机进口压力的变化,从80s~105s,H逐渐升高,Ma保持不变,发动机由巡航过渡到最大状态过程中,发动机进口流量逐渐增大,发动机进口压力逐渐减小,模拟腔的压力始终跟踪进口压力的变化,同样,其它时间段模拟腔的压力均能跟踪进口压力的变化,而且整个仿真过程中,压力的相对误差在3%以内,这说明μ综合控制器具有对发动机飞行包线内进口压力的伺服跟踪性能,同时具有对发动机稳态和过渡态变化时的抗干扰性能。高空模拟腔温度的仿真结果如图 13所示。由图 13可知,从10s~30s,H保持不变,Ma由0.4增大到0.8的过程中,模拟腔的温度始终跟踪进口压力的变化,由局部放大图可知,温度的最大偏差为1K,同样,其它时间段模拟腔的温度均能跟踪进口温度的变化,而且整个仿真过程中,温度的相对误差在1%以内,这说明μ综合控制器具有对发动机飞行包线内进口温度的伺服跟踪能力,同时也具有对发动机进口流量变化的抗干扰性能。
高空模拟腔两路空气进气量的控制作用曲线如图 14所示,由图可知,当发动机进口流量变化时,二路空气进气量做出了快速的变化,以快速调节闭环偏差,从而保证了模拟腔压力和温度对发动机进口压力和温度的伺服跟踪。模拟腔两路进气调节阀阀位的控制作用曲线如图 15所示,图中红色虚线为1号调节阀控制曲线,黑实线为2号调节阀的控制曲线。
本文提出的一种积分型μ综合控制方法实现了高空模拟腔的温度、压力伺服控制,得到如下结论:
(1) 对温度、压力回路采用不同的性能加权函数和对控制器输出幅值进行加权函数的设计,采用D-K迭代算法设计的μ综合控制器,具有双回路控制的伺服性能和抗干扰鲁棒性能。
(2) 以某涡扇发动机的仿真试验为例,对高空模拟腔积分型μ综合控制方法进行了仿真验证,仿真结果表明,高空模拟腔温度仿真的最大偏差为1K,压力的相对误差在3%以内。
[1] |
Peter A Montgomery, Rick Burdette, Brian Krupp. A Real-Time Turbine Engine Facility Model and Simulation for Test Operations Modernization and Integration [R]. ASME 2000-GT-0576.
(0) |
[2] |
Davis M W Jr. Analysis of Aeropropulsion Test Facility Aerodynamic and Structural Issues Using Compression System Numerical Simulations[R]. ASME 2014-GT-25068.
(0) |
[3] |
Peter A Montgomery, Rick Burdette, Larry Wilhite, et al. Modernization of a Turbine Engine Test Facility Utilizing a Real-Time Facility Model and Simulation[R]. ASME 2001-GT-0573.
(0) |
[4] |
Peter A Montgomery, Rick Burdette. Evolution of a Turbine Engine Test Facility to Meet the Test Needs of Future Aircraft Systems[R]. ASME GT 2002-30605.
(0) |
[5] |
Milt Davis, Peter Montgomery. A Flight Simulation Vision for Aeropropulsion Altitude Ground Test Facilities[J]. Transactions of the ASME, 2005, 127(1).
(0) |
[6] |
Schmidt K, Merten R, Menrath M, et al. Adaptation of the Stuttgart University Altitude Test Facility for BR700 Core Demonstrator Engine Tests[R]. ASME 98-GT-556.
(0) |
[7] |
Gary M Elfstrom. History of Test Facility Design Expertise at Aiolos Engineering Corporation[J]. AIAA Journal, 2007, 45.
(0) |
[8] |
Majid Borairi. Design and Commissioning of a Multivariable Control System for a Gas Turbine Engine Test Facility[R]. AIAA 2006-3151.
(0) |
[9] |
Peter M Pachlhofer, Joseph W, Dennis J Dicki. Advance in Engine Test Capabilities at the NASA Glenn Research Center's Propulsion System Laboratory[R]. ASME 2006-GT-90181.
(0) |
[10] |
王曦, 姚华. 弹用涡喷发动机鲁棒性能μ综合控制[J]. 推进技术, 2003, 24(3): 236-239. (WANG Xi, YAO Hua. μ Synthesis Control with Robust Performance for Missile Turbojet Engine[J]. Journal of Propulsion Technology, 2003, 24(3): 236-239.)
(0) |
[11] |
王曦, 韩乃湘, 李喜发, 等. 航空发动机鲁棒H∞/PI状态反馈控制[J]. 推进技术, 2003, 24(4): 364-367. (WANG Xi, HAN Nai-xiang, LI Xi-fa, et al. Robust H∞/PI State Feedback Control for Aeroengine[J]. Journal of Propulsion Technology, 2003, 24(4): 364-367.)
(0) |
[12] |
Doyle J C. Analysis of Feedback Systems with Structured Uncertainties[J]. IEE Proc.Part D, Control Theory Appl, 1982, 129: 242-250. DOI:10.1049/ip-d.1982.0053
(0) |
[13] |
Doyle J C. A Review of μ: for Case Studies in Robust Control[C]. Munich: Proceedings of 10th IFAC World Congress, 1987.
(0) |
[14] |
Enns D F. Structured Singular Value Synthesis Design Example: Rocket Stabilization[C]. NewYork: American Control Conference, IEEE, 1990.
(0) |
[15] |
侯敏杰. 高空模拟试验技术[M]. 北京: 航空工业出版社, 2014.
(0) |
[16] |
Da-Wei Gu D C, Petkov P H, Konstantinov M M. Robust Control Design with MATLAB[M]. London: Springer London, 2005.
(0) |
[17] |
章卫国, 李爱军, 刘小雄, 等. 鲁棒飞行控制系统设计[M]. 北京: 国防工业出版社, 2012.
(0) |
[18] |
Green M, Limebeer D J N. Linear Robust Control[M]. Englewood Cliffs: Prentice Hall, 1995.
(0) |