湍流分离流中颗粒的扩散机制

邱骁,丁珏,王忠杰,翁培奋

(上海大学上海市应用数学和力学研究所,上海　200072)

湍流分离流中颗粒的扩散机制

邱骁,丁珏,王忠杰,翁培奋

(上海大学上海市应用数学和力学研究所,上海200072)

后台阶流动包含分离流重要的流动特性,采取欧拉-拉格朗日耦合算法对后台阶分离流动中颗粒扩散运动进行数值研究.气相场采取大涡模拟方法,亚格子模式基于标准的Smagorinsky模式,颗粒相运动采取轨道法模拟.计算所得气相的流向平均速度和平均脉动速度与实验结果吻合较好,验证了模型和方法的正确性.基于此,数值分析后台阶两相流动的特性以及流场涡结构的发展和演化过程.结果表明：两相流中颗粒的扩散特性既受到颗粒粒径的影响,又与颗粒和涡结构的相互作用时间有关.后台阶流场中增加结构物时,流场涡结构发生变化,即与扰动源保持一定距离后,涡数量增多,流场中颗粒分布不均匀,较多颗粒聚集在涡的外缘.

后台阶流动；分离流；欧拉-拉格朗日方法；气粒两相流动；大涡模拟

在自然界和工程应用中,气粒两相流是一种普遍存在的湍流流动现象,如锅炉中煤粉颗粒的燃烧、发动机内燃料的扩散和运动、循环流化床内的流动、大气环境气溶胶的扩散等.尽管气粒两相流动广泛存在,但是人们对颗粒与气相的相互作用机制的认识较为缺乏.目前,仅有少数几种气粒两相流动的相互作用机制能够通过分析得到.现实中的气粒两相流往往是非均匀的和各向异性的,经常包含颗粒与壁面间的相互作用以及流动分离的作用.颗粒在湍流流场中的扩散机制十分复杂,因此研究两相湍流输运特性,讨论颗粒在湍流分离流动中的扩散机制,对于揭示颗粒运动规律、强化混合和提高燃烧效率以及大气环境的治理具有重要的学术价值和工程应用价值[1-4].

在分离流研究中,最常用的简化模型是后台阶流动.后台阶流动是一种既包含了分离流各种重要的流动特性,又具有简单几何形状的典型分离流动.Yu等[5],Mohammed等[6]和Chen等[7]用大涡模拟方法数值模拟了气粒两相的后台阶流动,分析了初始速度差对颗粒分布的影响.王兵等[8]数值研究了颗粒在大涡结构中的弥散效应,讨论了不同St数对颗粒分布的影响.Fessler等[9-10]实验研究了气体驱动下不同颗粒在后台阶中的流动和分布,其结果为理论和数值研究提供了参考.

为了进一步研究复杂分离流的特性,本工作基于数值模拟方法深入讨论了两相湍流中颗粒与涡的相互作用,比较结构扰动对颗粒扩散的影响,为分离流动两相湍流的调节与控制提供理论依据.

1　物理模型和数值方法

本工作基于欧拉-拉格朗日耦合算法对后台阶分离流动中的颗粒扩散进行数值研究.气相流场采取大涡模拟的方法,其基本思想是通过精确求解某个尺度以上所有湍流尺度的运动,从而能够捕捉到雷诺平均N-S(Reynolds-averaged Navier-Stokes,RANS)方程无法求解的许多非稳态、非平衡过程中出现的大尺度效应和拟序结构,同时又克服了直接数值模拟需要求解所有湍流尺度而带来的巨大计算开销的问题,因而被认为是最具有潜力的湍流数值模拟发展方向[11-12].

1.1气相控制方程

1.1.1过滤方程

大涡模拟求解的控制方程是经过滤波处理的N-S方程组.根据大涡模拟的思想,将流体中的变量通过空间过滤的方式分为可解尺度量和不可解尺度量.例如,在不可压缩流场中,将流场中变量f化为两部分,即

式中,D表示流体区域,G表示过滤函数.过滤函数G满足

式中,∆为过滤尺度.

过滤函数有很多种,但是对于OpenFOAM这类采取有限体积的计算方法,采用盒式过滤方法更合适.盒式滤波函数为

式中,xi为任意网格节点的坐标，∆xi为网格尺度.

1.1.2亚格子模型

大涡模拟的基本思想是对可解尺度湍流进行直接数值求解,将不可解尺度湍流对可解尺度湍流产生的影响通过亚格子模式计算.常用的亚格子模型有亚格子的涡黏和涡扩散模型及其改进模型、结构性亚格子模式、理性亚格子模式等.

本工作采取的大涡模拟(large eddy simulation,LES)亚格子模型为标准Smagorinsky模型.Smagorinsky模型是1963年由Smagorinsky提出的第一个亚格子模型,其中不可压缩湍流的亚格子应力的表达式为

式中,νT为亚格子涡黏系数,ij为可解尺度的变形率张量,δij为克罗内克尔函数.

1.1.3大涡模拟的气相控制方程

气相流动的大涡模拟中,对大于滤波宽度的大尺度涡运动进行直接数值模拟,而对小于滤波宽度的小尺度涡运动则引入亚格子模型进行模拟.亚格子模型采用Smagorinsky模型.考虑气相和颗粒相的双向耦合,过滤后得到直角坐标系下的不可压缩Navier-Stokes方程.

连续性方程：

亚格子应力τij采用涡黏性假设确定,有

式中,µT为亚格子涡黏性系数.本工作中采用Smagorinsky亚格子模式,即

1.2颗粒相控制方程

在气相流动达到稳定后,采用拉格朗日颗粒轨道跟踪法对颗粒的运动进行模拟.根据颗粒所受到的作用力,通过牛顿第二定律积分得到颗粒的速度,再通过对颗粒的运动方程积分,得到颗粒的位移.假设颗粒为球体,颗粒相的密度远大于气相的密度,颗粒间的相互作用忽略不计.考虑颗粒受到的流场的曳力,颗粒相在直角坐标系中的控制方程为

式中,xp为颗粒的位置,Up为颗粒速度,Uf为对应于颗粒当前位置的流场速度,τp为颗粒弛豫时间,Cd为流场曳力系数,ρf为流场密度,ρp为颗粒密度,dp为颗粒的直径.

1.3数值方法

本工作中计算区域与Fessler等[10]的实验基本相同,其中后台阶的扩张比为5∶3,后台阶高度为26.7 mm,气体入口速度为10.5 m/s.基于入口速度及入口槽道半宽的雷诺数为13 800,基于入口速度和后台阶的雷诺数为18 400.为了确保流动的充分发展,台阶后流向的尺度为34H,为实验尺寸的两倍.计算区域的几何尺寸如图1所示,台阶入口前段采用20×30×30网格进行划分,后台阶入口到出口段,流向采用拉伸网格,整个区域采用600×50×30网格进行划分.

图1　计算区域几何尺寸Fig.1　Geometrical dimension of computational domain

1.3.1边界条件

气相入口速度由文献[10]提供的实验数据确定,并加入随机扰动.气相的壁面速度采用无滑移边界条件,出口采用对流出口边界条件.入口及壁面处压力均采用Neumann边界条件,出口压力设定为101.325 kPa,前后均采取周期性边界条件.颗粒相的加入与文献[10]中的情形相同,在流场稳定之后喷入颗粒.求解的计算时间步长为10-4s.采用PISO(pressure implicit with splitting of operator)算法对速度压力进行迭代计算.

1.3.2数值方法的验证

图2为计算得到的流场的流线图.计算得到的回流区长度为7.6H,文献[10]中回流区长度为7.4H,对比可知,计算结果与实验数据吻合较好.

图2　流线图Fig.2　Streamline chart

图3和4给出了加入直径为70µm,密度为2 500 kg/m3,颗粒质量携带率为4%的固体颗粒后,在x/H=2,5,7,9,14时截面流场流向平均速度、流向脉动速度的计算值与实验数据的比较.结果显示：流场流向速度的数值计算结果与文献[10]中的实验数据基本吻合.这验证了本工作中采取的模型和数值方法是正确的.

图3　x/H=2,5,7,9,14时的流向平均速度Fig.3　Mean streamwise velocities when x/H=2,5,7,9,14

图4　x/H=2,5,7,9,14时的脉动速度Fig.4　Fluctuating velocities when x/H=2,5,7,9,14

2　后台阶两相湍流的数值计算结果与讨论

2.1两相湍流涡量的演化

通过计算可知,在后台阶流动中由于边界的变化,流体在台阶后方出现流动分离.涡卷起、运动、合并,在台阶下游形成不断变化的回流区,并于下壁面7.6H位置处再附.随后,流体继续向前发展,与上壁面产生的涡相互作用,合并和破碎之后流出计算区域.

图5为加入直径为70µm,密度为2 500 kg/m3,颗粒质量携带率为4%的颗粒前后,流场中分离流动涡结构的演化情况.图5(a)为未加入颗粒时流场发展稳定湍流涡量的瞬时图,可以看到,后台阶回流区涡量变化剧烈,上下壁面附近涡的数量基本相同.图5(c)和(e)分别为随时间发展流场涡量和颗粒的瞬时分布.当颗粒刚进入流场(t=1.1 s)时,颗粒对流场湍流的影响不明显；随着颗粒的扩散,颗粒对流场的作用越来越明显.当t=1.75 s时,上壁面附近聚集较多的颗粒,流场大涡的数量明显减少,涡破碎强烈；靠近下壁面附近,颗粒密度分布较小,涡破碎的剧烈程度相比上壁面附近要低.

图5　涡量场的瞬时演变和颗粒密度的法向分布Fig.5　Instantaneous evolution of vorticity field and particle wall-normal density distributions

2.2两相湍流颗粒扩散的机制

2.2.1不同St数(斯特罗哈数)颗粒的扩散特性

St数为颗粒弛豫时间τp与流动特征时间τf的比值,它反应了流体中颗粒的运动行为.当气相场动态稳定(t=1s)后,分别加入粒径dp为2,20,50,100,150和200µm的颗粒(入口速度与气相场速度相同),其中颗粒的密度为2 500 kg/m3,颗粒质量携带率为4%.图6为计算所得的t=1.8 s时气相的瞬时涡量图.图7为t=1.8 s时不同粒径颗粒在流场中的瞬时分布.

计算结果显示：颗粒St数不同,颗粒在流场中的分布有明显不同.St＜1时颗粒(dp= 20µm)的弛豫时间τp远小于流动特征时间τf,颗粒的扩散运动对湍流流场充分响应,部分颗粒被流场大涡带入回流区域.对于St～o(1)(dp=50µm)的颗粒,流场涡结构对颗粒分布的影响减弱,在流场大涡外缘及涡区域附近颗粒的浓度很高,回流区仍有部分颗粒存在.随着粒径的增加,颗粒受到流场涡结构的影响越来越小.对于St≫1的情况,大部分颗粒沿着初始流动方向运动,只有较少的粒子从主流中扩散出来.由此可知,流场中颗粒相的扩散运动与颗粒的St数密切相关,即受到大涡结构和颗粒惯性的共同作用.

图6　气相的瞬时涡量Fig.6　Instantaneous vorticity of gas-phase

图7　不同粒径颗粒在流场中的瞬时分布Fig.7　Instantaneous distributions of different particles in the flow field

2.2.2两相初始速度差对颗粒扩散的影响

α=Up/Uf表示气相与颗粒相的速度差,其中Up为颗粒相初始速度,Uf为气相初始速度.流场中加入颗粒的密度为2 500 kg/m3,颗粒质量携带率为4%.

计算结果如图8所示.可以看出：在不同的两相初始速度差条件下,颗粒的运动规律和分布情况有所不同；对于St≪1,St～o(1)的颗粒,当α=1.0,0.5,0时,颗粒的分布情况相似,这说明两相速度差对St数较小的颗粒影响不明显.

图8　x/H=25时颗粒的流向平均速度Fig.8　Mean streamwise velocity of particles when x/H=25

图9为t=1.8 s时不同初始滑移条件下颗粒的分布情况.对于St≫1(dp=200µm)的情况,当α=1时,颗粒的运动轨迹基本上不受流场大涡结构的影响,颗粒沿着初始运动方向穿透大涡；随着α的不断减小,流场对颗粒运动的作用越来越明显；在α=0时,颗粒集中在流场大涡外缘及相邻大涡附近的区域.这说明两相速度差对St数较大的颗粒影响显著.分析原因主要是较大粒径的颗粒初始速度较小,颗粒经过大涡的时间较长,其扩散和分布受到流场涡结构的影响亦较大.因此,流场涡对颗粒运动和颗粒分布的影响与颗粒经过大涡的时间有关.

图9　t=1.8 s时不同初始滑移条件下颗粒的分布情况Fig.9　Particle distributions under different initial slip conditions when t=1.8 s

2.2.3增加结构物对颗粒扩散运动的影响

基于前面的讨论,进一步分析后台阶上壁面处加入结构物对两相湍流颗粒扩散运动特性的影响,其中结构物的高度h为30 mm,宽度w为10 mm,与入口的距离d为10 mm.

当流场稳定(t=1 s)后,加入粒径为100µm的颗粒,颗粒密度为2 500 kg/m3,颗粒质量携带率为4%.气相和颗粒相的初始速度为10.5 m/s.当t=1.8 s时,流场涡量和颗粒分布如图10所示.

图10　上壁面有无结构物时流场中的涡量分布和颗粒分布Fig.10　Distributions of vorticity and particles in the flow field with and without upper surface structure

由图10可以看出,上壁面加入结构物后,后台阶的回流区长度较未加入结构物时略小,且涡的数量明显增加,流体扰动更加剧烈.

图11　t=1.8 s时流场涡量最大值Fig.11　Maximum values of vorticity in the flow field when t=1.8 s

图11给出t=1.8 s时两个不同截面流场沿垂向涡量最大值的变化情况.可以看出：当上壁面没有结构物时,在x/H=14截面处流场涡量的最大值(1 150 s-1)均出现在上下壁面附近,说明此区域存在大尺度的涡,而其他区域涡量值较小；当上壁面有结构物时,x/H=14截面处流场的涡量值均增大,最大值为700 s-1,且主要位于主流区域,说明流场结构发生了变化,流场涡数量增多.在x/H=25截面处,由于远离扰动源,两种流场涡量的变化趋势相似,即流场中无结构物时,流场的涡量最大值为3 948 s-1；加入结构物后,流场的涡量最大值变为3 576 s-1.由此可知,加入结构物后,流场中小尺度涡的数量明显增加,大尺度涡的数量减少,流场中颗粒分布不均匀.

3　结论

本工作通过建立气粒两相流模型,对两相湍流后台阶流动进行数值研究,并分析了不同St数颗粒、不同入口速度下的颗粒运动规律和分布情况.

(1)在气粒两相湍流中,颗粒对流场涡的演化有影响,流场大涡不断破碎成小尺度的涡；在局部区域,颗粒聚集,浓度较高,颗粒对流场的作用明显.

(2)在气粒两相湍流中,颗粒的扩散受到流场涡结构与颗粒相互作用的影响.当颗粒St≪1时,颗粒相的运动完全被气相场的运动控制；当颗粒的St～o(1)和St～o(10)时,在大涡结构的作用下,颗粒主要集中在流场大涡外缘及相邻大涡附近的区域；当颗粒的St≫1时,流场大涡结构对颗粒分布的影响减弱,颗粒在流场中的分布主要由初始速度决定.此外,流场大涡对颗粒运动和颗粒分布的影响与颗粒经过大涡的时间有关.颗粒经过大涡的时间越长,其运动和分布受到大涡的影响也越大.

(3)加入结构物后,流动分离现象更加显著.中小尺度涡的数量增加,大尺度涡的数量减少.颗粒聚集在涡的外缘,流场中颗粒分布不均匀.

[1]EATON J K.Two-way coupled turbulence simulations of gas-particle flows using point-particle tracking[J].International Journal of Multiphase Flow,2009,35(9)：792-800.

[2]HADINOTO K,JONES E N,YURTERI C,et al.Reynolds number dependence of gas-phase turbulence in gas-particle flows[J].International Journal of Multiphase Flow,2005,31(4)：416-434.

[3]YAN F,LIGHTSTONE M F,WOOD P E.Numerical study on turbulence modulation in gas-particle flows[J].Heat&Mass Transfer,2007,43(3)：243-253.

[4]LAITONE J A.A numerical solution for gas-particle flows at high Reynolds numbers[J].Journal of Applied Mechanics,1981,48(3)：465-471.

[5]YU K F,LAU K S,CHAN C K.Large eddy simulation of particle-laden turbulent flow over a backward-facing step[J].Communications in Nonlinear Science and Numerical Simulation,2004, 9(2)：251-262.

[6]MOHAMMED H A,ALAwI O A,WAHID M A.Mixed convective nano-fluid flow in a channel having backward-facing step with a baffle[J].Powder Technology,2015,275：329-343.

[7]CHEN Y T,NIE J H,ARMALy B F,et al.Turbulent separated convection flow adjacent to backward-facing step-effects of step height[J].International Journal of Heat and Mass Transfer, 2006,49(19)：3670-3680.

[8]王兵,张会强,王希麟.颗粒在大涡结构中的弥散[J].力学学报,2005,37(1)：105-109.

[9]KULIcK J D,FESSLER J R,EATON J K.Particle response and turbulence modification in fully developed channel flow[J].Journal of Fluid Mechanics,1994,277(1)：109-134.

[10]FESSLER J R,EATON J K.Turbulence modification by particles in a backward-facing step flow[J]. Journal of Fluid Mechanics,1999,394：97-117.

[11]SAGAUT P.Large eddy simulation for incompressible flows：an introduction[M].3rd ed.Berlin：Springer-Verlag,2006.

[12]张兆顺,崔桂香,许春晓.湍流理论与模拟[M].北京：清华大学出版社,2005：54-129.

本文彩色版可登陆本刊网站查询：http://www.journal.shu.edu.cn

Particles diffusion mechanism in turbulent separation flow

QIU Xiao,DING Jue,WANG Zhongjie,WENG Peifen
(Shanghai Institute of Applied Mathematics and Mechanics,Shanghai University, Shanghai 200072,China)

Back-facing step flow contains important characteristics of a separation flow. Based on this,a numerical simulation on particles diffusion mechanism in turbulent separation flow was conducted by Euler-Lagrangian coupled method.The gas phase governing equations were solved by large eddy simulation method and the subscale grid model was solved by Smagorinsky model.The particle phase was tracked with Lagrangian method.The accuracy of the solver was tested by comparing the simulation results with the experiment data.Based on this solver,the characteristics of the two phase flow,and the developing and evolutional process of vortex were numerically analyzed.The results showed that the diffusion of particles was related to particle’s diameter,the interaction time between particles and vortex structures.The structure of vortex changed when an obstacle was added to the flow field.The number of vortex increased at a certain distance away from the obstacle.The particles concentrate at the edge of the vortex and the distribution of particles were nonuniform in the flow field.

backward-facing step flow；separation flow；Euler-Lagrangian method；gasparticle flow；large-eddy simulation

O 359

A

1007-2861(2016)05-0586-11

10.3969/j.issn.1007-2861.2015.03.002

2015-06-12

国家自然科学基金资助项目(11472167)

丁珏(1973—),女,副研究员,博士,研究方向为多相流体力学等.E-mail:dingjue_lu@shu.edu.cn