含蠟天然氣孔板式節流閥內蠟沉積特性數值模擬

張永星, 顧鴻博, 陳曉玲

(1.中國石油大學(北京)機械與儲運工程學院, 北京 102249; 2.中國石油大學(北京)克拉瑪依校區工學院, 克拉瑪依 834000)

含蠟天然氣是指含有較高碳數的正烷烴,或含有由直鏈、支鏈烷烴所組成高分子碳氫化合物等復雜混合物的天然氣。在溫度降低到析蠟溫度以下時,會有蠟晶析出,蠟晶經過不斷運動,最終在一定條件下,氣相中的蠟顆粒會沉積于管道內壁。

目前國內外天然氣集輸系統均面臨著不同程度的蠟沉積問題。以新疆克深2氣田處理裝置為例,原料氣在經過J-T閥時,由于節流膨脹導致低溫分離器發生蠟堵,最終影響裝置的正常運行[1]。在克拉蘇氣田大北區塊天然氣處理廠、塔河九區天然氣處理廠、克拉2氣田、塔里木油田凝析氣藏等同樣發生了不同程度的蠟堵現象[2-5]。伊朗南部的南帕爾斯天然氣加工廠站場設備內也發現有蠟顆粒的沉積,引起站內設備的安全問題,影響生產進度[6]。因此,天然氣管道處理站管道蠟沉積不僅會減小管道的有效流通面積,降低管道的輸送能力,嚴重時還會造成蠟堵事故[7],引發安全問題。雖然含蠟天然氣在集輸系統各個環節都有著不同程度的蠟沉積現象,但由于節流閥前后會引起氣體的壓力、溫度大幅變化、從而使閥內蠟沉積所帶來的問題顯得更為突出。

目前針對含蠟天然氣在節流設備中的蠟沉積問題的研究較少。李俊逸等[8]采用數值模擬的方法研究得到了硫沉積速率隨氣流進口速度的增大而先增大后減小的規律。劉成川等[9]建立了元素硫在多場耦合作用下的沉積模型,認為沉積主要取決于壓力。溫度等參數變化。劉巖等[10]采用計算流體力學與離散單元法相耦合,認為顆粒在沉積過程中主要受二次流的影響,大部分沉積在壁面兩側附近。Sun等[11]通過Fluent進行氣固兩相模擬,并認為顆粒運動與斯托克斯數緊密相關。Wang等[12]借助DPM(deformable part model)模型對天然氣內水合物顆粒進行模擬,得到了管道底部顆粒物濃度較大、頂部顆粒物濃度較小的趨勢規律。

然而,目前針對傳統孔板式節流閥的研究主要集中在焦耳湯姆遜效應及其水合物生成等現象[13-16],圍繞節流閥內蠟沉積規律的研究相對較少,對孔板式節流閥內蠟沉積規律認識不明確,難以有效地對實際生產中蠟沉積進行主動防治。因此,開展節流閥設備蠟沉積問題的研究,揭示節流閥內高流速下的蠟沉積規律,對提高含蠟天然氣安全高效運行具有重要意義。為此,現借助計算流體力學軟件,通過UDF(user-defined function)編寫自定義函數,實現蠟的沉積與剝離,建立孔板式節流閥的蠟沉積仿真模型,分析孔板式節流閥內蠟沉積特性。

1 數值模型

1.1 模型假設

由于天然氣中含蠟量較低,故假設流動為稀疏的氣-固兩相流,不考慮由于溫降可能導致凝析油的產生;顆粒為球形,且符合正態分布;顆粒主要受到流體的曳力、重力、阻力、浮力等作用;由于在較低溫度的條件時,蠟的黏度變化不大,暫時忽略溫度對于黏度的影響。

1.2 數學模型

1.2.1 連續相湍流方程

k-ε(k為湍動能,ε為湍流耗散率)是目前使用最廣泛的湍流方程[17],對流場求解具有普遍的適用性。該模型通過分別求解湍動能及湍流耗散率方程來確定湍流長度和時間尺度,而湍流耗散率方程則是通過物理推理并類比相似方程得出[18],其具體控制方程見文獻[19]。

1.2.2 離散相運動方程

采用DPM模型對顆粒運動軌跡進行追蹤,以x方向為例,顆粒所受的平衡方程為

(1)

式(1)中:t為時間,s;u為天然氣流動速度, m/s;us為固體顆粒速度,m/s;ρ為天然氣密度,kg/m3;ρs為固體顆粒密度,kg/m3;gx為x方向的重力加速度分量,m/s2;μ為天然氣動力黏度,Pa·s;ds為蠟顆粒直徑,μm;CD為顆粒阻力系數;Fx為顆粒所受x方向的外力,N。

1.3 臨界速度模型

臨界速度模型[20]通過研究顆粒的彈性性質,以顆粒物與壁面撞擊時的法向速度作為判斷依據,當顆粒的法向速度大于臨界速度,顆粒就會反彈;反之,顆粒就會發生黏附。圖1所示為單個顆粒在流體中的受力示意圖。

FL為流體作用在顆粒上的吸引力;a為與壁面接觸面積的半徑;Fpo為黏附力;FD為曳力;Dp/2-α為顆粒圓心到壁面距離圖1 單個顆粒受力示意圖Fig.1 Schematic diagram of individual particle stress

具體計算公式為

(2)

(3)

(4)

(5)

式中:Vcr為臨界速度;E為混合的楊氏模量;Es為壁面的楊氏模量;Ep為顆粒的楊氏模量;νs為壁面的泊松比;νp為顆粒的泊松比;Dp為顆粒直徑;ρp為顆粒的密度;R為運動恢復系數,這里取0.7。

Soltani等[21]通過大量實驗對不同的剝離機制進行研究,發現圓形顆粒發生剝離的主要方式為滾動。Dasetal[22]同樣發現,滾動剝離是圓球形顆粒剝離的主要形式。故假設當蠟晶顆粒發生剝離時,只發生滾動剝離,且為單個顆粒剝離,直接剝離與滑動剝離不會發生。滾動剝離方程[23]為

(6)

式(6)中:α為顆粒與壁面接觸后顆粒直徑所減少的部分;FL為流體作用在顆粒上的剪切力;a為顆粒與壁面的接觸半徑;針對球形顆粒,FD為所受到的曳力,可以表示為

(7)

式(7)中:uc為顆粒所在位置的流速;Acp為顆粒截面積;f為修正系數;Cu為Cunningham修正系數,由Sundaram[24]提出,具體表達式為

(8)

式(8)中:Kn為克努森數,用來表征稀疏效應的影響,不可忽略[25]。

Fpo具體表達式為

(9)

式(9)中:Fpo為黏附力;WA為黏著功,J/m2。黏著功取決于顆粒與壁面材料,通常通過實驗得到,關于對蠟顆粒黏著功研究較少,在此假設與硅顆粒相同,取0.039 J/m2。

2 模型驗證

2.1 物理模型及網格劃分

以文獻[26]中的孔板式節流閥為研究對象,結構見圖2。其中直管段直徑D=76 mm,孔板的孔徑d=30 mm,孔板的厚度t=10 mm,左側為節流閥入口,右側為出口。為保證湍流穩定性,在節流閥進口和出口分別增加10D和15D的直管段。采用ICEM軟件對計算域進行六面體網格劃分,并插入膨脹層以加密壁面附近的網格,如圖3所示。

圖2 孔板式節流閥及路徑示意Fig.2 Orifice type throttle valve and path diagram

圖3 孔板式節流閥網格模型Fig.3 Orifice throttle valve mesh model

2.2 求解策略

計算過程中將二階中心差分格式用于離散擴散項,二階迎風差分格式用于其他項,時間步長為1×10-4,先對連續相進行穩態計算,收斂后加入顆粒相進行非穩態耦合計算。采用SIMPLE(semi-implicit method for pressure-linked equations)算法模擬流場特征[27],并通過UDF實現蠟晶顆粒的沉積、反彈、剝離等特征。

2.3 定解條件與工況設定

連續相采用壓力入口及壓力出口,壁面邊界為無滑移壁面。天然氣密度采用理想氣體模型,動力黏度為1.087×10-5kg/(m·s)。在距入口0.76 m處建立顆粒射流源平面,其中顆粒質量流量為1.3×10-6kg/s。采用58號石蠟材料,石蠟顆粒密度為1 980 kg/m3、比熱容為986 J/ (kg·K)。同時對顆粒加載虛擬質量力及壓力梯度力。離散相進、出口邊界處采用逃逸邊界,壁面采用UDF自定義函數。氣流入口壓力為5～8 MPa。設定節流閥孔徑為20～60 mm,喉道長度為10～40 mm。數值模擬工況設定見表1。

表1 數值模擬工況

2.4 網格無關性驗證

為了確保數值模擬結果獨立于網格,以喉道內中軸線處速度值為參考,分別對72萬、100萬、128萬、142萬以及173萬網格數進行無關性分析,結果見圖4。

圖4 網格無關性分析Fig.4 Mesh independence analysis

通過圖4分析可知,當沿流動方向距離大于0.76 m后,72萬網格數時的速度明顯低于其他網格數時的速度,而對于他網格數,速度曲線基本重合。綜合考慮計算精度及計算速度后,最終選用128萬網格數進行后續數值計算。

以文獻[26]中五處溫度監測點作為參考,將模擬結果與文獻中實驗節點的溫度值進行比較,具體結果見表2。

表2 文獻溫度與模擬溫度對比

通過表2可知,Fluent模擬結果與實驗結果最大誤差為4.7%,誤差均不超過5%,說明所選取的模型適用于本研究。

2.5 蠟沉積模型驗證

1995年,美國部分學者發現在孔板式流量計中存在固體顆粒,會沉積到管壁,影響測量準確性[28]。為驗證臨界速度模型的準確性,在2.1節孔板式節流閥模型的基礎上,從入口(x=0)處注入蠟晶顆粒,將數值模擬與現場中觀察沉積位置進行對比。圖5為基于臨界速度模型由數值模擬得到的蠟沉積分布,而圖6為文獻[28]中實驗觀察到的蠟晶沉積位置。

圖5 模擬蠟沉積位置Fig.5 Simulate wax deposition locations

圖6 現場固體顆粒沉積位置[28]Fig.6 On-site solid particle deposition location[28]

2020年,Peng等[29]學者對頁巖氣開采過程中孔板內顆粒沉積進行分析。為進一步驗證蠟沉積模型的準確性,將數值模擬與實驗中得到的沉積位置進行對比。圖7為基于臨界速度模型由數值模擬得到的沉積分布,圖8為文獻[29]中實驗得到的沉積位置。

圖7 模擬顆粒沉積位置Fig.7 Simulate particle deposition location

圖8 文獻中實驗沉積位置[29]Fig.8 Experimental deposition location in literature [29]

從數值模擬得到的圖5可以看出,管壁四周均有沉積的蠟顆粒,但大部分蠟顆粒沉積在管道下方,與文獻中的實驗結果圖6沉積趨勢相同。從模擬得到的圖7可以看出,孔板內均出現顆粒沉積,隨著氣相速度增大,沉積量逐漸增加,并且在孔板入口截面處出現少量沉積,與文獻[29]中實驗結果圖8沉積趨勢相同。通過數值模擬與實驗結果的比較,發現臨界速度模型可以較為準確的描述蠟顆粒沉積情況。因此,文中后續模擬分析將采用臨界速度模型對蠟沉積規律進行研究。

3 結果分析

3.1 壓力對沉積率的影響

圖9所示為最大沉積量隨壓力變化趨勢,從中可以看出,當管內入口壓力從5 MPa增加到6 MPa,最大沉積量由1.06×10-6kg增加至1.67×10-6kg,增加了36.52%,當壓力從6 MPa增加到7 MPa時,最大沉積量由1.67×10-6kg增加至1.92×10-6kg,增加了13.02%,而當壓力從7 MPa增加到8 MPa時,最大沉積量則由1.92×10-6kg增加至2.11×10-6kg,只增加了9%。發現隨著進口壓力的增加,管道內壁的最大沉積量逐漸增加,但沉積量增加的趨勢逐漸變緩。蠟沉積量隨入口壓力呈現上述變化趨勢可能主要由以下兩個方面的原因。一方面由于隨著壓力的增加,管道內的流速也隨之增加,從而使管內湍流強度增加,增大了蠟顆粒撞擊壁面形成沉積的幾率,所以蠟沉積量隨著入口壓力的增大而增大;另一方面,隨著壓力增大,顆粒法相速度隨之增大,達到了部分較小粒徑的沉積條件,但由于顆粒粒徑變小,且由于粒徑呈正態分布,所以小顆粒數量相對較少,因此沉積量增加的趨勢逐漸變緩。

圖9 沉積量變化趨勢Fig.9 Variation trend of deposition rate

圖10所示為管道內蠟沉積的位置。從圖10可以看出,蠟顆粒雖然在節流閥后的整個管道內均有沉積,但主要沉積位置發生在喉道出口(x=0.77～0.775 m)處以及距喉道出口0.17 m(x=0.94 m)左右的管壁處。第二處沉積位置隨著壓力增大逐漸向管道出口移動,壓力增加1 MPa,沉積位置向后移動約0.005 m,并且沉積寬度為0.01 m左右。當壓力為5 MPa時,沉積位置主要在0.935～0.945 m;當壓力為6 MPa時,沉積位置主要在0.938～0.948 m;當壓力為7 MPa時,沉積位置主要在0.945～0.955 m;當壓力為8 MPa時,沉積位置主要在0.948～0.958 m。管道內蠟顆粒主要沉積在上述兩個位置可能是湍動能強度和渦流共同作用的結果。

圖10 模擬蠟顆粒沉積位置示意圖Fig.10 Schematic diagram of simulated wax particle deposition location

圖11所示為沿流動方向中軸線處湍動能分布情況,從中可以看出湍動能在喉道出口(x=0.77 m)和距喉道出口0.17 m處附近,即兩處主要的沉積位置,出現了兩個明顯的峰值。說明在該處速度的波動幅度更大,并且這種劇烈程度隨著壓力的增加而增大,從而導致上述兩處位置蠟沉積量隨著壓力增加而增加。同時發現第二處湍動能峰值位置隨著壓力的增大逐漸后移,該趨勢與主要沉積位置隨著壓力的增大相應后移相一致。另外還發現在喉道出口(x=0.77 m)位置處,隨著壓力的增加,湍動能的增大幅度逐漸減小,并在7～8 MPa時十分接近,而在距喉道出口0.17 m處附近,隨著壓力的增加,湍動能的增大幅度差別不大,從而導致隨著壓力的增大后者的湍動能大于前者的湍動能。這也是導致在距喉道出口0.17 m處蠟沉積量大于喉道出口處沉積量的原因之一。

圖11 沿流動方向湍動能變化Fig.11 Variation of turbulent kinetic energy along the flow direction

圖12為兩個主要沉積截面處垂直于流動方向的湍動能分布。從圖12(b)中可以看到,湍動能出現中間下凹兩端凸起的駝峰形主要是因為軸線處速度較為穩定,隨著氣流向徑向方向流動,形成漩渦,不穩定程度加劇,湍動能逐漸增大,達到峰值。圖13所示為所研究最低壓力(5 MPa)和最高壓力(8 MPa)下管道內的流場圖。氣體經過節流閥后,速度明顯增加,并在沉積處產生兩處渦流。隨著壓力增大,渦流影響范圍沒有發生明顯變化,但變化程度更為劇烈。當顆粒隨流體運動過程中,在渦流位置由于流體方向發生突然改變,顆粒在慣性作用下撞擊到壁面,發生沉積。因此,隨著壓力增大,顆粒與壁面發生撞擊的幾率越大,最終導致管內沉積量增加。

圖12 沉積截面處垂直于流動方向的湍動能變化Fig.12 Variation of turbulent kinetic energy at the deposition section perpendicular to the flow direction

3.2 蠟顆粒直徑的影響

圖14所示為顆粒粒徑與最大沉積量之間的關系。可以看出,在節流直徑、喉道長度一致時,最大沉積量隨粒徑增大而增大,并在粒徑為50 μm時發生陡增。以入口壓力為5 MPa時為例,粒徑由30 μm增加到50 μm時,最大沉積量由0.067 7×10-6kg增加到0.142 6×10-6kg,增加了52.52%;粒徑由80 μm增加到100 μm時,最大沉積量由0.739 4×10-6kg增加到2.145 4×10-6kg,增加了65.53%。這是因為隨著粒徑的增大,蠟顆粒具有更大的慣性。因此,當渦流方向突然改變時,在慣性的作用下,較大粒徑的蠟顆粒更易保持原來的方向運動,脫離原流場,從而撞擊至節流閥管壁的幾率增大,導致更大的沉積量。

圖14 粒徑對沉積量的影響Fig.14 Effect of particle size on deposition amount

3.3 節流孔徑的影響

圖15所示為入口壓力為5 MPa時,節流孔直徑與最大沉積量的關系。分析可知,節流孔徑增大,管壁上蠟顆粒的最大沉積量逐漸下降。當孔徑由20 mm增加至30 mm時,管壁最大沉積量由1.293×10-6kg降至0.952 7×10-6kg,下降26.37%;當孔徑繼續增大至40 mm時,管壁沉最大沉積量由0.952 7×10-6kg降至0.493 9×10-6kg,下降48.15%;當管徑為50 mm和60 mm時,沉積量接近為0。從不同孔徑下流線示意圖16可以看出,由于管內為超音速流動,隨著孔徑的增大,速度也隨之增大。并且渦流明顯變小,顆粒大部分隨氣流向前流動,顆粒撞擊壁面幾率減小,節流孔出口壁面處對蠟顆粒的阻礙作用減小,因此在大孔徑下,蠟顆粒更容易被氣流帶走,而不發生沉積。

圖15 孔徑對沉積量的影響Fig.15 Effect of pore size on deposition amount

3.4 喉道長度的影響

圖17所示為管道內最大沉積量隨喉道長度的變化規律。入口壓力為5 MPa,隨著喉道長度增加,管內最大沉積量逐漸上升,但增加幅度逐漸減小。由10 mm增加至20 mm時,最大沉積量增加了20.83%,湍動能增加15.02%;當增加至30 mm時,沉積率增加了13.99%,湍動能增加3.38%;當增加至40 mm時,最大沉積量增加了6.43%,湍動能增加0.038%。

圖17 喉道長度對沉積量的影響Fig.17 Effect of throat length on deposition amount

圖18為沿流動方向中軸線處湍動能變化,隨著喉道長度的增加,喉道出口截面處湍動能略有增大,增幅逐漸減小,而距喉道出口0.17 m處左右,湍動能逐漸減小,并且峰值出現后移。這是因為隨著喉道長度的增加,氣流經過節流后更容易達到穩定的狀態,導致湍動能變化不大。這表明此時流動已經逐漸穩定,更接近于圓管內湍動能分布。由于湍動能變化幅度減小,最終導致最大沉積量增長率也隨之減小。

圖18 沿流動方向湍動能變化Fig.18 Variation of turbulent kinetic energy along the flow direction

4 結論

(1)傳統孔板式節流閥的蠟沉積主要出現在2個位置:節流閥節流出口(x=0.77～0.775 m)處壁面處以及距節流閥喉道出口0.17 m左右的管壁處,并且管道內顆粒最大沉積量隨壓力的增加而增加,但增加趨勢逐漸變緩。沉積寬度為0.01 m,當壓力增加1 MPa時,第二處沉積位置后移約0.005 m。

(2)在節流直徑、喉道長度一定時,管道內最大沉積量隨著粒徑的增大而增大,并且當粒徑大于50 μm后,最大沉積量明顯增加。

(3)管道內壁最大沉積量與節流孔徑呈負相關。其中,當節流孔徑增加到50 mm時,最大沉積量幾乎為0;最大沉積量與喉道長度呈正相關,但隨著喉道長度的增加,管內最大沉積量的增長趨勢逐漸變緩。