李震東 許 磊 包士毅

(浙江工業(yè)大學化工機械設計研究所 杭州 310023)

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

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

0 引言

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

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

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

1 數(shù)學模型

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

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

連續(xù)相方程為

動量守恒方程為

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

1.2 離散相控制方程

離散相運動可分為平移和旋轉(zhuǎn),可以用牛頓第二定律來計算,具體如下:

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

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

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

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

式中,L為從顆粒中心到接觸點的矢量。

1.3 顆粒形狀模型

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

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

1.4 流體和固體耦合模型

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

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

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

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

式中,ω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為單個顆粒體積。

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

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

1.5 磨損模型

采用SIEM 模型[12]來計算彎管磨損情況,磨損模型如下所示:

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

2 數(shù)值方法

2.1 幾何模型

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

圖1 Chen 等人[22]所采用的實驗裝置示意圖

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

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

2.2 具體參數(shù)設置

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

表1 模擬參數(shù)表

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

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

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

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

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

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

圖5 驗證算例的管道中顆粒分布

3.1 顆粒形狀對彎管沖蝕磨損的影響

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

圖6 驗證算例的管道中速度場與壓力場云圖分布

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

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

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

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

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

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

圖10 歐拉角示意圖[12]

圖11 不同形狀系數(shù)顆粒歐拉角3 個分量的均值和標準差

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

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

圖13 長橢球與扁橢球顆粒進動角與章動角的均值和標準差

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

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

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

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

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

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

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

4 結(jié)論

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

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

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

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