组合导向浮阀塔板多相流的数值模拟

张 朦，张海涛，张杰旭，张成芳

华东理工大学石油加工研究所，上海 200237

组合导向浮阀塔板多相流的数值模拟

张 朦，张海涛，张杰旭，张成芳

华东理工大学石油加工研究所，上海 200237

为了深入了解组合导向浮阀塔板的流体力学性能，在阀片全开状态下，利用计算流体力学（CFD）方法对700 mm中试规模的组合导向浮阀塔塔板上的气液两相三维流场进行模拟研究，采用欧拉-欧拉非稳态模型（Eulerian），考虑气液两相间的动量传递现象，将清液层高度作为衡量非稳态流场收敛的判别标准。通过关联实验数据，得出适合组合导向浮阀塔板的清液层高度关系式，并与模拟结果进行对比。模拟结果表明，组合导向浮阀塔板的液面梯度和弓形区域回流面积小，具有良好的流体力学性能。

计算流体力学 组合导向浮阀塔板 清液层高度

随着工业生产中对板效率要求的不断提高，在化工、石油化工和炼油等工业中广泛应用的F1型浮阀塔板呈现出各种弊端，使浮阀塔板的操作性能受到了限制。为了适应工业化的要求，路秀林等［1］发明的导向浮阀塔板具有优良的流体力学性能，如液面梯度小、板上液层分布薄而均匀和压降低等，因此得到了广泛应用［2］。王忠诚等［3］通过对导向浮阀塔板的改进，提出了导向作用更强的导向梯形浮阀塔板。张杰旭等［4］综合考虑了导向浮阀和导向梯形浮阀的优缺点，开发了组合导向浮阀塔板，通过两种浮阀的合理配比使操作性能更佳。

由于传统塔板结构设计中主要依据的实验方法和经验关联式逐渐呈现出各种弊端，新兴的计算流体力学（CFD）方法逐渐成为塔板设计开发中的重要工具。国内外学者采用仿真模拟的方法对此做了大量的研究：Yoshida［5］和Zhang等［6］分别对筛板建立了二维拟单相流模型，Yoshida忽略气相的影响，在层流状态下模拟液相流场，而Zhang等考虑了垂直的气相阻力作用。Krishna等［7，8］通过建立欧拉双流体三维模型模拟筛板上的流体力学行为，相间动量交换相采用Bennett［9］提出的液相体积分率关系式，先后对矩形和圆形筛板进行了模拟计算。Gesit等［10］用Bennett等［9，11］曳力系数在工业规模筛板塔的流体力学性能方面做了大量研究。Li等［12］模拟计算了F1型浮阀塔板的三维流动状态，但其采用矩形塔板，不能反映工业中常用的圆形塔器弓形区域的流动。崔觉剑［13］和Wang等［14］分别模拟计算了气相通过折边固定阀塔板和锯齿形导向浮阀塔板的单向流流场，并与实验结果进行了对比，结果表明模拟结果与实验结果相吻合。Zarei等［15］将MVG塔板与筛板的流体力学性能进行了比较，通过对两种塔板清液层高度、泡沫区高度和压降的比较说明MVG塔板流体力学性能的优越性。Jiang等［16］借鉴前人的研究方法，将CFD方法运用在对新型固阀塔板流场的研究中，并建立了适合三角固阀的气相分率式。

在公开的文献中对基础板型的CFD模拟研究较多，而对组合导向浮阀塔板多相流的仿真模拟还鲜有报道。对组合导向浮阀塔板流体力学性能的模拟结果有助于深入了解不同工况下塔内气液两相的流动和接触状态，这为工业生产中优化操作条件和改进塔器内部构件提供借鉴。因此，本工作利用FLUENT6.3.26对700 mm中试规模的组合导向浮阀塔板上气液两相复杂的三维流场进行数值模拟，模型中的气液两相采用空气-水为媒介。

1 气液两相流场的CFD模拟

通过建立欧拉-欧拉两相流三维模型对塔板上的流体力学性能进行研究，在欧拉多相流模型中，各相被看作是互相贯穿的连续介质。

1.1 塔板上气液两相流动的数学模型

1.1.1 控制方程

两相连续性方程（质量守恒方程）：

动量守恒方程（N-S）：

气液相体积分率之和为1，为了使方程组封闭，假设气液两相处在同一压力场下。湍流模型采用混合相κ-ε湍流模型来封闭欧拉模型。

1.1.2 塔板上气液两相动量传递过程分析

动量守恒方程是整个控制方程组的核心，既包括相内动量传递又包括相间的动量传递，式（3）和（4）中的动量交换源项MG L代表两相间的动量传递。气液两相之间的相互作用力主要包括曳力、虚拟质量力和升力等。Sokolichin等［17］曾证明虚拟质量力对结果影响不大，在分析中一般忽略不计；由于对分散气相确定升力系数比较困难，所以升力一般也忽略不计，而曳力是阻止气泡无限加速的主要作用力，伴随气泡的存在而存在，对塔板上气泡运动的影响也是最重要。为了简化模型，动量交换源项一般只考虑曳力的作用，其表达式如下：

其中，CD是相间动量交换系数或曳力系数。

由于气相是以气泡群的形式穿过液层的，所以采用Krishna等［18］提出的气泡群曳力模型：

其中，Vslip是气泡群相对于液相的滑移速率，一般由表观气速和气相平均分率进行估算，如式（7）：

将式（6）和（7）代入式（5）得到：

如果初始化时板上液相体积分率为零，板上只有气相存在，此时就不存在不同相之间的动量交换，而式（8）在此情况下并不符合实际情况，所以在式（8）的基础上加一项因子将其修正为式（9）。此因子的分子是局部液相分率，分母是平均液相分率，对此项在计算域内进行积分约等于1，也就是说此因子的加入对结果影响不大。Krishna［7］也曾表示在较宽的操作范围内，这种修正方法可以确保CFD模拟所得平均气相分率与实验数据保持吻合，说明这种修正是符合实际情况的。由于现阶段很难准确地预测每个气泡的大小，通过对比式（5）与（8）可以发现，在选用式（8）计算气液两相之间的相互作用时，不用考虑每个气泡的大小对曳力的影响，从而可提高计算的准确性。

气相穿过不同鼓泡元件的流体力学行为是不同的，平均气相分率也会不同。对于组合导向浮阀塔板参照前人［12，16］的做法，建立新的清液层关系式。清液层高度主要受表观气速、板上持液量和堰高的影响，其一般形式如下：

通过对组合导向浮阀塔板的实验数据［4］进行回归分析得到参数a1～a5的最佳值：a1为0.54，a2为 3.087 8，a3为130.1，a4为7.65，a5为0.6，从而得到组合导向浮阀塔板气相平均分率的表达式为：

1.2 几何模型及网格划分

1.2.1 几何模型及边界条件

图1是两种浮阀的结构及其参数，浮阀开度为8 mm，导向孔的高度为4 mm。图2为计算流动区域的几何模型及边界条件，塔板的浮阀排布方式采用工业中最常见的形式，在液相进口和弓形区域安排梯形浮阀。在建模时对浮阀的阀盖和阀脚以及塔板的厚度做了简化处理，由于计算区域的几何形状具有对称性，所以只计算半个塔盘。几何模型的结构参数：直径为700 mm，高度为300 mm，堰高为40 mm，降液管底隙高度为35 mm，塔板上浮阀的排布采用196 mm×110 mm的叉排方式。

图1　组合导向浮阀Fig.1 Combined guided floating valve

图2　几何模型及边界条件Fig.2 Geometric model and boundary

进口边界条件：

气液两相都采用“速度入口”边界条件，且在进口处速度分布均匀；气相进口只允许气相进入，液相进口只允许液相进入：

液相：

气相：

进口κ，ε边界条件采用文献［12，16］的经验关联式：液相：

气相：

出口边界条件：气液相出口取压力出口边界条件，液相出口只有液相流出，气相出口只有气相流出。

1.2.2 网格独立性验证

网格划分的合理与否对计算结果的正确性和计算成本的大小有重要影响，由于靠近塔板区的速度梯度较大，所以将计算域分成上下两个区进行网格划分，对靠近塔板区域进行加密处理。为了验证网格的独立性，采用5种网格划分方法，以清液层高度作为判别标准。图3是不同网格数目下清液层高度的变化趋势，可看出，随着网格数目的增多，模拟值与实验值越接近，当网格数目超过174 565时变化趋于平缓。考虑到网格数越多，计算成本越高，所以采用网格数为174 565的划分方式，靠近塔板区（0～30 mm）的网格尺寸采用5 mm，远离塔板的区域（30～300 mm）网格尺寸采用8 mm，图4是塔板网格划分示意图。

使用FLUENT 6.3.26在处理器为Intel（R） Core（TM） i7-4930K CPU @ 3.40GHz 3.70GHz的PC机上对模型进行计算。FLUENT是基于有限体积法对计算域进行离散化的，为保证计算的稳定性，采用一阶迎风格式对控制方程进行离散。选用非定常的压力基-隐式求解器，时间步长设为0.002 s，在默认的松弛因子下，选择SIMPLE算法对压力速度耦合方程求解，各控制方程中参数的收敛标准设为10-5。

图3　不同网格数目下的清液层高度模拟结果Fig.3 Clear liquid height under different grid counts

图4　塔板的网格划分Fig.4 The grid map of tray deck

2 模拟结果与讨论

2.1 模型验证

塔板上气液两相的流动是一个动态平衡过程，在模拟过程中实时监控塔板上的持液量，如果清液层高度在一定值附近小幅度上下波动，则可以认为计算收敛。清液层高度是板上液相平均分率和计算域高度的乘积，本工作中关于塔板的流体力学参数都是由计算平衡后一段时间内的平均值所得到。如图5在20 s后塔板上清液层高度基本不再有大幅度的改变，此时便可以认为计算达到了收敛。

图6为清液层高度的CFD计算值与实验拟合值的比较。由图6可以得出如下规律：（1）在相同的液流强度和气速下，清液层高度随堰高的增大而增大；（2）在相同的气速和堰高下，清液层高度随液流强度的增大而增大；（3）在相同的液流强度和堰高下，清液层高度随表观气速的增大而减小。通过将相同操作条件（气液流量）下CFD模拟得到的清液层高度与拟合式（10），（11）和（12）的计算值进行对比可以看出，模拟值与拟合式的计算值比较吻合，这说明通过对实验数据进行回归分析所建立的拟合式与CFD模型相匹配，适合用于CFD模拟计算，从而证明了所建CFD模型的正确性。

图5　清液层高度的瞬态监测Fig.5 Transient clear liquid height monitored as a function of time

图6　清液层高度的CFD计算值与实验拟合值的比较Fig.6 Comparison of hclof CFD simulations with correlation of this paper

2.2讨论分析

在LS为15 m3/（m·h），VS为0.6 m/s，hw为40 mm的条件下进行模拟，图7是y为93 mm截面的瞬态气相分率分布图，塔板上气液两相的流动是一个从非稳态到稳态的瞬态变化过程，在20 s以后气相分率分布基本趋于稳定，说明气液两相的流动进入动态平衡状态。

图7 不同时刻y为93 mm截面的瞬态模拟结果Fig.7 Snapshots of the front view for y at 93 mm

图8为塔板上不同高度的气相分率分布图，气相穿过塔板阀孔进入板上液层，夹带液滴向上运动，在向上运动的过程中气相体积分率逐渐变大，分散的气泡逐渐汇聚形成连续区域。从图（c）可以看出气相在穿过导向浮阀时的流动路径是从浮阀两侧和顶部导向孔进入板上液层。

图8　不同高度气相分率分布Fig.8 Gas hold-up profile at different elevations

泡沫区高度是塔板重要的流体力学参数，代表雾沫夹带程度的大小，根据Gesit等［10］的观点，一般将塔板上液相分率大于10%的气液两相接触区域定义为泡沫区。图9是x-y水平面的气相平均分率随高度的变化曲线，通过对比不同表观气速下的气相平均分率可以看出：当表观气速增大时更多液滴被夹带进入泡沫区，从而造成泡沫区高度的增加；在阀盖12 mm附近平均气相分率有一个下降的趋势。说明在阀盖上存在流动死区，与文献［16］中对三角固阀塔板的模拟相符。

图9　水平截面平均气相分率随高度的变化曲线Fig.9 Average gas hold-up at different elevations

图10 x为0截面液相平均水平速度沿y轴的变化曲线Fig.10 Average liquid velocity profile along y direction

图10是x为0截面的液相平均速度的水平分量从塔板中心到塔壁的变化曲线，可以看出从塔板中心到弓形区域速度整体逐渐变小，呈上下振荡形式。液体在流经导向浮阀时，受到从导向孔水平喷出的气体的推动作用，液速在此处达到峰值；而在相邻浮阀之间，由于受到气流的对冲作用，液速变小。

图11是不同水平面的气液相速度矢量图，可以看出，气体流经导向浮阀的流动路径以及气相对液相的作用。图（b）是浮阀附近的气相速度矢量放大图，可清晰看出气相从浮阀两侧和顶部导向孔喷出；图（c）是z为16 mm水平面的液相速度矢量图，由于受到气体的推动作用，沿气流方向的液速大幅增加。从导向孔流出的气体与液流方向一致，可推动液体向溢流堰方向流动，从而达到减小液面梯度、降低板压降的目的。

对于圆形塔器，由于弓形区域液相的流路较长而且流速较慢，容易出现回流现象。图（d）是 z 为16 mm截面上弓形区域局部矢量放大图，从图可以看出，组合导向浮阀塔板较传统塔板中的液相回流现象明显减弱。这是因为在液相进口和弓形区域安排了导向梯形浮阀，它的导向作用比矩形导向浮阀更强，从而加速了液体从进口到出口的流动，说明组合导向浮阀塔板具有良好的流体力学性能。

图11　不同横截面气液相速度矢量Fig.11 gas-liquid velocity profile at different elevations

由于浮阀的导向作用，导向浮阀塔板的鼓泡作用不像筛板样均匀，被气相带起的液滴整体趋势是向出口方向回落，还有少部分反向回落。从图12可以看出主要有三个漩涡：（1）液体进口上方存在回流现象，因为在液体进口上方存在低压区，在压强差的作用下会在竖直方向上发生回流；（2）阀盖附近存在漩涡，气体流经阀盖时受到阻挡而形成低压区，从而造成回流现象；（3）在溢流堰附近，溢流堰下的液层受到溢流堰的阻挡，形成回流来补充鼓泡区的液相空缺。

3 结 论

利用FLUENT6.3.26对组合导向浮阀塔板在阀片全开状态下进行仿真模拟。通过对实验数据的回归分析建立适合组合导向浮阀塔板的平均相体积分率关联式；将CFD模拟结果与所建立的相分率关系式进行对比，验证了所建CFD模型的正确性；通过改变操作条件探究了堰高、表观气速和液流强度对清液层高度的影响，表明组合导向浮阀塔板具有良好的流体力学性能，液面梯度小、板上流场分布较均匀和弓形区域流动死区小等。

符号说明

［1]路秀林， 赵 培. 导向浮阀塔板［J］. 化工装备技术， 1992， 13（1）: 1-5. Lu Xiulin， Zhao Pei. The directed floating valve tray［J］. Chemical Equipment Technology， 1992， 13（1）: 1-5.

［2]李成才. 导向矩形浮阀塔板在分离塔上的应用［J］. 科学之友， 2006， （4）: 14-15. Li Chengcai. Application of guide rectangular floating valve plate to methyl acetate separate tower［J］. Friend of Science Amateurs， 2006，（4）: 14-15.

［3]王忠诚， 姜 斌， 张雅芝， 等. 导向梯形浮阀塔板的特点及工业应用［J］. 石油化工， 1998， 27（4）: 271-275. Whang Zhongcheng， Jiang Bin， Zhang Yazhi et al. Characteristics of guided trapezoid valve tray and its industrial applications［J］. Petrochemical Technology， 1998， 27（4）: 271-275.

［4]张杰旭， 赵 培， 李玉安， 等. 组合导向浮阀塔板的开发研究［J］. 化工进展， 2005， 24（增）: 69-72. Zhang Jiexu， Zhao Pei， Li Yu'an， et al. Development and study of combination directed valve trays［J］. Chemical Industry and Engineering Progress， 2005， 24（sup）: 69-72.

［5]Yoshida H. Liquid flow over distillation column plates［J］. Chem Eng Comm， 1987， 51: 261- 275.

［6]Zhang M Q， Yu G C. Simulation of two dimensional liquid phase flow on a distillation tray［J］. Chinese Journal of Chemical Engineering，1994， 2（2）: 63-71.

［7]Krishna R. CFD simulations of sieve tray hydrodynamics［J］. Chem Eng Res， Design， 1999， 77（7）: 639-646.

［8]van Baten J M， Krishna R. Modeling sieve tray hydraulics using computational dynamics［J］. Chem Eng Journal， 2000， 77（3）: 143-151.

［9]Bennett D L. New pressure drop correlation for sieve tray distillation columns［J］. AIChE Journal， 1983， 29（3）: 434-442.

［10]Gesit G， Nandakumar K， Chuang K T. CFD modeling of flow patterns and hydraulics of commercial-scale sieve trays［J］. AIChE Journal，2003， 49（4）: 910-924.

［11]Colwell C J. Clear liquid height and froth density on sieve trays［J］. Ind Eng Chem Pro Des Dev， 1981， 20（2）: 298-307.

［12]Li Xingang， Liu Dexin， Xu Shimin. CFD simulation of hydrodynamics of valve tray［J］. Chem Eng Pro， 2009， 48（1）: 145-151.

［13]崔觉剑. 折边固定阀塔板流场的实验研究和数值模拟［D］. 杭州: 浙江工业大学， 2007.

［14]Wang Lianghua， Cui Juejian， Yao Kejian. Numerical simulation and analysis of gas flow field in serrated valve column［J］. Chinese Journal of Chemical Engineering， 2008， 16（4）: 541-546.

［15]Zarei T， Rahbar R， Zivdar M. Computational fluid dynamic simulation of MVG tray hydraulics［J］. Korean Journal of Chem Eng， 2009，26（5）: 1213-1219.

［16]Jiang S， Gao H， Sun J S， et al. Modeling fixed triangular valve tray hydraulics using computational fluid dynamics［J］. Chem Eng Prog，2012， 2（52）: 74-84.

［17]Sokolichin A， Eigenberber G. Gas-liquid flow in bubble columns and loop reactors: part I. detailed modeling and numerical simulation［J］. Chem Eng Sci， 1994， 49（24）: 5735-5746.

［18]Krishna R， Urseanu M I. Rise velocity of a swarm of large gas bubbles in liquids［J］. Chem Eng Sci， 1999， 54（2）: 171-183.

Numerical Simulation of Multiphase Flow in Combined Guided Valve Tray of Distillation Column

Zhang Meng， Zhang Haitao， Zhang Jiexu， Zhang Chengfang
Institute of Chemical Engineering， East China University of Science and Technology， Shanghai 200237， China

A three dimensional computational fluid dynamics （CFD） model was developed to predict the multiphase flow profile of the combined guided valve tray under the full open condition. Eulerian-Eulerian framework was chosen to solve the equations of both gas and liquid phases. The clear liquid height was monitored as the numerical convergence standard of unstable flow field. The inter-phase momentum exchanging term was based on a new correlation for the clear liquid height， which was obtained from experimental data. By compared with the experimental counterparts， the simulation results were verified and more simulations were carried out for various operational conditions. The results showed that the gradient of liquid level and reflux area were small， and the fluid mechanics performance was good.

computational fluid dynamics； combined guided valve tray； clear liquid height

TQ021.1；TQ018

A

1001—7631 （ 2015 ） 02—0106—09

2014-11-18；

2015-01-26。

张 朦（1989—），女，硕士研究生；张海涛（1970—），男，副教授，通讯联系人。E-mail: zht@ecust.edu.cn。