分子蒸馏技术是一种在高真空条件下利用不同物质分子运动平均自由程差异来实现物料分离的新型液-液分离技术[1, 2]。该技术具有蒸馏温度低、物质受热时间短、分离程度高等特点，有利于高沸点、热敏性和高黏度物料分离，现广泛应用于食品、医药、油脂加工、石油化工等行业[3~5]。分子蒸馏设备大致可分为静止式、降膜式、刮膜式及离心式四种[6]。其中，刮膜式分子蒸馏器是在降膜式分子蒸馏器基础上引入刮膜器，刮擦壁面上的液膜形成均匀厚度薄膜，增加传质传热效率，进而提高分离效率。但同时因为刮膜器的作用，液膜内流动、传质传热机理也更加复杂。近年来，国内外学者对刮膜式分子蒸馏设备中液膜的传质传热进行了大量的理论模型研究。Komori和Takata[7, 8]从理论和实验两方面研究了头波对刮膜式分子蒸馏过程中分离效率的影响。Micov[9]和 Nguyen[10]则分别对刮膜式分子蒸馏过程中的分离因数和液膜组成的变化进行了研究，并建立一维的数学模型。Nguyen等[10]则在完全理想状态下建立一维等温蒸馏模型，忽略传质阻力和再蒸发，确定了馏出液浓度、残留液浓度、进料速率、蒸发面积以及收率间的相互关系。李国兵等[11]忽略曲面效应，建立了传质传热数学模型，并用Rung-Kutta方法联合平均隐式差分法对其进行数值求解，并在此基础上对邻苯二甲酸二丁酯(DBP)-癸二酸二丁酯(DBS)混合物分子蒸馏过程进行数值模拟[12]，得到温度和浓度分布在三维上的变化规律，同时考察了温度、浓度、分离因数等的影响因素。李慧和王珂鑫[13]则依据进料量、液体浓度、刮膜电机转速以及蒸发器尺寸等主要因素之间的结构关系，建立刮膜式分子蒸馏的传质模型，并通过MATLAB 进行数值仿真，确定每个参数对蒸馏过程的影响。目前已有的数学模型虽然能为分子蒸馏装置的设计和理论研究提供参考，但不能形象直观地看到刮膜式分子蒸馏器内液膜厚度、温度和浓度的变化。随着计算流体力学的发展，采用数值法计算可为液膜传热传质过程的研究提供可靠的流场数据。许松林等[14, 15]利用CFD软件分析了蒸发器上液体负荷、刮膜器旋转速度对蒸发器上液膜厚度变化及液膜中速度矢量的影响。
本文在已有的刮膜式分子蒸馏流体力学的基础上，采用FLUENT软件对DBP-DBS混合物的刮膜分子蒸馏过程进行数值模拟，得到分子蒸馏器内温度和物料浓度变化，为进一步研究分子蒸馏器内传质传热研究提供参考。
2 模型与方法 2.1 物理模型及网格划分在刮膜式分子蒸馏设备中，物料从进料口滴到刮膜器顶部固定圆盘，从刮膜器顶圆盘周边与蒸发壁面之间的缝隙流下，然后被刮板刮动，均匀分布在蒸发壁面上。图 1为简化的刮膜式分子蒸馏三维计算模型，蒸发器外径是100 mm，内径70 mm，液相从入口处的圆环进料，内外两环相距1 mm，内部2个直径为8 mm转子处于同一平面，呈180度分布，转子距蒸发壁面0.3 mm，蒸发壁面高200 mm。为使模拟结果显示的更精确，在前处理网格划分时，尽量使用结构化网格。先对径向切面进行网格划分，如图 2所示，使用边界层细化网格划分外圆环，submap划分切除了两个转子圆面的内圆环面，最后用cooper法对整个三维结构进行网格划分。
|
图 1 刮膜式分子蒸馏器结构示意图 Fig.1 Sheme of a wiped film molecular distillator |
|
图 2 网格划分和转子附近入口局部网格划分 Fig.2 Grid scheme for flowing area and partial grid scheme for the inlet near the rotor |
控制方程中物理性质由系统各个相共同决定，是体积分数加权平均值，VOF多相流模型追踪气液相界面。在VOF模型中，由气态DBP和DBS组成的气相为主相，液态DBP和DBS组成的液相为第二相。纯组分物性参数见文献[16, 17]，混合物物性参数计算公式见文献[18, 19]。液相轻重组分黏度混合规则按质量比重混合，气相轻重组分黏度混合规则按理想气体混合规则。组分方程中的质量源相和能量方程的计算与文献[20, 21]一致，只考虑蒸馏过程中的蒸发现象，不考虑冷凝过程和冷凝面的反蒸发现象，忽略分子蒸馏器内气相分子对液膜表面分子蒸发速率的影响，即蒸发速率方程可用Langmuir-Knudsen方程计算。结合研究的二元物系DBP-DBS的控制方程如下：
连续性方程：
| $\frac{\partial }{\partial t}\left( \rho \right)+\nabla \cdot \left( \rho \overrightarrow{u} \right)=0$ | (1) |
动量方程：
| $\frac{\partial }{\partial t}\left( \rho \overrightarrow{u} \right)+\nabla \cdot \left( \rho \overrightarrow{u}\overrightarrow{v} \right)=-\nabla \rho +\nabla \cdot \left[ \mu \left( \nabla \overrightarrow{u}+{{\overrightarrow{u}}^{T}} \right) \right]+\rho \overrightarrow{g}+\overrightarrow{F}$ | (2) |
体积分数方程：
| $\frac{\partial {{\alpha }_{q}}}{\partial t}+\overrightarrow{u}\cdot \nabla {{\alpha }_{q}}=0$ | (3) |
物质属性方程：
| $\mu ={{\alpha }_{l}}{{\mu }_{l}}+(1-{{\alpha }_{l}}){{\mu }_{q}}$ | (4) |
组分传输方程：
| $\frac{\partial }{\partial t}\left( {{\rho }_{q}}{{\alpha }_{q}}{{y}_{1}} \right)+\nabla \cdot \left( {{\rho }_{q}}{{\alpha }_{q}}\overrightarrow{u}{{y}_{1}} \right)={{S}_{lq}}+\nabla \cdot \left( {{\rho }_{q}}{{\alpha }_{q}}D_{1}^{q}\nabla {{y}_{1}} \right)$ | (5) |
能量方程：
| $\frac{\partial }{\partial t}\left( \rho E \right)+\nabla \cdot \left( \overrightarrow{u}\left( \rho E+P \right) \right)=\nabla \cdot \left( {{k}_{eff}}\nabla T \right)+{{S}_{h}}$ | (6) |
此外，由于刮板的存在，壁面液体出现周期性的湍流，因此选择RNG k-ε湍流模型。
2.3 边界条件及数值求解为简化计算，可设初始蒸馏室内充满DBS蒸汽[20]，液相进料混合物中DBP的质量分数为0.5，入口和出口边界条件分别设置为速度入口、压力出口。在旋转坐标系下设置混合流体，其旋转速度与转子转速相同。蒸发面与冷凝面均设为移动墙，其相对速度设为坐标旋转速度的相反值；将转子设为固定墙，相对速度为零，右手定则确定转子的旋转方向。选择增强壁面方式处理近壁面区域的流动，同时考虑重力。采用隐式有限体积法求解控制方程，使用PISO方法求解压力速度耦合项。连续性方程、动量方程、湍动能方程和湍动能耗散方程均采用二阶迎风格式离散计算。在求解方程时，先计算流动方程和未加质量源相的传质方程，待液膜完全稳定后，再加入质量源相和能量方程，求解所有方程。
3 结果和讨论 3.1 液膜温度分布从转子头波前方开始，监测液膜表面上温度、浓度等数据。以距离转子中心的角度大小来说明距离头波的远近。角度越大，距离头波越远。且由于两个转子的对称性，只需要监测至180度。由图 3(a)可知，在流体流动方向，液膜表面温度沿轴向逐渐增大。且距离进料口越远的平面，其液膜表面温度越高，这是因为液膜在转子刮擦和重力的双重作用下，液膜厚度沿轴向会逐渐变薄，在相同的加热壁面作用下，其液膜表面温度也就越高。而距离进料口较近的平面，因不断有新物料流下更新液膜，其液膜表面温度较低。在同一平面上，因液膜暴露时间不一样，距离转子前头波越近的液膜其表面温度越高，温度变化也比较剧烈，距离头波越远，液膜受热蒸发的时间越久，因表面蒸发带走的热量越多，其温度也越低。从图 3(b)看出，由于转子的转动，液膜温度沿径向变化保持一致。在蒸发壁面附近因存在温度边界层，温度下降的较快，有较大的温度梯度，之后热量在热传导及对流作用下，从壁面传到液膜表面，温度逐渐降低，产生沿径向变化的温度梯度。在液膜蒸发表面因物料蒸发带走热量，也会产生较大的温度梯度，导致温度下降。本文数值模拟的温度沿三维方向的变化与文献[12]结论一致。
|
图 3 液膜温度分布 Fig.3 Liquid film temperature distributions (a) surface temperature distribution of liquid films along circumferential direction at different heights (b)liquid film temperature distribution at z = 0.1m at different circumferential positions |
由图 4(a)可知，在转子刮擦作用下，液膜上易挥发组分开始蒸发，易挥发组分浓度沿轴向逐渐减小。z=0.05m平面因距离进料口较远，液膜较薄，易挥发组分浓度降低较多。而在同一平面上，距离头波越远的液膜，其在刮膜器刮擦后暴露的时间越长，表面易挥发组分的浓度就越低。故易挥发组分浓度沿周向是逐渐减小的。从图 4(b)可知，浓度边界层主要存在液膜表面上，而在靠近蒸发壁面处易挥发组分浓度则基本保持不变。但随着扩散传质的进行，浓度边界层逐渐变厚，会扩散到液膜主体，故易挥发组分浓度沿径向逐渐减小，产生径向浓度梯度。此外，靠近头波处的液膜因液膜更新更快，其易挥发组分浓度减少的较少。而远离头波处的液膜，因液膜较薄且蒸发时间长，其浓度减少量更大。
|
图 4 液膜浓度分布 Fig.4 Liquid film concentration distributions (a) DBP surface concentration distribution along circumferential direction at different heights (b) liquid film DBP concentration distribution at z = 0.1m at different circumferential positions |
由Langmuir-Knudsen方程[22]可知，物料的蒸发速率是温度和易挥发组分浓度的函数。因液膜表面温度及易挥发组分浓度沿轴向、周向及径向呈三维变化，因此蒸发速率也呈三维变化。由图 5(a)可知，易挥发组分的蒸发速率沿蒸馏器轴向逐渐增大，距离进料口越远的平面，其蒸发速率越大。虽然易挥发组分浓度沿轴向是逐渐降低，但液膜表面温度沿轴向是逐渐增大，且温度增大的影响超过了浓度降低的影响。因此从整体上看来，易挥发组分的蒸发速率沿轴向是逐渐增大的。而在同一平面上，距离头波越远的液膜表面，其蒸发速率先迅速降低，之后降低幅度趋于平缓。这是因为液膜表面温度和易挥发组分浓度均沿周向是逐渐降低的，在距离头波越近的液膜，其降低幅度也越大，之后降低幅度趋于平缓，故蒸发速率沿周向变化也呈现相同趋势。从图 5(b)可以看出，易挥发组分的蒸发速率沿径向也是逐渐降低的，在蒸发壁面附近降低的幅度最大，这是因为蒸发壁面附近存在较大的温度梯度，且温度和易挥发组分浓度沿径向也是逐渐降低的。
|
图 5 蒸发速率分布 Fig.5 Liquid film evaporation rate distributions (a) DBP surface evaporation rate distribution along circumferential direction at different heights (b)liquid film DBP evaporation rate distribution of at z = 0.1m at different circumferential positions |
壁面加热温度固定时，刮膜设备中蒸馏物料进料量是影响液膜流动状况的主要因素。进料量较小时，转子刮擦作用下易出现干板，液膜分布不均出现局部过热。当进料量较大，转子来不及刮擦，物料就在重力作用下流出，不利于物料受热蒸发分离。因此固定其他工艺条件不变，研究进料量对液膜传质传热影响。由图 6(a)、(b)、(c)可知，同一平面上进料量越大，其液膜表面温度越低，沿周向变化幅度也越小；进料量越大，因液膜更新太快，蒸发的也就越慢，易挥发组分浓度降低的也越少，蒸发速率也就越小。在进料量为3.6 kg⋅h-1和9.0 kg⋅h-1时，液膜表面温度和蒸发速率沿周向是逐渐降低的，易挥发组分浓度降低幅度也较大；当进料量增大到12.7 kg⋅h-1和18 kg⋅h-1时，头波附近温度较低，且因进料量过大，导致停留时间过短，物料受热时间短，同一平面上温度都较小，易挥发组分浓度减少量也较小。另外，按文献[21]上方法，获取z=0.1 m平面上平均液相总传质系数，得到了不同进料量下液相总传质系数的变化规律(图 6(d))。进料量越大，液膜就越厚，大部分热量就用于加热液膜，而只有少部分热量用于液膜表面蒸发，因此易挥发组分减少的也越少。故进料量的增加增大了传质阻力，不利于传质。因此适当减少进料量，增加物料受热时间，能够提高物料蒸发速率，利于提高蒸馏效率。
|
图 6 进料量对z = 0.1 m平面液膜表面温度、浓度、蒸发速率和液相总传质系数的影响 Fig.6 Effects of flow rate on surface temperature, concentration, evaporation rate and total mass transfer coefficient of z = 0.1 m (a) surface temperature (b)concentration of DBP (c) evaporation rate (d) the total mass transfer coefficient |
固定其他工艺参数，研究壁面加热温度对液膜传质传热影响。由图 7(a)、(b)、 (c)可知，在同一平面，液膜表面温度随加热温度增大而增大，且液膜温度沿周向变化呈相同趋势，都是逐渐减小；易挥发组分浓度沿周向也是逐渐减小，且加热温度越高，减少的就越多，但当温度由393 K增加到413 K时，浓度减少的幅度也逐渐变小；物料的蒸发速率也随加热温度的增大而显著增大。由图 7(d)可知，液相总传质系数也随壁面加热温度的增大而增大。由此可知，增大加热温度能显著增加蒸发效率，但在现实生产中，较大的加热温度不仅消耗能源，增加生产成本，且过高的温度也会增加热敏性物料热分解的危险，因此应合理选择加热温度。
|
图 7 加热温度对z = 0.1 m平面液膜表面温度、浓度、蒸发速率和液相总传质系数的影响 Fig.7 Effects of heat temperature on surface temperature, concentration, evaporation rate and total mass transfer coefficient of z = 0.1m (a) surface temperature (b)concentration of DBP (c) evaporation rate (d) the total mass transfer coefficient |
本文通过CFD对刮膜式分子蒸馏器内传质传热过程进行数值模拟，得到液膜表面温度、易挥发组分浓度、蒸发速率等参数沿三维方向变化，并研究了进料量和壁面加热温度对其影响，得到以下结论：(1) 液膜表面温度沿周向逐渐降低，距离头波越远的液膜，其表面温度越低，沿液体流动方向逐渐增加，沿径向逐渐减少，在靠近壁面和蒸发表面温度梯度较大；
(2) 易挥发组分蒸发速率沿三维方向变化与温度相同；易挥发组分浓度沿周向和轴向均逐渐减小，在靠近壁面处浓度几乎保持不变，但总体上沿径向呈减小趋势；
(3) 液膜表面温度和易挥发组分的蒸发速率随进料量的增加而降低，随加热温度的增加而增加；易挥发组分浓度随进料量的减少而增加，随加热温度的增加而减少；液相总传质系数随进料量的增加而减少，随加热温度的增大而增大。
符号说明：
| D1 | — 扩散系数，m2⋅s-1 | Ts | — 液膜表面温度，K |
| E | — 总能量，J | TW | — 加热温度，K |
| F | — 外加体积力，N | u | — 流体速度，m⋅s-1 |
| G | — 易挥发组分蒸发速率，mol⋅m-2⋅s-1 | xDBP | — 易挥发组分质量分率 |
| g | — 重力加速度，m⋅s-2 | y1 | — 易挥发组分摩尔分率 |
| K | — 液相总传质系数，kg⋅m-3⋅s-1 | a | — 体积分率 |
| keff | — 传热系数，W×m2×K-1 | δ/δ0 | — 无因次液膜厚度 |
| P | — 压力，Pa | θ | — 距离转子中心的角度，° |
| Sh | — 能量源相，J×kg-1⋅mol-1 | μ | — 流体动力黏度，mPa×s |
| Slq | — 传质源相，kg⋅m-1⋅s-1 | ρ | — 流体密度，kg⋅m-3 |
| T | — 温度，K | l, q | — 下标液相，气相 |
| [1] | HANG Zhi-jun(杭智军), YING An-guo(应安国), WU Hao(武浩) . The progress of molecular distillation technology application and research in the pharmaceutical industry[J]. Chemical Production and Technology(化工生产与技术) , 2012, 19 (2) : 37-40 |
| [2] | Lutisan J, Cvengros J . Mean free path of molecular on molecular distillation[J]. The Chemical Engineering Journal and the Biochemical Engineering Journal , 1995, 56 (2) : 39-50 DOI:10.1016/0923-0467(94)02857-7 |
| [3] | Bhandarkar M, Ferron J R . Transport processes in thin liquid films during high-vacuum distillation[J]. Industrial & Engineering Chemistry Research , 1988, 27 (6) : 1016-1024 |
| [4] | Komori S, Takata K, Murakami Y . Flow structure and mixing mechanism in an agitated thin-film evaporator[J]. Journal of Chemical Engineering of Japan , 1988, 21 (6) : 639-644 DOI:10.1252/jcej.21.639 |
| [5] | Komori S, Takata K, Murakami Y . Flow and mixing characteristics in an agitated thin film evaporator with vertically aligned blades[J]. Journal of Chemical Engineering of Japan , 1989, 22 (4) : 346-350 DOI:10.1252/jcej.22.346 |
| [6] | MicovM , Lutisan J, Cvengros J . Balance equations for molecular distillation[J]. Separation Science and Technology , 1997, 32 (18) : 3051-3066 DOI:10.1080/01496399708000795 |
| [7] | Nguyen A D, Goffic F L . Limits of wiped film short-path distiller[J]. Chemical Engineering Science , 1997, 52 (16) : 2661-2666 DOI:10.1016/S0009-2509(97)00089-4 |
| [8] | LI Guo-bing(李国兵), XU Chun-jian(许春建), ZHOU Ming(周明) . Mathematic Model of Mass Transfer and Heat Transfer of Wiped Film Molecular Distillation(刮膜分子蒸馏传质和传热的数学模型)[J]. Journal of Tianjin University(天津大学学报) , 2003, 36 (3) : 261-266 |
| [9] | LI Guo-bing(李国兵), XU Chun-jian(许春建), ZHOU Ming(周明) . Numerical Simulation on Wiped Film Distillation of DBP-DBS(DBP-DBS刮膜分子蒸馏过程数值模拟)[J]. Journal of Tianjin University(天津大学学报) , 2004, 37 (3) : 1026-1034 |
| [10] | LI Hui(李慧), WANG Ke-xin(王珂鑫) . Mass transfer model and simulation of wiped film molecular distillation(刮膜式分子蒸馏传质模型及其仿真)[J]. Journal of Chemical Industry and Engineering(化工学报) , 2015, 66 (3) : 196-202 |
| [11] | ZHANG Ping(张平), LI Ya-xia(李雅侠), QIN Ran(秦然) . Three-dimension numeric simulation of liquid flow on oblique sieve tray(倾斜塔板液相三维流场的数值模拟)[J]. Journal of Chemical Engineering of Chinese Universities(高校化学工程学报) , 2014, 28 (6) : 1217-1222 |
| [12] | GU fang(谷芳), LIU Chun-jiang(刘春江), YU Li-ming(余黎明) . The CFD simulation of mass-transfer process in falling film with countercurrent two-phase flow(气-液两相降膜流动及传质过程的CFD研究)[J]. Journal of Chemical Engineering of Chinese Universities(高校化学工程学报) , 2005, 19 (4) : 438-444 |
| [13] | WANG Yan-fei(王燕飞), XU Song-lin(许松林) . Simulation of computational fluid dynamics for molecular distillation process[J]. Chemical Engineering (China)(化学工程) , 2010, 38 (1) : 30-33 |
| [14] | WANG Xiao-ming(王小明), XU Song-lin(许松林), XU Chang(徐畅) . Study on the effect of wiper structure on the liquid film flow condition in wiped film molecular distillation apparatus[J]. Journal of Chemical Engineering of Chinese Universities(高校化学工程学报) , 2012, 26 (1) : 13-18 |
| [15] | YAO Wen-jin(姚文瑾), LI Yang(李阳), XU Song-lin(许松林) . Residence time distribution of liquid in the wiped film molecular distillator and the wall structure optimization of the wiped film molecular distillator(刮膜式分子蒸馏器液体停留时间分布及壁面结构优化)[J]. Journal of Chemical Engineering of Chinese Universities(高校化学工程学报) , 2015, 29 (2) : 278-283 |
| [16] | Perry E S, Weber W H . Vapor pressures of Phlegmatic liquids(П).High weight esters and silicone oils[J]. Journal of the American Chemical Society , 1949, 74 (11) : 3729-3730 |
| [17] | Yaws C L . Chemical properties handbook[M]. New York: McGraw-Hill, 1998 . |
| [18] | Kawala Z, Stephan K . Evaporation rate and separation factor of molecular distillation[J]. Chemical Engineering & Technology , 1989, 12 (3) : 406-413 |
| [19] | Reid R C, Prausnitz J M, Poling B E . The properties of gases and liquids[M]. New York: McGraw-Hill, 1987 . |
| [20] | Lijun C, Haibo D, Aiwu Z . Hydrodynamics and transport processes during centrifugal short-path distillation[J]. Chemical Engineering & Technology , 2013, 36 (5) : 851-862 |
| [21] | Jiang Y, Lijun C, XigangY . CFD simulation of the heat and mass transfer process during centrifugal short-path distillation[J]. Chemical Engineering & Technology , 2015, 38 (1) : 173-180 |
| [22] | Kaplon J, Kawala Z, Skoczylas A . Evaporation rate of a liquid from the surface of a rotating disc in high vacuum[J]. Chemical Engineering Science , 1986, 41 (3) : 519-522 DOI:10.1016/0009-2509(86)87033-6 |