随着我国化石能源利用的不断增加，快速上升的氮氧化物排放量对环境造成了严重的破坏。控制和减少煤燃烧产生的NO_x，主要是控制燃料型NO_x的生成。由于烟气脱硝装置投资和运行费用都十分昂贵，国际上的主流研究方向之一是采用低NO_x燃烧技术。王鹏涛等^[1］分析了三次风比对NO_x排放的影响规律。王康等^[2］研究了不同煤粉细度及煤粉体积浓度对NO_x生成的影响，结果表明煤粉细度的减小明显降低了NO_x的体积分数。赵星海等^[3］结合空气分级及富氧燃烧技术，对锅炉炉内温度场以及 NO_x排放特性进行了数值模拟。刘建全等^[4］对不同结构的燃烧器NO_x生成特性进行了对比，分析了不同运行因素对NO_x排放特性的影响。张秀霞^[5]通过实验研究了锅炉分级燃烧中将一部分热风作为火下风（below fire air, BFA）与空气分级及燃料分级相结合的的燃烧技术，获得了使NO_x排放量最低的最佳配风量。Ribeirete等^[6］评估了空气分级对锅炉整体性能以及污染物排放的影响，通过实验得到了最佳过量空气系数。Houshfar等^[7］通过实验研究了不同燃料、不同分级燃烧的规律，获得了最低NO_x排放量时的工况。王文洁等^[8]以75 t·h⁻¹循环流化床锅炉为研究对象，采用数值模拟方法对其物料流动与燃烧时NO_x排放特性进行了研究。

炉膛内煤粉燃烧是极为复杂的过程，涉及诸如煤粉燃烧流动、传热传质等许多物理化学过程。场协同理论从流场、压力场和温度场相互配合的角度重新审视传热的物理机制，能够统一认识现有传热现象的物理本质，取得了较好的应用效果，但其在锅炉燃烧中的应用极少，因此有必要结合场协同理论分析炉膛整体流动及NO_x排放。

1 计算对象

某330 MW的四角切圆锅炉采用浓淡分离直流型燃烧器和同心反切燃烧技术组织燃烧，一次风在炉膛中心形成一个顺时针小切圆，二次风在外围形成一个逆时针大切圆。一次风喷口A、B内浓外淡，C、D、E下浓上淡；与一次风相间布置的二次风喷口共有7个。在主燃烧区上部有1层紧凑燃尽风（OFA）和4层分离燃尽风（SOFA）。基准工况下喷口无摆角。燃烧器横截面及喷口布置如图1所示。锅炉燃用煤的煤质特性如表1所示，表中：Q_{net, ar}为低位发热量；M_ar、A_ar、V_ar 、F_car分别为水分、挥发分、灰分、固定碳质量分数；C_ar、H_ar、O_ar、N_ar、S_ar分别为碳、氢、氧、氮、硫元素质量分数。

2 计算模型及条件

采用混合分数/概率密度函数模拟气相湍流燃烧，湍流模型为RNG k−ε模型，采用P1辐射模型计算辐射传热，P1辐射模型中发射率（黑度）为0.6，采用双平行反应模型模拟煤粉挥发分的析出，采用动力扩散控制燃烧模型模拟焦炭燃烧，采用随机轨道方法跟踪煤粉颗粒。压力−速度耦合的求解方法为SIMPLE方法，对流项为二阶中心差分，扩散项为二阶迎风差分。

基于多场协同理论，各协同角表达式分别为

$ \qquad \alpha = \arccos \dfrac{{\left| {U \cdot \nabla u} \right|}}{{\left| U \right|\left| {\nabla u} \right|}} =\arccos \left| {\dfrac{{u\dfrac{{\partial u}}{{\partial x}} + v\dfrac{{\partial u}}{{\partial y}} + w\dfrac{{\partial u}}{{\partial z}}}}{{\sqrt {{u^2} + {v^2} + {w^2}} \sqrt {{{\left(\dfrac{{\partial u}}{{\partial x}}\right)}^2} + {{\left(\dfrac{{\partial u}}{{\partial y}}\right)}^2} + {{\left(\dfrac{{\partial u}}{{\partial z}}\right)}^2}} }}} \right| $

(1)

$ \qquad \beta = \arccos \dfrac{{\left| {U \cdot \nabla T} \right|}}{{\left| U \right|\left| {\nabla T} \right|}} = \arccos \left| {\dfrac{{u\dfrac{{\partial T}}{{\partial x}} + v\dfrac{{\partial T}}{{\partial y}} + w\dfrac{{\partial T}}{{\partial z}}}}{{\sqrt {{u^2} + {v^2} + {w^2}} \sqrt {{{\left(\dfrac{{\partial T}}{{\partial x}}\right)}^2} + {{\left(\dfrac{{\partial T}}{{\partial y}}\right)}^2} + {{\left(\dfrac{{\partial T}}{{\partial z}}\right)}^2}} }}} \right| $

(2)

$ \qquad \gamma = \arccos \dfrac{{\left| {\nabla T \cdot \nabla u} \right|}}{{\left| {\nabla T} \right|\left| {\nabla u} \right|}} = \arccos \left| {\dfrac{{\dfrac{{\partial T}}{{\partial x}} \cdot \dfrac{{\partial u}}{{\partial x}} + \dfrac{{\partial T}}{{\partial y}} \cdot \dfrac{{\partial u}}{{\partial y}} +\dfrac{{\partial T}}{{\partial z}} \cdot \dfrac{{\partial u}}{{\partial z}}}}{{\sqrt {{{\left(\dfrac{{\partial T}}{{\partial x}}\right)}^2} + {{\left(\dfrac{{\partial T}}{{\partial y}}\right)}^2} + {{\left(\dfrac{{\partial T}}{{\partial z}}\right)}^2}} \sqrt {{{\left(\dfrac{{\partial u}}{{\partial x}}\right)}^2} + {{\left(\dfrac{{\partial u}}{{\partial y}}\right)}^2} + {{\left(\dfrac{{\partial u}}{{\partial z}}\right)}^2}} }}} \right| $

(3)

$ \qquad \varphi = \arccos \dfrac{{\left| {\nabla p \cdot \nabla u} \right|}}{{\left| {\nabla p} \right|\left| {\nabla u} \right|}} = \arccos \left| {\dfrac{{\dfrac{{\partial p}}{{\partial x}} \cdot \dfrac{{\partial u}}{{\partial x}} + \dfrac{{\partial p}}{{\partial y}} \cdot \dfrac{{\partial u}}{{\partial y}} + \dfrac{{\partial p}}{{\partial z}} \cdot \dfrac{{\partial u}}{{\partial z}}}}{{\sqrt {{{\left(\dfrac{{\partial p}}{{\partial x}}\right)}^2} + {{\left(\dfrac{{\partial p}}{{\partial y}}\right)}^2} + {{\left(\dfrac{{\partial p}}{{\partial z}}\right)}^2}} \sqrt {{{\left(\dfrac{{\partial u}}{{\partial x}}\right)}^2} + {{\left(\dfrac{{\partial u}}{{\partial y}}\right)}^2} + {{\left(\dfrac{{\partial u}}{{\partial z}}\right)}^2}} }}} \right| $

(4)

$ \qquad \theta = \arccos \dfrac{{\left| {U \cdot \nabla p} \right|}}{{\left| U \right|\left| {\nabla p} \right|}} =\arccos \left| {\dfrac{{u\dfrac{{\partial p}}{{\partial x}} + v\dfrac{{\partial p}}{{\partial y}} + w\dfrac{{\partial p}}{{\partial z}}}}{{\sqrt {{u^2} + {v^2} + {w^2}} \sqrt {{{\left(\dfrac{{\partial p}}{{\partial x}}\right)}^2} + {{\left(\dfrac{{\partial p}}{{\partial y}}\right)}^2} + {{\left(\dfrac{{\partial p}}{{\partial z}}\right)}^2}} }}} \right| $

(5)

式中：T为温度；p为压力；u、v、w分别为X、Y、Z三个方向的速度分量；α为速度矢量U与速度梯度 $ \nabla $ u之间的协同角，β为速度矢量U与温度梯度 $ \nabla $ T之间的协同角，γ为温度梯度 $ \nabla $ T与速度梯度 $ \nabla $ u之间的协同角，φ为压力梯度 $ \nabla $ p与速度梯度 $ \nabla $ u之间的协同角，θ为速度矢量U与压力梯度 $ \nabla $ p之间的协同角。

一次风从 A、B、C、D层燃烧器喷入，E层燃烧器停用。二次风喷口、紧凑燃尽风喷口、分离燃尽风喷口全开。采用均等配风方式，总风量为353.7 m³·s⁻¹，总煤量为36.8 kg·s⁻¹，其中燃烧器喷口的煤粉浓淡比为7:3。一次风率为21%，二次风率为49%，燃尽风率为30%。按照电厂实验数据给定锅炉边界条件，一次风进口温度为343 K，二次风和燃尽风进口温度为614 K，冷灰斗底部壁面平均温度为472 K，炉膛壁面平均温度为654 K，过热器区域烟道壁面平均温度为792 K，再热器区域烟道壁面平均温度为696 K，省煤器区域烟道壁面平均温度为549 K。炉膛整体为规则形状，采用结构化网格处理。为减小计算中的伪扩散，将炉膛主燃区横截面划分为辐射状，对燃烧器喷口附近区域进行加密处理，并进行了网格无关性验证，网格数目约为190万。

3 结果分析

将炉内同一位置处的烟气含氧量、烟气温度、NO_x质量浓度的模拟结果与文献[1]中的实验结果进行了比较，结果如表2所示，由表中可见，实验值和模拟值的相对误差均小于3.5%。该模拟结果可以较好地反映炉内的温度，表明本文采用的模型和方法可行。

图2为炉膛中心截面协同角分布。图2(a)为炉膛中心截面协同角α分布。在炉膛燃烧区壁面附近及中心部分区域，协同角α较小，表明速度与速度梯度的协同性较好，而介于这两部分之间的大块区域，协同角α较大。图2(b)为炉膛中心截面协同角β分布，其规律类似于图2(a)中协同角α的分布规律，协同角β大的区域占比相对较大，协同角β越小，速度与温度场协同性越好，流动换热能力增强^[9-11］。

图2(c)为炉膛中心截面协同角γ分布。在整个燃烧区域，尤其是喷口附近，协同角γ较大，且分布比较均匀。协同角γ越大，表明强化传热的综合性能越高。图2(d)为炉膛中心截面协同角φ分布，类似于γ分布，协同角φ分布也比较均匀，大角度φ主要集中在壁面附近。φ越大，表明压力梯度与速度梯度的协同性越好，流体的流动阻力越小。图2(e)为炉膛中心截面协同角θ分布，在炉膛中心大部分区域以及炉膛顶部及折焰角处，流线与压力梯度接近垂直，协同角θ较大，表明传热功耗大。但是在喷口附近，协同角θ较小。

图3为部分截面烟气温度分布。射流所在的截面温度分布能直观反映射流轨迹和火焰形状。通过分析一次风截面温度可以看出，距离燃烧器喷口较近的射流根部温度较低，这对燃烧器安全运行是有利的。炉内特别是燃烧器所在区域存在强烈的湍流，形成了高温切圆，且燃烧器所在区域为高温区，最高温度接近1 550 K。

表3为A/B层与Z轴不同竖直摆角工况，表4为DD/DE层不同摆角工况，表5为 EE/OFA层不同摆角工况。表3～5中Z轴正方向为从炉膛底部到顶部方向。分别改变A/B层、DD/DE层、EE/OFA层喷口摆角，结果发现，在改变的角度范围内，协同角β、α、γ、φ及θ整体变化不大，但是A/B层的γ值变化范围稍大些，结合图2（a)～(e)的情况来看，局部区域内协同角较优。

图4为不同工况下沿锅炉高度方向NO排放质量浓度分布。增大DD、DE层入射方向沿顺时针方向的偏转角度时，炉膛底部的NO质量浓度先增加后减少。沿锅炉高度方向上看，最高NO排放质量浓度均出现在主燃区。该区域内煤粉集中燃烧，是产生NO的主要区域。工况1的主燃区NO质量浓度最高，这是由于工况1中主燃区中心煤粉燃烧充分，产生更多的燃料型NO。而工况5、6、7与工况1（基准）相比，燃尽风层和折焰角区域的局部NO质量浓度升高。随着DD、DE层入射方向沿顺时针方向的偏转角度增大，煤粉在主燃区的停留时间缩短，使得NO质量浓度下降，更多的未燃尽煤粉在炉膛上方燃烧，从而使整体燃烧效率下降。

4 结　论

对不同条件下的炉膛燃烧流场进行了数值模拟，分析了多场协同角度、烟气温度以及NO排放质量浓度，主要结论为：在炉膛燃烧区壁面附近及中心部分区域，协同角α较小，而介于这两部分之间的大块区域，协同角α较大。协同角β大的区域占比相对较大。在整个燃烧区域，尤其是喷口附近，协同角γ较大，且分布比较均匀。协同角φ分布比较均匀，大角度φ主要集中在壁面附近。在炉膛中心大部分区域以及炉膛顶部及折焰角处，协同角θ较大，但是在喷口附近，协同角θ较小。