维普资讯 http://www.cqvip.com 第２７卷第３期 弹箭与制导学报 ・ １７９ ・ 基于燃面推移的内流场数值仿真 谢丽宽，马 拯，俞红博，于冀国 （海军驻西安地区军事代表室，西安 ７１００２５） ［摘要］基于处理移动边界的动网格技术．采用隐式的、非定常的耦合解法，选用 ￡双方程的Ｒｅａｌｉｚａｂｌｅ湍流 模型和壁面函数法．对某小型固体火箭发动机粘性的、可压缩的、高雷诺数的燃烧室一喷管湍流流场进行了一 体化的数值计算．得到其在移动边界条件下瞬态的流场分布．并分析了其完整的内弹道各参数随着装药燃面 推移而连续变化的规律．通过数据拟合出的三条内弹道特性参数曲线与试验数据符合较好．说明文中所采用 的内弹道数值仿真计算方法可靠。 ［关键词］固体火箭发动机；数值模拟；燃面 ［中图分类号］Ｖ４３５ ［文献标志码］Ａ Ｂａｓｅｄ ｏｎ ａ Ｍｏｖｉｎｇ Ｂｕｒｎｉｎｇ Ｓｕｒｆａｃｅ ｏｆ Ｎｕｍｅｒｉｃａｌ Ｓｉｍｕｌａｔｉｏｎ ｔｏ ｔｈｅ Ｗｈｏｌｅ Ｆｉｅｌｄ ｉｎ ｔｈｅ Ｃｈａｍｂｅｒ ａｎｄ Ｎｏｚｚｌｅ ＸＩＥ Ｌｉ—ｋｕａｎ．ＭＡ Ｚｈｅｎｇ．ＹＵ Ｈｏｎｇ—ｂｏ．ＹＵ Ｊ ｉ－ｇｕｏ （Ｍｉｌｉｔａｒｙ Ｒｅｐｒｅｓｅｎｔａｔｉｖｅｓ Ｏｆｆｉｃｅ ｏｆ Ｎａｖｙ ｉｎ Ｘｉ’ａｎ Ａｒｅａ，Ｘｉ’ａｎ ７１００２５．Ｃｈｉｎａ） Ａｂｓｔｒａｃｔ：Ｂａｓｅｄ ｏｎ ｔｈｅ ｔｅｃｈｎｉｑｕｅ ｆｏｒ ｄｙｎａｍｉｃ ｍｅｓｈｅｓ ｄｅａｌｉｎｇ ｗｉｔｈ ｍｏｖｉｎｇ ｂｏｕｎｄａｒｙ ｐｒｏｂｌｅｍｓ，ａ ｃｏｕｐｌｉｎｇ ａｎｄ ｉｍｐｌｉｃｉｔ ｓｏｌｕｔｉｏｎ ｏｆ ｕｎｓｔｅａｄｙ—ｓｔａｔｅ ｆｌｏｗｓ ｉｓ ｐｅｒｆｏｒｍｅｄ ｗｉｔｈ ｔｈｅ Ｒｅａｌｉｚａｂｌｅ ￡ｔｕｒｂｕｌｅｎｃｅ ｍｏｄｅｌｓ ａｎｄ ｗａｌｌ ｆｕｎｃｔｉｏｎｓ ｉｎ ｔｈｉｓ ｐａｐｅｒ． Ａ ｉｎｔｅｇｒａｔｉｖｅ ｎｕｍｅｒｉｃａｌ ｓｉｍｕｌａｔｉｏｎ ｈａｓ ｂｅｅｎ ｃａｒｒｉｅｄ ＯＵｔ ｔＯ ｉｎｖｅｓｔｉｇａｔｅ ｔｈｅ ｆｕｌｌ ｔｕｒｂｕｌｅｎｔ ｆｌｏｗｓ ｉｎ ｔｈｅ ｃｈａｍｂｅｒ－ｎｏｚｚｌｅ ｆｉｅｌｄ ｏｆ ｓｏｍｅ ａ ｓｍａｌｌ—ｓｃａｌｅ ｓｏｌｉｄ ｒｏｃｋｅｔ ｍｏｔｏｒ（ＳＲＭ）．ｗｈｉｃｈ ｉｎｔｅｒｎａｌ ｆｌｏｗ ｆｉｅｌｄｓ ｉｓ ｖｉｓｃｏｕｓ．ｃｏｍｐｒｅｓｓｉｂｌｅ ａｎｄ ｈｉｇｈ—Ｒｅｙｎｏｌｄｓ— ｎｕｍｂｅｒ．Ａｎｄ ｔｈｅｎ．ａ ｔｉｍｅ－ｄｅｐｅｎｄｅｎｔ ｆｌｏｗｓ ｉｓ ｏｂｔａｉｎｅｄ ｗｉｔｈ ｄｙｎａｍｉｃ ｂｏｕｎｄａｒｉｅｓ，ｒｅｓｕｌｔｉｎｇ ｉｎ ｔｈｅ ａｎａｌｙｚｉｎｇ ｏｆ ｉｎｔｅｇｒｉｔｙ ｒｅｌａｔｉｏｎｓ ａｍｏｎｇ ａｌｌ ｐａｒａｍｅｔｅｒｓ ｖａｒｙｉｎｇ ｃｏｎｔｉｎｕｏｕｓｌｙ ｗｉｔｈ ｔｈｅ ｃｏｍｂｕｓｔｉｏｎ ｏｆ ｔｈｅ ｐｒｏｐｅｌｌａｎｔ ｓｕｒｆａｃｅ．Ｔｈｅ ｔｈｒｅｅ ｉｎｔｅｒｎａｌ ｂａｌｌｉｓｔｉｃ ｔｒａｊｅｃｔｏｒｙ ｐｒｏｐｅｒｔｙ ｃｕｒｖｅｓ ｂｙ ｔｈｅ ｄａｔａ ｆｉｔｔｉｎｇ ａｇｒｅｅ ｗｅｌｌ ｗｉｔｈ ｔｅｓｔ ｄａｔａ ｃｕｒｖｅｓ，ｓｈｏｗｉｎｇ ｔｈａｔ ｔｈｅ ｎｕｍｅｒｉｃａｌ ｓｉｍｕ— ｌａｔｉｏｎ ｍｅｔｈｏｄ ｏｆ ｉｎｔｅｒｎａｌ ｂａｌｌｉｓｔｉｃ ｐｒｏｐｅｒｔｉｅｓ ｏｆ ｓｏｌｉｄ ｒｏｃｋｅｔ ｍｏｔｏｒ ｉｓ ｒｅｌｉａｂｌｅ． Ｋｅｙ ｗｏｒｄｓ：ｓｏｌｉｄ ｒｏｃｋｅｔ ｍｏｔｏｒ；ｎｕｍｅｒｉｃａｌ ｓｉｍｕｌａｔｉｏｎ；ｂｕｒｎｉｎｇ ｓｕｒｆａｃｅ Ｒｅａｌｉｚａｂｌｅ湍流模型和壁面函数法，对某小型固 １ 引言 体火箭发动机粘性的、可压缩的、高雷诺数的燃 固体火箭发动机内弹道性能参数的统计指 烧室一喷管湍流流场进行了一体化的数值计算， 标（期望值和方差），主要通过对发动机的试验数 得到其在移动边界条件下瞬态的流场分布，并分 据进行统计处理得到。固体火箭发动机的价格 析了其完整的内弹道各参数随着装药燃面推移 和试验费用昂贵，一个型号的固体火箭发动机在 而连续变化的规律。 研制过程中往往只能进行很少的几次试验，可供 发动机内弹道性能指标统计用的样本数较小，得 ２计算模型及网格划分 到的指标误差相应较大。通过固体火箭发动机 ２．１ 计算模型 内弹道性能数字仿真得到的内弹道性能指标，不 带旋流修正的 ｅ模型（即Ｒｅａｌｉｚａｂｌｅ ｅ模 仅花费较少，而且还可以提供较多的样本数。 ），旋流修正意味着模型在雷诺压力中有数学 文中基于处理移动边界的动网格技术，采用 约束，以保证湍流的连续性。带旋流修正的 ｅ 隐式的、非定常的耦合解法，选用 ｅ双方程的 模型与标准志一ｅ模型相比有两个主要的不同点： ＊收稿日期：２００６一Ｏ７—１７； 修回日期：２００６—１０—１８ 作者简介：谢丽宽（１９８０一）．男，浙江瞿州人，硕士研究生，研究方向：导弹动力系统仿真。 维普资讯 http://www.cqvip.com ・ ｌ８Ｏ ・ 弹箭与制导学报 为湍流粘性增加了一个公式；为耗散率增加了新 的传输方程（这个方程来源于一个为层流速度波 动而作的精确方程）。带旋流修正的肛ｅ模型比 标准ｋ－ｅ型在强流线弯曲、漩涡和旋转有更好的 表现。而且它对于旋转流动、强逆压梯度的边界 ２．２ 计算区域与计算网格 对上述计算模型，其计算区域为一具有运动 边界的、时变的流场区域，在计算过程中该变化 着的流场区域是通过编制燃面推移程序实现的， 任取四个时刻如图ｌ阴影所示（箭头表示边界移 动方向）。随着推进剂装药的燃面以某种特定的 燃速向后燃烧退移，所计算的流场区域不断增 层流动、流动分离和二次流有很好的计算结果。 壁面函数法是一种仅用来计算在流动中壁 面对其附近网格点影响的经验函数Ｌ２ｊ。它减少 了在湍流边界层中完全求解整个区域的必要，并 允许靠近壁面的网格点放置到离壁面相对较远 的地方。因此壁面函数的引用，大大提高了在湍 流模拟中的计算效率。 近壁处的速度分布可用： Ｕ 一 １ Ｉｎ（ ） （１） 式中 三 （２） 其中：是为冯・卡门（ｖｏｎ Ｋａｒｍａｎ）常数，其值为 ０．４２。Ｅ为经验常数，大小为９．８ｌ。【， 为流场中Ｐ 点的平均速度，志 为Ｐ点的湍流动能，Ｙ 为Ｐ点 距墙面的距离， 为流场的动粘度。上述法则适 用于Ｙ ＞ｌ１．２２５范围，对于其它区域，其速度 分布为Ｕ 一 。 温度分布由下式来表示： Ｔ＊三！！： 二！： 一 Ｐｒ ＋ １ ＰＰｒ ≯ ｛ ，［÷ ｃ ］＋ ｆ专ｌＤ ｛ ，Ｌ，刍＋（Ｐｒ—Ｐｒ，） ｝（ ＞ ） （３） 式中Ｐ由下式来计算： Ｐ一 （ （ 一 ）（鲁）“ （４） 各参数的意义为： ／为流体的热导率， 为 墙的热通量， 为靠墙面网格上Ｐ节点的温度， 丁 ．为墙壁温度，ｎ为摩尔普朗特数，ｎ，为湍流 普朗特数，在墙壁处为０．８５，Ａ为Ｖａｎ Ｄｒｉｅｓｔ常 数，【， 为在Ｙ 一 处的平均速度大小。 大，燃气的流动发生在变化着流场内。 二：…………………一……二二二 二二：＿Ｉ 二：二二。 一…… 乙 ＿一＿Ｉ 图ｌ 燃面推移乐意图 文中采用了区域分区、强制类型转换、自适 应、边界层设置、局部加密和区域扩充等技术。基 本做法是：首先根据计算区域几何形状采用区域 分区技术，如喷管部分共分为７个分区；依据网 格划分的具体要求更改节点类型；根据对流场特 性的预测，采用边界层设置技术和局部加密技术 重新布置网格分布；对流场参数变化剧烈的区域 如喉部，采用局部加密技术；最后再根据试算的 结果，用自适应方法修改网格。二维网格单元有 三角形、四边形、五边形或楔形网格等形状，考虑 到四边形网格几乎对所有算法都适用，故文中采 用结构化的四边形网格。在生成网格时，为了减 弱交叉扩散（非正交源项，要控制网格单元的倾 斜比，尽量使四边形单元的短边与流线垂直。 给定初始时刻的网格划分（见图２），随着燃 面向左推移时，下一时刻的网格划分由上一时刻 网格系统通过程序（ＵＤｔｅ）实现，所增加的网格 单元层厚度设为ｌ０ｍｍ。 图２ 初始时刻网格系统 维普资讯 http://www.cqvip.com 第２７卷第３期 基于燃面推移的内流场数值仿真访｝丽宽等 ・ ｌ８ｌ ・ ３ 边界条件与初始条件 ３．１ 边界条件 燃烧室一喷管一体化计算的边界条件主要 有：燃面边界条件、固体壁面边界条件、对称平面 或对称轴线边界条件、喷管出口边界条件。 ａ．燃面边界条件 燃面边界条件分为加质边界条件和移动边 界条件，燃面加质质量流率为： ｍ— Ｐ ｐ—Ｐｐａｐ （ｓ） 式中：ｒ＾为推进剂燃速：助推级４０．２２６ｍｍ／ｓ，续 航级２Ｏ．７５ｍｍ／ｓ； —１．７８ｇ／ｃｍ。。边界移动速 度与燃气加质量是耦合的，不同的燃速、压强决 定了燃气加质量的大小，此关系通过编制程序 ｍａｓｓ—ｆｌｕｘ．ｃ实现川三０．３３，在设置燃面边界条 件Ｈ，－ｉ￣ｎ载。 ｂ．固体壁面边界条件 固体壁面采用无滑移条件，且假设壁面绝 热，气体同固壁之间量及热量的交换，也无 组分的交换。其边界条件为： “一０， 一０， 一０， 一０ （６） ｄ ｄ７ Ｃ．对称平面或对称轴线边界条件 由于对称性，对称轴线上物理量的法向梯度 必须为零，或者其本身为零。即：对于所有标量 （１０、 、￡、 等）和平行于对称轴线矢性量 如速度分量、法向应力等，有 一０；而对于垂 ｏｙ 直于对称轴线的参量Ａ如速度分量、标量通量 、剪应力等，其值本身为零，即Ａ一０。 ｄ．喷管出口边界条件 喷管出口为超音速区，全部参数由邻近的上 游参数二阶外推插值而得 ］，如喷管出口压力由 外插求得；亚声速时需给定一个边界条件，一般 出口压力取为外界静压。 ３．２ 初始条件 对于非稳态、湍流流场计算来说，初场的给定 是一个比较关键的问题，若初场偏离收敛解过大， 则易造成迭代发散。文中先通过对稳态的、固定 边界的流场进行计算，所得结果作为非稳态的、移 动边界的流场的初值，给定总温为３４８８．５０Ｋ，燃 气气体分子量为２４．９６６，定压比热容为１００６．４３Ｊ／ ・Ｋ，导热率按照Ｋｉｎｅｔｉｃ－ｔｈｅｏｒｙ进行计算，粘性 系数用Ｓｕｔｈｅｒｌａｎｄ公式来计算 。 ４ ＵＤＦ求解流程 文中是用耦合法，在耦合求解器求解过程中 （见幽３），在迭代循环开始之前定义初始化函数 （同样用ＤＥＦＩＮＥ—ＩＮＩＴ定义）。然后迭代循环 开始执行用户定义的调整函数（ＤＥＦＩＮＥ—ＡＤ— ＪＵＳ－Ｉ、定义的）。ＦＩ ＵＥＮＴ同时求解连续、动量 和（适合的地方）能量的控制方程（同时解一套物 质输运或矢量方程），其余的求解步骤与分离式 求解相同。且计算模型必须基于如下假设： （１）流场按二维轴对称计算，燃气在燃烧室、 喷管中的流动为一体化计算； （２）装药燃速为双燃速（不同时间段燃速不 同），不考虑燃面退过程中的化学反应，燃烧等效 为燃面加质； （３）燃气等效为具有某种热力学性质与输运 性质的单一理想气体； （４）不考虑两相流，燃气为纯气相； （５）不考虑辐射换热、壳体烧蚀及导热。 图３ ＵＤＦ求解流程图 ５ 结果与讨论 燃烧室压力取燃烧室出口处的总压值（但是 整个截面上总压并不一致均匀），即初始燃面处 的总压值；推力由燃烧室压力按Ｆ—Ｃ Ｐ Ａ 求 维普资讯 http://www.cqvip.com ・ ｌ８２ ・ 弹箭与制导学报 得，其中由ｄ／ｄ，＝３．２１７ ４查表得ＣＦ—１．４３６， 流量取喷管人口截面处值。 ＿ ＾整个计算过程共２７０个时间步，每个时间步 芝 为０．２ｓ。在取时间点时借鉴“稳态非线性问题迭 代一个层次，相当于线性非稳态问题前进一个时 层”ｌｓ ，并仿照ＴＵＳ格式向下一时刻取一点（对 缝 应于空间层次上的上游）。０ｓ取流场初值下的收 蛩 敛值，０．２ｓ．０．４ｓ、０．６ｓ、０．８ｓ、１ｓ（这些为工作时间 而不是计算运行时间）为已得收敛值的下平均， 随后１～１１．９ｓ之间各时间点为前３点、本身、下 一时刻共５个时间步收敛值的平均，ｌ１．９ｓ、 ｌ２．１Ｓ、１２．３ｓ的值为前点和本身值的加权平均， 随后ｌ２．５～５１．７ｓ取点方法同ｌ～儿．９ｓ；５１．９ｓ 后为拖尾段，其压力和流量分别为Ｌ６ ： ｃ。一丁ｘ／Ｆ￣ｆ 式中：Ｔ，取燃烧室出口处值，ｒ取０．６６８０。 按以上方法每个参数共得到６Ｏ个数据点， 再拟合可得燃烧室压力一时间、发动机推力一时 间和质量流率一时间三条曲线，如图４中虚线所 示。图中实线为４５℃时某次发动机试车实测值。 通过对比可以看出： （１）在平衡段、收尾段三个特征参数随时问 变化趋势一致（启动段数值计算模型没有考虑， 而仅是给初场值计算所得）； （２）平衡段三个特征参数的计算值比实测值 略有偏低，但都可认为在正常的波动范围内，实测 值波动很可能是装药的燃烧特性的不稳定所致； （３）收尾段三个特征参数的计算值明显低于 实测值，主要原因是：试车发动机在收尾段残药 比较多； （４）助推级转续航级时，计算所得曲线比较 陡直，主要是因为：一是在文中的数值计算中，当 燃速改变时没有比较细致地考虑燃烧的理化模 型，使得计算值比实际值低，这合乎预判；另外也 跟数据点的抽样方式有关，虽然抽样时己考虑了 ．． 三ｌ２ 杂ｆ！ ‘２ ６ １ ０ ｌ ４ ７ ｌ０ Ｉ ３ Ｉ６ Ｉ ９ ２２ ２５ ２８３Ｉ ３４ ３７４０４３４６４９ ５２ ５５５８ ／／ｓ 图４ 内弹道参数曲线与试验数据的拟合 前后时刻的影响，但仅作了简单的加权平均处 理，导致结果有些偏差。 （５）通过数值计算拟合出的三条内弹道特性 参数与实际情况基本吻合，说明文中所采用的计 算方法可靠、计算结果可信。 ［参 考 文 献］ ［１］方丁酉，张为华，杨涛．固体火箭发动机内弹道学 ［Ｍ］．长沙：国防科技大学出版社．１９９７． ［２］苏铭德．黄素逸编著．计算流体力学基础［Ｍ］．北 京：清华大学出版社，１９９７． ［３］ 候晓．固体火箭发动机两相湍流的数值研究［Ｄ］． 西安：西北工业大学．１９９０． ［４］叶万举．常显奇．曹泰岳编．固体火箭发动机气体 工作过程理论基础［Ｍ］．长沙：国防科技大学出版 社，ｉ９８４． ［５］Ｓｅｕｎｇ Ｗｏｏｋ Ｂａｅｋ．Ｉｎ Ｈｙｕｎ．Ｎｕｍｅｒｉｃａｌ ａｎａｌｙｓｉｓ ｏｆ ｉｇｎｉｔｉｏｎ ｔｒａｎｓｉｅｎｔ ｉｎ ａｎ ａｘｉｓｙｍｍｅｔｒｉｃ ｓｏｌｉｄ ｒｏｃｋｅｔ ｍｏｔｏｒ ｅｑｕｉｐｐｅｄ ｗｉｔｈ ｒｅａｒ ｉｇｎｉｔｉｏｎ ｓｙｓｔｅｍ［Ｒ］． ＡＩＡＡ ９７—２７ｌ５． ［６］宜敏，张中钦，等编著．固体火箭发动机原理［Ｍ］． 北京：国防工业出版社．ｉ９８５．
