基于分离因数法的简捷炼油塔模型开发及应用

2023年6月周磊等：基于分离因数法的简捷炼油塔模型开发及应用
蜡油和减压渣油等产品。

减压塔共有19块塔板，设置3条侧线采出，减一线流量为23.5t/h ，流股压力为-0.0991MPa 。

减二线流量为110t/h ，流股压力为-0.0987MPa 。

减三线流量为195t/h ，流股压力为-0.0984MPa 。

减底流量为280t/h ，流股压力为-0.0981MPa 。

2.1 初馏塔计算
为验证模型的可行性，对上述常减压装置建模。

调整模型参数，默认切割温度点组分分配比ln(d r /b r )为0。

初馏塔、常压塔采用切割温度模式，
减压塔采用收率模式计算，将常减压计算结果与实际D86
蒸馏曲线数据和流量相比对。

图2　简捷炼油塔计算框图
·
·2823
化工进展， 2023， 42（6）
表3为初馏塔输入的模型参数。

从图4初馏塔计算结果可以看出，分离因数模型在模拟初馏塔塔顶蒸馏曲线在前后两端与实际值较好吻合，与实际偏差都在2%~4%之间。

计算初顶油采出量为153.31t/h ，与实际偏差为0.046%。

2.2 常压塔计算
表4为模拟常压塔输入的模型参数。

从图5常压塔模拟计算结果中可以看出，分离因数模型计算的蒸馏曲线与实际曲线贴近，各线产品采出平均偏差均在5%以内，说明可以很好模拟常压塔分离产品馏分。

表5为计算出的各线产品流量，与实际偏差均在1.5%以内。

2.3 减压塔计算
表6为模拟减压塔输入的模型参数。

图6为减
压塔计算模拟结果对比，由于减压塔整体馏分偏
图3　常减压装置流程实例
表3　初馏塔输入参数
产品名称初顶气初顶油初底油
采出类型气相顶油底油
切割温度/℃—5.24
180.03顶部
分离因数—55
底部分离因数—104
流股压力/MPa 0.0950.11
0.134表4　常压塔输入参数
产品名称常顶气常顶油常一线常二线常三线常四线常底油
采出类型气相顶油侧线1侧线2侧线3侧线4底油
切割温度/℃—-1.39
159.97239.32282.92401.69434.11顶部分离因数—5
3515303522底部分离因数—101010298
10流股压力/
MPa
0.030.1
0.110.110.120.120.12表2　轻端组分分析
序号1234567
组分H 2N 2O 2
CO
CO 2H 2S
CH 4质量分数0.0000920.0103
0.0005750.0004020.006280.01390.00611序号89
1011121314组分C 2H 6C 3H 8NC 4IC 4
NC 5C 2H 4C 3H 6
质量分数0.03220.2890.4860.149
0.00580.0000480.000031
D 86温度/℃
体积分数/%
图4　初馏塔D86蒸馏曲线结果对比
表1　原油TBP 蒸馏数据（质量基准，d 15.6/15.6重度为0.872）
蒸馏曲线百分比/%
05
103050709095100温度/℃0.38
47.7292.68
214.38334.89476.51671.08773.09863.71·
·2824
2023年6月
周磊等：基于分离因数法的简捷炼油塔模型开发及应用
重，初馏点的计算值与实际值偏差稍大。

但整体趋势上看，计算值与实际值基本一致，减一线、减二线、减三线计算值与实际数据平均偏差分别为0.41%、2.3%、2.8%。

表7为减压塔计算的各产品流量，与实际偏差均在0.5%以内。

图7为基于分离因数法模型模拟常减压三塔的误差分析，其中塔顶油产品使用馏程大于50%点
的数据，其余产品流股使用去除初馏点的数据。

从误差分析结果可以看出，平均相对误差在5%以内，与实际蒸馏曲线相近。

表8为使用简捷炼油塔模拟常减压三塔的计算模式时间与迭代次数。

模拟整个常减压流程仅需13.7s （包括原油预处理及结果输出时间，CPU 为AMD Ryzen 7 4800H ），分离因数模型计算三塔的时间不到0.1s 即可完成各产品组分流量计算。

采用切割温度模式下，可直接快速地完成组分分离计算，规定收率模式下，产品流股迭代3~5次即可完成收敛，顶油流股迭代次数在10次左右，收敛速度快，简化炼油塔分离原油的复杂性。

从模拟实例结果中可以看出，采用分离因数法简捷炼油塔模型模拟常减压过程，各产品采出量与实际值贴近。

在模拟常压塔的各产品质量与实际数据趋势一致，数值较吻合。

减压塔计算值在初馏点有差距，其余点与实际数据基本一致，并且都在工程允许偏差范围5%以内。

表明该简捷炼油塔模型较好地模拟炼油塔分离过程，并且计算得到的产品蒸馏曲线与实际数据基本一致。

D 86温度/℃
体积分数/%
图5　常压塔D86蒸馏曲线结果对比表5　常压塔流量计算结果对比
产品名称常一线流量常二线流量常三线流量常四线流量
实际值/t·h -1
110.0078.00
122.0032.00模型值/t·h -1
108.4477.41
120.6631.85相对偏差/%
1.420.751.100.46
表8　常减压三塔计算模式时间与迭代次数
项目名称初馏塔常压塔
减压塔
选用模式切割温度切割温度收率
计算时间/s 0.0780.0930.094
迭代次数——27
D 86温度/℃
体积分数/%
图6　减压塔D86蒸馏曲线结果对比
表7　减压塔流量计算结果对比
产品名称减一线流量减二线流量减三线流量减底流量
实际值/t·h -1
23.50
110.00195.00280.00模型值/t·h -1
23.53
109.57195.91280.98相对偏差/%-0.120.39
-0.47-0.35表6　减压塔输入参数
产品名称减顶气减顶油减一线减二线减三线减底油
采出类型气相顶油侧线1侧线2侧线3底油
收率/%—2.023.78
17.5831.4345.08顶部分离
因数
—5
30373055底部分离
因数
—10104
9030流股压力/MPa
-0.0992-0.0992-0.0991-0.0987-0.0984-0.0981
初
顶
油
常
顶
油常一线常二线常三线常四线减一线减二线减三线相对误差/%
图7　常减压三塔误差分析
·
·2825
化工进展， 2023， 42（6）
3 结论
开发了一种基于分离因数的简捷炼油塔模型，可以快速、简便地计算出各塔产品质量和流量，模拟炼油塔油品分离过程，并进行实例验证，结论如下。

（1）该模型将炼油塔原油分离过程近似于连续的分离单元，大大简化计算炼油塔的复杂程度，兼顾简捷性和准确性，模拟常减压装置的产品蒸馏曲线与实际相对比都在5%偏差以内，适合作为实际工厂生产调度模型。

（2）分离因数模型相对于其他模型具有更强的灵活性，调整分离因数可以模拟计算炼油塔不完全分割的馏分产品及性质，而且没有大多数简化模型认为切割点在相邻产品前后两端中点，并位于原油蒸馏曲线上的假设，模型更加符合实际过程。

符号说明
B
j
——j单元的底部流量，kmol/h
b
j,i
——j单元中i组分底部流量，kmol/h
b r——切割温度组分底部流量，kmol/h
C B,N，C D,N——分别为上段和下段直线公式的截距
D
j
——j单元的顶部流量，kmol/h
d
j,i
——j单元中i组分顶部流量，kmol/h
d r——切割温度组分顶部流量，kmol/h
F
j
——j单元的进料流量，kmol/h
F total——总进料流量(不计组分水)，kmol/h
f
j,i
——j单元中i组分进料流量，kmol/h
k1，k2——分别为上段和下段直线公式的斜率SI——分离因数
T——分离单元
T b——切割温度，K
T
i
——i组分的正常沸点，K
N——产品采出数目
x top,i，x top,o——分别为i组分和关键组分的顶部馏出物液相
组成，摩尔分数
x bot,i，x bot,o——分别为i组分和关键组分的底部采出物液相
组成，摩尔分数
Y
j
——j单元的产品收率，%
α
i,o
——i组分与关键组分之间相对挥发度
下角标
i——组分编号
j——分离单元编号
n——组分总个数
o——关键组分
r——沸点为切割温度的组分
参考文献
[1]黄德先, 江永亨, 金以慧. 炼油工业过程控制的研究现状、问题与展
望[J]. 自动化学报, 2017, 43(6): 902-916.
HUANG Dexian, JIANG Yongheng, JIN Yihui. Present research situation, major bottlenecks, and prospect of refinery industry process control[J]. Acta Automatica Sinica, 2017, 43(6): 902-916. [2]ZHANG Nan, SMITH Robin, BULATOV Igor, et al. Sustaining high energy efficiency in existing processes with advanced process integration technology[J]. Applied Energy, 2013, 101: 26-32. [3]刘广杰, 孙晓岩, 毕荣山, 等. 复杂炼油塔模拟的改进联立方程内外
层算法[J]. 化工学报, 2020, 71(5): 2118-2127.
LIU Guangjie, SUN Xiaoyan, BI Rongshan, et al. Improved equation oriented inside-out method of complex crude distillation column simulation[J]. CIESC Journal, 2020, 71(5): 2118-2127.
[4]SMITH Robin, OCHOA-ESTOPIER Lluvia M, JOBSON M. The use of reduced models in the optimization of energy integrated processes[J]. Chemical Engineering Transactions, 2013, 35: 139-144. [5]KHOR Cheng Seong, VARVAREZOS Dimitrios. Petroleum refinery optimization[J]. Optimization and Engineering, 2017, 18(4): 943-989.
[6]YANG Y, BARTON P I. Integrated crude selection and refinery optimization under uncertainty[J]. AIChE Journal, 2016, 62(4): 1038-1053.
[7]URIBE-RODRIGUEZ Ariel, CASTRO Pedro M, GONZALO Guillén Gosálbez, et al. Global optimization of large-scale MIQCQPs via cluster decomposition: Application to short-term planning of an integrated refinery-petrochemical complex[J]. Computers & Chemical Engineering, 2020, 140: 106883.
[8]JI X C, HUANG S M, GROSSMANN I. Integrated operational and financial hedging for risk management in crude oil procurement[J]. Industrial & Engineering Chemistry Research, 2015, 54(37): 9191-9201.
[9]KALLRATH Josef. Planning and scheduling in the process industry[J]. OR Spectrum, 2002, 24(3): 219-250.
[10]KETABCHI Elham, MECHLERI Evgenia, ARELLANO-GARCIA Harvey. Increasing operational efficiency through the integration of an oil refinery and an ethylene production plant[J]. Chemical Engineering Research and Design, 2019, 152: 85-94.
[11]ZHAO Hao, RONG Gang, FENG Yiping. Effective solution approach for integrated optimization models of refinery production and utility system[J]. Industrial & Engineering Chemistry Research, 2015, 54(37): 9238-9250.
[12]LI Jie, XIAO Xin, BOUKOUVALA Fani, et al. Data-driven mathematical modeling and global optimization framework for entire petrochemical planning operations[J]. AIChE Journal, 2016, 62(9): 3020-3040.
[13]BROOKS R W, WALSEM F D V, DRURY J. Choosing cutpoints to optimize product yields[J]. Hydrocarbon Processing, 1999, 78(11): 53-56， 58， 60.
[14]PINTO J M, JOLY M, MORO L F L. Planning and scheduling models for refinery operations[J]. Computers & Chemical Engineering, 2000, 24(9/10): 2259-2276.
[15]NEIRO Sérgio M S, PINTO José M. A general modeling framework for the operational planning of petroleum supply chains[J]. Computers & Chemical Engineering, 2004, 28(6/7): 871-896.
[16]GUERRA O J, ROUX G A C L. Improvements in petroleum refinery planning: 2. case studies[J]. Industrial & Engineering Chemistry Research, 2011, 50(23): 13419-13426.
[17]GUERRA O J, ROUX G A C L. Improvements in petroleum refinery planning: 1. Formulation of process models[J]. Industrial & Engineering
··2826
2023年6月周磊等：基于分离因数法的简捷炼油塔模型开发及应用
Chemistry Research, 2011, 50(23): 13403-13418.
[18]ZHANG J, ZHU X X, TOWLER G P. A simultaneous optimization strategy for overall integration in refinery planning[J]. Industrial & Engineering Chemistry Research, 2001, 40(12): 2640-2653. [19]LI Wenkai, HUI Chiwai, LI Anxue. Integrating CDU, FCC and product blending models into refinery planning[J]. Computers & Chemical Engineering, 2005, 29(9): 2010-2028.
[20]MENEZES Brenno C, KELLY Jeffrey D, GROSSMANN Ignacio E. Improved swing-cut modeling for planning and scheduling of oil-refinery distillation units[J]. Industrial & Engineering Chemistry Research, 2013, 52(51): 18324-18333.
[21]ALATTAS Abdulrahman M, GROSSMANN Ignacio E, PALOU-RIVERA Ignasi. Integration of nonlinear crude distillation unit models in refinery planning optimization[J]. Industrial & Engineering Chemistry Research, 2011, 50(11): 6860-6870.
[22]MAHALEC Vladimir, SANCHEZ Yoel. Inferential monitoring and optimization of crude separation units via hybrid models[J]. Computers & Chemical Engineering, 2012, 45: 15-26.
[23]FU G, SANCHEZ Y, MAHALEC V. Hybrid model for optimization of crude oil distillation units[J]. AIChE Journal, 2016, 62: 1065-1078.
[24]LI Fupei, QIAN Feng, YANG Minglei, et al. Product tri-section based crude distillation unit model for refinery production planning and
refinery optimization[J]. AIChE Journal, 2021, 67(2): e17115. [25]FU G, MAHALEC V. Comparison of methods for computing crude distillation product properties in production planning and scheduling[J]. Industrial & Engineering Chemistry Research, 2015, 54(45): 11371-11382.
[26]陈志奎. SHORTCUT——石油蒸馏塔简捷模拟块[J]. 石油炼制与化
工, 1986, 17(10): 1-6.
CHEN Zhikui. SHORTCUT—Simple simulation block of petroleum distillation column[J]. Petroleum Processing and Petrochemicals, 1986, 17(10): 1-6.
[27]GEDDES R L. A general index of fractional distillation power for hydrocarbon mixtures[J]. AIChE Journal, 1958, 4(4): 389-392. [28]GILBERT R J H, LEATHER J, ELLIS J F G. The application of the Geddes fractionation index to crude distillation units[J]. AIChE Journal, 1966, 12(3): 432-437.
[29]Aspen Plus. Getting started modeling processes with solids[R]. Cambridge: Aspen Technology, Inc, 2006.
[30]张海娟, 唐全红. 基于CrudeManager的常减压装置模拟研究[J]. 计
算机与应用化学, 2013, 30(8): 887-890.
ZHANG Haijuan, TANG Quanhong. Research on simulation of atmospheric and vacuum distillation unit basis of CrudeManager[J]. Computers and Applied Chemistry, 2013, 30(8): 887-890.
·
·2827
化
工进展
Chemical Industry and Engineering Progress
2023 年第 42 卷第 6 期
鼓泡流态化气泡间相互作用下相间传质过程的模拟
章凯1，金捍宇1，刘思宇2，王帅1
（1 哈尔滨工业大学能源科学与工程学院，黑龙江 哈尔滨 150001；2 华能国际电力股份有限公司，北京 100031）摘要：气泡-乳化相相间传质行为是流化床反应器设计与调控的关键问题，气泡间的相互作用增加了相间传质过程的复杂性。

本文利用CFD-DEM 方法对三维射流鼓泡床内垂直分布的双气泡传质过程展开模拟，基于建立的对流-扩散传质系数独立确定方法，分析了气泡存在相互作用时不同传质机制的演变过程，对比了双气泡与单气泡传质过程的差异性。

研究结果表明，与单气泡相比，上方气泡的存在会使得下方气泡轴向拉伸，减小了下方气泡体积并提高了气泡上升速度。

对相间传质过程而言，上方气泡浓度的释放增加了下方气泡的扩散传质阻力，且使其出现逆向对流传质现象。

不同入射时间间隔也会进一步影响传质过程，气泡复杂的传质过程将是气泡自身气流特性、形状演变与气泡间浓度干扰等多因素作用的结果。

关键词：流态化；介尺度；气泡；传质；CFD-DEM 方法
中图分类号：TF803.11+4 文献标志码：A 文章编号：1000-6613（2023）06-2828-08
Simulation of mass transfer process under the bubble interaction in
bubbling fluidization
ZHANG Kai 1，JIN Hanyu 1，LIU Siyu 2，WANG Shuai 1
(1 School of Energy Science and Engineering, Harbin Institute of Technology, Harbin 150001, Heilongjiang, China;
2
Huaneng International Power Co., Ltd., Beijing 100031, China)
Abstract: The study of mass transfer behaviors between bubble phase and emulsion phase is plays an important role in the design and regulation of fluidized bed reactor. The interaction between bubbles increases the complexity of mass transfer process. In this paper, computational fluid dynamics-discrete element method (CFD-DEM) was used to simulate the mass transfer process of two vertically distributed bubbles in a three-dimensional injected bubbling fluidized bed. On the basis of the independent determination method of mass transfer coefficients by convection and diffusion, the evolution of different mass transfer mechanisms with mutual interaction between bubbles was analyzed and the difference of mass transfer processes of two-bubbles and single-bubble was compared. The results showed that compared with the single bubble, the lower bubbles were axially stretched due to the existence of the upper bubbles. The volume of lower bubbles was reduced and the bubble rising velocity was increased. For the mass transfer process, the release of gas species from the upper bubbles increased the resistance of the diffusion-induced mass transfer of the lower bubbles and the reverse convective mass transfer was present. The interval of injection time further affected the mass transfer process. The complex mass
研究开发
DOI ：10.16085/j.issn.1000-6613.2022-1525
收稿日期：2022-08-18；修改稿日期：2022-10-26。

基金项目：国家自然科学基金（52076060）。

第一作者：章凯（1997—），男，博士研究生，研究方向为多相流体力学。

E-mail ：**************。

通信作者：王帅，教授，博士生导师，研究方向为多相流体力学。

E-mail ：*****************.cn 。

引用本文：章凯, 金捍宇, 刘思宇, 等. 鼓泡流态化气泡间相互作用下相间传质过程的模拟[J]. 化工进展, 2023, 42(6): 2828-2835.
Citation ：ZHANG Kai, JIN Hanyu, LIU Siyu, et al. Simulation of mass transfer process under the bubble interaction in bubbling fluidization[J]. Chemical Industry and Engineering Progress, 2023, 42(6): 2828-2835.
·
·2828
2023年6月章凯等：鼓泡流态化气泡间相互作用下相间传质过程的模拟
transfer process of bubbles depended on many factors including the gas flow characteristics, the evolution of bubble shape and the concentration interference between bubbles.
Keywords: fluidization; mesoscale; bubble; mass transfer; computational fluid dynamics-discrete element method
鼓泡流化床因其具有良好的颗粒混合效果和传热传质特性在各个领域应用广泛，对鼓泡床反应器内相间传递机理的认知是反应器优化设计的基础[1]。

气泡作为鼓泡流态化中最为典型的特征结构，其在床层内部的气流交换中扮演着重要角色，认知气泡相间传质机理对于反应器反应过程整体预测与调控具有重要意义[2]。

近年来，随着计算机技术的发展，数值模拟技术越发成为研究鼓泡流化床内流动与相间质量传递特性的重要手段[3-4]。

Vepsäläinen等[5]利用计算流体力学方法（CFD）研究了气泡相间传质行为，研究结果表明，在B类和D类颗粒床层内，增大气流速度会提高气泡对流通流，进而增大传质系数。

Tan 等[6]利用计算流体力学-离散单元模型耦合方法（CFD-DEM）模拟了B类和A类颗粒床层内单气泡传质过程，研究发现，Davidson等[7]提出的势流理论模型能够较好地预测B类颗粒床层中无晕气泡传质系数，然而对于A类颗粒床层中的含晕气泡，传统模型低估了传质过程中扩散系数的影响。

对流和扩散是气泡相间传质的两种主要方式，目前对气泡相间传质系数的计算大多基于整体统计方法[6, 8]，获得的整体传质系数可近似代表大颗粒床层中对流传质项的贡献（气泡类型为无晕气泡，对流传质作用远大于扩散传质）。

然而对于小颗粒床层中形成的含晕气泡，扩散传质作用增强，仅依赖整体法无法评估两种传质机制的分别贡献。

在本文作者课题组[9]先前的工作中，一种基于气泡边界积分的方法被用于分别统计含晕气泡对流与扩散传质系数，该方法对单气泡相间传质系数的统计与
K&L模型[2]预测结果取得了较好的一致性。

对于实际应用中的鼓泡床反应器来说，气泡存在着相互干扰、聚并、破碎等行为[10]，气泡间的相互作用进一步加深了反应器流动与传质过程的非均匀性与复杂性[1, 3]。

研究表明，与单气泡相比，自由鼓泡流化床中气泡传质行为会随着气泡通流的升高而显著增强[8, 11]。

为了深入理解气泡间相互作用对含晕气泡相间传质过程的影响机制，本文将基于CFD-DEM方法对双气泡射流鼓泡床内气泡传质行为展开模拟，并结合构建的对流与扩散传质系数统计方法，定量分析不同传质机制在气泡动态演变过程中的变化，明晰气泡传质过程的内在机理，为鼓泡床反应器设计优化与调控奠定基础。

1 数学模型
1.1 颗粒相控制方程
在CFD-DEM模型中，颗粒相被离散化处理，颗粒形状均假设为球形。

颗粒i的平动和转动方程表示如式(1)、式(2)所示[12]。

m
i
d u p,i
d t=m i g+F d,i+F c,i(1)
I
i
dωp,i
d t=T i(2)
式中，颗粒间作用力F c,i由软球模型中的线性弹簧-阻尼模型求得[13]；F d,i表示颗粒相与气相相间作用力，由压力梯度力和曳力两部分组成，如式(3)所示。

F d,i=-V i∇P g+βi V iεp,i(u g,i-u p,i)(3)
式中，β表示气固相间曳力系数，本文采用Gidaspow曳力模型[14]，如式(4)~式(6)所示。

β=
ì
í
î
ïïï
ï
ï
ï
ï
3
4C Dρgεgεp,i
|
|u g-u p,i
d pε
-2.65gεg≥0.8 150εp,i(1-εg)μg
εg d2p
+1.75ρgεp,i
|
|u g-u p,i
d pεg<0.8
(4)
C D=
ì
í
î
ï
ïï
24
Re p(1+0.15Re0.687p) Re p<1000
0.44 Re p≥1000
(5)
Re p=ρgεg
|
|u g-u p,i d p
μg(6) 1.2 气相控制方程
气相视为连续介质相，连续性方程和动量方程可分别如式(7)、式(8)所示[15]。

∂()εgρg
∂t+∇⋅
(ε
gρg u g)=0(7)∂()
εgρg u g
∂t+∇⋅
(ε
gρg u g u g)=-εg∇P g+∇⋅(εgτg)+
εgρg g-1V cell∑i=1n F d,i(8)其中，气相密度ρg通过求解理想气体状态方程获得。

气相剪切应力张量τg表示为式(9)。

·
·2829
化工进展， 2023， 42（6）
τg =μg {
}
éëùû∇u g +()∇u g T -23
(
)∇⋅u g I
(9)
气相组分输运方程如式(10)所示。

∂
∂t (
εg ρg X g,j
)
+∇⋅(
εg ρg X g,j u g )
=∇⋅(
εg ρg D g,j ∇X g,j )
(10)
式中，D g,j 表示气相组分j 的扩散系数。

1.3 气泡-乳化相传质系数计算方法
对流与扩散传质系数的统计方法如图1所示[9]。

假设气泡内示踪气体浓度C j ,b 高于床层内气体浓度C j ,e ，在气泡与床层组分交换中，一方面，床层内间隙气流从气泡底部流入并推动气泡内气体于顶部流出，其对应于对流传质贡献；另一方面，由于气泡相与床层内气体浓度差异的存在，在浓度梯度的驱动下，气相组分会由高浓度的气泡相扩散至低浓度的乳化相，这对应于扩散传质贡献。

在本文中，设定空隙率阈值为0.75以识别气泡
边界[11]，气泡边界被一系列微小控制体（CFD 网格）离散划分。

对于每个微小控制体，可以通过网格内数据计算获得对流与扩散传质分量，将三个方向上的传质分量沿气泡边界积分可获得对流与扩散传质质量流率，表示如式(11)、式(12)所示[9]。

V be,con =
∫Ωεg
[γx (u
g -u b )A x +γy (v g -v b )
A y +
γz
(w g
-w b
)A z
]⋅C
j
(11)
V be,dif =
∫
Ωεg éë
êê-γx D g
∂c j ∂x A x -γy D g ∂c j ∂y A y -γz D g ∂c j ∂z A z ù
ûúú⋅C tot (12)
式中，Ωb 表示气泡边界；u b 、v b 、w b 表示气泡
在x 、y 、z 三个方向上的速度分量，可通过计算气泡形心位移随时间的变化速率获得；A 表示控制体
在某一方向上的面积；C j 和c j 分别表示组分j 的物质的量浓度和物质的量分数；C tot 表示气相组分的总
物质的量浓度；γ是用于确定传质方向的单位向量，在本文中，γ的正和负分别为气体流入和流出气泡的方向。

传质系数可表示为负的质量流率与气泡体积和气泡相与乳化相浓度差之比，传质系数分别计算如式(13)~式(15)所示。

K be,con =-V be,con
()
C j ,b -C j ,e V b (13)K be ,dif =
-V be,dif
(
)
C j ,b -C j ,e V b
(14)K be =K be,con +K be,dif
(15)
2 模拟工况设置
本文采用的三维射流鼓泡床模型如图2所示。

鼓泡床为圆柱形几何结构，直径为14.4mm ，高度为28mm ，初始堆积颗粒高度为14.4mm 。

床底部中心设有0.8mm×0.8mm 的矩形入射孔口，用于注射高速示踪气体丙烷（C 3H 8），流化气体采用氮
气（N 2）。

对于单气泡射流，其操作步骤如下。

（1）床层底部全部通入最小流化速度N 2，床内形成初始流化状态；
（2）示踪气体C 3H 8在入射时间t inj1内以入射气
流速度u inj 从入射孔口注入床层，入口除孔口区域仍通以最小流化速度的N 2；
（3）入射时间t inj1结束，仅关闭孔口气流，形
成的气泡结构在初始流化床内自由上升。

图2
　三维射流鼓泡床模型示意图
图1　传质系数统计方法示意图
·
·2830
2023年6月章凯等：鼓泡流态化气泡间相互作用下相间传质过程的模拟
双气泡的形成是在单气泡操作基础上注入第二个气泡。

（a）第一个气泡在入射结束后Δt inj时间内自由
上升，随后在入射时间t inj2内以入射气流速度u inj从
入射孔口注入示踪气体C3H8，底部其余部分仍然通以最小流化速度的N2；
（b）入射时间t inj2结束，第二个气泡结构形成，
仅关闭孔口气流，形成垂直分布的双气泡结构在初始流化床内自由上升。

模拟基于开源软件平台MFIX开展[12,16]，采用速
度入口和压力出口边界条件，壁面设为无滑移壁面，具体参数如表1所示。

3 结果与讨论
3.1 模型验证
基于CFD-DEM方法对气泡相间传质行为的模拟已得到验证[11]。

在本文作者课题组先前的工作中利用CFD-DEM方法对单气泡流动行为展开了模拟[17]，并与Boyce等[18]的实验结果展开对比。

图3给出了气泡在不同位置处形状参数的定量结果，气泡形状参数定义为气泡的最大横向尺寸W b与最大轴向尺寸H b之比，从图中可以发现，气泡在上升过程中逐渐横向拉伸而趋于扁平，模拟与实验结果基本吻合。

这里为了进一步验证模型对三维气泡的预测结果，表2给出了不同入射时间下气泡上升速度模拟和模型计算结果的对比，从表中可以看出，模拟和模型预测结果的相对误差均在5%以内，说明本文构建的CFD-DEM模型能够较为准确地描述三维气泡的运动行为。

3.2 双气泡与单气泡传质结果对比
为了分析双气泡与单气泡上升过程中气泡结构的动态演变，图4给出了相同入射时间与气流入射速度下，不同时刻床内固相体积分数与三维气泡结构图。

从图中可以看出，单气泡射流后气泡在上升过程中更倾向于形成球形结构，但在双气泡射流中，下方气泡在上方气泡影响下被轴向拉伸，气泡在上升过程中逐渐演变为子弹状，而后顶部趋于扁平，整体接近圆柱形。

表1　模拟主要参数
参数
颗粒直径/mm
颗粒密度/kg·m-3
弹簧常数/N·m-1
摩擦系数
恢复系数
最小流化速度/m·s-1
入射气流速度/m·s-1 N2/C3H8（298K）密度/kg·m-3孔口尺寸/mm×mm
网格尺寸/mm×mm×mm
数值
0.125
2525
25
0.1
0.9
0.013
4.0
1.145/1.803
0.8×0.8
0.4×0.4×0.4
表2　不同入射时间下气泡上升速度模拟值与理论值对比[7]
入射时间/ms
16
20
24
模拟值/m∙s-1
0.137
0.139
0.143
理论值/m∙s-1
0.14
0.146
0.15
相对误差/%
2.1
4.8
4.7
图3　气泡形状参数模拟与实验结果[18]
对比图4　单气泡与双气泡结构动态演变过程
·
·2831。