基于動網格與UDF技術的內桿做行星運動的環空中流場數值模擬

龐博學，楊樹人，劉麗麗，劉 超，張 瑩，馬立輝

(東北石油大學，黑龍江 大慶 163000)

基于動網格與UDF技術的內桿做行星運動的環空中流場數值模擬

龐博學，楊樹人，劉麗麗，劉 超，張 瑩，馬立輝

(東北石油大學，黑龍江 大慶 163000)

應用動網格技術以及UDF程序成功實現了內桿在油井偏心環空中自轉同時公轉的行星運動狀態。完成了內桿在井筒偏心環空中不同偏心率位置轉動的情況下，流體在井筒環空中流動的數值模擬。計算結果表明，內桿偏心率增加，流體周向運動增強，湍流高黏度區域擴大，層流區域相應增加。存在一臨界偏心率εr，在其兩側，流體作用力分別將內桿離心外推和向心內拉。且內桿偏心率越大，流體內推力就越大，流體對內桿公轉運動方向上的阻力也越大。

動網格；UDF；行星運動；Fluent；偏心率

對于內桿靜止的偏心環空中的流動，前人做了大量的研究工作。與內桿靜止的偏心環空流動相比較，當內桿自轉與公轉同時進行時，流體在此種環空中的流動要更加復雜，人們對這種流動的研究少之又少[1]。在石油鉆井工程中，尤其是對于水平井和定向井，鉆桿因受自身重力的影響，在井筒中是彎曲的，再加上電動機驅動鉆桿的自轉，在井筒與鉆桿所形成的偏心環空中的鉆井液的流動就屬于內桿做行星運動的環空中的流動問題[2]。與此相似，在采油工程尤其是螺桿泵采油中，由于內桿偏心，抽油桿在自轉的同時還沿外管軸線公轉，在該環空中采出液的流動也屬于內桿做行星運動的環空中的流動問題[3]。由此可見，此種流動在工程實際中常常遇到。

本文旨在應用計算流體動力學的方法，通過編譯UDF程序實現內管行星運動狀態，應用Fluent動網格技術對內桿行星運動的環空中的流動問題進行數值求解。此舉可以為定向井和水平井等工況鉆井液的水力參數設計，以及螺桿泵采油生產工況的優化提供理論指導方向，具有學術價值和工程實踐意義[4]。

1 動網格技術在內桿行星運動中的應用

1.1 模型簡化

內桿行星運動環空流場簡化模型如圖1所示。將石油工程實際中所描述的問題做如下簡化：流體處于內、外2個偏心的圓管構成的環空中，外圓管靜止，內圓管以恒定的角速度ω自轉，又以恒定的角速度Ω沿外圓管軸線公轉，并有沿外圓管軸線方向的壓力P作用在流體上。

圖1 內桿行星運動環空流場簡化模型

1.2 動網格方法

為實現環空中內桿的行星運動，數值模擬應用動網格技術，如圖2所示。內桿運動使得流體幾何區域變化劇烈，動網格更新方法采用局部網格重構模型(remeshing)，動區域選擇指定剛體運動(rigid body)。模擬計算過程中，因動邊界導致局部網格超過定義的最大尺寸時網格一分為二，局部網格受擠壓緊縮低于定義的最小尺寸時合二為一，當網格扭曲程度超過定義的傾斜度時亦重新劃分，以適應運動區域的形態變化。被重新劃分的網格單元緊鄰動網格節點，運動邊界附近網格得到及時更新[5]。

圖2 內桿自轉及內桿自轉且公轉動網格瞬態圖

動網格任意控制體V，通用標量Φ的積分守恒方程為：

時間導數項可由一階向后差分項用如下公式表示：

式中，n和n+1是當前時間及下一層時間。

第n+1時間層上的體積Vn+1為：

式中，dV/dt是控制體的時間導數。

為滿足網格守恒律，控制體的體積時間導數為：

每個控制容積面上的點積為：

式中，δVj是整個時間步Δt上控制容積面j膨脹引起的體積改變。

2 UDF的編譯及求解模型的選擇

2.1 內桿行星運動的用戶自定義函數(UDF)的編譯

UDF是為滿足計算個性化需求，可以動態地連接到Fluent求解器上的用戶自編的程序。通過UDF調用動網格宏(DEFINE_CG_MOTION)可實現剛體運動邊界的控制。內桿運動UDF編譯的主體思想如模型簡化圖。自轉通過DEFINE_CG_MOTION宏中的omega[z]實現，公轉通過控制內桿旋轉速度矢量vel[x]、vel[y]以及動網格矢徑完成。

2.2 數值求解條件及模型選擇

定義流體密度ρ為0.8×103g/m3，流變指數n為0.6，外管內徑R為40 mm，內桿外徑r為20 mm，內桿自轉速度ω為30 rad/s，公轉速度Ω與ω相同，為30 rad/s，環空流體軸向壓力梯度為300 Pa/m。

此次內桿行星運動的環空流場數值計算對控制方程的離散應用有限體積法，采用分離式求解器求解。選擇k-epsilion 湍流模型，并開啟能量方程，用中心差分格式離散擴散項，選擇PRESTO!格式對壓力項離散，其余用一階迎風離散，壓力速度耦合方程采用PISO算法，使用編譯型UDF文件定義剛體速度，動網格采用局部網格重構方法，對內桿行星運動的環空流場進行非定常求解[6]。

3 偏心環空螺旋流與內桿行星運動的環空中的流動

偏心環空螺旋流是僅有內桿自轉的環空中的流動，而內桿行星運動的環空中的流動則是在偏心環空螺旋流的基礎上，增加內桿以角速度Ω繞外管軸線公轉的情況。分別模擬計算偏心環空螺旋流與內桿行星運動的環空中的流動有助于分析內桿公轉對流場的影響。

3.1 偏心環空螺旋流與內桿行星運動環空中流體的速度場計算

通過流場周向速度分布圖(見圖3)可以得出結論：在僅有內桿自轉的偏心環空螺旋流中，由于黏性力的作用，內桿的自轉促使其壁面附近區域流體產生周向速度，速度標量值以內桿軸線為圓心呈環形區域分布，且由內桿壁面沿徑向逐漸減小。而在內桿行星運動的環空中，由于內桿在環空中的偏心公轉，被轉動波及的流體區域更大，流體周向速度最大值保持在靠近內管壁環空窄間距側。

圖3 流場周向速度分布圖

3.2 偏心環空螺旋流與內桿行星運動環空中流體的壓力場計算

通過流體動態壓力分布圖(見圖4)可以得出結論：在偏心環空螺旋流中，流體作用力峰值出現在內桿壁表面，且內桿周向受力基本均勻；而在內桿行星運動的環空中，流體作用力峰值保持在內桿公轉迎流面窄間距處，即流體有阻止內桿公轉并將內桿向環空內側推回的趨勢。

圖4 流體動態壓力分布圖

4 內桿偏心程度對流場的影響

在內桿自轉且公轉的環空流動中，流體與內桿存在相互作用力，內桿在流體中轉動的同時伴有橫向擾動， 這種擾動使得內桿公轉半徑即相對于外管的偏心距在時刻變化，流體對內桿的徑向作用力與內桿偏心程度相互作用。由于內桿偏磨、失穩等情況與偏心率密切相關，因此有必要針對不同偏心率對內桿行星運動流場的影響計算分析。

定義內桿偏心率ε=e/emax，其中e是內桿相對于外管軸線的偏心距，emax是最大偏心距(即內桿外壁接觸外管內壁的完全偏心)，此模型中，emax=20 mm。分別模擬計算內桿偏心率ε=0.25、ε=0.5以及ε=0.75時內桿行星運動的環空流體場。

4.1 不同偏心率環空流體周向速度場分布

通過不同偏心率環空流場周向速度分布圖(見圖5)可以得出結論：由于內桿在環空中的偏心公轉，隨著偏心距的增加，被轉動波及的流體區域擴大，流體周向運動增強，速度標量值增大。

圖5 流場周向速度分布圖

4.2 不同偏心率環空流體湍流強度分布

通過不同偏心率環空流場湍流強度分布圖(見圖6)可以得出結論：內桿偏心距對流場湍流強度影響明顯，湍流強度峰值集中出現在環空寬間距側，且隨著偏心距的增加，湍流強度減弱，層流區域增加。

圖6 流場湍流強度分布圖

4.3 不同偏心率環空流體高黏度區域分布

通過不同偏心率環空流場湍流黏度分布圖(見圖7)可以得出結論：由于內桿在環空中旋轉導致漩渦，渦流擴散使得流體產生遠大于分子水平的表觀黏度，即湍流黏度，使組分黏度增加，湍流黏度峰值區域集中出現在寬間距側，但繞內桿轉動方向出現旋轉偏移。且偏心距越大，湍流高黏度區域越大。

圖7 流場湍流黏度區域分布圖

4.4 偏心率對內桿在環空中所受流體力的影響

內管在環空流體中轉動的過程中，流體作用力將影響內管受力平衡，此作用力由流體壓強和周向剪應力綜合作用產生[7]。將流體作用力F分解為沿公轉路線徑向上的流體力Fj和沿公轉路線切向上的流體力Fq。定義Fj正方向指向公轉軌道軸心，Fq正方向與內桿公轉運動速度矢量方向一致。分別模擬內桿偏心率ε為0.05、0.20、0.35、0.50、0.65、0.80及0.95情況下行星運動時的流體場，計算得出內管所受流體沿公轉路線徑向的作用力Fj和沿公轉路線切向的作用力Fq，繪制如圖8所示曲線。

圖8 內桿在不同偏心率ε下所受流體作用力曲線

通過不同偏心率環空流場內桿所受流體作用力曲線可以得出結論：內桿在公轉軌道徑向上(Fj方向)所受流體作用力一般情況下不平衡。存在臨界偏心率εr(約為0.65)，使得當ε<εr時，Fj始終小于零，流體作用力有將內桿向離心方向外推的趨勢；ε>εr時，Fj始終大于零，流體作用力有將內桿向公轉軌道軸心方向內推的趨勢，且偏心率ε越大，內推作用力越大。而在公轉軌道切向上(Fq方向)，內桿所受到的流體作用力始終小于零，說明流體力一直在阻止內桿沿公轉軌道的運動，內桿偏心率越大，流體阻力越大。

5 結語

應用Fluent動網格技術和UDF可以成功實現偏心環空中內桿的行星運動。對內桿在不同偏心率下行星運動的環空中的流動研究得出如下結論。

1)隨著內桿偏心率ε的增加，被轉動波及的流體區域擴大，流體周向運動增強。

2)內桿旋轉導致渦流擴散使得流體產生遠大于分子水平的湍流黏度，其峰值區域在環空寬間距側，但繞內桿轉動方向出現旋轉偏移。隨著偏心距的增加，湍流高黏度區域擴大，湍流強度減弱，層流區域增加。

3)存在臨界偏心率εr：當ε<εr時，流體作用力將內桿向離心方向外推；當ε>εr時，流體作用力將內桿向軸心內推，偏心率越大，內推作用力越大，流體對內桿公轉運動方向上的阻力越大。

[1] 崔海清,張淑云,修德艷,等. 流體在內管做行星運動的環空中流動的壓力梯度[J]. 中國石油大學學報：自然科學版,2008(6):76-78，87.

[2] 裴曉含. 冪律流體在內管做行星運動的環空中流動時內管壁的受力分析[D].大慶：大慶石油學院,2007.

[3] 季海軍. 冪律流體在內管做行星運動的環空中的流動[D].大慶：大慶石油學院,2005.

[4] 楊樹人. 粘彈性流體偏心環空非定常流的數值計算[D].大慶：大慶石油大學,2006.

[5] 張來平,鄧小剛,張涵信. 動網格生成技術及非定常計算方法進展綜述[J]. 力學進展,2010(4):424-447.

[6] 楊樹人,龐博學,劉麗麗. 天然氣井應用渦流工具排水采氣的流場分析[J]. 石油礦場機械,2014(10):13-16.

[7] 馮順新,李啟兵,符松. 內管自轉且公轉時環管冪律流的受力分析[J]. 計算力學學報,2008(2):212-217.

責任編輯彭光宇

NumericalSimulationofFlowFieldinEccentricAnnulusWhoseInnerTubehasPlanetaryMotionbasedonDynamicMeshandUDF

PANG Boxue，YANG Shuren，LIU Lili，LIU Chao，ZHANG Ying，MA Lihui

(Northeast Petroleum University, Daqing 163000，China)

Movement of the inner tube rotation and revolution in eccentric annulus was achieved successfully by means of dynamic mesh and UDF. The flow in the annulus whose inner tube at different eccentricity rotates and revolves simultaneously was simulated. The results showed that: with the eccentricity increase, circumferential movement of fluid enhance, aera of high turbulent viscosity and laminar region expanded. There was a critical eccentricity on each side of which fluid force pushed inner tube outward and pulled it inward. The larger the eccentricity was, the greater fluid force that pulled inner tube inward was. Furthermore, resistance to the motion of the inner tube in the direction of revolution by fluid becomed larger.

dynamic mesh, UDF, planetary motion, Fluent; eccentricity

TH 311

：A

龐博學(1991-)，男，碩士研究生，主要從事復雜流體流動與數值模擬等方面的研究。

2014-12-28