固体火箭发动机声涡耦合数值仿真*

韩磊，卞云龙，张卫平，孙再庸

(1.内蒙动力机械研究所，内蒙古 呼和浩特　010010； 2. 中国航天科工集团第二研究院，北京　100038)



仿真技术

固体火箭发动机声涡耦合数值仿真*

韩磊1,2，卞云龙1，张卫平1，孙再庸1

(1.内蒙动力机械研究所，内蒙古 呼和浩特010010； 2. 中国航天科工集团第二研究院，北京100038)

摘要：研究固体火箭发动机工作中发生的声涡耦合导致的不稳定燃烧，以及燃烧室空腔变化对压力振荡的影响。对VKI缩比发动机流场采用大涡模拟(LES)方法进行数值模拟，确定数值方法可行。采用有限元和大涡模拟结合的方法，对某固体火箭发动机声涡耦合进行仿真分析，对声场分别采用理论方法和有限元数值方法对声模态和声学特性计算，确定有限元方法精度可靠及大涡模拟方法对流场旋涡流动不稳定捕捉准确，计算结果和试验中压力振荡结果相符。表明研究的该固体火箭发动机发生了声涡耦合引起的纵向1阶声不稳定，且获得了自由容积对压力振荡频率及幅值的影响规律。

关键词：声涡耦合；不稳定燃烧；压力振荡；数值模拟；声模态；有限元；大涡模拟

0引言

发动机燃烧工作过程中，由于某种随机扰动影响燃烧速率，燃烧室内流场和燃烧过程、空腔的反射、声场等相互作用，引起的不规则的燃烧室压力变化，称为不稳定燃烧，也常称为振荡燃烧[1]。压力振荡会导致较严重的推力振荡，轻则发动机振动，中断燃烧，严重的情况下，和飞行器耦合共振，导致系统可靠性降低、失效乃至灾难性的爆炸事故发生[2-3]。

近年来，随着固体运载要求的提出以及导弹技术发展，对固体火箭发动机提出了高装填、大推力比、初始大推力、大长细比等要求。国内的发动机设计者为了满足高性能指标，提高设计标准，导致战术地空导弹用固体火箭发动机屡次出现燃烧不稳定现象。除此之外，国外研制的大型助推器，如美国的SpaceShuttle和欧洲的Ariane5[4-5]采用的大型分段式固体发动机相继出现了较为严重的不稳定燃烧，不稳定燃烧研究也越来越重要[6]。

涡脱落产生的声学振荡被认为是不稳定燃烧问题中一个关键问题，最初Flandro和Jacobs[7]提出由于涡脱落激发的声模态是引起不稳定燃烧的一个重要因素，早期证据表明复杂几何装药的战术导弹上面级发动机的不稳定燃烧由旋涡脱落驱动。后来有研究表明声涡耦合导致的声加强剪切层不稳定的现象会放大了压力震荡，当涡脱落频率与声振频率相同时会形成了不稳定燃烧[8]。

1问题描述

本文同时采用理论计算和有限元结合的方法，针对在工作末期发生了不稳定燃烧的某固体火箭发动机不同燃烧时刻燃烧室声腔进行分析和数值计算。并对该发动机翼柱烧尽、振荡前、振荡时、振荡后等关键工作时刻建模，对不同时刻燃烧室流场采用大涡模拟方法数值模拟，初步分析认为由于流动不稳定性产生了涡脱落，导致声涡耦合，从而产生不稳定燃烧，研究该不稳定燃烧发生时刻是否发生声涡耦合，对耦合情况下振荡特性分析，并研究工作过程中自由容积变化对振荡特性的影响。

2数值方法与算例校验

2.1声学共振频率理论方法

燃烧室固有声学频率计算公式[9]:

(1)

(2)

(3)

式中：n为模态阶数；c为当地声速。

2.2声学共振频率有限元方法

采用基于有限元的方法,直接离散三维波动方程,可得到复杂装药声腔的声学特性。使用迦辽金法对简化的亥姆霍兹方程进行离散,声腔模态方程的单元矩阵[10]形式为

(4)

式中：Kf为声刚度矩阵；Mf为声质量矩阵；p为声特征向量；ωa为特征根。模型表面定义0位约束，使用有限元法求得圆周频率ωa进而求得声振频率fa。

2.3大涡模拟方法

通过在傅里叶空间或构型空间将随时间变化的N-S方程进行滤波可得到控制方程。考虑到气体的可压缩性，利用Favre平均对控制方程按式(5)简化。

(5)

本文不考虑化学反应，仅计算单组分工质，滤波后连续方程、动量方程与能量方程分别为

(6)

(7)

(8)

式中：“-”表示Reynolds平均;“～”表示Favre平均。

(9)

(10)

(11)

(12)

式中：常数Cε为1.1。

2.4大涡模拟算例校验

为研究障碍物导致涡脱落，VKI实验室设计了含有绝热环和潜入式喷管空腔缩比冷流实验发动机，Anthoine[12]还对该发动机进行试验研究。流场计算区域尺寸结构如图1所示。VKI发动机尺寸已知，且经过冷流实验和数值模拟验证，可以对本文采用的大涡模拟方法进行校验，为此选用此发动机进行方法校验。

图1　VKI几何模型Fig.1　VKI model and domain

旋涡具有三维特性，产生、拉伸过程均为三维空间过程。二维模型计算对旋涡的捕捉和模拟总有欠缺，但三维模拟计算量太大，Mason[13]和张峤[14]等人采用二维和三维大涡模拟对比后表明二维模拟可以成功的对ETM-03，RSRM，VKI等发动机进行数值模拟，且对振荡频率等特性很好的预估，为此本文综合计算精度和效率采用二维模型。

采用结构化网格，结构化网格质量好，易于计算收敛；方便控制局部区域网格稀疏，易于对流动复杂区域和关注区域网格加密；减少其他区域网格，从而减少计算时间。此处对障碍区域、产生旋涡流动区域、旋涡和喷管头部碰撞区域进行网格加密，局部加密网格如图2所示。

图2　局部网格加密Fig.2　Refined mesh zone

在发动机流场内选取关键点监测压力波动，入口为质量流率入口，温度为285K，质量流率为66.3kg/m2的理想气体，采用轴向进气[15],比热比为1.4，Pr数取0.71，出口为压力出口，压力为101 325Pa，温度为285K。流场初始计算先用稳态计算得到一个初始流场，待流场稳定后，改湍流模型为大涡模拟，为了避免连续方程和动量方程中心差分格式产生数值振荡，采用BCD格式进行离散，并同时改为时间步长为4e-7的非稳态计算，待流场稳定开始监测头部压力。

图3为部分时间段内发动机头部监测点压力曲线，由图可以看出，发动机内流场存在明显的压力振荡，压力振荡区间为(170 000～210 000Pa)，对压力曲线进行快速傅里叶变换后频率分布图如图4，在0～2 000BZ频率内，其中振荡主频为424BZ，压力振荡主要振荡频率为424,848,1 272,1 697BZ，呈明显倍频，主频压力振幅是其他振荡频率下振幅的7倍。

图3　压力振荡曲线Fig.3　Pressure-time oscillation curve

图4　频谱分析Fig.4　Pressure spectrum analysis

Anthoine等人对此发动机分别进行过数值模拟和试验，本文得到的振荡频率与其试验和数值模拟结果对比如图5和表1所示，四阶主频与Anthoine实验和数值结果吻合，其中一阶模态误差最大，最大误差仅为3.4%，表明本文所采用的大涡模拟方法对研究旋涡脱落导致流场不稳定可行。

图5　不同研究方法下频率对比Fig.5　Comparison of frequencies in differentmethods

模态阶数LES结果/BZ实验结果/BZ误差(%)14244103.428488753.13127212800.64169717402.5

3数值方法与算例校验

3.1计算模型

本文研究的发动机为大长径比固体火箭发动机，装药为前内孔后翼，为此选取以下4个关键工作时刻进行研究，发动机在6s时后翼烧尽，燃烧室空腔变化趋于均匀，在15.5s左右发生振荡燃烧，后翼烧尽后很长一段时间均会产生压力振荡，压力振荡频率为173BZ。在6～15.5s工作时间段内燃烧室内自由容积随着燃面的退移产生显著的变化，会导致燃通比等参数明显变化。为此，本节主要以后翼燃烧结束时刻6s，振荡发生前的时刻14s，压力振荡发生开始时刻15.5s，振荡后期时刻16s，4个时刻的发动机建模，计算模型如图6。

图6　关键时刻几何模型Fig.6　Key models and domain

压力监测点在流场区域面积最小的16s时刻模型中选取，监测点位置分别为：燃烧室头部，突扩前，突扩后，突扩段中间，突扩段尾部，如图7所示。

图7　压力监测点分布Fig.7　Distribution of pressure monitoring points

3.2边界条件与计算参数

燃面在不同时刻略有变化，燃烧室工作压强在不同燃烧时刻也不同，这些因素在本文研究中不是主要因素，为此忽略燃面、压强的细微不同，所有时刻燃面入口统一设置为质量流率入口，质量流量均为15.5s时刻的预估质量流率，出口为压力出口101 325Pa，出口温度300K，燃气参数见表2。喷管出口为超声速，出口界面参数由外推得到，且由于气体在喷管加速至超声速，下游对上游无影响，壁面采用无滑移壁面。

表2　燃气参数

3.3燃烧室声振频率理论计算

应用式(1)计算燃烧室声振频率如表3所示。

表3　理论计算纵向频率

3.4燃烧室声学模态有限元仿真计算

理论计算声腔声模态和振型时要求声腔几何构型简单，且边界条件明确，鉴于许多燃烧室初始声腔极复杂，本文中燃烧室相对简单，但依然存在阶梯凸台等特征，为此对燃烧室三维声腔建立有限元模型，喷管部分超声速，相当于燃烧室两端封闭，只选取燃烧室主体部分，流体单元选取fluid30，模型表明定义零位移约束，取声介质密度为4.0kg/m3,平均声速为1 061m/s。

图8　各阶振型声压分布Fig.8　Acoustic pressure distribution of every order axial vibration modes

图8中为某固体火箭发动机发生不稳定燃烧时前4阶纵向振型声压分布云图，前4阶声振频率如表4所示，一阶纵向频率165BZ与试验发动机压力振荡频率相近，声压分布中，发动机头尾部为声压波腹，表明发生了与一阶纵向声频耦合的声涡耦合。图9为有限元方法和理论计算方法下声振频率对比，对比结果表明随着模态阶数增大，两者计算误差也有所增大，最大误差约为5%。

表4　有限元计算声振频率

图9　不同研究方法下频率对比Fig.9　Comparison of frequencies with differentmethods

3.5燃烧室各时刻压力振荡分析

在发动机工作6s后翼柱烧尽，此后发动机内腔变化趋于均匀，不同时刻计算入口条件一致的情况下，待监测压力波动稳定后，观测得知从6s到最后16s过程中均会产生一定的压力振荡。6s时压力在波动中下降，14，16s时刻压力振荡在波动中上升，分别对以下4个时刻压力波动周期平均计算后发现压力波动周期处于5.86～6.26ms之间，压力波动频率均处于159～161BZ之间。若对压力波动滤去波动量后， 14，15.5，16s时压力波动类似正弦波动，15.5s时的压力波动频率如图10，周期性明显。尽管都存在压力振荡，可以看出14，15.5s时刻压力振荡幅值较大，15.5s压力振荡相比其他时刻明显高一个数量级，6s时由于后阶凸台明显，会产生明显的转角涡脱落，可能会导致头部压力波动较大，14s时的压力波动幅值是前2个时刻的10倍左右。15.5s压力振荡幅值远大于其他时刻，是前两时刻压力波动20倍左右。模拟结果得到的15.5s的压力振荡最大，和发动机实际工作情况相符。

图10　15.5 s时刻压力振荡Fig.10　Pressure fluctuation at 15.5 s

监测点head处的平均压力随发动机工作时间变化如图11所示，随着发动机工作，自由容积的变大，平均压力先逐渐降低，在振荡发生时刻附近上升，15.5s时最大，然后下降，平均压力的最大波动范围在5%左右。

对压力时间曲线，进行开窗快速傅里叶变换后，得到频率、振幅曲线， 均出现压力振荡， 如图12所示，15.5s时振幅最大，波峰对应的前2阶频率分别为167，332BZ，其167BZ与一阶纵向声频161BZ和165BZ接近。随着时间变化，燃烧室自由容积增大，1阶压力振荡频率几乎不变，其他频率峰值逐渐平移，待发生压力和声振耦合共振后，随着时间推移，产生共振的压力振荡频率消失，到燃烧室工作末期，压力振荡几乎消失。在燃面药形设计中避开14～15.5s时刻的变化压力振荡明显时的空腔，即可只发生低频、小振幅微弱振荡，而避免发生共振。

图11　不同时刻平均压力Fig.11　Average pressure at different times

图12　不同时刻振荡频率和振幅Fig.12　Pressure spectrum at different times

4结束语

大涡模拟和有限元方法结合的方法可以对某火箭发动机数值仿真与分析，表明该发动机发生了纵向1阶频率主导的声涡耦合声不稳定燃烧，压力振荡的时刻特性与发动机试车工作情况相符。随着发动机工作，自由容积增大的过程中，平均压力先增大最然后下降，在振荡发生时达到峰值，平均压力的最大波动范围在5%左右。表明本文方法可以对声涡耦合引起的不稳定燃烧情况进行预计，可以为发动机设计和抑振提供参考校验。

参考文献：

[1]孙维申.固体火箭发动机不稳定燃烧[M].北京：北京工业学院出版社，1987.

SUNWei-shen.CombustionInstabilityinSolidRocketMotors[M].Beijing:BeijingInstituteofTechnologyPress, 1987.

[2]刘佩进，金秉宁，李强.战术导弹固体发动机燃烧不稳定研究概述[J].固体火箭技术，2012，35(4)：446-456.

LIUPei-jin,JINBing-ning,LIQiang.ASurveyofCombustionInstabilityinTacticalSRM[J].JournalofChinaOrdance，2012，35(4)：446-456.

[3]王宁飞，张峤，李军伟，等.固体火箭发动机不稳定燃烧研究进展[J].航空动力学报,2011,26(6)：1406-1412.

WANGNing-fei,ZHANGQiao,LIJun-wei,etal.ProgressofInvestigationonCombustionInstabilityofSolidRocketMotors[J].JournalofAerospacePower,2011,26(6):1406-1412.

[4]BLOMSHIELDFS.HistoricalPerspectiveofCombustionInstabilityinMotors:CaseStudies[R].AIAA2001-3875.

[5]胡大宁，何国强，刘佩进, 等. 翼柱型药柱固体火箭发动机不稳定燃烧研究[J].固体火箭技术，2010，33(5)：502-506.

HUDa-ning,HEGuo-qiang,LIUPei-jin,etal.StudyonInstableCombustionofSolidRocketMotorwithFinocylGrain[J].JournalofChinaOrdance,2010, 33(5)：502-506.

[6]BLOMSHIELDFS.SummaryofMulti-DisciplinaryUniversityResearchInitiativeinSolidPropellantCombustionInstability[R].AIAA2000-3172.

[7]FLANDROGA,JACOBSHR.Vortex-GeneratedSoundinCavities[R].AIAA1973-1014.

[8]熊永亮，杨丹，郜冶，等.固体火箭发动机燃烧室内声波振荡数值研究[J].推进技术,2006,27(6):484-487.

XIONGYong-liang,YANGDan,GAOYe,etal.NumericalStudyofAcousticOscillationinSolidRocketMotors[J].JournalofPropulsionTechnology,2006,27(6):484-487.

[9]唐金兰.固体火箭发动机原理[M].北京：国防科技出版社,2013.

TANGJin-lan.SolidRocketMotorMechanical[M].Beijing:NationalDefenseTechnologyPress,2015.

[10]杨向明,刘佩进,陈晓龙.翼柱型装药固体火箭发动机燃烧室声场分析[J].宇航学报，2008，29(5):1593-1597.

YANGXiang-ming,LIUPei-jing,CHENXiao-long.AnalysisofAcousticPropertywithComplicatedGrainShapeinCombustionChamberofSolidRocketMotor[J].JournalofAstronautics, 2008, 29(5): 1593-1597.

[11]NICOUDF,DUCROSF.Subgrid-ScaleStressModelingBasedontheSquareoftheVelocityGradientTensor[J].Flow,TurbulenceandCombustion,1999,62(3):183-200.

[12]ANTHOINEJ,BUNCHLINJM,GUERYJF.EffectofNozzleCavityonResonanceinLargeSRMNumericalSimulation[J].JournalofPropulsionandPower,2003,19(3);374-384.

[13]MASONDR,MORSTADTRA,CANNONSM.PressureOscillationandStructuralVibrationsinSpaceShuttleRSRMandETM-3Motors[R].AIAA2004-3898,2004

[14]张峤，李军伟，王伟臣，等.固体火箭发动机涡声耦合特性数值研究[J].推进技术，2011，32(3):349-354.

ZHANGQiao，LIJun-wei，WANGWei-chenetal.NumericalAnalysisonVortexAcousticCouplingCharacteristicsinSolidRocketMotors[J].JournalofPropulsionTechnology,2011,32(3):349-354.

[15]STELLAF,PAGLIAF.PressureOscillationinSolidRocketMotorsNumericalStudy[J].AerospaceScienceandTechnology,2011,15(1):53-59.

Numerical Simulation on Acoustic Vortex Coupling in Solid Rocket Motors

HAN Lei1,2, BIAN Yun-long1, ZHANG Wei-ping1, SUN Zai-yong1

(1. Dynamic Machinery Institute of Inner Mongolia, Inner Mongolia Hohhot 010010, China；2.TheSecondAcademyofChinaAerospaceandIndustryCorporation,Beijing100038,China)

Abstract:In order to study solid rocket motors (SRM) combustion instability induced by acoustic vortex coupling and the effect of combustion chamber free volume on pressure oscillation. Large-eddy-simulation (LES) is used to simulate flow field domain of VKI scale down motor to show that this numerical method is capable of application. A method combining finite-element-method (FEM) and LES is applied to simulate and analyze acoustic vortex coupling for a certain motor. Theoretical method and FEM are used to calculate acoustic modes and acoustic characteristics, and it shows that the accuracy of FEM is reliable. For catching of the flow instability detail, LES is accurate. Pressure oscillation of simulation is consistent with the experiment result. And the result shows that the first-order axial mode acoustic instability induced by acoustic vortex coupling occurs in the SRM. The effect of the free volume on pressure oscillation is also acquired.

Key words:acoustic vortex coupling；unstable combustion；pressure oscillation ；numerical simulation；acoustic mode；finite-element-method(FEM)；Large-eddy-simulation(LES)

*收稿日期：2015-05-20；修回日期：2015-07-06

作者简介：韩磊(1989-)，男，宁夏银川人。硕士生，主要从事固体火箭发动机数值仿真。

通信地址：010010内蒙古呼和浩特市新华东街65号航天科工六院41所E-mail：docarezone@qq.com

doi:10.3969/j.issn.1009-086x.2016.03.030

中图分类号：V435；TP391.9

文献标志码：A

文章编号：1009-086X(2016)-03-0194-07