首页/文章/ 详情

JFM: 多分散气体-固体颗粒多相流的统一气体动力学波-颗粒方法

6月前浏览4659








本文摘要:(由ai生成)本文提出了一种多尺度统一气体动力学波-颗粒(UGKWP)方法,用于模拟多分散气体颗粒多相流。该算法采用气体动力学格式(GKS)模拟气相,并针对多分散固体颗粒相提出了UGKWP方法。UGKWP方法可以根据局部Knudsen数在欧拉流体方法和拉格朗日颗粒方法之间自动切换,优化了模拟效率。通过模拟气固流化系统,验证了该方法的有效性,流场的时间平均变量与实验吻合良好。此外,还研究了冲击颗粒-床的相互作用,证明了该算法在高度可压缩情况下多分散气体-颗粒系统中的应用潜力。最终,文章指出UGKWP方法在模拟具有多尺度特性的多分散气体-颗粒多相流方面的优势。  


 

本文提出了一种统一的气体颗粒多相流算法。采用气体动力学格式(GKS)模拟气相,并针对多分散固体颗粒相提出了多尺度统一气体动力学波-颗粒(UGKWP)方法。对于每个分散的固体颗粒相,UGKWP中确定性波和统计颗粒的分解是基于局部单元的Knudsen数。固体颗粒相的方法可以在小单元的克努森数下变成欧拉流体方法,在大单元的克努森数下成为拉格朗日颗粒方法。由于单个颗粒相的不同物理性质,如颗粒直径、材料密度等,这成为模拟克努森数变化较大的分散颗粒相的优化算法。气体颗粒流的GKS-UGKWP方法统一了欧拉-欧拉和欧拉-拉格朗日方法。UGKWP中固体颗粒相的粒子和波分解及其耦合演化是为了平衡物理精度和数值效率。模拟了气固流化系统的两种情况,即一个循环流化床和一个湍流流化床。捕获了流化颗粒的典型流动结构,流场的时间平均变量与实验测量结果吻合良好。此外,通过所提出的方法研究了冲击颗粒-床的相互作用,验证了在高度可压缩情况下多分散气体-颗粒系统的算法,其中研究了颗粒云的动态演化过程。


SIMPOP
分散固体颗粒相的统一气体动力学波    
 

 

 

1.分散相的控制方程

分散相的演化受动力学方程控制

其中

 

并且注意,跨组分碰撞的质量守恒可以通过上述gik自动满足。阻力的一般形式可以通过阻力模型来评估

在数值模拟中,(2.4)中的τ将由为固相选择的阻力模型闭合,稍后将详细介绍。此外,考虑的另一个相互作用力是浮力,浮力可以建模为


2.统一气体动能波-粒子法
在本小节中,介绍了分散相演化的UGKWP。通常,分裂算子用于在固相ts的数值时间步长内通过以下程序求解(2.1):


为简洁起见,由Ld1、Ld2和Ld3更新的变量表示为


首先,我们关注Ld1。没有外力和固体颗粒跨物种碰撞的分散相动力学方程为

为简洁起见,本小节将忽略下标k。动力学方程的积分解可以写成

在UGKWP中,宏观保守变量和微观分布函数都将在有限体积框架下更新。单元i的单元平均宏观变量Wi由守恒定律定义

其中Wi定义为的单元平均宏观变量

Fij可以通过以下公式计算

Si由下式定义:


代表由于固体粒子的非弹性碰撞而损失的能量

 


其中e∈[0,1]是用于确定非弹性碰撞中损失能量百分比的恢复系数,e的值为0.8。

将随时间变化的分布函数(2.11)代入(2.14),通量可以重写为

在细胞界面获得局部平衡态g0的过程和g(t)的构造与GKS中的过程相同。对于二阶精度,单元界面周围的平衡状态g写成

空间导数的系数可以从宏观变量的相应导数中获得

时间导数的系数可以由相容条件确定

现在,平衡态g(x,t,u)中的所有系数都已确定,其积分变为

系数:

因此,平衡态的通量由下式给出

此外,粒子自由传输的通量贡献是通过跟踪从f0采样的粒子来计算的。因此,单元平均宏观变量的更新可以写成

 

粒子分布的演变可以写成

第一次与其他粒子碰撞之前的自由传输时间表示为tc,然后tc的累积分布函数为

每个粒子的自由流动时间tf分别由下式确定

其中t是时间步长。因此,在一个时间步长内,所有粒子可以分为两组:无碰撞粒子和碰撞粒子,它们由时间步长t和自由流动时间tf之间的关系决定。具体来说,如果tf=t,则该粒子是无碰撞的,并且在整个时间步长中完全跟踪其轨迹。相反,如果tf<t,则该粒子是碰撞粒子,并且其轨迹被跟踪到tf。在模拟中,碰撞粒子将在tf处被消除,通过计算其通过细胞界面通量的贡献,该粒子携带的相关质量、动量和能量将合并到相关细胞中的宏观量中。更具体地说,在时间间隔t∈[0,tf]内的自由流过程中的粒子轨迹由

确定,此外有:

现在,已经确定了(2.25)中的所有项,并且可以定义宏观变量Wi。  

所有粒子都已追踪到时间tf。tf=t的无碰撞粒子将在时间步长结束时存活;而tf<t的碰撞粒子在其第一次碰撞后将被删除,并且假设其在该单元中进入平衡状态。因此,在每个时间步长结束时,单元i中碰撞粒子的流体动力学宏观变量可以直接通过

Wp i是细胞中剩余的无碰撞粒子的质量、动量和能量。这里,宏观变量Whi说明了所有消除的碰撞粒子,这些碰撞粒子可以在下一个时间步长开始时基于麦克斯韦分布从Whi重新采样。现在,已经给出了宏观变量和微观粒子的定义。上述方法就是所谓的统一气体动力学粒子(UGKP)方法。  

上述UGKP可以进一步发展,以在效率和内存减少方面获得优化的UGKWP方法。在UGKP方法中,在每个时间步长中将所有粒子分为无碰撞粒子和碰撞粒子。碰撞粒子在第一次碰撞后被删除,并在下一个时间步长开始时从Whi重新采样。然而,只有重新采样的粒子中没有碰撞的部分才能在下一个时间步长中存活下来,并且所有重新采样的碰撞粒子都将再次被删除。幸运的是,这些碰撞粒子的输运通量可以在不使用粒子的情况下进行分析评估。因此,我们根本不需要从Whi重新采样这些碰撞粒子。根据累积分布(2.27),无碰撞粒子的比例为exp(-t/τ),因此,在UGKWP中,只有来自单元i中流体动力学变量Whi的无碰撞粒子才会被重新采样,其总质量、动量和能量

现在,与UGKP相同,在UGKWP中,粒子自由流动的净通量,包括前一时间步长的剩余粒子和重新采样的无碰撞粒子,可以通过以下公式计算

因此,UGKWP中的宏观流量变量定义如下

来自未采样碰撞粒子的通量函数(Zhu等,2019b;Liu等2020)定义如下:

其中

 

以第k个分散相为例,其与其他分散相的碰撞可以通过

取欧拉域中的矩ψ=u,fk=gk+O(τk),我们可以得到

在本文中,假设第i个和第k个分散相之间的辅助速度Uik为

现在,我们需要确定τik,即第i个和第k个分散相之间的碰撞时间。通常,多分散特定流中常用的参数是βik,称为所谓的固体间阻力模型,与τik有以下关系:

将(2.40)代入(2.41),可以明确地得到τik的表达式

另外,可以通过解析解获得

本文将使用Mathiesen基于颗粒流动力学理论(KTGF)提出的固体间阻力模型(Mathiesen等1999)

其中pc,ik是第i和第k分散相之间的碰撞压力

由于非弹性粒子碰撞中的颗粒温度值较低,因此忽略了密集粒子碰撞对颗粒温度的影响(Fox&Vedula,2010)。  

最后,在第三部分Ld3,有:

其中一个固体粒子的加速度可以分解为三部分

其中:

取Ld3方程上在欧拉域中的矩,我们得到

第k分散相与气流之间的气固阻力:

隐式离散化:

不含ap的加速度可以表示为:

第k个固相的宏观变量定义为

与MP-PIC方法的处理一样,ap在最后定义为(Snider 20012007)

在本文中,由(2.60)获得的ap进一步受到以下稳定性条件的约束:

其中Δcell是网格尺寸,kc是小于1的安全系数,如本文中使用的0.8。现在加速度可以完全确定为

第k个固相Un+1k的宏观速度更新为

此外,剩余自由传输粒子的速度和位置定义为

上述过程用于在一个时间步长ts中更新分散颗粒相。


3.固体颗粒相的Kn和流态
参数Knk代表第k个分散颗粒相的克努森数,由碰撞时间τk与宏观流动特征时间尺度的比值定义

特征时间tref取固相的时间步长ts,τk为固体粒子碰撞之间的时间间隔。在本文中,τk定义为(Passalacqua等,2010;Marchisio和Fox 2013)

g0是具有以下形式的径向分布函数:

其中c=t/s,max是总固体体积分数t与多分散固体混合物的允许最大值s,max的比率。第k个分散相的流动状态由Knk决定。通常,对于稀释流,固体颗粒之间的碰撞频率较低,导致Knk较大,UGKWP会对固体颗粒进行采样和跟踪,自动保持非平衡状态。相反,在高浓度区域,颗粒之间的高碰撞频率意味着固相处于平衡状态,UGKWP中不会对任何颗粒进行采样。在e=1的连续流态的极限下,上述(2.1)的UGKWP方法可以恢复以下流体动力学方程的解:

在(2.69)中,pk是第k个分散固相的压力,它是动压,碰撞压力和摩擦压力的总和。对于稠密的特定流(Dou等人,2023)已经开展了诸多研究。在本文中,碰撞压力由下式计算

第k分散相持久的颗粒间接触和摩擦在固相接近堆积时起着重要作用。本文采用Johnson–Jackson模型(Johnson&Jackson 1987; Houim&Oran 2016)

为了避免过度填充问题,在固相的UGKWP方法中采用了所提出的接近填充条件的通量限制模型(Yang等,2022c)。  

气相气体动力学方案  
1.气相控制方程
气相被视为连续流,控制方程是NS方程,其源项反映了相间的相互作用(Gidaspow 1994;Ishii和Hibiki 2006)

 

其中,应变率张量σ为

 

根据Toro(2013),(3.1)可以写成以下形式:

2.气体演化的气体动力学方案
体流动由具有相间相互作用的NS方程控制,其解将由相应的GKS获得,这是UGKWP在连续介质状态下的一种限制方案。通常,气相(3.4)在一个时间步长中的演变可以分为三部分

由Lg1、Lg2和Lg3定义的变量表示为

首先,没有nozzle的动力学方程由下式所示:

其中u为速度,τg为气相的弛豫时间,fg为气相分布函数,gg为相应的平衡态(麦克斯韦分布)。局部平衡状态gg可以写成

碰撞项满足相容条件

f在单元界面的积分解可以写成

(3.12)中的初始气体分布函数可以构造为

其中H(x)是Heaviside函数,fl0和fr0是一个单元界面左侧和右侧的初始气体分布函数。初始气体分布函数构造为

其中gl和gr是单元界面左侧和右侧的麦克斯韦分布函数,其可以完全由宏观保守流动变量Wl和Wr确定。

于Chapman–Enskog展开,分布函数的非平衡部分满足

界面周围的平衡状态g被建模为

g可以由相容性条件确定

在确定了初始气体分布函数f0和平衡状态g中的所有参数后,将(3.13)和(3.18)代入(3.12),单元界面处的时间相关分布函数f(x,t,u,ξ)可以表示为

计算一个时间步长上的通量传输

其中nij是细胞界面的法向量

在第二部分中,Lg2 nozzle项

其中,

在第三部分中,Lg3用于相间相互作用

第二个方程表示气相与多分散相之间的动量交换

1GKS-UGKWP方法求解多分散气-颗粒两相流的流程图


其中βt和Ut是整个固相的等效动量传递系数和速度

βt和Uk的计算是基于固相的n+1态的变量。对于上述方程,可以得到Ug的解析解

最后,气相的能量可以通过下式确定

在演化过程中,ts、k和tg将基于Courant–Friedrichs–Lewy(CFL)条件进行计算;固相将首先更新一个固体时间步长ts=min(ts,k);然后,将基于气体时间步长tg更新气相,直到i tg,i=ts,并且将完成ts中的气-颗粒两相流的演变。图1给出了GKS-UGKWP多分散气体颗粒流的流程图,与单分散对应物的主要区别在于固相的计算,用蓝色方框标记。

 

1循环流化床工况下固体颗粒的特性(Niemi 2012)  

数值模拟  

1.循环流化床

1.1. 案例描述

第一种情况是具有两个分散固相的循环流化床(CFB)(Niemi 2012;Wang等2015)。实验数据将用于验证GKS-UGKWP方法。与之前的研究一样(Wang et al.2015),本文采用了D×H=0.4m×3m的二维域。整个区域使用均匀的矩形网格,网格数为66×500,相应的网格尺寸Δx≈Δy=6×10^-3m。根据实验测量(Niemi 2012),固体颗粒的总存量为2.85kg,每个分散相的质量分数、直径和材料密度列于表1中。在这种情况下,最大固体体积分数取为t,max=0.55。最初,固相均匀分布在整个区域,根据表1所示的质量,初始固体体积分数为1=0.0498和2=0.0140,假设立管厚度为T=1.5 cm,这与Wang等人(2015)中使用的值相同。在模拟中,固体粒子可以自由地离开顶部边界处的域。为了简化模拟,左边界和右边界是固定壁,因此逃逸的固体颗粒将从底部边界在计算域中得到补充,而不是像Wang等所采用的那样打开右边界。(2015)。速度Ug=2.25 m s^-1的气体通过底部边界流入区域,使固体颗粒流化。在这种情况下,Gibilaro提出的广泛使用的阻力相关性用于两个分散相(Gibilaro等,1985;Mckeen&Pugsley,2003),可以写成

 

1.2.结果  

在这种情况下,模拟时间为10.0秒,并使用6.0到10.0秒的结果进行平均。物理上,为了研究立管中不同垂直位置的流动特性,实验中在h=0.32、0.40、0.80、1.20m处设置了四个压力计。将上述四个高度的总固体体积分数t和固相的总垂直速度Us取平均值,并与图2中的实验测量值进行比较。注意,整个固相的总垂直速度Us是通过固体体积分数加权的单个速度获得的。网格进一步细化到80×600个网格点,t和Us的时间平均结果如图3所示。与图2中的结果相比,精细网格也获得了类似的结果,因此以下分析仅基于66×500点的较粗网格。

 

2:采用GKS-UGKWP方法,共66×500个网格点,得到了不同立管高度下固相Us的时间平均总固体体积分数t和时间平均总垂直速度分布,并与实验结果进行了比较  


图2显示,GKS-UGKWP方法的数值结果与实验测量结果基本一致,验证了所提出方法的可靠性和准确性。在h=1.20 m时,计算结果t=2%–3%在某种程度上低于实验值6%,这可能是由于边界处理,使得从顶部边界逸出的颗粒在底部边界得到补充,而不是从侧壁得到补充。阻力模型的选择也会影响解决方案。如引言中所述,考虑细观亚网格流结构的阻力模型比传统的均匀阻力模型表现更好(王和李,2007;朱等人2021)。如果采用中尺度阻力模型,目前的方法可以在h=1.20m左右给出更好的预测,这将在未来进行研究。在实验中,通常假设固体浓度与垂直立管上的压力梯度相同。这里,由时间平均压力梯度评估的预测t也如图4所示。与图2(a)相比,这两种方法的固体浓度值在底部密集区域(即h=0.23、0.40 m)显示出一些偏差,但在其他区域(即h=0.80、1.20 m)一致。


3采用GKS-UGKWP方法,在总共80×600个网格点上,获得了不同立管高度下固相Us的时间平均总固体体积分数t和时间平均总垂直速度的分布,并与实验测量结果进行了比较  

 

4通过GKS-UGKWP方法和与实验测量值的比较,根据不同立管高度下的时间平均压力梯度  


此外,固体颗粒在不同时间的快照如图5所示。通常,固体颗粒倾向于在提升管底部和壁附近积聚,导致这些区域的浓度相对较高。此外,瞬时结果清楚地显示了瞬时共存和动态干预的稀/密流区域。空间演化的固体体积分数很难被混合EE/EL方法平稳地捕获。在对单分散CFB情况的研究中也发现了上述特征。  

对于每个分散固相,Knk,由Knk=τk/ts与第k个分散相的局部碰撞时间τk定义,如图6所示。图6显示了波分量、波轮廓和粒子分量,即通过其垂直速度着色的采样粒子集。每个固相的垂直速度也在图6中给出。对于这两种固相,由于颗粒在提升管的近底部和近壁区域的积聚和碰撞,Kn通常在这些区域较小。两个分散的固体颗粒相根据各自的Kn调整其重量以适应UGKWP中的波和颗粒分量。GKS-UGKWP用于多分散流的一个明显优点是,每个分散相都可以进行独立的波和粒子分解。


5总固体体积分数在(a)t=6.0 s、(b)t=7.0 s、(c)t=8.0 s、(d)t=9.0 s和(e)t=10.0 s时的瞬时快照


6Kn、固体体积分数、UGKWP波中逐波固体体积分数的瞬时快照、UGKWP波中采样粒子的集  合以及t=8.0 s时固相的垂直速度Us:(a)Kn1、(b)1、(c)波1、(d)P1、(e)Us、1、(f)Kn2、(g)2、(h)波2、(i)P2、(j)Us,2。下标1的图是第一个固相(小颗粒)的结果,下标2的图是第二个固相(大粒子)的结果。Knk的价值是由Kn的传说所决定的。固体体积分数k和相应的波分量波k由eps图例着色。粒子集Pk中的离散粒子和固相的垂直速度Us,k由Us图例着色,其中k=1,2。Kn的传说是指数分布的。请注意,在t=8.0 s时,1和2的总和正好等于t,如图5所示。  


2.湍流流化床
2.1.案例描述  
致密湍流流化床(TFB)被应用于石油精炼行业,并进行了实验和数值研究(Gao等,2009;Li等,2009)。在这个问题中,涉及两种颗粒,如流化催化裂化(FCC)催化剂(细颗粒)和millet(粗颗粒),其详细性质如表2所示。表2显示,两种类型的颗粒的密度非常接近,而颗粒尺寸非常不同。这种流动条件给单独具有单个固相的气体-颗粒系统的数值方法带来了挑战,其中GKS-UGKWP中多个固相的跟踪似乎适合于这个问题。本文研究了初始床高H0=1.155m,气体流速Ug=0.53ms-1的情况。计算域为D×H=0.5 m×4 m,由40×300的均匀矩形网格覆盖。在本研究中,最大固体体积分数取s,max=0.65。在模拟开始时,所有固体颗粒都均匀分布在整个提升管中,小颗粒相和大颗粒相的初始体积分数分别为1=0.118和2=0.070。从顶部边界逸出的固体粒子将通过底部边界再循环回计算域。对于气相,在顶部边界采用标准大气条件,气体通过底部边界吹入立管,固体颗粒的密度为1463.3 kg m^-3,按其初始体积分数加权。与上述CFB情况相同,颗粒相和气相的侧壁分别采用混合边界条件和无滑移壁边界条件。


2TFB中细颗粒和粗颗粒的性质(Gao等,2009)


在GKS-UGKWP中,每个固相都可以为多分散系统选择最准确、最合适的阻力模型。如表2所示,细颗粒和粗颗粒分别为Geldart A和D基团,并采用不同的阻力模型来评估FCC催化剂和millet颗粒相中的气固相互作用。对不同的阻力模型及其修正进行了研究和比较(Gao等人,2009)。更具体地说,对于粗颗粒,使用Gidaspow模型(Gidaspow1994)

而对于精细的FCC催化剂颗粒,则采用四区阻力模型

其中dk是固体颗粒的直径,(4.3)中的d*k是FCC催化剂相的有效直径,取300μm以更好地与实验测量结果一致(Gao等人,2009;李等人,2009)。(4.3)中提出的阻力模型充分考虑了颗粒团簇的影响,本文也采用了这种模型。此外,在(4.2)和(4.3)中,定义如下


2.2. 结论  

10.0到15.0 s的时间平均结果如图7所示。提升管中的固相在底部区域显示出较高的浓度,在顶部区域显示出较低的密度,并且在1.7–2.0 m左右的非常小的区域中发生急剧转变。总体而言,预测的固相表观密度与实验测量结果一致。值得一提的是,在数值模拟中,选择一个包含子网格信息的合理准确的阻力模型对于准确预测非常重要(Wang&Li 2007;Zhu等2021)。此外,图7还显示了两个固相的剖面图,并显示了沿立管高度的相似趋势。

图8给出了10.0–15.0 s范围内的固相密度的瞬时快照,这表明了TFB的典型特征。底部致密/中间过渡/上稀释区域的共存模式被很好地捕捉到。在顶部稀释区域中的所有固体颗粒都是FCC催化剂(Geldart A)类型。图8显示了粒子团现象,这在这些区域的阻力建模中存在困难。在底部致密区,内部有可变固体颗粒的气泡呈现出复杂的模式,并与致密固相缠结。

 

7(a) 用GKS-UGKWP方法得到了整个固相的时间平均表观密度kρk沿立管高度的分布,并与实验结果进行了比较。(b) 每个分散相的时间平均固体体积分数的分布k,k=1,2和它们的总和t  

 

8整个固相表观密度的瞬时快照kρk在(a)t=10.0 s,(b)t=11.0s,(c)t=12.0s,(d)t=13.0 s,(e)t=14.0 s和(f)t=15.0 s  

 

9Knk、固体表观密度kρk、UGKWP波中的固体波表观密度kρk和t=11.0s时UGKWP Pk中的采样粒子集的瞬时快照:(a)Kn1、(b)1ρ1、(c)波1ρ1,(d)P1、(e)Kn2、(f)2ρ2、(g)波2ρ2和(h)P2。下标1和2分别代表第一(FCC催化剂颗粒)固相和第二(millet颗粒)固相。Knk由Kn图例着色,固体表观密度kρk和波分量波kρk由表观密度图例着色,粒子集Pk中的离散粒子由Us图例着色(固体粒子的垂直速度),其中k=1,2。Kn的色谱是指数分布的  


由于两个固相(60和930μm)的直径相差15倍,因此呈现它们各自的分布是很有趣的。图9显示了t=11.0s时各相的表观密度,恰好是图8(b)中11.0s时的表观浓度。明显的特征是,在上稀释区域,只有较小的(FCC催化剂)颗粒不存在millet颗粒,如图7的时间平均曲线所示。此外,对于h≈1.7m的两个固相,存在一个分离区,上面有一个大的Kn,下面有一个小的Kn。与FFC催化剂相相比,millet相表现出强烈的非平衡态,在底部致密区具有较大的Kn。由局部Kn确定的波和粒子的分解如图9所示,通过轮廓表观密度和散射粒子集。对于FCC催化剂相,拉格朗日粒子完全决定了其在上稀释区域的演化,而波分量在底部区域占主导地位,有大量实际粒子及其碰撞。对于millet大颗粒相,即使在底部密集区域,也会对其演化过程中的大量颗粒进行采样和跟踪。  

有必要了解与气体-颗粒系统中的稀/密和非平衡/平衡相关的流动特性。稀流或浓流通常由固体体积分数决定,而非平衡/平衡由固体颗粒相的Kn决定。稀流和稠密流区域可以与非平衡和平衡状态相关联,特别是对于稠密颗粒流,其颗粒直径或材料密度有许多差异。在UGKWP中,Kn用于根据局部非平衡的程度来确定波和粒子的分解。每个分散相在不同时间的Kn快照如图10所示。对于两个分散相,在h<1.7m的底部致密区域,非平衡/平衡状态在空间上是不固定的,并且变得动态地可相互转换。所提出的适用于多分散流的GKS-UGKWP适用于通过波和粒子的最优分解来处理每个单独的固相。

 

10:催化裂化催化剂颗粒相Kn在不同时间的瞬时快照:(a)t=10s,(b)t=12s,(c)t=14s和(d)t=15s    

millet颗粒相在不同时间Kn的瞬态快照:(e)t=10S,(f)t=12S,(g)t=14S和(h)t=15S。Kn的图例呈指数分布  


3. 激波与稠密粒子幕的相互作用

3.1. 案例描述  

冲击波与固体颗粒床的相互作用是一个极具挑战性的问题,对于一种能够捕捉超音速冲击波传播并计算中等/稠密情况下气固相互作用和颗粒碰撞的数值方法来说(Ling等人,2012;Tian等人,2020)。在这里,通过GKS-UGKWP方法在二维空间中研究了平面激波与粒子幕的相互作用,并将模拟解与实验测量结果进行了比较(Ling等人,2012)。如图11(a)所示,气管中马赫数(Ma)=1.66的平面冲击从左向右(x方向)移动,并遇到宽度为L=2 mm的初始静止粒子幕。从冲击对固体粒子床的冲击开始,将发生向左移动的反射冲击和向右移动的传输冲击。同时,在高速气流的驱动下,固体颗粒将跟随传递的冲击前沿向右移动。计算域X×Y为[-0.375m,0.125 m]×[-0.04 m,0.04 m],由具有1000×160个单元的均匀矩形网格覆盖。在实验中,固体颗粒的直径分布为106-125μm。而在计算中,模拟是基于固体颗粒直径为d1=110μm和d2=120μm的两个分散相。所有固体颗粒的材料密度为ρs=2520 kg m^-3,恢复系数e=0.95。最初,固体颗粒均匀分布。气相在该区域的初始状态与实验相同:pg=8.27×10^4Pa,Ug=0ms^-1和Tg=296.40K,气体常数R=287.05 J(kg K)^-1,比热比γ=1.4。在左边界处,Ma=1.66时的预冲击条件由pg=2.52×10^5 Pa,Ug=304.16 m s^-1和Tg=423.79 K给出。在右边界处采用自由边界条件。此外,对于顶部和底部边界,分别对气相和固体颗粒相采用无滑移和滑移边界条件。为了监测气流的压力,在初始颗粒幕的左前方上游68.6mm和下游64.2mm处设置了两个测量位置。在这种情况下,气相的内部自由度由K(~t)=K0+1.0×(~t/~tend)^3建模,以模拟由于与分散固体颗粒的相互作用而增加的气流湍流强度。这里,~t=t/(L/us)是粒子幕初始宽度L和冲击速度us=572.00 m s^-1的归一化时间,~tend是归一化模拟时间,在这种情况下取350.0。

 

11:激波Ma=1.66与稠密粒子幕相互作用的示意图    

(a)没有近壁间隙的整个区域;(b)具有近壁间隙的局部颗粒分布(87%的展向颗粒帘)  


第k个分散相的固体颗粒上的阻力可以写成一般形式

利用(2.4)中D的定义,可以得到:

其中根据Clift(1970)提出的标准阻力相关性,有:

12:数值结果(用“Num”表示)和实验测量(用“Exp”表示):(a)上游和下游计量点的时间相关气体压力;(b)粒子云锋的轨迹  

3.2. 结果  

测量位置的时间相关压力如图12(a)所示,并与实验测量结果进行了比较。一般地,反射冲击和透射冲击都可以很好地被捕获。此外,固体粒子云上游锋和下游锋的轨迹如图12(b)所示,与实验数据吻合良好。~t=193.0时固体颗粒分布的瞬时快照如图13所示。固体颗粒在现有区域分布不均匀:中心区域(约x=0.045m)的浓度高于上游和下游锋附近的区域。此外,可以观察到轻微的粒子团现象。对第一和第二分散相的进一步比较表明,在上游前沿和下游前沿,具有较小颗粒d1=110μm的第一固相比具有d2=120μm的第二相移动得更快(约3 mm)。

在实验中,固体颗粒帷幕是由颗粒从储层自由下落到试验段中产生的(Ling等,2012)。如实验中所指出的,粒子幕在翼展方向(图11中的y方向)上约占87%。因此,粒子幕和墙之间存在间隙,本文也对此进行了研究。根据颗粒幕在管中占据87%的实验,将考虑靠近壁的13%的间隙。使用新的设置,新计算的压力和粒子轨迹如图14所示。有趣的是,在时间~t=50–80时,下游表压位置的压力变化在传递的冲击通过后,气体压力出现早期下降和后期上升,与之前100%颗粒幕占据的计算相比,与实验数据的一致性更好,如图14中的蓝色圆圈所示。

 

13:归一化时间~t=193.0时颗粒相表观密度的瞬时快照:(a)整个固体颗粒相(来自相1和相2的颗粒的总和),(b)d1=110μm的第一颗粒相和(c)d2=120μm的第二颗粒相。  

 

图14:数值结果(用'Num'表示)和实验测量(用'Exp'表示)
(a)上游和下游标定点的随时间变化的气体压力;(b)粒子云锋的轨迹  
结论  

本文发展了一种适用于多分散气体-颗粒系统的多尺度GKS-乌格奎普方法。特别地,用于具有耦合的分散固体颗粒相和气相的系统的依赖于单元分辨率的动态模型被直接用于设计相应的多尺度方法。为了有效地捕捉平衡和非平衡颗粒演化过程,根据每个固体颗粒相各自的局部Kn,在UGKWP中将颗粒分布函数分解为波动和离散颗粒分量。UGKWP将选择一种最佳方式来描述固体粒子的动力学,方法是分别演化欧拉确定性波,并在物理精度和数值效率之间取得平衡的情况下跟踪单个拉格朗日粒子。在确定固体颗粒动力学时,已经考虑了固体颗粒的性质,例如颗粒尺寸、浓度和混凝土材料。UGKWP的一个显著特点是,波-粒子分解可以在固体粒子相的相应平衡和非平衡状态下,自动将方案引向气体-固体粒子系统的EE和EL方法。研究了流化床中稠密多分散流动的两种情况。具体地,在FCC催化剂反应器中,大的millet颗粒和小的FCC催化剂颗粒的直径在颗粒尺寸上相差15倍。数值实验表明,只有细小的FCC催化剂颗粒出现在提升管的顶部区域,而millet颗粒只存在于底部密集区。此外,对于FCC催化剂相,离散颗粒描述在顶部稀释区域中起关键作用,而波动描述在底部密集区中贡献最大,即使存在大量真实颗粒。对于大颗粒相,millet颗粒只存在于底部密集区,需要一个离散的颗粒描述来捕捉局部的非平衡态。上述观察表明了波-粒分解在描述分散固相方面的灵活性和在动态跟踪平衡和非平衡流动演化方面的适应性。同时,用当前的方法研究了Ma = 1.66的激波与多分散固体颗粒相的相互作用。数值解,如气体流在标准点的压力和固体粒子云前沿的轨迹,与实验测量值进行了比较。尽管目前的研究仅涉及具有两个分散相的颗粒系统,但是算法本身可以容易地扩展到具有多个分散固体颗粒相的气体-颗粒系统。总的来说,固体颗粒相的波和颗粒演化的分解在模拟具有多尺度性质的多分散气体-颗粒多相流方面具有优势。



翻译转载自《Journal of Fluid Mechanics》"Unified gas-kinetic wave-partical method for polydisperse gas-solid particle multiphrase flow"



来源:多相流在线
Marc碰撞多相流湍流UGUM理论GID材料多尺度控制
著作权归作者所有,欢迎分享,未经许可,不得转载
首次发布时间:2024-05-26
最近编辑:6月前
积鼎科技
联系我们13162025768
获赞 108粉丝 110文章 302课程 0
点赞
收藏
作者推荐

CFDPro水动力仿真 │ NS方程组的衍生派系— 水动力方程组

本文摘要:(由ai生成)本文回顾了水动力方程的发展,特别是圣维南方程的贡献。水动力方程在水利、海洋和气象等领域广泛应用,如洪水预测、海岸工程等。文章介绍了水动力方程的一维、二维和三维表达形式。中国水利部正推进数字孪生流域建设,水动力模型作为关键模型,对水利现代化和智慧化至关重要。水动力方程的发展史说起NS方程,大家都知道方程名称中的两位科学家先驱——纳维和斯托克斯。可是还有一名科学家应该被大家知晓,那就是法国科学家圣维南。圣维南全名为AdhémarJeanClaudeBarrédeSaint-Venant,一般用Saint-Venant来代表他。他是十九世纪法国著名的机械师、数学家,主要从事力学、弹性、静水学和流体动力学方面的研究,推导出了非恒定明渠流浅水方程。这个方程也被称为圣维南方程,是现代水利工程中的一组重要方程。圣维南方程是水动力方程一维化形态。圣维南不仅独立的推导出了非恒定明渠流浅水方程,他还在纳维之后,重新推导了关于粘性流动的方程,考虑到了内部的粘性应力,并完全放弃了纳维采用的分子方法,并于1843年发表了相关论文。该论文首次正确地识别了粘性系数及其作为流动中速度梯度的乘数的作用。他还进一步将这些乘积识别为由于摩擦作用从而在流体内部产生的粘性应力。随后,在基于不可压条件和静水压假设等条件下,二维水动力方程组被推导出来,也就是被大众熟悉的浅水方程。浅水方程的假设前提为水体的深度方向尺度远远小于水平方向尺度,因此可以将深度方向的速度分布等效为均匀流。连续方程中的密度项也通过积分转化为水深的表达式。浅水方程的出现极大的方便了科学家通过数学物理方法寻找水体运动规律。但是在水库或是深度方向需要精细考虑的场景,二维的浅水方程就不够准确了。因此科学家将深度方向再进行分层细化,每层之间再通过约束方程进行质量和动量的传递,从而产生了三维浅水方程。水动力方程的表达形式先看圣维南方程组(一维水动力)的表达形式:连续性方程:动量方程:其中,Z表示自由水面高度,S表示过流断面面积,g表示重力加速度,J表示水头损失因子,β表示断面流速均匀性系数,qa表示质量源项。该方程组看起来跟我们认知中的一维NS方程组相差甚远,那是因为在对河(渠)道进行一维化概化时,描述水体运动规律的参数维流量Q和过流面积S。而过流面积又可以表达为水深h的函数S=f(h),流量可以表达为流速v和过流面积的函数Q=f(v,S)。图1河(渠)道过流断面示意图图2河(渠)道沿程坡降示意图在经过一些数学上的变换后,圣维南方程组就可以变为向量形式:其中,守恒量项,通量项,源项。浅水方程组(二维水动力)的表达形式沿用向量形式:其中,守恒量项,x方向通量项,y方向通量项,底坡源项,摩擦源项。图3河道二维流动示意图三维浅水方程与二维较为类似,由于篇幅原因,就不在此赘述,其垂向分层思想由下图所示。图4三维水动力垂向分层示意图从一、二维方程的形式我们可以看出,在NS方程中存在密度项不在了,取而代之的是表达几何属性的面积S或者水深h,其实这就是基于不可压缩假设做了深度方向的均匀化与积分转换。水动力方程的应用领域和价值水动力方程在海洋、气象以及水利领域发挥了重要作用。圣维南方程常用于明渠及简单河道的流量与水位关系的规律性计算;二维浅水方程常用于模拟河流中水体的运动状态,如水位、流速等。在洪水的预测和控制中,二维浅水方程能帮助水文部门预测洪水的演进路径、淹没范围及过程历时,为防洪减灾提供定量分析的重要依据。在海岸工程中,利用波浪模型,浅水方程可用于模拟波浪在近岸浅水区的传播和变形,可以分析波浪对海岸线的侵蚀作用,对港口建设、海岸防护以及海滩保护的方案设计提供重要科学支撑。在湖泊和水库的管理中,浅水方程结合水质模型,能够帮助管理机构了解水体循环、混合和水质演变规律,为水资源的调配、水质保护以及水生态的修护提供决策支持。对于环境保护领域,浅水方程可以用于模拟水流对污染物的输运和扩散过程,评估水质污染风险,为水环境治理和生态保护提供科学依据。在水利工程设计与运行项目中,浅水方程可以用于分析水流条件、评估水工程结构物的稳定性和安全性,帮助设计单位优化设计方案,确保工程实施的安全可靠。目前,中国水利部正在大力推进基于水利专业模型的数字孪生流域建设,而水动力模型正是重要的水利专业模型,所以说数值计算技术不仅在工业界是重要的设计和校验手段,在水利行业同样也是备受关注的水利现代化和智慧化的科学方法。来源:多相流在线

未登录
还没有评论
课程
培训
服务
行家
VIP会员 学习 福利任务 兑换礼品
下载APP
联系我们
帮助与反馈