液态食品连续欧姆加热的数值模拟

徐 巍，姜 欣，沈五雄，田 婷，周家华,*

液态食品连续欧姆加热的数值模拟

徐 巍1，姜 欣2，沈五雄，田 婷1，周家华1,*

(1.华南农业大学食品学院，广东 广州 510642；2.华南理工大学轻工与食品学院， 广东 广州 510641)

为研究液态食品欧姆加热过程的温度场和速度场分布，实验采用有限元数值计算的方法，对圆柱状加热室中自来水和豆浆的连续欧姆加热过程进行三维稳态数值模拟，分析不同截面加热室速度场和温度场，以及不同流速对加热室的温度分布的影响。加热室中心温度模拟计算结果与参考文献值对比，自来水和豆浆的误差分别为2.32℃和2.65℃，数值模拟很好地预测了自来水和豆浆加热过程中的速度场和温度场。

欧姆加热；有限元；数值模拟；电导率

食品热加工单元是食品加工中的重要过程，它不仅影响食品的最终品质，而且还关系到食品生产的能耗。欧姆加热也称为通电加热、电阻加热，是利用物料本身的导电性直接加热物料的新型热加工技术[1-3]。由于欧姆加热是不存在传统意义加热面的直接加热，具有加热均匀、升温快、效率高的特点[4]，特别对热敏性食品如牛奶、豆浆、果汁等的加热具有很大优势[5-6]。因此，国内外对液态食品欧姆加热进行了大量研究，Tao等[7]研究了新鲜橙汁的欧姆加热过程，结果表明温度能够影响新鲜橙汁的电导率。李法德等[8]自行设计制作了能够实现温度程序控制的欧姆加热装置，对豆浆的电导率随温度变化规律进行了实验分析。周亚军等[9]研究了含水果颗粒液态食品物料的欧姆加热速率，以及物料的温度场分布。耿敬章[10]研究了苹果汁的欧姆加热杀菌，分析了电压、pH值、时间及加热体积等对苹果汁中的嗜酸耐热菌杀菌效果的影响。

液态食品连续欧姆加热虽然得到一些应用，但由于电场和热对流的存在使得对温度场监测存在一定的难度，特别缺乏对流速较低区域和出口温度的监测[11]。李陆星等[12]研究了自来水和豆浆在连续通电加热过程中加热室的温度场，并采用热电偶测定加热室中心温度。由于热电偶仅能测定某些固定点温度，而不能有效监测加热室整体的温度分布，同时热电偶还会影响加热室速度场和温度场分布，具有一定的局限性。林向阳等[13]利用磁共振成像(magnetic resonance imaging，MRI)技术连续无间断地测量食品在欧姆加热过程中的温度变化，绘制物料在加热过程中的温度分布图。但是MRI方法不仅复杂而且成本高，仅限于少数研究机构应用，难以应用于工业化欧姆加热。

有限元法(finite element method)是随着电子计算机的发展而迅速发展起来的一种现代计算方法[14]。有限元数值模拟不仅可以有效地监测欧姆加热过程中温度和速度场分布情况，而且节省大量的时间和成本，是一种简单、高效的方法。本研究利用有限元数值模拟软件COMSOL Multiphysics建立液态连续欧姆加热的物理模型，以自来水和豆浆为对象研究三维稳态欧姆加热温度场和速度场分布，并分析不同入口流速对加热室温差的影响。

1 物理及数学模型

1.1 物理模型

本研究采用管道式加热室，电极板位于加热室的两侧，两侧电极板分别开8个直径为0.006m的小孔，物料从左侧开孔流入，右侧开孔流出，加热室的长度为0.09m，电极板的直径为0.0736m，厚度为0.006m，物理模型见图1。

图1 加热室物理模型Fig.1 Physical model of ohmic heater

1.2 数学模型

考虑到装置的复杂性及模拟的可行性，为简化计算做以下假设：

1)流动体系为牛顿不可压缩流体。2)液态食品各向同性，其密度、热导率、比热容和黏度在研究的温度范围内变化不大，忽略对研究结果的影响；电导率则采用温度的函数表示。3)加热室绝热性好，加热管道外壁与外界环境无热交换。4)电极板的加热功率损失为零。

1.3 数学控制方程

1.3.1 质量守恒方程

由于自来水为不可压缩流体，故可以得到自来水的连续方程[15]为：

式中：ux、uy、uz分别为流体在x、y、z坐标上的速度分量/(m/s)。

1.3.2 动量守恒方程[15]

式2～5中：η为动力黏度/(Pa·s)；T为温度/ K；k为湍流动量/(m2/s2)；ε 为动量扩散比率/(m2s3)；ρ为密度/(kg/m3)； 为哈密顿算子；各个经验常数为：Cu=0.09、Cε1=1.44、Cε2=1.92、σε=1.3、σx=1.0。

1.3.3 能量守恒方程[16]

式中：Cp为比热容/(J/(kg·K))；K为热传导率/ (W/(m·K))；Q为物料的内源热即电流产生的热量/J；

电流产生的热量：

2 有限元分析

2.1 边界条件设置

湍流边界条件：左侧入口(z轴正方向)自来水、豆浆的速度分别为0.05、0.62；出口的压强都为0Pa。

温度边界条件：加热室外壁边界条件采用第一类边界条件，即热绝缘假定物料与环境的热交换量为零；左侧入口温度T=293.15K；右侧对流出口=0。

电场边界条件：左侧电极板控制U=220V；右侧电极板接地U=0V；模型外表面绝缘。

2.2 物性参数

影响食品欧姆加热温度的关键因素是电导率,而电导率为食品本身的物性参数，主要为温度的函数。为准确模拟欧姆加热过程，本研究考虑温度对自来水和豆浆电导率的影响。同时由于在实验温度范围内，自来水和豆浆的密度、热导率、比热容和黏度变化较小，对结果影响不大，为简化计算，忽略温度对其影响。自来水的密度、热导率、比热容和黏度分别为1000kg/m3、0.6W/(m·K)、4180J/(kg·K)和0.001Pa·s，固形物含量为5.01%的豆浆[17-20]密度、热导率、比热容和黏度分别为1010kg/m3、0.58W/(m·K)、3980J/(kg·K)、0.002Pa·s。并根据参考文献[21]提供的电导率函数，设置自来水的电导率为0.00408T－1.106992，豆浆的电导率为0.00728T+0.00331。

3 数值模拟及结果分析

网格划分的好坏直接影响模型计算结果，网格过密会造成计算量巨大，求解时间过长，过疏则会导致计算结果不准确。网格的划分方法采用软件自由网格划分，其划分的过程是在物理模型的各条边上插入节点，节点间隔为边长的十分之一，然后根据节点形成四面体网格。为了更好地研究加热器进出口处的速度场和温度场分布，本研究对物理模拟中开孔处进行网格加密，四面体网格数为21292个。

3.1 自来水加热速度场分布

图2 加热自来水距左侧电极板0.01m处截面速度分布Fig.2 Velocity distribution at the tap water cross section at a distance of 0.01 m from the left electrode plate

图3 加热自来水中间截面的速度分布Fig.3 Velocity distribution at the middle tap water cross section

研究选取Z轴(加热室轴向方向)方向上各距左、右侧电极板0.01m处截面和中间截面，通过研究这3个截面速度场分布来分析加热室中整体速度场分布。

图4 加热自来水距右侧电极板0.01m处截面速度分布Fig.4 Velocity distribution at the tap water cross section at a distance of 0.01 m from the right electrode plate

由图2、4可知，入口截面和出口截面速度分布类似以十字形分布；从图2入口截面速度分布可以看出，入口流速以各个圆形开孔为中心向四周逐渐降低，最高流速为0.02m/s，并在x轴方向贴近圆筒处出现流速最低点，为0.002m/s；由图4可知，出口截面处速度以十字形向周围降低，中心最高为0.0055m/s，开孔夹角45°方向贴近壁面处流速最低，为0.001m/s；由图3可知，在加热室轴向中心横截面处物料的流速为圆形分布，流速由中心到壁面逐渐降低，流速最高为0.007m/s，最低为0.0005m/s。物料从入口截面处的十字分布发展为中间截面的圆形分布；最大流速从入口截面0.02m/s逐渐降低至出口截面0.0055m/s；速度较低区域出现在加热室出口壁面处。

3.2 自来水温度场分布

图5 加热自来水距左侧电极板0.01m处截面温度分布Fig.5 Temperature distribution at the tap water cross section at adistance of 0.01 m from the left electrode plate

图6 加热自来水中间截面的温度分布Fig.6 Temperature distribution at the middle tap water cross section

图7 加热自来水距右侧电极板0.01m处截面温度分布Fig.7 Temperature distribution at the tap water cross section at a distance of 0.01 m from the right electrode plate

模拟自来水连续欧姆加热过程温度场分布，研究选取Z轴(加热室轴向方向)方向上各距左、右侧电极板0.01m处截面和中间截面，通过研究这3个截面温度场分布来分析加热室中整体温度场分布。

图5～7表明，加热室温度分布大致以加热室中心为对称，在接近壁面的45°方向处温度较高，而在开孔沿线方向温度较低。当流体流经加热室时，流体质点除沿加热室轴向运动外，还在加热室的径向上有速度分量，质点间彼此碰撞、相互混和[21]。由于加热室采用圆筒状并且物料入口为十字开孔，这就使得加热室极板间开孔沿线方向速度较大，而贴近加热室壁面45°方向区域的流速较低，形成滞留区域，滞留区的物料与加热室主体的物料对流较慢，加热时间长，因此形成了加热室的高温区。

加热室最高温度从入口截面处的54℃逐渐降低到中间截面的44℃然后又升高到出口截面的50℃。对比速度场图可以看出，加热室两端电极板壁面45°处出现速度最低点，此处物料遇到电极板流速较低，加热时间相对较长，由于物料温度和加热时间成正比，因此温度就会比其他区域高，而加热室开孔方向速度较高，物料温度就会较低。

[21]中加热室中心温度测定值为28.44℃，而本研究模拟计算值为30.76℃，二者相差2.32℃。出现差别的原因可能是：1)模拟假设加热室管壁为热绝缘，而实际存在管壁与空气的热传导，有一定的能量损失。2)模拟假定电极板的加热功率损失为零，而在实际生产中电极板存在一定的热损失。3)对比文献中采用热电偶监测温度值，但热电偶外表面无论是导电还是不导电，都会干扰电场分布和流体流动，其测定的只是一个近似值。考虑到以上因素，模拟计算结果和文献测定值相差很小，模拟计算可以很好地反映欧姆加热的温度场分布。

3.3 豆浆加热的温度场分布

为了进一步验证数值模拟的有效性，本研究同时还模拟了豆浆的连续欧姆加热，模拟选取的加热室温度截面和自来水相同，分别为距左、右极板0.01m处和中间截面。

图8 加热豆浆距左侧电极板0.01m处截面温度分布Fig.8 Temperature distribution at the soymilk cross section at a distance of 0.01 m from the left electrode plate

图9 加热豆浆中间截面的温度分布Fig.9 Temperature distribution at the middle soymilk cross section

图10 加热豆浆距右侧电极板0.01m处截面温度分布Fig.10 Temperature distribution at the soymilk cross section at a distance of 0.01 m from the right electrode plate

从图8～10可以看出，豆浆的温度分布和自来水温度分布整体相似，但其温度值普遍高于自来水。豆浆中含有可导电的固形物如蛋白质胶体颗粒以及无机盐离子[17]，使得豆浆的电导率要高于自来水的电导率，导致豆浆吸收的热量要稍高于自来水。加热室中心圆点处的温度数值模拟值为40.98℃，而文献[21]中的实验值为38.33℃，二者相差2.65℃。导致豆浆温度模拟值和实验温度值的误差原因，除了上述自来水产生温差的原因外，还有实际生产中由于豆浆中含有固形物，在加热的过程中固形物容易在电极板上形成污垢，影响电极板中加热效果的原因。

3.4 不同流速下的温度分布

为了考虑欧姆加热过程中的温度均匀性和整体温度的高低，优化液态连续欧姆加热过程，研究流速和温度场的关系，实验模拟了自来水在不同进口速度下加热室的最大温度差(即过热点与过冷点的温度差)和中心温度分布。

图11 不同进口流速下自来水加热室最大温度差Fig.11 Change in the largest temperature difference with inlet velocity

图11 表明，进口流速由0.02m/s增大到0.09m/s时，最大温度差从143℃逐渐下降到22℃。进口流速较小时，整体最大温度差较大，易发生局部过热现象，速度大于0.05m/s时，加热室内最大温度差较小，且随进口流速增大变化不明显，此时加热室内整体温度分布较均匀。

图12 不同进口流速下自来水加热室中心温度Fig.12 Change in the center temperature with inlet velocity

由图12可知，在流速不断增大的同时，加热室的中心温度从46.5℃下降至27℃，在流速从0.02m/s升至0.05 m/s阶段，中心温度下降速率较大，而在接近0.05m/s左右时温度值变化趋于平缓。说明中心温度可反应物料的整体温度，流速增大会造成整体温度变小。

从以上入口流速和加热室最大温度差和中心温度场的关系可知，随着进口流速增大，加热室的最大温度差、中心温度逐渐降低。这是因为当流速增大时，物料在加热室的停留时间减短并且湍动混合程度增大，从而降低了加热室最大温度差和中心温度。进口流速与加热室最大温度差和中心温度密切相关，进口流速的选择要从加热室的最大温度差和中心温度方面考虑。

4 结 论

利用有限元数值模拟的方法研究液态食品连续欧姆加热的速度和温度分布。流体流过微型圆孔进入加热室，在接近加热室两端电极板45°处速度较低，形成滞留区域，导致流体在此处停留的时间长，所以在流速最低处出现过热现象。模拟结果表明，液态食品的连续欧姆加热过程中，不同进口流速对加热室的温度分布有很大影响，自来水和固形物含量为5.1%豆浆的中心温度模拟计算值与文献测定值的差别仅为2.32℃和2.65℃。研究表明，有限元计算模拟能够很好地预测液态食品连续欧姆加热过程的速度分布与温度分布，精确的材料参数和边界条件也利于提高模拟预测的准确性。

参考文献：

[1]李里特. 食品物性学[M]. 北京: 中国农业出版社, 1998: 382-397.

[2]徐军, 陈永红. 食品欧姆加热原理及其影响因素[J]. 食品工业, 1999 (5): 42-45.

[3]王锋, 孙骊, 薛文通. 食品通电加热技术的研究进展[J]. 包装与食品机械, 2002, 20(6): 12-13.

[4]周亚军, 殷涌光, 王淑杰, 等. 食品欧姆加热技术的原理及研究进展[J]. 吉林大学学报: 工学版, 2004, 34(2): 324-329.

[5]PARROT D L. Use of ohmic heating for aseptic processing of food particulates[J]. Food Technology, 1992, 46(12): 68-72.

[6]el MOCTAR A O, PEERHOSSAINI H, le PEURIAN P, et al. Ohmic heating of complex fluids[J]. International Journal of Heat and Mass Transfer, 1993, 36(12): 3143-3152.

[7]TAO Qihua, JINDAL V K, van WINDEN J. Design and performance evaluation of an ohmic heating unit for liquid foods[J]. Computers and Electronics in Agriculture, 1993, 9(3): 243-253.

[8]李法德, 孙玉利, 李陆星. 连续通电加热条件下豆浆的电导率[J]. 农业工程学报, 2008, 24(12): 275-278.

[9]周亚军, 殷涌光, 刘微, 等. 含水果颗粒液态食品物料通电加热速度与温度场[J]. 农业机械学报, 2004, 35(5): 136-139.

[10]耿敬章. 欧姆加热对苹果汁中嗜酸耐热菌的杀灭作用研究[D]. 西安:陕西师范大学, 2006.

[11]SALENGKE S, SASTEY S K. Models for ohmic heating of solidliquid mixtures under worst-case heating scenarios[J]. Journal of Food Engineering, 2007, 83: 337-355.

[12]李陆星, 李法德, 耿爱军, 等. 液态食品连续通电加热的温度场研究[J]. 食品科技, 2007(6): 42-45.

[13]林向阳, 阮榕生, 何承云, 等. 食品中欧姆加热温度分布监测方法的研究[J]. 食品科学, 2005, 26(6): 272-275.

[14]江见鲸, 何放龙, 何益斌. 有限元法及其应用[M]. 北京: 机械工业出版社, 2006: 1-5.

[15]谢舜韶, 谷和平, 肖人卓. 化工传递过程[M]. 北京: 化学工业出版社, 2007: 39-40.

[16]de ALWIS A A P, FRYER P J. A finite-element analysis of heat generation and transfer during ohmic heating of food[J]. Chemical Engineering Science, 1990, 45(6): 1547-1559.

[17]李修渠, 李里特, 辰巳英三. 豆浆的电导率[J]. 中国农业大学学报, 1999, 4(2): 103-106.

[18]张爱萍, 孟瑞锋, 刘东红. 豆浆浓度在线检测方法与控制装置的研究[J]. 包装与食品机械, 2009, 27(4): 8-10.

[19]鲍鲁生, 刘兰芝, 倪国娟. 豆浆黏度的研究[J]. 食品科学, 1982, 3(8): 23-25.

[20]李云飞, 葛克山. 食品工程原理[M]. 北京: 中国农业大学出版社, 2002: 844-845.

[21]李陆星. 连续通电加热系统的温度分布与极板污染的实验研究[D].济南: 山东农业大学, 2007.

Numerical Simulation of Liquid Food during Continuous Ohmic Heating

XU Wei1，JIANG Xin1， SHEN Wu-xiong1，TIAN Ting1，ZHOU Jia-hua1,*
(1.College of Food Science, South China Agricultural University, Guangzhou 510642, China；2.College of Light Industry ang Food Science, South China University of Technology, Guangzhou 510641, China)

Ohmic heating is a novel heating method characterized by uniform heating, fast temperature rising and high efficiency. But its application is hindered by difficulty in instant monitoring of temperature during heating. In this study, finite element method was employed to simulate 3-D stationary continuous heating tap water and soymilk in cylindrical ohmic heater. The influence of inlet velocity on temperature and velocity fields of different cross sections was analyzed. The results were validated against the experimental results from published literature, and the differences were 2.32 ℃ for tap water and 2.65 ℃ for soymilk. Numerical simulation is effective in predicting velocity and temperature fields in liquid food continuous ohmic heating.

ohmic heating；finite element method；numerical simulation；electrical conductivity

TS201.1

A

1002-6630(2010)19-0074-05

2010-01-07

徐巍(1986—)，男，硕士研究生，研究方向为食品工程。E-mail：xuweiscau@yahoo.com.cn

*通信作者：周家华(1968—)，男，教授，博士，研究方向为食品工程和化学工程。E-mail：jiahuazhou@scau.edu.cn