李 博,田 丰,喻勇涛,郝彬彬
(1. 沈阳航空航天大学自动化学院,辽宁 沈阳 110136;2. 中国航发沈阳发动机研究所,辽宁 沈阳 110015)
随着当今科技发展水平的进步,在军事领域,武器装备的科技含量越来越高,尤其体现在航空领域,如机动性强,航程远的固定翼战机;亦如可以垂直起降,可在空中悬停作业的旋翼飞机[1]。随着各国短距起飞/垂直降落类型战斗机的问世,已将上述两种类型飞机的优点融为一体,即能发挥出很强的机动性,又能降低对起飞降落场地的需求,使军事上的战略战术更加灵活多变[2]。
战斗机科技含量提高的同时,也对战斗机设计和制造提出了更高的要求。综合考虑飞机与发动机的性能、架构等方面影响,将飞机与发动机进行一体化设计是一种有效的选择,相比于传统的飞机与发动机通常采用分离建模,再进行组合测试的设计研发体系,飞/发一体化设计能够充分考虑飞机与发动机之间的耦合特性,将飞机与发动机之间的相互影响降至最低,以至于将飞机的整体性能达到最高[3-4]。
本文以F-35B战机配备的F-135发动机为背景,根据飞机与发动机之间的耦合关系,建立基于Simulink的飞/发一体化模型。
F-35B战机既拥有普通固定翼战机机动性强,航程远的强大战斗力,同时也具备能够短距起飞/垂直降落的特点,这种特点主要源自于F-135发动机独特的结构,如图1。
图1 F-135发动机结构示意图
F-135发动机主要由升力风扇、三轴承喷管以及滚转喷管组成。其中三轴承喷管位于发动机后方,可向下偏转;滚转喷管位于发动机两侧,机翼下方,另外从主发动机低压轴引出驱动轴,用来驱动升力风扇。
F-35B战机拥有三种飞行状态,短距起飞/垂直降落状态、巡航飞行状态和悬停飞行状态。当飞机处于巡航飞行状态时,升力风扇停止工作,三轴承喷管向后偏转,即图1中巡航位置,提供巡航飞行过程中的全部推力;短距起飞/垂直降落状态以及悬停飞行状态时,三轴承喷管的三段筒体通过相互交汇处的机械结构进行旋转,进而改变推力矢量的方向,向下偏折90°,即图1中垂直推力位置,与升力风扇一起提供飞行所需的大部分升力,滚转喷管从主发动机中引气,为飞行提供一部分升力,并且负责调节飞机的横向姿态,保持机体稳定[5]。
飞/发一体化建模能使模型模拟其动态、静态特性,并且能为飞机整体综合控制提供便利[6]。本模型将分为两个部分,分别是发动机动力模块和机体模块。
发动机动力模块包括三轴承喷管、升力风扇和滚转喷管,各部件是否开启及其推力矢量大小和方向均可独立控制。
发动机模块输入量包括:主风扇物理转速nzon、升力风扇物理转速nlift、高度H、速度V0(地面坐标系);输出量包括:升力风扇推力Flift、三轴承喷管推力Fzon、左滚转喷管推力FgzL、右滚转喷管推力FgzR。
3.1.1 三轴承喷管模型设计
当飞机飞行高度处于11000米以下时,可由飞机当前高度H通过式(1)、(2)分别计算出此时的Ts0和Ps0
Ts0=288.15-0.0065H
(1)
Ps0=101325(1-H/44331)5.25588
(2)
结合当前飞机所处高度的大气压力Ps0和大气温度Ts0,根据绝热空气指数k1,以及当前飞机马赫数M0,可通过式(4)、(5)计算出进口气流总温T1、进口气流总压P1,其中马赫数M0由飞机速度V0通过式(3)计算得出
M0=V0/340
(3)
(4)
(5)
由主风扇物理转速nzon,进口气流静温T1和悬停设计点进口总温T1d,利用式(6)得到主风扇相对换算转速nzon,c
(6)
文献[7]中提及进口导流叶片角度θ1,主风扇进出口增压比πzon和主风扇换算转速nzon,c可以通过三维插值法计算求得主风扇的效率ηzon以及主风扇换算流量qm,a,c。由主风扇换算流量qm,a,c,进口气流总温T1以及进口气流总压P1通过式(7)可以计算得出主风扇实际流量qm,a,再加上主燃油流量qm,f,即为发动机全部输出流量,如式(8)所示,其中一部分为滚转喷管出口流量qm,gz,余下部分为三轴承喷管实际流量qm,z。
(7)
qm,z=qm,a+qm,f-qm,gz
(8)
为了方便计算,简化发动机内部散失的热量与压力,令三轴承喷管的总压、总温近似等于主风扇出口。由上文得到的进口气流总温T1以及进口气流总压P1可以通过式(9)式(10)计算得到主风扇出口总温T2和总压P2,k1为绝热空气指数,πzon为主风扇进出口增压比,ηzon为主风扇的效率。
(9)
P2=P1πzon
(10)
通过式(11)计算得出主风扇出口静压Ps2,利用所得Ps2通过式(12)得到喷管出口马赫数Mazon。马赫数Mazon再经过式(13)计算得到出口气流速度Vzon,其中Cp2为比定压热容比。
(11)
(12)
(13)
将上述求得的数值带入式(14)中,即可求得三轴承喷管推力Fzon,其中Azon为三轴承喷口横截面积。
Fzon=qm,z(Vzon-V0)+(Ps2-Ps0)Azon
(14)
3.1.2 升力风扇模型设计
升力风扇是实现垂直起降的重要升力来源,在悬停或起飞降落过程中可为飞机提供最大约80KN的升力。当机体产生俯仰运动时,可通过叶栅调节推力方向,平衡机体。
(15)
(16)
(17)
(18)
(19)
(20)
(21)
(22)
(23)
将上述求得的各数值带入式(24)中,即可求得升力风扇推力Flift,其中Alift为升力风扇出口横截面积。
(24)
3.1.3 滚转喷管模型设计
滚转喷管位于机体两侧,机翼正下方,当机体产生滚转运动时,从发动机外涵道引气,为滚转喷管提供动力,通过改变喷管角度,调节机体平衡[10]。其中涉及的中间变量包括主风扇实际流量qm,a、出口总温T2、出口总压P2、静压Ps0。Cgz为滚转喷管引气比,与上文计算出的主风扇实际流量即内外涵道总流量qm,a通过式(25)可以计算出滚转喷管出口流量[11]。
qm,gz=qm,aCgz
(25)
由于滚转喷管与三轴承喷管高度几乎相同,由式(2)可知,二者所处位置的大气静压也同为Ps0,取Vgz近似等于Vzon。由于滚转喷管需要从外涵道引气,因此滚转喷管的进口总压就等于主风扇出口总压P2。文献[12]提及σ15可由滚转喷管进口总压P2和动压Pd2通过神经网络映射函数计算得到,具体过程不做赘述,通过式(26)即可求出滚转喷管出口总压P3。
P3=σ15P2
(26)
将上述求得的各个数值带入式(27)中,即可求得滚转喷管推力FgzL和FgzR,其中Agz滚转喷管出口横截面积。
Fgz=0.5qm,gz(Vgz-V0)+0.5(P3-Ps0)Agz
(27)
机体模块是飞发一体化建模的重要组成部分,用于模拟在不同大小和方向的矢量推力的作用下,飞机所处的位置及姿态。通过与发动机模块的数据交流,可以计算出机体在各个方向的力和绕各轴旋转的力矩,并采用simulink软件中的6DOF工具计算出机体的运动状态与飞行姿态。
机体模块输入量包括:升力风扇推力Flift、三轴承喷管推力Fzon、滚转喷管推力FgzL、俯仰角θ、滚转角φ;输出量包括:速度Vex、Vey、Vez,位置Xex、Xey、Xez,机体姿态(滚转角φ、俯仰角θ、偏航角ψ),机体坐标系与地面坐标系变换矩阵DCMbe,机体速度Vbx、Vby、Vbz,机体角速度ωbp、ωbq、ωbr,机体角加速度dωbp/dt、dωbq/dt、dωbr/dt,机体加速度Abx、Aby、Abz;中间量包括沿三个轴的力Fx、Fy、Fz,以及绕三个轴旋转的力矩L、M、N,其中机体坐标轴以机体质心为原点,质心到机体前方为X轴正方向,质心到机体右侧为Y轴正方向,过质心且垂直于X轴与Y轴向下为Z轴正方向;地面坐标系与最初时刻机体坐标系完全重合,当机体运动随机体姿态改变时,地面坐标系不变。
由于整个飞机受力情况复杂,采用综合分析将十分繁琐,因此需要先将各部件推力分解在地面坐标系中分解,再通过simulink中的6DOF工具计算出机体的运动状态与运行姿态。
3.2.1 三轴承喷管推力分解
三轴承喷管在运行中,既可以上下偏转,调节俯仰姿态,也可以左右偏转,调节偏航姿态,因此可以在地面坐标系中将三轴承喷管推力分解为Fx,zon、Fy,zon、Fz,zon三个方向的推力与绕Y、Z轴旋转的力矩Mzon、Nzon。
当三轴承喷管左右偏折角度为λzon,ψ时,将推力Fzon在机体坐标轴中分解为Fy,zon、Fxz,zon。
Fy,zon=Fzonsinλzon,ψ
(28)
Fxz,zon=Fzoncosλzon,ψ
(29)
(30)
(31)
(32)
(33)
设三轴承喷管喷口至质心距离为xzon,则绕Z轴力矩Nzon与绕Y轴力矩Mzon如下:
Nzon=Fy,zonxzon
(34)
(35)
3.2.2 升力风扇推力分解
升力风扇推力Flift方向只受叶栅角度λlift和机体俯仰角θ影响,只在俯仰面内,即XOZ面内运动,因此可以将推力在地面坐标系中分解为X,Z轴方向的力Fx,lift、Fz,lift,与绕Y轴的力矩Mlift。
(36)
(37)
(38)
(39)
设升力风扇至质心距离为xlift,则绕Y轴力矩Mlift如下
(40)
3.2.3 滚转喷管推力分解
由于滚转喷管只能调节滚转姿态,因此可以将滚转喷管推力在地面坐标轴中分解为绕X轴的力矩Lgz和X、Y、Z轴方向的力Fx,gz、Fy,gz、Fz,gz。
假设机体此时左右滚转喷管的偏转角度分别为λgzL、λgzR,将FgzL与FgzR在机体坐标轴中分解。
Fy,gz=FgzLsinλgzL+FgzRsinλgzR
(41)
(42)
(43)
(44)
设左右滚转喷管至质心距离为xgzL和xgzR,则绕X轴力矩Lgz如下
LgzL=FgzLcosλgzLxgzL
(45)
LgzR=FgzRcosλgzRxgzR
(46)
Lgz=LgzR-LgzL
(47)
3.2.4 机体姿态解算
设飞机整体质量为m,重力加速度g取9.8,由式(48)可知机体所受重力为Fg,且重力方向不随机体运动而变化,始终沿地面坐标轴Z轴向下。
Fg=mg
(48)
由各部件推力分解可通过式(49)-(54)求出机体整体受力及绕各轴旋转的力矩。
Fx=Fx,zon-Fx,lift-Fx,gz
(49)
Fy=Fy,gz-Fy,zon
(50)
Fz=Fz,g-Fz,zon-Fz,lift-Fz,gz
(51)
L=Lgz
(52)
M=Mlift+Mzon
(53)
N=Nzon
(54)
将上述Fx、Fy、Fz与L、M、N作为输入,通过simulink中的6DOF工具,可以计算出机体的运动状态与飞行姿态,6DOF工具如图2。
图2 6DOF输入与输出
根据文献[12],以该模型处于H=9.7km处,M0=1.2时的特性作为本模型的设计点参数,其对应的模型仿真结果为三轴承喷管推力Fzon=78100N,升力风扇推力Flift=81100N,左右滚转喷管推力FgzL=FgzL=8940N,具体部件设计点参数见表1。
表1 模型设计点参数
为了确保本模型的正确性与有效性,将对模型进行数值仿真。将模型输入调整为高度H=9.7km,马赫数M0=1.2,即机体悬停在9.7千米高空时的状态,且姿态角以及各喷管角度全部置零,令初始转速n=7000r/min,在1000秒内,转速匀速提升至10000r/min,观察三轴承喷管、升力风扇和滚转喷管的推力、机体模块中各轴的分力以及绕各轴旋转的力矩的变化,如图3至图7。
图3 三轴承喷管推力仿真曲线
图4 升力风扇推力仿真曲线
图5 滚转喷管推力仿真曲线
图6 分解推力仿真曲线
图7 力矩仿真曲线
随着发动机转速从7000r/min匀速提升至10000r/min,如图3至图5,三个部件的推力也明显提升。各推力经过机体模块分解后,如图6,由于各部件偏转角度全部为0°,因此升力风扇和滚转喷管推力垂直向下,产生沿Z轴负方向的升力Fz,三轴承喷管推力水平向后,产生沿X轴正方向的推力Fx,其中升力风扇的升力Flift产生如图7的绕Y轴旋转的力矩M。
令马赫数M0=1.2,转速n=8500r/min固定不变,将高度H在1500秒内从0km匀速提升至18km,将三轴承喷管推力Fzon仿真结果与外文文献[12]的研究数据结果进行比较,变化趋势基本相同,确保了动力模型的正确性与有效性,具体数据见表2。
表2 三轴承喷管推力数据对比
为验证一体化模型的有效性,将对飞机从0m垂直爬升至1000m并保持稳定的过程进行仿真验证。首先将转速设置为8500 r/min,其它各控制量、飞机位置与运行姿态等输入量的初始值全部设置为0,并开始运行1500s,并观察飞行高度、俯仰角、滚转角以及偏航角来判断飞机是否成功垂直爬升至1000m,并保持机体稳定,如图8至图11。
图8 高度仿真曲线
图9 俯仰角仿真曲线
图10 滚转角仿真曲线
图11 偏航角仿真曲线
飞机高度在250s左右达到1000m,并在750s左右趋于稳定高度;俯仰姿态在起飞时产生微小偏差,仅为-0.3°,并且迅速调整进入平稳状态;滚转姿态由于受左右喷管交替作用,产生持续等幅振荡,但幅度极小,仅为0.75°,并不影响机体整体运行状态;偏航角起飞阶段处于减幅振荡,偏差迅速削减,很快达到平稳状态,完成了使机体平稳上升至1000m高度并保持稳定的目标。验证了本文飞/发一体化建模的有效性。
本文基于F-135发动机,进行了针对于短距起飞/垂直降落的飞/发一体化建模研究。建立了由升力风扇、三轴承喷管及滚转喷管组成的发动机模块模型,为飞行提供动力来源;建立了机体模块模型,将发动机模块产生的力进行分解并计算力矩,利用simulink中的6DOF工具,计算机体的运动状态和运行姿态;对模型输出进行仿真,将仿真结果与文献数据进行对比,确保了模型的正确性和有效性。