燃烧室热-声激励及响应的模拟研究

臧 也，田 晶，张凤玲，艾延廷，王 志(沈阳航空航天大学 辽宁省航空推进系统先进测试技术重点实验室，沈阳 110136)

燃烧室是航空发动机的重要组成部件，其火焰筒为薄壁结构，内部流场极为复杂。燃料在燃烧室中燃烧，经由涡轮和尾喷口将内能转化为机械能从而产生推力[1]。燃烧会产生温度与压力的波动，甚至产生热-流-固耦合振荡，引起结构疲劳破坏。研究燃烧室内部声压激励及结构响应情况，对分析航空发动机燃烧室结构疲劳破坏产生机理,改进燃烧室结构设计具有重要意义。

近年来，国内外学者针对燃烧室内部燃烧及其影响做了大量研究。Valter Bellucci[2]等人对燃烧室放热过程进行数值模拟，分析了火焰与声压之间的耦合关系；Jaap F.Van Kampen等人[3]建立了一维声学模型，研究热声不稳定性；Sattelma-yer等人[4]对环形燃烧室和单个燃烧筒内的旋流流场进行研究，分析了不同流场对燃烧放热率的影响，并建立了火焰模型;Prakash等人[5]研究了不同燃气当量比下火焰的动态特性，分析了声压与火焰放热率之间的相互影响;K.D.Artur等[6]进行预混燃烧下的热-声不稳定性分析，研究了温度场与声压力场的相互影响关系;Huls R.等[7]基于声-弹性有限元模型和声振试验研究火焰筒壁面振动，并较好地预测了燃烧室内的声振水平;舒歌群等[8-9]应用可视化高速摄影技术和声学仿真软件开展了单点和多点激励下的燃烧室压力波瞬态响应计算，研究了燃烧室内的热声现象及压力场变化规律。在对燃烧室进行的众多研究中，研究重点主要集中在火焰波动与燃烧噪声之间的互相影响，准确获取燃烧室内燃烧的压力曲线对分析燃烧噪声有着重要的意义，利用试验手段测试得到发动机燃烧室内的压力曲线来研究和分析燃烧噪声是目前常用的方法。而一般使用各种模拟软件计算得到的燃烧室内压力曲线都是光滑曲线，如果使用模拟软件计算后的光滑曲线评估燃烧噪声，则会带来很大的失真性。燃烧室实际噪声频率范围宽，可视为多种频率声波的叠加，所以本文取单种频率的声波进行瞬态模拟研究分析。通过将入口条件设置为脉冲式入口，通过CFD软件模拟燃烧过程及声波传播过程，从而获得单一频率下的高频振荡的燃烧室内压力曲线，从而分析不同工况下的热声激励变化及燃烧室壁面响应结果。

1 结构建模及计算原理

1.1 燃烧室结构建模

燃烧室结构如图1所示，燃烧室整体长为2 000 mm。燃烧室入口分为空气入口与燃气入口，其中空气入口直径为47 mm，两个燃气入口直径均为12 mm。综合考虑运行时间与计算精度，在距离入口1 000 mm处取一段长为400 mm的结构段作为研究对象。本文模拟燃气进气质量流量不变，通过给定脉冲燃料入口创造不稳定燃烧条件，产生压力波动模拟声波与燃烧室壁面的双向耦合。燃烧室实际噪声频率范围宽，可视为多种频率声波的叠加，本文取一种频率的声波进行研究分析。

图1 结构示意图

数值研究中常采用的网格类型主要分为结构化网格与非结构化网格。结构化网格允许使用较少的网格单元并且计算量小，而非结构化网格具有更高的数值耗散率，对复杂几何模型流动更有利[10]。综合以上因素，燃烧室结构段计算采用结构化网格，其计算网格数为4800，计算节点数为34200；内部流场使用非结构化网格，其计算网格数为47248，计算节点数为52756。本文采用理论模型计算出结构段各阶固有频率，并针对该结构制作实际模型并进行模态实验，如图2所示。

图2 燃烧室结构段试验件和试验整体布局

计算结果与实验结果对比如表1所示，各阶固有频率的计算结果与模态实验结果基本吻合，说明选取的模型及网格划分合理，可以用于计算。

表1 结构段前六阶固有模态频率对比

1.2 数值模拟流程

在本文的数值模拟过程中，结构温度场响应速度远低于流场温度变化速度，故将稳态温度结果作为固体域初始边界条件，计算中边界面采用绝热计算条件。本文通过ANSYS Workbench数据交换平台提供的System Coupling模块实现燃烧室结构段与内部流场的双向耦合[11]。燃烧室结构段变形分析采用ANSYS软件中的瞬态结构动力学分析方法，流场特性采用CFX中SST湍流模型进行分析。燃烧室与燃烧流场双向流固耦合求解采用双重循环迭代方法[12]。耦合计算流程如图3所示。将所处工况下的稳态温度场作为初始条件导入模型，根据给定的进出口压力对动量方程进行求解，确定流场中各节点的初始压力和速度，将过程的总时间分成若干个时间步，在每个时间步中做稳态计算。

图3 计算流程图

单步耦合计算原理过程如图4所示。在某一时刻循环开始，以上一时刻流场压力p、速度v的分布和燃烧室结构段变形位移x作为初始条件进行计算，在流体域进行若干步计算收敛后，再通过网格插值将流场压力p、速度v的分布信息传递到流固耦合面。燃烧室固体耦合面以其为边界条件计算得到燃烧室结构瞬态动力响应，然后结构变形的位移x等信息再通过网格插值传递给流场耦合面作为流场耦合面的边界条件，当流体域与固体域的位移、载荷都达到收敛状态时，则完成一次双向耦合迭代计算，然后继续进入下一时刻循环。通过双向流固耦合计算，可获得任一时刻燃烧流场压力、速度及温度分布特性和燃烧室结构运动变形特性。

图4 单步耦合计算原理图

1.3 声压波模拟方法及验证

在燃烧室中，压力振荡的根源是燃烧的不稳定性，可将火焰视为声源。如果火焰在诸如管或燃烧室的密闭环境中燃烧，则产生的声波从边界反射。火焰对声波有很大影响，既可以削弱声波也可以进一步对声波进行放大。火焰放大某处声波的条件公式[13]为：

(1)

公式(1)表明，当声压振荡和放热振荡之间的相移小于900时，相关的反馈机构产生声能，反之，声波被削弱。因此，当评估燃烧的不稳定性时，关于燃烧室中声场的火焰相位行为信息非常重要。在实际情况下，声能在入口、出口，通过壁的振动和边界层中被吸收。这意味着关于火焰振幅特性的信息也很重要。当能量的产生大于耗散时，在每个振荡周期期间，系统中声能被增强，并且由于非线性效应，压力振幅增加直到饱和。燃烧噪声可由波动方程来表达[14]

(2)

其中，t是时间，p是声压，ρ是平均流量密度，c是声速,γ是比热比，Q′为单位质量工资的热释放量。对于具有刚性壁的长棱柱管，声波传播方程可简化到一维形式

p(x,t)=A1ej(ωt-kx)+A2ej(ωt-kx)

(3)

设波沿x方向传播如图5所示，声压表示为p(x,t)，

图5 声波传播示意图

(4)

波动方程为

(5)

将声波传播方程简化为一维形式

p(x,t)=A1ej(ωt-kx)+A2ej(ωt-kx)

(6)

其中，A1和A2分别是在正负x方向上以平均声速c0行进的声压波的振幅，假设不存在反射波，则A2=0。k为角频率和声速之比，即k=ω/c0=2π/λ。再设x=0处的声源振动时，在毗邻介质中产生pmejωt的声压,则声压的振幅值为pm，声场中的声压为p=pmej(ωt-kx)。

本文将由燃烧不稳定产生的压力波动视做声压波动，通过压力波与燃烧室壁面进行互相作用实现燃烧室内部流场与结构的双向耦合计算。应用CFD软件模拟声波传递过程，在燃烧室流场入口处施加正弦压力波动激励，设置流场中流体域为理想空气。模拟结果与通过公式计算得到理论结果对比，如图6所示。

由图6可知，火焰波前缘声波传递规律符合理论计算结果，不同位置的声波运动曲线的计算结果均与理论计算结果相符，说明所选择的网格密度及单步计算时间符合声波传递计算需要，此外也证明采用CFD技术分析声波传递是可行的。

图6 声波传递理论结果与模拟结果

1.4 双向耦合计算边界条件

图7给出了双向耦合建模的边界条件示意图。本文对燃烧室的燃烧过程进行模拟，各计算参数如表2所示。模型流体域入口采用速度入口，燃料为甲烷[15]。为了模拟持续的声压波动过程，将燃气入口分别设置为不同频率的脉冲入口[16]。

图7 边界条件示意图

表2 燃烧室工作状态参数

出口采用压力出口，相对压力为0。湍流模型采用SST模型。燃烧室壁面与内部流场耦合面如图7中标签处所示，流体域耦合面为与固体域耦合面接触的对应表面。耦合面在瞬态多场求解器中设置为可移动边界条件。为了便于计算收敛，非定常求解计算以定常结果作为初始值。

2 模拟结果及分析

2.1 声压模拟结果

不同时刻下燃烧室内压力分布云图如图8所示。图片显示了声波产生及传播的过程，传递方式可近似为平面波。

图8 不同时刻声压云图

取燃烧室中轴线上的点进行压力分析，燃烧室内各处声压幅值在燃烧室内不同位置变化情况如图9所示。

图9 燃烧室内各处声压幅值对比图

由图9可知，声压在入口处具有较大幅值，到出口端声压幅值逐步减少，燃烧室前半段声压力在模拟工况范围内，随频率升高而降低，后半段声压力随频率升高而升高。

选取入口端一点A和出口端一点B，分析不同位置的声压情况，结果如图10所示。

由图10可知，燃烧室点火时，产生的声压波是多种波动合成，声压在接近入口处幅值较大，远离出口处声压幅值显著降低。结合声压增强公式分析其原因，是由于入口处的声压振荡与火焰放热振荡相位差较小，声波得到加强，远离入口时，声压振荡与放热振荡的相位差逐渐增大，声波增强效果减弱。此外，燃烧室内声压变化还与燃烧室形状有关。

图10 燃烧室内A、B两点声压-时间曲线

声压作用在耦合面上，燃烧室产生的声压应力结果如图11所示。

图11 不同时刻声压应力云图

不同时刻声压应力分布不断改变，但燃烧室边缘一直处于应力集中状态，这些区域产生的局部应力集中现象，会导致疲劳裂纹的萌生和扩展，危及航空发动机和燃气轮机的安全运行[17]。提取一点的应力时间曲线如图12所示。将应力曲线进行傅里叶变换得到频域图如图13所示。由图13可知，声压应力主要由声压波的前5倍频合成，其中声压1倍频的影响最大。

图12 燃烧室上一点声压应力-时间曲线

图13 应力频域图

2.2 燃烧室温度模拟结果

取55 Hz工况下的不同时刻流场温度云图如图14所示。

图14 不同时刻温度场云图

云图显示了火焰放热振荡的半个周期，当火焰放热振荡与声波振荡相移较小，就会产生热声共振现象。燃烧室内温度在前半段燃烧区变化较大，在后半段分布较为平均，导致了前半段热声共振现象比较明显，声压幅值较大，后半段声压幅值降低，这与声压分布的趋势一致。温度场作用下的结构段热应力云图如图15所示。

热应力分布的规律与声压应力略有不同，但是有部分应力集中区域重合，如燃烧室的棱角处，声压应力与热应力较大，此处为易产生疲劳破坏的区域。此外热应力与声压应力量级相同，声压应力为交变应力，易引起结构共振。

图15 不同时刻热应力云图

3 结论

(1)本文用CFX软件模拟声波传递的过程，计算结果与理论计算结果基本吻合，为双向耦合过程中声波传递特性模拟提供了一个有效的方法；

(2)燃烧室内声压幅值受燃烧放热振荡影响，入口处声压由于热声耦合效应而得到增强，远离入口处的声压增强效应显著降低；

(3)温度引起的热应力与声压引起的应力量级大小相同，且声压应力为交变应力，可引起结构共振；

(4)声压产生的应力波主要由声压的1～5倍频合成，其中1倍频对应力影响最大。

参考文献(References)：

[1] 高金海,王建军,马艳红,等.基于热-流-固耦合方法火焰筒壁温三维数值模拟[J].北京航空航天大学学报,2011,37(3):300-304.

[2] NOWAK D,BELLUCCI V,CERNY J,et al.Numerical modeling of thermoacoustic oscillations in a gas turbine combustion chamber[C]// ASME Turbo Expo,2006.

[3] VAN J F K,HULS R A,KOK J B W,et al.One-dimensional acoustic modeling of thermoacoustic instabilities[C]// The Tenth International Congress on Sound and Vibration.KTH,2003.

[4] FANACA D,SATTELMAYER T,Hirsch C,et al.Comparison of the flow field of a swirl stabilized premixed burner in an annular and a single burner combustion chamber[R].ASME 2009-GT-59884.

[5] PRAKASH S,LIEUWEN T,ZINN B T,et al.Acoustic based rapid blowout mitigation in a swirl stabilized combustor[R].ASME 2005- GT-68589.

[6] POZARLIK A.Vibro-acoustical instabilities induced by combustion dynamics in gas turbine combustors[DB/OL].http://core.ac.uk/download/pdf/1478150.pdf,2017-08-20.

[7] HULS R A.Acousto-elastic interaction in combustion chambers[D].Enschede:University of Twente,2006.

[8] 陈昕巍.内燃机燃烧过程热声耦合机理的现象学研究[D].天津:天津大学,2006.

[9] 韦静思.内燃机燃烧过程中热声耦合机理的研究[D].天津:天津大学,2010.

[10]MINA S,JIM B W K,ARTUR K P.Sensitivity of the numerical prediction of flow in the limousine combustor on the chosen mesh and turbulent combustion model[C]//San Antoio:Proceedings ASME Turbo Expo 2013,2013.

[11]宋学官，蔡林，张华.ANSYS流固耦合分析与工程实例[M].北京：中国水利水电出版社，2012.

[12]钱若军,董石麟.流固耦合理论研究进展[J].间结构,2008,14(1)：4-1.

[13]HATHOUT J P.Thermoacoustic Instability[DB/OL].http://www.researchgate.net/publication/2578539_Thermoacoustic_Instability,2017-08-20.

[14]BLITZ J,MAYER W G.Temkin.Elements of acoustics[M].Elements of acousfics/.Wiley,1981.

[15]马洪安,解茂昭,曾文,等.航空发动机燃烧室燃烧过程与排放物生成的反应动力学数值模拟[J].航空动力学报,2013,28(2):297-306.

[16]薛帅杰,洪流,杨伟东.声学扰动对燃烧室声学特性的影响研究[J].推进技术,2016,37(2):201-208.

[17]袁善虎,蒋洪德,陈海燕,等.一种基于能量参数的非局部缺口疲劳寿命预测方法[J].推进技术,2017,38(3):653-658.