超燃冲压发动机喷管推力性能理论预测

2022-02-16张子健马凯夫刘云峰

气体物理 2022年1期

韩信, 张子健, 马凯夫, 刘云峰

(1.中国科学院力学研究所, 北京100190; 2.中国科学院大学工程科学学院, 北京 100049; 3.香港理工大学航空工程跨领域学部, 中国香港 999077)

引言

超燃冲压发动机概念自提出以来, 有60多年研究历史(1958—2020年), 取得了很大的进展[1-9].但是, 如何评估超燃冲压发动机的推力性能是一个难题.美国X-43A氢燃料超燃冲压发动机Ma=7的飞行试验虽然只进行了11 s, 但测得有加速度, 然而Ma=9.8的飞行试验则未能测得明显的加速度[10], 说明Mach数越高, 净推力越难获得.美国空军X-51A碳氢燃料冲压发动机的第1次飞行试验虽然运行时间达200 s, 但是没有达到从Ma=4.75 加速到Ma>5的预计目标.第4次飞行试验虽然获得成功, 但是最终只加速到Ma=5.1.飞行试验结果说明, 碳氢燃料超燃冲压发动机比氢燃料超燃冲压发动机更难获得净推力.

除了数值模拟和试验测量手段以外, 对超燃冲压发动机的推力性能进行理论分析和评估, 找出影响推力的关键参数和物理规律就显得格外重要.目前主要有两类分析方法, 第一类方法是基于Tsien等在1949年提出的一维加热管流理论[11].Birzer等发展了氢燃料超燃冲压发动机燃烧室准一维模型[12], O′Brien等发展了有限化学反应速率的准一维模型, 可以预测燃烧室内的压力分布[13], Vanyai等开展了超燃冲压发动机理想燃烧过程的理论分析[14].第二类方法是基于Riggins等提出的推力势分析方法[15-16].Sislian等利用该方法分析了激波诱导燃烧冲压发动机的推力性能[17-20], 但是该方法也只能给出定性的结果.

超燃冲压发动机的燃烧室一般是等截面的, 燃烧后的流动仍然保持超声速流动.喷管采用扩张型喷管, 推力主要由高温高压的燃烧气体在尾喷管中膨胀做功产生.在本文中, 利用气体动力学等熵膨胀理论, 对理想气体在喷管中的膨胀做功过程进行了分析, 通过对喷管壁面的压力积分得到发动机的推力, 提出无量纲推力的理论表达式.通过理论表达式可以明显地发现影响发动机推力的关键参数和物理规律.与以前准一维理论研究方法不同的是, 本研究的特点是将超燃冲压发动机的燃烧室与喷管分离, 把燃烧室出口参数作为理论分析的边界条件, 得到一般化的物理规律.

1 理论分析

图1给出了亚燃冲压发动机和超燃冲压发动机的工作原理图[21].对于亚燃冲压发动机, 超声速气流经过进气道减速后, 变成亚声速气流, 在燃烧室内燃烧, 产生高温高压燃烧产物.高温高压燃烧产物经过Laval喷管膨胀, 产生推力.亚燃冲压发动机有两个推力部件, 一是燃烧室, 二是尾喷管.对于超燃冲压发动机, 超声速气流经过进气道减速后, 仍然保持超声速流动, 在燃烧室内燃烧产生高温高压气体, 经过喷管膨胀产生推力.由于高温高压燃烧产物仍然是超声速的, 超燃冲压发动机的喷管是扩张型喷管.隔离段和燃烧室几乎是等截面的, 因此, 超燃冲压发动机只有尾喷管一个推力部件.

图1 亚燃冲压发动机和超燃冲压发动机的示意图[21]

与亚燃冲压发动机相比, 超燃冲压发动机只有尾喷管一个推力部件, 因此, 只须研究尾喷管的流动特征即可得到发动机的推力性能.二者的不同在于, 亚燃冲压发动机的燃烧产物是从Ma=1开始膨胀, 而超燃是从Ma>1开始膨胀.因此, 超燃冲压发动机燃烧产物的Mach数是影响推力性能的一个关键参数.由于推力主要是由喷管壁面的压力积分产生, 因此, 喷管入口的静压是影响推力的另一个关键参数.接下来, 利用理想气体等熵膨胀关系式推导得到超燃冲压发动机尾喷管推力的无量纲关系式.

对于Laval喷管中的理想气体等熵膨胀过程, 气体压力和喷管面积随喷管出口Mach数的变化关系式分别见式(1)和(2)

(1)

(2)

其中,Ma是喷管出口Mach数,p和p0分别是静压和总压,A和A*分别是喷管出口截面积和喉道截面面积,γ是比热比.

对喷管壁面压力进行积分, 得到无量纲推力

(3)

其中,F*是特指气体从喉道(Ma=1)开始膨胀产生的推力,p*A*是喷管喉道的静压和面积的乘积,Ma是喷管出口Mach数,γ是气体比热比.

式(3)的计算结果见图2, 从图2可以看出, 无量纲推力随喷管出口Mach数的增加呈先增加后趋于平缓的趋势.当喷管出口Mach数较高时, 即Ma>5时, 无量纲推力超过80%, 继续增加喷管出口Mach数, 无量纲推力变化缓慢.此时, 继续增大喷管尺度, 不但无量纲推力增加不大, 反而会急剧增加喷管的体积和重量.

图2 无量纲推力与喷管出口Mach数的关系曲线

对于超燃冲压发动机, 因为高温高压的燃烧气体是从Ma>1开始膨胀, 因此, 根据式(3)推导得到的无量纲推力见式(4)

(4)

其中,F是超燃冲压发动机喷管产生的推力,Ma1和Ma2分别是超燃冲压发动机喷管入口和出口的Mach数.

对于一个给定的喷管, 求解式(4)须知道Ma2.假设是等熵膨胀, 那么Ma2是Ma1和喷管面积比A2/A1的函数, 见式(5)

(5)

因为推导过程中存在一个隐函数, 式(5)无法给出显式的解析解, 但是可以通过数值方法求解.

式(4)给出的推力是通过喷管喉道参数p*A*进行无量纲化的.在实际工程应用中, 须用燃烧室出口参数或喷管入口参数p1A1进行无量纲化, 二者关系见式(6)

(6)

在实际的工程应用中, 发动机的喷管尺度是固定不变的, 而燃烧室出口或喷管入口的Mach数是变化的.用数值方法对A2/A1=5, 10, 15, 20, 25的无量纲推力进行了求解, 求解结果见图3.上述面积比大约对应理想气体Laval喷管Ma=3, 4, 4.5, 4.7, 5的状态, 基本覆盖了超燃冲压发动机喷管的尺度范围.注意: 图3中的推力是用喷管入口参数p1A1进行无量纲化的.

从图3中可以看出几个明显的流动规律.首先, 对于固定的喷管面积比, 无量纲推力随着Ma1的增加, 先增加而后趋于平缓, 到Ma1=2时, 变化非常缓慢.这是因为Ma1不同, 导致喷管Ma2不同, 无量纲推力发生变化.其次, 随着喷管面积比的增大, 无量纲推力增大, 但是趋于收敛.当A2/A1>15以后, 变化很小, 此时对应的大约是Laval喷管Ma=4.5的状态, 从图1已知, 此时已经获得了大于80%的无量纲推力.继续增大面积比, 无量纲推力增加不大.

图3 超声速燃烧喷管产生的无量纲推力与入口Mach数的关系

为了得到更直观的理论公式, 对式(3)和(6)的理论求解结果进行拟合(假设γ=1.3), 得到拟合公式, 拟合曲线见图4和5, 拟合公式见式(7)和(8), 最终得到拟合后的无量纲推力式(9)和(10).

图4 喷管无量纲推力的拟合结果

(7)

(8)

(9)

(10)

从式(10)可以明显看出, 超燃冲压发动机的无量纲推力主要受两个参数影响, 一是喷管入口Mach数, 随着Mach数的增大, 对发动机的无量纲推力是有利的, 说明高Mach数超燃冲压发动机和斜爆轰冲压发动机在推力性能方面是有一定优势的.二是喷管出口和入口Mach数的比值Ma2/Ma1, 可以看出, Mach数比值越大, 相当于喷管尺度越大, 喷管产生的推力越大.但是, 在工程应用中, 发动机的喷管尺度是固定不变的, 因此这一项的变化不是很大.最后, 从式(10)还可以看出, 为了提高超燃冲压发动机的绝对推力, 主要途径就是提高燃烧产物的静压p1.

图5 喷管喉道参数与出口参数的关系

2 数值算例验证

开展了斜爆轰冲压发动机的数值模拟研究, 用数值模拟结果对上述理论分析结果进行验证, 这里只简要介绍数值模拟, 具体结果详见文献[22-26, 29].图6给出的是发动机模型图, 图7是发动机的三维立体图.发动机采用二维结构, 总长度为2.2 m, 高度为0.55 m, 宽度为0.5 m.进气道为一级进气道, 压缩角度为15°, 长度为1.6 m.燃烧室是等截面的, 长度为410 mm, 高度为76.5 mm, 向下倾斜15°.尾喷管长度是400 mm, 单边膨胀15°, 面积比为2.4.在进气道前缘处用3个小支板实现氢气燃料的喷射和混合过程, 在燃烧室入口处留有边界层抽吸缝.

图6 斜爆轰发动机示意图(单位: mm)

图7 斜爆轰发动机的三维设计图

采用CFD++软件进行数值模拟研究, 控制方程采用带化学反应的Reynolds平均的N-S方程(RANS), 湍流模型采用SSTk-ω全湍流模型[27], 化学反应模型采用9组分19步反应的基元反应模型[28].数值方法、湍流模型以及化学反应模型的合理性在之前的论文中已有详细的说明[22,24,29], 本文不再赘述.为了节省运算时间, 将整个计算域划分为两部分.对于进气道中的燃料喷射和混合问题, 采用三维数值模拟.由于燃烧室是二维矩形燃烧室, 因此, 对于燃烧室内的燃烧过程, 采用二维数值模拟, 其网格如图8所示, 网格总数约3×106, 入口参数从三维数值模拟结果中提取.Ma=7.0，总温是3 372 K，总压是2.26 MPa.氢气喷口Mach数为声速，总温300 K，总压是2.4 MPa，平均当量比为1.2.壁面采用无滑移等温壁面边界条件, 壁面温度是300 K.

图8 燃烧室计算域和网格图

如图6所示, 自由来流经过前体进气道压缩后, 形成一道斜激波(oblique shock wave, OSW), 斜激波在燃烧室上壁面反射形成斜爆轰波(oblique detonation wave, ODW), 斜爆轰波正好打在燃烧室下壁面唇口处.进气道形成的湍流边界层在燃烧室入口前面被抽吸缝抽吸, 避免对燃烧流场的影响.高温高压的爆轰产物经过等截面燃烧室后, 进入尾喷管膨胀做功并产生推力.燃烧室的上壁面是可以前后移动的.如果斜爆轰波不是打在下壁面唇口处, 而是打在下壁面内部, 激波/边界层相互作用就会产生分离泡.分离泡前面产生斜激波, 激波/激波相互作用会产生过驱动的正爆轰波.在数值模拟中, 通过调节燃烧室上壁面的位置, 可以产生斜爆轰波(ODW)和过驱动正爆轰波(normal detonation wave, NDW)两种典型的流场结构.

图9给出了斜爆轰模态(ODW)的燃烧室和喷管流场压力云图和沿一条流线上的Mach数和压力分布曲线(云图方向旋转为水平方向).喷管出口和入口的面积比为2.4.从图上可以看出, 一道斜爆轰波打在下壁面唇口处.从流线图可以看出, 燃烧室入口处Ma=3.5, 经过斜爆轰后, Mach数降低为Ma=1.5, 但是仍然保持超声速流动.受复杂波系的影响, 燃烧室内的Mach数有微小波动.高温高压爆轰产物经过喷管膨胀, Mach数升高, 喷管入口处Ma=1.5(x=1.6 m处), 喷管出口处约为Ma=2.3.

(a)Parameter distribution along the streamline

图10给出了正爆轰模态(NDW)的流场压力云图和沿流线的Mach数和压力分布.在正爆轰模态中, 燃烧室的上壁面后移, 斜爆轰打在燃烧室下壁面上, 引起边界层分离, 形成分离泡.分离泡前面产生斜激波, 激波/激波相互作用形成过驱动正爆轰波.从图中流线可以看出, 流线与正爆轰波面是正交的.经过正爆轰波后, 燃烧室内流动变为亚声速流动, 经过非定常流动后, 又变为超声速流动.喷管入口处约为Ma=1.2(x=1.6 m处), 喷管出口处约为Ma=2.1.

(a)Parameter distribution along the streamline

图11，12分别给出了斜爆轰模态和正爆轰模态的喷管入口截面和出口截面参数分布, 从图中可以看出, 截面上的参数分布是非常不均匀的.为了用理论方法计算喷管产生的无量纲推力, 须给出喷管入口和出口截面上的静压和Mach数.因此, 采用平均方法计算得到各个参数的平均值, 利用平均值再计算得到喷管沿流动方向产生的推力.同时, 将数值模拟结果进行积分得到推力.将二者进行比较, 比较结果见表1.注意: 表1中的推力是沿着喷管上壁面方向的, 不是水平方向.推力计算时, 喷管宽度假定为1 m.

(a)Inlet

表1 理论推力与数值模拟壁面压力积分推力比较

从表1可以看出, 尽管燃烧室和喷管中的流场非常不均匀, 但是理论预测结果与数值积分后的推力吻合得较好, 斜爆轰模态的偏差为-1.78%, 正爆轰模态的偏差为3.28%.在该算例中, 斜爆轰模态的推力约为正爆轰模态推力的2.12倍, 主要是因为前者喷管入口压力约为后者压力的2.28倍.这是因为, 在正爆轰模态中, 由于分离泡的影响, 燃烧室入口有效流道面积减少, 燃烧气体经过分离泡膨胀后在喷管入口处静压降低.图13给出了正爆轰模态的燃烧室流场数值纹影局部放大图, 可以看出, 对于正爆轰模态, 受分离泡影响, 燃烧室流道面积变小, 形成了气动喉道, 而斜爆轰则不受影响.喷管入口和出口Mach数的不同也会产生一定的影响, 但是影响很小.上述结果证明了理论分析的可靠性.