风速对毛豆颗粒流态化行为影响的仿真分析

杨晓燕，钱 丹，楼晓华，李林建，许 彤，刘冬瑞

（1.四方科技集团股份有限公司，江苏南通 226371；2.南通四方罐式储运设备制造有限公司，江苏南通 226371）

0 引言

流态化技术具有混合效率高、传热传质快等优点，在生物质能源、矿业和食品等行业应用广泛［1-2］。目前在食品行业，流态化技术主要用于蔬果类食品的快速干燥、冻结等。毛豆速冻常采用流态式速冻设备，将毛豆颗粒置于水平网筛上，在高速低温气流作用下，毛豆颗粒层产生“悬浮”，从而实现快速冻结效果［3］。

毛豆颗粒较大，属于Geldart D颗粒范畴，在食品技术领域，对于毛豆流态化的研究报道不多见。在生物质、能源领域中，有许多关于Geldart D颗粒流态化行为的研究［4-5］。在食品颗粒流态化速冻领域，常见有通过布风板的结构［6-7］、不同的风机频率和控温方式［8］、脉动气流［9］等方式改善食品颗粒流态化品质的举措。对于C类（粒径小于50 µm）和D类（粒径大于2 mm）颗粒，流化速度高，混合质量受到形状影响较大［10］。

在毛豆颗粒速冻过程中，为解决颗粒边缘部分难以流化、易粘连等问题，本文利用CFD-DEM气固耦合方法，探索毛豆颗粒流态化的行为特征，分析在不同进口风速下，毛豆颗粒的床层膨胀比和混合效率，为其速冻选择合理的进口风速，保障其能够充分混合均匀，并进行相关试验验证。为食品行业流态化设备改进和参数优化提供一定的技术支持。

1 数学模型

在CFD-DEM求解过程中，气相使用Navier-Stokes方程进行求解［11］；颗粒作为离散相，使用力学平衡方程求解［12］。对于1组具有相同大小的颗粒，颗粒的平动和转动约束方程如下：

式中 mp——颗粒质量；

vp——颗粒平动速度；

Fc—— 颗粒之间以及颗粒和边界间的接触力；

Ff→p——流体对颗粒的作用力；

g——重力加速度；

Jp——颗粒转动惯量；

ωp——颗粒转动速度；

Mc——颗粒旋转产生的扭矩；

Mf→p—— 流体相的速度梯度对颗粒产生的额外扭矩。

相互作用力方程：

式中 FD——阻力项；

FN→D—— 非阻力项，包括压力梯度力F∇P、升力FL、附加的虚质量力FVM以及其他力Fothers。

式中VP——颗粒体积；

∇P——局部压力梯度。

使用阻力系数计算流体对颗粒的作用力［13］：

式中 CD——阻力系数；

ρf——气体密度；

A'——颗粒在流体流动方向上的投影面积；

u-vp——颗粒和流体之间的相对速度。

根据不同的颗粒形状以及颗粒相的密度，阻力系数有多种计算方法。针对流化床稠密相颗粒，使用Huilin & Gidaspow阻力系数模型，计算式如下：

其中，ψ为混合参数，是与流体体积分数αf相关的函数，具体表示为：

CDErgun和CDWen&Yu为另外2种计算阻尼系数的方法，分别适用于较高粒子浓度（αs≥0.2）和较低粒子浓度（αs＜0.2），且CDErgun考虑了压力损失的情况，其具体公式为［14］：

式中φ——颗粒的球形度。

2 建模和方法

2.1 参数和仿真条件设置

流态化范围设置为500 mm×200 mm×1 000 mm（长×宽×高），毛豆初始堆积高度为100 mm，颗粒填充率为20%，见图1。

图1 毛豆颗粒初始堆积Fig.1 Initial stacking status of soybean particles

划分流体网格，网格大小为20 mm，气体密度为1.225 kg/m3，黏度为1.789 4×10-5kg/（m·s），环境温度为293.15 K；风力进口采用速度进口，方向垂直于截面；出口与大气相通，采用压力出口；CFD时间步长为0.002 s。

给定毛豆的最小过筛尺寸为10 mm，为能使颗粒通过的最小方形孔，孔尺寸纵向长宽比y/z=1.8，横向长宽比z/x=0.85，颗粒总数为7 175。在颗粒仿真软件中，设置颗粒与颗粒之间的静摩擦系数为0.2，动摩擦系数为0.3，恢复系数为0.9；颗粒与边界之间的静摩擦系数为0.3，动摩擦系数为0.3，恢复系数为0.3；仿真总时长为10 s，时间步长为0.02 s。不考虑颗粒的破碎、侵蚀等。

在食品流态冻结过程中，食品颗粒被放置在网孔均匀分布的网带上，如图2所示。

图2 网带进风示意图Fig.2 Schematic diagram of net belt air inlet

在仿真的过程中，如考虑网孔的情况，则网格的数量将呈指数级的增加，大大延长计算时间。FALAH［15］将9个内径为2 mm的圆孔等效为1个相同面积的矩形孔，进行相应的耦合计算，得到较为准确的仿真结果。因此，将有网孔的情况简化为无网孔进行讨论。

根据最小流态化速度vmf预测公式［16］：

式中 μf——气体黏度；

Ar——阿基米德数，

根据WEN等［17］的研究，C1=33.67；C2=0.036 5，计算得到毛豆颗粒vmf=1.584 m/s。

考虑壁面、颗粒间摩擦以及预测公式偏差等因素，结合相关食品流态化测试数据，设定无网带工况下进风速度为2.4 m/s。参考FALAH［15］简化网孔的做法，已知网带孔隙率为17.3%，则在有网带情况下，为保证流量相等，网孔的进口风速应为2.4/17.3%=13.8 m/s。其余条件保持相同，进行耦合仿真。

床层高度和压降是评价流态化过程的重要参数。图3为有、无网孔时压降、床层高度随时间变化曲线。在10 s的流态化过程中，床层高度呈周期性的波动。在仿真中，无网孔的平均压降峰值和吹起高度峰值均要大于有网孔的情况，但是整体偏差较小。因此，可以使用无网孔的工况代替有网孔的工况进行相似的仿真计算。

图3 压降、床层高度随时间变化曲线Fig.3 Curve of pressure drop and bed height over time

2.2 仿真结果

2.2.1 颗粒流动分布特征

0～2 s内毛豆颗粒的分布和流速分布（+z方向）如图4所示。该时间段刚好为床层高度波动的第1个周期。由于初始堆积较密集，毛豆颗粒在0.1 s呈现出整体向上悬浮的状态。随后在0.2～0.3 s出现颗粒松散分层的现象。随着床层的继续膨胀，颗粒内部中心区域出现椭圆形的气泡，0.4～0.7 s气泡逐渐增大。伴随着气泡内部翻腾的颗粒个数减少，和气泡边缘的颗粒不断碰撞，椭圆形气泡逐渐缩减，颗粒也逐渐从气泡边缘回到床层底部。在1～2 s期间，毛豆颗粒出现分层涌动的现象，在两侧出现小型气泡，气泡的边缘颗粒数较多，气泡较小，出现时间短，床层的上表面不断波动。

图4 毛豆颗粒分布（+z方向）Fig.4 Soybean particle distribution （+z direction）

2.2.2 最大床层高度和压降

设置气体流速区间为0.2～3.2 m/s，速度增量为0.2 m/s，进行16组不同进口风速下气固耦合仿真，得到压降均值和最大床层高度变化曲线，如图5所示。在固定床阶段，压降值与进口风速呈正相关；在流化床阶段，压降值增长幅度较小。拟合曲线的转折点在2.0 m/s，即毛豆颗粒的最小流态化速度为2.0 m/s，略大于预测公式计算值（1.584 m/s）。

图5 毛豆颗粒在不同风速下的最大床层高度和压降变化Fig.5 Maximum bed height and pressure drop changes of soybean particles under different wind speeds

在小于2.0 m/s风速工况下，床层最大高度在100 mm左右，近似于颗粒的初始堆积高度，此时未发生流化。从2.0 m/s开始，床层高度随着进口风速的增大而增大，进入流化状态。

2.2.3 膨胀比和混合率

在流态化的过程中，当气体不断注入床层，颗粒间的间隙增大，出现床层膨胀的现象。为了寻找毛豆颗粒的最佳流态化速度，引入膨胀比与混合效率。采用CAXA软件对图像进行处理，计算得到毛豆在流态化过程中的膨胀比：

以+z方向为正方向选取图像，在CAXA软件中使用样条曲线的方法，对颗粒床层外缘轮廓进行描绘，并测量计算出堆积面积。

在风速2.0 m/s之前，床层高度未发生明显变化，因此以此速度值为起始，计算不同风速下的膨胀比，如图6所示。随着进口风速增大，毛豆颗粒的膨胀比增大，即床层高度增高，与最大床层高度的变化情况相似。尽管风速越大，床层膨胀越明显，然而进口风速越大，对风机要求越高，能耗越高。

图6 不同风速下膨胀比Fig.6 Expansion ratio under different wind speeds

将初始时刻的毛豆颗粒均分为上、下2层，在最后时刻记录上、下2层颗粒在测量区域的占比，如图7（a）所示。从图7（b）毛豆颗粒的混合情况看，在2.0 m/s时，上层与下层颗粒基本没有混合；在2.6 m/s时，上层与下层颗粒的混合率均在50%左右，说明在床层高度方向上，颗粒混合达到比较均匀的状态。随着风速的增大，混合率基本保持不变。因此，综合考虑膨胀比、风机能耗以及混合效果，在初始堆积高度为100 mm的工况下，2.6 m/s为毛豆颗粒的最佳流态化速度。

图7 不同风速下混合效果Fig.7 Mixing effect under different wind speeds

3 流态化试验设置

流化床实验平台如图8所示。风机将蒸发器提供的冷量通过空气输送至风道中，经过网带传递到食品颗粒上，使食品颗粒产生翻滚、沸腾甚至悬浮的状态，实现食品颗粒的冷冻，最后热空气从网带上方流回蒸发器进行热交换。蒸发器外部连接冷凝器、压缩机等设备，完成制冷循环。在网带上用塑料板隔出500 mm×200 mm×400 mm的范围进行试验，将毛豆平铺在网带上。该网带的孔隙率为17.3%，毛豆的堆积厚度为100 mm。通过工业相机进行拍摄，以50帧/s的帧数拍摄10 s。

图8 流化床实验平台Fig.8 Fluidized bed experimental platform

图9是毛豆的流态仿真和试验工况的对比结果，截取5个时刻的图像。试验拍摄图片显示，毛豆床层表面出现连续的波动起伏，在风力作用下局部区域突起，随后伴随着少量毛豆颗粒的喷射行为。这与仿真现象较为一致，证明仿真结果较为准确。

图9 仿真和试验结果对比Fig.9 Comparison of simulation and experiment results

4 结语

本文以毛豆颗粒为研究对象，基于欧拉-拉格朗日方法对毛豆流态化行为进行描述，通过CFD-DEM仿真手段探究不同风速下毛豆的最大床层高度、压降变化，得到毛豆的初始流态化速度为2.0 m/s。从膨胀比和混合率指标出发，确定在初始堆积高度为100 mm时，毛豆颗粒的最佳流态化速度为2.6 m/s。研究对食品流化床设备的优化设计，以及进一步提高食品速冻过程中的混合传热效率具有一定的参考价值。