有限導(dǎo)流壓裂定向井耦合流動模型

賈 品， 程林松， 黃世軍， 方思東

（中國石油大學(xué)（北京）石油工程學(xué)院，北京 102249）

有限導(dǎo)流壓裂定向井耦合流動模型

賈 品， 程林松， 黃世軍， 方思東

（中國石油大學(xué)（北京）石油工程學(xué)院，北京 102249）

考慮人工裂縫、斜井筒和地層中的耦合流動，將壓裂定向井的流動分為地層向裂縫滲流、裂縫內(nèi)流動和斜井筒變質(zhì)量管流.采用鏡像反映和勢疊加原理，建立地層中勢分布和滲流數(shù)學(xué)模型.應(yīng)用邊值理論，將人工裂縫面離散為裂縫網(wǎng)格，并與斜井筒變質(zhì)量管流模型進(jìn)行耦合，建立有限導(dǎo)流壓裂定向井耦合流動數(shù)學(xué)模型，并形成相應(yīng)的迭代算法.應(yīng)用分析表明：該模型與Prats圖版法相比，對于垂直縫吻合度較高；裂縫導(dǎo)流能力和井斜角對產(chǎn)量及裂縫面壓力分布影響較大.

壓裂定向井；有限導(dǎo)流；耦合模型；網(wǎng)格離散；產(chǎn)能

0 引言

隨著灘海低滲透油田的大規(guī)模開發(fā)，為節(jié)約鉆井成本，鉆井平臺大多建在陸地上.開發(fā)此類油藏普遍采用定向井，為獲得較高產(chǎn)量，定向井通常進(jìn)行壓裂［1－3］.壓裂定向井對于開發(fā)灘海油藏、多層系油藏和鉆井平臺不易建設(shè)的油田尤其有效.因此，研究壓裂定向井的產(chǎn)能亟待進(jìn)行.從研究現(xiàn)狀看，國內(nèi)外學(xué)者［4－6］只有對壓裂定向井的壓力動態(tài)做過研究，對裂縫有限導(dǎo)流時的產(chǎn)能研究幾乎為空白.

基于前人對無限導(dǎo)流壓裂定向井壓力動態(tài)的研究［4－5］，將流動分為油層向裂縫滲流，裂縫內(nèi)流動和井筒變質(zhì)量管流.利用疊加原理、微源匯、邊值理論和管流理論，同時考慮裂縫有限導(dǎo)流和井筒壓力損耗，推導(dǎo)壓裂定向井滲流耦合模型，建立壓裂定向井的產(chǎn)能公式，可以預(yù)測壓裂定向井的產(chǎn)能和壓力分布.

1 壓裂定向井及裂縫三維空間描述

將傾斜裂縫面離散為若干裂縫網(wǎng)格，并將井筒劃分為若干井筒微元段.各裂縫網(wǎng)格有不同的網(wǎng)格壓力和油層流量，各井筒微元有不同的流壓和管流流量.

油層為上下封閉，均質(zhì)等厚水平無限大油藏，油層厚度為h.假設(shè)定向井（見圖1）沿X-Z平面，井斜角為θ，在油層中的長度為Ls，將其劃分為M個微元段（見圖2），每個井筒微元長度為ΔLs．壓裂裂縫沿著Y方向延伸，縱向上完全穿透油層，裂縫半長為Lf，寬度為w，劃分整個裂縫面為2N ×M個裂縫網(wǎng)格（見圖3），裂縫網(wǎng)格沿裂縫延伸方向（Y方向）長度為ΔLf，沿井筒方向長度為ΔLs，與井筒微元長度一致.

設(shè)傾斜井筒上端點的坐標(biāo)M0（x0，y0，z0），式1給出了井筒上第i微元段上任意一點坐標(biāo)

圖1 壓裂定向井示意圖Fig.1 Fractured deviated well

其中，1≤i≤M，t取值為0≤t≤1表示在第i微元段上的任意一點，例

圖2 離散井筒示意圖Fig.2 Discreted wellbore

圖3 裂縫面網(wǎng)格示意圖Fig.3 Fracture grid

如t＝0或t＝1表示第i微元段上的兩個端點．由于井筒半徑相對于裂縫半長很小，因此在計算裂縫面上任意一點坐標(biāo)時，井筒半徑可以忽略不計，則裂縫面上沿井筒方向第i網(wǎng)格，沿裂縫延伸方向第j網(wǎng)格上任意一點坐標(biāo)

其中，1≤i≤M，1≤j≤2N；0≤t≤1，0≤s≤1表示在第（i，j）裂縫網(wǎng)格上X和Y方向上的任意一點.

2 耦合流動數(shù)學(xué)模型

2.1 油層向裂縫流動

根據(jù)劉想平［8］，黃世軍［9］等人的研究，油層中任意一長為L，產(chǎn)量為Q的線匯在無限大油層M（x，y，z）點產(chǎn)生的勢分布，將此結(jié)論引入到本文中．基于鏡像反映原理，將上下封閉水平無限大油層中裂縫面的第（i，j）裂縫網(wǎng)格映射成無界地層中的傾斜板源，其z方向中點坐標(biāo)Z為

其中zm，（i，j）為第（i，j）裂縫網(wǎng)格z方向的中點坐標(biāo).

設(shè)裂縫面的第（i，j）裂縫網(wǎng)格中與井筒平行的某條線匯（XZ平面上的線匯）y方向坐標(biāo)為Y．利用疊加原理，可得該線匯在上下封閉油藏中M（x，y，z）點產(chǎn)生的勢

等式中qr，（i，j）為油層到第（i，j）裂縫網(wǎng)格的流量，r1、r2分別為線匯的兩個端點到M（x，y，z）的距離．εi，j定義為εi，j（Z，x，y，z；Y）＝ln［（r1＋r2＋ΔLs）／（r1＋r2－ΔLs）］．

對與井筒平行的線匯（XZ平面線匯）沿著裂縫延伸方向（Y軸正方向）積分，即對等式（3）的變量Y沿y方向積分，可得裂縫面的第（i，j）整個裂縫網(wǎng)格在上下封閉油層中M（x，y，z）點產(chǎn)生的勢

積分上下限為y2＝y(tǒng)0＋（j－N）ΔLf，y1＝y(tǒng)0＋（j－1－N）ΔLf．由于以上積分形式復(fù)雜，本文采用數(shù)值積分進(jìn)行計算.

將整個裂縫面共劃成了2N×M個裂縫網(wǎng)格．從左到右，從上到下依次對裂縫網(wǎng)格編號，序號為I＝j(luò)＋（i－1）2N．設(shè)供給邊界處的勢為Φe，則有

設(shè)供給邊界處的壓力為pe，第I裂縫網(wǎng)格的壓力為pfI，油層流入第I裂縫網(wǎng)格的地層流量為qrI，則每個裂縫網(wǎng)格壓力和油層流量的關(guān)系矩陣為

其中φi，j表示第J裂縫網(wǎng)格在第I裂縫網(wǎng)格中心處產(chǎn)生的勢，φe，J表示第J裂縫網(wǎng)格在油層邊界處產(chǎn)生的勢，Km為基質(zhì)滲透率，μ為流體粘度.由以上流動方程可以看出，給定裂縫網(wǎng)格壓力，可以得到油層到裂縫網(wǎng)格的流量，反之可得到裂縫網(wǎng)格壓力.

2.2 裂縫有限導(dǎo)流

裂縫面流動為穩(wěn)態(tài)流動，即滲流過程為邊值問題.以裂縫網(wǎng)格為研究對象，則流入流出該網(wǎng)格的流量包括兩部分（見圖4）：①相鄰網(wǎng)格的流入流出流量；②油層到該網(wǎng)格的流量.

利用微元法對第（i，j）裂縫網(wǎng)格的控制方程進(jìn)行推導(dǎo)．裂縫網(wǎng)格為均勻網(wǎng)格，i方向網(wǎng)格大小為Δx，j方向網(wǎng)格大小為Δy，油層流入裂縫網(wǎng)格的流量為qr，（i，j），裂縫寬度為w．根據(jù)物質(zhì)平衡原理及達(dá)西滲流方程，以裂縫網(wǎng)格為微元，則有以下離散形式的控制方程

圖4 裂縫網(wǎng)格流量示意圖Fig.4 Scheme of flux to fracture grid

其中pf，（i，j）為第（i，j）裂縫網(wǎng)格壓力，pwf，i為第i井筒微元流壓，Kf·w為裂縫導(dǎo)流能力.

因此對于裂縫面流動，網(wǎng)格壓力可由等式（7）和等式（8）計算.裂縫面共有2NM個網(wǎng)格，則共有2NM個流動方程，若給定油層到裂縫網(wǎng)格的流量qr，（i，j）（1≤i≤M，1≤j≤2N）和井筒微元壓力pwf，i（1≤i≤M），通過求解2N×M個線性方程組，可以得到每個裂縫網(wǎng)格壓力pf，（i，j）（1≤i≤M，1≤j≤2N）.

2.3 井筒變質(zhì)量管流

壓裂定向井的裂縫完全穿透油層時，油井生產(chǎn)過程中沿定向井井筒各處均有流體從裂縫流入井筒.從定向井底部到頂部，井筒內(nèi)流體質(zhì)量不斷增加，井筒流動為變質(zhì)量管流流動.國內(nèi)外學(xué)者［10－12］對變質(zhì)量流動做了大量研究，模型計算結(jié)果準(zhǔn)確，且利用應(yīng)用.引用于樂香等人［12］的結(jié)論，應(yīng)用到壓裂定向井微元井筒中（見圖5），可以得到壓裂定向井沿井筒的壓力損失計算模型：

與井筒微元相鄰的裂縫網(wǎng)格，控制方程為

等式（9）右端第一項為摩擦阻力項，第二項為重力項，第三項為上下游端面壓力項.其中，Δpwf，i為第i井筒微元與第i＋1井筒微元井底流壓壓差；φ為油層孔隙度；f1為生產(chǎn)井筒管壁摩擦系數(shù)；f2為由于從裂縫網(wǎng)格流入生產(chǎn)井筒的流體所造成的摩擦系數(shù)；qwl，i為第i井筒微元的管流流量；qwr，i為與第i井筒微元相鄰裂縫網(wǎng)格流入該井筒微元的流量；ρ為流體密度；D為井筒直徑，i的取值范圍為1≤i≤M－1．

其中qwr，i可以由以下等式計算

圖5 井筒微元及與其相鄰裂縫網(wǎng)格示意圖Fig.5 Wellbore grid and neighboring fracture grid

等式（10）右邊第一項為與第i井筒微元左端鄰的裂縫網(wǎng)格流入井筒微元的流量，第二項為與右端相鄰的裂縫網(wǎng)格流入井筒微元的流量，i的取值范圍為1≤i≤M．

3 耦合流動計算模型

若壓裂定向井為定井底流壓生產(chǎn)，則pwf，1為已知量，即油井生產(chǎn)時的定井底流壓.假設(shè)計算出了裂縫網(wǎng)格壓力和井筒微元壓力，則由等式（10）和等式（11）可得出裂縫網(wǎng)格流入每個井筒微元的流量qwr和管流流量qwl．剩下未知量Δpwf，i和pwf，i共有2M－1個，聯(lián)立等式（9）和等式（12）共有2M－2個方程，加上pwf，1為已知量，共有2M－1個方程組，則2M－1個未知量可以求解，pwf，i（2≤i≤M）可以計算出來.

生產(chǎn)井井筒微元管流流量qwl和井底流壓pwf滿足以下約束關(guān)系式

三部分流動模型有很強(qiáng)的耦合性，每一部分流動模型的計算都需要其他流動模型的計算結(jié)果.因此，為簡化計算過程，并提高模型收斂性和計算效率，本文基于以下的耦合方法和迭代算法求解生產(chǎn)井產(chǎn)量和油層壓力分布，對于矩陣方程，采用預(yù)處理共軛梯度法求解.具體計算流程圖如下（見圖6）.

具體耦合計算步驟如下

步驟1：井筒微元流壓pwf，i（1≤i≤M）賦初值；

步驟2：裂縫網(wǎng)格壓力pf，（i，j）（1≤i≤M，1≤j≤2N）賦初值；

步驟3：由等式（6）計算油層流入裂縫網(wǎng)格流量qr，（i，j）（1 ≤i≤M，1≤j≤2N）；

步驟4：由等式（7）和（8）計算裂縫網(wǎng)格壓力pf，（i，j）（1≤i ≤M，1≤j≤2N）；

步驟5：判斷步驟4計算的裂縫網(wǎng)格壓力與步驟2的網(wǎng)格壓力初值是否達(dá)到收斂條件.若收斂進(jìn)行步驟6，否則返回步驟2；

步驟6：由等式（9）計算井筒微元流壓pwf，i（1≤i≤M）；

圖6 耦合計算流程圖Fig.6 Flowchart of coupling calculation

步驟7：判斷步驟6計算出的井筒微元流壓與步驟1的井筒微元流壓初值是否達(dá)到收斂條件.若收斂則計算結(jié)束，否則返回步驟1，繼續(xù)計算，直到流壓pwf，i收斂.

其中，步驟5和步驟7的收斂條件分別為以下形式：

步驟5收斂條件

4 應(yīng)用分析

4.1 模型檢驗

為檢驗本文模型的準(zhǔn)確性，根據(jù)我國東部某灘海油田壓裂定向井的地層參數(shù)和完井?dāng)?shù)據(jù)（見表1），利用耦合流動模型進(jìn)行計算，與Prats圖版法［13］進(jìn)行對比.Prats于1961年研究了壓裂直井穩(wěn)態(tài)滲流問題，給出了裂縫不同無因次導(dǎo)流能下的有效井筒半徑，此方法只適用于計算垂直裂縫產(chǎn)能，而不能用來計算傾斜裂縫產(chǎn)能.因此為與Prats方法對比，設(shè)模型中的井斜角為0°.表1給出了主要計算參數(shù)，對比結(jié)果見圖7.

表1 地層及完井參數(shù)Table 1 Reservoir and completeion parameters

從圖7中可以看出，本文井筒無限導(dǎo)流時的產(chǎn)量與Prats圖版法產(chǎn)量相差很小，井筒無限導(dǎo)流產(chǎn)量與Prats圖版法產(chǎn)量幾乎重合，吻合度較高.

4.2 裂縫導(dǎo)流能力的影響

裂縫為有限導(dǎo)流時，每個裂縫網(wǎng)格的壓力Pf和油層到網(wǎng)格的流量qr不同.取壓裂定向井井斜角為60°，不同裂縫半長，不同裂縫導(dǎo)流能力下的單井產(chǎn)量和裂縫面壓力分布見圖8.

油田現(xiàn)場壓裂時，人工裂縫導(dǎo)流能力一般在300 mD·m～3 000 mD·m.從圖8中可以看出，在此范圍內(nèi)產(chǎn)量隨裂縫導(dǎo)流能力不斷增加；當(dāng)裂縫導(dǎo)流能力接近10 000 mD·m，產(chǎn)量不再由明顯變化，裂縫接近無限導(dǎo)流.裂縫半長為120 m時，產(chǎn)量隨著裂縫導(dǎo)流能力增加從12 m3·d－1增加到20 m3·d－1，說明裂縫導(dǎo)流能力對產(chǎn)量影響較大.

在不同裂縫半長下，裂縫導(dǎo)流能力的增加對單井產(chǎn)量的貢獻(xiàn)不同.裂縫半長為200 m時，產(chǎn)量最大增加11 m3·d－1；而裂縫半長為40 m時，產(chǎn)量最大增加3 m3·d－1，因此在定向井壓裂時，對于裂縫半長較大的生產(chǎn)井，增加裂縫導(dǎo)流能力對提高單井產(chǎn)能尤為重要.

圖9為裂縫傾角60°，半長80 m，導(dǎo)流能力300 mD·m，裂縫面壓力分布，圖中縱坐標(biāo)為生產(chǎn)井井筒方向，橫坐標(biāo)為裂縫延伸方向.

圖7 本文模型與Prats圖版法對比Fig.7 Results of coupled model and Prats method

圖8 不同裂縫半長不同裂縫導(dǎo)流能力下產(chǎn)量Fig.8 Productivity with different fracture length and conductivity

圖9 裂縫面壓力分布Fig.9 Fracture pressure distribution

根據(jù)圖9，裂縫導(dǎo)流能力較低時，裂縫面沿橫軸端部與井筒（橫軸中間部位）壓差較大，為2 MPa，井筒處上下兩端也存在壓差，最大為0.2 MPa，整個裂縫面壓差較大.

4.3 井斜角和裂縫半長的影響

從圖10可以得出，裂縫半長對壓裂定向井的產(chǎn)能有較大影響，而當(dāng)井斜角較大時，裂縫面與油層的接觸面積增加，其對產(chǎn)能也有一定影響.因此，以下對井斜角和裂縫半長對產(chǎn)能的影響做以討論.

裂縫導(dǎo)流能力為2 000 mD·m，井斜角為20°～80°，裂縫半長為40 m～200 m，不同裂縫半長和不同井斜角下的產(chǎn)量對比如下圖.

從圖10中不難得出，同一裂縫半長下，單井產(chǎn)量隨著井斜角的增加而增大，當(dāng)井斜角小于60°時，產(chǎn)量增加緩慢，井斜角大于60°時，產(chǎn)量增加很快，以裂縫半長40 m為例，井斜角小于60°時，產(chǎn)量增加1 m3·d－1，井斜角大于60°時，產(chǎn)量增加5 m3·d－1.井斜角增加，相同裂縫半長下，裂縫與油藏的接觸面積不斷增加，近井地帶滲流阻力減小，產(chǎn)量增加.不同裂縫半長下，井斜角的增加對單井產(chǎn)量的貢獻(xiàn)幾乎相同.因此，在定向井開發(fā)油藏時，較大的井斜角可以顯著提高生產(chǎn)井的產(chǎn)油能力.

圖11中得到，井斜角一定時，裂縫半長越大，生產(chǎn)井產(chǎn)量越大，壓裂定向井的產(chǎn)量與裂縫半長呈近似直線關(guān)系，說明由于裂縫半長增加而增大的等效井徑與裂縫半長為指數(shù)關(guān)系.同樣，不同井斜角下裂縫半長的增加對生產(chǎn)井的產(chǎn)量貢獻(xiàn)幾乎相同.

圖10 不同井斜角下產(chǎn)量Fig.10 Productivity with different inclined angle

圖11 不同裂縫半長下產(chǎn)量Fig.11 Productivity with different fracture length

5 結(jié)論

1）井筒存在壓力損失時的產(chǎn)量比井筒無限導(dǎo)時產(chǎn)量差距較大，且沿井筒產(chǎn)量呈上大下小分布.

2）裂縫半長較大時，裂縫導(dǎo)流能力的增大會顯著增加生產(chǎn)井的產(chǎn)量；裂縫半長較小時，產(chǎn)量增加不明顯；對于裂縫半長較大的裂縫，提高裂縫導(dǎo)流能力尤為重要.無論裂縫半長大小，生產(chǎn)井產(chǎn)量隨裂縫導(dǎo)流能力的增加有最大值，達(dá)到最大值后，產(chǎn)量不再增加.

3）井斜角小于60°時，井產(chǎn)量隨井斜角產(chǎn)量增大而緩慢增加；井斜角大于60°時產(chǎn)量增加變快.單井產(chǎn)量隨裂縫半長的增加而增大，且與裂縫半長呈直線關(guān)系.

［1］Niu Zengqian，Sui Xiangyun，Zhang Ping，et al.Research of fracturing technics for high angle deviated hole［J］.Oil Drilling ＆Production Technology，2005，27（5）：61－63.

［2］Yan Jinchuan，Diao Su，Zhu Liping，et al.Optimization and field application of hydraulic fracturing design of directional wells ［J］.Petroleum Geology and Recovery Efficiency，2008，15（5）：102－104.

［3］Diao Su，Yan Jinchuan，Ren Shan，et al.Research and field application on hydraulic fracturing technology of directional wells in western region of Sichuan［J］.Journal of Southwest Petroleum University（Science＆Technology Edition），2009，31（1）：111－115.

［4］Cinco H，Ramey H J，Miller F G.Unsteady-state pressure distribution created by a well with an inclined fracture［J］.SPE 5591，1975.

［5］Habte A D，Dinh A V，Tiab D.Pressure analysis of a well with an inclined asymmetric hydraulic fracture using type curves ［J］.SPE 140638，2010.

［6］Dinh A V，Tiab D.Pressure analysis of well with an inclined hydraulic fracture［J］.SPE 120540，2008.

［7］Zhang Yi.Research on productivity of horizontal and high deviated hydraulic fractured well［D］.Chengdu：Southwest Petroleum University，2003.

［8］Liu Xiangping，Zhang Zhaoshun，Liu Xianger，et al.A model to calculate pressure drops of horizontal wellbore variable mass flow coupled with flow in a reservoir［J］.Journal of Southwest Petroleum Institute，2000，22（2）：36－39.

［9］Huang Shijun，Cheng Linsong，Zhao Fenglan，et al.The flow model coupling reservoir percolation and variable mass pipe flow in production section of the stepped horizontal well［J］.Journal of Hydrodynamics，2005，20（4）：463－471.

［10］Novvy R A.Pressure drop in horizontal wells：When can they be ignored？［J］.SPE 24941，1992.

［11］Dikken B J.Pressure drop in the horizontal wells and its effects on their production performance［J］.JPT，1990：1426－1433.

［12］Yu Lexiang，Zhou Letian，Zhang Qi，et al.Pressure gradient model for variable mass fluid flow in horizontal wellbore［J］. Journal of University of Petroleum，China，2001，25（4）：47－50.

［13］Prats M，Hazebrock P，Strickler W R.Effect of vertical fractures on reservoir behavior-imcompressible-fuild case［J］.SPEJ，1962，22（5）：87－94.

［14］Lian Peiqing，Tong Dengke，Cheng Linsong，et al.A coupling model low permeability reservoir and fractured horizontal wellbore in nonsteady state［J］.Chinese J Comput Phys，2010，27（2）：203－210.

［15］Cheng Linsong，Pi Jian，Lian Peiqing，et al.A computational method of productivity of horizontal well in naturally fractured reservoirs［J］.Chinese J Comput Phys，2011，28（2）：230－236.

［16］Ye Shuangjiang，Jiang Hanqiao，Li Junjian，et al.Productivity calculation of infill horizontal wells in mixed well pattern［J］. Chinese J Comput Phys，2011，28（2）：693－697.

［17］Du Dianfa，Li Dongdong，Shi Dayou，et al.A study on heavy oil well test［J］.Chinese J Comput Phys，2011，28（3）：385－396.

A Coupling Flow Model of Finite-conductivity Fractured Directional Well

JIA Pin，CHENG Linsong，HUANG Shijun，F(xiàn)ANG Sidong
（Faculty of Petroleum Engineering，China University of Petroleum（Beijing），Beijing 102249，China）

Considering artificial fracture，inclined wellbore and formation coupling flow，flow to directional well is divided into three parts：Flow from reservoir to fracture，flow in fracutre and variable mass flow in production pipe.Principles of potential superposition and mirror reflection as well as the concept of infinitesimal line congruence are used to model flow from reservoir to fracture and potential distribution in reservoir.By discreting fracuture to 2D grid，boundary theory is used to couple flow dynamics in fracture with variable mass flow in inclined pipe.A comprehensive coupling model for finite-conductivity fractured directional well is shown. Iteration method is used to solve the model.A practical case shows that for vertical fracutre result of coupling model agrees with that of Prats method.Fracture conductivity and inclined angle have great effects on productivity and pressure.

fractured directional well；finite-conductivity fracutre flow；coupling flow model；discreted grid；productivity

date：2013－09－17；Revised date：2013－11－27

TE312

A

2013－09－17；

2013－11－27

國家973（2013CB228000）及國家自然科學(xué)基金（51174215／E0403）資助項目

賈品（1990－），男，碩士研究生，從事油藏工程和滲流理論研究，E-mail：jiapin1990＠163.com

1001-246X（2014）05-0559-08