基于离散元方法的高碱煤灰沉积过程数值模拟研究

2021-05-15金默刘道银陈晓平

化工学报 2021年4期

金默，刘道银，陈晓平

（东南大学能源转换及其过程测控教育部重点实验室，江苏南京210096）

引言

准东煤是一种低燃点、高热值的优质动力用煤，已探明储量约为390 Gt，在我国煤炭储存上占据重要地位[1-2]。然而，由于准东煤中碱金属（Na 与K）含量较高，在准东煤的燃烧过程中，碱金属经历蒸发-冷凝过程并在灰颗粒和受热面表面形成黏性表层，使得煤灰颗粒在对流受热面更易积灰，从而降低了锅炉热效率，严重威胁锅炉机组的正常安全运行[3-4]。近年来研究者们多采用数值模拟方法研究锅炉相关运行参数对受热面积灰的影响规律，但对高碱煤积灰问题研究尚不透彻。因此，开展高碱煤灰颗粒在对流受热面上积灰过程的模拟研究对解决高碱煤积灰问题有重要指导意义。

为模拟灰颗粒的积灰过程，研究者们基于计算流体动力学(CFD)方法开发了多种积灰子模型。现有的积灰子模型主要包括临界黏度模型[5-7]、物理平衡模型[8-10]以及熔化相占比模型[11-12]。以上模型通过比较颗粒发生碰撞时颗粒的黏度、速度和温度等参数与对应临界值的大小决定颗粒发生黏附的概率。上述模型中的相关临界值是根据实验经验公式计算得到的，受颗粒成分等多种因素的影响[13]，不同文献中相关临界值的取值有较大差别（如临界黏度的取值范围在8～108Pa·s）。Christoph 等[14]通过对比以上积灰子模型的模拟结果后认为，熔化相占比模型的准确度更高。

在上述积灰子模型的基础上，研究者们对积灰机理开展了深入探索。现阶段明确的颗粒积灰机理主要包括惯性碰撞、热泳力、沉积以及湍流扩散[15-18]。对于粗颗粒（粒径大于10 μm），惯性碰撞为主导积灰机理；对于细颗粒（粒径小于10 μm），热泳力为主导积灰机理；Xin等[18]结合实验数据与模拟结果，分别考虑了惯性碰撞、热泳力以及碱金属蒸气冷凝对准东煤灰颗粒沉积的影响，研究发现，直接由碱金属蒸气冷凝形成的积灰质量仅占初始积灰层总质量2%；Lee 等[19]结合300 kW 煤粉炉积灰实验数据，建立了管壁初始积灰层生长模型，碱金属蒸气冷凝在温度较低的灰颗粒和受热面表面形成黏性表层，导致颗粒捕集率明显提升；Walsh 等[20]通过对受热面沉积物的成分分析，指出高碱煤灰颗粒以及受热面表面黏性表层混合物的主要成分为Na2SO4，其熔点在1200 K 左右；同时，Walsh 结合Wibberley 等[21]和Wall 等[22]的分析，提出了包覆黏性表层颗粒的碰撞黏附模型，模型中颗粒黏附概率正比于碰撞时颗粒表面能与颗粒动能的比值；Liu等[23]在Walsh 碰撞黏附模型基础上，优化了已有的积灰磨损模型，充分考虑到不同粒径颗粒对已形成沉积物的磨损作用，其模拟结果与实验结果能较好吻合。以上模型虽然能有效预测积灰生长过程，但都在一定程度上简化了颗粒的碰撞/黏附过程，同时忽略了颗粒间的相互作用。本文将基于DEM 方法，采用JKR 软球模型[24]计算每个颗粒的受力情况，获取每个时间步长内颗粒与周围颗粒、壁面的距离，从而判断颗粒的黏附状态，预测受热面的积灰生长过程。

本文以流化床锅炉对流过热器为研究对象，利用开源流体力学计算软件MFIX 进行对流过热器区域包覆有黏性表层的高碱煤灰积灰过程数值模拟。重点考察灰颗粒粒径、烟气流速、烟气温度及壁面温度对管壁积灰特性的影响规律，为解决对流受热面的积灰问题提供理论指导。

1 数学模型

1.1 基本假设与边界条件

本文主要研究包覆有黏性外层的高碱煤灰颗粒在流化床锅炉对流过热器受热面上的积灰特性，所选择模拟对象为唐智等[25]所搭建的小型流化床积灰实验台。由于实验台积灰区域体积较大且本文搭建的模型对计算能力的要求较高，为加快计算速度本文对积灰区域进行了简化。在保证积灰管附近流场、温度场特征分布的条件下，本文取积灰管区域一定轴向厚度的矩形截面（尺寸为40 mm×40 mm×4 mm）为模拟区域，模拟区域中心位置为直径24 mm的恒温圆管。烟气入口设为速度入口边界条件；烟气出口设为压力出口边界条件；x 方向两边界面设为绝热无滑移壁面条件；z 方向两边界面设为对称面条件；圆柱壁面设为恒温无滑移壁面条件，具体计算域及边界条件设置如图1所示。为简化计算，进行如下假设：（1）灰颗粒和受热面表面覆盖有一层均匀的Na2SO4黏性表层；（2）Na2SO4黏性表层的厚度远小于颗粒粒径[19]，忽略其对粒径的影响；（3）除Na2SO4黏性表层表面能外，颗粒的其他性质不受温度影响。

1.2 DEM模型

颗粒相使用DEM 模型跟踪体系内每个颗粒的运动状态[26-27]，每隔固定的时间步长计算一次颗粒位置，若两个颗粒间存在重叠则计算颗粒间作用力，从而确定下一步长内颗粒的速度与位移。由于模拟中所使用的颗粒均包覆有Na2SO4黏性表层且粒径均大于10 μm，热泳力和碱金属蒸气冷凝机理对积灰结果影响不大，因此模型中仅考虑颗粒重力、颗粒间作用力以及流体对颗粒的曳力。颗粒的运动方程如下

图1 计算域与边界条件Fig.1 Simulation domain and boundary condition

阶段1：靠近阶段，随着颗粒间距的减小，非接触范德华力逐渐增加。靠近阶段Fn(δ)的计算公式[29]如下

图2 颗粒碰撞时法向接触力/位移关系(颗粒粒径dp=80 μm)Fig.2 Calculated normal force-displacement diagram of particles

1.3 计算方法

在实际模拟过程中，受限于计算能力，若按实际飞灰浓度开展数值模拟，单个工况实际耗时预计将超过600 d。为缩短模拟时间，通过增加入口颗粒质量流量人为加快积灰生长速度。模拟时，入口飞灰颗粒质量流量取3.8×10-3kg/s，此时，不同烟气入口速度下积灰壁面附近的飞灰颗粒浓度为3.9 ～7.9 kg/m3，约为燃煤流化床过热器区域实际飞灰浓度的260～520 倍。需要说明的是，虽然模拟中飞灰浓度远高于实际浓度，但模拟中飞灰体积分数最高仅有0.05%，飞灰颗粒对流场的影响颗粒可以忽略[31]。模拟工况参数根据流化床高温过热器区域环境参数设定[32]，颗粒性质参数参考Liu 等[23]的实验数据，具体数值见表1、表2。

表1 模拟工况参数Table 1 Parameters of simulated working conditions

表2 颗粒性质参数Table 2 Particle property parameter

采用开源流体力学计算软件MFIX 的有限体积法求解流动和传热控制方程。黏性颗粒碰撞模型利用Fortran 语言编写并耦合至MFIX 软件自带的DEM 模型中。颗粒间、颗粒与壁面间传热按MFIX中默认的DEM 配套传热模型求解，即在颗粒发生碰撞时。计算域网格为四边形网格，综合考虑计算精度与计算时长，网格总数选取为12800 个，其中x 方向40 个，y 方向40 个，z 方向8 个；连续性方程、动量方程和能量方程的收敛标准取10-6；流体计算时间步长设为10-5s，颗粒相计算时间步长依据工况条件进行动态调整，其量级在10-7s。

2 计算结果与讨论

2.1 模型验证

为定量对比不同工况下壁面的积灰特性，对管壁的灰颗粒碰撞率和黏附率进行了统计，碰撞率和捕集率定义为

式中，mimpact为碰撞壁面的灰颗粒总质量，kg；mstick为发生黏附的灰颗粒总质量，kg；minject为注入计算域的总颗粒质量，kg。

在已有的研究中[9,11,13]，惯性碰撞被认为是积灰中的主导机理，主要体现为较大的灰颗粒（粒径大于10 μm）由于其自身惯性无法跟随障碍物周围流线而与障碍物发生碰撞，这一现象通常用颗粒的Stokes 数描述[33]，Stokes 数为颗粒松弛时间与流体特征时间之比

式中，ρp为颗粒密度，kg/m3；u∞为流体无穷远处速度，m/s；μg为流体黏度，Pa·s；D 为管壁直径，m。当颗粒粒径较大，速度较快，即St＞1 时，颗粒运动轨迹不易受流场影响，此时颗粒的碰撞率较高。为了能快速直观获取颗粒碰撞率，Wessel 等[34]依据势流计算式推导出颗粒碰撞率Pimpact与Stokes数的计算式

式中，Steff为颗粒等效Stokes 数；Rep为颗粒Reynolds数；β1、β2、β3分别取1.54、-0.548、0.202[34]。

在不同烟气流速、不同粒径下飞灰颗粒碰撞率的模拟值和采用式（26）计算得到的计算值分布如图3 所示。由图可见模拟值与计算值基本一致，在一定程度上验证了本文模型的准确性。

在以往的研究中，研究者通过实验，采集了不同实验台中高碱煤灰颗粒积灰捕集率[23,35-36]。根据实验中给出的工况参数，采用本文建构的DEM 模型进行数值模拟，捕集率实验值与模拟值对比见表3。由表中数据可见实验值与模拟值接近而实验值均高于模拟值。可能原因是模拟中颗粒取单一粒径，取值为实验中颗粒的平均粒径，而本文建构的DEM模型忽略了热泳力等作用于细颗粒的积灰机理，从而模拟得到的颗粒捕集率偏低。

图3 不同风速、粒径下碰撞率模拟值与计算值对比Fig.3 Comparison of impact efficiency between simulation value and calculation value under different gas velocity and diameters

表3 捕集率实验值与模拟值对比Table 3 Comparison between experimental results and simulated results

2.2 烟气流速和颗粒粒径对受热面积灰特性的影响

不同粒径和入口烟气流速下的颗粒碰撞率和捕集率的计算结果如图4所示。不同粒径飞灰颗粒与壁面的碰撞率和捕集率随入口烟气流速的变化趋势一致，碰撞率随入口烟气流速的增加而升高，捕集率随入口烟气速度的增加而降低。

当颗粒粒径小于40 μm 时，计算域网格划分要求更高，颗粒数量也快速提升，现有计算能力难以满足完整积灰过程的模拟。因此颗粒粒径在40 μm以下的工况只模拟了0.1 s 的积灰时间，其碰撞率与捕集率仍被记录于图4。比较发现，颗粒的碰撞率和捕集率颗粒的碰撞率随粒径增长而增长，捕集率则随粒径增长而下降；粒径较大的颗粒积灰特性不易受烟气流速影响；相较于速度，粒径对于碰撞率的影响更为明显。在复现Liu等[23]实验的模拟中，烟气流速高(12 m/s)而颗粒粒径小(平均粒径17 μm)，颗粒碰撞率(3.28%)远低于其余工况下的碰撞率，该工况下颗粒捕集率也仅有0.22%。

图4 不同粒径、烟气流速下颗粒捕集率和碰撞率Fig.4 Capture and impaction rate of particles with different diameters and gas velocity

图5 积灰特性Fig.5 Schematic diagram of deposition characteristics

为进一步比较不同工况下的积灰特性，本文采用最大积灰厚度H 和积灰扩展角α（图5）描述积灰形状。

计算结果表明，积灰主要形成于管壁迎风面，管壁背风面基本没有颗粒沉积，这是由于管壁背风面的积灰主要由细颗粒组成且主导积灰机理为热泳力[23]，而模拟中颗粒粒径较大且模型中未考虑热泳力机理，因此本文仅讨论管壁迎风面的积灰特性。典型工况下3 s 模拟时长内壁面积灰生长过程如图6所示，积灰生长过程中，迎风面积灰整体形状大体保持不变，沿中轴线对称，最大积灰厚度达到2.67 mm，积灰扩展角为134°。

不同风速下管壁积灰特性见表4，最大积灰厚度与积灰扩展角均随风速提升而下降，但下降的幅值较小，分别为0.18 mm 和4°。管壁积灰厚度分布由图7所示，在不同风速下，积灰的形状大体保持一致且沿中线对称，由于风速提升而减少的积灰主要发生在30°～70°和110°～150°这两个区域。

图6 管壁积灰生长过程Fig.6 Deposition growing process on probe surface

表4 不同风速下颗粒积灰特性Table 4 Deposition characteristics of particles under different gas velocity

2.3 烟气温度和壁面温度对受热面积灰特性的影响

壁面温度对积灰特性的影响见表5。由表5 可知，壁面温度变化对颗粒碰撞率无明显影响，而随着壁面温度的上升，颗粒的捕集率、最大积灰厚度及积灰扩展角均有提升。这是因为壁面温度上升导致沉积于壁面的颗粒温度上升，从而使颗粒表面Na2SO4黏性边界层的表面能提升，颗粒间、颗粒与壁面间的最大黏性力增大，导致颗粒更容易在壁面发生沉积。

图7 不同风速下管壁积灰厚度分布Fig.7 Distribution of deposition thickness on probe surface under different gas velocity

表5 不同壁面温度下颗粒积灰特性（dp=60 μm，v=12 m/s，Tgas=900℃）Table 5 Deposition character of particles under different wall temperature

表6 不同烟气温度下颗粒积灰特性（dp=60 μm，v=6 m/s，Twall=550℃）Table 6 Deposition character of particles under different gas temperature

不同烟气温度对应的不同积灰特性参数见表6。与壁面温度变化对积灰特性的影响相似，不同壁面温度下灰颗粒的碰撞率基本一致，其他积灰特性参数随烟气温度的升高而提升。相较于壁面温度，烟气温度升高对积灰特性的影响较小，可能的原因是颗粒/烟气间对流传热系数低于颗粒/壁面、颗粒/颗粒间的热导率，因此烟气温度对积灰平均温度影响较小，对积灰特性的影响不明显。