射流环流反应器中颗粒悬浮的数值模拟

庄黎伟，戴干策

(华东理工大学化学工程联合国家重点实验室，上海 200237)

射流环流反应器是由射流驱动形成环流的多相反应器，具备较好的气液分散和传热传质性能。根据不同的工艺要求，可以有不同的结构型式[1-4]，许多研究者进行了相关研究工作。Blenke 等[8]通过实验，作出了环流雷诺数和射流雷诺数之比、导流筒和反应器直径比以及总阻力系数三者的关系曲线。Tebel 等[6,7]提出了代表性粘度、代表性雷诺数概念并将其应用到假塑性流体中，用来描述环流反应器中流体在低雷诺数下的环流速度和射流速度之比随雷诺数的变化规律，补充了Blenke 等的研究工作。Tebel 等[2,8]研究了气液固三相射流环流反应器的流体力学性质。Mathpati 等[1]通过模拟和粒子图像测速法比较了不同湍流模型在湍流圆射流模拟中的准确性，指出标准κ-ε以及已有的改进κ-ε模型在预测产生项和输运项方面的缺陷，并提出了混合κ-ε模型，该模型在湍流圆射流的模拟中与雷诺应力模型（RSM）和大涡模拟（LES）具有可比拟的准确性，但在运算时间上未体现出优势。Szafran 等[9]运用标准κ-ε结合双流体模型研究了气液固三相射流环流反应器的流体力学特性，并与实验结果作比较，在气液两相和气液固三相体系中模拟预测的液体环流速度的误差分别为4%和15%，但是气含率的预测值与实验值相差最高达270%，说明多相流模型的选择有缺陷。

已有文献中鲜有关于射流环流反应器内颗粒悬浮模拟研究的报道，为确定射流所能驱动的环流速度的大小，本工作将采用离散相模型（DPM）模拟射流环流反应器中的颗粒悬浮。结合液相流场和离散相模型，通过预测颗粒在反应器内运动轨迹来得到颗粒的悬浮情况，并根据催化剂颗粒的悬浮情况实现催化剂粒径和射流速度的选择，并进一步优化设备结构。

1 物理模型和数值计算方法

1.1 控制方程

对于定态不可压缩流动，连续性方程和动量方程为

式中雷诺应力项需要借助模型来使方程组封闭。双方程k-ε模型已成功应用于平面射流，但是对于圆射流的发展速率和湍流参数的预测值偏高[1]。RSM 湍流模型基于各项异性假设，克服了涡粘性模型的缺点，具有更强的描述复杂湍流的能力。雷诺应力输运方程为

左端两项分别为雷诺应力随时间的变化率和对流项，右端分别为扩散项、应力产生项、浮力产生项、压力应变项、粘性耗散项、系统旋转产生项和自定义源项，经模型化后封闭雷诺应力方程组。

1.2 离散相模型

催化剂体积分数小于12%，忽略颗粒与颗粒之间的相互作用以及颗粒体积分数对连续相的的影响。在拉格朗日框架下模拟颗粒运动及其轨迹，颗粒在已求解的连续相流场中受力运动方程为

其中下标x代表笛卡尔坐标系下的x方向，右端第二项为单位质量催化剂颗粒所受重力与浮力之差，FD（u-up）为单位质量催化剂颗粒所受曳力

其中相对雷诺数

假设催化剂是直径为dp的球形颗粒，曳力系数CD可以表示为

其中a1，a2，a3为模型常数。由于ρp大于ρ无需考虑虚拟质量力，同时dp远大于1 μm 布朗力等微颗粒受力可忽略，所以额外加速项Fx为零。

在湍流场中，通过连续相的平均速度u来预测催化剂颗粒的运动轨道，同时用脉动速度u′来描述湍流对催化剂颗粒的分散作用。

1.3 模型建立及数值计算

射流环流反应器的结构型式如图1所示，液体射流从底部中心进入导流筒内，固体催化剂依靠环流在体系内悬浮，实现连续操作。设备筒体高度H为3.32 m，直径D为1.2 m，导流筒高度Hd为2.5 m，直径De为0.8 m，导流筒两端离设备距离Ht和Hb分别为0.42 m 和0.4 m，下端距离射流喷嘴口Hn为0.072 m，入口处直径dn为0.022 m。催化剂从液体入口处加料。根据设备的几何和安装特点，用二维轴对称模型对几何控制体进行简化。划分网格时，对射流出口以及射流发展区域进行局部加密，在设备壁面处保证第一层网格中心处于幂律区（即30＜Y+＜300），导流筒和环隙区域采用结构化网格，两端采用非结构化网格。为验证网格无关性，进行了不同数量网格的划分。

采用商业软件Fluent，RSM 湍流模型进行单相模拟，模拟物料为20 ℃纯水，进出口分别采用速度入口和压力出口边界条件，湍流参数根据水力直径和湍流强度给出，近壁面区采用标准壁面函数，固体壁面无滑移条件，中心线采用轴线条件，压力速度为PISO 耦合方式，压力项采用二阶离散格式，其余采用Quick 离散格式，动量、湍流动能、耗散率的亚松弛因子分别为0.3、0.5、0.5，其余为默认值。连续性方程残差为10-3，其余为10-4，同时监控环隙循环流量大小以判断收敛。模拟了不同射流速度Vn下反应器内循环流量的变化，以及导流筒直径，设备底部形式（见图2）和射流喷头长度的影响。

在得到单相流场后，用非稳态DPM 模型模拟不同粒径dp催化剂颗粒的悬浮。进出口均采用逃逸条件，轴线及壁面均为反射，催化剂进口质量流率为1.67 kg/s，进料时间持续6.5 s，追踪约3 900 个催化剂颗粒，进料结束颗粒体积分数约为0.1%，颗粒密度为3 000 kg/m3。

图1 射流环流反应器Fig.1 Schematic diagram of jet loop reactor

图2 锥形底（左）和流线底（右）Fig.2 Cone bottom（left）and streamline bottom（right）

2 结果与讨论

2.1 模型验证

为了消除模拟结果对于网格数目的依赖性，分别在3 种网格数下进行了参数灵敏性的检验。在离射流口轴向距离yn为322 mm 的位置上进行数据采集，对比结果如图3所示。由图可知，3 种网格密度下，湍流动能与速度大小差别均在10%以下，根据检验结果，采用1.2×105的网格划分方式进行模拟，以兼顾计算效率和准确性。

为了考察射流速度径向分布情况，模拟了射流速度4.38 m/s 时，不同轴向距离x处速度的径向分布，以及不同射流速度下，922 mm 处速度的径向分布，并分别与理论曲线进行对比，结果如图4和图5所示。在本研究的射流环流反应器中，不同轴向距离yn（yn＞30dn）上速度的径向分布具有自相似特性[10,11]。在湍流圆射流（本工作所有工况中最低射流雷诺数约为96 500）下有以下经验公式[12]：

图3 湍流动能和速度的径向分布Fig.3 Radial distribution of turbulent kinetic energy(a) and velocity(b)

其中，U为某yn位置上的速度大小，Umax为射流中心速度，r1/2为速度为射流中心速度一半时的径向位置。

由图4、图5可知，不同yn位置上以及不同射流速度的相同yn位置上，无量纲速度U/Umax和无量纲径向距离r/r1/2均落在同一条曲线上，除了在较大径向距离上出现了偏差，原因可能是射流处于导流筒限制的有限空间内，与自由射流有所区别（De＞50dn），同时，从图4中r/R在0.1～0.35 出现的平台区可以看出射流环流反应器导流筒内射流的特点，射流外缘的环流与射流在边界处进行动量交换，当yn值较小时，动量交换还未完成，随着射流的发展（yn＞30dn），射流才进入自相似阶段。由以上两个原因造成了在较大径向位置处模拟结果与理论值的偏差，在若干文献中[12-14]也出现了类似情况。

图4 不同yn 处速度的径向分布Fig.4 Radial distribution of velocity at different x

图5 不同射流速度下的速度径向分布Fig.5 Radial distribution of velocity at different Vn

由以上模型验证的结果可以得出，数值计算合理。

2.2 循环流量与射流流量之比

在射流环流反应器中，环流沿整个流道的阻力损失由射流补偿，以维持环流的运转，根据动量守恒原理[8]可得循环流量与射流流量比：

Tebel 等[6,7]的研究认为，随着雷诺数增加到一定值，阻力系数ζs 仅与环流反应器的几何结构有关，由式（10）可知，在一定的几何结构下，RQ为常数。RQ的大小反映了在特定设备中射流驱动环流的能力。将Blenke 等[5]得出的实验拟合曲线用于本研究的设备，根据De/D得出ζs 的值，由此作出RQ与Vn的关系曲线，并与模拟结果进行对比。

2.2.1 导流筒与设备直径比对RQ的影响

不同导流筒与设备直径比下的总阻力系数ζs 会有不同，根据文献[5]将两者的对应关系列于表1中。

表1 De/D 与ζs 的对应关系Table 1 The corresponding relation of De/D and ζs

由（10）式得出不同De/D下的RQ，与6 种射流速度下的模拟值进行比较，结果如图6所示。由图6可知，模拟值与理论值相差在10%以内，且当De/D为0.67 时，相同射流流量所能驱动的循环流量最大，此时的总阻力系数比De/D为0.5 时低了44.4%。由此可知，在保持其他结构参数不变的情况下，De/D为0.67 左右时，射流驱动环流的能力相对较强，可见合理选取De/D的值有助于降低能耗。以下均采用De/D为0.67 的情况进行研究。

图6 不同射流速度下循环流量与射流流量比Fig.6 Flow rate ratio of circulation to jet at different jet velocity

2.2.2 设备底部型式对RQ的影响

由于流体在设备底部发生流向的转变，且催化剂颗粒的沉积发生在底部，因此，有必要选择合适的底部型式以降低能耗和防止催化剂沉积。原有设备采用碟形底（见图1），考虑到流线型设计可能有助于降低局部阻力，于是设计了一种流线底（见图2右），比较两种型式的底部结构对RQ的影响来探究合理的底部型式。根据两种底部结构的模拟情况（见图6(d)），所述两种底部型式，射流驱动环流的能力无显著差别。

2.2.3 射流喷口与导流筒下缘的距离对RQ的影响

由2.2.2 可知，在一定结构尺寸下，设备底部型式对于RQ的影响很小，所以，采用图2(b)的底部型式，Hn为400 mm，同Hn为72 mm 的碟形底进行比较，研究Hn对RQ的影响，模拟所得结果如图6(e)所示。由图6(e)可知，两种Hn值下模拟得出的RQ值相差低于3.6%，可认为在给定范围内，Hn对于射流驱动环流的能力的影响可忽略。

2.3 催化剂颗粒的悬浮

2.3.1 液相流场

催化剂颗粒的悬浮决定于射流环流反应器内的流场。图7(a)是Vn为10 m/s 时反应器内的速度矢量图，箭头表示速度方向。反应器底部高速射流驱动流体环流运动，导流筒内流体向上流动，到达顶部，一部分流体从顶部出口离开，大部分绕过导流筒上缘偏折进入环隙向下运动，重新到达底部与射流汇合，形成一个环流周期。射流卷吸使导流筒内壁附近形成了漩涡，又由于导流筒的限制，漩涡呈现狭长的外形。催化剂借助射流和环流在导流筒内完全悬浮。图7(b)是Vn为10 m/s 时反应器内的湍流动能分布云图。为了体现湍流动能的区分度，显示范围最大值为0.1 m2/s2，空白区域表示湍流动能超过0.1 m2/s2。高速射流从喷咀处向上喷射，形成强剪切，湍流动能最大，随着射流的发展，其与周边流体发生动量交换，湍流动能沿径向分布逐渐均匀。湍流脉动促进催化剂颗粒的混合和分散。

图7 速度矢量和湍流动能的分布Fig.7 Contours of velocity vectors (a) and turbulent kinetic energy (b)

2.3.2 底部几何构型对催化剂悬浮的影响

选取图1和图2左两种底部几何构型（包括底部型式和喷咀高度），在6 种Vn下，考察了16 种粒径的催化剂颗粒在设备内的悬浮情况。图8～11 为不同Vn下，各种粒径催化剂颗粒在设备内的质量浓度分布。为了体现浓度的区分度，质量浓度显示为0～10 kg/m3，空白区域表示质量浓度超过10 kg/m3，即该区域颗粒浓度相对较高。所有颗粒浓度分布云图均对应某一时刻（颗粒投放结束之后1 000 s，至少17 个环流周期），此时颗粒浓度分布趋于稳定，且反应器内剩余颗粒数目均在2 000 以上。

由图8可以看出，催化剂从底部喷射出约1 000 s 之后，反应器内部颗粒浓度分布趋于稳定。在非稳态追踪过程中，大部分0.05 mm 的颗粒随着流体进行环流运动，由于粒径较小，有较好的流动跟随性，但底部出现催化剂颗粒小范围的聚集，虽未沉降，但悬浮高度较低；0.1 mm 的催化剂的悬浮情况与0.05 mm 的差别主要在于颗粒分布的均匀性，从图中可以看出，在射流出口两侧的局部环流和反应器底部，颗粒浓度明显高于其他位置；当dp为0.3 mm 时，大量催化剂颗粒在反应器底部沉降，只有小部分颗粒由于射流的喷动和卷吸作用，悬浮于导流筒内部，且悬浮高度低于导流筒上缘；大部分0.5 mm 的催化剂颗粒已经沉降于反应器底部，无法再次悬浮。

图8 不同粒径催化剂颗粒的质量浓度分布Fig.8 Concentration distribution of particles with different diameter (Vn=4.38 m/s, plate bottom)

图9 不同粒径催化剂颗粒的质量浓度分布Fig.9 Concentration distribution of particles with different diameter (Vn=4.38 m/s, cone bottom)

通过比较图8和9 可以看出，相同粒径催化剂颗粒的悬浮情况，低喷咀高度的锥形底明显好于高喷头高度碟形底。以dp为0.5 mm 为例，图9中催化剂颗粒虽然悬浮高度较低，且悬浮均匀性较差，但在实际生产过程中，可以通过调节射流速度来提高悬浮高度和均匀性，而图8中催化剂颗粒过早积底，悬浮终止。由于在完全悬浮状态下，部分催化剂颗粒进入局部低速区或者与壁面发生碰撞，不可避免的发生沉降，图8中的底部型式不利于沉降催化剂的再次悬浮，而图9的底部型式却克服了这一缺点，当催化剂颗粒沉降至底部，在锥形底部的坡度诱导和环流冲刷的双重作用下，催化剂颗粒再次汇聚到射流喷咀处，高速射流再次将催化剂喷向导流筒内，实现悬浮。可见，合理的底部几何构型有助于消除催化剂的积底，改善催化剂的悬浮状况。所以接下来设备将采用低喷头高度的锥形底的底部几何构型，在此基础上探究催化剂粒径dp和环流速度Ve对催化剂颗粒悬浮的影响。

2.3.3 粒径和环流速度对催化剂悬浮的影响

在一定的设备结构下，RQ为定值，射流速度保持不变，环流速度也不变。由图9可以看出，在射流速度为4.38 m/s 时，催化剂随着粒径dp的减小，悬浮高度逐渐增加，从0.5 mm 的催化剂悬浮至导流筒2/3 高度处到0.05 mm 的催化剂均匀悬浮于整个设备。相对于大粒径的催化剂颗粒，较小粒径的催化剂充分提高了设备体积的利用率，增大了发生反应的区域，在实际生产中有积极意义。所以，为了提高设备利用率，根据图9，当射流速度为4.38 m/s 时，应采用粒径在0.3 mm 以下的催化剂。

环流速度Ve和催化剂粒径dp对催化剂颗粒悬浮的影响在本质上是一样的，提高环流速度和降低催化剂粒径都提高了催化剂对流体作用力的响应。通过图10和图11的对比可知，在高环流速度下，相同粒径的催化剂颗粒悬浮高度越高或者越能充满整个设备，设备的体积利用率越高，同时悬浮的均匀性也相应提高。当射流速度从4.38 m/s 提高至16 m/s 时，dp为0.5 mm 的催化剂已经几乎可以充满整个设备，继续提高环流速度对于催化剂悬浮状况的改善意义不大。

图10 不同粒径催化剂颗粒的质量浓度分布Fig.10 Concentration distribution of particles with different diameter (Vn=10 m/s, cone bottom)

图11 不同粒径催化剂颗粒的质量浓度分布Fig.11 Concentration distribution of particles with different diameter (Vn=16 m/s, cone bottom)

所以，在实际生产过程中，由于催化剂颗粒不可避免的存在粒径上的差异，环流速度越高能使更多催化剂悬浮，提高设备的体积利用率，但环流速度的提高则会造成较高的能耗。本研究中催化剂粒径不大于0.5 mm 时，10 m/s 的射流速度则可以基本满足悬浮要求。

2.3.4 临界射流速度

由2.3.3 可知，为了提高设备的体积利用率，应使催化剂颗粒完全悬浮于整个设备，在此悬浮状态下，存在Vn（或Ve）和dp的对应关系，即使一定粒径催化剂完全悬浮的最低Vnc。Tebel 等[2]通过动量守恒定律和导流筒内外物料守恒得出了射流环流反应器中液固两相操作射流速度的表达式，他们认为，存在一个最大的射流速度，此时催化剂颗粒均匀悬浮于导流筒内部，环隙无颗粒，继续提高射流速度，催化剂颗粒将越过导流筒进入环隙，作环流运动。他们通过实验验证了该表达式

其中uss(ce)为ce浓度下的颗粒群最终沉降速度，ce和c分别为导流筒内颗粒体积浓度和整个设备颗粒平均体积浓度，Δρ为催化剂颗粒和连续相的密度差。在本研究中，c约为0.001，ρm约为ρ，uss(ce)近似为ups，根据无量纲数K判别颗粒周围流体的流型（层流，过渡，紊流）[15]来得到最终沉降速度

在本工作中，当颗粒直径dp介于0.15 mm 和2 mm 时，流型处于过渡区，此时最终沉降速度

层流和紊流流型下最终沉降速度表达式见文献[15]。

理论上Vnc即为Vnmax。通过模拟6 种射流速度下16 种粒径催化剂的悬浮，并比较了使同一粒径催化剂均匀悬浮于导流筒内临界速度的模拟值和理论值，结果见图12。由图12可知，本工作的射流环流反应器，在低颗粒体积浓度（约为0.1%）下，密度ρp为3 000 kg/m3的催化剂在20 ℃纯水中完全悬浮于导流筒内部，所需的临界射流速度模拟值和理论值偏差在（dp＞0.6 mm）10%以内，而当颗粒较小时，偏差较大（最高达64%）。

对于图（12）中粒径较小时模拟值和理论值的偏差较大（最高达64%）的原因，作者认为可能是式（10）不适用于预测粒径较小的催化剂的悬浮，由式（12）可知，颗粒最终沉降速度对颗粒直径敏感，而射流环流反应器中颗粒悬浮的本质就是通过环流克服颗粒的沉降，临界射流速度也应该对颗粒敏感，这与本研究的模拟结果一致（见图12），而且，通过计算，当粒径足够小时，保证颗粒体积浓度不变，由式（10）得出临界射流速度维持在12 m/s，这有悖于颗粒悬浮的机理。因此，根据模拟得到dp和Vnc的对应关系，同时考虑到颗粒悬浮需要克服沉降以及射流环流反应器内颗粒需要完全悬浮于导流筒内，将dp和Vnc的关系以无量纲数的形式表示，经过数据拟合得到新的关联式（13）。该关联式适用于与本研究结构类似的射流环流反应器，催化剂颗粒体积浓度不大于0.1%，粒径为0.2～1 mm。在设计此类射流环流反应器时，可以参照式（13）进行催化剂粒径和射流速度的选择。

图12 临界射流速度Fig.12 The critical jet velocity Vnc

3 结 论

a）对给定的设备结构，在环流雷诺数Rem为8×104～3×105时，RQ为常数；在保持其他结构参数不变情况下，导流筒与设备的直径比De/D取0.67 左右时，RQ最大；Hn对RQ的影响可忽略。

b）随颗粒粒径增大，为使颗粒悬浮于整个设备，所需环流速度几乎线性增大。根据模拟结果拟合出关联式（13），可作为此类射流环流反应器设计的参考。

c）具有一定坡度的设备底部设计，如锥形底，有助于催化剂颗粒在底部向中心汇聚，同时降低Hn的值，能显著改善催化剂颗粒的悬浮状况，消除催化剂积底的影响。高射流喷咀高度的碟形底的设备中，增大颗粒粒径（dp=0.5 mm），催化剂完全积底，无法再次悬浮，而在喷咀高度Hn较低的锥形底设备中，考察的所有粒径催化剂均无积底影响。

综上所述，选用雷诺应力湍流模型计算流场，结合离散相模型，模拟估算射流环流反应器中颗粒的悬浮是可行的。

符号说明

c——颗粒体积浓度

CD——曳力系数

dn——射流喷头直径，m

dp——颗粒直径，m

D——反应器直径，m

De——导流筒直径，m

gx——当地重力加速度，m/s2

H——反应器高度，m

Hb——导流筒下端与设备底端距离，m

Hd——导流筒高度，m

Hn——射流入口与导流筒下端的距离，m

Ht——导流筒上端与设备顶端距离，m

K——流型判别无因次数

P——时均压力，Pa

Qe——环流流量，m3/s

Qn——射流流量，m3/s

u——流体速度，m/s

u′——脉动速度

ui——时均速度分量

up——催化剂颗粒运动速度，m/s

ups——单颗粒最终沉降速度，m/s

uss——颗粒群最终沉降速度，m/s

Ve——环流速度，m/s

Vn——射流速度，m/s

Vnc——临界射流速度，m/s

xi——笛卡尔坐标分量

Y+——无量纲距离

yn——距喷咀出口轴向距离，mm

ε——湍流动能耗散率，m2/s3

ζs——总阻力系数

η——最终沉降速度纠正系数

κ——湍流动能，m2/s2

μ——流体分子粘度，Pa·S

ρ——流体密度，kg/m3

ρm——悬浮液密度，kg/m3

ρp——催化剂颗粒密度，kg/m3

[1]Mathpati C S, Deshpande S S, Joshi J B.Computational and experimental fluid dynamics of jet loop reactor [J].AIChE Journal, 2009,55(10):2526-2544.

[2]Tebel K H, Zehner P.Fluid dynamic description of jet-loop reactors in multiphase operation [J].Chemical Engineering &Technology,1989, 12(1):274-280.

[3]Tinge J T, Rodriguez Casado A J.Influence of pressure on the gas hold-up of aqueous activated carbon slurries in a down flow jet loop reactor [J].Chemical Engineering Science, 2002, 57(17):3575-3580.

[4]Havelka P, Lineket V, Sinkule J, et al.Hydrodynamic and mass transfer characteristics of ejector loop reactors [J].Chemical Engineering Science, 2000, 55(3):535-549.

[5]Blenke H, Bohner K, Schuster S.Beitrag zur optimalen Gestaltung chemischer Reaktoren [J].Chemie Ingenieur Technik, 1965,37(3):289-294.

[6]Tebel K H, Zehner P, Langer G, et al.Homogenisieren strukturviskoser Flüssigkeiten in Schlaufenreaktoren und Rührkesseln [J].Chemie Ingenieur Technik, 1986, 58(10):820-821.

[7]Tebel K H, Zehner P.Ein Konzept der repräsentativen Viskosität für Mischprozesse [J].Chemie Ingenieur Technik, 1985, 57(1):49-51.

[8]Zehner P.Suspendieren von Feststoffen im Strahlschlaufenreaktor [J].Chemie Ingenieur Technik, 1980, 52(11):910-911.

[9]Szafran R G, Kmiec A.Application of CFD modelling technique in engineering calculations of three-phase flow hydrodynamics in a jet-loop reactor [J].International Journal of Chemical Reactor Engineering, 2004, 2(2):1179.

[10]Hussein H J, Capp S P, George W K.Velocity measurements in a high-Reynolds-number, momentum-conserving, axisymmetric,turbulent jet [J].Journal of Fluid Mechanics, 1994, 258(1):31-75.

[11]Panchapakesan N R, Lumley J L.Turbulence measurements in axisymmetric jets of air and helium.Part 1.Air jet [J].Journal of Fluid Mechanics, 1993, 246(1):197-223.

[12]Schlichting H, Kestin J.Boundary-layer theory [M].McGraw-Hill New York, 1968:1-817.

[13]Wasewar K L, Sarathi J V.CFD modelling and simulation of jet mixed tanks [J].Engineering Applications of Computational Fluid Mechanics, 2008, 2(2):155-171.

[14]Jayanti S.Hydrodynamics of jet mixing in vessels [J].Chemical Engineering Science, 2001, 56(1):193-210.

[15]时 均, 汪家鼎, 余国琮, 等.化学工程手册 [M].北京:化学工业出版社, 1996:8-10.