刘畅,李树森,张立峰
北京科技大学冶金与生态工程学院 北京 100083
摘要
利用水模型与数值模拟相结合的方法,对Ruhrstahl-Heraeus (RH)精炼过程中的多相流体流动及混匀现象进行模拟研究。根据相似原理建立与实际210 t RH精炼装置几何相似比为1：5的水模型,采用粒子图像测速(PIV)技术获取水模型中心截面处流场分布。数值模拟采用多相流模型(VOF)和离散相模型(DPM)相耦合的计算方法,湍流模型分别选用k-ε模型和大涡模拟(LES)模型。对比测量值与计算值,结果表明2种湍流模型均能较好地预测RH内流场分布;而采用LES模拟能够获得RH内瞬态的速度分布及漩涡的产生和耗散过程。测量并计算了水模型钢包内整体的混匀时间分布,结果表明上升管附近区域的混匀时间大于下降管附近区域的混匀时间。开发了气泡膨胀的数学模型,并将其用于钢液-氩气体系的模拟计算,结果表明气泡膨胀过程对钢液流动的影响显著。
关键词：
二次精炼过程作为冶炼过程的重要环节,对改善后续工序操作和提高成品钢的质量有显著作用。常见的精炼方法有钢包吹氩、喂线、真空脱碳法(VD)、钢包炉精炼法(LF)和循环脱气法(RH)等,一般具有搅拌、加合金、调节成分等作用。RH精炼装置具有脱气、脱碳、均匀温度与成分、促进夹杂物上浮等作用[1],因此被广泛用于超低碳钢、硅钢等钢种的生产[2]。其内部流体流动现象对于精炼效率有显著影响,因此对RH内部流场的研究具有重要意义[3,4,5]。但由于工业实验具有成本高、危险性高等特点,许多学者采用物理模拟或数值模拟的方法研究冶金反应器内的流动、传热、反应及夹杂物行为等现象。
对于RH真空精炼过程中多相流动及传输现象的模拟研究,以往学者多采用将钢液处理为连续相而气泡为离散项[6,7,8],或钢液相与气相均处理为连续相[9~11]2种形式。RH精炼过程中,Ar气作为驱动气体,与钢液的交互作用对整个RH反应器内的流体流动具有重要影响[11,12,13]。然而,在实际冶炼的钢液-氩气体系中,气泡在上升管内受热膨胀以及碰撞破碎聚合过程极其复杂,目前对吹入的气相在上升管内膨胀过程的模拟鲜见报道[14,15]。Szekely和Fang[16]首次给出了真空条件下气泡在液相中上浮过程的气泡直径变化的微分方程。之后Chen等[7]研究表明,气泡膨胀过程对流场有较大影响。
对于湍流模型的选择,以往研究中多采用Reynolds平均N-S方程的方法(RANS)进行求解[8,17,18];例如k-ε模型,由于模型假设湍流为各向同性,因而难以对强湍流状态下的流体的涡流现象进行描述。大涡模拟(large eddy simulation,LES)只对大尺寸脉动运用数值计算方法直接求解,通过过滤方法消除湍流中小尺度脉动,并将其对可解尺度的脉动的作用通过亚格子应力模型来求解,具有较好的适用性[19]。
本研究采用物理模拟与数值模拟相结合的方法,研究了RH内部流场及混匀现象。根据Fr准数相似原理建立水模型[20],利用粒子图像测速(particle image velocimetry,PIV)系统实现对水模型内瞬时流场的非接触式测量,获取速度矢量及湍流特性的分布。建立多相流(volume of fluid,VOF)模型和离散相(discrete phase model,DPM)模型相耦合的三维数学模型,对RH精炼过程进行模拟,在此基础上对比了2种湍流模型。在计算稳态流场的基础上,增加示踪剂浓度方程的求解,获得钢包内整体的混匀时间分布。在经验证数学模型的基础上,提出了一种新的数学模型,考虑Ar气泡在钢液中的膨胀过程,在计算过程中对气泡直径进行修正。
VOF模型通过求解一组动量方程来模拟2种及2种以上互不相溶的流体,并在整个计算区域内追踪各相的体积分数,用于获得相间界面轮廓。
连续性方程：
式中,下标q为相数,αq为q相的体积分数;ρq为q相的密度,kg/m3;vq为q相的速度,m/s。
动量方程：
式中,v为混合相的速度,m/s;P为压力,Pa;F为力的源项,N/m3;α为体积分数;ρ为密度,kg/m3;μ为黏度,kg/(ms);g为重力加速度,m/s2;下标l、g分别表示液相和气相。
气泡在流体中运动受力包括曳力、重力、浮力、虚拟质量力、压力梯度力以及提升力[21],其运动方程可以表示为：
式中,FD、FG、FB、FVM、FP和FL分别为曳力、重力、浮力、虚拟质量力、压力梯度力和提升力,m/s2;Re为Reynolds数;db为气泡直径,m;CD、CVM、CL分别为曳力系数、虚拟质量力系数、升力系数;下标b表示气泡相。
气泡初始直径为与单个吹气孔的吹气流量有关的参量[22,23]：
式中,db,0为气泡初始直径,m;Qn为单个吹气孔的吹气流量,m3/s;σ为气泡与流体间的表面张力,N/m;dn为单个吹气孔的直径,m。
考虑气泡在钢液中被加热,且上升过程中所受静压力逐渐降低,因此气泡会迅速膨胀。根据温度和静压力的变化,修正后的气泡直径db计算公式为：
式中,Tl为钢液温度,K;T0为气体初始温度,K;P0为大气压,Pa;PV为真空室压力,Pa;Δh为气泡所处位置距离钢液面的高度,m。
气泡在上升管内存在复杂的破碎和聚合现象,本研究中未考虑破碎及聚合过程发生的详细机理[24,25],只关注气泡存在的最终状态,考虑流场中湍动能耗散率对气泡可能存在的最大直径的限制,计算公式[26]如下：
式中,db, max为气泡可能存在的最大直径,m;Wecr为临界Weber数;ε为湍动能耗散率,m2/s3。
采用标准k-ε方程描述气液两相流中的湍流特性,可直接求解出湍动能及其耗散率的分布,对应表达式分别为：
式中,k为湍流动能,m2/s2;μt为湍流黏度,kg/(ms);Gk为平均速度梯度引起的湍流动能,kg/(ms3);Sk为湍流动能的源项,kg/(ms3);Sε为湍动能耗散率的源项,kg/(ms3);常数C1ε=1.44,C2ε=1.92,σk=1.0,σε=1.3。
湍流模型选用大涡模型计算湍流涡黏性,为构造亚格子应力的封闭模式,本研究采用Smagorinsky-Lilly模型求解涡黏度,计算公式为：
式中,LS为亚格子的混合长度,m;
式中,κ为Karman常数;d为与最近壁面的距离,m;CS为Smagorinsky常数,本研究中取值0.1;Δ为所处位置的网格尺寸,Δ=V1/3;i、j为张量的自由标。
在计算稳态流场的基础上,加入一定量的示踪剂并求解其传输方程,表达式如下：
式中,C为示踪剂浓度;Sct为湍流Schmidt数,计算中取值为1[27];SC为源项。定义示踪剂浓度变化稳定在±5%[28]所需的时间为混匀时间tm。
钢包及真空室上表面边界条件均为压力出口,其中对于水-空气体系和钢液-氩气体系,真空室上表面压力分别为96523和100 Pa。采用k-ε模型计算时,离散相模型中选用随机游走模型。考虑到水模型体系压力变化较小,在空气-水体系的计算中未考虑气泡在上升管内的膨胀过程,即气泡直径在计算过程中为常数。图1给出了大涡模拟计算用网格,在拟定气液交界面处对网格进行了加密处理。考虑到大涡模拟计算结果与网格尺寸关系密切,综合考虑计算效率和准确性,最终选用网格总数约217万。
图1计算网格尺寸及边界条件
Fig.1Geometric model and mesh configuration
本研究根据几何相似与动力学相似,建立与原型为210 t的RH装置相对应的水模型,几何相似比为1∶5,模型尺寸如图2所示。
为保证水模型与原型之间传输现象的一致性,使钢液与水的Froude准数Fr相等,Reynolds数Re在同一自模化区：
式中,V为速度,m/s;L为特征长度,m。
实际体系与水模型的相关物性参数如表1所示。
图2水模型尺寸及监测点布置示意图
Fig.2Geometrical dimensions of experimental water model and monitors distribution in it (unit: mm)
表1RH水模型与原型物性参数对照表
Table 1Physial parameters of materials for water model and prototype
本研究通过PIV技术测量水模型内部流场。测量时,需在流场中散布适当密度且跟随性良好的示踪粒子,2束激光以固定的时间间隔照亮拍摄平面,同时利用63062型CCD (charge coupled device)相机等摄像设备获取示踪粒子的图像,对得到的图像序列进行互相关分析,得到二维速度矢量分布。实验所用激光发射器为Vlite 380脉冲固体激光器系统,图像分析软件采用InSight 4G。示踪粒子采用空心SiO2微球,粒径约10 μm,密度约0.99 kg/m3。
采用应激-响应方法来测量混匀时间,在钢包内共选取20个监测点,如图2所示,分布于钢包内中心纵截面上,上升管入口、下降管出口、上升管和下降管之间以及上升管与钢包壁之间不同深度位置。通过加入一定量的饱和NaCl溶液,同时监测钢包内不同位置的电导率变化,定义最终电导率变化不超过稳定值的±5%[28]所需的时间为混匀时间tm。实验用数据采集装置为DJ-800型多功能监测系统,采集频率10 s-1,采集时间6 min,5次平行实验求平均值。
对比相同工况下的水模型实验与数值计算结果,如图3和4所示真空室与钢包内的速度矢量,标尺中箭头长度与颜色均表示速度大小。三者的结果中,真空室内的流场稍有不同,采用LES模拟的情况下,液面处的速度值与测量值更为接近;三者的钢包内流体的流动形态结果大致相同,流体从下降管进入钢包中向底面运动的过程中速度逐渐衰减,撞击底部之后向四周扩散。
对比中心截面上2条特征线上的速度值,即钢包内下降管中心竖直线和钢包内距离底部500 mm水平线,如图5所示。2种计算模型与测量值吻合均较好,表明当前模型可以用于对RH精炼过程中流体流动的模拟。
图3水模型真空室内流场分布测量值与计算值对比
Fig.3Comparison of velocity distributions measured (a) and calculated by usingk-εmodel (b) and LES model (c) in vacuum chamber
图4水模型钢包内流场分布测量值与计算值对比
Fig.4Comparison of velocity distributions measured (a) and calculated by usingk-εmodel (b) and LES model (c) in ladle
k-ε模型采用Reynolds平均的方法,在求解过程中,假设3个方向的脉动速度相等,而实验测量中不同方向的脉动速度并不相等,测量结果如图6所示。因此,采用k-ε模型难以获取流体的瞬态湍流信息,其模型本身是具有一定局限性的。
图5中心截面特征线上计算与测量速度值对比
Fig.5Comparison of measured and calculated velocity magnitudes(a) along the horizontal line 500 mm from ladle bottom(b) along the center line of down-leg snorkel
图6水模型测量不同位置处脉动速度
Fig.6Fluctuating velocities for different positions in physical model (U’—Xcomponent of fluctuating velocity,V’—Ycomponent of fluctuating velocity)
进一步对比2种湍流模型,图7给出了速度为0.1 m/s的等值面图,由于k-ε模型对流场进行了平均处理,等值面分布较为平滑,而LES模型计算的瞬态分布较不均匀。两者均表明在真空室和钢包下降管区域附近的速度较大。
图8给出了LES模型的瞬态速度分布,矢量图中给出了流体在RH反应器中循环的过程,除在下降管与壁面间形成的完成整体循环的大涡流外,在流体运动过程中也存在小漩涡的生成和消失,整体速度分布较为不规则。从图8b可以看出,气泡在上升管内上浮速度大于所处位置钢液流速。
图7不同湍流模型速度等值面对比
Fig.7Iso-surfaces ofv=0.1 m/s fork-εmodel (a) and LES model (b)
图8大涡模拟瞬态速度分布
Fig.8Instantaneous velocity distributions for LES model(a)Y=0, velocity vector(b)Y=0, volume fraciton of gasαg=0.5, velocity contour
图9时均速度云图及气泡分布
Fig.9Averaged velocity distribution and bubble distribution calculated by usingk-εmodel (a) and LES model (b)
图9给出了2种计算模型下气泡的分布情况,如图9a所示,采用k-ε模型时气泡的分布在各方向上较为均匀,气泡在液相中的停留时间最长约1 s左右;采用LES模型时(图9b),不同位置吹气孔吹入的气泡运动轨迹不尽相同,且气泡在液相中的停留时间最长可达2 s,在真空室内的散布范围也较大。从液面波动情况的对比来看,k-ε模型的计算结果中液面平稳,难以观察到波动情况;LES模型的计算结果中液面波动较为明显。
混匀时间的测量实验中,示踪剂从真空室侧壁上的加料口加入,因此在数学模型计算中标记真空室内液相中心一定区域为示踪剂添加位置,如图10所示,2种湍流模型中的示踪剂加入位置保持一致;监测点布置与水模型保持一致,如图2所示。
图11a给出了实验测量4个监测点的电导率随时间变化的曲线。可以看出,示踪剂循环几个周期后,逐渐达到混匀状态。其中,位于下降管出口位置的监测点(Point 2-4)的电导率变化最为明显,几次出现峰值,之后达到稳定。图11中3种情况下得到的曲线变化趋势大致相同,但基于LES模型的计算结果稍滞后于k-ε模型。
将水模型实验测量的钢包中心截面上共20个点的混匀时间绘制成云图,如图12a所示。图12b和c分别为2个模型计算结果,计算值与测量结果的分布趋势大致相同。三者的结果均表明,钢包内不同位置处的混匀时间不尽相同;钢包内中心截面上,下降管附近的混匀时间整体上低于上升管附近的区域。
图10初始示踪剂分布
Fig.10Initial distribution of tracer in vacuum chamber (C0—initial concentration of tracer)
图11水模型监测点电导率或浓度变化曲线
Fig.11Variation of conductivities or tracer concentrations (C) of monitors measured (a) and calculated byk-εmodel (b) and LES model (c) (C∞—final trace concentration)
图12水模型混匀时间分布测量值与计算值对比
Fig.12Comparison of mixing time distributions measured (a) and calculated by usingk-εmodel (b) and LES model (c)
将上述经验证的模型用于计算实际尺寸的RH反应器中钢液-氩气体系的流动特性。Ar气与高温钢液接触后,其性质会迅速变化。如图13a和b所示分别为有无考虑气泡膨胀过程的情况下,上升管内气泡的分布。不考虑气泡膨胀过程时,气泡逐步汇聚在上升管中心。考虑膨胀过程,气泡上升过程中所处位置的静压力逐渐降低,气泡直径逐渐增大;当气泡上浮进入到真空室内,流体中的湍动能耗散率处于较高的水平,气泡直径受到限制;继续上浮直至脱离钢液进入气相。因此气泡在整个上浮过程中,气泡直径呈现先增大后减小的趋势。
图14给出了上述2种气泡直径分布下对应的速度分布。加入气泡膨胀模型后,RH装置内的整体流场明显强于未考虑气泡膨胀过程的情况,表明气泡膨胀现象对流场的影响较大。图14a中真空室内气液界面稍高于图14b,这是由于后者在真空室上表面的气泡直径较小造成的。气泡从液面脱离时的直径大小会影响液面的波动情况,直径越大,波动越大。
图13上升管内气泡直径分布
Fig.13Bubble diameter distributions without (a) and with (b) bubble expansion in up-leg snorkel (db—bubble diameter)
图14中心截面上速度分布
Fig.14Velocity distributions on the sectionY=0 without (a) and with (b) bubble expansion in up-leg snorkel
(1) LES模型能够更好地反应RH反应器内漩涡的产生和耗散。在流体循环的过程中,钢包中心截面上,下降管与钢包壁面间存在2个大的涡流。
(2) 获得了钢包内整体的混匀时间分布,表明混匀时间的数值与监测点的位置有关,上升管附近区域的混匀时间高于下降管附近的混匀时间。
(3) 气泡膨胀过程对RH反应器内的流体流动具有重要影响,考虑气泡膨胀过程的情况下,钢液循环速率更快。气泡从钢液钢液面处脱离时的直径越大,引起的钢液面的波动越明显。
致谢 感谢稀贵金属绿色回收与提取北京市重点实验室(GREM)和北京科技大学高品质钢研究中心(HQSC)及绿色冶金与冶金过程模拟仿真实验室(GPM2)的资助。
1 数学模型
1.1 多相流VOF模型
1.2 气泡运动DPM模型
1.3 湍流k-ε模型
1.4 湍流LES模型
1.5 混匀时间计算模型
1.6 计算参数与边界条件
2 物理模型
2.1 模型参数
2.2 实验条件
3 结果与讨论
3.1 模型验证
3.2 湍流模型对比
3.3 混匀时间计算模型
3.4 气泡膨胀模型
4 结论
来源--金属学报