基于沖蝕-動網格耦合的繞絲篩管沖蝕過程數值模擬

孫巖 樓一珊 曹硯峰 文敏 翟曉鵬 趙軒康

1.長江大學石油工程學院；2.中海油研究總院有限責任公司

繞絲篩管是一種重要的機械防砂工具。繞絲篩管在防砂過程中，擋砂介質受地層砂連續不均勻的沖擊，壁面從點痕出現到沖蝕痕跡深度增加，并最終造成繞絲篩管破損的過程被稱為繞絲篩管的沖蝕進化(Evolution of the surface profile)［1］。篩管的沖蝕破壞將使防砂失效、油氣井大量出砂，導致油氣井減產、停產。因此，研究防砂篩管的沖蝕現象將對油氣井防砂及篩管設計具有重要指導意義。

為研究防砂篩管的沖蝕磨損現象，Mckeon等［2］通過實驗研究了繞絲與基管間距大小對繞絲篩管沖蝕磨損程度的影響；Gilledpie等［3］對不同類型篩管進行沖蝕分析并建立了壽命預測模型；劉永紅等［4］通過實驗討論了流速、砂粒直徑、砂粒濃度及沖蝕角對割縫篩管沖蝕磨損率的影響；劉新峰等［5］通過研究各因素不同類型篩管的沖蝕影響程度對如何延長篩管使用壽命提出建議；李效波等［6］通過SEM、EDS對篩管沖蝕腐蝕的形貌進行了宏微觀分析，根據實驗結果對迭加機理進行驗證；廖華林等［7］通過實驗建立了考慮堵塞問題的礫石充填篩管防砂裝置并分析篩管失效原因及堵塞的內在機制。部分學者通過數值模擬手段進一步研究了繞絲篩管的沖蝕破壞問題。Deng等［8］通過流場計算及侵蝕模型的引入研究不同生產壓差與含砂量對篩網的侵蝕情況；鄧自強［9］采用Fluent軟件對割縫篩管進行了液固兩相流數值模擬，探究篩管割縫處最大沖蝕率與流速、砂粒濃度和含砂水平之間的關系。目前，數值模擬方法主要用于計算平均沖蝕率以研究沖蝕破壞現象，然而，篩網的破壞主要是篩網受砂礫作用而變形所致。目前研究未考慮到篩網在沖蝕過程中的變形問題。

建立了篩管井下工況的固-液兩相流模型，并結合沖蝕動網格耦合方法，在考慮篩管受固體顆粒沖蝕變形、篩管內部過流流量變化的情況下研究繞絲篩管沖蝕進化過程。

1 計算理論

為解決繞絲篩管沖蝕破壞問題，需要解決以下幾方面問題，得到不同參數對繞絲篩管篩網變化的影響規律：(1)模型基本理論；(2)有限元模型；(3)模型計算條件和邊界條件；(4)模型計算設定參數。

1.1 沖蝕-動網格耦合理論

在實際工況下，設備材料因受沖蝕作用導致質量損失引起幾何變形，觀測到沖蝕率隨時間變化，導致固體顆粒撞擊速度發生變化，從而加速或減緩沖蝕的過程，這將影響顆粒流動和沖蝕過程并最終導致材料失效。沖蝕率本質上是瞬態的，因為在繞絲篩管發生質量損失后，流動邊界條件發生了變化，然而目前的計算模型僅提供一個初始侵蝕速率，這就使得試驗結果有顯著的誤差。因此筆者利用ANSYS Fluent中基于準穩態方法運行的沖蝕—動網格模型來研究繞絲篩管沖蝕進化的過程。該模型在沖蝕過程的每一段時間間隔捕獲因沖蝕作用導致的材料表面幾何變形的沖蝕質量損失。并通過該方法重新劃分網格,網格位置由動態網格子模型使用物理時間步長更新，并計算該區域的局部流動條件，預測每一時間段的沖蝕率，實現了流-固-沖蝕相互作用的完全耦合［10］。

沖蝕動網格耦合，單個面的網格變形計算如式(1)［10］

式中，ERf為單位時間內單位面積的侵蝕率，kg · (m2· s)-1；Δt為網格更新的時間步長，s；ρwm為每單位體積質量的設備材料密度，kg · m-3；Δxf為單個面網格變形量，m。

1.2 沖蝕模型

對篩管壁面沖蝕速率采用McLaury模型來計算，該模型主要被用來預測含砂流體對材料壁面的沖蝕磨損［11-13］，形式如式(2)

式中，E為表面磨損率，kg/(m2· s)；V為顆粒速度，m/s；A、k為經驗系數，A=F · HBk；HB為材料布氏硬度，N/mm2；F為經驗常數；γ為顆粒與壁面沖擊角度，°。

沖擊角函數可用分段多項式表達

式中，γlim為過度角，15°；b,c,x,y,z,w為經驗系數。

1.3 篩管模型

繞絲篩管由基管(帶孔中心管)、縱筋和不銹鋼繞絲組成。不銹鋼繞絲保持一定縫隙纏繞在縱筋上，形成繞絲縫隙。繞絲縫隙在繞絲篩管中起到擋砂作用，繞絲縫隙的寬度可以決定繞絲篩管的擋砂精度。含砂流體穿過繞絲縫隙，通過基管上具有一定密度和孔徑的圓孔(基管孔)流入井筒內(圖1)。

圖1 繞絲篩管模型Fig.1 Model of wire wrapped screen

以Mckeon等［2］試驗條件為基礎，(表1)建立了繞絲篩管CFD流體域幾何模型，如圖2所示。設置出口、入口為邊界條件，模型參數：繞絲間縫寬0.371 mm；繞絲尺寸2.286×2.667 mm；基孔直徑9.525 mm；縱筋尺寸2.286×2.667 mm；縱筋間距9.525 mm；入口截面積0.004 558 1 m2。

在計算域用ANSYS Meshing進行網格劃分，采用對復雜幾何結構有很強適應性的非結構性網格，對幾何模型進行網格獨立性驗證，不同網格節點數下最大速度大小如圖3所示，由圖3可知當網格節點數達到350 000時，最大速度基本不變，最終網格劃分結果如圖4所示；劃分了節點數為389 632，單元數為1 993 663個的非結構化網格，網格質量SkewnessMax=0.81，網格質量較好，滿足計算要求。

1.4 邊界條件及數值計算方法

湍流模型采用Realizablek-ε模型，離散相模型入口面為射流源，入口邊界條件設置為速度入口，出口邊界條件采用outflow，壁面為無滑移壁面，采用連續相模型和離散相模型(Simple算法)［14］求解。動量、湍流能與湍流耗散率均采用二階迎風差分格式離散。沖蝕-動網格耦合設置平滑步為5，最大節點移動比例設置為20%，采用可變時間步法設置初始時間步長為10，每個流動迭代次數設置為700步。根據Mckeon等［2］試驗給出的試驗條件及數據設置入口速度為0.664 m/s、粒徑為0.24 mm、砂粒質量流量為0.016 3 kg/s、砂粒密度為2 300 kg/m3、含砂濃度為4.5×103m3、入口流量為3.359 6×10-3m3/s、水力直徑為0.071 m，流體相為液態水。

表1 篩管單位時間內的比沖蝕率與沖蝕質量損失計算結果Table 1 Calculation results of specific erosion rate and erosion mass loss of screen per unit time

圖2 流體域幾何模型Fig.2 Geometric model of fluid domain

圖3 壁面剪切應力隨節點數變化曲線Fig.3 Variation of wall shear stress with the number of nodes

圖4 網格劃分結果Fig.4 Grid division result

2 CFD模型驗證

通過CFD后處理模塊得到繞絲篩管沖蝕率分布云圖，可看出沖蝕磨損嚴重區域主要集中在基管孔周圍，以及繞絲交于縱筋且靠近基管孔一側區域，將數值模擬結果與實驗結果進行比較(圖5)可知數值模擬結果較為準確地預測了繞絲篩管發生沖蝕磨損的區域。通過函數計算計算出篩管的比沖蝕率，結果表明單位時間內沖蝕造成的質量損失整體趨勢下降，與Mckeon等［2］實驗結果對比，整體趨勢基本一致(圖6)。

圖5 繞絲及縱筋實驗結果與數值模擬結果對比Fig.5 Comparison between experimental results of numerical simulation results of wrapping wire and vertical member

圖6 計算結果與實驗結果對比Fig.6 Comparison between calculation results and experimental results

數值計算結果與實驗計算的比沖蝕率結果較為相近(表1)，經計算平均誤差為1.976%,沖蝕質量損失的結果平均誤差為2.73%，單位時間內繞絲篩管的沖蝕質量損失均在0.07 g左右，計算值與實驗值接近，進一步佐證了該模型的可靠性。因未考慮到沖蝕進化現象而導致計算結果誤差較大的問題在引用沖蝕-動網格耦合計算后得到了有效的解決。

數值計算與實驗結果誤差原因：數值模擬時未考慮砂粒間相互作用，計算結果較大；由于實際砂粒間存在相互作用致使砂粒能量消耗，降低了砂粒動能，速度隨之衰減，因此實際比沖蝕率較低。

3 繞絲篩管動網格沖蝕進化過程分析

3.1 流場隨沖蝕時間變化規律

如圖7所示，流體流經繞絲篩管后壓力逐漸降低；隨時間延長進口壓力不變，出口壓力增大。經CFD后處理函數計算壓差可知隨沖蝕時間延長，攜砂流體流經篩管后壓差逐漸降低，主要原因是繞絲篩管逐漸發生沖蝕破壞，節流效應下降，壓力向下擴展。

圖7 不同沖蝕時間下流場截面壓力分布云圖Fig.7 Cloud chart of pressure distribution on flow field section at different erosion time

過流流速最高的位置出現在基管孔內及繞絲縱筋靠近基管孔的兩側(圖8)，主要原因是攜砂流體流經繞絲縫隙，受節流效應，使過流流速增加，通過基管孔再次節流增速噴射入井筒內。隨時間的延長，繞絲篩管壁面區域流場流速變化不大。經CFD后處理函數計算流速，由計算結果可知，區域流速隨時間延長整體呈下降趨勢。這是由于繞絲篩管受沖蝕破壞，節流效應下降，區域整體流速逐漸降低。

圖8 不同沖蝕時間下流場截面流速分布云圖Fig.8 Cloud chart of flow velocity distribution on flow field section at different erosion time

3.2 繞絲篩管沖蝕隨時間變化規律

繞絲篩管沖蝕變形量如圖9云圖所示。

圖9 繞絲篩管沖蝕變形量云圖(a: 繞絲、縱筋壁面；b: 基管入口端面)Fig.9 Cloud chart of erosion deformation of wire wrapped screen (a: wall surface of wrapping wire and vertical member; b: end face of base pipe inlet)

圖9云圖中紅色區域為發生明顯沖蝕變形的區域，主要集中在繞絲下游一側、繞絲正上方與基管孔周圍，更直觀地展示了繞絲篩管沖蝕破壞開始于何處。隨時間延長，攜砂流體不斷的沖刷篩管，導致繞絲篩管沖蝕變形量不斷增加，沖蝕變形區域逐漸擴大。

攜砂流體以0.664 m/s的初始流速流入篩管，繞絲篩管沖蝕進化過程中單位時間內的累計沖蝕變形量及過流流速變化趨勢如圖10所示。篩管內過流流速平均增量為7.83 m/s，隨時間延長，繞絲篩管沖蝕變形量增加，繞絲篩管內過流流速呈遞減趨勢下降1%左右。這可能是由于繞絲篩管在攜砂流體集中流動區域被沖蝕破壞，導致攜砂流體流動面積擴大，比沖蝕率降低，從而導致過流流速下降。

圖10 繞絲篩管沖蝕變形量與流速趨勢圖Fig.10 Trend map of erosion deformation and flow velocity

由圖10可看出，繞絲累計沖蝕變形量最多，隨時間延長累計沖蝕變形量由3.40×10-9m增長至1.30×10-8m，縱筋累計沖蝕變形量由3.21×10-9m增長至1.30×10-8m，而基管孔累計沖蝕變形量相對較少。這是由于繞絲與縱筋相交處形成暫時的封閉空間，攜砂流體首先撞擊在繞絲壁面，部分砂粒發生一次撞擊彈開后再次或多次撞擊在繞絲壁面直至從繞絲與縱筋縫隙中流入井筒(圖11)，導致其沖蝕破壞較嚴重。

通過CFD后處理函數計算繞絲篩管撞擊砂粒質量流量，計算結果證明繞絲受攜砂流體撞擊頻次最多，單位時間內撞擊砂粒質量流量在0.055～0.040 3 kg/s。其次是縱筋，但縱筋單位時間內砂粒撞擊質量流量變化較繞絲較大，由繞絲篩管結構可知(圖1)，縱筋之間間隔比繞絲間間隔較大，隨時間增加，繞絲篩管受沖蝕破壞，導致攜砂流體流動區域發生變化，從而使縱筋受撞擊頻次變化較大。因此，在篩管設計時可結合實際工況模擬分析，設計合理的繞絲、縱筋、基管孔間隔尺寸來延長篩管壽命。

圖11 砂粒運動軌跡流速云圖Fig.11 Cloud chart of sand trajectory and velocity

采用數值模擬手段考慮了繞絲篩管在沖蝕過程中的變形問題。以上闡述了繞絲累計沖蝕變形量的變化規律，隨時間延長可以看出，繞絲及縱筋的幾何形狀因受沖蝕磨損作用而發生明顯變化(圖12)，縫隙逐漸變寬。變寬的縫隙在砂粒流動的過程中起到優勢通道的作用，節流作用隨流動面積增大而減小，導致過流流速下降，如圖13所示。

圖12 不同時間繞絲縱筋受沖蝕后幾何形狀變化云圖Fig.12 Cloud chart of geometric erosion deformation of wrapping wire and vertical member at different time

繞絲篩管累計沖蝕變形量與沖蝕總質量損失成正比，且均與比沖蝕率成反比，這也證明了繞絲篩管受沖蝕磨損破壞導致流體流動區域增大，擋砂精度變差，從而使得比沖蝕率下降，最終會導致防砂失效。

圖13 繞絲篩管比沖蝕率與累計沖蝕變形量及沖蝕總質量損失趨勢Fig.13 Trend map of specific erosion rate, cumulative erosion deformation and erosion induced total mass loss of wire wrapped screen

4 各因素對繞絲篩管的沖蝕影響

影響繞絲篩管沖蝕率的主要影響因素包括流體流速、含砂質量分數、砂粒粒徑、篩管材料性、砂粒硬度等。本文數值模擬研究考慮的主要因素為流速、含砂質量分數、砂粒粒徑。Abrams［15］在通過對地層顆粒運移的研究提出1/3橋堵原理，即固相顆粒在地層孔喉處形成“砂橋”的條件是顆粒尺寸為地層孔喉尺寸的1/3［16］。Gillespi等［17］通過實驗也得出砂粒直徑小于篩網擋砂精度3倍以下時，不存在砂粒堵塞現象的結論。因此，本文在各因素對繞絲篩管沖蝕影響分析的數值模擬研究中不考慮砂粒堵塞情況，采用直徑不大于123 μm的砂粒，僅探究細粉砂對繞絲篩管的沖蝕影響。

4.1 流速

沖擊角90°、砂粒直徑120 μm、含砂質量分數0.5%，入口流速分別取0.664、0.833 2、1.062 4、1.261 6 m/s進行數值模擬計算，繞絲篩管沖蝕率隨流速變化曲線及沖蝕分布云圖如圖14、圖15所示。繞絲篩管沖蝕率與攜砂流體的流速呈指數關系，這與API中標準的金屬材料過流磨損情況相同。隨著含砂流體流速增加，砂粒所具有的動能增大，對篩管的沖蝕破壞能力增強。同時，攜砂流體流速增加，質量流量隨之增加導致單位時間內壁面遭受更多的砂粒撞擊，沖蝕率增加，沖蝕質量損失增大。

圖14 篩網沖蝕率隨流速變化曲線Fig.14 Variation of erosion rate with flow velocity

圖15 不同流速下篩網沖蝕分布云圖(a: 繞絲、縱筋壁面；b: 基管出口端面)Fig.15 Cloud chart of erosion distribution of screen mesh at different flow velocities(a: wall surface of wrapping wire and vertical member; b: end face of base pipe inlet)

4.2 含砂質量分數

沖擊角90°、砂粒直徑120 μm、流速0.664 m/s，含砂質量分數分別取0.3%、0.5%、0.7%、0.9%進行數值模擬計算，繞絲篩管沖蝕率隨含砂質量分數變化曲線及沖蝕分布云圖如圖16、圖17所示。

圖16 篩網沖蝕率隨含砂質量分數變化曲線Fig.16 Variation of erosion rate with sand mass fraction

圖17 不同含砂質量分數下篩網沖蝕分布云圖(a: 繞絲、縱筋壁面；b: 基管出口端面)Fig.17 Cloud chart of erosion distribution of screen mesh at different sand mass fractions (a: wall surface of wrapping wire and vertical member;b: end face of base pipe inlet)

繞絲篩管沖蝕率與含砂質量分數呈線性關系，繞絲篩管壁面沖蝕破壞區域擴大。原因是增加含砂質量分數即增加砂粒濃度，就代表增加了砂粒與繞絲篩管入口碰撞的次數，單位時間內砂粒撞擊篩管壁面的頻率增加，導致篩管壁面沖蝕率增加，沖蝕質量損失增加。

4.3 砂粒直徑

沖擊角90°、含砂質量分數0.5%、流速0.664 m/s，砂粒直徑分別取120、100、80、60 μm進行數值模擬計算，繞絲篩管沖蝕率隨砂粒直徑變化曲線及沖蝕分布云圖如圖18、圖19所示。繞絲篩管沖蝕率與砂粒直徑基本呈線性關系，沖蝕破壞區域擴大。原因是攜砂流體流速一定，其使砂粒具有的最大動能一定，隨砂粒直徑增加導致慣性增大，砂粒對繞絲篩管壁面的切削作用增加，從而加劇了繞絲篩管壁面的沖蝕磨損。這一結果符合Misra和Finnie［18］的研究結論之一：材料沖蝕量隨砂粒直徑的增大而增加。

圖18 篩網沖蝕率隨砂粒直徑變化曲線Fig.18 Variation of erosion rate with sand particle diameter

圖19 不同砂粒直徑條件下篩網沖蝕分布云圖(a: 繞絲、縱筋壁面；b: 基管出口端面)Fig.19 Cloud chart of erosion distribution of screen mesh at different sand particle diameters(a: wall surface of wrapping wire and vertical member; b: end face of base pipe inlet)

5 未堵塞工況下沖蝕模型

5.1 沖蝕模型的建立

考慮了流速、含砂質量分數、砂粒直徑3個因素，同時考慮篩管材料硬度等，提出一個未堵塞工況下的沖蝕模型

式中，Ev為沖蝕速率，kg/(m2· s)；Vg為篩網過流流速，m/s；Dm為砂粒直徑，m；Cp為含砂弄濃度，kg/m3；m為攜砂流體質量流量，kg/s；A為篩網面積，m2，本數模面積為0.004 558 1 m2；Δ為沖蝕時間，s；K、a、b、c、d為各參數系數，通過對數模結果進行擬合獲得。

5.2 沖蝕模型的擬合

根據不同條件下的數值模擬數據(表2)，運用1stOpt軟件對沖蝕模型進行非線性擬合(表3～表4)，優化算法采用麥夸特法(Levenberg-Marquardt)+通用全局優化法。

表2 各方案數值模擬條件Table 2 Numerical simulation conditions of each plan

表3 擬合結果相關性Table 3 Correlation of fitting results

表4 各系數擬合值Table 4 Fitting value of each coefficient

擬合結果與數模結果對比如圖20所示。最終獲得了未堵塞工況下沖蝕模型

圖20 擬合結果與數值模擬數據對比Fig.20 Comparison between fitting results and numerical simulation data

基于上述沖蝕模型，結合Gillespi等［17］通過實驗得出當繞絲質量損失約為其質量的2%時，可能會導致防砂失效的結論，在繞絲篩管的使用壽命預測時，將這一條件作為繞絲篩管防砂失效的一臨界條件，更準確的預測繞絲篩管使用壽命。

6 結論及展望

(1)研究表明，繞絲篩管沖蝕變形區域主要集中在繞絲下游一側、繞絲正上方與基管孔周圍，直觀的表示出繞絲篩管沖蝕破壞開始于何處。隨時間增加，由于繞絲篩管在攜砂流體集中流動區域不斷受到沖蝕破壞，導致攜砂流體流動面積擴大，比沖蝕率降低。這對使用沖蝕測試數據來推斷井的總壽命或防砂失效時間具有重要參考意義。

(2)數值模擬結果表明，采用沖蝕-動網格耦合模型，在篩網發生沖蝕變形后再次進行網格劃分更好的模擬了實際篩網的沖蝕進化現象，較準確的預測繞絲篩管受攜砂流體沖蝕篩管表面所造成的的幾何變形及失效破壞現象，因此使用此數值手段模擬結果對繞絲篩管的壽命預測及結構設計有著重要參考意義。

(3)在不考慮砂粒堵塞工況(粒徑小于123 μm)下，繞絲篩管沖蝕率隨流體流速增加呈指數增長，隨含砂質量分數及砂粒直徑的增加呈線性增長，沖蝕磨損區域擴展現象顯著。

(4)由于實際工況較為復雜，涉及因素較多，難以用一個公式或多個公式來精確的描述，因此試驗和數值模擬都存在一定的誤差，但采用沖蝕-動網格耦合技術進行的數模研究結果誤差較小，通過對其數模結果擬合，建立沖蝕預測模型來計算沖蝕率在一定意義上是可取的。