傾斜氣井臨界攜液流速預(yù)測(cè)新模型

王武杰，崔國民，魏耀奇，潘杰

（1.上海理工大學(xué)能源與動(dòng)力工程學(xué)院，上海市動(dòng)力工程多相流動(dòng)與傳熱重點(diǎn)實(shí)驗(yàn)室，上海 200093；2.國家石油天然氣管材工程技術(shù)研究中心，陜西寶雞 721008；3.西安石油大學(xué)石油工程學(xué)院，西安 710065）

0 引言

準(zhǔn)確預(yù)測(cè)臨界攜液流速對(duì)診斷氣井積液狀況、解決含水氣井穩(wěn)產(chǎn)問題、挖掘煤炭地下氣化開發(fā)潛力至關(guān)重要[1-3]，同時(shí)，其對(duì)促進(jìn)致密油氣快速規(guī)模發(fā)展，早日實(shí)現(xiàn)中國“能源獨(dú)立”戰(zhàn)略具有積極意義[4-5]。越來越多的實(shí)驗(yàn)數(shù)據(jù)表明，傾斜氣井的攜液能力最差，預(yù)測(cè)也最為復(fù)雜和困難[6-7]。因此，有必要對(duì)傾斜氣井臨界攜液流速的預(yù)測(cè)展開深入研究。分層流作為傾斜管道氣液兩相流動(dòng)的主要流型被普遍認(rèn)可。分層流流型下，井筒內(nèi)氣液兩相的周向分布對(duì)流動(dòng)特性有至關(guān)重要的影響，尤其對(duì)液相流量較小而氣相流量較大的情形更甚[8-9]。目前對(duì)臨界攜液狀態(tài)下的氣液相分布缺乏量化研究，且尚無傾斜管道氣液兩相臨界流動(dòng)條件下的相界面摩擦因子計(jì)算模型?；诖?，本文分別考察單位管長(zhǎng)氣液兩相系統(tǒng)重力勢(shì)能、表面自由能和氣相動(dòng)能的改變量，并根據(jù)能量最小原理確定臨界攜液狀態(tài)下氣液相分布的結(jié)構(gòu)參數(shù)；同時(shí)采用實(shí)驗(yàn)數(shù)據(jù)修正相界面摩擦因子計(jì)算模型，使預(yù)測(cè)模型封閉；并將模型應(yīng)用于傾斜氣井臨界攜液流速的預(yù)測(cè)，以評(píng)估模型的應(yīng)用效果。

1 實(shí)驗(yàn)系統(tǒng)

本文通過可視化實(shí)驗(yàn)，模擬氣井?dāng)y液的動(dòng)態(tài)過程，實(shí)驗(yàn)裝置如圖1所示。

圖1 實(shí)驗(yàn)裝置圖

分別采用水和甘油、水的混合溶液（二者體積分?jǐn)?shù)均為 50%）模擬井下積液，利用空氣模擬天然氣。實(shí)驗(yàn)條件下空氣密度ρa(bǔ)=1.184 kg/m3，空氣動(dòng)力黏度μa=1.849×10-5Pa·s；水密度ρw=997.05 kg/m3，水動(dòng)力黏度μw=8.900 8×10-4Pa·s；水、甘油混合溶液密度ρgw=1 139.56 kg/m3，水、甘油混合溶液動(dòng)力黏度μgw=6.882×10-3Pa·s；空氣、水表面張力σaw=0.071 97 N/m，空氣與水、甘油混合溶液表面張力σagw=0.067 11 N/m。

分別選用外徑50 mm、壁厚5 mm和外徑70 mm、壁厚5 mm的有機(jī)玻璃管進(jìn)行實(shí)驗(yàn)，步驟為：①在管內(nèi)低洼處注入一定體積的液體，通過風(fēng)機(jī)增壓向玻璃管內(nèi)注入壓縮空氣，計(jì)量壓力和流量；②采用變頻器實(shí)現(xiàn)氣體流速的連續(xù)可調(diào)，直至低洼處的積液完全鋪展于上傾段并與管壁相對(duì)靜止（此時(shí)達(dá)臨界氣體流速）。

2 模型建立與求解方法

2.1 受力分析與控制方程

臨界氣體流速條件下，傾斜管內(nèi)流型為波狀分層流，液膜在重力、相界面和壁面切應(yīng)力的作用下處于平衡狀態(tài)。此時(shí)，近氣相側(cè)的液膜在相界面切應(yīng)力的作用下向下游運(yùn)動(dòng)，而近壁面?zhèn)鹊囊耗ぴ谥亓ψ饔孟禄亓?，即液膜?nèi)部存在循環(huán)流動(dòng)。傾斜井筒微元段內(nèi)液膜受力如圖2所示。上述流動(dòng)過程可簡(jiǎn)化為穩(wěn)態(tài)一維流動(dòng)，并忽略加速壓降，可對(duì)氣、液相分別建立動(dòng)量方程：

圖2 傾斜井筒內(nèi)液膜流動(dòng)與受力分析

流動(dòng)過程中近似認(rèn)為氣相與液相的平均壓力梯度相等，聯(lián)立（1）式和（2）式可得控制方程：

上式中，幾何參數(shù)與氣液相分布有關(guān)，應(yīng)力參數(shù)與持液率、各相速度和摩擦因子有關(guān)。通過閉合幾何參數(shù)和應(yīng)力參數(shù)可使控制方程閉合。

2.2 氣液相分布

根據(jù)Brauner等[10]的研究，液膜的周向分布如圖3所示。其中，O和O1分別表示井筒流道的截面軸心和虛擬圓圓心?？刂品匠讨械膸缀螀?shù)采用下式計(jì)算：

圖3 井筒周向的液膜分布

不同相界面假設(shè)下的幾何參數(shù)的計(jì)算式見表1。

表1 不同相界面下的幾何參數(shù)計(jì)算關(guān)系式

2.3 能量最小原理

在相界面形狀由水平變?yōu)閺澢倪^程中，重力、切應(yīng)力和表面張力將引起系統(tǒng)勢(shì)能、動(dòng)能和表面自由能的變化。臨界條件下，單位管長(zhǎng)氣液兩相系統(tǒng)總能的改變量可表示為：

基于能量最小原理，單位管長(zhǎng)系統(tǒng)總能量最小時(shí)，相界面結(jié)構(gòu)最穩(wěn)定，即：

①單位管長(zhǎng)重力勢(shì)能的改變量[10]

②單位管長(zhǎng)表面自由能的改變量

③單位管長(zhǎng)氣相動(dòng)能改變量

根據(jù)七分之一次方定律，圓管內(nèi)湍流氣相的速度分布表示如下：

單位管長(zhǎng)氣相動(dòng)能改變量采用下式計(jì)算：

2.4 剪切應(yīng)力計(jì)算

Taitel等[11]給出了相關(guān)應(yīng)力計(jì)算表達(dá)式：

Banafi等[12]實(shí)驗(yàn)研究指出 Taitel模型可準(zhǔn)確預(yù)測(cè)氣液兩相低持液率流動(dòng)條件下的氣壁面摩擦因子，表達(dá)式如下：

層流時(shí)，Cg=16，m=1；湍流時(shí)，Cg=0.046，m=0.2。

Biberg[13]指出，傾斜管道光滑分層流條件下，平均液壁面切應(yīng)力與平均相界面切應(yīng)力滿足下式：

上式中f（φ0P）是 難以獲得解 析 解 的積分項(xiàng)。Biberg[13]給出以下近似計(jì)算關(guān)系式：

式中，an和bn的取值見表 2[13]。該近似值的最大絕對(duì)誤差不超過2×10-6。

表2 近似關(guān)系式（（26）式）中參數(shù)an和bn的取值

準(zhǔn)確計(jì)算相界面摩擦因子是預(yù)測(cè)相界面切應(yīng)力的關(guān)鍵。傾斜管道中，界面波和持液率對(duì)相界面摩擦因子的影響較大。界面波動(dòng)特性主要取決于機(jī)械能從氣相向液相的傳遞速率以及液體中黏性耗散速率，而這一過程受傾角的影響[14]。隨著管道傾角的增加，重力對(duì)界面波生長(zhǎng)的約束減弱，使得傾斜管中的液層表面更易于觀察到界面波。本文認(rèn)為在傾斜管道中，重力對(duì)傾斜管內(nèi)臨界狀態(tài)下的界面波動(dòng)作用更加顯著。因此，一方面將界面波等效為管壁的相對(duì)粗糙度，另一方面，引入傾角對(duì)相界面摩擦因子進(jìn)行修正[14]，具體過程如下：

其中，無因次相界面濕周定義為：

再考慮傾角的影響對(duì)其進(jìn)行修正：

式中，A和B與流體物性和管道條件有關(guān)[15]。

2.5 模型求解方法

控制方程（（3）式）可以視為持液率的函數(shù)，即：f（ε）=0。零始終是方程f（ε）=0的一個(gè)解。的函數(shù)圖像如圖4所示。當(dāng)usg=usg,cr時(shí)，f（ε）=0只有一個(gè)非零解，且該非零解在（0，0.1）區(qū)間求得，求解步驟如圖5所示。

圖 4 函數(shù) f（ε） 1/3在不同流動(dòng)條件下的變化曲線

圖5 模型求解過程示意圖

3 結(jié)果分析與實(shí)例驗(yàn)證

3.1 氣液相分布預(yù)測(cè)結(jié)果與分析

對(duì)控制方程無因次化處理后可得：

其中

由上式可知，臨界狀態(tài)下的氣液相分布與邦德數(shù)、氣相表觀韋伯?dāng)?shù)、持液率和井斜角有關(guān)?；谀芰孔钚≡?，量化氣液兩相系統(tǒng)總能量隨相界面結(jié)構(gòu)的變化，并以空氣和水為例，分別考察上述 4個(gè)因素對(duì)管內(nèi)氣液相分布的影響，計(jì)算結(jié)果見圖6（圖中紅色數(shù)據(jù)點(diǎn)為各工況下的能量最小值點(diǎn)）?？梢钥吹剑孩匐S著管徑的增加，最穩(wěn)定的氣液相分布（對(duì)應(yīng)能量最小值）所對(duì)應(yīng)的相界面彎曲角從125°（微凸）逐漸趨于180°（水平）；反之，管徑越小，在潤(rùn)濕性和表面張力的影響下相界面彎曲角越小，呈現(xiàn)凸型結(jié)構(gòu)的液相分布（見圖6a）；②隨著井斜角的增加，重力抵消潤(rùn)濕性和表面張力的作用效果更加顯著，使最穩(wěn)定的氣液相分布所對(duì)應(yīng)的相界面彎曲角從120°（微凸）趨近于180°（水平）（見圖 6b）；③結(jié)合圖 6c、圖 6d可知，氣體流速和氣相密度增加，使氣相動(dòng)能也增加，在潤(rùn)濕性和表面張力的影響下相界面彎曲角減?。ㄗ兺梗┯欣诮档拖到y(tǒng)表面自由能；④隨著持液率的增加，液相重力勢(shì)能的影響增大，有利于穩(wěn)定原來水平相界面的形狀，同時(shí)界面微凸也可減小持液率增加導(dǎo)致的系統(tǒng)表面自由能的增加。但是在液相密度遠(yuǎn)大于氣相密度的低壓系統(tǒng)下，液相重力勢(shì)能起主導(dǎo)作用（見圖6e），隨著持液率的增加，最穩(wěn)定的氣液相分布所對(duì)應(yīng)的相界面彎曲角從105°逐漸趨于165°（近水平）。綜上所述，低持液率流動(dòng)條件下，受表面張力和潤(rùn)濕性的影響，相界面呈現(xiàn)微凸，符合如圖7所示的實(shí)驗(yàn)觀察結(jié)果。

圖6 液相分布的模型預(yù)測(cè)結(jié)果

圖7 實(shí)驗(yàn)觀察的臨界狀態(tài)下的氣液相分布示意圖

3.2 實(shí)驗(yàn)結(jié)果與分析

圖 8為臨界攜液氣體流速和臨界壓力梯度的實(shí)驗(yàn)結(jié)果。由圖可知，臨界攜液流速和臨界壓力梯度隨著傾角的增加呈非線性增加。在相同的管道傾角下，液相密度或者黏度越大，臨界攜液流速越大。這與Birvalski等[16]的實(shí)驗(yàn)結(jié)果一致。同時(shí)，管徑增加使臨界攜液氣體流速增加。這是因?yàn)樵谙嗤忠郝氏?，界面形狀的變化不及Bo和Wesg數(shù)對(duì)臨界攜液流速的影響顯著，要使控制方程依然成立，必須提高氣體流速，以抵消管徑增加引起的Bo/Wesg的增量。

圖8 臨界攜液流速和臨界壓力梯度的實(shí)驗(yàn)結(jié)果

3.3 傾斜管內(nèi)相界面摩擦因子計(jì)算模型

采用Birvalski等[16]和本文的實(shí)驗(yàn)數(shù)據(jù)，基于簡(jiǎn)化后的水平相界面形狀封閉幾何參數(shù)，對(duì)不同傾角下的相界面摩擦因子修正參數(shù)進(jìn)行反算，并通過數(shù)據(jù)擬合得到相界面摩擦因子修正參數(shù)（f（α）），計(jì)算結(jié)果如表3所示。采用表 3中相界面摩擦因子修正參數(shù)對(duì)相界面摩擦因子計(jì)算公式（（31）式）進(jìn)行修正后，模型的預(yù)測(cè)值達(dá)到了非常高的精度。皮爾遜相關(guān)系數(shù)相對(duì)較高。

表3 不同實(shí)驗(yàn)條件下的相界面摩擦因子修正參數(shù)

由表3可知，f（α）的斜率為正值，即表明相界面摩擦因子隨著傾角的增大整體呈現(xiàn)遞增的趨勢(shì)。主要是因?yàn)閮A角的增加使得約束界面波動(dòng)的重力分量減小，導(dǎo)致界面波動(dòng)更加劇烈和明顯。同時(shí)，f（α）小于1，表明在微傾管道內(nèi)低持液率流動(dòng)條件下，相界面摩擦因子要小于相同氣體流速下水平管道相界面摩擦因子預(yù)測(cè)值，主要原因是當(dāng)管道接近水平時(shí)，傾角對(duì)界面波的影響不如持液率的影響顯著。本文基于實(shí)驗(yàn)數(shù)據(jù)，結(jié)合無因次化分析，給出相界面摩擦因子修正參數(shù)的計(jì)算式：

3.4 模型對(duì)比與工程實(shí)例驗(yàn)證

將模型預(yù)測(cè)結(jié)果與實(shí)驗(yàn)數(shù)據(jù)進(jìn)行對(duì)比（見圖9）。本文模型預(yù)測(cè)的臨界攜液流速隨管徑的增加而增加，隨液相密度和黏度的增加而增加（同管徑下，甘油溶液的臨界攜液流速大于水的臨界攜液流速），這一規(guī)律與符合 Birvalski等[16]和 Rastogi等[17]的實(shí)驗(yàn)結(jié)果；臨界壓力梯度的預(yù)測(cè)結(jié)果隨管徑的增加而減小，符合實(shí)驗(yàn)結(jié)果，管徑的增加使液相水力直徑增加是壓力梯度減小的重要原因。對(duì)比本文實(shí)驗(yàn)數(shù)據(jù)，新模型預(yù)測(cè)臨界攜液流速和臨界壓力梯度的平均相對(duì)誤差分別為1.19%和3.02%。同時(shí)，臨界攜液流速隨傾角的增加先增大后減小，最大臨界攜液流速對(duì)應(yīng)的管道傾角隨管徑的增加而增加，同時(shí)也受流動(dòng)介質(zhì)物性的影響。本文模型模擬計(jì)算發(fā)現(xiàn)，標(biāo)況下天然氣攜礦化水時(shí)，最大臨界攜液流速對(duì)應(yīng)的管道傾角分別為 30°（D=50.67 mm）、37°（D=62.00 mm）、40°（D=75.9 mm）和 45°（D=100.53 mm）。

圖9 本文模型預(yù)測(cè)結(jié)果與實(shí)驗(yàn)結(jié)果對(duì)比

本文模型應(yīng)用于高壓系統(tǒng)時(shí)需要對(duì)預(yù)測(cè)結(jié)果進(jìn)行轉(zhuǎn)換。本文從能量角度出發(fā)，認(rèn)為高、低壓管道系統(tǒng)內(nèi)，臨界條件下單位體積氣液兩相間傳遞的機(jī)械能應(yīng)相等，即：

為便于比較，將臨界攜液流速轉(zhuǎn)換為標(biāo)況下的臨界攜液流量：

當(dāng)天然氣的實(shí)際流速大于傾斜氣井中臨界攜液流速最大值時(shí)，氣井才能連續(xù)攜液。氣田現(xiàn)場(chǎng)需結(jié)合方位測(cè)井資料來確定傾斜氣井的最大臨界攜液流速。選取中國普光氣田和延長(zhǎng)氣田42口氣井的生產(chǎn)數(shù)據(jù)，并與實(shí)際產(chǎn)氣量繪制交會(huì)圖（見圖10），可以看到本文模型僅誤判1口氣井，誤判率2.38%，說明本文模型的預(yù)測(cè)結(jié)果具有較高的可信度，可對(duì)傾斜氣井積液進(jìn)行有效判斷。

圖10 本文模型預(yù)測(cè)結(jié)果與對(duì)比

4 結(jié)論

低持液率流動(dòng)條件下，相界面形狀更容易受到壁面潤(rùn)濕性和表面張力的影響呈現(xiàn)彎曲，管徑越小、持液率越小、井斜角越小、氣體流速越大、氣相密度越大時(shí)，相界面彎曲越明顯。

氣井臨界攜液流速隨井斜角的增加呈現(xiàn)非線性先增后減的趨勢(shì)；最大臨界攜液流速對(duì)應(yīng)的井斜角隨井筒直徑的增加而增大，同時(shí)也受氣液兩相物性的影響。

本文模型預(yù)測(cè)臨界攜液流速與臨界壓力梯度的平均相對(duì)誤差分別為1.19%和3.02%，現(xiàn)場(chǎng)傾斜氣井積液誤判率2.38%，可對(duì)傾斜氣井積液進(jìn)行有效判斷。

符號(hào)注釋：

an，bn，A，B——常系數(shù)，無因次；AT——油管截面面積，m2；Ag——?dú)庀嗔鲃?dòng)通道面積，m2；ΔAgl——單位管長(zhǎng)氣液相接觸面面積的變化，m2/m；Al——液相流動(dòng)通道面積，m2；ΔAlw——單位管長(zhǎng)液壁接觸面面積的變化，m2/m；Bo——邦德數(shù)，無因次；Cg，m——與流態(tài)有關(guān)的系數(shù)，無因次；D——油管內(nèi)徑，m；DH,g——?dú)庀嗨χ睆剑琺；Dl（φ0P）——不同液相分布角下的等效層流直徑，m；fg——?dú)庀嗄Σ烈蜃?，無因次；fi——相界面摩擦因子，無因次；fl——液相摩擦因子，無因次；f（α）——相界面摩擦因子修正參數(shù)，無因次；f（φ0P）——與液相分布角有關(guān)的函數(shù)，無因次；F——中間變量，無因次；g——重力加速度，m/s2；k——相界面絕對(duì)粗糙度，m；n——常系數(shù)編號(hào)；p——壓力，Pa；Qcr——臨界攜液流量，m3/d；r——距離井筒軸心的徑向位置，m；R——油管半徑，m；Reg——?dú)庀嗬字Z數(shù)，無因次；R1——虛擬圓半徑，m；Sg——?dú)庀酀裰?，m；Si——相界面濕周，m；si——無因次相界面濕周，無因次；Sl——液相濕周，m；T——熱力學(xué)溫度，K；ug——?dú)庀嗾鎸?shí)流速，m/s；g——平均氣相真實(shí)速度，m/s；ug,max——最大氣相真實(shí)速度，m/s；ui——相界面真實(shí)流速，m/s；ul——液相真實(shí)流速，m/s；usg——?dú)庀啾碛^流速，m/s；Δusg——?dú)庀啾碛^流速迭代步長(zhǎng)，m/s；usg,cr——臨界攜液氣相表觀流速，m/s；usg,cr,hp——高壓系統(tǒng)下的臨界攜液氣相表觀流速，m/s；Wesg——?dú)庀啾碛^韋伯?dāng)?shù)；x——流動(dòng)方向坐標(biāo)，m；YGPC——彎曲相界面的重心位置，m；YGP——水平相界面的重心位置，m；Z——壓縮因子；α——管道傾角，（°）；β——井斜角，（°）；ΔEk,g——單位管長(zhǎng)氣相的動(dòng)量改變量，J/m；ΔEp——單位管長(zhǎng)氣液兩相系統(tǒng)的重力勢(shì)能改變量，J/m；ΔEs——單位管長(zhǎng)氣液兩相系統(tǒng)的表面自由能改變量，J/m；ΔEtot——單位管長(zhǎng)氣液兩相系統(tǒng)的總能量改變量，J/m；δl——液膜厚度，m；δl（φ）——液膜厚度分布，m；ε——持液率；ε0——持液率迭代初值；Δε——持液率迭代步長(zhǎng)；θ——潤(rùn)濕角，rad；μa——空氣動(dòng)力黏度，Pa·s；μg——?dú)庀鄤?dòng)力黏度，Pa·s；μgw——水、甘油混合溶液動(dòng)力黏度，Pa·s；μl——液相動(dòng)力黏度，Pa·s；μw——水動(dòng)力黏度，Pa·s；ρa(bǔ)——空氣密度，kg/m3；ρg——?dú)庀嗝芏?，kg/m3；ρg,hp——高壓系統(tǒng)下的氣相密度，kg/m3；ρgw——水、甘油混合溶液密度，kg/m3；ρl——液相密度，kg/m3；ρw——水密度，kg/m3；σaw——空氣、水界面的表面張力，N/m；σagw——空氣與水、甘油混合溶液的表面張力，N/m；σgl——?dú)庖航缑娴谋砻鎻埩Γ琋/m；τi——相界面剪切應(yīng)力，Pa；——平均相界面切應(yīng)力，Pa；τig——相界面反作用于氣相的相界面切應(yīng)力，Pa；τgi——?dú)庀嘧饔糜谙嘟缑娴南嘟缑媲袘?yīng)力，Pa；τwg——?dú)獗诿婕羟袘?yīng)力，Pa；τwl——液壁面剪切應(yīng)力，Pa；τwl——平均液壁面切應(yīng)力，Pa；φ——任意液相分布角，rad；φ0——彎曲相界面的液相分布角，rad；φ0P——水平相界面的液相分布角，rad；φPA——相界面彎曲角，rad。下標(biāo)：g—?dú)庀?；l—液相。