基于动态增厚模型的甲烷/空气预混旋流燃烧回火大涡模拟

夏 昊，王金华，卫旭涛，韩 旺，吴凤荣，黄佐华

(1.西安交通大学动力工程多相流国家重点实验室，西安 710049；2.爱丁堡大学工程学院，爱丁堡EH8 3JL，苏格兰，英国；3.中国石油天然气管道工程有限公司，廊坊131000)

回火是一种存在于燃气轮机中常见的燃烧不稳定问题.火焰向上游传播会损伤燃烧器头部，降低燃烧器寿命.主要的回火机理主要有自点火、边界层回火、燃烧诱导涡破碎导致的回火、燃烧不稳定性引起的回火等.燃气轮机中燃烧器结构以钝体旋流燃烧器为主.钝体旋流燃烧器回火机理以边界层回火和燃烧诱导涡破碎机理为主.边界层回火主要是由于在边界层内，流场速度小于在壁面熄火距离以上范围内的火焰传播速度，从而导致回火.这个模型被称为临界速度梯度模型.但是，在湍流条件下，临界速度梯度无法准确预测边界层回火.Putnam[1]在临界速度梯度模型的基础上提出了无量纲的 Peclet 数模型.Peclet 数模型考虑燃料组分以及燃烧室内压力等其他因素对回火的影响.Kurdyumov 等[2]以及Eichler 等[3]提出了Damköhler 数模型，在临界速度梯度模型的基础上，考虑了火焰速度和火焰厚度变化对回火预测的影响.Gruber 等[4]通过直接数值模拟(DNS)对低雷诺数湍流槽道流中的边界层回火进行了模拟，结果表明在湍流边界层中，火焰诱发了逆流区的产生，而逆流区又促进火焰的上游传播，当逆流区消失时，火焰的上游传播速度减小.边界层回火是流场与火焰相互作用的结果，而Peclet 数模型和Damköhler 数模型都忽略了火焰与流场的相互作用，不能够真实地描述边界层回火中发生的物理过程.燃烧诱导涡破碎是旋流火焰发生回火的主要原因.燃烧诱导涡破碎最先由Fritz 等[5]在研究旋流燃烧器预混管中的回火现象时提出.在实验中，由于火焰始终沿着中心线向上游传播，而稳态时中心流的轴向速度比火焰传播速度高得多，同时燃烧室内并未出现燃烧不稳定现象.因此，他们提出一种新的机理，即燃烧诱导涡破碎促使旋流火焰发生回火.反应区域的斜压转矩在回流泡中产生了负涡量，进而产生负的轴向速度，回流泡到达了一个更上游的位置，同时整个流场与反应达到一个稳定状态；随着当量比的增加，反应区会变得更加紧凑，火焰跟随着回流区运动.此时，由斜压转矩产生的负的周向涡量也会增加，当量比超过临界值时，就可能引起回火.但是对于钝体旋流火焰，由于中心钝体的存在，中心区域流场受钝体壁面影响，仅仅用燃烧诱导涡破碎无法完全描述回火过程.

目前，国内外已开展一些对回火的数值模拟研究.Lietz 等[6]通过大涡模拟对纯甲烷和甲烷/氢气钝体旋流火焰回火过程进行预测.虽然在宏观结构上能够较好预测，但是在回火速度等方面仍存在差异.虽然钝体旋流火焰回火机理方面已有一些研究，但对回火过程中边界层回火和燃烧诱导涡破碎机理的作用机理认识仍不够清楚.因此需要进一步实验和模拟研究钝体旋流火焰回火过程中火焰-流场-壁面相互作用机理.笔者通过大涡模拟结合增厚火焰模型对得克萨斯大学旋流燃烧器回火过程进行模拟.国内关于钝体旋流火焰回火机理的研究相对较少.Xu 等[7]采用计算和实验相结合的方法研究甲烷与富氢及纯氢合成气贫预混燃烧的回火极限，并分析了不同燃料、不同稳定方式的回火特性.数值和实验结果均表明，回火极限可以关联为Peclet 数模型.Tian 等[8]选用雷诺应力模型耦合有限速率/涡黏模型对二维氢气轴对称旋流预混火焰进行模拟，研究了预混区出口周向涡量对燃烧诱导涡破碎回火的影响.综合来看目前国内针对富氢及纯氢燃料在燃气轮机背景钝体旋流燃烧回火方面系统性实验研究和高精度模拟研究非常缺乏.对于回火的非定常性，大涡模拟由于具有较低的计算成本和较高的湍流解析度而被广泛使用.本文利用增厚火焰模型大涡模拟耦合甲烷两步机理对甲烷预混火焰在得克萨斯大学旋流燃烧器中的回火过程进行模拟.分别采用绝热边界条件及热损失边界条件对回火过程进行模拟并对比热损失对火焰结构及回火速度的影响，并对预混钝体旋流火焰回火过程进行分析.

1 数值模拟方法与计算设置

1.1 数值模拟方法

本文采用大涡模拟的方法对甲烷预混旋流火焰回火进行模拟.湍流模型采用Dynamic Smagorinsky 模型[9].燃烧模型采用动态增厚模型.增厚F 倍的火焰面可以通过求解方程(1)获得：

式中：xk为组分k 的摩尔分数；Dkj为组分j 和 k 的二元扩散系数，通过Chapman-Enskog 表达式计算.对式(1)滤波后得到式(3)：

式中：n 为火焰面增厚之后包含的网格数，本文n＝5；Δmesh为网格尺寸；δL为层流火焰厚度；f(Ω) ＝tanh(αΩ)，α控制过渡区域厚度[12].最终，组分方程重整为

1.2 计算设置

本文模拟的对象为德克萨斯大学旋流燃烧器，其主要结构如图1 所示.混合段长度为150 mm，出口直径为52 mm.混合段中心为直径25.4 mm 的圆柱形中心体，中心体底部安装轴向旋流器.燃烧室结构为圆柱形，直径为100 mm，高度为150 mm.轴向旋流器由8 个与轴向呈60°夹角的叶片构成，旋流数约为0.9[13].

图1 德克萨斯大学旋流燃烧器[13]Fig.1 UT swirl combustor[13]

计算域设置如图2 所示，计算域入口为燃烧器预混段入口，在出口上方设置延伸段以防止出口回流.入口处速度边界为预先计算的湍流流场.初始时刻预混段内为预混气，在出口上方设置高温区域.壁面为无滑移边界.本文计算的回火工况为当量比0.8，平均速度为2.5 m/s，化学反应机理采用包含6 组分、两步反应的2S_CH4_BFER[14].在计算回火工况之前，首先计算当量比为0.6 的稳定工况.之后将上游当量比改至回火当量比.

图2 计算域及边界设置Fig.2 Computational domain and boundary settings

本文采用的求解器基于开源 CFD 平台OpenFOAM 编写.压力-速度耦合通过PIMPLE 算法计算.采用隐式欧拉格式对时间推进进行求解.整个计算域网格数约为1 千万.对混合段内尤其是钝体壁面附近的网格进行加密，壁面附近网格尺度约为0.4 mm.通过网格敏感性评估，发现当前网格可以对流场统计预测.

2 结果与分析

2.1 冷态流场对比

由于回火主要发生在预混段内，本文主要对比预混段内不同高度的冷态流场.分别提取预混段出口上游5 mm、29 mm、54 mm、73 mm 处的轴向平均速度与周向平均速度与实验值进行对比，如图3 所示.图中虚线为模拟值，散点为实验值.由图3 可知，周向速度峰值出现在中心体附近，而轴向速度峰值出现在外壁面附近.除了在预混段出口上游5 mm 处靠近中心体处计算得到的平均轴向速度和周向速度略微低于实验值以及周向速度在上游73 mm 处整体略低于实验值，整体冷态流场可以较好地反映实验流场.由于下文火焰结构分析区域主要集中在远离出口处以及旋流器的中间位置，因此可以认为这些偏差对后续结果分析影响不大.

图3 冷态流场对比Fig.3 Cold flow comparison

2.2 温度边界对回火速度的影响

中心体壁面分别采用绝热以及考虑热损失两种热边界条件对回火过程进行计算.考虑热损失通过在温度边界设置壁面导热系数，外部环境换热系数及外界环境温度实现.本文中考虑热损失工况环境温度设置为350 K，环境换热系数为10 W/(m2·K)，壁面导热系数为10 W/(m·K).图4 为两种不同热边界条件下计算得到的回火过程中瞬时热释放.在绝热边界条件下，火焰在中心体壁面附近热释放较强，而在考虑热损失情况下，火焰在壁面处热释放减弱.这主要是由于额外的热损失使壁面温度低于绝热边界条件下的壁面温度，从而导致反应减弱.而在绝热边界条件下，由于边界层内速度低，未燃气驻留时间更长，且没有热损失，因此热释放率更高.此外，在考虑热损失后，火焰整体热释放率低于绝热边界条件下得到的火焰热释放率.

图4 不同热边界条件下火焰热释放分布Fig.4 Flame heat release rate under different thermal boundary conditions

本文提取温度1 400 K 等势面以表征进展变量为0.7 的火焰面，如图5 所示.在预混火焰回火过程中，火焰面出现火舌结构.由于存在热损失，火舌顶点径向位置与中心体表面有一定距离，而在绝热边界条件下火舌完全附着在中心体表面.根据一定时间内火舌顶点传播的轴向距离，可以计算得到平均回火速度.不同边界条件下计算得到的平均回火速度及实验得到的平均回火速度如表1 所示，在绝热边界条件下，平均回火速度约为 0.71 m/s，远高于实验值0.58 m/s[15].在考虑热损失后，平均回火速度为0.57 m/s，接近于实验值.考虑热损失后，回火速度准确度相比于绝热边界条件有明显改善.

图5 瞬时温度1 400 K 等势面Fig.5 Instantaneous isosurface of 1 400 K

表1 不同热边界条件下回火速度Tab.1 Flashback speed under different thermal boundary conditions

2.3 回火过程分析

基于上一部分的分析，本部分根据考虑热损失得到的计算结果对回火过程进行进一步的分析.图6为回火过程中火焰面瞬时结构的变化过程.在回火过程后期，由于燃烧器出口下游为高温已燃气，温度较高，仍然能得到1 400 K 等势面，但已没有热释放区域.由图可知，在回火过程中，火焰出现火舌结构，火舌随预混段内部流场作周向运动，并沿中心体向下传播.图7 为回火过程中某一瞬时t 时刻截面轴向速度及二维涡量分布.图中黑线表示进展变量0.7 的火焰面.图7(b)中涡量场分布为当地涡量在截面法向上的分量.在回火过程中，上游受火焰面的作用形成的涡量场促进了边界层内的回流，使火焰更容易向上游传播.而当火焰面传播到上游区域，当地流场受到火焰影响而产生回流，进一步促进了火焰向上游传播.因此，火焰不断向上游传播.

图6 回火过程Fig.6 Flashback process

图7 轴向速度场与二维涡量场Fig.7 Axial flow velocity and two-dimensional vorticity

提取每一时刻的火焰面顶点，可以得到图8 所示的火舌顶点轨迹.由图8(a)可得，火舌周向运动方向与预混段内部流场方向相同，沿中心体壁面向下传播.在火舌传播一定距离后，在某一位置会出现新的火舌，因此火舌顶点轨迹会出现断点.火舌轴向运动轨迹如图8(b)所示.在回火过程中，火舌轴向运动并不是匀速的，从出口至上游150 mm 处，火舌轴向传播速度逐渐减小.蓝线为对轨迹线性拟合得到的平均速度.当火焰传播到旋流器附近时，火舌顶点几乎不再向下传播，仅沿着钝体周围作圆周运动.这主要由于当火焰传播到旋流器附近时，截面积减小，局部流速高于火焰传播速度.钝体旋流火焰主要的回火机理为边界层回火与旋流燃烧诱导涡破碎回火.火舌作圆周运动表明此时旋流燃烧诱导涡破碎占主导.由图8(a)可得，火舌的轨迹并非完全连续.当火舌传播一定距离后，在某一位置形成新的火舌继续传播.而在火舌传播过程中，有些时刻火舌顶点圆周运动并不明显，表明此时火舌顶点相对来流的运动较为明显，受边界层回火机理影响较大.

图8 回火火舌顶点三维轨迹及轴向轨迹Fig.8 Three dimensional trajectory and axial trajectory of flame tongue

3 结论

本文通过动态增厚火焰大涡模拟结合两步机理对德克萨斯大学旋流燃烧器中甲烷预混火焰回火进行模拟，分别计算了绝热条件下及考虑热损失下的回火过程.主要结论如下：

(1) 动态增厚模型耦合两步反应机理能够对甲烷预混旋流火焰边界层回火宏观结构及回火速度进行有效预测.热损失对回火速度预测有很大影响.相比于绝热边界条件，考虑热损失后回火速度的预测更接近实验值.而采用绝热边界会使回火速度提高约22%.

(2) 预混钝体旋流火焰回火过程中会形成火舌结构.火舌顺着预混段内流场方向作周向运动，同时沿着中心体壁面向上游传播.

(3) 对于预混钝体旋流火焰，火舌结构的传播并非完全连续，而是会出现新的火舌.当边界层回火机理占主导时，火舌顶点相对来流传播变明显，周向运动会减弱.