随着畜牧业和养殖业技术的高速发展，饲料产量快速增长。饲料制成颗粒具有提高饲料的消化吸收率、改善蛋白质的消化、延长保质期等优点，其生产备受关注和重视[1-2]。经挤压膨化或搅拌混合后，饲料颗粒湿含量约为20%，不利于产品储藏，将饲料颗粒干燥至安全湿含量可以避免其发生霉变[3]。振动流化床干燥是一种强化干燥技术，具有提质、增效、降耗等优点[4-5]。振动流化床结构简单、易于维修、安全性好，应用于饲料工业潜力巨大。振动流化床干燥是在流化床干燥的基础上加装振动电机，使干燥室内物料在振动和气流的双重作用下实现流态化，使热空气与物料充分混合，解决了物料结块、难以流化和干燥不均匀等问题。振动流化床干燥在农业、食品、化工、医药、矿冶等工业行业的应用逐渐增多[6]。Zamani等[7]通过正交实验和响应曲面法，探讨甘菊花振动流化床干燥的最佳工艺，发现最佳工艺参数为振动频率10.88 Hz、进口热空气温度64.08 ℃、流速20.63 m3/h。Gulzar等[8]根据振动流化床干燥的传热传质系数高、干燥时间短特点，进行甜瓜种子干燥的试验，发现振动频率6~11 Hz、振幅10 mm时，进口热空气温度与干燥时间呈反比；振动流化床干燥温度为80 °C时，只需干燥25 min。Sinha等[9]进行茶叶振动流化床烘干机干燥热力学㶲分析，发现采取适当措施可以避免热传导损失，大幅提高能耗的㶲效率。Stakić等[10]针对振动流化床对种子（细颗粒物料）对流干燥情况进行数值模拟和试验，发现较高的气体温度和速度能够实现快速干燥，且对于较深的床层物料更明显。Cano-Pleite等[11]运用双流体模型对立式振动流化床中孤立气泡的尺寸和运动进行了数值研究，发现振动幅度会影响气泡直径和速度振荡的大小，气泡特性的延迟主要受床振动频率的影响。张健平等[12]借助振动流化床干燥试验装置，分析油菜籽流化床干燥工艺对流传热、传质系数的影响，建立回归模型。结果发现，对流传热，传质系数两个模型预测值与试验值的最大相对误差仅为4.83%和4.79%，表明该两个回归模型拟合度较好、可靠性较高。叶飞等[13]发现，使用组合电磁加热耦合热风、热管余热回收和流化床干燥等技术优化茶叶加工关键工艺，可以提高茶叶的滋味和色泽品质，使产能提高37.5%、热效率提高100%、用工成本减少50%。伊宁等[14]根据振动流化床干燥工艺的特点，结合生产经验探讨聚丙烯酰胺振动流化床干燥法的最佳工艺条件。骆振福等[15]通过分析气泡在振动流化床中的受力情况，推导出气泡在振动流化床中生成的时间计算关联式，发现振动对气泡的形成有抑制作用。李秀芹等[16]认为，振动及气流的双重作用可以驱使整个干燥床层物料松散、流化均匀，减少沟流、腾涌、附壁等现象的发生，振动的引入有助于床层中传热传质过程的强化。刘安源等[17]针对振动流化床干燥装置干燥特性，重点进行计算方法研究。上述研究表明，通过振动辅助实现气粒两相流态化，能够大幅提高干燥过程传热传质系数，实现干燥过程提质、增效、降耗目标。因此，干燥室气粒两相流动特性对于振动流化床干燥传热传质特性有至关重要的影响。颗粒饲料通常采用固定床干燥机、带式干燥机、滚筒干燥机或气流式干燥机进行干燥。这些干燥机普遍存在饲料产品质量不高、干燥效率较低、干燥能耗大等问题[6,18]。研究人员尝试将振动流化床干燥技术应用于饲料工业。苗帅等[19]介绍1种应用于大颗粒膨化饲料干燥的大型双质体振动流化床，通过大量现场调试，实现优化振动频率、热空气流速等重要运行参数，旨在实现节能高效目标。饲料颗粒的振动流化床干燥过程是一个复杂的稠密气粒两相流传热传质问题，现场调试费时费力、占地大，具有很大的局限性。因此，文章基于计算流体动力学方法，通过数值模拟手段，以振动流化床干燥室为研究对象，应用Eulerian-Eulerian双流体模型和颗粒动力学理论构建符合实际的饲料工业振动流化床干燥数学模型，研究饲料颗粒振动流化床干燥过程的气粒两相流特性。以干燥器振动频率、进口热空气速度为主要影响因素，数值模拟分析饲料颗粒在振动流化床中的体积分数分布情况以及颗粒流态化现象，寻找最佳工艺参数，为饲料颗粒振动流化床干燥传热传质过程的定量描述和优化控制提供参考。1数值模型的建立1.1物理模型饲料颗粒振动流化床干燥装置见图1。饲料颗粒振动流化床干燥装置主要由振动流化床干燥室、加热器、测控系统、隔板、鼓风机构、旋风分离器和引风机等组成。10.13557/j.cnki.issn1002-2813.2021.10.023.F001图1饲料颗粒振动流化床干燥装置饲料干燥主要在分布板上部进行，即振动流化床的上箱体内部，建立模型时以上箱体的内部为主要研究对象。为简化模型，减少计算机运算量，认为上箱体是1个具有对称性结构的长方体，热空气均匀进入干燥室内。只以上箱体的横截面作为振动流化床的物理模型便能较好的研究气固两相的传热传质过程，同时亦可实现降维。振动流化床干燥室横向截面见图2。干燥室物理模型的尺寸为300 mm×1 200 mm。其中，h0为静止时的床层高度；顶部为干燥后的废气出口；底部为热空气入口。10.13557/j.cnki.issn1002-2813.2021.10.023.F002图2振动流化床干燥室横向截面1.2数学模型结合振动流化床的干燥流程及饲料颗粒的物性，在建立数学模型时做如下假设：忽略干燥过程中饲料颗粒的直径变化；颗粒与颗粒之间由于温差会发生热传递不予考虑，以及干燥室壁的热损耗也不考虑；颗粒相与热空气相相互作用的曳力服从Gidaspow模型[20-21]；基于连续介质假设，将饲料颗粒相采用拟流体相处理，可列出各相的能量守恒方程、质量守恒方程和动量守恒方程来求解。颗粒相和热空气相的湍流模型采用Renormalization-group (RNG) k-ε湍流模型。1.3控制方程组1.3.1质量守恒方程热空气相连续性方程。∂∂u(atρt)+1v∂∂v(atρtvhtv)+∂∂w(atρthtw)=m˙bt (1)式中：at为热空气相的体积分数；ρt为热空气相的密度（kg/m3）；htv、htw分别为热空气相的速度矢量ht在笛卡尔坐标体系下的轴向和径向分量（m/s）；m˙bt为热空气相的增加速率[kg/(m3·s)]。饲料颗粒相的连续性方程。单位体积内颗粒相由于蒸发，会传递到热空气中。∂∂u(arρr)+1v∂∂v(arρrvhrv)+∂∂w(arρrhrw)=m˙br (2)式中：ar为颗粒相体积分数；ρr为颗粒相的密度（kg/m3）；hrv、hrw分别为热空气相的速度矢量hr在笛卡尔坐标体系下的轴向和径向分量（m/s）；m˙br为单位体积单位时间内从热空气相传递到饲料颗粒相的质量[kg/(m3·s)]。上述式子根据体积分数守恒满足at+ar=1而颗粒相蒸发的水分和热空气增加的水分满足式子m˙bt=-m˙br。1.3.2能量守恒方程在欧拉多相流模型中，针对颗粒和热空气相的能量守恒，采取用焓表示的能量方程。∂∂u(atρtHt)+∇⋅(atρtu⃗tHt)=-at∂P∂u+τ¯¯t:∇u⃗t+Qbt+m˙btHevp (3)同理，饲料颗粒相的能量方程描述如下：∂∂u(arρrHr)+∇⋅(arρru⃗rHr)=-ar∂P∂u+τ¯¯r:∇u⃗r+Qbr+m˙brHevp (4)式中：Ht、Hr为热空气和饲料颗粒相比焓（J/kg）；Qbt、Qbr为对流传热速率（W/m3），Qbr=-Qbt；Hevp为水的蒸发（或冷凝）潜热（J/kg）。1.3.3动量守恒方程热空气相的动量守恒方程如下：∂∂u(αtρtu⃗t)+∇⋅(αtρtu⃗tu⃗t)=-αt∇P+τ¯¯t⋅∇+atρtg⃗+R⃗br+m˙btu⃗bt (5)同理，可得到饲料颗粒相的动量守恒方程式。∂∂u(αrρru⃗r)+∇⋅(αrρru⃗ru⃗r)=-αr∇P-∇Pr+τ¯¯r⋅∇+arρrg⃗+R⃗bt+m˙bru⃗br (6)式中：u⃗t为热空气相速度矢量（m/s）；u⃗r为颗粒相的速度矢量（m/s）；P为操作压力（Pa）；Pr为饲料颗粒相压力（Pa）；u⃗bt、u⃗br为气固两相间的相对速度（m/s）。R⃗bt、R⃗br为气固两相间的曳力（N）；曳力计算如下：R⃗br=Kbr(u⃗t-u⃗r) (7)Rbr=-Rbt (8)式中：Kbr为动量传递系数（kg/s）。公式（1）~（8）组成了描述振动流化床干燥过程的控制方程组，公式（6）中颗粒相的压力梯度项∇Pr、公式（7）的动量传递系数Kbr都选用Gidaspow模型[20-21]表示；湍流模型中热空气相选用Renormalization-group k-ε模型预测，饲料颗粒相则采用Techen理论模型[22]描述。1.4饲料颗粒物性为模拟研究热空气速度和干燥器振动频率对饲料颗粒的影响，模型中饲料颗粒的物性参数须测出，模型参数列表见表1。10.13557/j.cnki.issn1002-2813.2021.10.023.T001表1模型参数列表参数名称数值绝干颗粒密度/(kg/m)1.17绝干颗粒比热/[kJ/(kg·K)]1.33初始干基湿含量/(kg/kg)0.21初始温度/K300.00空气的密度/(kg/m)1.25操作压力（绝对）/Pa1.01×105干燥室高度/m1.20干燥室宽度/m0.30静止床层高/m0.10静止床层空隙率0.35振幅/mm3.002模型的求解2.1网格的划分针对振动流化床干燥的干燥室情况，不考虑上箱体出气口影响及下箱体进气口影响，可将干燥室简化为矩形，并进行网格划分尺寸为5 mm×10 mm，网格的尺寸应大于颗粒的尺寸。振动流化床干燥室二维计算网格见图3。10.13557/j.cnki.issn1002-2813.2021.10.023.F003图3振动流化床干燥室二维计算网格2.2边界条件2.2.1入口边界入口边界为y=0 & -150 mm≤x≤150 mm。热空气相由下向上进入，径向速度为0 m/s。饲料颗粒相在入口边界的流入速度为0 m/s。2.2.2出口边界出口边界为y=1 200 mm & -150 mm≤x≤150 mm。出口压力设置为常压，采用压力出口，解决气体回流问题对模拟结果产生影响，保证模拟收敛。2.2.3壁面边界壁面边界为x=±150 mm & 0 mm≤y≤1 200 mm。假设壁面光滑，且为绝热材料，壁面传热忽略不计。2.3数值求解根据前文Eulerian-Eulerian模型对振动流化床气粒两相流传热传质的数学描述，建立干燥过程控制方程组。方程组因由较复杂的非线性偏微方程组成，须采用有限体积法对饲料颗粒干燥过程的方程组进行离散化处理，再数值求解。干燥过程的模拟还需考虑曳力模型、饲料颗粒参数模型、振动流化床的振动模型等，将以上模型通过用户自定义函数（UDF）嵌入到CFD求解器中，最终完成干燥过程的数值模拟。3结果与分析3.1振动频率对干燥室气粒两相流的影响已知振动流化床床层空隙率为0.35，振幅为3 mm。为使模拟结果有较明显的差异，同时为避免因振动频率选择跨度过大而错失最佳值，选取的模拟振动频率f分别为0、8、12、16 Hz，不同频率下的颗粒体积分数分布见图4。图4不同频率下的颗粒体积分数分布10.13557/j.cnki.issn1002-2813.2021.10.023.F4a1（a）f=0 Hz时，颗粒体积分数分布10.13557/j.cnki.issn1002-2813.2021.10.023.F4a2（b）f=8 Hz时，颗粒体积分数分布10.13557/j.cnki.issn1002-2813.2021.10.023.F4a3（c）f=12 Hz时，颗粒体积分数分布10.13557/j.cnki.issn1002-2813.2021.10.023.F4a4（d）f=16 Hz时，颗粒体积分数分布由图4（a）可知，当振动频率f为0 Hz时，振动流化床处于静止状态。热空气不断进入，饲料颗粒没有产生明显的流化现象，随着时间的推移出现了沟流现象。这是因为没有振动辅助，热空气流动很难使饲料颗粒悬浮，只能随着时间推移从颗粒团缝中穿出，会形成沟流。由图4（b）可知，当振动频率f为8 Hz时，较f=0 Hz时流化状态好，但前期也出现沟流现象。90 s后，饲料颗粒的翻滚程度不够，颗粒分布不均匀，无法和热空气进行较好的干燥。由图4（c）可知，当振动频率f为12 Hz时，床层无明显的沟流现象，整个流化过程较均匀，没有出现大幅度的波动和回流。饲料颗粒上下翻滚程度适宜，颗粒分布均匀，干燥效果好。由图4（d）可知，当振动频率f为16 Hz时，发现有大气泡的发生，床层右边的气泡破裂导致床层左边的颗粒进行补充，出现回流现象，使颗粒分布不均匀。随着时间的推移，腾涌逐渐发生。190 s时，部分颗粒被气流夹带上升到干燥室出口附近，不利于颗粒饲料的干燥。振动频率过大，容易导致颗粒还没经充分干燥，就被夹带出干燥室，也使得干燥过程热空气流速的调节范围受到限制。因此，为保证颗粒饲料在振动流化床干燥过程时，干燥室热空气与颗粒形成良好的气粒两相流态化，尽量避免沟流、腾涌和夹带现象，提高热质传递强度，选择振动频率f=12 HZ最适宜。3.2进口热空气速度对干燥室气粒两相流的影响进口热空气速度对干燥过程及干燥质量亦可产生重要影响。因此，为研究进口热空气速度对干燥室气粒两相流的影响，在上述确定振动频率为12 Hz的基础上，分别模拟热空气速度为0.5、0.8、1.1 m/s时饲料颗粒流化状态及体积分数分布情况。不同进口热空气速度下的颗粒体积分数分布见图5。图5不同进口热空气速度下的颗粒体积分数分布10.13557/j.cnki.issn1002-2813.2021.10.023.F5a1（a）热空气速度v=0.5 m/s时，颗粒体积分数分布10.13557/j.cnki.issn1002-2813.2021.10.023.F5a2（b）热空气速度v=0.8 m/s时，颗粒体积分数分布10.13557/j.cnki.issn1002-2813.2021.10.023.F5a3（c）热空气速度v=1.1 m/s时，颗粒体积分数分布由图5（a）可知，当进口热空气速度v为0.5 m/s时，床层流化程度不明显，随着干燥过程的进行，沟流现象逐渐出现，且越来越明显。190 s时，沟流导致颗粒局部团聚，饲料颗粒的体积分数约为45%。这是因为在振动频率为12 Hz的情况下，热空气的流速较小，不足以使床层饲料颗粒上下翻滚，不能形成良好的流化状态。由图5（b）可知，当进口热空气速度v为0.8 m/s时，干燥室床层无明显的沟流，流化质量好，颗粒体积分数分布较均匀。由图5（c）可知，当进口热空气速度v为1.1 m/s时，可以清晰地观察到床层中产生了气泡。随着干燥过程的进行，气泡逐渐上升变大、破裂，饲料颗粒水分不断蒸发，颗粒质量减小变轻。190 s时，饲料颗粒直接从干燥室顶端溢出。这是因为床层压降随床层高度的增大而降低，当压降越来越小，气泡便越来越大，直至到达床层顶部破裂，使得颗粒返混严重，这些现象都不利于干燥的进行。因此，当振动频率为12 Hz时，最佳的进口热空气速度为0.8 m/s。4结论饲料颗粒振动流化床干燥过程是一个复杂的稠密气粒两相流传热传质问题。振动流化床干燥器的振动频率、进口热空气速度都对颗粒饲料的体积分数分布、流动性有显著影响；当振动频率12 Hz、进口热空气速度0.8 m/s时，干燥室饲料颗粒流态化效果最佳，饲料颗粒相体积分数分布均匀。后续研究将此流动特性数值分析的基础上，结合实验进行干燥过程传热传质的定量描述和优化控制，为饲料工业振动流化床干燥技术创新发展提供参考。