顆粒形狀對(duì)彎管沖蝕磨損影響的數(shù)值模擬研究①

李震東 許 磊 包士毅

(浙江工業(yè)大學(xué)化工機(jī)械設(shè)計(jì)研究所 杭州 310023)

(過(guò)程裝備及其再制造教育部工程研究中心 杭州 310023)

(浙江工業(yè)大學(xué)機(jī)械工程學(xué)院 杭州 310032)

0 引言

管道是天然氣集輸管中重要的輸運(yùn)設(shè)備,其中的流體會(huì)攜帶固體顆粒撞擊管道壁面,進(jìn)而造成壁面質(zhì)量損失,產(chǎn)生沖蝕磨損問(wèn)題[1]。管道沖蝕磨損過(guò)程較為復(fù)雜,與管內(nèi)介質(zhì)、管道材料、管道所處環(huán)境等因素息息相關(guān)[2]。有研究指出,管道中彎頭處的沖蝕速率是直管段的50 倍左右[3]。因此,較為準(zhǔn)確地預(yù)測(cè)不同條件下管道彎頭處的磨損速率,對(duì)于預(yù)防輸氣管泄漏、減少經(jīng)濟(jì)損失等問(wèn)題有較大意義。

隨著計(jì)算流體力學(xué)的發(fā)展,由Tsuji 等人[4]提出了計(jì)算流體動(dòng)力學(xué)-離散單元法(computational fluid method and discrete element method,CFD-DEM)。該方法被應(yīng)用于二維流化床研究,且在隨后的多相流研究中得到廣泛應(yīng)用[5]。為了模擬更為真實(shí)的工況,一些顆粒形狀模型被開(kāi)發(fā)出來(lái),例如:超橢球模型、組合顆粒模型以及多面體模型[6]被用來(lái)對(duì)非球形顆粒進(jìn)行建模。對(duì)于磨損理論的研究早在二十世紀(jì)五六十年代就開(kāi)始了。1960 年Finnie[7]提出了脆性材料沖蝕模型,并以一彈性球體垂直撞擊脆性材料為例分析脆性材料的沖蝕機(jī)理,隨后多種沖蝕模型被學(xué)者們提出,這些模型對(duì)后來(lái)學(xué)者們運(yùn)用數(shù)值模擬方法來(lái)研究磨損問(wèn)題發(fā)揮了重要作用?；谏鲜隼碚?一些學(xué)者對(duì)氣固兩相流中管內(nèi)的非球形顆粒流動(dòng)進(jìn)行研究。Hilton 等人[8]采用超橢球模型研究了在水平氣力輸送管道中,顆粒形狀對(duì)稀釋流到段塞流這一轉(zhuǎn)換過(guò)程中的影響。Kruggel-Emden 等人[9]除了基本的球形顆粒外還建立了4 種多面體模型,討論了這5 種顆粒在氣力帶動(dòng)過(guò)程中流動(dòng)特性的差異。周甲偉[10]研究了顆粒形狀和旋流強(qiáng)度對(duì)彎管磨損的影響,其顆粒形狀采用組合顆粒法進(jìn)行建模,此外,他們還研究了煤塊顆粒破碎過(guò)程。Zeng 等人[11]則分析了彎管中V 型磨損帶形成的原因,并且比較了不同多面體形式顆粒對(duì)彎管磨損的影響。

本研究采用CFD-DEM 方法耦合切向撞擊能量模型[12](shear impact energy model,SIEM)來(lái)研究不同形狀顆粒對(duì)90 °彎管沖蝕磨損特性的影響。首先將標(biāo)準(zhǔn)算例與實(shí)驗(yàn)數(shù)據(jù)進(jìn)行比較驗(yàn)證上述方法的可靠性,其次,對(duì)比不同形狀顆粒對(duì)彎管的磨損情況。此外,將氣速這一重要因素考慮在內(nèi),對(duì)比不同氣速條件下顆粒形狀對(duì)彎管磨損的影響程度。

1 數(shù)學(xué)模型

1.1 連續(xù)相控制方程

由于本論文算例均為稀氣-固兩相流,為了節(jié)省計(jì)算時(shí)間,故均采用單向耦合方法,流體控制方程如下所示:

連續(xù)相方程為

動(dòng)量守恒方程為

式中,ρf為流體密度(假設(shè)流體不可壓縮,所以為常數(shù)),u為流體流速,p為流體壓力,μeff為流體有效粘度,g為重力加速度。

1.2 離散相控制方程

離散相運(yùn)動(dòng)可分為平移和旋轉(zhuǎn),可以用牛頓第二定律來(lái)計(jì)算,具體如下:

式中,m為顆粒質(zhì)量,dv/dt為加速度,Fc為接觸力,Fd為流體對(duì)固體的曳力,FlS為Saffman 升力,FlM為Magnus 升力,Fb為浮力,I為慣性矩,dω/dt為角加速度,Tc為接觸扭矩,Tf為流體作用下的扭矩。采用Cundall 等[13]提出的線彈性-阻尼模型計(jì)算Fc:

式中,Fc,n為接觸力的法向分量,Fc,t為接觸力的切向分量,kn為法向彈性系數(shù),kt為切向彈性系數(shù),δn為法向位移,δt為切向位移,ηn為法向阻尼系數(shù),ηt為切向阻尼系數(shù),vn為法向速度,vt為切向速度。阻尼系數(shù)可根據(jù)Ting 等人[14]提出的方程得到,當(dāng)接觸力2 個(gè)分量滿足式(8),則接觸力切向分量按照式(9)求解。

式中,fs為滑動(dòng)摩擦系數(shù)。

由于非球形顆粒在運(yùn)動(dòng)中存在各向異性,引入局部坐標(biāo)系來(lái)獲得非球形的轉(zhuǎn)動(dòng)慣量I,但接觸扭矩Tc仍按全局坐標(biāo)系計(jì)算求解,具體如式(10)所示。

式中,L為從顆粒中心到接觸點(diǎn)的矢量。

1.3 顆粒形狀模型

除了常規(guī)的球形模型以外,超橢球模型、組合顆粒模型以及多面體模型[6]也被廣泛應(yīng)用于DEM 建模方法中。本文采用由Barr[15]提出的超橢球方程來(lái)對(duì)非球形顆粒進(jìn)行建模,具體公式如下:

式中,s1和s2為顆粒形狀系數(shù),a、b、c為軸線方向的半軸長(zhǎng)。

1.4 流體和固體耦合模型

在單向耦合中,只考慮流體對(duì)固體所產(chǎn)生的力,其中包括曳力Fd、Saffman 升力FlS、Magnus 升力FlM和扭矩Tf。

曳力Fd采用Di Felice[16]所提出的曳力模型,具體如下:

式中,Fd0為顆粒無(wú)其他粒子作用時(shí)流體對(duì)其產(chǎn)生的阻力,αf為流體體積分?jǐn)?shù),CD為流體曳力系數(shù),dp為顆粒直徑,vf為氣體速度,vp為顆粒速度,Rep,α為流體雷諾數(shù),μf為流體粘度。

顆粒在流體作用下的線速度和角速度都較大,不可忽略,故Saffman 升力[17]必須考慮,其表達(dá)式如下:

式中,ωf為流體角速度,ClS為升力系數(shù),可根據(jù)Mei[18]的研究成果得到,具體如式(18)所示,Rep為顆粒雷諾數(shù)。

Magnus 升力是由顆粒旋轉(zhuǎn)而引起的,具體如下:

式中,Rep為Magnus 升力作用下的雷諾數(shù);Rer為顆粒旋轉(zhuǎn)的雷諾數(shù);ωp為顆粒角速度;ClM為升力系數(shù),可根據(jù)Oesterlé 等人[19]的研究成果得到,具體如式(23)所示。

顆粒浮力Fb如下:

式中,Vp為單個(gè)顆粒體積。

顆粒的轉(zhuǎn)動(dòng)可根據(jù)Rubinow 等人[20]的研究成果得到,具體如下:

式中,CR為旋轉(zhuǎn)扭矩系數(shù),可根據(jù)Rubinow 等人[20]以及Dennis 等人[21]的研究?jī)?nèi)容中獲得,具體如式(27)所示。

1.5 磨損模型

采用SIEM 模型[12]來(lái)計(jì)算彎管磨損情況,磨損模型如下所示:

式中,W為材料被碰撞后所損失的體積,EShear為顆粒碰撞壁面后顆粒損失的切向動(dòng)能,p為壁面塑性流動(dòng)壓力,t0為碰撞開(kāi)始時(shí)間,t為接觸持續(xù)時(shí)間,Fc,t為顆粒與壁面接觸時(shí)對(duì)壁面的切向力,vt為切向速度。

2 數(shù)值方法

2.1 幾何模型

本文所采用的幾何模型以及實(shí)驗(yàn)數(shù)據(jù)均來(lái)自Chen 等人[22]的研究?jī)?nèi)容及成果,實(shí)驗(yàn)裝置示意圖如圖1 所示。實(shí)驗(yàn)管材選用鋁,質(zhì)量流量為40 lb/d(2.08 ×10-4kg/s)的沙子從顆粒入口注入后在空氣壓縮機(jī)的帶動(dòng)下,沖擊測(cè)試單元處的彎管,并從出口流出,氣體速度為150 ft/s(45.72 m/s),彎管特定位置的磨損率通過(guò)輪廓測(cè)量?jī)x來(lái)獲取。

圖1 Chen 等人[22]所采用的實(shí)驗(yàn)裝置示意圖

計(jì)算中取實(shí)驗(yàn)示意圖中顆粒經(jīng)過(guò)彎管豎直向上沖擊彎管的直管以及后續(xù)彎管部分作為模擬模型,具體尺寸如圖2 所示。

圖2 管道幾何參數(shù)示意圖

2.2 具體參數(shù)設(shè)置

連續(xù)相采用商用Fluent 軟件進(jìn)行模擬,流體運(yùn)動(dòng)基于壓力求解器,用SIMPLEC 來(lái)求解氣體控制方程,對(duì)流相采用QUICK 方法。離散相固體顆粒采用商用DEM 軟件DEMSLab 來(lái)進(jìn)行模擬,具體參數(shù)如表1 所示,彎管網(wǎng)格劃分如圖3 所示。

表1 模擬參數(shù)表

圖3 模擬采用的網(wǎng)格

2.3 網(wǎng)格無(wú)關(guān)性分析

由于本文只針對(duì)彎管磨損進(jìn)行分析,因此只對(duì)彎管處網(wǎng)格無(wú)關(guān)性進(jìn)行分析。將彎管部分網(wǎng)格進(jìn)行調(diào)整,分別為21 008、23 028、25 048、27 068、29 088個(gè)網(wǎng)格單元,在相同邊界條件下得到彎管外拱處磨損率,并與實(shí)驗(yàn)數(shù)據(jù)進(jìn)行比較,結(jié)果如圖4 所示。從模擬結(jié)果上看,采用25 048 個(gè)網(wǎng)格的計(jì)算精度已經(jīng)達(dá)到要求,因此選擇該規(guī)格網(wǎng)格繼續(xù)進(jìn)行本文研究。

圖4 不同網(wǎng)格數(shù)量對(duì)應(yīng)彎管外拱磨損情況與Chen 等人[22]的實(shí)驗(yàn)數(shù)據(jù)對(duì)比結(jié)果

3 結(jié)果分析與討論

模擬中顆粒分布如圖5 所示(顆粒被放大10倍),從圖中可以看出,2 條黑色箭頭為顆粒在沖擊彎管、經(jīng)過(guò)壁面反彈后顆粒流的主要流向,這時(shí)顆粒在氣力的帶動(dòng)下,二次沖擊彎管、進(jìn)而造成V 型磨損帶,這與Zeng 等人[11]的研究結(jié)論一致。

圖5 驗(yàn)證算例的管道中顆粒分布

3.1 顆粒形狀對(duì)彎管沖蝕磨損的影響

為了研究顆粒形狀對(duì)彎管沖蝕磨損的影響,本文基于Chen 等人[22]的實(shí)驗(yàn)條件,在45.72 m/s 的氣速條件下(管道流場(chǎng)如圖6 所示),采用超橢球模型建立6種顆粒模型,具體形狀及尺寸如圖7(a)～ (f)所示。其中a、b、c半軸長(zhǎng)單位為m,后文圖中的(a)～(f)對(duì)應(yīng)6 種不同形狀的顆粒模型,并以上述6 種顆粒模型為研究對(duì)象探究不同顆粒形狀對(duì)彎管磨損的影響程度,單個(gè)不同形狀顆粒的質(zhì)量均相同。

圖6 驗(yàn)證算例的管道中速度場(chǎng)與壓力場(chǎng)云圖分布

圖7 顆粒形狀示意圖及其具體尺寸

圖8 上圖為不同形狀顆粒對(duì)彎管拱背的磨損率,下圖則表示不同形狀顆粒在氣力的帶動(dòng)下多數(shù)顆粒在與壁面碰撞前其速度能達(dá)到的最大值以及顆粒對(duì)彎管內(nèi)壁造成多次碰撞后的最低速度。總體來(lái)看,不同形狀顆粒造成的拱背處磨損率變化趨勢(shì)并無(wú)明顯差異,磨損率會(huì)隨著角度的增加,在40 °～50 °之間達(dá)峰值,隨后會(huì)迅速下降。

圖8 不同形狀顆粒造成的彎管拱背處磨損率以及對(duì)應(yīng)情況下顆粒的最大速度以及最小速度

對(duì)比圖8(a)～(d)中4 種顆粒模型下彎管拱背的磨損率,可以看出,隨著形狀系數(shù)的增加(顆粒球形度為球形顆粒表面積與與顆粒相同體積的顆粒表面積之比,正方體約為0.806;顆粒越接近方形,球形度越低),彎管的最大磨損率會(huì)先小幅下降,再逐漸上升;圖8 還顯示出彎管的最大磨損位置會(huì)朝大角度方向偏移,最終會(huì)穩(wěn)定在47 °左右。最大磨損率的變化趨勢(shì)與顆粒最大速度的變化趨勢(shì)相同,即隨著形狀系數(shù)的增大顆粒的最大速度先有所下降,然后逐漸增大,這也直接說(shuō)明彎管最大磨損率與顆粒最大速度有關(guān),而顆粒的最大速度與顆粒受到的曳力有關(guān),曳力越大,沖擊速度也越大。此外,圖8中不同形狀顆粒的速度極差也能反映顆粒經(jīng)過(guò)彎管對(duì)其產(chǎn)生的切向撞擊能的相對(duì)大小?？傮w看,隨著形狀系數(shù)的增加,顆粒的速度極差越大,管壁受到的切向撞擊能也越大,這也正好符合圖9 中(a)～(d)總磨損率逐漸上升這一規(guī)律。盡管(b)對(duì)彎管拱背最大磨損率略低于(a),但(b)的總磨損率卻較大,這也恰好說(shuō)明隨著形狀系數(shù)的增大,顆粒流對(duì)彎管的磨損范圍也會(huì)相應(yīng)增大。

圖9 不同形狀顆粒造成的彎管總磨損率

通過(guò)顆粒在管壁附近的姿態(tài)可判斷其對(duì)管壁的沖蝕磨損方式。歐拉角便是描述顆粒姿態(tài)的方法,XYZ坐標(biāo)系先后繞著Z軸、X′軸以及Z″軸旋轉(zhuǎn)Ψ度、θ度以及φ度得到局部坐標(biāo)系X?Y?Z″,如圖10所示,Ψ、θ、φ分別為歐拉角3 個(gè)分量進(jìn)動(dòng)角、章動(dòng)角和自轉(zhuǎn)角,這樣顆粒在某一位置的姿態(tài)便可通過(guò)歐拉角表示。通過(guò)統(tǒng)計(jì)多組數(shù)據(jù)中顆粒在彎管拱背附近的平均姿態(tài)以及歐拉角3 個(gè)分量的標(biāo)準(zhǔn)差來(lái)表示大多數(shù)顆粒在靠近管壁時(shí)的狀態(tài)。將90 °彎管分成18 份,分別對(duì)每段位置的顆粒歐拉角進(jìn)行統(tǒng)計(jì),由于球形顆粒的各向同性,故不需考慮,統(tǒng)計(jì)結(jié)果如圖11 所示。通過(guò)顆粒3 個(gè)分量的平均值,可以畫(huà)出3 種不同形狀系數(shù)顆粒在彎管拱背附近XY面的投影情況,如圖12 所示。不同形狀顆粒的姿態(tài)在5 °～20 °、65 °～85 °之間差異明顯,當(dāng)顆粒不為球形時(shí),彎管的磨損是顆?；瑒?dòng)摩擦和滾動(dòng)摩擦共同作用的結(jié)果,隨著顆粒逐漸趨于方形,顆粒發(fā)生滑動(dòng)摩擦的造成的磨損量增多,彎管的磨損量也會(huì)升高?？傊w粒形狀的改變直接影響著顆粒作用于管道壁面的撞擊能大小、作用范圍以及通過(guò)影響顆粒所受曳力間接影響彎管磨損。

圖10 歐拉角示意圖[12]

圖11 不同形狀系數(shù)顆粒歐拉角3 個(gè)分量的均值和標(biāo)準(zhǔn)差

圖12 不同形狀系數(shù)顆粒參數(shù)獲取位置以及顆粒在XY 面的投影

進(jìn)一步分析長(zhǎng)橢球(e)和扁橢球(f)對(duì)彎管磨損的結(jié)果,由于長(zhǎng)橢球顆粒與扁橢球顆粒繞Z軸旋轉(zhuǎn)任意角度都不會(huì)改變姿態(tài),因此僅需考慮顆粒歐拉角前2 個(gè)分量的平均值和標(biāo)準(zhǔn)差,結(jié)果如圖13 所示。兩種顆粒在彎管拱背附近XY面的投影情況,如圖14 所示。從圖13 可以看出,拱背不同位置的顆粒,盡管歐拉角2 個(gè)分量的標(biāo)準(zhǔn)差都較大,但由于顆粒關(guān)于Z軸與X軸中心對(duì)稱,圖14 中顆粒的投影仍可以代表大多數(shù)顆粒在XY面的投影。從圖13可以看出,拱背附近顆粒,在氣力帶動(dòng)下,長(zhǎng)橢球顆粒更易發(fā)生滑動(dòng)摩擦,扁橢球更易發(fā)生滾動(dòng)摩擦,顆粒發(fā)生滑動(dòng)摩擦損失的切向動(dòng)能多于滾動(dòng)摩擦,造成彎管的磨損也更嚴(yán)重。從圖8中(e)和(f)的速度極差也可說(shuō)明長(zhǎng)橢球顆粒對(duì)管壁的切向撞擊能要高于扁橢球,因此在38 °～55 °位置處,長(zhǎng)橢球造成的磨損率高于扁橢球,而在其他位置略低。這可能與顆粒-管壁的接觸次數(shù)和接觸時(shí)間有關(guān)。圖9(d)的總磨損率明顯高于圖9(e),仍可以由2種顆粒的速度極差來(lái)判斷,即由于單個(gè)(d)顆粒損失的動(dòng)能高于(e),經(jīng)過(guò)多次累加,就會(huì)造成形狀接近方形的顆粒對(duì)彎管的磨損多于長(zhǎng)橢球顆粒。

圖13 長(zhǎng)橢球與扁橢球顆粒進(jìn)動(dòng)角與章動(dòng)角的均值和標(biāo)準(zhǔn)差

圖14 顆粒參數(shù)獲取位置以及顆粒在XY 面的投影

圖15 為不同形狀顆粒對(duì)整個(gè)彎管磨損率影響的情況。當(dāng)顆粒流為球形顆粒時(shí),磨損嚴(yán)重區(qū)域呈橢圓形,管壁兩側(cè)的V 形磨損帶還不明顯;在形狀系數(shù)增大過(guò)程中,橢圓形變大,隨后磨損嚴(yán)重區(qū)域沿著管壁兩側(cè)擴(kuò)展,在形狀系數(shù)為8 時(shí)擴(kuò)展區(qū)域最為明顯。圖15 結(jié)果也表明長(zhǎng)橢球?qū)澒茉斐傻哪p低于形狀系數(shù)為8 的顆粒,但高于扁橢球。模擬云圖也可以說(shuō)明,隨著形狀系數(shù)的增加,彎管磨損嚴(yán)重區(qū)域逐步擴(kuò)大,V 形磨損帶也越發(fā)明顯,長(zhǎng)橢球?qū)澒墚a(chǎn)生的沖蝕磨損高于扁橢球,這也印證了之前的結(jié)論。圖16 為不同顆粒在管內(nèi)的位置分布,不同顆粒流豎直向上沖擊彎管,隨后沿著彎管外拱以及管壁兩側(cè)流出,造成了不同的沖蝕形貌,同樣可以觀察到上述現(xiàn)象。

圖15 不同形狀顆粒對(duì)彎管沖蝕磨損云圖

3.2 不同氣速下顆粒形狀對(duì)彎管磨損的影響

以圖16(a)、(d)、(e)和(f)4 種顆粒形狀為研究對(duì)象,以v=45.72 m/s 和v=30.48 m/s 這2 種氣速為入口邊界條件來(lái)對(duì)比不同氣速下顆粒形狀對(duì)彎管的磨損程度。

圖17 為不同顆粒在2 種氣速條件下對(duì)彎管拱背的磨損率,左圖為v=45.72 m/s 氣速條件,右為v=30.48 m/s 氣速條件。從圖中可以看出,無(wú)論改變氣速還是顆粒形狀,彎管拱背的磨損趨勢(shì)都并無(wú)太大差別,且氣速變動(dòng)并未引起最大磨損率位置發(fā)生改變。

圖18 為不同形狀顆粒對(duì)彎管總磨損率以及對(duì)應(yīng)彎管的磨損云圖。對(duì)比同種顆粒形狀在不同速度下的彎管總磨損率的差值可以看出,氣速的增大使得更接近方形顆粒的算例(d)中彎管磨損率的升高最為明顯。綜合來(lái)看,氣速的增大會(huì)使彎管各部分磨損量均勻上升,彎管磨損分布并無(wú)較大改變,磨損分布的不同主要是顆粒形狀不同所導(dǎo)致。

4 結(jié)論

(1)相同質(zhì)量流量以及氣速條件下,同質(zhì)量的非球形顆粒較球形顆粒造成的彎管內(nèi)壁磨損率更大,且顆粒形狀越接近方形,彎管磨損越嚴(yán)重,磨損范圍越廣。

(2)相比于扁橢球,長(zhǎng)橢球顆粒造成的彎管磨損更大,這主要由于長(zhǎng)橢球顆粒在管壁面更易發(fā)生滑動(dòng)摩擦,扁橢球易發(fā)生滾動(dòng)摩擦。

(3)形狀近似為方形的顆粒進(jìn)出彎管損失的動(dòng)能要多于長(zhǎng)橢球,這就導(dǎo)致長(zhǎng)橢球顆粒對(duì)彎管沖擊造成的彎管質(zhì)量損失較少。

(4)氣速的提升僅會(huì)均勻增大不同形狀顆粒對(duì)彎管的磨損率,而磨損分布以及最大磨損率位置并無(wú)較大差異。