戚菁菁，初同龍，孫伶

1.中國石油天然氣管道工程有限公司（河北 廊坊 065000）

2.中石油管道局工程有限公司 燃氣分公司（河北廊坊 065000）

3.中國石油管道分公司 科技處（河北 廊坊 065000）

0 引言

我國所產(chǎn)原油80%以上為高含蠟、高凝點、高黏度的蠟基質(zhì)原油[1]。對于易凝高黏的含蠟原油，由于管道周邊為自然溫度場，環(huán)境溫度相對較低，因此若管道在常溫下運行容易因沿線溫降過快導致沿程摩阻急劇升高而被迫停輸，甚至可能發(fā)生凝管事故。由于管道溫度場的變化會顯著影響熱油管道的最小輸量及安全運行[2]，因此準確模擬出含蠟原油管道周圍溫度場分布至關重要。

在計算管道內(nèi)油流溫度場時，常選用一維、二維或三維溫度場計算。管線結蠟層、管壁、防腐層及土壤之間傳熱在物理上是復雜的三維問題，然而軸向溫度梯度相對徑向溫度梯度非常小，因此工程上可以忽略軸向溫降影響，這樣三維導熱問題就簡化為二維土壤溫度場模擬問題。文獻[3]將土壤溫度場假設為一維分布即溫度變化只發(fā)生在縱向上，由于這種方法過于簡單，計算結果同實際測量的結果相比有非常大的誤差，因此工程上不予使用。

對于一維、二維和三維溫度場，其計算方法通常包括解析法、實驗方法以及數(shù)值計算等[1-3]。解析法按Fourier 的熱傳導方程進行求解，能夠獲得問題的精確解，并且可以為數(shù)值計算和實驗法提供比較依據(jù)。但是該方法存在很大的局限性，即使求解簡單的導熱問題也相當復雜，對于復雜幾何形狀的物體和非線性邊界條件下的導熱問題，應用解析法幾乎是不可能的。實驗方法是傳熱學的基本求解方法，通過實驗可以得出相對準確的結果，但是此種方法的適應性不好，針對某個特定的問題總是要進行新的實驗，耗費的人力、物力和財力巨大，并且實驗周期長。數(shù)值計算在很大程度上彌補了解析法和實驗方法的缺點，其適應性強，特別是對于復雜問題更顯其優(yōu)越性[2-4]。

針對非網(wǎng)格劃分對不規(guī)則區(qū)域具有良好適應性及容易實現(xiàn)網(wǎng)格自動化生成等特點[5-6]，本文采用拉普拉斯變換將土壤的半無限大區(qū)域轉換為有限矩形區(qū)域，考慮到土壤的物性變化，采用基于非結構化網(wǎng)格離散管線周圍的土壤溫度場。

該類型網(wǎng)格常用的基本單元是三角形和四面體，生成方法包括規(guī)則劃分法、Delaunay 三角剖分法以及陣面推進法等[5-6]。對二維而言，三角形網(wǎng)格可以充滿任意幾何形狀的區(qū)域。雖然其生成過程較為繁瑣，數(shù)據(jù)結構相對復雜，但具有自適應的優(yōu)點。因此選用二維非結構化網(wǎng)格的數(shù)值計算方法模擬計算凍土區(qū)含蠟原油管道溫度場。

1 數(shù)學模型

埋地管道正常運行過程中，傳熱過程包含3 部分，分別為管內(nèi)油流以熱對流方式傳熱至結蠟層，再通過結蠟層、管壁、防腐層等將熱量交換至周圍土壤，最后經(jīng)地面與大氣換熱。內(nèi)部因素和外界環(huán)境因素共同構成了影響熱油管道散熱的影響因素。其中內(nèi)部因素有油品熱物理性質(zhì)、管道輸量、管徑以及結蠟層、管壁、防腐層等。外界因素有埋地管段土壤物性、大氣溫度、風速等。其物理結構模型如圖1、圖2所示。

根據(jù)簡化，基于圖1、圖2，綜合考慮界面上各層之間的相互影響，建立熱油管道正常運行過程水力-熱力模型。

1.1 流體方程

質(zhì)量守恒方程：

式中：ρ 為密度，kg/m3；A 為管流截面面積，mm2；V 為流體平均速度，m/s；τ為時間，s；z為管道軸向位置，m。

圖1 埋地管道剖面

圖2 埋地管道示意圖

動量守恒方程：

式中：α 為管道軸向與水平方向的夾角，（°）；g 為重力加速度，m/s2；p 為油流截面平均壓力，MPa；f 為摩阻系數(shù)；D為管道內(nèi)直徑，mm。

考慮到油流的速度方向，V2應該改為V|V|，不失一般性，本文中予以簡化。

能量方程：

式中：s 為原油比熵，J/（kg·℃）；u 為原油比內(nèi)能，J/kg；h為原油比焓，J/kg；q為單位時間內(nèi)原油在單位管壁面積上的散熱量，J/（mm2·s）。

由質(zhì)量守恒方程、動量守恒方程和能量守恒方程得油流的換熱方程：

式中：Cp為原油定壓比熱容，J/（kg·℃）；T 為原油溫度，℃。

1.2 結蠟層、管壁和防腐層的傳熱方程

式中：Cj為第j 層（結蠟層、管壁和防腐層）比熱容，J/K；λj為第j 層（結蠟層、管壁和防腐層）導熱系數(shù)，kcal/（m·h·℃），其中，j=1、2、3，分別表示結蠟層、管壁和防腐層，其中1 kcal=4.184 kJ；r 為徑向位置，mm；θ為環(huán)向弧度，rad。

1.3 土壤導熱方程

式中：下標為s 的參數(shù)代表土壤的參數(shù)；x 為垂直于軸向的水平位置，m；y為深度，m。

1.4 關聯(lián)方程

管內(nèi)流體、蠟沉積、管壁、防腐層以及土壤的熱交換過程相互關聯(lián)，滿足如下方程：

式中：α0為流體對管內(nèi)壁的放熱系數(shù)，W/(m2·℃)；T0為管內(nèi)壁溫度，℃；Rj為各層厚度，mm。

1.5 滿足的邊界條件

埋地管道的界面在物理上是對稱的，其計算區(qū)域也具有對稱性，僅取管道截面右半部分研究，應滿足的邊界條件如下：

式中：αa為地表向大氣的放熱系數(shù)，W/(m2·℃)；ho為管道埋深，m；H 為縱向熱力影響區(qū)，m；Ta為大氣溫度，℃；Tn為恒溫層溫度，℃；De為管道外徑，mm。

2 計算區(qū)域離散化及數(shù)值計算方法

2.1 計算區(qū)域離散化

土壤導熱計算區(qū)域應為半無限大土壤介質(zhì)區(qū)域。但解析求解需對問題作簡化，因而所獲得的結果與真實值偏差較大[7]。因此可以通過控制網(wǎng)格劃分，獲得較高精確度的解。本研究的計算區(qū)域限制為原油管道的熱力影響區(qū)，采用數(shù)值計算方法。在數(shù)值計算之前要對求解區(qū)域離散（網(wǎng)格劃分），如圖3—圖5所示。

圖3 土壤區(qū)域密集網(wǎng)格劃分

圖4 土壤區(qū)域稀疏網(wǎng)格劃分

圖5 結蠟層、管壁和防腐層網(wǎng)格劃分及周圍的土壤網(wǎng)格

采用有限差分法（FDA）對含蠟原油管道進行離散化，如圖6 所示。計算節(jié)點從前一個泵站的出口（點1）到后一個泵站的入口（點n）。

圖6 管道離散圖

2.2 控制方程離散化

對結蠟層、管壁、防腐層及土壤采用控制容積法研究。

2.2.1 油流方程離散

采用有限差分法，對流體非穩(wěn)態(tài)水力-熱力耦合方程進行離散求解。如圖6所示，沿z方向將管段劃分為若干離散的小管段，在ΔZ時間段內(nèi)，對方程（4）兩邊進行離散，得到：

2.2.2 對結蠟層、管壁及防腐層的控制方程離散

結蠟層、管壁、防腐層中的導熱方程可以寫成如下統(tǒng)一形式：

圖7 極坐標中的網(wǎng)格系統(tǒng)

結合圖7，將式（19）整理成通用的離散化方程形式：

式中：

式中：下標E、S、W、N分別代表東、南、西、北4個方位。

2.2.3 土壤控制方程的離散

將計算節(jié)點置于三角形的重心，如圖7所示，節(jié)點P0可視為陰影三角區(qū)域的代表，稱該三角形為P0的控制容積。離散化導熱方程建立起計算節(jié)點P0與相鄰點P1,、P2、P3的溫度的代數(shù)關系式。土壤導熱方程（6）可以針對任意的控制容積寫成如下形式：

3 算例

選取我國東部地區(qū)某含蠟原油管道，該管道管體大部分處在凍土區(qū)，存在凍脹融沉風險，因此需要了解其溫度場變化，從而確定管道安全輸量，減小凍脹融沉風險。該管道參數(shù)如下：管頂埋深2.2 m；通過對管道運營的土壤進行試驗分析得到，土壤分層（0～3 m為含水20%亞砂土，3～10 m為含水25%粉質(zhì)黏土，10～20 m為基巖），熱力影響范圍為20 m×15 m；無保溫層。凍土導熱系數(shù)λ1=5.77 kJ/（m·h·℃），融土導熱系數(shù)λ2=5.77 kJ/（m·h·℃），土壤比熱容為C= 1 900 J/( kg·℃) ，土壤表面對大氣的放熱系數(shù)aK= 17. 5 W/ (m2·℃) ，土壤初始溫度T0=-12.5 ℃,下邊界溫度為-5 ℃，地面空氣溫度TK=-15 ℃。假定凍土厚度6 m，根據(jù)管道冬季各站場原油進出站溫度（最低2 ℃左右，最高5 ℃左右）分別假定油溫為2 ℃和5 ℃，預測油溫下土壤溫度場及凍融情況。

根據(jù)1.1及2.1可得到計算網(wǎng)格劃分如圖8所示。

圖8 網(wǎng)格劃分

采用自主開發(fā)的網(wǎng)格計算軟件，具體參數(shù)見2.1，計算結果如圖9—圖11所示。管道周圍綠色線所圍為融化區(qū)，即不存在凍土，下方接近水平的線是多年凍土下限，同時給出油溫5 ℃和2 ℃時，0 ℃等值線，從而確定凍土區(qū)域。

圖9 2 ℃油溫條件下土壤溫度場云圖和等值線

圖10 5 ℃油溫條件下土壤溫度場云圖和等值線

圖11 0 ℃等值線

1）油溫5 ℃時，管道熱力影響較大，不存在凍融情況。

2）油溫2 ℃，管道只對周圍10 m 內(nèi)土壤產(chǎn)生熱影響，溫度場零度線在管道外側，會造成凍土融化。

本文的土壤參數(shù)、地表溫度邊界、現(xiàn)場凍土情況和油溫歷史等條件采用了經(jīng)驗數(shù)據(jù)，在進行數(shù)值模擬后得出的結果對管道安全運行具有指導作用。

4 結束語

基于非結構化網(wǎng)格數(shù)值計算方法，針對含蠟原油管道外溫度場建立了二維非穩(wěn)態(tài)熱力計算的數(shù)學模型。將該模型應用于某凍土區(qū)含蠟原油管道溫度數(shù)值模擬，通過模擬得到該管道2 ℃及5 ℃油溫條件下土壤溫度場云圖和等值線圖，了解管道凍融沉降區(qū)域，對于指導管道安全運行和節(jié)能降耗具有重要意義。