基于残值函数法的双燃速药柱内弹道分析*

2020-05-13褚佑彪赵天泉董新刚

固体火箭技术 2020年6期

褚佑彪，赵天泉，董新刚，任萍

(中国航天科技集团有限公司四院四十一所，固体火箭发动机燃烧、热结构与内流场国防科技重点实验室,西安 710025)

0 引言

燃面计算用于确定装药在燃烧过程中燃烧表面积随燃烧时间的变化规律，直接影响发动机内弹道性能预示精度，是发动机内弹道设计的基础，在固体火箭发动机的设计中一直占有重要地位[1]。为满足先进固体动力技术的发展要求，燃烧室的药型设计越来越复杂。在单室双推发动机的设计过程中，为保证助推段与续航段的推力比，常采用高低燃速搭配的思路进行装药设计。当高低燃速推进剂在交界面处燃烧时，燃面会出现分离、交汇等复杂的拓扑结构变化，是燃面推移计算的一个难点。

目前，可实现燃面不等速推移的有实体造型法[2]、网格推移法[3-4]、最小距离函数法[5-6]、Level Set方法[7-8]和残值函数法[9]等。其中，实体造型法通过特征造型或驱动尺寸实现燃面推移过程的模拟。但是，对于结构复杂的药型，实体造型法的推移过程十分繁琐，并可能出现奇点，导致计算无法继续；对于多燃速药柱，燃面的分离需要在推移过程中构建新的几何特征，较难实现。网格推移法在初始药型上生成非结构网格，利用当地燃速进行推移，再通过网格重构获得新的燃面，该方法通用性较好，但其处理燃面交汇、分离等拓扑结构变化的稳定性较差。最小距离函数法通过计算药柱内部各点到初始燃面的距离，即最小距离函数(MDF)，进行燃面推移分析，可以实现不等速推移，但是其无法自然捕捉到由于燃速间断的存在而导致的燃面拓扑结构的变化。Level Set方法采用初值形式的偏微分方程将一个纯几何问题转变为用偏微分方程描述的数学问题，对于推进剂内部存在燃速间断的工况，由于燃气与推进剂之间也存在间断，在数值求解微分方程组时，会遇到间断相交的情况，进而导致燃面计算精度及稳定性较差。残值函数法[9]在笛卡尔网格上，采用残值函数记录燃面位置，通过惠更斯原理模拟燃面的推移过程。此算法可以准确捕捉由于燃速间断的存在导致燃面的交汇、分离、消失等复杂拓扑结构变化，且不需要进行网格重构，稳定性好。

本文应用残值函数法模拟了多燃速药柱的燃面推移过程，并与试验结果进行对比分析。在此基础上，开展参数研究，分析了高低燃速推挤剂交界面的位置对药柱燃面演化过程和内弹道的影响。

1 内弹道计算方法

为准确分析双燃速药柱在界面处的耦合燃速过程，本文采用残值函数法[9]进行燃面计算分析。该方法基于燃面推移的一般性原理——惠更斯球面波传播原理，不仅可以处理燃面的等速推移，还可以处理燃面的不等速推移，包括燃面上存在燃速间断的工况。由于残值函数法采用笛卡尔网格离散药柱，对燃面的推移区域可以有效地提前预估，因此可以将三维体循环缩减为具有一定厚度的曲面循环，大大降低了计算量。此外，在燃面计算过程中引入残值函数，不仅有效地控制燃面推移的累计误差，还可准确捕捉燃面推移过程中形状变化，自然处理界面拓扑结构变化，避免燃面的网格重构，提高算法的稳定性。

残值函数法中燃面的推移流程可简要描述为：

第一步：如图1所示，Pi为当前燃面上的网格点。预估燃面上的离散点Pi的影响区域，并向四周推移，在推移的区域内，推进剂转化为燃气，推移距离Δi由式(1)给出。

Δi=dt×ri+δi

(1)

式中 dt为时间步长；ri为网格点Pi处推进剂的燃速；δi为Pi处的残值函数。

第二步：计算并更新残值函数。Pj为药柱内的网格点且处于Pi的影响范围内。Pj处的残值函数δj由式(2)获得。

(2)

其中，m为当前燃面上能够影响到Pj的网格点数目，εij的值根据Pi与Pj之间是否存在燃速间断进行分类计算：

(1)如图1(a)所示，若Pi与Pj之间没有燃速间断，当Pi与Pj之间的距离dij趋于0时，相应的燃速ri和rj趋于相等，则推进剂沿直线燃烧，所以Pj处的εij可以由

εij=Δi-dij

(3)

获得。

(2)如图1(b)所示，若Pi与Pj之间存在燃速间断，类比折射定律可知，推进剂燃烧路径满足：

(a) No discontinuity in burning rate (b) Discontinuity in burning rate

risinθj=rjsinθi

(4)

式中θi和θj分别为两种推进剂燃烧路径与推进剂交界面法向量的夹角。

O点为推进剂燃烧路径与交界面的交点。因此，Pj处的εij可由式(5)获得：

εij=(Δi-diO)rj/ri-dOj

(5)

式中diO为Pi到O点的直线距离；dOj为O点到Pj点的直线距离。

第三步：确定新的燃面，即更新后的燃气区域与药柱区域的交界面。

结合内弹道基本方程，可以计算燃烧室内压强，获得发动机内弹道性能。

2 某双燃速药柱内弹道分析

单室双推发动机的设计过程中，为保证助推级与续航级的推力比，常采用高低燃速搭配的思路进行助推级和续航级装药设计。图2为某发动机药柱结构，绿色表示助推级药柱，蓝色表示续航级药柱。助推级采用低燃速大燃面药型，药柱结构为盲孔和翼柱相结合结构，续航级采用小燃面高燃速端燃药柱结构。推进剂性能参数见表1。

表1 推进剂性能参数

图2 药柱结构图Fig.2 Structure diagram of the grain

图3给出采用残值函数法获得的发动机内弹道数据，与试验结果进行对比可知，残值函数法准确地模拟了双燃速药柱燃烧过程燃面的变化过程。

图3 数值仿真与试验结果对比Fig.3 Comparison of results of numerical simulation and experiment

由试验结果可知，在助推级与续航级转级过程中，发动机压强由7.8 MPa急剧下降至1.5 MPa，然后逐渐增加至局部极大值6.4 MPa，最终逐渐减少至2.6 MPa，构成一个非常明显的压强起伏。图4给出典型时刻药柱结构图，图中给出药柱过轴线的剖面图，红色区域代表燃气填充区域，绿色代表助推级低燃速推进剂，蓝色代表续航级高燃速推进剂。如图4(a)所示，燃面与高低燃速推进剂交界面刚刚接触，此时，低燃速药柱侧面尚有3.5 mm厚推进剂。随着燃面的推移，在t=7.97 s时刻，燃面推移至高燃速药柱内，燃面迅速“膨大”，但是由于低燃速药柱燃面减少量显著，导致燃烧室压强迅速降低。在t=7.97～10.28 s范围内，燃面在高燃速药柱中不断“膨大”，且低燃速药柱的两侧同时燃烧，燃烧室压强逐渐回升。随着燃面的不断推移，药柱燃面逐渐趋于平面，面积减少，低燃速药柱燃尽，导致燃烧室压强趋于平稳。由以上分析可知，此双燃速药柱发动机的内弹道呈现显著起伏现象，主要是由于燃面由低燃速药柱向高燃速药柱过渡时，明显“膨大”所致。

(a) t=6.82 s (b) t=7.97 s

3 双燃速药柱交界面位置影响分析

双燃速药柱成型过程中，首先进行续航级药柱浇注，待续航级药柱预固化后，再浇注助推级药柱。两级药柱的交界面一般通过续航级装药量进行控制，控制精度不高，本节通过参数研究，进行交界面位置对发动机内弹道性能的影响评估，进而指导药型优化设计及药柱成型要求。

图5基于图2所示药柱将高低燃速交界面向两侧分别移动20 mm进行对比分析，推进剂性能参数保持不变，算例A高低燃速交界面位置为x=320 mm，算例B高低燃速交界面位置为x=280 mm。结果表明，算例A和算例B内弹道特性与原药柱存在显著差异。算例A在转级过程中，燃烧室压强陡升至14.7 MPa，超出助推级初始高压段的最大压强。算例B在转级过程中，燃烧室压强大幅降低至0.58 MPa。

图5 不同交界面位置下的压强-时间曲线Fig.5 Pressure-time curves at different interface positions

图6给出算例A在t=6.3 s时的燃面拓扑结构，此时低燃速药柱推移肉厚为70 mm，与图4(b)基本一致。对比可知，与原药柱相比，算例A在从低燃速药柱向高燃速药柱过渡时，高燃速燃面 “膨大”现象更加明显。与低燃速燃面叠加后燃面增加显著，使得发动机压强出现陡增，在6.3 s时发动机压强达到14.7 MPa。若高低燃速的交界面继续向右移动，发动机转级过程产生的压强峰会进一步增大，最终导致发动机壳体结构失效。

图6 算例A在6.3 s时的燃面拓扑结构图7 算例B在10 s时的燃面拓扑结构Fig.6 Topological structure of burning surface at t=6.3 s for Case A Fig.7 Topological structure of burning surface at t=10 s for Case B

如图7所示，当算例B燃面推移至高低燃速交界面时，已趋于平面，此时燃面面积最小且燃速较低，使得燃烧室内的压强降低至0.58 MPa，在此压强下，推进剂不能维持正常燃烧，影响发动机工作可靠性。由以上分析可知，双燃速药柱交界面的位置直接影响着发动机内弹道性能的稳定性，因此需要在药型参数设计及药柱成型过程中严格控制。