旋进涡核与火焰面耦合作用对燃烧稳定性影响的数值研究
2018-07-25王振林李祥晟庄士超杨诏
王振林, 李祥晟, 庄士超, 杨诏
(西安交通大学能源与动力工程学院, 710049, 西安)
现代低排放燃气轮机燃烧室中广泛采用贫预混燃烧方式,其所固有的燃烧不稳定性对装置的安全运行危害极大,对不稳定燃烧的产生机理进行研究进而获得有效的控制方法具有十分重要的意义。以往的文献研究表明,燃烧室中燃烧不稳定来自于当量比的扰动和气动力学引起的涡旋形成及破碎过程,其中当量比扰动对燃烧不稳定的影响已经开展了大量的实验及数值研究工作。在一个螺旋形的旋进涡核(PVC)引起的旋涡与燃烧不稳定间的关系方面,目前所开展的研究较少,尤其是数值研究与问题的复杂性密切相关。要捕获涡旋与不稳定燃烧过程,需要能够对微小的流场结构进行描述的精细模型,计算量巨大,尤其在湍流耦合燃烧时更是如此。
内回流区的形成通常会有PVC伴随着发生,其特征是一个离轴的周期性绕中轴旋转的旋进运动[1]。大量研究表明,PVC会通过多种方式影响旋流火焰的稳定。例如,PVC会加强燃料/空气的混合[2],强化未燃和已燃气体的混合[3],卷曲、拉伸以及局部熄灭化学反应区[4],PVC也有可能与火焰的热声耦合振荡相互作用[5]。在某些燃气轮机燃烧室中,冷态和热态情况下均存在PVC[1-5],而在另一些燃烧室中,冷态时存在PVC,热态情况下PVC却被抑制消失。Stohr等使用高速激光诊断法,实验研究了部分预混旋流抬举火焰的动力学特性[6],发现火焰根部的熄灭与火焰和PVC的相互作用紧密相关[7]。Oberleithner等通过实验和线性稳定性分析的对比研究,认为PVC是一个总的流体动力学不稳定模态的表现,PVC通常发生在冷态流动和M型抬举火焰中,但不存在于V型的依附火焰中,同时,喷嘴出口上方的密度场在PVC形成的过程中起重要作用[8]。Manoharn等的另一个研究[9]也证实了PVC会在由火焰引起的密度梯度较大的工况下受到抑制。Gorbunova等则采用热源主动加热的方法研究了热源在不同加热功率下PVC的频率和幅度的变化[10]。
国内学者对PVC流场结构也展开了相应的实验和数值研究。樊艳娜等通过实验方法对不同实验条件下燃烧室中的冷态流场进行了详细分析,并证实了PVC结构的存在[11]。张宏达等对悉尼旋流燃烧器的冷态流场进行了大涡模拟,研究了不同旋流数下的流场结构、旋进频率和旋进涡核[12],并采用本征正交分解的方法重构了有PVC存在的冷态湍流脉动速度场[13]。在PVC与火焰耦合作用的研究方面,目前尚未有相关公开文献报导。
本文以文献[14]采用的贫油预混燃烧室为研究对象,利用商业CFD计算软件ANSYS FLUENT结合UDF,对该燃烧室分别在热功率为10、35 kW时的燃烧过程进行大涡模拟,重点研究了燃烧工况下热声耦合不稳定燃烧发生时瞬态流场结构PVC与火焰及燃烧放热之间的关系。通过对数值模拟获得的瞬态冷热态流场和热态火焰特性的分析,确定了引起35 kW工况下不稳定燃烧的重要因素,为深入理解湍流燃烧中流场结构与燃烧放热的相互作用,从而为制定有效的不稳定燃烧控制方法提供了参照。
1 计算模型与数值模拟方法
1.1 计算模型
本文采用的模型燃烧室为Turbomeca设计的一个工业燃烧室,图1给出了计算流体动力学(CFD)计算区域的三维建模图,除了充气室、旋流器和燃烧室外,为了模拟无反射出口,计算域还包含了一部分大气环境。图2a给出了旋流器和燃烧室部分的二维结构示意图。大气温度及压力下的空气通过充气室进入一个由12个径向分布的旋流叶片和一个收缩钝体组成的旋流燃烧器中,与通过旋流叶片上1 mm直径的小孔进入旋流器中的燃料进行混合,以确保混合气体在进入到燃烧室前得到充分预混。燃烧室轴向截面为85 mm×85 mm,且每一面都装备了1.5 mm厚的石英来对整个火焰区域进行光学测量。更多关于模型燃烧室的细节可通过文献[14]获得。
图1 CFD计算域的三维建模图
1.2 网格划分及数值方法
计算域网格采用了四面体非结构网格和六面体结构网格相结合的混合网格,其中除燃烧室部分采用了结构化六面体网格以便于调整燃烧区域的网格数外,其他计算区域均采用了非结构化四面体网格。通过调整燃烧室核心区域内的结构化网格,对大涡模拟(LES)的网格无关性进行了验证,在该区域内采用71万和159万两种不同数量的网格,如表1所示。结果显示,两种计算网格之间偏差不超过5%,因此本研究在燃烧核心区采用了71万的网格,计算的总网格数约为347万。网格在轴向、径向、周向上均为非均匀网格,图2b给出了旋流器和燃烧室部分在中心截面上的网格划分信息。
(a)旋流器和燃烧室结构二维示意图
(b)旋流器和燃烧室网格划分示意图图2 旋流器和燃烧室结构及网格划分示意图
燃烧反应采用文献[14]给出的甲烷与空气的两步反应机理2S_CH4_BFER,该机理包含了6种组分(CH4,O2,N2,CO,CO2和H2O)。Navier-Stokes方程采用可压缩的LES方法进行求解。亚网格模型采用Wall-modeled LES模型,该模型仅在对数层内的区域激活RANS,而在边界层外采用大涡模拟。湍流-化学反应间的相互作用采用限速率/涡耗散模型来描述。在数值模拟中,空气和甲烷入口均采用质量入口边界条件;协流入口和大气环境均采用压力远场边界条件;其他壁面均采用绝热无滑移边界条件。空气和甲烷入口在不同热功率下的质量流量见表2。非定常计算中时间步长的确定取决于燃烧振荡的频率。由于本文研究的不稳定燃烧振荡的频率数量级为102Hz,对应振荡的一个周期约为1 ms,因此时间步长取1~10 μs。本文时间步长取5×10-6s,数值模拟中所有变量在时间及空间上均为二阶精度。热态结果通过对冷态结果进行点燃后计算获得。
表1 燃烧室中的网格信息
表2 不同热功率下的空气和甲烷质量流量
2 计算结果分析与比较
2.1 冷态结果
图3给出了两种工况下的瞬态速度场并通过等压面的方式给出了在冷态流场中存在的大尺度PVC流场结构。如图所示,两种工况的冷态流场中均出现了PVC结构,其中10、35 kW工况下分别使用35、250 Pa的等压面显示了燃烧室中存在的三维PVC结构。图4给出了10 kW工况下燃烧室轴向截面上的涡旋分布。由此可见,PVC在燃烧室轴向截面上产生的旋涡是“之”字型反向旋转分布的。
(a)10 kW (b)35 kW图3 两种工况下冷态流场中以等压面显示的PVC结构
图4 10 kW工况下燃烧室轴向截面上的涡旋分布
2.2 25 kW热态结果与实验数据的比较
本文将25 kW的热态结果与实验数据[14]进行了对比,该工况下的燃烧情况见文献[15]。图5给出了大涡模拟获得的时均速度和温度与实验数值的对比。总体来看,数值模拟结果与实验结果吻合较好,但在喷嘴附近的速度上数值模拟结果与实验结果存在差异。其原因在于数值模拟难于精确捕获喷嘴附近存在的强旋流流动,并且缺少旋流器中的钝体结构准确的型线数据,导致数值模拟的三维建模和实际形体间存在差异。此外,从图5中可以看到,温度的数值模拟结果与实验结果在燃烧室的上游部分亦存在差异,这是由于25 kW工况下存在的两种不同燃烧状态导致的,数值模拟获得的时均及RMS温度包含了PVC存在时的数据,但数值模拟结果与实验结果总体趋势吻合良好。
2.3 10 kW与35 kW热态时间平均流场分析
(a)时均轴向速度
(b)时均径向速度
(c)时均温度
(d)RMS温度图5 数值模拟与实验结果对比
图6给出了两种工况下的时均轴向速度场。从图中可以看到,流场是典型的收缩旋流火焰,包括了一个从喷嘴进入到燃烧室中的圆锥形射流、一个内回流区(IRZ)和一个外回流区(ORZ)3个主要的流场结构。速度梯度值高的地方发生在位于入流和IRZ之间的内剪切层(ISL)以及位于入流和ORZ的外剪切层(OSL)。图7给出了两种工况下燃烧室轴向截面中心线上的时均轴向速度曲线,二者最大的不同是存在于燃烧室上游区域的IRZ以及燃烧室下游区域速度的增长率。在10 kW的工况下,内回流区中的回流速率很小且沿着轴向逐渐增大;在35 kW的工况下,内回流区的回流速率很高,且轴向回流速度呈现出先增大后迅速减小的趋势。回流速度在轴向高度20 mm处达到最大值,在60 mm高度处降低为0,随后轴向速度呈正值。通过对比可发现,35 kW工况下IRZ中的回流现象比10 kW工况下更加剧烈。
(a)10 kW (b)35 kW图6 两种工况下时间平均速度流场
图7 轴向中心线上的速度变化
图8给出了两种工况下时均温度分布云图。图9给出了两种工况下燃烧室轴向截面中心上的时均温度曲线。两种工况下的时均温度场差异极大,与时均轴向速度场类似,巨大的差异依然存在于内回流区中。首先,两种工况代表两种不同类型的火焰:10 kW工况下,与喷嘴接触的是V型火焰,在内流区中从喷嘴出口到燃烧室下游均充满高温燃气;35 kW工况下,火焰未接触到喷嘴,该种类型的火焰被称为M型抬举火焰,内回流区中喷嘴出口附近依然是低温未燃的空气甲烷混合气体。对V型接触火焰观察可见,低温主要存在于射流中,内回流区中的高温已燃气体充斥着从喷嘴出口到燃烧室下游的整个区域。对M型火焰观察可见,燃烧室上游的回流区在轴向高度20 mm以下区域均为低温未燃的新鲜混合气体,而在20 mm以上区域的温度则在中心线两侧向下游逐渐升高,具体原因将在下文瞬态流场分析中给出。由图8可见,两种不同热功率下的时均温度场展示出了两种截然不同的火焰燃烧形态。
(a)10 kW (b)35 kW图8 两种工况下时间平均的温度分布
图9 两种工况下轴向中心线上的时均温度
(a)10 kW (b)35 kW图10 两种工况下时均放热分布
图10给出了两种工况下时均放热量q分布的云图。图11分别给出了两种工况下燃烧室轴向截面中心上的时均放热曲线。10 kW的放热区域主要集中在射流两侧与内回流区及外回流区之间的ISL及OSL中,而35 kW的放热区域范围则要大得多,主要存在于轴向高度40~60 mm处燃烧室径向两侧的位置。特别需要注意的是,在IRZ的轴向高度20 mm处也有放热较强的一小块区域,该局部区域发生较强放热的原因将在下文的瞬态场分析中给出。由图11可见,M型火焰在内回流区轴向高度20 mm处放热量达到最大。综合图7、图9和图11可以看到,在35 kW工况中,内回流区中轴向高度20 mm处燃烧室轴向中心线上的回流速度最大,是低温未燃气体与高温气体的临界点以及局部范围内放热最剧烈的位置。该位置处放热最剧烈意味着化学反应放热对气体加速作用最大,导致此处的轴向回流速度最高。
图11 两种工况下轴向中心线上的时均放热
2.4 10 kW与35 kW热态瞬态流场分析
上文分析表明,两种工况在流场及放热等场量的时均分布存在着巨大差异,本部分对瞬态流场做进一步的深入研究,详细分析瞬时温度场、放热场以及流场结构与放热之间的相互作用。
(a)10 kW (b)35 kW图12 两种工况下的瞬态流场
2.4.1 流场结构及其特性 图12给出了两种工况下的瞬态流场图。由图可见,两种工况的瞬态流场在结构上存在巨大差异。在10 kW工况下,燃烧室瞬态流场中依然存在一个稳定的射流区域,在射流的两侧分布着对称的涡结构。由于PVC在流场中存在的典型特征会在燃烧室轴向界面上产生反向旋转交错布置的涡旋,因此在该工况燃烧发生时燃烧室中在冷态情况下存在的PVC消失了,其中射流外侧与燃烧室壁面之间的旋涡由外回流区产生,而内剪切层中存在的旋涡则是由向燃烧室下游喷射的射流与内回流区中反方向流动的回流摩擦产生的。在35 kW工况下,热态流场同样存在着射流与燃烧室壁面间的外回流区,但是在射流与内回流区间的ISL中旋涡的布置与10 kW工况下有很明显的差异,此时的ISL中旋涡分布不再是中心线对称分布,而是呈现了“之”字型分布。这种分布的旋涡就是流场中PVC结构存在的典型特征,证明了在35 kW工况下PVC仍然是存在的。
2.4.2 火焰燃烧状态 图13给出了两种工况下燃烧室中压力监测点上的压力振荡信号。在10 kW工况下,计算初期压力波动比较明显,但随着计算的进行压力波动幅度急剧衰减。计算初期压力波动来自于定常计算中的数值波动,而计算后期仍然继续存在的小振幅压力波动,则是正常稳定燃烧时由于湍流脉动所造成的燃烧室压力脉动的起伏,不认为其发生了不稳定燃烧。在35 kW工况下,随数值计算的进行,燃烧室中的压力振幅逐渐增大并最终发展成了定型振荡,说明在该工况下燃烧室中发生了不稳定燃烧现象。两种工况下所发生的不同燃烧状态将在下文中从瞬态流场结构和放热之间的相互作用这一角度进行分析。
(a)10 kW (b)35 kW图13 两种工况下燃烧室内压力振荡曲线
2.4.3 流场结构与放热的相互作用 图14展示了35 kW工况下的瞬态流场与温度场云图。在驻线位置(图中虚线)处,来自内回流区中的高温燃气与由喷嘴进入燃烧室中的低温混合气正面对冲,并为低温混合气提供高温和化学基,从而促进了后者的点燃过程。驻点(图中黑点)和驻线位置的改变与PVC引起的旋涡的运动密切相关。首先,在t=0~0.5 ms期间,旋涡向左拉伸火焰并引起火焰面的卷曲,会导致火焰面的增大并加强低温未燃混合气体与高温燃气的混合。该局部区域的放热在旋涡增大火焰面和加强混合的双重作用下得到增强,而在这一过程中驻线的位置也随之转移到了喷嘴出口左侧。当t=0.75 ms时,燃烧室轴向截面在喷嘴出口的右侧出现了一个新的旋涡,该旋涡与火焰和放热的相互作用方式与t=0~0.5 ms时二者之间的相互作用方式相同。同时,驻线的位置随着右侧新旋涡的出现而转移到右侧。综上分析,发现由PVC引起的旋涡会影响火焰面的卷曲,并通过周期性地改变驻点和驻线的位置,引起局部范围内放热位置和强度周期性地变化。
(a)0 ms (b)0.25 ms (c)0.5 ms
(d)0.75 ms (e)1 ms (f)1.25 ms图14 5 kW工况下瞬态流场与温度场云图
对于M型火焰,在PVC与火焰放热相互作用的周期性运动中,虽然驻线和驻点的位置也随之周期性变化,但回流的高温燃气与来流的低温混合气的对冲主要发生在燃烧室轴向20 mm左右的高度处。这会导致回流速度在20 mm处遇到与之流向相反的低温来流,从而回流速度达到最大值,对应于图7b中燃烧室轴向中心线上的时均轴向速度在20 mm处达到最大值。驻线位置(高温燃气与低温来流对冲位置)也决定了M型火焰在内回流区中的温度场以20 mm的高度为高低温分界,对应于图9b中燃烧室轴向中心线上在高度20 mm左右处温度的急剧攀升。由于在高温燃气与低温来流发生对冲的驻线位置处发生了局部的化学反应放热,因此图10b中燃烧室轴向高度20 mm处是除了燃烧室中部以外的另一个放热较强的位置。
图15给出了两种工况下三维时均放热等值面,其中10 kW工况下的放热功率为0.04 W,35 kW工况下的放热功率为0.09 W。由图可见:在V型火焰中,放热区域位于射流的内外两侧的剪切层中,呈现出圆锥状;在M型火焰中,放热区域主要位于燃烧室轴向高度的中部位置,该位置对应于PVC的尾部区域。对图14分析可知,PVC通过自身的周期性运动对火焰面的拉伸以及对驻点驻线的影响主要存在于燃烧室轴向高度约20 mm处。可见,燃烧室中总放热波动是由PVC在喷嘴附近20 mm处的小范围放热而引起燃烧室中部大范围剧烈放热造成的。该过程为PVC在喷嘴附近引起火焰面拉伸及驻点驻线的位置发生了周期性变化,使未燃和已燃气体混合并随后在内剪切层中被引燃,从而发生反应放热。部分反应区域周期性地被从驻点附近拖拽到内剪切层中,作为内剪切层中化学反应的点火源,该化学反应随着PVC的运动在燃烧室中部位置,即在PVC耗散殆尽的位置最剧烈,从而导致燃烧室中总放热量的周期变化。这些过程与周期性的涡运动直接相关,表明PVC在该不稳定燃烧的机制中起重要作用。图16给出了燃烧中监测点上的压力及总放热量波动曲线。由图可见,在燃烧室中压力达到最大值时放热也最剧烈,且二者的波动频率一致。这说明在该情况下,燃烧室中发生了热声耦合现象且符合瑞利准则。因此,在M型火焰中,由PVC导致的燃烧室中总放热量的波动与压力波动的热声耦合是该情况下燃烧不稳定发生的一个确定的重要因素。
(a)10 kW (b)35 kW图15 两种工况下三维时均放热等值面
3 结 论
本文通过大涡模拟对10、35 kW两种工况下燃烧的冷热态流场进行了数值模拟,主要结论如下:
(1)两种热功率下的冷态流场中均存在PVC流场结构,该结构产生的旋涡在燃烧室的轴向截面上呈“之”字型分布。
(2)两种热功率下发生了不同的燃烧状态。10 kW工况下燃烧稳定进行且PVC结构消失,在内回流区中高温燃气与旋流器喷嘴相接触,火焰呈现V型;35 kW工况下发生了不稳定燃烧现象且PVC结构依然存在,火焰抬举于旋流器喷嘴上方,呈现M型。内回流区中PVC结构的存在与否是导致不同热功率下发生不同燃烧状态的原因。
(3)35 kW工况下发生的不稳定燃烧现象与PVC密切相关,PVC通过周期性地拉伸火焰面,导致旋流器喷嘴上方的驻点驻线发生周期性变化,使得燃烧室内发生了热声耦合现象。
(4)35 kW工况下燃烧室内总放热量变化的原因是,由PVC引起的局部周期性放热为点火源,在内剪切层中发生了化学反应,最终在燃烧室中游位置的PVC尾部发生了大范围且更剧烈的放热。