管道內壁粗糙度對沿程阻力影響的FLUENT數值模擬分析

姚雪蕾，袁成清，付宜風，白秀琴

(武漢理工大學a.能源與動力工程學院可靠性工程研究所；b.船舶動力工程技術交通行業重點實驗室,武漢 430063)

管道內壁粗糙度對沿程阻力影響的FLUENT數值模擬分析

姚雪蕾，袁成清，付宜風，白秀琴

(武漢理工大學a.能源與動力工程學院可靠性工程研究所；b.船舶動力工程技術交通行業重點實驗室,武漢 430063)

針對在海底管道的設計階段準確獲取沿程阻力的問題，采用FLUENT軟件運用數值模擬的方法，分析粘性流體在管道中的摩阻損失與管道內壁粗糙度的關系，將計算結果與理論結果進行對比。結果顯示，除流態從水力光滑區到混合摩擦區間有一個過渡階段存在誤差外，其余都吻合良好，證實輸油管道的內壁粗糙度只在混合摩擦區對摩阻有較大影響，改善管壁面粗糙度可以減小阻力。若流態處于水力光滑區到混合摩擦區之間，計算摩阻系數僅靠經驗公式是不合適的，應考慮采用數值模擬計算的方法。

輸油管道；粗糙度；沿程阻力；FLUENT；數值模擬

隨著計算流體動力學的發展，利用FLUENT分析粘性流體在管道中的阻力特性已經得到了很多實用性的成果，如針對粗糙的管壁[1-4]或者是產生局部阻力損失的典型閥件(突擴管、三通管、彎管等)[5-6]，通過數值模擬的方法可以模擬出常規實驗無法得到的不同流場的各種信息。然而，目前工程上在管道的設計階段，沿程阻力的確定一般只采用經驗公式數學計算的方法，故考慮運用數值模擬的方法對其準確性進行驗證，為管道設計中沿程阻力的獲取提供依據。

1 沿程阻力與管道內壁粗糙度

管道的沿程阻力損失與管長、管流的速度分布，以及管壁的粗糙度等因素有關，通常用hf表示，按達西(Darcy-Weisbach)公式計算。

(1)

式中:L——管道長度;D——管道內徑;v——流體流速;λ——摩阻系數。

相對當量粗糙度ε=2e/D。式中e為管壁的絕對當量粗糙度，是指管內壁凸起高度的統計平均值。研究表明，新的、清潔的管壁的絕對當量粗糙度僅取決于管材及制管方法，與管徑無關；使用后的管路則隨運行情況，如所輸流體性質、腐蝕程度、運行年限、清管方法等不同而有顯著變化[7]。

根據莫迪理論，一般可以把圖分為5個區域：層流區、層流轉變為紊流的過渡區、紊流水力光滑區、紊流混合摩擦區、紊流粗糙區[8-9]。在實際運行中，熱原油管道上最常見的流態是水力光滑區，輕油管道也多在水力光滑區，輸送低粘油品的較小直徑管道可能進入混合摩擦區。所以本研究重點在紊流水力光滑區和混合摩擦區，摩阻系數λ與相對粗糙度ε的關系。

2 經驗公式分析

從莫迪圖可以定性看出，水力光滑區中相對粗糙度ε變化，摩阻系數λ不變；而混合摩擦區的摩阻系數λ會隨相對粗糙度ε的減小而減少。然而，不能明確地觀察出λ隨ε定量變化的形態和趨勢。因此，以工程上常應用的米勒公式(式2)和Colebrook-White公式(式3)為基礎(見表1)，以相對粗糙度ε為變量，雷諾數Re為參量進行數學計算，研究粗糙度對輸油管道摩阻系數的影響，得出摩阻系數隨相對粗糙度變化曲線[10]。

計算結果見圖1。其中由于公式(2)為λ的隱函數，故利用埃特金法[11]進行迭代求得近似解?？傮w上可發現摩阻系數隨相對粗糙度的減小而減小，與莫迪圖結果一致。并且在相對粗糙度較大的區域，摩阻系數減小的更明顯，而在相對粗糙度較小的區域，摩阻系數減小幅度較小，并趨于穩定。高雷諾數管內流體要在相對粗糙度較小下才處于水力光滑區，相比低雷諾數下的流態更容易處于水力光滑區，穩定區域更大。但水力光滑區的摩阻系數隨雷諾數減小而增大。

表1 λ的經驗公式

圖1 摩阻系數隨相對粗糙度變化

在上述計算過程中，判斷流態為紊流水力光滑區還是混合摩擦區是以臨界雷諾數Re1為標志分界線，其值由管道層流邊界厚度與絕對當量粗糙度的比值推導而來。

(4)

對于某一直徑的管子，其臨界雷諾數的數值取決于管壁絕對當量粗糙度。故當管內油流的雷諾數接近臨界值時，合理確定粗糙度數值，成為準確判別流態區和計算λ的關鍵因素。從另一方面來看，式(4)中Re1值是根據勃拉休斯(Blasius)公式等經驗公式得來的，若采用不同的公式推導，將會得出不同的計算式。因此，根據臨界值來劃分流態區只是近似的，實際管道中流態從水力光滑區到混合摩擦區變化有一個過渡階段，使用經驗公式對過渡階段的λ值計算不夠準確，需查詢莫迪圖或進行試驗來修正。

雷諾數為104～105情況下，分別查莫迪圖[9]及用公式計算得到摩阻系數λ，繪制臨近過渡階段的摩阻系數隨相對粗糙度變化曲線見圖2。

圖2 過渡階段摩阻系數隨相對粗糙度變化

由圖2可見，若根據理論公式計算，粗糙度較小時，流體處于在水力光滑區摩阻系數不隨粗糙度變化而變化，曲線是一條直線；隨著粗糙度增大，當進入混合摩擦區時，λ才開始隨相對粗糙度增大而增大。而實際上摩阻系數λ并不是在臨界雷諾數時突變，而是在Re接近Re1前就開始逐漸增大。同時還可以發現，在這一過渡階段之外，公式計算得出的λ值與莫迪圖查詢的值都很相近，也證明了米勒公式和Colebrook-White公式在多數工業管道流體中是適用的。

3 基于Fluent的數值模擬

研究流體流動問題的完整體系由傳統的理論分析方法、實驗測得方法及CFD數值方法組成，各有優缺點，但又有很多聯系，在實際工程中將三者進行有機的結合能得到更合理的結果。用經驗公式對管道內流體進行摩阻計算實際上還不夠完善，特別是對于黏性流體在管中的流動，某些情況下得到的λ值并不準確，還需要通過實驗進行模擬和驗證。但實驗常受到諸多因素的限制，且經費投入較大，相形之下，CFD更有優勢。分析流體流動的基本方程主要是質量守恒、動量守恒、能量守恒三種方程。如果流動涉及不同組分的混合或是相互作用，還要遵守組分守恒定律。在實際計算時，還要考慮不同的流態，對于湍流，目前工程中一般是對瞬態N-S方程做時間平均處理，同時補充反映湍流特性的其他方程，如湍動能方程和湍流耗散率方程等[12]。

本研究中基于FLUENT軟件進行三維數值模擬。模擬對象為一段直管內油品的湍流流動。通過改變入口端流速v來控制雷諾數Re，并分別在不同管壁粗糙度下進行數值模擬，得出處于管中試驗段某兩端面的壓力的差值，根據伯努利方程和達西公式聯立計算即可得到沿程損失阻力系數λ，再將所得的值與莫迪圖曲線或經驗公式計算得出的數值對比，判斷其是否正確。

3.1 物理模型

在GAMBIT中建立直徑300 mm，長4 m的圓截面直管，其中前3 m是前置段，用來讓湍流充分發展，后1 m為試驗段。管內流體選取我國中原地區原油，基本參數：50 ℃密度846.6 kg/m3，運動粘度10.32 mm2/s。

3.2 劃分網格

管內流體的阻力主要是流體與壁面間的作用產生的，在壁面附近會有較大的速度梯度。故在端面上創建邊界層網格加密。具體設置為：first per centage(第一層邊界層網格的高度關于寬度的百分比)為20,Growth factor取1.2，rows取5層。端面網格劃分采用默認方式，單元格為四邊形，網格個數count設為50。得到mesh faces為473，體網格mesh volumes為100 276。網格劃分情況見圖3，圖中坐標X方向為流向，Y、Z方向為徑向。

圖3 模型網格劃分

后續模擬計算中，經過測試，在此網格密度上即使再提高其密度，對計算的精確性也無法進一步明顯的提高，而計算量反而大大增加。因此選用此密度進行網格劃分可達到計算要求。

3.3 邊界條件與初始化

將在GAMBIT中劃分好的網格文件導入FLUENT中。采用標準的k-ε湍流模型與標準壁面函數法(Standard wall function)，變化參數主要有兩個:①入口速度及相關的湍流參數;②壁面粗糙度。具體設置如下。

1)模型求解器(Solver)選擇壓力基、顯式、定常流。

2)粘性模型(viscous model)選擇雙方程模型(k-epsilon[2 eqn]Standard)。

3)流體材料設置為原油。

4)速度進口(velocity-inlet)邊界條件設置為入口端速度，變化范圍為2～6 m/s。湍流方法指定為湍流強度I與水力直徑DH,充分發展的管流核心湍流強度I可按經驗公式來估算。

(5)

式中:ReDH——以水力直徑為特征長度求出，水力直徑定義為流通面積與濕周比值的4倍。

5)出流(outflow)邊界條件保持默認值。

6)固壁(wall)邊界條件設置壁面粗糙度(wall roughness)，在此默認粗糙度厚度值KS(roughness height)等同于絕對當量粗糙度e。變化范圍為0.05～0.30 mm；粗糙度常數CS(roughness constant)保持默認值0.5，認為管壁粗糙度均勻分布。

7)其他條件如求解方法(Solutions)等默認FLUENT的初始設定，使用進口的量對全場完成初始化。

3.4 迭代計算

計算何時收斂對結果的準確性有重要影響，本研究中，模型的收斂采用兩個判斷標準，一是監視每個守恒值的殘差，X、Y及Z方向速度的殘差分別要求低于1×10-7，連續性的殘差低于1×10-2，湍動能k和耗散率ε的殘差低于1×10-6。另一個是出口流量值，只有當出口流量值穩定，且入口流量與出口流量近似相等時，才可以判定為收斂。完成以上所有設定后，開始迭代計算。

4 FLUENT計算結果

4.1 試算情況

以入口速度設置為2 m/s，粗糙度為0.05 mm的邊界條件為例，使用標準的k-ε湍流模式對該模型進行試算。迭代停止后各方向速度及連續性殘差監測見圖4。

圖4 殘差監控

由圖4可見，其收斂條件已經符合上述要求。另外，質量流量通量報告表明，出口與入口的質量流量誤差為4.577 636 7×10-5，也表明計算收斂，應用FLUENT軟件對本例的計算是可行的。

管道出口速度在徑向上的變化見圖5，軸線壓降見圖6。

圖5 出口速度徑向分布

圖6 軸線壓降

由圖5、6可見，管道后端速度呈充分發展狀態，壓降呈線性?？芍四Ｐ陀嬎阌蜻x取足夠長，可應用FLUENT進行模擬計算。

4.2 計算結果分析

通過FLUENT的Report菜單中的surface Integrals命令可以獲得圓管任一端面的面積加權平均壓強值。模擬選擇x=3.4 m和x=3.9 m兩端面的計算值，分別記為pa1和pa2，通過其壓強差得到沿程阻力損失進而算得摩阻系數λ。

(6)

式中:D——管徑0.3 m;L——距離0.5 m;ρ——密度,846.6 kg/m3;v——流速,m/s。

為了驗證FLUENT仿真的準確度，將模擬計算值λ1與通過按表1中經驗公式計算得出的數值λ2進行對比。結果對比見到表2。

由表2可見，若流態為水力光滑區，基于FLUENT仿真的模擬計算結果與經驗公式計算值基本接近，最大誤差為2.33%。在此狀態下對三維管道的數值模擬取得了合理的結果，并具有較高的計算精度。而流態為混合摩擦區時，基于FLUENT的模擬計算結果與經驗公式計算值誤差較大，最大誤差絕對值達15.65%，數值模擬沒有達到較好的效果，需要分析誤差原因。

1)流態為水力光滑區時，相對粗糙度ε變化，摩阻系數λ幾乎不變。模擬計算結果符合理論，證明了摩阻系數λ只和雷諾數有關，與相對粗糙度ε無關。

2)當雷諾數Re很接近臨近雷諾數Re1的值時，模擬計算的摩阻系數λ開始逐漸增大，但不顯著。例如流速為3 m/s，粗糙度為0.25 mm時λ值較之前不隨粗糙度變化的λ的值增大率為1.8%，與公式計算值的誤差為2.23%，此時Re比Re1只近似小于1 800。分析認為FLUENT的數值模擬結果可以證明，出現此變化的原因應該是管道中的流態變化有一個過渡階段，實際上摩阻系數λ并不是在臨界雷諾數時突變，而是在Re接近Re1前的過渡階段就開始逐漸增大。

圖7 數值模擬計算的λ隨ε變化

3)流態為混合摩擦區時，數值模擬的值與公式計算值有較大的誤差，使用標準k-ε湍流模型進行數值模擬的值小于經驗公式計算值。但變化趨勢相同，如圖7，和前面的分析一致，符合莫迪圖。證明λ與Re和ε均有關，在雷諾數相同時，λ隨ε增大而增大。

4.3 誤差原因

1)經驗公式使用的是科爾布魯克-懷特(Colebrook-White)公式，通過試算迭代得到的近似解是不精確的。

2)由于湍流模式和數值計算中都有不可避免的誤差，同時模擬的管道中流體可能處于水力光滑區至混合摩擦區間的過渡階段，不完全適用于科爾布魯克-懷特公式或數值模擬。從數據上來看，隨著雷諾數增大，λ誤差整體上從15%下降為14%。由此可以推測，流體雷諾數比臨界雷諾數Re1大得越多，較完全地處于混合摩擦區時，模擬的計算值越有符合經驗公式值的趨勢。

表2 模擬計算與公式計算結果對比

注：流態Ⅱ表示為水力光滑區，Ⅲ表示為混合摩擦區。

3)標準k-ε湍流模型在壁面區使用了不夠精確的近壁函數的半經驗公式，需要修正模型或精細網格。

為了驗證猜想2)，選取幾組參數進行數值模擬，讓Re值適中地處于Re1和Re2值之間。其中Re2為判斷混合摩擦區和粗糙區的臨界雷諾數，計算式見表1。數據處理同前。結果見表3。

結果表明，數值模擬得到的λ與經驗公式值誤差很小，符合猜想2)。通過FLUENT的數值模擬發現，水力光滑區至混合摩擦區間的過渡階段(Re值接近于Re1的值)的流態，不完全適用于科爾布魯克-懷特公式，或是使用標準k-ε湍流模型進行數值模擬。當流態較完全地處于混合摩擦區時的結果更吻合，能適用于工程。

為了驗證猜想3)，對模型劃分了更精細的網格。具體設置為：first percentage取18,Growth factor取1.1，rows取6層。端面網格劃分采用默認方式，單元格為四邊形，網格個數count設為5。得到mesh faces為616，體網格mesh volumes為146 608。使用標準k-ε湍流模式中非平衡壁面函數法(Non-Equilibrium Wall Function)。以流速為6 m/s，粗糙度為0.3 mm為例，結果見表4。

表3 計算結果

表4 計算結果比較

結果表明，計算精確性沒有明顯提高，而計算量卻大大增加。由此判斷猜想3)不是產生誤差的最主要原因。使用更精細的模型要付出更大的模擬計算量和時間，不適用于工程實際應用。

5 結論

1)數值模擬證實了輸油管道的內壁粗糙度只在流態為混合摩擦區對摩阻有較大影響，即摩阻系數λ隨相對粗糙度ε的增大而增大。且摩阻系數的變化幅度與粗糙度成正比，在相對粗糙度較大的區域，摩阻系數增大的更明顯，在相對粗糙度較小的區域，摩阻系數增大的幅度較??；

2)經驗公式推知流態區的劃分只是近似的，當流態處于變化臨界時摩阻系數用經驗公式計算的準確性有待驗證。仿真分析表明，管道中的流態從水力光滑區到混合摩擦區間有一個過渡階段，在這一階段摩阻系數λ隨與相對粗糙度ε增大而小幅度增大。因此，若想獲取準確的沿程阻力，在此階段計算摩阻系數僅靠經驗公式是不合適的，數值模擬計算是應考慮的方法；

3)在過渡階段之外，流態可以認為是完全處于水力光滑區或混合摩擦區。在此種情況下，經驗公式計算、莫地圖查詢與數值模擬得出的λ值都很相近，證明了米勒公式和Colebrook-White公式在多數工業管道流體沿程阻力計算中具有一定的適用性和準確性。

[1] VIJIAPURAPU，S， JIE C. Numerical Simulation of Turbulent Flow in a Ribbed Pipe Using Large Eddy Simulation[J]. Numerical Heat Transfer, Part A: Application(ISSN:1040-7782), 2007,51(12):1137-1165.

[2] VIJIAPURAPU，S，.& JIE C. Performance of turbulence models for flows through rough pipes[J], Applied Mathematical Modelling(ISSN: 0307-904X), 2010,34(6):1458-1466.

[3] ZHANG,Y B & Che,D F. Numerical simulation of turbulent flow and wall mass transfer in a rectangular channel roughened by v-shaped grooves[J],Numerical Heat Transfer,Part A:Application(ISSN:1040-7782),2014,66(5):551-581.

[4] 張建卓, 李康康, 高 猛,等. 粗糙度對流體泄漏和摩擦阻力的仿真研究[J].微計算機信息,2011,27(6):238-240.

[5] 孫 琳.圓管突擴(縮)局部阻力試驗及流場數值模擬[D].哈爾濱:東北農業大學,2010.

[6] 趙 月.基于CFD的管道局部阻力的數值模擬[D].大慶:東北石油大學,2011.

[7] 楊筱蘅.輸油管道設計與管理[M].北京:中國石油大學出版社,2001.

[8] 楊樹人,汪志明,何光渝,等.工程流體力學[M].北京:石油工業出版,2006.

[9] 詹德新,王家楣.工程流體力學[M].武漢:湖北科學技術出版社,2001.

[10] 楊曉鴻,朱薇玲.粗糙度對輸氣管道摩阻系數的影響[J].石油化工設備,2005,34(1):25-28.

[11] 羅祖帖,Excel在管路阻力計算中的應用[J].化肥設計,2004,42(2):15-17.

[12] 王福軍.計算流體動力學分析:CFD軟件原理與應用[M].北京：清華大學出版社,2004.

Numerical Simulation of Effect of Wall Roughness on the On-way Resistance Based on FLUENT

YAO Xue-lei, YUAN Cheng-qing, FU Yi-feng, BAI Xiu-qin

(a. Reliability Engineering Institute, School of Energy and Power Engineering; b. Key Laboratory of Marine Power Engineering and Technology (Ministry of Transport), Wuhan University of Technology, Wuhan 430063, China)

Aiming at the problem of obtaining accurately the on-way resistance of submarine pipeline at the design stage, numerical simulation method is adopted to analyze the relationship between loss of viscous fluid friction and wall roughness of the pipeline based on FLUENT. Comparing the calculated results with theoretical results and the actual situations, the results show there is a good agreement, except for the flow pattern transiting from the hydraulic smooth region to the mixed friction region. It is confirmed that the friction coefficient is strongly influenced by pipeline wall roughness only in the turbulent mixing friction region. In this case, reducing the surface roughness of the wall can cut down the resistance. If the flow pattern transits from the hydraulic smooth region to the mixed friction region, it is inappropriate to calculate friction coefficient only by the empirical formulas, and numerical simulation is the method that should be considered.

oil pipeline; roughness; on-way resistance; FLUENT; numerical simulation

10.3963/j.issn.1671-7953.2015.06.024

2015-07-17

國家自然科學基金優秀青年基金項目 (51422507)

姚雪蕾(1992-)，女，碩士生

TE832；P756.2

A

1671-7953(2015)06-0101-07

修回日期：2015-09-10

研究方向:摩擦學系統及表面工程

E-mail: yaoxuelei2013@outlook.com