发射箱后盖构型对箱内初始冲击波形成的影响

2023-05-31梁晓扬陆辰昱乐贵高

兵工学报 2023年5期

梁晓扬, 陆辰昱, 乐贵高

(南京理工大学机械工程学院, 江苏南京 210094)

0 引言

箱式发射技术可有效改善导弹的运输和贮运环境,大幅缩短战时装填准备时间,在目前的各型陆基和舰载武器装备中得到了广泛应用[1-3]。发射箱盖作为该技术的重要组成部分,其开启方式直接关系到作战时的发射响应速度,因此一直是发射领域的研究热点。早期,发射箱盖大多采用机械结构或爆破等形式打开,但是这些开盖方式均有不同程度的缺点,如机械盖开盖机构复杂、反应时间长;爆破盖维护成本高、火工品爆炸易对弹头处的精密元件造成损伤等。为满足新的发射要求,简单的开盖方式和便捷的易碎盖技术成为当前的重点研究方向。对于使用该技术的发射箱,主要有3种开盖方式[4],即弹头直接顶破、燃气胀破式(常见于同心筒发射装置)和燃气压力波冲破。其中,利用发射箱后盖反射压力波(由燃气排出时产生)冲破前盖,可有效避免弹头撞击前盖时发生损伤,或者燃气流高温烧蚀弹体和发射装置。多种贮运发射箱装备采用该开盖方式,如美国爱国者防空导弹发射系统。因此,形成提高箱内冲击波强度的设计方法,对发射装置的性能优化具有重要的指导意义。

近年来,国外研究人员采用数值方法和实验方法对开放式激波管[5-7]在短时间内出现的高瞬态流动现象和超声速射流撞击物面形成的流动结构[8-11]进行了深入研究。这些研究结果对于理解瞬态超声速射流的动态演变过程和冲击流场中出现的物理现象有很大的帮助,但其仍不能深刻认识发射箱内冲击波的产生过程。而目前,在国内有很多学者对燃气冲击波开盖技术开展了大量相关研究工作。苗佩云等[12]讨论了燃气开盖技术的必要性,并通过数值方法分析对比了4种燃气开盖方案的可行性,研究结果表明,4种方案均可通过后盖反射的燃气冲击波将发射筒前盖开启。在随后的研究中,姜毅等[13]、傅德彬等[14]采用数值模拟和试验测试方法对发射筒燃气开盖技术进行研究,得出了类似的结论。党海燕等[15]探讨了发射箱前后盖压力匹配的问题。郭锦炎等[16]分析了发动机建压速率对传递到前盖压力波能量的影响。张宝振等[17]采用CEL方法对压力波开盖问题进行研究,得到了前易碎盖的受载特性。以上研究成果在一定程度上为燃气冲击波开盖技术的研究提拱了参考依据,但很少有工作直接分析比较不同后盖构型对箱内反射压力波强度的影响。

本文采用数值模拟方法,结合网格动态自适应技术,对国内某型号发射箱发射过程进行了数值计算,详细地分析了箱内反射压力波的形成机制和传播特性以及其出现的物理现象。对比了后盖构型对箱内反射冲击波的影响,以期获得对于提高冲击波反射强度有意义的设计参数。

1 物理模型

1.1 几何模型

针对冲击波开盖技术发射箱内复杂流动特性的研究,为确保问题研究的真实性,本文采用的发射装置几何模型来自国内某出口型号装备,具体结构布局如图1所示。其中,图1(e)为平板型和穹顶型两种后盖。二者均由环氧泡沫材料制成,并在表面设置弱化槽以实现正反两个方向上不同的承压要求。导弹发动机结构如图1(c)所示,喷管喉部直径Dt=120 mm,出口直径De=350 mm,收敛半角为30°,扩张半角为22°。发动机堵盖材料为轻质泡沫(材料硬度较低,破碎后不会对后盖产生影响),通过胶粘的方式固定在喉部位置,并依据技术条件,满足破膜压强为2 MPa的打开要求。

图1 某发射箱模型及监测面分布示意图Fig.1 Schematic diagram of structure and monitoring surface distribution of a launch canister model

为定量分析箱内冲击波传播过程中强度的变化特性和作用于前盖上的开盖能量,在发射箱内,沿弹体轴线方向布置了多道监测面,其名称和位置分布如图1(a)所示。需要指出的是,两种结构的监测面位置均一致,其具体相对位置(距弹体尾端面)在表1中给出。

表1 发射箱内监测面位置分布

1.2 网格划分及边界条件

发射箱内形成的流动区域具有非常复杂的几何特征,很难对其直接生成高质量的网格,为此本文在几何模型处理方面做了相关简化,删除了对流场影响不大的细节部分,如减重孔、工艺圆角等,但保留了尾翼、导轨等较大的结构特征,以减小几何边界对冲击波传播特性的影响。图2(a)给出了整个计算区域的网格模型。在对计算空间进行离散时,考虑到燃气射流结构属于高度可压缩流动,非结构网格不仅会对计算精度带来一定的耗散,也会影响2阶计算精度收敛的稳定性。因此,本文采用分块网格划分方法,把整个计算空间分成若干子区域,然后对每个子块均采用高质量的结构化贴体网格进行离散,并对流动物理量变化剧烈区域,如边界层区域、喷管壁面附近、射流核心区域,生成较密的网格来更好地捕捉流动细节,具体物面网格分布细节如图2(b)和图2(c)所示。同时,为满足湍流模型壁面函数对物面附近网格的要求,划分10层网格模拟边界层,其中第1层网格距壁面高度约3 mm,相邻网格之间以等比例方式增长,以确保壁面附近的典型y+值小于30。

图2 计算区域网格和边界条件Fig.2 Computational domain and boundary conditions

边界条件的定义如图2中所示,喷管入口为燃气进口,采用滞止总压入口条件,总温为3 400 K,总压值采用样条曲线插值的方法实现,具体燃烧室总压-时间历程曲线如图3所示。数值计算中,发动机燃气和空气均视为热力学性质不同的两种理想气体,且相互之间不发生化学反应,具体燃气组分热力学参数见表2。计算区域外边界定义为压力出口,压强和温度值与周围环境条件相同。模型中发射箱、导弹等壁面边界均定义为绝热无滑移固壁条件,并采用非平衡壁面函数求解近壁面湍流。对于箱体后盖的处理,其开启前,后盖属性定义为固体介质,并记录后盖内表面压强变化情况;当满足开启后盖的设计条件时,即内表面平均压强为80 kPa,则将后盖由固体属性定义为流体介质,即气流可自由通过。

图3 发动机燃烧室总压-时间变化曲线Fig.3 Curve of total pressure varying with time of rocket engine combustion chamber

表2 发动机热力学参数

1.3 网格无关性验证

在数值计算过程中,网格的质量、疏密程度、单元面法线与流动方向是否对齐等因素,关系到控制方程中物理量的空间离散形式和收敛状态,这是直接影响方程组迭代误差的一个较为重要的方面。因此在数值研究过程中,有必要确定模型网格满足高质量和高分辨率的要求。为进行网格收敛性研究,本节采用上述网格划分方法,选取895万(Grid 1)、1 216万(Grid 2)、1 453万(Grid 3)三套网格来评估网格数量对计算结果的影响。图4给出了上述三套网格得到的发射箱后盖平均压强随时间变化曲线,由图4中的曲线变化趋势可以看到,三套网格均能较好地反映后盖受到冲击波作用后的响应历程。所不同的是当网格数量较密(Grid2和Grid3)时,曲线基本重合在一起,数值结果之间差异小于5%,可见模型网格数量在1 200万时,基本可以满足计算误差要求。而网格较为稀疏时(Grid 1),曲线与另外两条出现一定误差,其数值百分误差大于20%。综合考虑硬件资源和时间成本,最终计算模型网格数量控制在1 216万比较适宜。

图4 网格无关性验证Fig.4 Grid independence verification

2 数值方法

2.1 控制方程

基于火箭燃气超声速流动满足连续介质和理想气体假设,在不考虑化学反应和热传导情况下,建立三维多组分可压缩Navier-Stokes(简称N-S)方程,包含质量、动量和能量的控制方程[18-19]形式如下:

多组分输运方程为

(1)

式中:ρ为混合气体密度;Yl为组分l的质量分数;v为速度矢量;Rl为组分l在化学反应后的净生成率;Jl为组分扩散通量;Sl为自定义源项引起的生成率。

Jl的展开形式可写为

(2)

式中:Dl,m为混合物介质m中组分l的质量耗散系数;DT,l为组分l在温度T的热扩散系数;T为环境温度。

在笛卡尔坐标系Oxyz中,针对单一组分l,其守恒形式的三维可压缩N-S方程可表示为

(3)

(4)

(5)

(6)

式中:U为流动守恒变量;F、G、H为无黏通量;Fv、Gv、Hv为黏性通量;K为热传导系数;t为时间;p、e、τ、μe分别为压力、比动能、黏性应力张量、有效黏性系数;u、v、w分别为速度在x轴、y轴、z轴方向上的分量。

对于湍流模型,本文采用RNGk-ε(k和ε分别为湍动能和耗散率)两方程模型[20]。Lijo等[21]研究指出,与其他湍流模型相比,RNGk-ε模型在预测气体射流空间演变的时变过程有更佳的精度,因此可以更精确地求解射流发展过程中的流动细节。其湍动能方程(k方程)为

(7)

式中:μ为混合物黏性;μt为涡粘系数;Gk为平均速度梯度引起的湍动能的产生项;常系数σk=1.39。

湍动耗散率方程为

(8)

式中:常系数σε=1.39;C1ε=1.42;C2ε=1.68;Rε为耗散率方程ε的附加项,

(9)

Cμ=0.084 5,η=Sk/ε,S为平均特征应变率,η0=4.38,β=0.012。

基于有限体积法对上述主控方程离散和非定常流场求解时,空间格式上采用具有3阶精度的MUSCL中心差分格式离散,对流无黏通量求解采用2阶精度的AUSM+格式,以克服高马赫数流场中强间断位置出现的数值不稳定现象[22];时间推进采用LU-SGS(Lower-Upper Symmetric)隐式方法,且每个时间步的物理时间固定为1×10-5,每个时间步内有20个迭代子步,以确保全局流场物理量残差可以下降3个数量级。离散方程组的求解采用Gauss-Seidel隐式方法迭代。

2.2 网格动态自适应

在数值分析过程中,当开始创建计算域网格时,流场特征是未知的,因此生成的网格很难保证计算区域的不同位置都具备合理的网格分布,尤其是涉及到随时间演变的物理问题。由于计算效率和硬件资源的限制,通常不可能对流动区域使用均一的非常密的网格。为更好地解决上述问题,又不大幅增加计算量,在工程应用领域中,自适应网格技术是一种很好且经济的选择。

根据网格加密和粗化策略的不同,网格自适应方法的实现途径有多种,本文求解器采用当前常用的h-type自适应策略。该方法原理简单,通过对局部单元剖分或合并的方式来改变网格的疏密,可有效避免其他方法在求解三维复杂流场时计算稳定性较差的问题。其实现过程和数据结构如图5所示,为完整起见,更具体的细节可见参考文献[23-24]。

图5 二维网格自适应过程及数据结构示意图Fig.5 2D mesh adaptation process and data structure

为有效进行网格自适应,需要确定合适的自适应判据,以标识哪些网格单元需要加密或粗化。根据流场特点,自适应判据可以有多种形式,过去许多研究人员提出并讨论了许多判据指标来判断网格需要变化的区域。在本文中,针对燃气流场包含复杂波系结构的特点,基于Lovely等[25]的研究,综合考虑马赫数和压力梯度的影响,定义自适应判据如下:

(10)

(11)

在具体操作过程中,应用以下条件判断网格单元是否需要加密或粗化:

3 模型验证

为验证本文中数值方法的有效性和网格自适应方法对流场不连续位置的捕捉能力,以Lamont和Hunt等进行的超声速射流撞击不同角度的斜板实验作为验证算例;该实验给数值计算工作提供了高质量的纹影照片,因此是数值求解器的理想测试案例。试验中产生射流的喷管喉部半径r=10.7 mm,出口直径Dn=15 mm,收敛段型面由5.12r喉部半径的圆形与2r喉部半径的圆形相切而成,扩张段为15°半角锥形。喷管在PR为1.2和2的操作条件下运行,其中PR定义为喷管出口处压强与环境压强的比值。具体的实验描述和测试细节请见参考文献[26-27]。图6给出了本文计算得到的数值纹影结果与实验纹影图的对比情况。从图6中可以看到,模拟结果清晰地捕捉到了冲击流场的各种特征,并且激波结构与实验结果一致性较好,证明了本文算法的可靠性。此外,也可以发现网格自适应方法能够在保证数值稳定性的前提下,使接触不连续位置的网格分辨率得到明显提高。

图6 试验结果与本文计算结果流场对比图Fig.6 Flow field comparison between test results and calculated results

图7给出了标准k-ε湍流模型和RNGk-ε湍流模型得到的平板上的压力分布曲线,主要目的是为了评估本文湍流模型对复杂冲击流场的模拟能力。曲线横坐标参量以喷管出口半径进行无量纲化处理(s为喷管轴线与平板交点位置到测量点的距离,Rn为喷管出口半径),纵坐标采用喷管入口总压进行无量纲化处理(pc为喷管入口总压,p为相应位置的压力值)。从压强分布曲线中可以观察到,两种湍流模型预测到的结果与实验测试结果在趋势上保持一致,其压强峰值点也在一定程度上正确地反映了滞止点的位置,但是RNGk-ε湍流模型得到的物理量结果更加接近于实验数据,总体数值差异小于5%。因此,采用RNGk-ε湍流模型对欠膨胀射流进行数值求解比较合理。

图7 试验结果与数值结果平板上压强分布的比较Fig.7 Comparison of test results with numerical results of pressure distribution on the plate

4 计算结果及分析

4.1 发射箱内流动特性

非定常燃气流场动态演化过程中,流场特征(激波、接触间断、冲击波等)位置会随时间的推进在空间上发生迁移。因此,需要进一步确认本文自适应方法能够可靠地对相关网格进行动态调整。在本节中,首先依据前文的自适应方法,在初始网格上使用细化级别为3,最小尺寸为2 mm的细化参数,获得如图8所示的射流发展过程动态自适应网格图。由图8并结合马赫数场云图可以观察到:t=1 330 μs 时刻,当前的自适应方法对流场发展中特别重要的区域,如冲击波阵面、喷管附近初始涡环位置等,进行了有效的局部加密,取得了令人满意的效果;t=1 600 μs时刻,随着流动的发展,当前的自适应方法可以根据需要,将上一时刻网格进行了有效的动态调整,使其更符合物理问题对网格分辨率的要求。

图8 射流发展过程中网格动态自适应Fig.8 Dynamic mesh adaptation during jet development

图9详细地给出了平板型后盖情况下,起始冲击波与箱体后盖相互作用的早期阶段,以及向前盖方向传播冲击波的产生机制。此处的起始冲击波(见图9(a))最初是由发动机堵盖破膜瞬间产生的,其传播速度会领先于未膨胀的燃气射流在喷管内向下游传播。图9中t=0 μs时刻表示为发动机堵盖打开时间。由图9可见:在t=1 200 μs时刻,起始冲击波从左向右传播到箱盖端面,并在喷管边缘附近出现由于剪切层和激波衍射而形成的主涡环;t取值为1 400～1 600 μs时刻,起始冲击波与箱盖发生碰撞,同时在壁面上冲击点附近形成一个滞止区域,且随着时间的推移,此滞止区域会逐渐扩大。由于此区域附近的物理参量变化剧烈,流动结构在这里会发生极大的改变,导致流线不断沿径向向外偏转,会使跨过此处的气流逐渐转折成与壁面平行;t取值为1 800～2 900 μs时刻,沿径向传播的气流与箱体侧壁发生相互作用,并在箱体侧面拐角处附近发生另一个改变气流方向的梯度区域,引起流动偏转,从而出现向前盖方向传播的压力波。

图9 不同时刻箱内冲击波发展过程Fig.9 Developing process of the shock wave in the launch canister at different times

图10显示了平板型后盖情况下,冲击波在箱内传播及与前盖相互作用的整个过程。由图10可见:由于发射箱与弹体之间形成的主通道具有一定的导向性,在t分别为6.23 ms、8.73 ms和9.73 ms时刻,最初经箱体后盖与发动机起始冲击波相互作用产生的压力波在发射箱前半部分以稳定的平面冲击波形状向前发展,并且其波阵面强度随传播距离的变化出现一定程度的衰减。此后,随着箱内通道截面积的增加,超压峰值明显下降;在t=17.73 ms时刻,运动的冲击波到达前盖并与之发生碰撞,在前盖的边缘位置处形成规则反射;由于发射箱前盖具有一个凹壁轮廓,在t取值为18.23～18.73 ms时刻可以发现冲击波在前盖中心点位置汇聚,并出现局部冲击强度的迅速增加,显然这对前易碎盖的顺利吹破是有利的。

图10 冲击波在箱内传播及与前盖相互作用过程Fig.10 Shock wave propagation in the launch canister and interaction with the front cover

4.2 箱盖几何构型对反射冲击波的影响

为对比不同后盖构型对箱内冲击波形成的影响,图11定性地比较了平板型箱盖结构和穹顶型箱盖结构下不同时刻箱内的流动特性。在这两种情况下,两种箱盖端面与喷管出口平面之间有相同的距离,并且后盖打开压强相同,以便可以更好地对这两种结构进行评估。由图11可以明显地观察到两种流动之间的差异:两种情况下起始冲击波在相同时刻抵达箱盖并与之发生相互作用(t=1 330 μs);随后气流发生偏折,并沿径向方向发展,且二者呈现出相似的流动特征(t=1 830 μs)。这两种情况的主要区别在于,平板型箱盖与起始冲击波相互作用后,气流方向会折转成与箱体侧壁垂直;因此,当冲击波移动至箱体侧壁时将与之发生垂直碰撞;然后在壁面前方出现一道强烈的反射激波(t=1 930 μs)。紧接着,反射激波与迎面而来的气流产生相互作用,并在箱体拐角位置形成了一个高度不稳定的湍流区域。由于反射激波的强烈扰动耗散了流动中所包含的能量,最终会显著地影响箱内通道中的压力波强度,见图中t取值为2 030～2 630 μs时刻。对于穹顶型箱盖,由于其结构边缘有一定程度的抬升,偏转后的气流在沿箱盖表面发展过程中,受结构表面的限制,使得其流动方向逐渐发生改变,最终影响其撞击箱体侧壁时的角度。意味着冲击波与箱体侧壁相互作用后将明显地减缓箱体拐角位置处湍流区域的形成,从而可以有效地提高箱内初始压力波强度。

图11 不同后盖构型箱内冲击波形成对比(左为平板型,右为穹顶型)Fig.11 Comparison of shock wave formation in the caniser with different cover configurations (Left is plate type, right is dome type)

为了对冲击波在扫过箱内通道抵达前盖过程中的传播特性进行定量评估,图12分别显示了两种后盖构型下箱内各监测面和前盖上平均压强随时间变化曲线。由图12可以清晰地观察到:发射箱内冲击波在通过各监测面后,监测面上的超压峰值迅速上升,然后随着冲击波向下游的传播,各监测面上的超压峰值会出现一定程度的下降;根据曲线形态和变化趋势可以推测两种箱盖构型在与发动机起始冲击波相互作用后,均能在箱内形成稳定的向前盖方向传播的冲击波,从侧面表明了利用燃气压力波将前盖冲破在原理上是可行性。通过比较两种箱盖构型下各监测面上的压力峰值可以发现,穹顶型箱盖构型中,箱内压力波在弹体尾端面位置(S-0)上的超压峰值从50.4 kPa提升到了74.8 kPa,而最终前盖上的超压峰值从29.3 kPa提升到了43.6 kPa,综合分析,两种情况之间的差异达到了48%左右,这与前面的分析具有很好的一致性,表明了穹顶型箱盖可以有效地提高箱内压力波初始峰值的强度,促进在前易碎盖处可以形成更强的压力峰值。另外,根据发射箱前盖上的压力-时程曲线可知,当冲击波作用于前盖时,其压力峰值相较于前方的监测面S-5有所上升,这种现象主要是由于组成冲击波的一道道压缩波在与前盖相遇后得到进一步叠加,导致压力显著增大。