頁(yè)巖氣水平井全井筒臨界攜液流量模型

周朝 何祖清 付道明 羅璇 劉歡樂 孫志揚(yáng)

1. 中國(guó)石化石油工程技術(shù)研究院；2. 頁(yè)巖油氣富集機(jī)理與有效開發(fā)國(guó)家重點(diǎn)實(shí)驗(yàn)室；3. 中國(guó)石化江漢油田石油工程技術(shù)研究院

0 引言

多段壓裂頁(yè)巖氣水平井在生產(chǎn)早期有大量的壓裂液返排到井筒，造成頁(yè)巖氣井長(zhǎng)期帶液生產(chǎn)，井筒積液?jiǎn)栴}不容忽視［1］。臨界攜液流量法是現(xiàn)場(chǎng)普遍采用的積液預(yù)測(cè)方法，然而，以Turner 模型、Coleman 模型為代表的球狀液滴模型對(duì)于國(guó)內(nèi)氣井的積液預(yù)測(cè)精度并不理想，考慮液滴變形的李閩模型對(duì)于垂直井筒的積液預(yù)測(cè)精度較高，但是并未考慮井斜變化和產(chǎn)液量對(duì)井筒積液的影響［2］。之后，王志彬等［3］、潘杰等［4］考慮液滴變形和產(chǎn)液量，建立了垂直井筒臨界攜液流量模型，Belfroid等［5］、李麗等［6］、Chen 等［7］考慮井斜角變化，建立了傾斜(水平) 井筒臨界攜液流量模型。Wang 等［8］根據(jù)室內(nèi)實(shí)驗(yàn)數(shù)據(jù)修正Belfroid 模型，得到考慮井斜和液相流速的連續(xù)油管臨界攜液相關(guān)式，但應(yīng)用條件受限。綜合分析，目前的臨界攜液流量模型均未完整反映頁(yè)巖氣井的復(fù)雜井身結(jié)構(gòu)和返排液量變化特征，為此，以液滴模型為基礎(chǔ)，綜合考慮井筒產(chǎn)液量、液滴變形和造斜率變化，根據(jù)變形液滴能量平衡關(guān)系，建立了頁(yè)巖氣井全井筒臨界攜液流量模型。

1 模型建立

1.1 物理模型

以橢球形變形液滴為基礎(chǔ)(圖1)，推導(dǎo)建立了臨界攜液流量新模型。液滴從球形變形為橢球形，以單個(gè)變形液滴為研究對(duì)象，主要進(jìn)行如下假設(shè)：

圖 1 液滴變形示意圖Fig. 1 Schematic droplet deformation

(1)球形液滴以球心所在平面發(fā)生對(duì)稱變形；

(2)球形液滴變形后為光滑的標(biāo)準(zhǔn)橢球體，橢球體的兩個(gè)長(zhǎng)軸長(zhǎng)度相等；

(3)氣流方向垂直于橢球形液滴兩個(gè)長(zhǎng)軸所在平面；

(4)忽略液滴與氣流之間的質(zhì)量交換和熱量交換。

1.2 臨界攜液流量模型

根據(jù)液滴動(dòng)力學(xué)理論，當(dāng)液滴在氣流中保持滯止?fàn)顟B(tài)時(shí)，作用在液滴上的作用力相互平衡，在垂直井筒中，液滴所受重力等于液滴所受浮力與氣體曳力之和，即

代入各個(gè)力的表達(dá)式，有

式中，F(xiàn)G為液滴重力，N；Fb為液滴所受浮力，N；Fd為氣體曳力，N；ρl為液相密度，kg/m3；Vd為液滴體積，m3；g為重力加速度，m/s2；ρg為氣相密度，kg/m3；Cd為曳力系數(shù)，無因次；Ad為液滴迎風(fēng)面積，m2；ug為氣相流速，m/s。

液滴從球形變形為橢球形后，式(2)可改寫為

式中，de為橢球形液滴長(zhǎng)軸長(zhǎng)度，m；h為橢球形液滴短軸長(zhǎng)度，m。

根據(jù)液滴攜液理論，如果氣流能將井筒中的最大液滴攜帶到井口，則氣井能夠連續(xù)攜液，此時(shí)的氣體流速稱為臨界攜液流速［9］。定義橢球形液滴的軸比α為h/de，并考慮臨界攜液時(shí)的最大穩(wěn)定橢球形液滴長(zhǎng)軸長(zhǎng)度為dec，則有

式中，ucv為垂直井筒臨界攜液流速，m/s；α為橢球形液滴的軸比，小數(shù)；dec為最大穩(wěn)定橢球形液滴長(zhǎng)軸長(zhǎng)度，m。

根據(jù)橢球形液滴軸比實(shí)驗(yàn)擬合得到的相關(guān)式為［10］

將式(5)代入式(4)整理，得到垂直井筒臨界攜液流速

在井筒傾斜段和水平段，由于造斜率和井身結(jié)構(gòu)的變化，會(huì)影響井筒流體的流動(dòng)。在造斜率發(fā)生變化的地方，液滴與井筒碰撞并反彈，氣體及其夾帶的液滴的流動(dòng)方向發(fā)生改變，以匹配新的造斜率和井身結(jié)構(gòu)。考慮傾斜段和水平段中液滴與管壁碰撞和反彈造成的能量損失，對(duì)垂直井筒臨界攜液流速進(jìn)行修正［11］，得到全井筒臨界攜液流速

式中，ucw為全井筒臨界攜液流速，m/s；R為最大造斜率(最大全角變化率)，(°)/30 m。

全井筒臨界攜液流量為

式中，qc為臨界攜液流量(標(biāo)況)，m3/d；A為管道橫截面積，m2；p為壓力，MPa；Z為氣相偏差系數(shù)，無因次；T為溫度，K。

為求取臨界攜液流量，需首先確定關(guān)鍵參數(shù)dec、Cd和表面張力σ，并選取適合的井筒兩相流模型。根據(jù)式(8)計(jì)算出頁(yè)巖氣井井筒各個(gè)深度處的臨界攜液流量后，取各臨界攜液流量的最大值作為井筒積液判斷標(biāo)準(zhǔn)。

2 關(guān)鍵參數(shù)確定和模型求解

2.1 最大穩(wěn)定橢球形液滴尺寸

在環(huán)霧流中，液滴在氣相湍流作用下向上運(yùn)動(dòng)，或者下降形成積液。湍流動(dòng)能表征了湍流的劇烈程度［12］，氣相湍流動(dòng)能與液滴表面自由能間的平衡關(guān)系決定了最大穩(wěn)定橢球形液滴尺寸。最大穩(wěn)定橢球形液滴的總表面自由能最小，達(dá)到臨界積液狀態(tài)。

2.1.1 最大穩(wěn)定橢球形液滴總表面自由能

單個(gè)最大穩(wěn)定橢球形液滴的表面積為

式中，S為最大穩(wěn)定橢球形液滴的表面積，m2。

表面自由能等于產(chǎn)生表面面積所做的功［13］，則最大穩(wěn)定橢球形液滴的總表面自由能為

式中，ESmin為最大穩(wěn)定橢球形液滴總表面自由能，W；σ為氣液表面張力，N/m；vsl為液相表觀流速，m/s。

2.1.2 氣相湍流動(dòng)能

井筒中單位體積氣相湍流動(dòng)能為

式中，eT為單位體積氣相湍流動(dòng)能，W/m3；v′r為氣相湍流平均徑向速度脈動(dòng)，m/s；v′θ為氣相湍流平均切向速度脈動(dòng)，m/s；v′z為氣相湍流平均軸向速度脈動(dòng)，m/s。

由于平均徑向速度脈動(dòng)的均方根近似等于摩擦速度，則有［14］

式中，vsg為氣相表觀流速，m/s；fsg為以氣相表觀流速流經(jīng)光滑管道的氣相摩阻系數(shù)，無因次。

在環(huán)霧流中通常為充分發(fā)展的氣相湍流流動(dòng)，可假設(shè)氣相湍流為各向同性，有v′r=v′θ=v′z，將式(12)代入式(11)，則井筒中氣相總湍流動(dòng)能為［15］

式中，ET為氣相總湍流動(dòng)能，W。

2.1.3 最大穩(wěn)定橢球形液滴長(zhǎng)軸長(zhǎng)度

當(dāng)井筒處于積液的臨界狀態(tài)時(shí)，最大穩(wěn)定橢球形液滴總表面自由能與氣相湍流動(dòng)能相互平衡，并且成正比［15］

其中

式中，Ce為比例系數(shù)，無因次；θd為井斜角，°。

將式(10)和式(13)代入式(14)，可得最大穩(wěn)定橢球形液滴長(zhǎng)軸長(zhǎng)度

氣相摩阻系數(shù)可用Blasius 公式計(jì)算［16］，對(duì)于充分發(fā)展的氣相湍流流動(dòng)，雷諾數(shù)Re＞2 100，則有

式中，μg為氣相黏度，Pa·s；D為管道內(nèi)徑，m。

對(duì)于垂直井筒，有Ce=0.75，并考慮到達(dá)到最大穩(wěn)定橢球形液滴時(shí)有vsg=ucv，將式(17)代入式(16)整理得

式中，Ql為液相流量，m3/s。

2.2 曳力系數(shù)

曳力系數(shù)的計(jì)算比較復(fù)雜，對(duì)于球形液滴的曳力系數(shù)，可采用Stokes 定律、Newton 公式和標(biāo)準(zhǔn)阻力曲線近似計(jì)算［17］。但是，對(duì)于橢球形變形液滴來說，其迎風(fēng)面積與球形液滴相比發(fā)生變化，液滴顯著變形的情況下，變形液滴的曳力能達(dá)到相同條件下剛性球體曳力的兩倍以上［18］，因此，球形液滴的曳力系數(shù)計(jì)算方法對(duì)于橢球形液滴不再適用［19］。在井筒高雷諾數(shù)條件下，考慮液滴內(nèi)部流動(dòng)，對(duì)剛性橢球體的曳力系數(shù)進(jìn)行修正，得到橢球形液滴的曳力系數(shù)［3］

2.3 表面張力

目前的臨界攜液流量模型中，常取表面張力為常數(shù)，但是，表面張力的取值受到液相含量、溫度和壓力的影響，沿井深是差異分布的，如果取表面張力為常數(shù)，會(huì)造成臨界攜液流量計(jì)算結(jié)果出現(xiàn)偏差，最大偏差可達(dá)40%［20］。頁(yè)巖氣井中的液相為返排壓裂液，Sutton 公式［21］適用于天然氣-鹽水體系，能夠反映頁(yè)巖氣井產(chǎn)出頁(yè)巖氣和返排壓裂液的生產(chǎn)特點(diǎn)，因此選用Sutton 公式計(jì)算氣-液表面張力

式中，Cs為礦化度，mg/L。

2.4 臨界攜液流量模型求解

式(7)、式(8)構(gòu)成了頁(yè)巖氣井全井筒臨界攜液流量模型，模型中的關(guān)鍵參數(shù)dec、Cd和σ分別由式(18)、式(19)、式(20)確定，其中dec和Cd反映了井筒產(chǎn)液量與液滴變形的影響。由于臨界攜液流速ucw、橢球形液滴軸比α和曳力系數(shù)Cd均與最大穩(wěn)定橢球形液滴長(zhǎng)軸長(zhǎng)度dec相關(guān)，可消去ucw，采用Newton-Raphson 迭代法求解

其中

根據(jù)式(21)求解出dec后，即可求解出ucw和qc。

2.5 井筒兩相流模型優(yōu)選

臨界攜液流量受到井筒溫壓分布的影響，為了準(zhǔn)確計(jì)算臨界攜液流量，需要優(yōu)選適合的兩相流模型。由于頁(yè)巖氣水平井井斜角變化較大，因此選取可以計(jì)算水平和傾斜井筒壓降的兩相流模型：Beggs-Brill 模 型［22］、Baker Jardine 模 型［23］、Dulker-Flanigan 模型［24-25］和Mukherjee-Brill 模型［26］。根據(jù)涪陵現(xiàn)場(chǎng)39 口頁(yè)巖氣井的實(shí)測(cè)油(套)壓數(shù)據(jù)和生產(chǎn)數(shù)據(jù)，分別利用4 種兩相流模型計(jì)算井底流壓，并與實(shí)測(cè)井底流壓結(jié)果比較。4 種兩相流模型的計(jì)算結(jié)果如圖2 所示，誤差分析如表1 所示?？梢钥闯?，Mukherjee-Brill 模型的平均相對(duì)誤差的絕對(duì)值和均方根誤差均最小，壓力計(jì)算精度最高，可用于計(jì)算頁(yè)巖氣水平井的井筒壓力分布，滿足現(xiàn)場(chǎng)壓力計(jì)算精度要求。

表 1 兩相流模型誤差分析Table 1 Error analysis of two-phase flow models

圖 2 兩相流模型壓力計(jì)算結(jié)果Fig. 2 Pressure calculation result of two-phase flow model

3 應(yīng)用實(shí)例

3.1 現(xiàn)場(chǎng)參數(shù)

根據(jù)涪陵現(xiàn)場(chǎng)39 口頁(yè)巖氣井的井筒溫壓測(cè)試數(shù)據(jù)(以最大井筒流壓梯度0.45 MPa/100 m 作為積液界限參考［27］)和生產(chǎn)動(dòng)態(tài)數(shù)據(jù)，得到氣井實(shí)際積液狀態(tài)，并選用Belfroid 模型［5］、修正李閩模型［2］、李麗模型［6］、Chen 模型［7］、Wang 模型［8］和本文新模型分別計(jì)算現(xiàn)場(chǎng)氣井的臨界攜液流量，根據(jù)積液預(yù)測(cè)精度對(duì)新模型進(jìn)行驗(yàn)證。其中的修正李閩模型是指應(yīng)用Fiedler 形狀函數(shù)修正［28］，使模型可應(yīng)用于傾斜(水平)井筒的積液預(yù)測(cè)。

Belfroid 模型和修正李閩模型具有類似的表達(dá)式

式中，uc為臨界攜液流速，m/s；C為系數(shù)，Belfroid 模型中C取6.6，修正李閩模型中C取2.5；fw為壁面摩阻系數(shù)，無因次；NB為Bond 數(shù)，無因次。

39 口頁(yè)巖氣井生產(chǎn)數(shù)據(jù)和實(shí)際積液狀態(tài)如表2所示。天然氣相對(duì)密度0.56~0.57，返排液密度1.01~1.05 g/cm3，油管內(nèi)徑為62 mm 或50.6 mm，套管內(nèi)徑為115 mm。采用Yarborough-Hall 方法［29］計(jì)算氣相偏差系數(shù)。

表 2 頁(yè)巖氣井生產(chǎn)數(shù)據(jù)和積液狀態(tài)Table 2 Field parameters of shale gas wells

3.2 積液預(yù)測(cè)結(jié)果分析

根據(jù)表2 的頁(yè)巖氣井現(xiàn)場(chǎng)參數(shù)，以及各井的井筒溫壓分布數(shù)據(jù)，分別用6 種臨界攜液流量模型計(jì)算全井筒臨界攜液流量，并取所有臨界攜液流量中的最大值與實(shí)際產(chǎn)氣量進(jìn)行比較，判斷井筒積液狀態(tài)。6 種模型的積液預(yù)測(cè)精度和預(yù)測(cè)結(jié)果如表3 和圖3 所示，圖3 中的對(duì)角線左上方區(qū)域表示氣井不積液，右下方區(qū)域表示氣井積液，而處于接近積液狀態(tài)的氣井，其數(shù)據(jù)點(diǎn)應(yīng)位于對(duì)角線±15% 偏差范圍內(nèi)［4］(兩條虛線之間的區(qū)域)。將6 種臨界攜液流量模型的積液預(yù)測(cè)結(jié)果與井筒實(shí)際積液狀態(tài)比較，得到各模型的積液預(yù)測(cè)精度。

由表3 和圖3 可知，對(duì)于頁(yè)巖氣水平井的積液預(yù)測(cè)，Belfroid 模型的預(yù)測(cè)結(jié)果偏差最大，Belfroid 模型和Wang 模型均高估了臨界攜液流量，因此這兩種模型可以準(zhǔn)確預(yù)測(cè)積液井，但是對(duì)不積液井的預(yù)測(cè)誤差較大；修正李閩模型低估了臨界攜液流量，因此可以準(zhǔn)確預(yù)測(cè)不積液井，但是對(duì)積液井的預(yù)測(cè)誤差較大；李麗模型和Chen 模型對(duì)于各種積液狀態(tài)井的預(yù)測(cè)結(jié)果均有偏差；新模型的積液預(yù)測(cè)結(jié)果與實(shí)際積液情況的符合率最高，總體積液預(yù)測(cè)精度達(dá)92.3%，新模型可以準(zhǔn)確預(yù)測(cè)積液井和接近積液井，只是對(duì)3 口不積液井的預(yù)測(cè)結(jié)果出現(xiàn)偏差，盡管如此，新模型對(duì)不積液井的預(yù)測(cè)精度仍然能夠滿足現(xiàn)場(chǎng)應(yīng)用需求。綜合來看，新模型對(duì)于頁(yè)巖氣水平井的積液預(yù)測(cè)精度最高，可以有效指導(dǎo)現(xiàn)場(chǎng)積液判斷與排采工藝選擇。

4 結(jié)論

(1)依據(jù)液滴動(dòng)力學(xué)分析和能量分析，綜合考慮井筒產(chǎn)液量、液滴變形和造斜率變化引起的液滴能量損失，根據(jù)變形液滴能量平衡關(guān)系，建立了頁(yè)巖氣井全井筒臨界攜液流量模型，彌補(bǔ)了現(xiàn)有模型沒有同時(shí)考慮頁(yè)巖氣水平井復(fù)雜井身結(jié)構(gòu)和返排液量變化的缺陷。

(2)根據(jù)變形液滴能量平衡關(guān)系，確定了最大穩(wěn)定變形液滴長(zhǎng)軸長(zhǎng)度；選取了適用于頁(yè)巖氣水平井的曳力系數(shù)公式和表面張力公式；優(yōu)選Mukherjee-Brill 兩相流模型用于計(jì)算頁(yè)巖氣水平井的井筒壓力分布，均方根誤差為8.96%，能滿足現(xiàn)場(chǎng)壓力計(jì)算精度要求。

(3) 39 口頁(yè)巖氣井的實(shí)例分析表明，與現(xiàn)有臨界攜液流量模型相比，新模型對(duì)于頁(yè)巖氣水平井的積液預(yù)測(cè)符合率最高，積液預(yù)測(cè)精度達(dá)92.3%，新模型可以有效指導(dǎo)頁(yè)巖氣水平井積液判斷。