1. 引言
生物质热解技术是实现农林废弃物资源化利用的重要途径,流化床炉因传热效率高、原料适应性强成为生物质热解的主流反应器[1]。进气速度作为流化床炉的核心操作参数,直接影响炉内颗粒流化状态、气固混合效率及热量传递过程,进而决定热解反应的效率与稳定性[2]。目前关于流化床炉进气速度的研究多聚焦于气化过程,针对秸秆热解专用流化床炉的进气速度与温度场耦合规律研究仍需完善。
本文以秸秆为热解原料,构建流化床炉二维气固多相流仿真模型,对比不同进气速度下炉内颗粒固含率及温度场分布特征,确定最优进气速度参数,为秸秆热解流化床炉的工艺优化与高效运行提供理论依据[3]-[6]。
2. 仿真模型建立
2.1. 物理模型设置
为了秸秆颗粒高效热解,热解流化床内底部通常充满石英砂颗粒,即床层,从而实现其与新鲜的生物质进料的迅速热交换[7]-[9]。为简化计算流程,节省计算资源,采用二维模型建立炉腔部分,该气–固–固多相流态化体系的热解流化床设计及其二维模型几何示意图如图1所示。加热的氮气通过底部空气分配器进入流化床,促使石英砂颗粒处于膨胀流态化状态;秸秆颗粒通过上部进料口进入热解流化床。二维热解流化床炉腔的关键几何结构参数包括炉腔直径D、炉腔高度H、石英砂床层高度h0,如表1所示。
Table 1. Basic parameters of the furnace chamber
表1. 炉腔基本参数
参数 |
数值 |
炉腔直径D/mm |
1000 |
炉腔高度H/mm |
4000 |
石英砂床层高度h0/mm |
1000 |
Figure 1. Pyrolysis fluidized bed design and schematic diagram of its two-dimensional furnace chamber model
图1. 热解流化床设计及其二维炉腔模型示意图
2.2. 控制方程及假设
在本次模拟中,考虑了流体和固体间的流动和传热。应用基于循环流体力学框架的欧拉–多相流模型模拟了流化床中气动输送的流体。在此模型中,气相为流化气体,固相为秸秆与石英砂颗粒,所有相都被视为彼此相互渗透的连续介质,具有相同的压力。动量守恒方程通过引入动量交换系数和固体压力实现封闭,而固体的物理性质则根据颗粒流动力学理论进行确定[10]-[12]。
模拟流化床的基本假设主要包括:
1) 采用欧拉–多相流模型来模拟瞬态过程,其中气相和固相被视为相互穿透的连续体,且气相被视为不可压缩流体。
2) 为了封闭和求解多相流中的连续性和动量守恒方程,利用颗粒流动力学理论。
3) 由于颗粒间的距离相对较大,颗粒间的相互作用力主要体现为颗粒碰撞,因此在流态化过程中,可以忽略颗粒间的摩擦作用。
4) 粘度模型采用了层流假设,主要考虑了颗粒碰撞和气固阻力对床层流型的影响,而气流的湍流效应则被忽略。
5) 流化床在0.5 s后达到稳定的流化状态,其中流化气体与床层材料的混合和流化分布均匀,接着秸秆颗粒进入反应器,持续1.5 s。
气相为流化气体,固相包括石英砂和生物质颗粒两相。因此,基于假设控制连续性方程可以描述为:
(1)
(2)
(3)
(4)
式中:
——体积分数,%;
——密度,kg/m3;
——速度矢量,m/s;
——气相;
——石英砂固相;
——生物质颗粒固相。
动量守恒方程由Navier-Stokes方程计算,其中包括相之间的动量交换系数。三个相的动量方程表示为:
(5)
(6)
(7)
(8)
(9)
式中:
——压力,Pa;
——动量交换系数,kg/(m3∙s);
——粘性应力张量,Pa;
——重力加速度,m/s2;
——剪切粘度,Pa∙s;
——体积粘度,Pa∙s。
能量方程可表示为:
(10)
(11)
式中:
——定压流体比热容,J/(kg∙K);
——流体温度,K;
——流体导热系数,W/(m∙K);
——流体域源项,W/m3;
——固体密度,kg/m3;
——固体定压比热容,J/(kg∙K);
——固体温度,K;
——固体导热系数,W/(m∙K);
——固体域源项,W/m3。
流体与固体交界处,流体域和固体域的能量相同,但方向相反:
(12)
(13)
2.3. 工况设计
为分析热解流化床炉的内部温度分布,通过前期文献调研分析,提出氮气进气速度v为工况变量,具体工况参数如表2所示,其中,工况1为基础工况。
Table 2. Computational model operating conditions
表2. 计算模型工况
工况 |
炉腔高度H/mm |
氮气进气速度v/m∙s−1 |
秸秆颗粒粒径d/mm |
1 |
4000 |
1 |
1 |
2 |
4000 |
0.5 |
1 |
3 |
4000 |
1.5 |
1 |
2.4. 网格划分及独立性验证
使用Ansys Meshing对热解流化床炉腔二维模型进行网格划分,网格如图2所示。
Figure 2. Overall furnace chamber meshing and local meshing division diagram
图2. 整体炉腔网格及局部网格划分图
为了确保网格划分的精确度,本章进行了网格独立性验证,通过选取几种不同的网格尺寸来进行比较炉腔竖直中心轴方向最大风速和最小风速的变化。从图3中可以看出,随着网格数量的增加,炉腔竖直中心轴方向最大风速呈现出一定的增长趋势。具体而言,在秸秆入口质量流量为8 kg/s,氮气进口温度为773 K时进行网格无关性验证,当网格数量从92,061增加至125,729时,竖直中心轴方向最大风速的增幅为0.7%,最小风速的变化范围在1.5%以内,两者的变化范围相对较小,表明网格数量的增加的影响变小。
(a)
(b)
Figure 3. Grid independence verification. (a) The maximum wind speed at the vertical central axis of the furnace chamber as a function of the number of grid elements; (b) The minimum wind speed at the vertical central axis of the furnace chamber as a function of the number of grid elements
图3. 网格无关性验证。(a) 炉腔竖直中心轴最大风速随网格数量变化;(b) 炉腔竖直中心轴最小风速随网格数量变化
在进行模拟时,需要确保所选的网格尺寸满足模拟的精度要求。通过对比不同网格尺寸下的炉腔竖直中心轴方向最大风速和最小风速变化趋势,可以判断网格尺寸对模拟结果的影响。在本研究中,当网格数量为92,061时,模拟的精度已经达到要求。另外,考虑到计算资源的有限性,既需要在满足模拟精度的前提下,同时合理选择网格尺寸以降低计算成本。因此,在综合权衡计算资源和模拟精度后,本章最终选择网格数量为92,061的网格划分方式进行后续的模拟分析。该网格划分方式既能保证模拟结果的准确性,又能有效利用计算资源,提高计算效率。
2.5. 曳力模型及数值方法
气体与固体之间的阻力是动量的交换,这对欧拉–多相流模型的预测结果有很大的影响。为了节省时间,使用了压力–速度耦合的SIMPLE方案。对于空间离散化,梯度项选择the least squares cell-based scheme,压力项和能量项选择second-order upwind interpolation。松弛因子的设置保持默认值。收敛标准也保持默认值:能量方程为10−6,其他方程为10−3。时间步长∆t = 0.01 dp/u。对Re < 1000模拟采用层流模型,对Re > 1000模拟采用LES模型。通过监测
和Nu的曲线以及颗粒的平均温度,判断流动和传热是否达到稳定的。
用曳力系数
来表征颗粒与流体间的流动及传热。曳力系数
定义为:
(14)
2.6. 物性及边界条件
流体氮气的温度为773 K,床层颗粒预热温度为600 K,颗粒与流体的温差不大,为了简化研究,忽略温度变化引起流体物性的变化,因此将流体密度设置为773 K时的密度。考虑到颗粒为致密固体,且无相变和反应发生,因此在此范围内,石英砂密度和导热系数均取773 K时的值。物性参数见表3。
Table 3. Material physical property parameters
表3. 材料物性参数
物性 |
氮气(流动相) |
石英砂(固相) |
/kg∙m−3 |
0.432 |
2700 |
/J∙(kg∙K)−1 |
1.077 × 10−3 |
850 |
/W∙(m∙K)−1 |
0.0477 |
2.7 |
/Pa∙s−1 |
3.26 × 10−5 |
/ |
为简化计算,对模型的边界条件进行以下假设:1) 流动状态为二维瞬态流动;2) 流体是连续的、不可压缩的;3) 忽略流体的粘性热耗散及物性变化;4) 流体在壁面处均为无滑移边界。数值模拟过程中,设定固体炉壁的材料为钛合金,钛合金材料具有较高的强度和耐腐蚀性,广泛应用于各种高要求环境。流动介质设定为气相,气相具有复杂的物性参数和传热特性,能够更好地模拟实际工况。
入口条件设置为质量流量入口边界条件,这意味着入口处的质量流量是可控的,通过调整入口质量流量参数,可以在不同的雷诺数状态下进行流动模拟。入口处的温度设定为恒温773 K,出口条件设置为压力出口,以模拟流体在出口处的压力变化。
在模拟过程中,为了更准确地反映实际情况,在炉壁上施加了绝热条件。在实际热解过程中,炉壁的温度是相对均匀的,且外侧覆盖有保温层,因此可以等效为绝热壁面条件。此外,炉壁表面的边界条件设置为无滑移壁面,这是因为在流动过程中,炉壁表面不会发生滑移现象,这有助于确保模拟结果的准确性。
3. 不同进气速度对炉腔温度的影响研究
进气速度是影响温度均匀性的重要模型参数,通过分析不同进气速度下的体积分数分布和炉腔温度分布均匀性,确定炉腔的最佳操作工况参数。
3.1. 体积分数分布
为了进一步揭示床层内固相颗粒随时间演化的流动特征,图4与图5分别展示了工况1下的石英砂颗粒和秸秆颗粒在0~2 s内的瞬时固含率分布演化过程。图中颜色由红至蓝分别对应颗粒浓度从高到低的变化,其中红色表示固相颗粒高度富集区域,蓝色代表气相主导区域。
Figure 4. Instantaneous solid content distribution map for quartz sand
图4. 瞬时石英砂固含率分布图
Figure 5. Instantaneous solid content distribution map for straw pellet
图5. 瞬时秸秆颗粒固含率分布图
如图4所示,从0 s图像可见,初始状态下石英砂颗粒稳定堆积于床层底部,形成结构紧密的高固含区。随着氮气从底部注入,至0.5 s时,气体开始在床层内部形成扰动,局部区域出现体积分数下降的迹象,表明床层开始发生初步膨胀。进入1 s后,大量气泡在床底生成并上升,气固两相之间发生强烈相互作用,床层出现明显的流化结构,同时颗粒开始呈现非均匀分布。至1.5 s时,床层整体呈现鼓泡流化特征,固相颗粒在氮气气流的推动下形成多处低固含率的气泡通道,高固含颗粒团聚区呈块状分布,边界处颗粒剧烈扰动。最终在2 s时,整个床层已完全流化,气泡上升及破裂导致石英砂与秸秆颗粒在垂直方向形成错落交错的混合区,石英砂颗粒在气泡上升与破裂驱动下呈现明显的周期性抛洒与回落特征,床层中部形成动态流动主通道,固含率呈现非均匀分布状态,系统进入了稳定流化阶段。
与此同时,图5反映了秸秆颗粒在进料过程中的运动特征。0 s时秸秆尚未进入床层,整体固含率接近于零。至0.5 s时,秸秆从进料口注入反应器,并逐渐在床层中部和底部区域累积,局部固含率显著上升。1 s为秸秆分布的高峰期,此阶段秸秆颗粒在气固拖曳力和碰撞作用下充分进入主流区,与床料发生强烈混合,大量颗粒被均匀输送至床层中心区域,并在与石英砂碰撞和混合的作用下形成相对集中的混合区。1.5 s以后,随着时间推移,流化过程由局部扰动逐步演化为大尺度气泡流,固相颗粒的体积分布发生明显拓展,颗粒混合效果显著增强,为后续的热解反应提供了良好的微观传热与相间接触基础,秸秆颗粒逐渐向上迁移并稀释分布,表明其在石英砂流化体系中的迁移路径受到气泡扰动和流体夹带的共同作用。至2 s时,大部分秸秆颗粒已扩散至整个床层中上部区域,但仍保留有局部团聚现象,说明混合尚未完全均匀化。
可以看出,石英砂与秸秆颗粒在流化床内表现出典型的异种颗粒耦合流动特征。石英砂构成主要床料框架,为秸秆提供了受热与分散基础,而气体引发的流动扰动和气泡通道作用则促进了异种颗粒之间的相对运动和空间混合。这一动态过程为后续热解反应的均匀性和热效率提供了物理保障,也为多颗粒系统建模与验证提供了支持。
为探究进气速度对床层内气固两相行为的影响,图6与图7分别展示了在工况2 (进风速度0.5 m/s)与工况3 (进风速度1.5 m/s)条件下的石英砂颗粒与秸秆颗粒瞬时固含率分布图。两组图对比反映了在不同进气风速作用下,床层内颗粒扰动强度、流化行为及颗粒分布特征的显著差异。
Figure 6. Instantaneous solid content distribution map for quartz sand. (a) Operating Condition 2: Inlet air velocity 0.5 m/s; (b) Operating Condition 3: Inlet air velocity 1.5 m/s
图6. 瞬时石英砂固含率分布图。(a) 工况2:进风速度0.5 m/s;(b) 工况3:进风速度1.5 m/s
Figure 7. Instantaneous solid content distribution map for straw pellet. (a) Operating Condition 2: Inlet air velocity 0.5 m/s; (b) Operating Condition 3: Inlet air velocity 1.5 m/s
图7. 瞬时秸秆颗粒固含率分布图。(a) 工况2:进风速度0.5 m/s;(b) 工况3:进风速度1.5 m/s
从图6可观察到,在进风速度为0.5 m/s的工况2中,石英砂在0.5~1.5 s内逐步进入流化状态,气泡由底部生成并沿竖直方向上升,颗粒整体呈现稳定的鼓泡流化结构,床层中部形成多个低固含通道,颗粒分布相对均匀。而在工况3 (1.5 m/s)中,进风速度大幅提升导致气泡穿透力增强,床层扰动剧烈,出现了较强的颗粒团聚与局部空斥现象。至2 s时,石英砂已出现严重冲床与局部流空,部分区域固含率接近0,说明过高气速已破坏正常流化状态,呈现近似“掏空”趋势。
图7显示了秸秆颗粒在两种风速下的瞬时分布情况。在工况2中,秸秆颗粒在0.5~1.5 s内逐步扩散至床层中部区域,颗粒分布较为均匀,固含率维持在中低水平。此种风速条件有助于秸秆在石英砂床层中的稳定混合,有利于后续热解反应。在工况3中,虽然初期秸秆颗粒被快速带入床层,但受强气流扰动与颗粒间剧烈碰撞影响,0.5~2.0 s内颗粒分布极不稳定,局部团聚严重,尤其在上部区域出现秸秆“冲顶”与排斥回落现象,混合状态明显劣于低风速工况。
为进一步探讨中等进风速度条件下床层内气固两相行为的稳定性与混合特性,结合工况1 (进风速度1 m/s)下的石英砂与秸秆颗粒瞬时固含率分布情况。通过对比分析可清晰判断,进风速度为1 m/s时的气固流化状态表现出良好的稳定性与混合均匀性。适度的进风速度有利于形成稳定流化状态与颗粒均匀混合,而风速过大虽可短时增强颗粒扰动,但易引发床层结构破坏与气固相脱耦,导致颗粒返混不足,严重影响混合与传热效率。研究结果表明,进风速度是影响床层稳定性与流化质量的关键参数,其选取需综合考虑床料密度、颗粒粒径及系统热解需求,以实现最佳工况调控。
综合上述对比分析可以得出结论:1 m/s作为一种中等风速在维持床层稳定流态、保障颗粒充分混合与防止结构紊乱方面具备明显优势。其不仅确保了石英砂与秸秆颗粒在床层内的协同流化,也为后续热解反应提供了理想的传热基础与物理支持。由此可知,在不考虑极端冲刷或局部沉积的前提下,进风速度1 m/s为当前体系下的最优操作区间。
3.2. 温度分布
为了进一步评估热解过程后期秸秆颗粒在流化床中的温度演化行为,以及炉腔整体温度分布的动态变化特征,图8与图9分别展示了模拟时长为10 s内不同时间点的秸秆颗粒温度分布图与炉腔整体温度分布图。图中颜色由蓝到红依次对应温度从300 K至773 K的变化,直观反映了颗粒受热及热场扩散过程。
Figure 8. Straw pellet temperature distribution diagram
图8. 秸秆颗粒温度分布图
Figure 9. Temperature distribution diagram of the furnace chamber
图9. 炉腔整体温度分布图
在流化床炉腔中,当氮气进气速度维持在1 m/s时(工况1),气固混合物达到了稳定的流态化状态,且核心区域的温度保持均匀。为了进一步揭示秸秆颗粒在热解过程中所经历的温度变化过程,图8展示了不同时间节点下秸秆颗粒温度的二维分布情况。初始时刻(0 s),由于秸秆颗粒尚未进料,图中炉腔内温度场几乎保持在室温水平(约300 K),秸秆颗粒温度接近最低值。随着秸秆颗粒逐步进入反应器,并与高温床料(石英砂)发生碰撞与传热。至2.5 s时,部分区域出现温度上升,表明秸秆颗粒开始受热。到了5 s左右,大量秸秆颗粒已在床层内分布,温度显著升高至500~600 K区间,呈现出中部升温、边界滞后的典型结构。进入7.5 s后,床层中部的生物质温度普遍达到650 K以上,靠近核心热解区。在10 s时,生物质温度场趋于均匀,绝大多数颗粒温度已稳定在700~750 K区间,这一温度范围与生物质热解主反应温区高度重合。说明模型有效捕捉了颗粒在动态混合与加热过程中的传热特性,并能准确反映生物质由常温向热解温区跃迁的演化过程。
如图8所示,在0 s时刻,秸秆颗粒尚未进料或刚进入床层,颗粒温度普遍较低,仅分布于300~400 K区间。至2.5 s时,受床层中高温石英砂与流化气传热作用影响,秸秆颗粒温度快速上升,温度区整体推进至500~650 K范围,图像中部与下部开始出现大面积橙黄色区域。5 s时,大部分秸秆颗粒温度均已达到650~700 K,说明床层核心区域传热已相对充分,温度场趋于稳定。随后至7.5 s与10 s,部分边缘区域温度略有波动,可能与颗粒迁移有关,但整体仍保持在热解反应有效区间,验证了床层具备持续热供能力,能满足热解反应热需求。
图9进一步从宏观尺度上反映了炉腔整体温度场的传播过程,由于氮气预热气体从底部持续送入,炉膛温度呈现由下向上的扩展态势。在0 s阶段,仅底部区域保持相对较高温度,而上部区域尚处于初始冷态。至2.5 s,热量从底部迅速向上扩散,中下部区域温度提升显著。5 s时,高温区基本覆盖炉腔中部,表明炉腔加热过程高效,热传导能力强。至7.5 s与10 s,炉腔整体温度达到相对均衡状态,尤其在反应器中部区域,温度稳定在700~750 K范围,说明热量传输均匀,具备良好的热解反应空间。
从图8与图9的联合分析可知,在持续进料和稳定流化条件下,秸秆颗粒能够迅速吸热升温,进入热解所需温度区间。同时,炉腔整体温度场维持均衡上升,热量传递过程连续稳定,未出现局部过热或温度失衡现象。此外,基础工况的模拟结果充分反映了秸秆颗粒进料、受热、升温以及热解所经历的全过程,秸秆颗粒温度与炉腔气相温度之间具有较高的一致性,验证了模型在耦合气固传热与气力输送中的准确性。这表明数值模型对颗粒热行为的捕捉具有较高准确性,为进一步研究不同操作参数对床层温度分布与热解效率的影响提供了可靠依据。
为研究进气速度对生物质热解过程温度响应的影响,图10与图11分别展示了在工况2与工况3下的秸秆颗粒温度分布图与炉腔整体温度分布图。从图10可见,秸秆颗粒在进风速度为0.5 m/s时,升温过程较为缓慢,至5 s才进入600~650 K温区,颗粒整体温升过程显现出滞后性,且空间分布上升温均匀但传热效率较低。相比之下,在工况3 (图b)中,较高风速显著增强了颗粒扰动与热对流效果,秸秆颗粒在2.5 s内已达到较高温度,5 s时大部分颗粒温度稳定于680~730 K区间,说明较强的气流强化了热传输过程。然而,图中亦可观察到高风速带来的负面影响:部分区域颗粒温度出现非均匀分布,且边界处升温不充分,存在扰动强度过大导致的颗粒脱离高温主区现象。
(a) (b)
Figure 10. Straw pellet temperature distribution diagram. (a) Operating Condition 2: Inlet air velocity 0.5 m/s; (b) Operating Condition 3: Inlet air velocity 1.5 m/s
图10. 秸秆颗粒温度分布图。(a) 工况2:进风速度0.5 m/s;(b) 工况3:进风速度1.5 m/s
(a) (b)
Figure 11. Temperature distribution diagram of the furnace chamber. (a) Operating Condition 2: Inlet air velocity 0.5 m/s; (b) Operating Condition 3: Inlet air velocity 1.5 m/s
图11. 炉腔整体温度分布图。(a) 工况2:进风速度0.5 m/s;(b) 工况3:进风速度1.5 m/s
图11则进一步反映了不同进气速度下炉腔整体热场分布特征。在工况2 (图a)中,热量由底部向上逐步扩展,至10 s时炉腔中下部温度可达700 K以上,但热区推进速度较慢,上部区域温升受限。而工况3 (图b)表现出更强的热场响应能力,热量在5 s之前已显著扩展至中上部区域,至10 s时炉腔内温度分布更加饱满,体现出热传导效率的增强。不过,由于高气速造成气体滞留时间缩短,热量分布仍存在局部不均现象。
综合工况1可得,适度提升风速有助于强化炉腔热对流,提升秸秆颗粒升温速度,缩短热解反应起效时间,但当风速过高时,会带来颗粒迁移不稳定、局部过热与传热非均匀等问题。因此,在气固热解系统中,进风速度应控制在可实现有效扰动但不破坏床层结构的区间。因此,可以得出1 m/s工况下颗粒温升响应快、温区分布均衡、炉腔热场稳定性优,综合性能最佳。
4. 结论
针对进风速度参数开展了对比分析。模拟结果表明,在低风速(0.5 m/s)下,床层扰动有限,固相混合效率偏低,温升过程滞后;而高风速(1.5 m/s)虽可增强对流传热效率,但亦导致局部颗粒脱离、冲刷严重、温度分布不均等不良现象。对比分析进一步确认中等风速1 m/s在颗粒混合、热响应速度与热场稳定性之间实现了良好平衡,适合作为优化设计参数。
基金项目
项目来源:2026年度辽宁省教育厅高校基本科研项目–创新团队项目;项目名称:油页岩低碳干馏原创技术团队;项目编号:LJ222610148002。
NOTES
*第一作者。
#通讯作者。