焦云鹏, 周晓庆, 陈建华丨基于OpenFOAM的精馏塔内气液两相流传热传质模拟

文摘 2024-05-17 13:52 北京

【点击上方蓝字关注我们】

基于OpenFOAM的精馏塔内气液两相流传热传质模拟

焦云鹏 ^1,² 周晓庆 ^1,³陈建华 ¹

1. 中国科学院过程工程研究所多相复杂系统国家重点实验室，北京 100190
2. 中国科学院大学化学工程学院，北京 100049
3. 沈阳化工大学化学工程学院，辽宁沈阳 110142

DOI：10.12034/j.issn.1009-606X.223180

摘要对于精馏塔内的多相传质过程，计算流体力学(Computational Fluid Dynamics, CFD)不仅可用于模拟流动现象，而且可以用于研究传热传质过程及其与流动的相互作用规律。基于之前精馏塔板气液流动模拟的研究工作，本研究应用OpenFOAM平台的多相流求解器，考虑了能量方程与组分输运方程，构建了精馏体系的传热传质模型，分别对环己烷-正庚烷理想体系与乙醇-水非理想体系进行了模拟，并分析了塔板上气液两相流动、组分浓度及温度等参数的分布规律。对于理想体系，理想溶液传质模型能够较准确地预测多层塔板上的温度及浓度场；对于非理想体系，需要在理想溶液传质模型的基础上引入活度系数模型，为此对比了UNIQUAC与NRTL两种活度系数模型的效果。在当前模拟框架下，引入活度系数模型有助于提升非理想体系温度场与浓度场的模拟精度，两种活度系数模型预测的整体趋势与文献结果基本一致，但在数值上UNIQUAC模型与文献结果吻合更好。此外，气液两相流场与浓度场分布的对比分析表明，液相的循环流动能够提升塔板的局部传质效率，导致液相入口与塔板堰处的效率高于塔板中心，但是循环区域内气液两相无法及时更新，将导致塔板整体传质效率下降。以上研究有望用于精馏塔板的设计与优化，为提高精馏塔传质效率、降低精馏能耗提供指导，对于其他气液两相流的CFD传热传质模拟也具有参考价值。

关键词气液两相流；精馏；热态模拟；传热与传质；OpenFOAM

1 前言

精馏过程是一种常见的物质分离和纯化方法，广泛应用于石油化工、制药、食品等工业领域^[1,2]。在精馏过程中，传热和传质性能是影响分离效率和产品质量的关键因素。传统的精馏塔设计主要基于经验和实验数据拟合的经验关联式。例如，Bell等^[3]利用实验对矩形塔板的四种流速分布，即反抛物线、直线、抛物线以及均匀分布进行了板效率的计算，结果表明流速分布不均匀会降低塔板效率。Porter等^[4]将塔板上的流体区域进行划分，并在两个区域建立各自的传质方程。利用实验研究可以对精馏塔进行相对可靠的测试，为研究人员提供基础数据，进而开发与验证计算模型，以此优化和改进精馏塔的设计和操作；但是实验方法存在获取数据困难、经验关联式适用性差以及成本高等问题。

近年来，计算流体力学(Computational Fluid Dynamics, CFD)获得快速发展与应用，数值模拟可以提供详细的流场分布和物理参数，包括各种场的动态演化数据，已经成为揭示塔内流体流动、传热与传质等物理过程、优化设计塔板结构的有力工具。王晓玲^[5]考虑了湍流对传质的影响，利用商业CFD软件Star-CD模拟了全回流条件下10层塔板的环己烷-正庚烷分离，模拟结果与Yanagi等^[6]的实验结果匹配较好。郭舒^[7]利用商业软件Fluent对单气泡上升运动进行实验和数值模拟研究，探究了气泡对传质的影响。赵丹^[8]建立了多溢流塔板上液相流动和传质的计算模型，利用Fluent模拟并分析了塔板上液相流动结构以及气液相混合状态对塔板效率的影响。栗秀萍等^[9]搭建了超重力减压精馏装置，将乙醇-水作为模拟物系，研究了间歇全回流操作对传质过程的影响。陈光辉等^[10]利用Fluent分别对新型鼓泡破泡一体化精馏塔盘和筛板的气液流动进行数值模拟，模拟结果表明，破泡装置能够提高传质效率。Park等^[11]利用Fluent模拟直接接触膜蒸馏系统的传热和传质，并通过实验进行验证，蒸馏水产量的预测结果与实验值吻合较好。Zhao等^[12]利用Fluent研究了气相速度对真空膜蒸馏真空侧水蒸气质量分数的影响，结果表明气流可以实现更大的水蒸气通量，从而有效增强传质。Tan等^[13]将三维欧拉多孔介质模型与传质模型耦合，考虑了冷凝相变和传质过程，利用Fluent对大型空分装置随机填料空气冷却塔中的传质和传热进行研究。

随着CFD开源软件不断发展，以OpenFOAM为主的开源平台凭借着透明化、可自定义和低成本的特点逐渐成为了研究气液相间传热传质问题的有效工具。吕志超^[14]建立了新的两相界面传质模型，并基于OpenFOAM对单液滴传质展开数值模拟，研究了分散相尺寸、连续相表面张力以及黏度对传质过程的影响。Godino等^[15]基于OpenFOAM平台中实现了一个针对气液沸腾传热传质体系的新型求解器，该求解器将共轭传热模型植入双流体框架，并在数值验证中表现出较高的准确性，有望用于求解核工程中液体系统的热耦合问题。Lee等^[16]在reactingTwoPhaseEulerFoam求解器的基础上结合Maxwell方程，设计并植入了新的界面热和质量传递源项，开发了一个新的多物理场模型，成功模拟了微波加热二元液体混合物中气液相互作用、气泡形成机制与气液相平衡等过程。Zhao等^[17]基于OpenFOAM建立了一个通用的两相流和传质计算流体动力学模型，研究了多个界面力模型和多尺度气泡模型对鼓泡塔气液流动与传质过程模拟的影响，并通过可视化实验验证其精度和有效性。Khan等^[18]利用OpenFOAM中的双流体模型研究了不同边界条件下鼓泡塔内的传质过程，包括吸收和解吸以及并流和逆流过程，模拟结果与经典实验数据呈现较好的一致性。

相较气液流体动力学模拟，包含传热与传质的气液两相流建模和仿真研究仍然较少，特别是利用OpenFOAM对气液精馏体系传质传热的模拟研究。为了理解精馏塔内气液两相传热传质特性，本研究在前期工作^[19]的基础上，采用双流体模型研究精馏传热与传质过程，耦合能量方程与组分输运方程，构建完整的气液两相热态模拟方法。首先对环己烷-正庚烷理想体系进行模拟验证；其次，考虑到实际生产中大多数物系为非理想体系，并且组分间额外的相互作用导致其模型复杂度更高，因此将验证后的理想体系热模型进一步耦合活度系数模型，对文献[20]中的乙醇-水非理想体系进行模拟，研究了曳力模型与活度系数模型等对气液流动与传质传热模拟结果的影响，分析了塔板内速度场、温度场、浓度场与点效率的分布特征。模拟结果对理解精馏塔内的流动与传质过程具有实际意义，对其他化工多相流反应器内的热态模拟也具有参考价值。

2 模型与方法

采用OpenFOAM-8平台中的multiphaseEulerFoam求解器对精馏塔内的传热传质过程进行模拟。该求解器采用了多流体方法，将各相视为独立的流体并依次求解各相的守恒方程。为了模拟二元精馏体系中的传质传热过程，本研究在多流体框架的基础上耦合了能量方程与组分输运方程，通过求解各相的质量分数、速度、压力等物理量来描述多相流动状态。求解器的求解流程如图1所示。

图1 multiphaseEulerFoam求解器的计算流程图Fig.1 Flow chart of the multiphaseEulerFoam solver

2.1　动量方程

求解器multiphaseEulerFoam采用了速度-压力耦合的PIMPLE算法求解动量方程与连续性方程。首先使用隐式格式求解Navier-Stokes方程组以得到新的速度场，在此基础上再使用显式格式求解连续方程组，进而获得新的压力场。接下来，使用PISO算法修正压力场，并更新速度场和压力场，之后再次进行上述过程求解更新后的速度场和压力场，迭代直至收敛。由于采用了多流体框架，假设气液两相互相独立并互相渗透，气液两相的守恒方程如下：

(1)

(2)

式中，角标φ表示气相G或液相L，u为流体速度矢量(m/s)，为重力加速度，p为压力梯度(Pa/m)，ρ与α分别表示密度(kg/m³)与相含率，S_mt表示传质引起的质量源项[kg/(m³·s)]，M_D与M_mt分别表示流动和传质引起的相间动量源项[kg/(m²·s²)]。在本研究中，传质源项通过后文介绍的传质模型进行计算，而流动引起的相间动量源项采用了与前期流动模拟^[19]一致的方法，即仅考虑气液相间作用中起主导作用的曳力模型，其中详细的曳力模型^[21,22]公式如下：

(3)

(4)

在精馏塔板模拟过程中，湍流对流动、传热与传质过程均存在重要影响。本研究采用了混合相k-ε湍流模型^[23]，该模型被广泛用于具有明显体积分率变化的气液两相流模拟^[24,25]，其中湍动能k_m与湍流耗散率ε_m的计算方程如下：

(5)

(6)

2.2　能量方程

在OpenFOAM中，基本的能量方程具有内能(Internal Energy)与焓(Enthalpy)两种形式[式(7)]，两者可以通过式(8)进行转化，两者的主要区别在于热力学状态量的选择^[26]。内能是描述物质系统内部能量的状态量，表示各分子固有的能量之和；焓用于描述物质所具有的能量状态量，相较于内能增加了运动压力项。对于一般的CFD模拟，求解焓相较求解内能通常更加稳定和准确，因为焓包含全部能量变化，相较内能更为准确；而内能形式的能量方程更适用于处理湍流引起的涡旋效应与湍流热扩散等现象。由于精馏塔板上气液两相均具有明显的湍流特性，本研究选用了内能形式的能量方程。

(7)

(8)

(9)

(10)

式中，e和h分别表示内能(J/kg)和焓(J/kg)，T表示温度(K)。本研究采用了涡扩散模型量化湍流对能量扩散的影响，其中α_Eff表示混合物有效湍流热扩散系数[kg/(m·s)]，α_h为混合物焓的热扩散系数[kg/(m·s)]，α_t为湍流热扩散系数[kg/(m·s)]，由湍流黏度(m²/s)与湍流普朗特数进行计算。

在精馏塔板上的气液传热过程中，传质产生的热源项Q_mt是决定温度场的关键参数。本工作采用内能形式对传质引起的传热源项进行了量化，将其分为三部分：界面处传入内部的能量、潜热导致的能量变化以及动能变化。

(11)

(12)

(13)

式中，dm_f/dt为界面处的传质速率；K为传热系数，利用Ranz-Marshall模型进行量化；e_f表示生成内能(J/kg)；h_v为相变潜热(J/kg)，等于相间的总内能之差；e_K表示该相的动能(J/kg)；posPart与negPart是OpenFOAM内置的符号函数，分别表示其参数的正值与负值部分。在求解能量方程后，OpenFOAM将内能作为参数更新流体的热物理性质并输出，其中包括了温度、流体密度、黏度、热扩散系数及热导率等。

2.3 组分输运方程

气液两相的组分输运方程(Species Transport Equation)如下：

(14)

式中，Y_i表示各组分的质量分数；在组分输运方程中，采用基于涡流扩散系数的湍流梯度热流密度模型引入湍流对传热过程的影响，其中D_Eff表示有效质量扩散系数[kg/(m·s)]；dm_i/dt为组分i引起的传质量，前文的传质源项S_mt通过对各组分传质量求和得到；而动量源项M_mt则是通过各组分传质量和相应的速度，由求解器直接计算得到。本研究关注的重点是精馏过程中的传质现象，因此仅考虑液相的蒸发传质过程，蒸发速率采用界面组分模型进行描述，即：

(15)

(16)

(17)

(18)

式中，K_i_φ表示质量传递系数(m^-2)，由Higbie扩散质量传递模型进行计算；D为质量扩散系数，利用热扩散系数α_h,_i、相密度与Lewis数进行计算；Y_f为饱和状态下界面处的各组分质量分数，是影响传质过程的关键参数，其值受气液平衡关系约束，包括拉乌尔定律(Raoult's Law)和道尔顿分压定律(Dalton's Law of Partial Pressure)，即：

(19)

(20)

(21)

式中，p_i表示组分i的分压(Pa)，p_i^*表示组分i的饱和蒸气压(Pa)，p表示当前压力(Pa)，x_i与y_i分别为该组分在液相与气相中的摩尔分数，γ_i表示液相组分i的活度系数，当模拟理想体系时为1，模拟非理想体系时由活度系数模型进行计算，M_i表示组分i的摩尔质量(kg/kmol)，M表示气相混合物的平均摩尔质量(kg/kmol)，Y_f,_i表示组分i在界面处气相的质量分数。

由于乙醇-水精馏体系属于非理想溶液体系，需额外引入活度系数模型(Activity Coefficient Model)对上述计算得出的Y_f,_i进行修正。为此选用了两种常用的活度系数模型进行对比计算，分别是NRTL模型与UNIQUAC模型。

NRTL模型^[27]是由Renon和Prausnitz于1968年提出的一种非随机双液相模型。在该模型中，混合物的活度系数通过一组二元交互作用参数计算，以反映混合物中分子间相互作用的强度，其公式如下：

(22)

(23)

(24)

式中，x_i与x_j分别为两组分在液相中的摩尔分数，G_ij表示混合物中两种分子之间的平均相互作用强度，τ_ij为无量纲作用系数，a_ij, b_ij和c_ij为模型的经验参数。

UNIQUAC模型^[28]由Prausnitz和Guggenheim于1975年提出，该模型基于分子体积和官能团理论，通过计算不同分子间的相互作用能以及各组分在混合物中的相对比例得到混合物的活度系数，其公式如下：

(25)

(26)

(27)

(28)

(29)

(30)

式中，Φ_i表示平均体积分数，θ_i表示平均表面积分数，r_i为纯组分下的体积分数，q_i表示面积参数，配位数z通常取10。

此外，还需计算组分i的饱和蒸汽压p_i^*，常用的饱和蒸汽压模型包括ArdenBuck方程、Antoine方程与扩展Antoine方程等。本工作统一采用标准Antoine方程，分别对两种组分的饱和蒸汽压进行计算。Antoine方程^[29]是一种用于描述液体饱和蒸气压的经验模型，包含三个由实验测定的温度压力数据拟合得到经验参数，形式如下：

(31)

3 模拟设置与参数

本研究分别对环己烷-正庚烷理想体系与乙醇-水非理想体系进行模拟。理想体系的模拟对象为直径1.2 m的工业规模精馏体系^[6]，模拟物系为环己烷与正庚烷，塔板开孔率为14%，操作方式为全回流。非理想体系的模拟对象为包含180个筛孔的Noriler精馏塔板^[20]，模拟物系为乙醇与水，塔板开孔率为6.5%，堰高为0.08 m。在精馏体系的传质过程中，各组分的热物理性质决定了质量与热交换程度。为了与文献中实验数据和模拟结果进行对比验证，两种体系的结构参数、操作条件与物性参数分别参考自Yanagi等^[6]与Noriler等^[20]的数据以及NIST Chemistry Webbook^[30]，具体参数设置如表1所示。

表1 两种精馏体系^[6,20]的结构参数、操作条件及物性参数Table 1 The structural parameters, operating conditions and physical properties of two distillation systems^[6,20]

根据塔板上流场的对称性，两种体系均采用结构化网格对半块塔板进行模拟，气孔均被简化为单个网格面。由于本研究中理想体系的设备与前期研究^[19]中的塔板具有相同的尺寸，因此采用了经过前期研究验证的网格尺寸，网格总数为32 480，并将开孔率由5%修改为与Yanagi等^[6]的实验体系一致的14%，如图2(a)所示。非理想体系仅对塔板上方0.12 m的高度进行模拟。该体系网格采用了3种不同分辨率网格进行网格无关性验证，网格总数分别为Mesh 1 (10 440), Mesh 2 (36 000)和Mesh 3 (53 280)，其中Mesh 2如图2(b)所示。

图2 两种精馏体系单层塔板网格结构图：(a) 环己烷-正庚烷体系；(b) 乙醇-水体系Fig.2 Single tray grid structures of two distillation systems: (a) cyclohexane-n-heptane system; (b) ethanol-water system

边界条件通过指定计算域边界上的流体性质、流动行为或物理量的变化，对流体流动和传热过程的准确模拟起着重要作用。本研究中，两种体系的边界条件均为液体从液相进口处进入，从侧面矩形出口流出；气体通过筛孔进入，从顶面气相出口流出。非理想体系的边界条件均保持与文献[20]一致，以保证后续模拟结果可以与文献值进行对比。具体边界条件如表2所示，其中壁面处的边界条件在两种体系中完全相同。在传质过程的模拟中，气泡直径对传质系数与湍流扩散的量化有直接影响。实验研究^[31,32]表明，精馏塔内的气泡尺寸一般呈现大小气泡共存的双峰分布，其中大气泡与小气泡的直径分别约为10~20与2~5 mm，并且大气泡的占比远大于小气泡。由于本工作模拟的两种精馏体系均具有较高气速，所以假设气泡均为大气泡，选取20 mm作为气泡直径。

表2 两种精馏体系^[6,20]模拟的边界条件Table 2 The boundary conditions of two distillation systems^[6,20]

两种体系各组分的饱和蒸气压由Antoine方程进行计算，Antoine系数从NIST Chemistry Webbook^[30]获取，具体数值如表3所示。NRTL和UNIQUAC模型相关的经验参数来自于文献[33]，其中Water/Ethanol和Ethanol/Water分别表示式(22)~(30)中角标为ij和ji的二元相互作用参数。

表3 活度系数模型^[33]与Antoine方程中的经验参数Table 3 Parameters of the activity coefficient models^[33] and Antoine equation

4 结果与讨论

4.1　环己烷-正庚烷体系的模拟与验证

本研究对环己烷-正庚烷理想溶液体系进行模拟，并将浓度和温度的模拟结果与实验值进行对比以验证模型的准确性。采用与前期研究一致的气液两相流动模型，包括修正后的Bennett曳力模型与混合相湍流模型，分别对Yanagi和Sakata实验体系^[24]10层精馏塔的第3, 6和9块塔板进行模拟，模拟时长为60 s，平均时间窗口为30 s。图3分别展示了第3, 6和9块塔板出口处液相环己烷浓度与气相出口温度的模拟值与实验值。可以看出，实验测量与模拟预测的液相中环己烷浓度均随着塔板层数增加而上升。这是由于塔板层数降低时，作为轻组分的环己烷在塔内停留时间增加，气液接触的反应面积增大，导致传质效率增加，使更多轻组分从液相转移至气相中。浓度场与温度场模拟结果的平均相对误差分别为3.98%与0.21%，以上结果表明，当前热模模型能够较准确地描述精馏塔板上理想体系的传热传质过程。

图3 理想体系下不同塔板的模拟值与实验值：(a) 液相环己烷在液相出口处的摩尔分数；(b) 塔板出口温度Fig.3 The simulated and experimental values of different trays in the ideal system: (a) molar fraction of liquid cyclohexane; (b) outlet temperature of trays

4.2　乙醇-水体系的网格无关性与曳力模型修正

由于在实际精馏生产过程中大多数物系为非理想体系，因此利用构建的CFD热模模型对非理想体系的乙醇-水精馏塔板开展进一步模拟与分析。首先应用前述三种不同分辨率的网格进行网格无关性验证。为了保证不同网格分辨率下气孔分布保持一致，网格采用了局部加密方法，即对塔板底部区域(z=0~30 mm)与顶部区域(z=30~120 mm)分别采用了不同的网格分辨率。Mesh 1在底部区域处网格尺寸为3 mm×5 mm×5 mm，顶部区域为6 mm×10 mm×10 mm；Mesh 2采用了统一的3 mm×5 mm×5 mm网格尺寸；Mesh 3在底部区域处网格尺寸为1.5 mm×1.5 mm×5 mm，顶部为3 mm×5 mm×5 mm。三种网格预测的清液层高度如图4所示，网格数为10 440时预测的清液层高度远高于网格数量为36 000与53 280的预测结果，而后两种网格的预测结果基本一致，即可以近似认为36 000网格可满足网格无关性要求，所以后续计算中均采用该网格进行模拟。

图4 不同网格分辨率下清液层高度的模拟结果Fig.4 Simulation results of clear liquid height under different grid resolutions

图4中显示达到网格无关性后的清液层高度约为0.055 m，而文献[20]的结果约为0.047 m，两者存在较大差距。该结果表明初步模拟的气液相间分布存在误差，而曳力作为气液两相流中主要的相间作用力，对模拟的准确性有显著影响。因此，采用前期研究^[19]提出的引入额外气泡群影响的修正曳力模型与原始的Bennett曳力模型进行对比验证。如表4所示，采用修正曳力模型后清液层高度的预测值与文献[20]一致，而未修正的曳力模型高估了清液层高度。因此，在当前模拟框架下，采用原始的Bennett曳力模型会低估气泡群对液相的阻碍作用，而修正曳力模型能够有效改善模拟结果，并且该修正曳力模型具有较好的适用性。

表4 采用修正与未修正曳力模型预测的清液层高度Table 4 Clear liquid height predicted by using modified and unmodified drag models

图5给出了利用Bennett曳力模型和修正曳力模型预测的平均液相体积分数沿垂直塔板方向的分布曲线。将两种曳力模型的模拟结果与文献[20]的结果进行对比，可以看出三种液相体积分数具有相似的分布。由于大部分液相沉积于塔板底部，在塔板顶部区域液相体积分数接近于0。在塔板中部区域，气相与液相混合程度增强，平均液相体积分数约为0.5。在该区域本研究中采用的两种曳力模型均高估了液相体积分数，但修正曳力模型的预测结果与文献[20]更加接近，这是由于修正曳力模型增强了气液相互作用，导致气含率增加，从而使液相体积分数降低。在塔板底部靠近筛孔处，修正曳力模型的预测结果显示出与文献[20]相似的非均匀分布，进一步证明了修正曳力模型的可靠性。因此，在后续模拟中均采用修正曳力模型描述气液相间动量传递。

图5 平均液相体积分数沿塔板垂直方向的分布Fig.5 Distributions of average liquid volume fraction along the vertical direction of the tray

4.3　乙醇-水体系中物理量随时间的变化规律

分别对清液层高度、液相温度和液体中乙醇的质量分数进行了体积平均统计，图6展示了三者随时间的变化关系。从图6(a)可以看出，清液层高度在0~10 s呈上升趋势，10 s后围绕着一个平均值振荡，液相温度和液相中乙醇质量分数在30 s后也会在一定数值范围内振荡。这表明系统在30 s后达到了准稳态过程，因此本研究采用60~90 s的时均结果与文献[20]中的结果进行对比，如表5所示。可以看出，本工作预测的各物理量时均统计结果与文献[20]结果基本一致，说明本研究构建的传热传质模型能够较好地描述精馏塔内的传热传质过程。

图6 塔板上不同物理量随时间的变化趋势：(a) 清液层高度；(b) 液相温度；(c) 液相中乙醇质量分数Fig.6 The temporal tendencies of different variables on the tray: (a) clear liquid height; (b) liquid temperature; (c) mass fraction of liquid ethanol

表5 塔板上不同物理量的时均值Table 5 The time-averaged values of different variables on the tray

4.4　活度系数模型对气液传质与传热过程的影响

图7给出了未采用活度系数的模型与两种活度系数模型预测的液相乙醇质量分数沿液相主流方向的变化趋势，三种模型预测的质量分数变化与文献[20]的结果趋势相同。在液相入口处液体快速涌入并与气相发生混合，局部气液相间具有较大的温度与浓度梯度，促使液相中作为轻组分的乙醇快速转化为气相，造成液相中乙醇质量分数先迅速降低，然后沿液体主流方向逐渐减少。在数值上，文献[20]中报道的液相乙醇质量分数随主流方向的变化范围为0.785至0.760，未采用活度系数模型修正的结果相对低估了塔板上方的乙醇浓度，其原因是乙醇-水溶液是典型的正偏差溶液，乙醇与水分子间产生的额外相互作用导致混合溶液中乙醇溶解度升高，而采用理想模型忽略了异种分子间相互作用和偏离行为，因而使预测结果低于文献值。所以，对于非理想体系的传热传质模拟，进行活度系数的模型修正是必要的。

图7 沿液体主流方向的液相乙醇平均质量分数变化Fig.7 Variation of average mass fraction of liquid ethanol along the liquid main flow direction

图8显示了平均液相温度从液体入口到出口的变化，未采用活度系数模型的计算结果由于与实验值差距较大而未在图中展示。可以看出，两种活度系数模型的模拟结果与文献[20]的结果变化趋势一致。沿液相主流方向，液体初始温度较低，在液相入口附近气液两相迅速发生热量交换，导致液相温度迅速上升。而随着液相的流入，液相温度的上升速率减慢，一直保持到出口堰附近。在数值上，文献[20]的液相温度随液相主流方向的变化范围为355.07~355.47 K，UNIQUAC模型预测的结果变化范围为355.13~355.35 K，NRTL模型预测的结果变化范围为355.64~355.94 K，未采用活度系数模型的预测结果变化范围为352.9~352.3 K。这说明在采用理想溶液传质模型时，由于低估了气液相平衡时的液相乙醇饱和浓度导致了传质量偏高，因此在计算能量方程时传质源项的高估引起了过量的热量交换，使液相温度的整体预测值低于文献值。在采用当前模型及参数设置下，UNIQUAC模型的预测结果与文献[20]更接近，NRTL模型整体高估了当前模拟条件下的液相温度，其原因可能是UNIQUAC模型描述的乙醇与水分子间相互作用更符合现实。此外，由于气泡直径决定着传热与传质系数以及湍流扩散过程，因此气泡直径的设置也可能导致不同模型的预测值出现偏差，为了获取更加准确的传质与传热模拟结果，将来可以考虑引入更为恰当的气泡尺寸分布模型。

图8 沿液体主流方向的液相温度变化Fig.8 Variation of liquid temperature along the liquid main flow direction

4.5　气液两相流场对传质过程的影响

上文分别讨论了气液流动与传质传热过程的定量结果，而在实际情况中流场分布对相间传质过程存在直接影响。图9给出了高度为10 mm处的液相乙醇浓度分布与液速矢量分布图。根据图9(a)所示，乙醇质量分数在液相入口附近快速下降，而在塔板后半程则逐渐趋于平稳。同时，对比图9(a)与9(b)可知，传质作用剧烈的位置与液相涡旋位置几乎重合。这说明液相入口处气液相互作用引起了液相回流并形成了涡旋，这增加了液相的停留时间和气液混合程度。显著的温度梯度使气液相间热量发生剧烈交换，在浓度梯度及温度影响下乙醇由液相快速转化为气相并溢出。此外，图9(b)的矢量图显示出在液相出口处由于堰的阻碍作用使液相形成了涡旋，该现象可能导致液相出口处的液相停留时间延长，进而使液相乙醇的汽化量增加。

图9 在z=10 mm截面上的乙醇浓度与速度分布：(a) 液相乙醇质量分数分布；(b) 液相速度矢量图Fig.9 Ethanol concentration and velocity distribution at the cross section of z=10 mm: (a) distribution of mass fraction of liquid ethanol; (b) liquid velocity vector diagram

图10给出了塔板上气相和液相速度的二维流线图。由图10(a)可知，气相基本沿塔板垂直方向快速通过了泡沫区域，没有形成明显的气相循环。图10(b)和10(c)显示，液相通过泡沫区域时在液相进口和出口堰附近轴向和径向均形成了两个液相循环区域，进而导致液相的停留时间增加。

图10 不同区域的二维流线图：(a) 塔板的中心剖面(y-z)的气相速度流线图；(b) 塔板的中心剖面(y-z)的液相速度流线图；(c) 塔板上方10 mm横截面(x-y)的液相速度流线图Fig.10 Two-dimensional streamlines of different regions: (a) gas velocity in the central plane (y-z); (b) liquid velocity in the central plane (y-z); (c) liquid velocity in cross section (x-y) at 10 mm above the tray bottom

图11统计了从液相入口到出口的Murphree点效率，在塔板的液相进口和出口处局部传质效率较高，而中间部分的局部传质效率最小，这是因为液相在该区域具有稳定的流动路径，从而导致液相与气相的混合不够充分。以上结果表明气液两相流场对传质过程具有直接影响，液相循环流动通过提高气液接触时间来提高精馏塔板局部区域的传质效率。但需要注意的是，在涡旋区域内，气液两相无法及时更新，可能会对塔板整体的气液传质效率产生负面影响。因此，在精馏塔板的设计和优化过程中，需要同时考虑上述两种效应。综上所述，CFD模拟可以获得详细的内部流场和传热传质信息，为精馏塔的优化设计提供理论指导。

图11 沿液相入口至出口的Murphree点效率分布Fig.11 Murphree point efficiency distribution from liquid inlet to liquid outlet

5 结论

利用OpenFOAM平台中的multiphaseEulerFoam多流体求解器，建立了基于气液相平衡的传热传质模型，分别对环己烷-正庚烷与乙醇-水体系进行了验证与模拟，给出了塔板上气液两相质量分数、温度及速度等参数的分布规律。

(1) 模拟发现前期研究中提出的修正曳力模型同样适用于环己烷-正庚烷和乙醇-水精馏塔热态模拟体系，具有良好的准确性与适用性。

(2) 对于环己烷-正庚烷理想体系，采用理想溶液传质模型进行模拟，预测的出口浓度及温度平均相对误差分别为3.98%与0.21%，验证了当前CFD传热传质模型的可靠性。

(3) 对于乙醇-水非理想体系，统计结果表明，气液两相流场、浓度场与温度场均会随着模拟时间增加达到准稳态过程，之后清液层高度、温度以及组分质量分数分别在其平均值附近波动；对比了UNIQUAC模型、NRTL模型与未采用活度系数模型的模拟结果，采用活度系数模型能够显著提升对浓度场与温度场的预测精度；在当前模拟框架与操作条件下，两种活度系数模型均可以准确预测浓度场分布，而对于温度场的预测UNIQUAC模型相较于NRTL模型与文献[20]结果更为接近。

(4) 考察了气液两相流场对传质过程的影响，发现液相在入口与出口附近的循环流动具有强化局部传质的作用，但循环区域中气液两相无法及时更新，可能导致塔板整体传质效率下降。

在下一步工作中，计划构建适合精馏体系的气泡尺寸分布模型，以更好地反映气泡直径对流动传质过程的影响，从而降低当前模拟中温度场与浓度场的预测误差。综上所述，基于开源平台OpenFOAM的热态模拟能够获取精馏塔的详细内部流场及传热传质信息，有望为化工多相流设备的优化设计提供新的理论指导和解决方案。

略

Heat and mass transfer simulation of gas-liquid two-phase flow in a distillation column based on OpenFOAM

Yunpeng JIAO ^1,² Xiaoqing ZHOU ^1,³Jianhua CHEN ¹

1. State Key Laboratory of Multiphase Complex Systems, Institute of Process Engineering, Chinese Academy of Sciences, Beijing 100190, China
2. School of Chemical Engineering, University of Chinese Academy of Sciences, Beijing 100049, China
3. School of Chemical Engineering and Technology, Shenyang University of Chemical Technology, Shenyang, Liaoning 110142, China ）

Abstract: For multiphase mass transfer processes in distillation columns, computational fluid dynamics can be used not only to simulate the flow phenomena, but also to study the heat and mass transfer processes and their interactions with the flow. Based on previous work on gas-liquid flow simulation in distillation trays, this work applied the multiphase flow solver of the OpenFOAM platform. The energy and species transport equations were considered to construct the heat and mass transfer models for the distillation system. The ideal system of cyclohexane-n-heptane and the non-ideal system of ethanol-water distillation were simulated respectively. The distribution of the gas-liquid two-phase flow, component and temperature on the column trays were analyzed. For the ideal system, the ideal solution mass transfer model can provide accurate predictions for the temperature and concentration fields on the trays. However, for non-ideal systems, it was necessary to introduce the activity coefficient model on the basis of the ideal model. To this end, the effects of two activity coefficient models, UNIQUAC and NRTL, were introduced and compared. In the current simulation framework, the activity coefficient models were able to improve the simulation accuracy of temperature and concentration field. The overall trends predicted by the two models were generally consistent with the results in the literature, and the UNIQUAC model agreed better with the literature. In addition, the comparative analysis of gas-liquid two-phase flow field and concentration field distribution showed that the circulating flow of liquid phase can enhance the local mass transfer efficiency of the column tray, resulting in higher efficiency at the liquid inlet and the weir than the tray center. However, the gas-liquid renewal in the circulating region renewed promptly, which led to a reduction in the overall mass transfer efficiency of the column tray. This study can be used for the design and optimization of distillation columns, and it is also valuable for simulations of heat and mass transfer in other gas-liquid two-phase flow systems.

Keywords: gas-liquid two-phase flow；distillation；hot state simulation；heat and mass transfer；OpenFOAM

引用本文：焦云鹏, 周晓庆, 陈建华. 基于OpenFOAM的精馏塔内气液两相流传热传质模拟. 过程工程学报, 2024, 24(4): 391-402. (Jiao Y P, Zhou X Q, Chen J H. Heat and mass transfer simulation of gas-liquid two-phase flow in a distillation column based on OpenFOAM (in Chinese). Chin. J. Process Eng., 2024, 24(4): 391-402, DOI: 10.12034/j.issn.1009‑606X.223180.)

作者简介：焦云鹏，博士研究生，化学工程专业，E-mail: jiaoyunpeng@ipe.ac.cn；

作者简介：陈建华，副研究员，研究领域为介尺度模型与多相流数值模拟，E-mail: jhchen@ipe.ac.cn

基金信息：国家自然科学基金资助项目(编号：22078327)；多相复杂系统国家重点实验室资助项目(编号：MPCS-2021-A-11)

中图分类号： TQ021.3

文章编号：1009-606X(2024)04-0391-12

文献标识码： A

收稿日期：2023-06-27

修回日期：2023-08-25

出版日期：2024-04-28

网刊发布日期：2024-05-08

http://mp.weixin.qq.com/s?__biz=MzI0MzM2MjgyOQ==&mid=2247529526&idx=2&sn=579b53c52dcdddddb9c36aa9560afb2e

过程工程学报

《过程工程学报》(月刊)创刊于1976年，由中国科学院过程工程研究所主办、科学出版社出版。《过程工程学报》以过程工程科学为学科基础，重点刊登材料、化工、生物、能源、冶金、石油、食品、医药、资源及环境保护等领域中涉及过程工程的原创论文。

焦云鹏, 周晓庆, 陈建华丨基于OpenFOAM的精馏塔内气液两相流传热传质模拟

1 前 言

2 模型与方法

2.1 动量方程

2.2 能量方程

2.3 组分输运方程

3 模拟设置与参数

4 结果与讨论

4.1 环己烷-正庚烷体系的模拟与验证

4.2 乙醇-水体系的网格无关性与曳力模型修正

4.3 乙醇-水体系中物理量随时间的变化规律

4.4 活度系数模型对气液传质与传热过程的影响

4.5 气液两相流场对传质过程的影响

5 结 论