李國豪，鄧道明，宮敬

（中國石油大學(xué)（北京）油氣管道輸送安全國家工程實驗室，石油工程教育部重點實驗室，城市油氣輸配技術(shù)北京市重點實驗室，北京102249）

引 言

國內(nèi)外天然氣集輸經(jīng)常采用濕氣集輸?shù)姆绞?。Meng[1]指出標(biāo)況（101.3 kPa, 15℃）下液相體積流量與氣相體積流量之比小于1100 m3/106m3時為低液相負荷流動。濕氣管道一般為低液相負荷流動，例如，我國某頁巖氣田濕氣管道氣液分離之前的平均液相負荷一般為200 m3/106m3左右，最大的平均液相負荷為1000 m3/106m3。但是即便液相負荷如此低和井場集氣站進行了氣液分離，該頁巖氣田的一些氣井和早期集輸管線仍存在較嚴重積液問題。

早期的一些學(xué)者認為，低液相負荷下天然氣兩相管流為彌散霧狀流，即認為低液相負荷流動不會產(chǎn)生嚴重積液現(xiàn)象。例如Adewumi 等[2]在美國國家自然科學(xué)基金、工程基金和氣體研究所(GRI)的資助下，對低液相負荷天然氣-凝析液管流進行了一系列研究，其模擬計算理論主要基于漂移的彌散霧狀流模型。但是其模型計算結(jié)果沒有得到實驗或生產(chǎn)數(shù)據(jù)的廣泛支持。顯然濕氣管道的彌散霧狀流模型難以解釋經(jīng)過分離器分離的濕氣在下游集輸管道里仍然積液嚴重的原因。

正常情況下，濕氣管道內(nèi)低液相負荷的流動為氣液兩相分層流。當(dāng)氣田開發(fā)初期開井?dāng)?shù)較少或氣田開發(fā)后期氣井產(chǎn)量較低時，集輸管內(nèi)氣體流速較低，如果氣相對液相的拖曳力不足以拖動液相向前流動，氣液分層流便開始變得不穩(wěn)定，最終液體會在上傾管道低洼處積集，積液現(xiàn)象因此產(chǎn)生。另外過大的集輸管道直徑也會引起積液。本文將上傾管道液膜不穩(wěn)定的最大氣相表觀速度（或者說液膜穩(wěn)定的最小氣相表觀速度）作為臨界氣速。工程實踐表明，起伏的濕氣集輸管道積液是普遍現(xiàn)象。

濕氣集輸管道積液的存在會導(dǎo)致諸多風(fēng)險：（1）誘使管內(nèi)產(chǎn)生內(nèi)腐蝕，進而影響生產(chǎn)運行[3]；（2）集輸管道積液意味著存在更大的相間滑脫損失，這會導(dǎo)致井口回壓增加，從而降低氣井產(chǎn)量，甚至可能造成氣井積液；（3）集輸管道積液還伴隨著管道集輸效率降低、超壓安全風(fēng)險、下游處理量不穩(wěn)定等問題。因此，研究濕氣管道的積液對集輸管道的設(shè)計和運行都具有重要的指導(dǎo)意義。

1 積液預(yù)測模型評述

一些學(xué)者采用商業(yè)多相流軟件OLGA 對濕氣管道進行積液研究。例如劉建武等[4]、張愛娟等[5]建立起伏的集輸管道的OLGA 計算模型，經(jīng)過OLGA 計算獲得了管道全線的積液量等相對宏觀的參數(shù)。但是，OLGA 軟件作為一個黑箱，使用者并不知曉其具體的計算模型和方法，難以利用OLGA 建立積液機理計算模型。

目前主流的積液機理模型有最小壓降模型、持液率突變模型、零液壁剪切應(yīng)力模型和最小界面剪切應(yīng)力模型4種。

1.1 最小壓降模型

最小壓降模型是指：逐漸減小氣相流量，當(dāng)氣液兩相管流的總壓降最小時，管內(nèi)由摩阻壓降主導(dǎo)轉(zhuǎn)換為重力壓降主導(dǎo)，管內(nèi)發(fā)生積液。Belt[6]、Yuan等[7]均支持該模型。Fan 等[8]的實驗也表明，當(dāng)水平傾角大于5°時，最小壓降與臨界氣速相對應(yīng)。

但是，最小壓降并不能完全解釋積液的機理。盡管臨界氣速測量方法相同，Alsaadi 等[9]和Wang等[10]的實驗結(jié)果卻截然不同：Alsaadi等[9]的實驗管徑是76.2 mm，在實驗中觀察到最小壓降對應(yīng)的氣速要小于臨界氣速；Wang 等[10]的實驗管徑是30 mm 和40 mm，實驗結(jié)果卻顯示最小壓降對應(yīng)的氣速要略微大于臨界氣速。這表明最小壓降可能與臨界氣速并不完全對應(yīng)。Skopich 等[11]也觀察到最小壓降有時存在于穩(wěn)定流中，有時又存在于不穩(wěn)定流中。

1.2 持液率突變模型

對于近水平氣液兩相流流動，僅改變氣相流量，由動量方程計算出的持液率可能會如圖1所示。圖1 中陰影區(qū)域?qū)?yīng)的持液率有多個，學(xué)術(shù)上將該區(qū)域稱為持液率多解區(qū)。許多學(xué)者從持液率解穩(wěn)定性的角度對多解區(qū)進行研究：Landman[12]利用KH 波理論對界面波穩(wěn)定性進行分析，得出持液率的高解是不存在的，中間解可能存在，低解是始終存在的；Barnea 等[13]在Landman[12]的基礎(chǔ)上考慮了分層流的結(jié)構(gòu)穩(wěn)定性，指出中間解不滿足結(jié)構(gòu)穩(wěn)定性，只有最低解是穩(wěn)定的。這表明在圖1多解區(qū)的左邊界處，持液率由A點突變到B點。

Espedal[14]、Langsholt 等[15]、Kjolaas 等[16]、Fan[17]均在實驗中觀察到了持液率突變現(xiàn)象。Biberg 等[18]將分層流持液率突變點作為臨界積液點。沈偉偉[19]綜合采用“最小滑脫”準(zhǔn)則和分層流多解方法判斷持液率突變，從而計算相應(yīng)的臨界氣速。該模型的缺點是適用范圍可能較窄，因為當(dāng)管道傾角和液相負荷較大時，持液率可能不突變。

圖1 持液率多解區(qū)Fig.1 Regional schematic diagram for multiple liquid holdup solutions

1.3 零液壁剪切應(yīng)力模型

零液壁剪切應(yīng)力模型最早由Turner 等[20]提出，其認為環(huán)狀流液膜轉(zhuǎn)向的臨界狀態(tài)是液壁剪切應(yīng)力為零，將液壁剪切應(yīng)力為零的狀態(tài)作為氣井積液的臨界點；然而他們認為零液壁剪切應(yīng)力模型的預(yù)測值與現(xiàn)場數(shù)據(jù)不符，采用液滴模型預(yù)測臨界氣速。

Fan[17]將該模型引入近水平管的分層流模型，將液膜最底部的壁面剪切應(yīng)力為零作為積液的臨界判據(jù)。Alsaadi[21]在Fan[17]的基礎(chǔ)上進行簡化，將液膜的平均壁面剪切應(yīng)力為零作為積液的臨界判據(jù)，進而把混合動量方程中的氣速和持液率解耦，節(jié)約了計算時間。若將分層流的液膜單獨作為分析對象，氣液界面剪切應(yīng)力是主動力，液壁剪切應(yīng)力是被動力，液壁剪切應(yīng)力為零是液膜轉(zhuǎn)向的現(xiàn)象，其轉(zhuǎn)向的本質(zhì)原因是氣體攜帶能力不足，從界面剪切應(yīng)力的角度出發(fā)能夠進一步揭示積液的原因。

1.4 最小界面剪切應(yīng)力模型

Barnea[22]提出了垂直管環(huán)狀流的最小界面剪切應(yīng)力模型，其假設(shè)環(huán)狀流液膜沿管道周向均勻分布，得出對于一定的液相表觀速度，界面剪切應(yīng)力僅為液膜厚度/持液率的函數(shù)，如圖2 所示。圖2 中界面剪切應(yīng)力曲線存在最小值，Barnea[22]認為在最小界面剪切應(yīng)力左側(cè)的環(huán)狀流是穩(wěn)定的，其右側(cè)的環(huán)狀流是不穩(wěn)定的。Luo 等[23]和Shekhar 等[24]均認為液膜沿管道周向均勻分布的假設(shè)并不合理，并做了相應(yīng)修正：Luo等[23]用管道底部的最大液膜厚度進行計算；Shekhar 等[24]根據(jù)Paz 等[25]的實驗數(shù)據(jù)，考慮了最大液膜厚度與平均液膜厚度的關(guān)系，最終通過平均液膜厚度求出臨界氣速。沈偉偉等[26]在Shekhar等[24]的研究基礎(chǔ)上考慮了液滴夾帶的影響。

圖2 最小界面剪切應(yīng)力示意圖Fig.2 Schematic diagram of minimum interfacial shear stress

Brito[27]將最小界面剪切應(yīng)力模型從環(huán)狀流擴展到了分層流，并認同最小界面剪切應(yīng)力左側(cè)對應(yīng)穩(wěn)定流，右側(cè)對應(yīng)不穩(wěn)定流的觀點。本文對穩(wěn)定區(qū)與不穩(wěn)定區(qū)的劃分也是如此。Brito[27]通過窮舉的算法尋找最小界面剪切應(yīng)力對應(yīng)的臨界氣速。 Brito 模型的不足是通過Blasius 經(jīng)驗關(guān)系式來計算液壁剪切應(yīng)力，因為Kowalski[28]測量氣液兩相分層流壁面剪切應(yīng)力的實驗顯示：氣壁剪切應(yīng)力與Blasius 關(guān)系式吻合得較好，但液壁剪切應(yīng)力與Blasius 關(guān)系式明顯不符。

2 模型建立

Langsholt 等[15]氣液兩相流的實驗顯示：氣液界面波在高氣速下是3D小尺度波（滾波），隨著氣速減小，界面波的波幅逐漸變大，依次出現(xiàn)2D 短波和2D長波，最終流型由分層流變?yōu)榧俣稳鳌an 等[8]將實驗中出現(xiàn)的分層流流型分為孤波（光滑分層流）和滾波兩種，孤波在高氣速下出現(xiàn)，滾波在接近流型轉(zhuǎn)變時出現(xiàn)。Alsaadi[21]的實驗中所有的分層流均是滾波分層流和長波分層流，長波分層流在接近流型轉(zhuǎn)變時出現(xiàn)。以上實驗現(xiàn)象均表明分層流在失穩(wěn)過程中，界面波的波幅會變大。

可以推測，積液和界面波的變化有關(guān)，界面波的變化又直接與界面摩擦因子和界面剪切應(yīng)力相聯(lián)系?？梢哉J為積液是由于氣相給予液相的拖曳力不夠，即認為氣液界面剪切應(yīng)力為最小值時為積液臨界狀態(tài)。

石油工業(yè)濕氣管道的管徑一般在2"以上，例如我國某頁巖氣開發(fā)采氣管線公稱直徑為DN65，濕氣集輸管道的最大管徑為DN550。Fan[29]在3"管道的低液相負荷氣液兩相流動的實驗中觀察到氣液界面是接近水平的。鄧道明等[30]曾建立高壓大管徑天然氣兩相流動計算模型，通過將模型計算結(jié)果與生產(chǎn)數(shù)據(jù)比較，認為平界面分層流模型與高壓大直徑天然氣管道工藝計算更為貼合。這里假設(shè)分層流的界面為水平界面，即管道橫截面液膜分布如圖3 所示。利用Biberg[31]基于雙極坐標(biāo)系導(dǎo)出的液相速度分布解析式，將分層流界面剪切應(yīng)力最小作為積液判據(jù)，建立新的濕氣管道積液預(yù)測模型。

圖3 濕氣管道內(nèi)分層流液膜分布Fig.3 Film distribution of stratified flow in wet gas pipeline

2.1 動量方程

由于分層流流型是典型的分離流，故采用經(jīng)典的雙流體模型；該模型的主要特點是考慮了相間的相互作用力。氣、液相的動量方程為：

式中，下角標(biāo)g 和l分別代表氣相和液相，i表示界面；幾何參數(shù)A 代表各相所占管道的橫截面積，S代表各相的濕周；剪切應(yīng)力τw表示各相壁面剪切應(yīng)力，τi表示界面剪切應(yīng)力；dp/dz表示各相沿流動方向（管道軸向）壓力梯度；ρgsinθ 表示各相重力沿軸向的分力。

2.2 分層流液膜區(qū)液相流量公式化

由于濕氣管道氣液兩相流動是低液相負荷的流動，其液相Reynolds 數(shù)一般要比轉(zhuǎn)捩Reynolds 數(shù)（Re=2000）低，因此其液膜流動可認為是層流。Biberg[31]曾基于雙極坐標(biāo)系，假設(shè)氣液兩相分層流的氣相為湍流、液相為層流、氣液界面為平界面，并且假設(shè)氣液界面上的剪切應(yīng)力均勻分布，利用N-S 方程推導(dǎo)出了雙極坐標(biāo)系下液膜區(qū)的速度分布。對液膜區(qū)的液相速度求面積分，可以得到液相流量：

φB(φ)和φi(φ)是液面角φ 的函數(shù)。如圖4 所示，當(dāng)φ∈[0,π/2]時，φB(φ)和φi(φ)均單調(diào)遞增，且當(dāng)φ→0時，φB(φ)和φi(φ)→0。

圖4 權(quán)函數(shù)φB(φ)和φi(φ)Fig.4 Weight functions φB(φ)and φi(φ)

2.3 界面剪切應(yīng)力

忽略液相斷面上由位置高度的變化而引起的壓力梯度，則氣、液兩相的壓力梯度dp/dz應(yīng)相等，將式(3)代入氣相動量方程式(1)，可以得到界面剪切應(yīng)力τi與液面角φ和氣壁剪切應(yīng)力τwg之間的關(guān)系：

式中，氣壁剪切應(yīng)力τwg的計算式為式(8)，其中氣壁摩擦因子fg的計算采用Blasius關(guān)系式[32]。

將式(8)代入式(7)，并利用Biberg[31]推出的液面角φ與持液率Hl的顯式表達式，得式(9)

圖5 τi與Hl和usg的關(guān)系[式(9)]Fig.5 Interfacial shear stress τivs holdup Hl and superficial gas velocity usg[Eq.(9)]

表1 3″上傾空氣-水管道積液實驗基本參數(shù)Table 1 Air-water experimental parameters in 3″upwardly inclined pipeline

由式(9)可知，當(dāng)D、θ、ρg、ρl、usl已知時，τi僅是Hl和usg的函數(shù)；圖5 為在表1 參數(shù)下τi與Hl和usg關(guān)系的曲面圖，工程上濕氣集輸管道因其較低的液相負荷，氣液界面剪切應(yīng)力與持液率和氣相折算速度關(guān)系也基本如圖5 所示。由圖5 可知，當(dāng)usg不變、Hl∈[0,0.5]（即φ∈[0,π/2]）時，由于式(9)等號右邊第一項隨Hl的減小而遞增，第二項隨Hl的減小而遞減，使得圖5 呈現(xiàn)下凹的特性，在極值點?τi(Hl)/?Hl=0 處取得最小界面剪切應(yīng)力τi_min，對應(yīng)的持液率稱為臨界持液率Hl_c。另外，從圖5 中可知，Hl_c對usg的變化不敏感（圖中與usg軸基本平行的深色帶區(qū)域），這為后續(xù)快速求解臨界氣速創(chuàng)造了條件。

τi還可通過封閉關(guān)系式(10)得出，式中的界面剪切因子fi計算值對臨界氣速的預(yù)測至關(guān)重要；但是目前的實驗手段并不能直接測出τi的大小，fi的作用機理仍然是一個未完全解決的問題，目前似乎還沒有一種普適性的fi經(jīng)驗關(guān)系式。本文模型中的fi采用Alsaadi[21]提出的關(guān)系式(11)。

圖6 與τi的交點Fig.6 Intersection of τi from Eq.(9)and from Eq.(10))。

2.4 臨界氣速預(yù)測

3 模型評估

圖7 臨界氣速計算流程Fig.7 Flow chart for calculating critical gas velocity

為了對新模型進行評估，本文收集了現(xiàn)有文獻中可以獲得的88 組實驗數(shù)據(jù)用于評估模型。數(shù)據(jù)來源于表2 中的文獻，其中Alsaadi[21]和Rodrigues[33]提供了分層流向假段塞流轉(zhuǎn)變的數(shù)據(jù)，其余文獻給出了臨界氣速的實驗數(shù)據(jù)。但是各研究者對于臨界氣速的識別方法有所不同。Langsholt 等[15]將持液率突變對應(yīng)的氣速作為臨界氣速；Brito[27]、Fan[17]、Nair[34]、Alsaadi等[9]將液膜回流對應(yīng)的氣速作為臨界氣速。上述幾位研究者對液膜回流的檢測方法又有所不同，Brito[27]和Nair[34]往液膜中注入鹽水，然后通過探針來感知液膜有無回流；Fan[17]則是往液膜注入顏料，通過注入點上游顏色變化來判斷液膜有無回流；Alsaadi 等[9]僅通過高速攝像機對液膜回流進行感知，由于攝像機只能捕捉到相對宏觀的現(xiàn)象，因此其探測液膜回流對應(yīng)的氣速精度可能比前兩種方法低。

如圖8所示，對于接近水平管道來說，不論是實驗數(shù)據(jù)還是模型預(yù)測值都存在如下趨勢：隨著傾角和液相負荷的增加，臨界氣速也隨之增加。

就臨界氣速而言，新模型預(yù)測值普遍高于Alsaadi 等[9]的實驗值，但略微低于Fan[17]的實驗值；整體而言新模型預(yù)測值更接近Fan[17]的實驗結(jié)果。對此有以下解釋。

（1）Alsaadi 等[9]與Fan[17]實驗的管徑均為3″，壓力為常壓，流體介質(zhì)均為空氣/水，當(dāng)液相表觀流速為0.01 m/s 時，兩者的實驗條件一致，但是兩者的臨界氣速實驗結(jié)果有差異。這可能是由于檢測液膜回流的方法不同所致，前者僅采用高速攝像機，后者注入了顏料，前者感知力要弱于后者，故前者的臨界氣速實驗值要小于后者的實驗值。

表2 積液實驗條件參數(shù)匯總Table 2 Parameters of liquid accumulation experiments

圖8 新模型與文獻[9,17]的實驗數(shù)據(jù)對比Fig.8 Comparison of predictions from new model with experimental data Ref.[9,17]studies

（2）界面摩擦因子對于臨界氣速的預(yù)測至關(guān)重要，但目前沒有一種界面摩擦因子的經(jīng)驗關(guān)系式能夠適用所有工況。所以模型的預(yù)測值與實驗值的偏差也可能是由于界面摩擦因子的不準(zhǔn)確所致。

圖9展示了各模型預(yù)測的臨界氣速與來源于表2 中88 個相應(yīng)實驗數(shù)據(jù)點之間的偏差。其中Alsaadi[21]和Rodrigues[33]的實驗臨界氣速取的是分層流最小氣速和假段塞流最大氣速的平均值。所選取的模型為各判據(jù)中最具有代表性的模型，其中新模型和Brito[27]模型均是以分層流的最小界面剪切應(yīng)力為積液判據(jù)，Alsaadi[21]模型以分層流的零液壁剪切應(yīng)力為積液判據(jù)，而沈偉偉[19]模型基于持液率突變準(zhǔn)則。圖9的主圖表示各模型預(yù)測值與全部實驗數(shù)據(jù)的偏差。整體上看，圖9(b)Brito[27]模型的預(yù)測效果最差，Brito[27]模型的預(yù)測值較多位于±25%相對偏差線之外；另外3 組模型的預(yù)測值則大部分在±25%偏差線之內(nèi)。圖9的縮略圖為各模型預(yù)測值與來源于Rodrigues[33]、Langsholt 等[15]、Espedal[14]的比較接近高壓力、大管徑實驗值的對比；其中Langsholt 等、Espedal 用較重的SF6 代替空氣，用于模擬高壓力下氣相。從縮略圖中可明顯看出只有新模型的預(yù)測值基本在±10%相對偏差線之內(nèi)。

表3所示，與全部實驗數(shù)據(jù)點相比，新模型預(yù)測值的平均相對偏差為5%，標(biāo)準(zhǔn)差為16%；而與接近高壓力、大管徑的數(shù)據(jù)對比，新模型的平均相對偏差為-1%，標(biāo)準(zhǔn)差為8%。綜合來說，本文模型最優(yōu)，文獻[19]模型次之，文獻[27]模型的預(yù)測結(jié)果與實驗數(shù)據(jù)偏差最大。

表3 各模型預(yù)測偏差Table 3 Prediction deviations for new and other models

4 結(jié) 論

積液現(xiàn)象發(fā)生的本質(zhì)是氣相對液相的攜帶能力不足?；谧钚〗缑婕羟袘?yīng)力的積液判據(jù)，引入平氣液界面分層流液膜區(qū)的速度場描述，建立了近水平濕氣管道積液臨界氣速預(yù)測新模型，提出了求解臨界氣速的簡捷算法。本研究結(jié)論如下。

（1）對于近水平管道，模型預(yù)測的臨界氣速隨著傾角和液相負荷的增加而增加，這與實驗的變化趨勢一致。

圖9 各模型預(yù)測的臨界氣速與表2中實驗值的偏差Fig.9 Comparison of predicted critical gas velocities from new and other typical models with experimental data

（2）利用已有的實驗測量臨界氣速數(shù)據(jù)，對新模型和Brito、Alsaadi、沈偉偉三種近水平管積液模型進行評估，結(jié)果顯示新模型預(yù)測的臨界氣速與已有的實驗數(shù)據(jù)最為吻合，沈偉偉的模型次之，Brito模型的預(yù)測結(jié)果與實驗數(shù)據(jù)偏差最大。新模型與全部實驗數(shù)據(jù)的整體平均相對偏差為5%，標(biāo)準(zhǔn)差為16%。

（3）與其中的接近高壓力、大管徑數(shù)據(jù)相比，新模型預(yù)測精度更加優(yōu)于其他模型，此時新模型平均相對偏差為-1%，標(biāo)準(zhǔn)差為8%。

符 號 說 明

A——面積，m2

d——管內(nèi)徑，m

Fr——Froude數(shù)

f——摩擦因子

g——重力加速度，m/s2

Hl——持液率

hl——液膜厚度，m

p——壓力，Pa

dp/dz——壓降梯度，Pa/m

Ql——液相體積流量，m3/s

R——管道半徑，m

Re——Reynolds數(shù)

S——濕周，m

s——標(biāo)準(zhǔn)差

ug——氣相實際流速，m/s

ul——液相實際流速，m/s

usg——氣相表觀流速，m/s

usg_c——臨界氣相表觀流速，m/s

usl——液相表觀流速，m/s

ε——絕對偏差

θ——管道軸向與水平方向的夾角，rad

μ——動力黏度，Pa·s

ρ——密度，kg/m3

τi——界面剪切應(yīng)力，Pa

τw——壁面剪切應(yīng)力，Pa

φ——液面角，即氣液界面周長對應(yīng)的半圓心角，rad

下角標(biāo)

g——氣相

i——界面

l——液相

w——壁面