徐 榕，程 明，赵坚行，刘 勇

（1.南京航空航天大学能源与动力学院，南京 210016；2.沈阳发动机设计研究所，沈阳 110015）

1 引言

旋流器是航空发动机燃烧室的主要部件之一，其性能好坏直接影响整个燃烧室的综合性能。斜切径向旋流器（又称旋流杯）是当代航空发动机燃烧室中经常使用的双级旋流器，因具有工作性能稳定可靠、流量及旋流特性相容性好等优点，已被CFM56、F101等涡扇发动机和T700等涡轴发动机采用。

为了提高对该类旋流器设计及优化的能力，使其在现代发动机燃烧室中更好地应用，国内外学者进行了大量研究。文献[1]利用PDPA测量旋向相反的CFM56发动机旋流杯下游液雾流场，试验研究了不同液体性质对液雾特性的影响。但因旋流器尺寸小、内部流场由多通道流域组成，结构复杂，受测量手段限制，很难通过试验研究深入了解其内部气流的流动情况。为此，有些学者采用数值方法模拟其内部流场，通过数值分析进一步掌握其工作原理。如文献 [2]采用Fluent软件对CFM56发动机旋流杯火焰筒头部流场进行计算，所得结果与PDPA的测量数据符合较好。在此基础上，提出可把数值模拟与试验研究相结合的方法发展成航空发动机燃烧室设计的工具。文献[3]采用PIV测量旋流杯环形燃烧室冷态和2相燃烧流场，同时还利用Fluent软件进行相应的数值计算。文献[4]在贴体坐标系下，研究旋流杯速度场及其湍流特性。

随着计算机技术和计算燃烧学的迅速发展，数值分析方法在了解燃烧室内部工作过程、指导燃烧室优化设计中的重要性日益增强。因此，文献[5]提出了高保真度的概念，即把发动机整体或燃烧室整体作为数值研究对象，在计算燃烧室流场时，尽可能保留燃烧室实际的复杂结构，以使模拟结果更真实地反映实际燃烧室内气流的流动情况。

本文把实际的斜切径向旋流器环形燃烧室作为研究对象，采用自编的燃烧室三维2相燃烧流场数值仿真专用程序，计算燃烧室全流程流场；在任意曲线坐标系下，采用多区域耦合法[6]，对其2相湍流燃烧全流程流场进行数值模拟。

2 三维网格生成

因研究对象结构复杂，且其内部流通区域为多连通域，故本文采用微分方程和分区相结合的方法编制网格生成程序，生成计算所用的三维贴体网格，供流场计算用。

微分方程法网格生成方程的一般形式为[7]

式中：Pm为控制源项，用来调整区域内网格的分布，直接影响网格生成质量（i,j =1,2,3 ）为逆变度量张量；ξi为任意曲线坐标；xk为圆柱坐标。

为了保证网格合理分布，满足流场计算的要求，本文采用分区法，把燃烧室分为前置扩压器、突扩段、旋流杯、帽罩、火焰筒和燃烧室内、外环冷却通道等6个部分，然后再将其组合成燃烧室整体网格，网格数为208×119×72。

图1 环形燃烧室结构

图2 斜切径向旋流器网格

由于旋流器周向均布，为了计算方便，本文选取包含1个旋流杯的环形燃烧室1/20扇形部分作为计算区域。图1、2分别为按上法生成的环形燃烧室整体网格及旋流器网格。由图可知，它由突扩扩压器、旋流杯、带帽罩的火焰筒和燃烧室内、外环冷却通道等部分组成。在火焰筒上、下壁面上，还分别开有主燃孔、掺混孔以及11排由导流板和许多小孔组成的气膜冷却槽；气膜孔出口的导流板和旋流器出口的挡溅盘因结构复杂，很难生成结构化网格，因此一般都将其忽略。但考虑到前者直接影响气膜冷却效果，后者对回流区形成以及主燃区流场有决定性作用，为了更真实地模拟实际发动机燃烧室流场，本文生成了挡溅盘和气膜孔出口导流板的网格。

3 基本控制方程及数学模型

3.1 气相基本控制方程

本文采用标准k-ε模型描述湍流特性，采用2阶矩-EBU湍流燃烧模型估算化学反应速率，采用六通量热辐射模型考虑热辐射对气流温度和壁温的影响，采用颗粒轨道模型考察液滴运动轨迹及其沿轨道变化的过程。气相采用Euler方法处理，液相采用Lagrange方法处理，气、液2相之间耦合采用PSIC法。在三维任意曲线坐标系(ξ,η,ζ)下，控制方程的通用形式为

式中：变量φ可分别表示速度u、v、w，湍流动能k及其耗散率ε，焓h，混合分数f，燃油质量分数mfu，燃油质量分数脉动均方值g，辐射通量 Rx、Rr和 Rθ；Γφ为各变量的运输系数，Sφ为气相场自身源项，为油珠蒸发产生的源项；U、V和W为任意曲线坐标系下的速度；gij为协变度量张量；J为坐标转换雅可比行列式，具体含义详见文献[7]。

3.2 液相基本方程

液相采用Lagrange方法处理。在曲线坐标(ξ,η,ζ)下，油珠运动方程可写为

应用4阶Runge-Kutta方程求解式（3），得到任意曲线坐标系(ξ,η,ζ)下的油珠运动速度，然后求得任意曲线坐标系下计算区域内油珠的运动轨迹，再由逆变换确定油珠在物理平面上的位置(x,r,θ)。

3.3 湍流燃烧模型

2阶矩-概率密度模型的基本思想是假定浓度脉动用2阶关联矩方程封闭，而对温度脉动和浓度脉动相关项采用简化Pdf模拟，并近似认为温度与浓度脉动的概率密度函数相互独立，则其时均反应速率表达式为

式中各相关项可由输运方程求得

因为式（4）虽考虑了温度与浓度脉动对化学反应速率的影响，但没有充分考虑湍流对化学反应的作用，而实际上湍流流动对燃烧过程的影响较大。为了弥补上述不足，本文采用了2阶矩-模型与EBU湍流燃烧模型相结合的形式,即2阶矩-EBU(SOM-EBU)湍流燃烧模型[8，9]，该模型取二者中的较小值来计算化学反应速率

4 结果与分析

根据上述计算方法与数学模型，对斜切径向旋流器环形燃烧室全流程流场进行了预测，分析了试验与最大等3种不同工况对燃烧室全流程2相燃烧流场的影响；其部分结果如图3～12所示。

本文为了验证数值模拟的可靠性，先在2种试验工况（case1,case2）下进行全流程2相燃烧流场计算，把所得的计算结果与试验数据进行比较。在此基础上，再在最大工况(case3)下进行相应的计算。各工况下的进口气流温度和油气比见表1。

表1 进口气流参数

图3为通过旋流器中心截面(K=42)的流线，图4为图3中帽罩内速度矢量局部放大。

图3 燃烧室K=42截面流线（case1,计算）

图4 K=42截面帽罩内速度矢量放大(case1,计算）

从图3中可见，在扩压器上、下突扩段有明显的旋涡作为气动壁面随进口状态的变化自动调整，削弱了燃烧室对进口气流畸变的敏感性，使得进入火焰筒的气流保持稳定。旋流器出口的高速旋转射流与主燃孔进来的气流相互作用，形成强逆向压力梯度，在火焰筒头部产生了上、下2个较强的旋涡，形成稳定的中心回流区，可作为点火源用；此回流区的大小与火焰的稳定性直接相关。内、外环通道冷却空气从火焰筒内、外环壁面上的气膜孔进入火焰筒，以便进行冷却。

由图4中可见，帽罩内在斜切孔旋流器的上、下存在旋涡，与突扩段处的旋涡有些雷同，起到提高头部压降、削弱对进口气流畸变的敏感的作用。由此可知，帽罩内气流流动状况对进入旋流器的气流速度分布有一定影响，预测该处气流速度分布可更真实反映旋流器斜切孔进口速度分布，以及为分析该速度分布对燃烧室内2相燃烧流场影响提供依据。

图 5（a）、（b）分别表示在试验工况1条件下，计算与试验获得的通过旋流器中心截面的火焰筒头部局部区域的热态流场流线分布。由图中可知，计算与试验所得的热态流场及其中心回流区基本相似。

在采用PIV测量火焰筒头部流场时，由于受到实际火焰筒尾部收缩以及内、外环弯曲壁面的影响，从火焰筒出口处进入的片状激光能够照亮的流动区域以及CCD相机从观察窗捕捉粒子图像的范围都受到一定限制，因此，PIV能获得的速度信息的区域相对较小，仅包含火焰筒头部的局部(图 5（b）)。另外，在旋流器出口附近，试验所得的热态气流速度分布(图 5（b）)与计算得到的相应的速度分布(图 5（a）)不完全相同，这是因为受到油雾锥的干扰，在该处CCD捕捉到的大多是油滴而非示踪粒子，因此，在图 5（b）中所显示的旋流器出口附近处的速度场为油滴速度分布，而不是气流速度分布；但随着轴向距离增加，油滴迅速蒸发，并与气流混合，形成混气进行燃烧，此时，CCD才能捕捉到示踪粒子，测得热态速度场。在图中可看到油滴和示踪粒子的影响区域间存在明显的交界面。

图5 K=42截面火焰筒头部局部放大流线

图6、7分别为在最大工况(case3)下，火焰筒内燃油轨迹和气流温度分布以及燃油质量分数分布。

图6由K=42纵截面温度分布与燃油轨迹，以及分别通过主燃孔(I=120)、掺混孔(I=155)与燃烧室出口 (I=205)3横截面温度分布组成。从图6中可知，燃油在主燃孔前蒸发成油蒸气，并与空气混合，形成可燃混气，在火焰筒头部主燃区进行燃烧，大部分燃油在该区烧完(图7)。故主燃区为火焰筒内气流温度最高的区域。

图6 火焰筒内燃油轨迹以及温度分布(case3,计算)

图7 燃油浓度分布(case3,计算)

图 8（a）、（b）分别为在最大工况与试验工况（case1）下，计算所得的旋流器中心截面(K=42)温度分布。由图中可知，2种工况下所得的温度场基本相似，但是，在最大工况下，因进口气流温度与油气比都比试验工况下的高，燃烧又较完全，故出口温度明显升高，分布也更为合理。

图8 K=42截面温度分布

图9为在最大工况下的旋流器中心截面 (K=42)CO2质量分数分布。由图中可见，因化学反应主要发生在火焰筒主燃区，大部分燃油在此烧完，故该区CO2的质量分数最大，后随着轴向距离增加，掺混孔射流以及气膜孔冷空气流进入，与高温燃气流掺混，使CO2质量分数逐渐减小。

图10为最大工况下通过掺混孔横截面(I=155)的温度分布。从图中可见掺混孔空气射流与主流混合的情况，以及在2个掺混孔之间的近壁处有个较高的温度区。这是由于有少量燃油在主燃区未完全燃烧，在主燃孔气柱扰流作用下形成的低速区域进行燃烧。

图9 K=42截面CO2质量分数分布(case3，计算)

图10 掺混孔I=155截面温度分布(case3，计算)

图 11（a）、（b）分别为在最大工况与试验工况（case1）下所得的燃烧室出口截面 (I=205)温度分布。由图可知，与试验工况下的相比，在最大工况下所得的出口截面温度较高，而且分布较均匀、合理，这是由于最大工况为设计工况，其进口条件(如温度及油气比等)更有利于燃烧。

图12为出口温度径向分布。从图中可见，在试验工况case 1和case2下所得的温度值与试验数据基本相符，可见本计算所用的数学模型与计算方法较为合理。

图11 出口截面I=205温度分布

图12 出口温度径向分布曲线

case2下的出口温度比case1下的略高，这是因其进口油气比稍大，加入燃油也略多些。在最大工况(case3）下所获得的分布曲线形状，虽与试验工况的略有区别，但符合燃烧室出口温度分布要求；此外，因其进口温度与油气比都较大，故出口温度也较高。由此可见，进口工况对燃烧室出口温度径向分布有着重要影响。

5 结论

（1）本文在任意曲线坐标系下，对包括突扩压器、帽罩、斜切径向旋流器、火焰筒以及内、外环冷却通道在内的环形燃烧室全流程流场进行了计算。计算结果表明：以数值模拟实际燃烧室形状，可提高数值仿真的保真度；计算燃烧室全流程流场能更真实地反映实际燃烧室内的气流流动、传热与燃烧情况，有助于进一步了解燃烧室内的各种复杂现象。

（2）分析了不同进口工况对全流程流场的影响，结果表明进口气流参数变化对燃烧流场的影响较大。

（3）计算所得的流场以及出口温度分布结果与试验数据符合得较好，表明本文计算方法合理，计算程序可靠，所得研究结果可为某型燃烧室优化设计提供可靠依据。

[1]Jeng S M，Flohre N M.Fluid Property Effects on Non Reacting Spray Characteristics Issued from a Counter-Rotating Swirler[R].AIAA 2005-356,2005.

[2]Giridharan M G,Mongia H C,Jeng S M.Swirl Cup Modeling-Part III:Spray Combustion in CFM-56 Single Cup Flame Tube[R].AIAA 2003-0319,2003.

[3]党新宪.双旋流环形燃烧室试验研究与数值模拟[D].南京：南京航空航天大学，2009.

[4]赵坚行，伍艳玲，周琳.先进燃烧室涡流杯流场计算[J].工程热物理学报，2000,21(1)：110-114.

[5]Turner M G,Norris A,Veres J P.High-Fidelity Three Dimensional Simulation of the GE90[R].NASA TM 2004-212981,2004.

[6]雷雨冰.数值模拟环形燃烧室整体流场 [D].南京：南京航空航天大学，2000.

[7]赵坚行.燃烧的数值模拟 [M].北京:科学出版社，2002.

[8]陈兴隆，周力行，张健.甲烷-空气同轴突扩湍流燃烧的新二阶矩模型的数值模拟[J].工程热物理学报，1999,20(6)：759-763.

[9]蔡文祥.环形燃烧室两相燃烧流场与燃烧性能数值研究[D].南京：南京航空航天大学，2007.