一般而言, 吸气式脉冲爆震发动机(air-breathing pulse detonation engine, APDE)由进气道、燃烧室与尾喷管组成。典型的脉冲爆震发动机循环包括爆震波起始、爆震波传播、燃烧产物的排出以及混气的填充过程^[1], 其非稳态的工作过程导致进排气系统的设计与现有的稳态设计存在很大差别。当工作频率很高时, 可以近似认为APDE的工作过程是连续的, 基于此假设可以借鉴现有的稳态设计方法来进行APDE的进气系统设计。即便如此, 亦须考虑进气道与爆震燃烧室相互作用的问题。例如, 在无阀工作模式下, 反传高压气体可能从进气道进口溢出, 产生负推力, 在超声速来流工况下, 甚至会直接破坏进气道进口激波系, 导致进气道无法起动; 同时, 进气道内流场结构的变化会影响燃烧室的性能, 进气道与燃烧室间形成的上述反馈机制会导致两者之间产生严重的耦合作用^[2-3]。目前, 已有不少学者采用数值模拟与实验的方法对APDE超声速进气道进行了研究。Oh等^[4]对轴对称超声速进气道扩张段内激波与声波的相互作用进行了详细的非定常数值模拟, 将稳态计算结果的出口背压作为基准量, 在此基础上施加了正弦压力脉动分量, 分析了由声波振荡引起的激波/边界层相互作用。Mullagiri等^[5-6]在进气道出口安装了4个具有90°相位差的活塞, 模拟了多管APDE进口工作状态, 发现压力振荡的幅度随着阻塞比的增大而增大, 随着激励频率的增加而减小。上述研究对反传压缩波进行了大幅简化, 仅考虑了压力对进气道的影响, 且在压力波形相似性上也有很大提升空间。后来, Segal等^[7]通过在扩压器出口横截面拐角处向主流周期性喷注二次流模拟脉冲爆震燃烧室反压, 研究了振荡背压对脉冲爆震发动机超声速进气道的影响。结果表明, 楔板上的分离泡会破坏流动, 形成复杂的压缩波与膨胀波系。关祥东等^[8]在瞬态数值模拟过程中将进气道出口边界条件设定为压力出口/固体壁面, 模拟机械阀门开闭形式的推力壁对进气道内部流场的影响。虽然该假设略显粗糙, 但计算结果仍有助于获得流场变化规律的初步认识。王丁喜等^[9]通过数值模拟在某轴对称超声速进气道出口处添加了3种不同的扰动, 得到了进气道内结尾正激波随出口条件变化的运动规律。以上研究对反传压缩波的模拟更加逼真, 但未考虑脉动燃气温度等因素的影响, 与实际情形存在一定差异。

基于上述问题, 出现了一批将进气道与燃烧室相接以保证反压真实性的研究。Wang等^[10]通过实验研究了无阀式APDE反压传播的规律, 结果表明, 进气道中的扰动强度随填充度的变化呈非线性波动, 未进一步深入分析进气道内流动机理。Chen等^[11]设计了进气道中心锥与气门结构, 通过实验研究了其对APDE推力与压力/燃气反传的作用。结果表明, 中心锥对进气稳定性有积极影响, 气动阀会增加进气道阻力, 但对抑制反压有明显效果。上述研究虽获得了真实的反压, 但对进气道结构考虑不周, 均采用皮托式进气道。除此之外, 中心锥对流道的堵塞作用会使得燃烧室进口流场分布不均匀, 总压畸变程度增大。

近年来, 涌现出若干真实进气道与燃烧室匹配工作的研究。温玉芬等^[12]采用非定常仿真方法, 在混压式进气道与燃烧室、尾喷管耦合模型上组织了单次爆震, 得到了超声速进气道在反传压缩波作用下的内部流场及工作特性, 对比了不同扩张段进气道与燃烧室匹配工作的差异。然而, 其研究的重心是进气道工作过程, 并未分析推力与比冲等性能参数。Zangiev等^[13]设计了带有往复式机械阀的燃烧室与超声速进气道模型, 模拟了马赫数为3.0与5.0条件下的多循环爆震。结果表明, 在28 km高空、马赫数为5.0条件下, APDE的比冲(1 600 s)和比推力(1.03 kN/(kg·s))均大于冲压发动机的比冲(900~1 200 s)和比推力(0.5 kN/(kg·s)), 而耗油率与冲压发动机相当, 均为0.21 kg/(N·h)。该研究结果证明了APDE采用冲压进气方式工作的潜在优势, 但未对进气道工作过程进行分析。Qiu等^[14]一维数值模拟的结果表明, 与孔板型气动阀相比, 超声速内压式进气道使APDE有可能获得与理想Brayton循环相当的推力增益, 能否与Brayton循环效率相当取决于结构优化是否能够最大程度地减小总压损失, 指出了APDE超声速进气道优化的必要性。

综上所述, 国内外在APDE超声速进气道的研究方面积累了一些经验, 但仍存在以下问题: ①由于对燃烧室复杂的进口条件简化程度过大, 得到的结果与真实进气道的工作情况差异较大; ②目前大多数研究均采用无阀方案, 侧重分析反压对进气道的影响, 很少有学者考虑研究机械阀对进气道与燃烧室的影响; ③兼顾进气道与燃烧室的研究或侧重进气道的工作, 或侧重燃烧室的性能, 未对二者相互作用的影响做出详尽分析。为此, 有必要研究二者的匹配工作特性。本研究建立了以马赫数2.5为设计点的三波系混压式进气道与脉冲爆震燃烧室模型, 并在此基础上设计了3种工作方案。采用带化学反应的非定常仿真方法, 模拟通过直接起爆方式形成的自持传播爆震波, 并分析不同方案中的反传效应对进气道内流场的非定常干扰过程及对燃烧室内流动特性与推进特性的影响。

1 物理模型

本文的仿真对象为由进气道与燃烧室构成的内流域及部分远场流域。进气道外压段根据Oswatitch最佳波系理论计算获得。针对该模型设计了3种工作方案, 如图 1所示, 燃料选择氢气, 氧化剂为来流空气。待进气道建立起稳定的波系后, 3种方案均按图 2所示的循环工作, 其中, 填充混气指自加入氢气质量源项起, 至混气填充至燃烧室出口止; 点火与爆震波传播阶段指自设置高温高压区起, 至混气燃烧完止; 泄压阶段指自燃气完全排出燃烧室起, 至反传压力对上游的影响消失时止。方案A与方案B均为无阀方案, 燃烧室与进气道内流域连通, 方案B较方案A的燃烧室在出口壁面上增加了5个泄压小孔, 旨在加快泄压过程, 提高工作频率; 泄压孔在进气、起爆、爆震波传播与排气阶段关闭, 仅在泄压阶段打开。方案C为机械阀方案, 在进排气与泄压阶段, 燃烧室头部阀门打开, 与进气道内流域连通, 进气道出口泄压孔关闭; 点火时燃烧室头部阀门关闭, 作为燃烧室头部推力壁, 同时为了减小进气道冲压阻力, 将其出口泄压孔打开。3种方案中, H₂质量源项与高温高压区的起始位置均相同。

三波系进气道设计马赫数为2.5, 设计点飞行高度为11 km, 其几何参数见表 1, 一维等熵计算所得的气动参数见表 2。

3种方案的燃烧室长高比L_c/H_c均为9.7, 方案B与方案C中燃烧室出口的泄压小孔均为5个, 中心位于距燃烧室出口25%L_c处, 单个泄压孔宽4 mm。方案C中, 阀门关闭时其顶端距燃烧室出口398 mm, 与燃烧室中轴线夹角30°, 进气道泄压孔宽38 mm, 参数如图 3所示。

2 数值模拟方法 2.1 网格无关性验证

本研究在上述模型的基础上, 将进气道进口分别向上游和外侧延长2H和H, 燃烧室出口分别向下游和外侧延长10H和H之后的流体域作为计算域。在整个计算域中填充结构化网格。网格单元总数约12万。对流体域中所有的近壁面网格沿垂直于流动的方向进行加密, 内流域大部分近壁面网格的y⁺保持在35左右, 外流域大部分近壁面网格的y⁺为40左右, 方案A的计算网格如图 4所示。

为排除网格对仿真结果的影响, 分别采用网格单元总数为15万(细网格)、12万(中网格)与8万(粗网格)的3套网格在如2.2节中所述的计算设置与边界条件下进行了仿真对比, 选取对网格疏密程度较敏感的壁面静压进行分析。图 5给出了通过3种网格计算的进气道建立起激波系后压缩面的静压分布, 可以看出细网格与中网格的计算结果在量级与分布规律上一致性较高, 粗网格的计算结果较其余两者偏差较大, 为降低仿真过程的总体计算量, 同时不致对结果带来较大误差, 选用中网格进行计算。

2.2 控制方程、湍流模型及边界条件

本研究采用商业CFD软件进行仿真, 对流项采用隐式格式进行计算, N-S方程的求解采用有限体积法, 其中无黏对流通量采用Roe平均方式进行MUSCL(monotonic upwind scheme for conservation laws)插值, 黏性通量采用2阶中心差分方式进行离散, 湍流模型选用k-ε模型。

仿真时, 设置自由来流马赫数为2.5, 静参数按11 km高度标准大气(p₀=22 670 Pa)条件给定。进气道建立起波系后, 监测燃烧室进口空气流量, 按化学恰当比将氢气质量源项添加在燃烧室进口。当氢气基本充满燃烧室时, 在燃烧室头部添加压力为2 MPa, 温度为3 000 K, 宽度约为燃烧室长度5%的高温高压区, 模拟直接起爆过程。化学反应机理采用FLUENT中的氢气-空气单步反应机理, 组分输运模型中的化学反应项选择体积反应模型, 湍流-化学相干模型选择有限速率模型, 化学反应求解器选择刚性求解器。根据阿累尼乌斯定律计算得到的化学反应特征时间约为1×10^-6 s, 选取计算时间步长为化学反应特征时间的10%, 即1×10^-7 s。在一个循环周期里, 方案B、C中燃烧室出口泄压小孔处边界条件在燃烧产物基本排出燃烧室, 即排气过程完成时由壁面改为压力出口, 方案C中燃烧室进口阀门在点火时刻由内部面改为壁面, 同时进气道出口侧壁上泄压孔由壁面改为压力出口。在有效推力由正变为负时, 燃烧室进口阀门由壁面改为内部面, 同时进气道出口侧壁上泄压孔由压力出口改为壁面。

3 算例检验

为检验数值模拟方法对进气道内非定常流动计算的准确性, 采用上述方法对文献[15]中的二元混压式超声速进气道进行了仿真, 并将仿真结果与实验结果进行对比。在文献[15]的研究结果中, 改变了进气道出口背压, 获得了进气道从不起动到再起动的工作过程, 对进气道内非稳态流动特性也有较详细研究, 故该文献中的实验结果可以用来检验本文第2节中的研究方法。

图 6为自由来流马赫数2.51, 攻角2°, 进气道处于亚临界状态时, 使用本文仿真方法获得的压缩面静压分布规律与文献[15]中实验测得的进气道压缩面静压分布规律一致。

图中, 纵坐标表示压缩面静压, 通过来流总压对其进行无量纲化(p/p₀); 横坐标表示进气道沿流向的位置, 通过进气道长度对其进行无量纲化(x/L), 将外压面顶点作为横坐标零点。可以看出通过数值仿真方法获得的压缩面静压分布与实验测量值的分布规律相符, 最大相对误差存在于x/L=0.3, 约为13.6%, 15个测点的平均相对误差为4.3%。综上, 可以认为本文仿真方法能够较准确地获得进气道内非定常流动特性。

4 结果与分析 4.1 爆震波传播特性

图 7是爆震波与反传压缩波无量纲轴向位置及两者波锋面无量纲压强(通过来流静压进行无量纲化)随时间变化的曲线。可看出, 自点火开始, 反传压缩波与爆震波锋面均具有稳定的移动速度。爆震波相对于地面参考系的绝对速度为2 602.61 m/s, 减去当地气流速度857.5 m/s后, 爆震波相对于未燃混气的移动速度为1 745.11 m/s。数值模拟计算得到的未燃混气压强为56 514.3 Pa, 对应条件下使用CEA估算的波速为1 686.4 m/s, 相对误差为3.5%。根据文献[16], 数值模拟所得波速较CEA计算的波速高出5%以内是可以接受的。采用2.2节中的仿真方法无法获得爆震波胞格结构, 但对于研究所关注的爆震波物性参数, 可以认为本文模拟方法获得的爆震波特性可靠。从爆震波形成至传出燃烧室出口, 爆震波锋面与反应区都未发生解耦, 波后已燃气体相对于波前未燃气体的静压比值p/p₁约为15.87, 静温比T/T₁约为10.81, 如表 3所示, 与CEA计算结果基本一致。综上, 可以认为数值研究可以较准确地模拟爆震波的直接起爆与传播过程。

4.2 进气道内流动特性分析

图 8为点火前进气道喉道截面(x=0.17 m)与燃烧室进口截面(x=0.59 m)处气流沿水平方向的速度。可以看出, 在进气道喉道处, 上下边界层厚度不一, 且主流区流速也有较大的空间不均匀度, 这是进气道内流道的几何不对称性导致的。气流经过0.42 m长的扩张段, 到达燃烧室进口时, 主流区速度基本保持在668.3 m/s, 且上下壁面边界层厚度也基本一致, 故由进气道几何结构不对称性对爆震燃烧带来的影响可以忽略不计。

方案A与方案B中进气道内流动特征相似, 仅反传距离有所不同; 而方案C中机械阀的关闭会导致机械阀前产生斜激波, 与方案A与方案B进气道内流动特征有较大差别。故在本小节主要对方案A与方案C进气道内流动特性进行分析。

在方案A中, 自点火开始, 高压燃气形成的反传压缩波与爆震波便分别向上下游传播。波锋面轴向位置及压强的时间历程如图 9所示, 自起爆后爆震波波峰压强始终保持在1.3 kPa左右, 反传压缩波波峰压强持续减小, 且该过程在爆震波传出燃烧室后仍在持续进行。受燃气高压的影响, 反传压缩波形成初期以正激波形态存在, 在逐渐向上游传播过程中与边界层发生干扰, 以分叉激波形态存在, 后随着进气道出口背压的降低, 形成了“莲藕状”的激波串结构。图 10为反传压缩波到达最大位置处时扩张段内的数值纹影, 可以看出分叉激波的2个分叉点已远离壁面, 相交于中轴线附近, 在其后形成了一系列斜激波串, 诱导附面层发生分离, 使附面层厚度增加, 主流区流道变窄, 流速增大, 增加了进气道出口截面流场的不均匀度。此外, 激波串中各级激波由后至前依次向后移动, 激波强度减小, 间距增大, 整个流场表现出较强的湍流特征。点火初期, 反传压缩波波峰压力为0.52 kPa, 当爆震波传播至燃烧室出口时, 反传压缩波波峰压力降低至0.39 kPa, 持续向上游传播的同时压强减小, 波锋面压强的时间历程见图 7。

在方案C中, 当阀门关闭时, 进气道出口空气冲击在阀门上, 在阀门前形成一道21°的斜激波, 波后马赫数降低至1.03且流动方向改变, 从进气道侧壁上的泄压孔流出。在斜激波后靠近阀门顶端的位置还形成了一道弓形激波, 气流在弓形激波后被滞止, 形成流动分离范围较小的回流区, 如图 11所示。在整个阀门关闭的过程中, 斜激波未向进气道上游移动, 未破坏进气道内流场结构。

3种方案中, 反传压缩波能达到的最大位置分别为62.6%L、57.4%L与26.2%L, 均位于扩张段内, 进气道喉道及上游的流动未受到反传效应的影响, 见图 12。

4.3 燃烧室内流动特性分析

3种方案中, 燃烧室内高压区的压力波均会在上下壁面间反射振荡。由于燃烧室几何结构具有对称性, 且进口不存在遮挡, 方案A与方案B燃气中的压力振荡仅存在于竖直方向; 而在方案C中, 由于燃烧室头部与水平方向有一夹角, 压缩波在头部壁面被反射, 同时具有水平与垂直于爆震波方向的速度, 在爆震波传出燃烧室前共形成了3级反射压缩波, 如图 13所示。随着时间的推移, 一系列反射压缩波逐渐追赶上爆震波, 在爆震波后方形成弓形高压区, 但自始至终反射压缩波均未与爆震波波面重合, 且爆震波波速与方案A、B中没有反射压缩波作用的波速相当。爆震波与各级反射压缩波轴向位置及与水平方向夹角如图 14所示, 在反射压缩波形成初期, 燃气温度较高, 故波速也较大, 追赶上爆震波后与爆震波保持同速, 同时各级反射压缩波与水平方向夹角也逐渐增大, 最后保持在90°左右。

在阀门打开前, 燃烧室头部边界层附近形成了顺时针回流区, 阀门打开后回流区中心移动至主流区, 且影响范围扩大, 同时将下方燃气卷吸至上游, 进入扩压段内, 如图 15所示。

回流区上方燃气从阀门打开瞬间开始始终具有与主流方向相同的速度, 故可以认为方案C中的燃气未因压力驱动而反传, 下壁面附近的燃气向上游移动是由湍流脉动导致的。然而, 方案C中阀门关闭时, 燃烧室内燃气无法迅速排出, 主流区燃气在点火后3.89×10^-3 s时才被完全吹除。而在方案A、B中, 燃气反传的最大距离为0.03%L, 主流区的燃气在点火后2.85×10^-3 s时就被完全吹除。

4.4 燃烧室推力特性分析

为衡量3种方案性能, 需要对有效推力进行分析, 有效推力F_eff由(1)式计算

(1)

式中: q表示流量; u表示气流轴向速度; A表示截面积; 下标out表示燃烧室出口截面; in表示进气道进口截面; outer表示进气道外壁面, 故表示外流作用于进气道外壁面上的压力在发动机轴线上的投影。

图 16为通过(1)式计算出的方案A、B与方案C 2个循环内的有效推力, 可以看出, 方案B在工作过程中几乎未获得正向的瞬态有效推力, 仅在排气与泄压过程中受到的阻力有所减小, 平均有效推力约-1.446 kN, 方案C的平均有效推力约为1.125 kN, 见表 4。方案B的有效推力较方案C小的原因在于, 方案B中点火后形成局部高压区, 其上下游均为低压区, 高压燃气会同时向上下游膨胀。反传燃气中一部分压力势能通过逆流向上游膨胀被耗散, 将燃气中沿x轴负方向的动量传递给来流空气, 使燃气与来流空气形成的混气沿水平轴正方向的速度降低, 但其并未作用在固体壁面上, 不仅无法产生推力, 还会影响进气道正常工作。由于进气道喉道之后是扩张型通道, 反传燃气中仅很少一部分能量可以作用在壁面沿y方向的投影上, 产生部分推力, 这也是方案B在点火时阻力稍稍减小的原因。

在方案C中, 燃烧室头部阀门打开时, 反传压缩波仍会向上游移动, 但移动距离较方案B短, 反传压缩波能达到的最大位置也较方案B靠后, 且平均移动速度也更低, 如图 12所示。表明反传燃气传递给来流空气的动能与压力势能更少, 正是由于可以作用在固体壁面上, 燃气才会更少地向上游膨胀做负功, 将绝大部分能量用来向燃烧室出口膨胀做正功产生推力。

在超声速工况下, 机械阀不仅能最大程度地利用燃气做功, 还能加快循环速度。如表 4所示, 方案C中单个循环的时间约为12.77 ms, 较方案B的46.25 ms减少了约72.4%。原因在于, 方案B中, 循环频率受限于泄压过程, 而方案C中的机械阀能对高压燃气起到遮挡作用, 可以使更多的压力势能通过在燃烧室出口膨胀的方式转化为动能, 削弱了高压燃气向上游缓慢的膨胀作用, 加快了泄压过程持续的时间。

综上所述, 在超声速来流工况下, 尽管机械阀带来了一定阻力, 但在增推与提高工作频率方面具有一定优势。

4.5 阀门开闭的影响

由图 16b)可以看出, 阀门关闭时推力由-2 kN骤降为-6 kN。因起爆位置距阀门有一定距离, 阀门关闭与点火同时进行, 阀门关闭时反传压缩波无法有效作用在推力壁上, 增加的4 kN阻力是由阀门关闭导致的, 对应的瞬态流场见图 12。从现有结果看, 增加的额外阻力并不大, 但其在时间上的积累很可观。减小该阻力影响的途径有2种: 一是通过合理的阀门型面设计等方法实现; 二是减小该瞬时阻力持续的时间, 即在反传压缩波正好移动到阀门位置时再关闭阀门。在实际应用中, 阀门的精确控制很难实现, 所以应尽量在设计过程中缩短起爆区域到阀门的距离, 缩短阻力持续的时间。

此外, 阀门打开瞬间, 推力由1 kN骤降为0.5 kN。在本方案中, 阀门打开的依据是有效推力由正变负, 实际上在该时刻阀门两侧尚存在压力差, 虽然固体结构不再获得正向推力, 但燃气仍对固壁做功, 起减小阻力的作用。判断阀门开启时间的正确方法应该是, 待阀门后高压燃气的压力降至与阀门前冲压作用产生的压力相当时打开阀门, 充分利用燃气中的压力势能。

5 结论

针对设计点工况马赫数为2.5的三波系二元混压式超声速进气道与脉冲爆震燃烧室进行建模, 并在此基础上设计了3种工作方案, 通过数值仿真方法模拟了其多循环工作过程, 并分析了反传效应对进气道与燃烧室的非定常干扰过程及其对燃烧室内流动特性与推进特性的影响。主要结论如下:

1) 燃烧室出口泄压与燃烧室进口机械阀均能缩短反传压缩波向上游传播的距离, 且能加快工作频率。与基准方案A相比, 燃烧室出口泄压措施会使时均推力损失64 N, 比冲减小0.2 s; 机械阀遮挡措施会使时均推力增大2 507 N, 比冲增大178.5 s, 且频率也会由21 Hz增加至78 Hz。

2) 机械阀会对燃气中的反传压缩波产生反射作用, 增大燃烧室内的湍流度, 在开阀瞬间使燃气被卷吸至上游进气道扩张段内。燃气中的反射压缩波会追赶上爆震波, 但不会对爆震波形态与波速产生影响。

3) 在3种方案中, 推力大小不同是因为向上游膨胀损耗的程度不同, 反传距离不同是因为燃烧室前方存在不同程度的固体阻塞, 固体阻塞能够抑制压力反传, 同时加快循环频率。故为了增大推力、缩短反传距离与获得更高的循环频率, 在点火及爆震波传播过程中, 应在燃烧室进口尽可能多地扩大物理遮挡物阻塞比, 但不应对进气道出口流动造成堵塞, 增加流动阻力。这样一方面爆震波可以更少地向上游膨胀损耗, 另一方面可以削弱反传对进气道的影响。此外, 在不对推力产生可观损耗的前提下, 泄压过程应尽可能地快, 以提高频率。