食品药品冷冻干燥数值模拟的现状与展望

詹 博，周国燕，曹斌宏，王巍悦

（上海理工大学低温生物与食品冷冻研究所，上海 200093)

食品药品冷冻干燥数值模拟的现状与展望

詹 博，周国燕，曹斌宏，王巍悦

（上海理工大学低温生物与食品冷冻研究所，上海 200093)

数值模拟用于食品药品的冷冻干燥，可以对其进行指导和预测，以提高冷冻干燥效率。本文介绍了国内外冷冻干燥数值模拟的研究现状并着重介绍了冷冻干燥模型一些新的改进，分类归纳了各种模型的假设条件及优缺点，回顾并总结了针对各种数学物理模型的求解方法并分析了各种方法的优缺点。对冷冻干燥数值模拟的发展方向进行了展望，指出：结合各种测量途径建立专门的模型是冷冻干燥数值模拟的一个发展方向。

冷冻干燥，数值模拟，求解方法，综述

真空冷冻干燥技术是一项用于食品、药品、生物制品长期保存的方法。经过冻干后的制品脱水彻底，营养成分损失少、质量轻、复水性好、便于运输和长期保存。不过真空冷冻干燥技术也存在一些不足之处：由于冷冻干燥在低温和真空的条件下对物料进行加热干燥，所以能耗较大，且所需时间较长。同时，由于冻干过程及冻干制品结构的复杂性，冻干过程中所需的很多参数都难以测量得到，冻干工艺也往往需要凭借经验进行大量实验来确定。若采用数值模拟来研究和预测冻干过程，可以对冻干工艺进行指导，以便更有效地控制冻干过程，提高生产效率。

1 冷冻干燥数学物理模型

目前应用较多的冷冻干燥模型有三类：URIF模型、非稳态模型和准稳态模型。

1.1 URIF模型

URIF模型全称为冰峰均匀退却模型，最早由Sandall[1]等人于1967年提出。该模型是建立在假设升华界面温度、物料表面温度恒定以及干燥层内水蒸气分压力不变的稳态过程基础上的。URIF模型的主要假设条件是：冰区和干区之间存在着一个明显的分界面，其厚度可以忽略;冰峰均匀退却;所形成的多孔层是绝干物质。该模型优点是简单、所需参数少、求解容易、能较好地模拟形状单一和组织结构均匀的物料的升华干燥过程[2]。但此模型在计算升华干燥阶段即初始60%～90%水分脱除时比较准确，而解析干燥阶段则较差。故此模型只能用于对质量要求不高的食品冻干，或仅用于对产品的冻干时间进行预测。

URIF模型在我国应用较广。2002年，谢国山和王立业[3]在URIF模型的基础上，假设结合水的脱除随物料温升的变化率为常数，建立了用于升华干燥和解析干燥的一维模型。2003年，谢国山[4]又对URIF模型进行了优化，将升华界面温度表示为随时间变化的线性分布，代替了URIF模型采用的恒定界面温度。2005年，王博[5]等在URIF模型的基础上建立了二维传热传质模型。2007年，苏庆勇和贺广兴[6]建立了针对圆柱型物料的URIF模型。

1.2 非稳态模型

非稳态模型最早由Liapis和Litchfield[7-8]于1979年提出。1998年，Sheehan和Liapis[9]提出了瓶装药品冷冻干燥的多维非稳态模型，其实质是二维轴对称模型。该模型的主要假设条件是：升华界面厚度无限小;由水蒸气和非凝性气体构成的二元混合气体充满整个干燥层的空隙;在升华界面，水蒸气分压力和冻结层的冰相平衡;已干层中，固相和气相处于热力学平衡状态;冻结层被认为是均质的，其导热系数、密度及比热均不随时间变化，且溶解气体部分忽略不计。该模型的特点是：简化条件相对来说较少，能较好地模拟冻干过程，与实际情况比较接近，但求解较困难，所需物性参数较多。近年来有不少学者在此基础上又做了进一步改进：2004年，Zhai[10]等人提出了考虑瓶底弯曲影响的二维轴对称非稳态模型，用来预测和描述冰界面的运动;2005年，Gan[11]等建立了考虑瓶塞和瓶壁温度影响的二维轴对称非稳态模型;2007年，Salvatore和Antonello[12]在Sheehan和Liapis模型的基础上建立了一维非稳态模型，并建立了更为精确的瓶壁非稳态传热模型;2007年，彭润玲[13]考虑了多孔介质的内部结构，给出了螺旋藻冻干过程热质传递的分形模拟。

2007年，Jin Hyun Nam和Chi Sung Song[14]在非稳态方程的基础上建立了针对板状牛肉的离散后的冷冻干燥模型。虽然没有考虑小瓶对换热的影响，但与之前所提到的非稳态模型不同的是，Jin Hyun Nam和Chi Sung Song考虑了样品侧面的传热传质对整个冷冻干燥过程的影响，并且将升华区域也纳入模拟，而不是仅仅计算升华界面的移动速度。经过计算和实验验证，发现在板状样品侧面参与换热与传质的情况下，即使在样品宽度10mm，厚度5mm的情况下，侧面的传质就已占很大一部分。研究发现，升华界面的运动轨迹，是由纵向和径向同时向中心移动。在一次升华进行到最后阶段时，升华界面已基本脱离底面。说明对于板状样品，若侧面参与传热传质，传统的升华界面移动速度公式将导致模型误差较大。

1.3 准稳态模型

由于冷冻干燥过程非常缓慢，在较短的时间内，温度、压力和升华界面位置等物理量变化非常小，可以近似用常量处理，对非稳态模型或URIF模型进行处理，即可得到准稳态模型[15]。

准稳态模型最早由Dyer和Sunderland[16]于1968年提出。1996年，李云飞和王成芝[17]在URIF模型的基础上，假设冻结层和干燥层温度以及水蒸气和非凝性气体的密度不随时间变化，建立了准稳态模型。2006年，姚志华和陈杰平[15]在一维非稳态模型的基础上，假设冻结层和干燥层温度以及水蒸气和非凝性气体的密度不随时间变化，建立了准稳态方程。

2007年，Salvatore和Antonello[12]在一维非稳态模型的基础上，对瓶装物料建立了两组用于配合在线软测量的一维准稳态模型。其中，软测量[18]的工作原理，就是在常规检测的基础上，利用辅助变量与主导变量的关系，通过软件计算，得到主导变量的测量值。其核心是建立用来预测主导变量的可靠的软测量模型。

第一组准稳态模型假定干燥层在极短时间内水蒸气通量沿径向方向保持不变，冻结层温度不随时间变化，忽略瓶壁的传热及吸附水的解析且所有热量都用于升华。第二组模型较为精确且复杂，冻结层、干燥层和通过瓶壁的导热均被考虑在内，质量流与升华界面移动速度的公式与第一组模型中的一样，能量平衡考虑瓶壁换热。

Salvatore和Antonello设计两组模型的目的是用于软测量监控法，第一组模型的精度对于大部分情况都适用，而第二组较复杂的模型在瓶壁温度的测量上更有优势。同时，Salvator和Antonello在文中提到，在结合软测量的情况下，其设计的准稳态模型比非稳态模型容易计算。由于软测量提供的对参数的校正，使准稳态模型精度也达到了模拟的要求，故而在这种情况下，用准稳态模型来模拟冷冻干燥问题，优势更大。

2 求解方法

目前求解数学物理模型的方法，主要分为三类，解析解法、数值解法、模拟及实验法[19]。

2.1 解析解法

解析解法[19]包括直接积分法、分离变量法、叠加法、Duhamel积分法、映像法、小扰动法、保角映射法、拉普拉斯变换法及近似积分法等。解析解法的优点在于，不仅推导严格、表达清晰，而且为其他解法提供了比较标准。但解析解法仅适于一些较简单的模型。

目前在冷冻干燥领域，解析解法应用较少，而且由于一般都是用于简单的模型，故求出的解与实测值相差较大。不过解析解法对于冷冻干燥模型的求解依然有一定的指导意义：2002年，谢国山和王立业[3]用解析解法对一维URIF模型进行了计算，用以研究干燥室压力和物料厚度对冷冻干燥过程的影响;2005年，王博[5]用解析解法对其建立的较简单的二维URIF模型进行了求解，计算值与实测值基本吻合，但依然存在较大误差;2009年，罗瑞明等[20]用解析解法对干切牛肉一维解析干燥过程的模拟进行了计算，经过实验验证得出其模型能较好地模拟特定条件下的解析干燥过程。

2.2 数值解法

数值解法[21]是在近几十年来发展起来的一种求解方法，它在求解许多较复杂的传热传质及流动问题时比较有效，故而在冷冻干燥数值模拟上应用较广。数值解法的适应性强，只要所采用的离散化分析与求解方法得当，其结果是相当精确的。由于数值解法本身有许多相互关联和约束原则，因此只要按照其规定的程序和原则进行分析与计算，其结果逼近于实际问题是可以得到一定保证的。数值解法的缺点是，所得结果是一系列的离散数据，而不是以公式出现的。另外，数值解会有误差，而且解的稳定性也有一定的前提和条件。数值解法主要分为三类：界面追踪法、固定界面法和固定区域法。

2.2.1 界面追踪法 界面追踪法[21]是直接对原控制方程及其边界条件离散化并进行求解。这类方法包括固定步长法、变空间步长法和变时间步长法[22]。这类解法适用于一维问题，但推广到二、三维情况会使方程形式变得极为复杂，而且对存在非单调、多个界面的情形也是不适用的。

2007年，贺广兴[23]用固定步长法对一维平板状食品的冷冻干燥过程进行了模拟计算。之后，苏庆勇和贺广兴[6]用变时间步长法对基于URIF模型的圆柱形物料进行了模拟计算，并用胡萝卜进行了验证实验。所得的模拟值与实验值吻合良好。

2.2.2 固定界面法 固定界面法[21]是通过某种转换将移动界面变为固定界面，从而将移动区域问题化为固定区域问题求解。固定界面法包括等温面移动法和自变量变换法。

在求解多维非稳态模型时，Sheehan和Liapis[9]应用了自变量变换法，在模型中加入3个无量纲变量来确定升华界面的位置，用双正交配置将径向参数用无量纲曲线坐标代替，将半径用无量纲径向坐标代替。之后使用Gear法（一种向后差分的数值积分方法，适用于求解刚性常微分方程组)对模型进行求解，并得出了温度分布和干燥时间。

2.2.3 固定区域法 固定区域法[21]不需要跟踪两相界面的位置，将分区求解的导热问题化成整个区域上的非线性导热问题来处理，因而具有很大灵活性，很容易推广到多维多界面的情况。这类方法有显热容法、热焓法和等效热容法。

2.2.3.1 显热容法 显热容法[24]是把物质的相变换热过程看做是在一个很小的温度范围内有一个很大的显热容，从而把分区描述的相变问题转变为单一区域上的非线性导热问题，达到整体求解的目的。其中显热容可用公式求得或通过DSC测得[25]。

2.2.3.2 焓法 焓法是将热焓和温度一起作为待求函数，在整个统一区域建立一个统一的能量方程，利用数值方法求出焓分布[25]。

2.2.3.3 等效热容法 等效热容法[26-27]是在划分微元或节点后，将各微元或节点的焓值进行求和平均，从而得到单元等效热容。

2005年，Jin Hyun Nam和Chi Sung Song[14]在升华-解析模型的基础上，采用等效热容法并对其进行了改进，建立了离散后的非稳态模型。为了求解方便，Jin Hyun Nam和Chi Sung Song引入体积平均参数的概念，采用冰饱和度来判断每个微元的状态，用以区分干燥区、升华区和冻结区。其中，冰饱和度即冰体积占空隙体积的分数;体积平均参数有体积平均密度、体积平均定压比热和体积平均热导率，这些参数都是微元内关于冰饱和度的平均热物性参数。通过体积平均参数，Jin Hyun Nam和Chi Sung Song建立了适用于模型所有区域的基于冰饱和度S的能量守恒及质量守恒模型，并用TDMA矩阵求解器求解各方程。

2.3 模拟法及实验法

模拟法[19]包括：热质模拟、动量模拟、电模拟、水力模拟等。当数学模型过于复杂，数值计算工作量过大，以致难以实施时，模拟法是最有效、最可靠的办法，是观察和探讨各种复杂的多孔介质传热传质问题的必经途径。不过模拟所得的结果需用实验进行验证，并进行修正，才能满足精确度的需求。

针对冷冻干燥过程的模拟中，使用较多的是计算流体力学（CFD)方法。与传统的模型相比，CFD技术具有极大的优势，经费开支更小，需要较少的实验模型，但可以在较短的时间得到模拟结果，为缩短开发周期奠定了基础[28]。CFD应用优点为：可以更细致地分析、研究流体的流动、物质和能量的传递等过程;可以容易地改变实验条件、参数，以获取大量在传统实验中很难得到的信息资料，使整个研究、设计所花的时间大大减少;可以方便地用于那些无法实现具体测量的场合，如高温、危险的环境。根据模拟数据，可以全方位的控制过程和优化设计。

在众多CFD软件中，Fluent[29]作为一款专业的CFD计算软件，功能全面，适用性广，在国内使用广泛。它以近代流体力学与数值方法为理论基础，以计算机为工具，可对流体力学的各类问题进行数值实验、计算机模拟和分析研究，从而可以在很大程度上替代耗资巨大的流体动力学实验。针对不同的领域，Fluent开发了某些专用模块，其最大的优势是强大的前处理程序模块Gambit，使得建立模型、网格化处理等非常方便，特别适用于外形复杂的CFD模拟过程。同时Fluent自带的用户自定义函数（UDF)功能允许用户自定义多种边界条件，配合Fluent多种离散格式和数值解法，可以计算各领域的复杂流体问题。

由于国外的技术封锁，各类CFD软件进入我国较晚，而对于冷冻干燥内部机理的模拟，更是近几年才开始进行。2008年，彭润玲[13]等人对螺旋藻建立了非稳态模型并用Fluent结合Matlab对其进行了求解。在求解时，升华界面移动速度根据升华界面能量平衡方程利用Matlab软件算得，作为非定常边界利用UDF输入Fluent。能量控制方程和质量传输方程则用Fluent软件中的耦合求解器进行求解。计算所得的结果与实测值接近。

3 结论与展望

本文主要介绍了最近几年冷冻干燥数值模拟中新的模型、针对旧模型的修正以及相对应的各种求解方法，得到以下结论：

3.1 目前对冷冻干燥问题的数值模拟已开始向非稳态方向发展。非稳态模型精度较高，不过考虑因素也较多，导致求解较困难。在求解时，可根据实际情况进行适当简化以减少计算量。同时，在特定条件下，用准稳态模型来代替非稳态模型也是不错的选择。

3.2 对于非瓶装物料的数值模拟问题，由于侧面参与传热传质，导致升华界面不再是均匀退却，传统的升华界面移动速度公式已不再适用，需要加以改进。

3.3 在进行数值模拟的计算时，由于计算量较大，需要大量的编程，若运用专业CFD软件来进行计算，可以节省时间，提高效率，同时计算也更精确。不过，由于各种CFD软件专门针对冷冻干燥所建立的模型较少，在运用软件进行数值模拟时还是需要使用者对软件有较高的熟练程度，能够根据冷冻干燥的特点对软件自带模型进行修改和自定义。同时，由于所有模型都是进行多种假设后得出的，在进行数值计算的同时，也需要不断地进行实验验证，修正数学模型中的假设参数，以更精确地指导和控制冷冻干燥过程。

综上所述，近年来计算机及各种编程软件和计算软件的发展，使得对冷冻干燥进行更精确的模拟成为可能。结合各种专业计算软件对冷冻干燥建立更精确的非稳态模型并进行求解，或结合各种测量途径建立专门的简化模型将会是冷冻干燥数值模拟的一个发展方向。

[1]Sandall C O，King C J，Wilke C R.The relationship between transport properties and rate of freeze-drying of poultry meat[J]. AIChE Journal，1967，13：428-438.

[2]彭润玲，徐成海，张世伟，等.真空冷冻干燥过程传热传质理论研究的动态[J].真空与低温，2005，12（1)：51-61.

[3]谢国山，王立业.真空冷冻干燥过程的模拟研究[J].化工装备技术，2003，24（1)：6-9.

[4]谢国山.真空冷冻干燥的实验研究和数值模拟[D].大连：大连理工大学，2003.

[5]王博，李成华，赵青松.基于URIF模型的冷冻干燥过程二维描述[J].农机化研究，2006（3)：71-75.

[6]苏庆勇，贺广兴.基于URIF模型的圆柱型物料冷冻干燥过程的二维模拟及验证[J].食品科技，2008（1)：81-84.

[7]Liapis A I，Litchfield R J.Numerical solution of moving boundary transport problems in finite media by orthogonal collocation[J].Computers and Chemical Engineering，1979（3)：615-621.

[8]Litchfield R J，Liapis A I.An adsorption-sublimation model for freeze dryer[J].Chemical Engineering Science，1979，34（9)：1085-1090.

[9]Sheehan P，Liapis A I.Modeling of the primary and secondary drying stages of the freeze drying of pharmaceutical products in vials：numerical results obtained from the solution of a dynamic and spatially multidimensional lyophilization model for different operational policies[J].Biotechnology and Bioengineering，1998，60（6)：712-728.

[10]Zhai SL，Su HY，Taylor R，et al.Pure ice sublimation within vials in a laboratory lyophiliser;comparison of theory with experiment[J].ChemicalEngineeringScience，2005，60：1167-1176.

[11]Gan K H，Bruttini R，Crosser O K，et al.Freeze-drying of pharmaceuticals in vials on trays：Effects of drying chamber wall temperature and tray side on lyophilization performance[J]. International Journal of Heat and Mass Transfer，2005，48（9)：1675-1687.

[12]Velardi S A，Barresi A A.Development of simplified models for the freeze-drying process and investigation of the optimal operating conditions[J].Chemical Engineering Research and Design，2008，86（1)：9-22.

[13]彭润玲，刘长勇，徐成海，等.生物材料冻干过程分形多孔介质传热传质模拟[J].农业工程学报，2009，25（9)：318-322.

[14]Nama J H，Song C S.Numerical simulation of conjugate heat and mass transfer during multidimensional freeze drying of slab-shaped food products[J].International Journal of Heat and Mass Transfer，2007，50，4891-4900.

[15]姚志华，陈杰平.食品真空冷冻干燥准稳态模型的研究[J].安徽科技学院学报，2007，21（1)：37-40.

[16]Dyer D F，Sunderland J E.Heat and mass transfer mechanisms in sublimation dehydration[J].Journal of Heat Transfer，1968（11)：379-384.

[17]李云飞，王成芝.真空冷冻干燥中能量控制下的准稳态模型研究[J].华东工业大学学报，1997，19（1)：34-38.

[18]侯岩松，颉刚，张敏，等.软测量建模软件的设计与开发[J].电脑知识与技术，2011，7（4)：753-756.

[19]刘伟，范爱武，黄晓明.多孔介质传热传质理论与应用[M].北京：科学出版社，2006：26-27.

[20]罗瑞明，董平，李亚蕾，等.干切牛肉冷冻干燥中解析干燥过程的动态模拟及优化[J].农业工程学报，2009，25（7)：271-277.

[21]Crank J.Free and moving boundary problems[M].Oxford：Clarerdon Press，1984：163-281.

[22]刘中良，潘延龄.一维动边界问题的变时间步长方法[J].华东石油学院学报，1986，10（3)：28-39.

[23]贺广兴.平板状食品冷冻干燥过程的理论计算[J].桂林航天工业高等专科学校学报，2008（4)：58-60.

[24]杨卫波，施明恒，刘光远，等.基于显热容法的地源热泵地埋管换热器周围土壤冻结特性研究[J].暖通空调HV&AC（专题研讨)，2008，38（4)：6-10.

[25]叶宏，何汉峰，葛新石，等.利用焓法和有效热容法对定形相变材料熔解过程分析的比较研究[J].太阳能学报，2004，25（4)：488-491.

[26]李海梅，顾元宪，申长雨.平面相变热传导问题等效热容法的有限元解[J].大连理工大学学报，2000，40（1)：45-48.

[27]周业涛，关振群，顾元宪.求解相变传热问题的等效热容法[J].化工学报，2004，55（9)：1428-1433.

[28]谢晶，瞿晓华，施骏业.计算流体力学在食品工业中的应用[J].上海水产大学学报，2004，13（2)：170-175.

[29]刘妍玲，王世清，张岩，等.计算流体力学在食品加工贮藏中的应用[J].保鲜与加工，2008（2)：15-18.

Present status and prediction in numerical simulation of pabular and pharmaceutical freeze drying

ZHAN Bo，ZHOU Guo-yan，CAO Bin-hong，WANG Wei-yue
（Institute of Cryo-medicine and Food Refrigeration，Shanghai University of Science and Technology，Shanghai 200093，China)

Numerical simulation could predict and guide the pabular and pharmaceutical freeze-drying to improve the efficiency of the process.The present status of research in numerical simulation of freeze-drying was introduced and put emphasis on new improvements.Sum up the assumed conditions and relative merits of models.The method for each model was reviewed and summarized，and relative merits of each method were analyzed.The trend of numerical simulation of freeze-drying was predicted and pointed out that：establishing special model combined with varied measuring methods was a trend of numerical simulation of freeze-drying.

freeze-drying;numerical simulation;solution method;review

TS205.7

A

1002-0306（2012)07-0414-04

2011－04－25

詹博（1985-)，男，硕士研究生，研究方向：食品冷冻干燥。

国家自然科学基金青年基金项目（50206013);上海市重点学科建设项目（S30503);上海市教委科研创新项目（09YZ230)。