单个焦炭颗粒在燃烧环境中的曳力研究

屈思博，张 纬，由长福,2

(1.清华大学能源与动力工程系热科学与动力工程教育部重点实验室，北京 100084；2.清华大学山西清洁能源研究院，太原030000)

气-固两相流广泛地存在于颗粒输送、煤粉燃烧等领域，准确地构建颗粒在不同环境中的动力学模型对于指导工业设计具有重要意义.燃烧环境中的颗粒不仅与周围气体发生动量交换，还发生相间的能量交换与气体组分的变化，引起多物理场的变化.因此，气-固两相间的相互作用机制更加复杂.对于离散颗粒系统，研究者通过实验或者数值模拟等方法构建了单个颗粒的曳力模型，并根据不同的颗粒形态和工况提出了相应的修正，使得相关模型得到了广泛的应用[1-4].但是，多数研究主要集中于冷态环境，而对燃烧环境的研究较少.

近年来，研究者较多地关注于颗粒表面在特定情况下存在的斯蒂芬流对颗粒受力及其表面附近的流场结构的影响.Kurose 等[5]与Chen 等[6]采用数值模拟方法，通过在颗粒表面设置速度源分别研究了斯蒂芬流对单个颗粒以及颗粒群的曳力的影响，但是，他们均将颗粒置于常温环境，没有考虑在燃烧环境下气体组分和物性的变化对颗粒曳力的影响.Zhang 等[7]同样采用数值模拟方法，通过设置不同的化学反应速率研究了燃烧强度对单个焦炭颗粒曳力的影响，所考虑的颗粒雷诺数 pRe 均较小(小于20)，颗粒表面附近的流场结构没有发生较大的变化，限制了修正模型的适用范围.对于煤粉颗粒燃烧，斯蒂芬流的影响主要发生在燃烧前期，即大量水分和挥发分的释放在颗粒表面产生质量流，而在焦炭燃烧阶段，质量流是较小的，对颗粒曳力的影响有限.因此，有必要对在包含多种流场结构的较宽的颗粒雷诺数范围下研究燃烧颗粒曳力的变化.

焦炭燃烧是煤粉燃烧的主要阶段，提供80%以上热量，并且一般占燃烧总时间的90%，决定了煤粉燃烧的主要特性[8].目前，效率较高的固定网格数值模拟方法包括浸入式边界法[9]和虚拟区域方法等[10]，后者由于可以较好地实现颗粒刚体约束，因而对大密度比的气-固两相流的模拟具有较好的适用性[11].因此，本研究以燃烧的焦炭颗粒为对象，虚拟区域法为数值模拟工具，对半径为100～500 mm 的固定焦炭颗粒在303 K(非燃烧)与1 173 K(燃烧)温度下的受力及其流场进行定量研究.考虑到煤粉颗粒在实际的工业应用中以悬浮燃烧的形式为主，本研究设置了线性增大的入口气体速度边界条件，主要关注颗粒所受曳力与重力相等时的动力学特性，并且通过与颗粒在常温下的各物理量作对比提出修正公式，为未来研究处于燃烧状态下颗粒群的动力学特性以及指导设计高效清洁的燃烧技术提供基础.

1 数值方法

气-固两相反应流中存在复杂的相间动量交换、传质和传热，由此形成多物理场.本研究基于动量守恒定律，采用经典的纳维-斯托克斯方程来描述流体运动，以获得流场的速度与压力分布；基于质量守恒和能量守恒定律，采用流体力学经典的传质与传热方程来获得流场的组分与温度分布.对于颗粒能量守恒，颗粒的温度变化同样可由传热基本理论来描述.虚拟区域法采用结构化网格统一描述颗粒与气体.如图1 所示，颗粒(蓝色圆点)与气体(红色三角点)区域分别适用固体(刚体)以及流体控制方程及其理论解，颗粒与气体的相界面(黑色圆点)通过其所在网格的控制方程内设置相应的源相以描述气-固两相的相互作用.根据虚拟区域法思想，本研究将流体的控制方程形式应用于全计算域，并在方程中设置“虚拟源项”，以使颗粒区域的方程的解满足颗粒约束，即：在数值上保持其由颗粒控制方程求得的结果相一致.全计算域的流体各控制方程中“虚拟源项”的具体形式、方程离散与数值求解算法，以及虚拟区域法在动力学与化学反应方面的验证详见参考文献[11-12].

图1 虚拟区域法计算域Fig.1 Computational domains of the fictitious domain method

值得说明的是，对于气体连续相介质(空气)，本研究采用理想气体假设，由公式(1)来描述气体的密度、压强与温度关系.

式中：ρg是流体密度，kg/m3；p 是压强，N/m2；R是理想气体常数，其值为8.314×103J/(K·mol)；T 是温度，K；Mg是流体摩尔质量，kg/mol.

对于颗粒受力，可由牛顿经典力学定律来描述，如公式(2)所示.

式中：Fg是颗粒所受流体的曳力(Fd)，N；n是颗粒中心指向表面作用点的单位法向量；τ是流体作用在颗粒表面应力张量，N/m2；I是压强在颗粒表面作用点指向中心的单位向量；ds 是颗粒表面作用点的微元面，m2.该式表明，本研究所计算的颗粒曳力由两部分组成，即颗粒表面的压强积分：与黏性应力积分：.

本研究采用基本的化学反应动力学来描述物质变化.假设化学反应n 的反应常数kn可由阿伦尼乌斯定律描述，如公式(3)所示：

式中：An是反应n 的指前因子，其单位根据具体的反应略有差别；En是反应n 的活化能，J/mol.

受具体的煤种特性、热解与燃烧环境的影响，焦炭在实际燃烧过程中发生的化学反应较为复杂.为简化问题，本研究采用表1[11-12]中所列适用范围较广的3 类单步不可逆反应来描述焦炭的燃烧.表中：是反应n 的反应热，J/(mol·s).反应(Ⅰ)和(Ⅱ)用于描述发生在焦炭颗粒表面的气-固非均相反应，反应(Ⅲ)用于描述发生在流场空间内的气体均相反应.

表1 化学反应动力学参数Tab.1 Kinetic parameters for chemical reactions

2 模拟设置与流场分析

2.1 颗粒位置与边界条件

如图2 所示，本研究设置了一个矩形的计算域，按照颗粒半径 rp，将计算域的长与宽分别设置为40 rp和 20 rp，将颗粒设置于入口中线，距离入口10 rp处.入口设置为线性增大(加速度设置为10 m/s2)的变速度边界条件，出入口两侧设置为周期性边界条件.气体 O2和N2的质量分数之比设置为0.23∶0.77，出口压力 p0设置为 105Pa.颗粒密度设置为 103kg/m3.以颗粒迎风点为起始点(θ=0°)，按照顺时针方向在颗粒表面均匀地设置11 个采样点，以研究各物理量在颗粒表面的分布情况.本研究对于处于燃烧状态的气体的物性参数参考自NIST 的数据[13].

图2 模拟设置Fig.2 Simulation settings

2.2 流场分析

本研究以100 μm 为间隔，设计了半径为100～500 μm 的固定焦炭颗粒分别在303 K 与1 173 K 环境温度下的受力与化学反应的计算.图3～图5 分别展示了rp分别为100 μm、300 μm 和500 μm 时在某一给定气体来流速度(均大于颗粒的终端速度)且流场充分发展时，焦炭颗粒分别在常温环境下和燃烧环境下所形成的速度场和压力场以及在燃烧环境下所形成的温度场.可知，随着颗粒尺寸和气体来流速度增大，即：Re 增大，流场从稳定的层流逐渐向不稳定转变，迎风与背风区的压差也增大.特别是在常温环境下，在颗粒尺寸和气体来流速度较大的情况下，颗粒两侧出现周期性的旋转方向相反的涡脱落，形成了卡门涡街.但是，在燃烧环境下，颗粒所形成的流场并没有因为颗粒尺寸和气体来流速度增大而发生明显的转变，一直保持了比较稳定的层流形态.这是因为气体在高温下拥有不同的物性参数所导致的，尤其是气体的黏度随温度升高而增大，即：在相同速度条件下，气体的动量减小，内摩擦增大.尽管保持了相同的颗粒尺寸和来流速度，但是燃烧状态的颗粒雷诺数更小.因此，颗粒在常温与在高温环境下形成的流场形态不同，并且在 ReT较高时，两者的差别更加明显.同时，焦炭与空气发生化学反应会导致的气体组分的改变，影响气体物性，进而影响流场形态以及颗粒的受力.图5 展示了颗粒在燃烧环境下所形成的温度场.可知，随着颗粒尺寸和气 体来流速度增大，可发生气-固反应的表面面积以及气流的传质与传热速率增大，颗粒中心的温度增大，且相对于整个观察区的颗粒燃烧区更集中.

图3 焦炭颗粒所生成的速度场与流线分布Fig.3 Velocity field and streamlines from flow past char particle

图4 焦炭颗粒所生成的压力场分布Fig.4 Pressure field from flow past char particle

图5 燃烧环境下焦炭颗粒所生成的温度场分布Fig.5 Temperature field from flow past char particle under combustion conditions

3 结果讨论

3.1 颗粒表面应力分布

本研究通过分别计算颗粒表面的黏性应力和压强的积分，并将两者求和来计算颗粒曳力.在气体来流速度和环境温度的综合作用以及颗粒尺寸的影响下，颗粒的黏性力与压力对颗粒曳力贡献程度是不同的.图6 给出了rp=100μm 和ug为0.5 m/s、2.5 m/s时，颗粒表面黏性应力在来流方向的分量和压强的分布情况.可知，整体上，在相同的气体来流速度条件下，相比于常温环境，颗粒在燃烧时所受黏性应力更大.当气体来流速度较小时，温度环境对压强的影响并不明显；而当气体来流速度较大时，相比于常温环境，颗粒在燃烧时所受压强更小.一方面，这符合在燃烧环境下气体黏性较大，从而导致颗粒所受黏性力较大的推测；另一方面，也体现了在燃烧环境下流场更加稳定，颗粒所受压力较小的趋势.表2 给出了该条件下颗粒所受黏性力与压力的积分，可知，压力对颗粒曳力的贡献较高，特别是在颗粒雷诺数较大的情况下.

表2 rp＝100 μm 时颗粒所受黏性力与压力的积分Tab.2 Integral of viscous force and pressure on particle(rp＝100 μm)

图6 rp＝100 μm 时颗粒表面应力分布Fig.6 Stress distribution on particle surface(rp＝100 μm)

图7 给出了ug=2 m/s 和rp为300 μm、500 μm时，颗粒表面黏性应力在来流方向的分量和压强的分布情况.可知，在相同的温度环境下，颗粒尺寸越大，颗粒所受黏性应力越小，但是颗粒尺寸对压强的影响并不明显；而在相同的颗粒尺寸条件下，不同的温度环境对压强的影响更大，相比于常温环境，颗粒在燃烧环境下所受压强更小.表3 给出了该条件下颗粒所受黏性力与压力的积分(ug= 2 m/s)，可知，整体上在相同颗粒尺寸和气体来流速度的条件下，相比于常温环境，颗粒保持了在燃烧环境下所受黏性力更大，所受压力更小的趋势.

图7 ug＝2 m/s 时颗粒表面应力分布Fig.7 Stress distribution on particle surface(ug＝2 m/s)

表3 ug＝2 m/s 时颗粒所受黏性力与压力的积分Tab.3 Integral of viscous force and pressure on particle(ug＝2 m/s)

3.2 颗粒终端速度变化

本研究设置了线性增大的入口气体速度边界条件，当颗粒所受曳力与重力相等时颗粒处于悬浮状态，并把相应的来流气体速度的作为颗粒的终端速度.图8(a)～(e)分别给出了半径为100～500 μm 的焦炭颗粒在变入口气体速度条件下的曳力变化情况及其与重力的相对大小关系.可知，当气体来流速度和颗粒尺寸较小，即颗粒雷诺数较小时，相比于常温环境，颗粒在燃烧时的曳力更大，终端速度更小.但是，在整体上，当颗粒雷诺数较大时，相比于常温环境，颗粒在燃烧时的曳力更小，终端速度更大，且随着颗粒雷诺数增大，两者差别明显.

图8 焦炭颗粒所受曳力变化Fig.8 Variations of drag force on char particle

根据前文的分析可知，当颗粒雷诺数较小时，颗粒在常温和燃烧环境下所形成的流场形态均比较稳定，颗粒所受压力的差别并不明显，但是高温以及化学反应使得气体黏性增大，导致颗粒所受的黏性力增大.因此，相比于常温环境，颗粒此时在燃烧环境下的曳力更大；而当颗粒雷诺数较大时，由于燃烧环境中的气体具有更大的黏性，从而相比于常温环境，颗粒在燃烧环境下所形成的流场更加稳定，所受压力更小，且两者压力的差别更加明显.同时，相比于黏性力，压力对颗粒曳力贡献较高.因此，相比于常温环境，颗粒此时在燃烧环境下的曳力更小.

3.3 颗粒曳力特性修正

本研究分别统计了半径为100～500 μm 的焦炭颗粒在非燃烧与燃烧状态下的终端速度，即：ut,0与ut，并将两者的比值作为终端速度的经验修正值，颗粒终端速度修正随颗粒直径的变化趋势如图9(a)所示.根据其变化特性，本研究将拟合曲线的函数形式设为 y=a · e[-b/(x+c)]+d，通过 ut/ut,0与 dp若干组经验值，本研究得到曲线的4 个参数a、b、c、d 分别为2.01、5.01×10-4、1.14×10-5、0.61，燃烧焦炭颗粒终端速度的修正公式如(4)所示.拟合的决定系数R2=0.999 5.

类似地，本研究通过公式(5)分别了计算颗粒处于悬浮状态时在常温与燃烧环境下曳力系数，即Cd与Cd0.

式中：A 为颗粒的迎风面积，m2；由于本研究在模拟中将颗粒固定，因此，在数值上，颗粒相对气体的滑移速度 (us-ug)等效于气体来流速度(取绝对值)，在 本研究中，具体为 ut.本研究将颗粒在常温与燃烧环境下的曳力系数的比值作为经验修正值，颗粒的曳力系数修正随颗粒的终端速度修正的变化趋势如图9(b)所示.根据其变化特性，本研究将拟合曲线的函数形式设为 y=ae(-bx)+c，通过 Cd/Cd0与 ut/ut,0若干组经验值，本研究得到曲线的3 个参数a、b、c 分别为39.71、2.51、0.79，即焦炭颗粒处于悬浮燃烧状态时曳力系数的修正公式如(6)所示.拟合的决定系数R2=0.999 8.

图9 修正曲线Fig.9 Correction curves

4 结论与展望

为明确单个焦炭颗粒在非燃烧与燃烧状态下的曳力变化，本研究采用虚拟区域方法对半径为100～500 μm 的固定焦炭颗粒分别在303 K(非燃烧)与1 173 K(燃烧)温度下的受力及其流场进行数值模拟.通过流场分析可知，相比于常温环境，颗粒在燃烧环境下所形成的流场更加稳定，这是因为在相同的颗粒尺寸和气体来流速度条件下，燃烧状态的颗粒雷诺数更小.通过颗粒表面的应力分布分析可知，整体上，在相同的颗粒尺寸和气体来流速度条件下，相比于常温环境，颗粒在燃烧时所受黏性力更大，所受压力更小.通过颗粒终端速度分析可知，在低颗粒雷诺数情况下，燃烧环境中气体黏性更大，颗粒黏性力更大，从而使得颗粒曳力更大；而在高颗粒雷诺数情况下，燃烧环境中气体黏性更大使得流场更加稳定，颗粒压力更小，从而使得颗粒曳力更小.最后，本研究分别统计了半径为100～500 μm 的焦炭颗粒在非燃烧与燃烧状态下的终端速度以及相应的曳力系数，并分别将这两个物理量在不同温度环境下的比值，即ut/ut,0与 Cd/Cd0，作为经验修正值，提出了修正公式，以期为未来研究处于燃烧状态下颗粒(群)的动力学特性以及指导设计高效清洁的燃烧技术提供基础.值得说明的是，本研究对颗粒描述采用简化的实心球体假设，与实际工业应用中颗粒的多形态多孔结构有所差别.燃烧环境温度的设置参考自平均温度较低、温度范围较窄的循环流化床.因此，对于具有更加复杂结构的颗粒以及以煤粉锅炉为代表的具有更高温度的燃烧环境[14-15]，需要未来作进一步的研究.