姚 軍, 李志豪, 孫 海

(中國(guó)石油大學(xué)(華東)石油工程學(xué)院,山東青島 266580)

頁(yè)巖氣藏屬于致密性氣藏,孔滲性質(zhì)差,常規(guī)開采難度大,必須采用水力壓裂措施進(jìn)行儲(chǔ)層改造,其兩個(gè)關(guān)鍵的技術(shù)是水平井和水力壓裂。壓裂技術(shù)在實(shí)際油氣藏實(shí)施很成功,然而在其設(shè)計(jì)優(yōu)化上做的工作還不足。據(jù)統(tǒng)計(jì)大約有三分之一的壓裂都沒(méi)有實(shí)現(xiàn)預(yù)期的產(chǎn)能[1],主要原因是裂縫和油藏特征的復(fù)雜性以及大量設(shè)計(jì)變量。研究問(wèn)題主要包括水平井位置選擇、井間距調(diào)整、裂縫條數(shù)與裂縫半長(zhǎng)的設(shè)計(jì)問(wèn)題。對(duì)于壓裂優(yōu)化這個(gè)“黑箱”問(wèn)題,亟需建立一套完整的優(yōu)化設(shè)計(jì)體系來(lái)指導(dǎo)現(xiàn)場(chǎng)施工作業(yè)。從井參數(shù)到裂縫參數(shù)設(shè)計(jì)會(huì)耗費(fèi)大量的時(shí)間,選擇合適的優(yōu)化方法尤為重要。目前已有較多學(xué)者對(duì)裂縫參數(shù)優(yōu)化進(jìn)行了研究,主要是基于無(wú)梯度算法和代理模型的使用。Ma等[2]對(duì)比分析了SPSA、GA和CMA-ES算法優(yōu)化井位置和壓裂段數(shù),充分展現(xiàn)了無(wú)梯度算法的實(shí)用性。代理模型在傳統(tǒng)油氣藏優(yōu)化里應(yīng)用廣泛,如產(chǎn)能優(yōu)化[3]、水驅(qū)[4]、氣驅(qū)[5]、化學(xué)驅(qū)[6]、蒸汽輔助重力采油[7]和歷史擬合[8]等;非常規(guī)油氣藏里代理模型應(yīng)用也很多,Gorucu等[9]基于ANN構(gòu)建了致密氣藏代理模型,前過(guò)程預(yù)測(cè)生產(chǎn)結(jié)果,逆過(guò)程反求最優(yōu)參數(shù);Wilson等[10]將原模型簡(jiǎn)化為內(nèi)外擁有不同孔滲性質(zhì)的簡(jiǎn)化物理模型,在此基礎(chǔ)上進(jìn)行優(yōu)化模擬,大大提高了優(yōu)化效率;Al-mudhafar等[11]基于二次多項(xiàng)式模型和徑向基模型在固定裂縫條數(shù)前提下優(yōu)化了單井時(shí)裂縫的其他參數(shù),結(jié)論是徑向基網(wǎng)絡(luò)模型得到的結(jié)果更加準(zhǔn)確。當(dāng)前研究存在的不足在于,多數(shù)學(xué)者都是在保持其他參數(shù)不變的情況下,對(duì)某幾個(gè)裂縫參數(shù)進(jìn)行敏感性分析,而不是從整體的角度去優(yōu)化,得到的結(jié)果往往不是最優(yōu)的;無(wú)梯度算法固然靈活,但在保證全局搜索時(shí)一定程度上降低了尋優(yōu)速度。筆者建立頁(yè)巖氣藏水平井流動(dòng)模型,將代理輔助分層粒子群算法(SHPSO)應(yīng)用到井工廠模式下,從水平井空間位置到裂縫的條數(shù)、間距和半長(zhǎng),以凈現(xiàn)值最大化為目標(biāo)優(yōu)選壓裂參數(shù)。

1 頁(yè)巖氣藏模型與求解函數(shù)

1.1 頁(yè)巖氣藏模型

基于商業(yè)數(shù)值模擬軟件ECLIPSE構(gòu)建頁(yè)巖氣藏水平井流動(dòng)數(shù)學(xué)模型。將氣水兩相、多重介質(zhì)模型、非達(dá)西流動(dòng)、水平井摩阻和吸附解吸因素考慮到模型中。使用多重介質(zhì)模型表征從基質(zhì)到裂縫的暫態(tài)流動(dòng),如圖1所示,Vp表示孔隙體積。

圖1 多重介質(zhì)模型示意圖

不同于致密氣藏,頁(yè)巖氣藏富含有機(jī)質(zhì),以吸附態(tài)形式存在于干酪根和黏土顆粒表面的氣體比例可多達(dá)60%[12]。本文中采用Langmuir等溫吸附公式描述吸附解吸效應(yīng),對(duì)應(yīng)的方程為

(1)

式中,θ為尺度系數(shù);ps為標(biāo)況壓力,1.013×105Pa;R為通用氣體常數(shù);Ts為標(biāo)況溫度,273 K;pc為當(dāng)前壓力,Pa;VL為L(zhǎng)angmuir體積,m3/kg;pL為L(zhǎng)angmuir壓力,Pa。

此外,壓裂后的頁(yè)巖氣藏,氣體在高導(dǎo)流水力裂縫中的流動(dòng)為高速非達(dá)西流,采用Forchheimer方程進(jìn)行描述,即

(2)

其中

式中,q為體積流量,m3/s;k為巖石滲透率,10-3μm2;Kr為相對(duì)滲透率;A為流動(dòng)區(qū)域面積,m3;μ為流體黏度,Pa·s;ρ為流體密度,kg/m3;β為Forchheimer參數(shù);φ為孔隙度;S為飽和度;a、b、c、d為常數(shù)。

頁(yè)巖氣藏壓裂后,水力裂縫激活天然裂縫,從而形成儲(chǔ)層改造體積(stimulated reservoir volume,SRV),極大改善了水力裂縫周邊區(qū)域的孔滲性質(zhì)。本文中SRV定義為圍繞水力裂縫外圈網(wǎng)格部分,如圖2灰色區(qū)域所示。

圖2 SRV示意圖

1.2 求解函數(shù)

本文中采用凈現(xiàn)值(net present value,NPV)為適應(yīng)度評(píng)價(jià)函數(shù)。模型涉及到的變量主要有12個(gè),分別為井1的坐標(biāo)x和y、橫向井間距L1、縱向井間距L2、裂縫條數(shù)Ni(i=1～4)、裂縫半長(zhǎng)Lxfi,

u=argmaxJ.

(3)

其中

lmin≤ui≤lmax,i=1,2,…,12.

式中,J為凈現(xiàn)值;Qg為年產(chǎn)氣量,m3;Qw為年產(chǎn)水量,m3;rg為天然氣價(jià)格,元/m3;rw為水處理費(fèi)用,元/m3;b為折扣率;m為總井?dāng)?shù);n為生產(chǎn)總時(shí)間段數(shù);Δt為折扣時(shí)間段,a;O為每口井每天的操作費(fèi)用,元;Cw為鉆豎直井段成本,元;N為裂縫條數(shù);Cf為單條裂縫造價(jià),元;Lw為水平井長(zhǎng)度,m;Cp為單位長(zhǎng)度水平井造價(jià),元/m;ui為第i個(gè)變量;lmin和lmax分別為第i個(gè)變量的下限和上限值。

水平井分布示意圖如圖3所示。

圖3 水平井分布示意圖

2 基于代理輔助的分層粒子群算法

代理輔助分層粒子群(SHPSO)算法由Yu等[13]首次提出,在解決多維問(wèn)題快速尋優(yōu)時(shí)表現(xiàn)良好。本文中將該算法運(yùn)用到頁(yè)巖氣產(chǎn)能優(yōu)化。該算法將代理模型與粒子群算法結(jié)合,既保留了粒子群算法的全局尋優(yōu)能力,又融入了代理模型的預(yù)測(cè)功能,相比普通智能算法擁有更高的尋優(yōu)速度和精度。

2.1 徑向基網(wǎng)絡(luò)模型

代理模型是指通過(guò)訓(xùn)練一復(fù)雜度低的數(shù)學(xué)模型,替代原有模型進(jìn)行設(shè)計(jì)優(yōu)化[14]。常用的代理模型有多項(xiàng)式響應(yīng)面法、克里金方法、徑向基函數(shù)法和神經(jīng)網(wǎng)絡(luò)方法等。本文中使用徑向基網(wǎng)絡(luò)方法(RBFNN)構(gòu)建代理模型。RBFNN是一種3層神經(jīng)網(wǎng)絡(luò),包括輸入層、隱含層和輸出層。其基本思想是用徑向基作為隱單元的“基”構(gòu)成隱含層空間,這樣就可以將輸入矢量直接映射到隱空間,而不需要像人工神經(jīng)網(wǎng)絡(luò)那樣進(jìn)行權(quán)連接,結(jié)構(gòu)示意圖如圖4所示。本文中輸入層為1.2節(jié)中的12個(gè)參數(shù)變量,隱含層為轉(zhuǎn)化后的高維線性可分空間,輸出層為凈現(xiàn)值。RBFNN具有唯一最佳逼近特性,且無(wú)局部極小問(wèn)題存在,具有很強(qiáng)的輸入輸出映射功能,學(xué)習(xí)算法簡(jiǎn)單方便。

圖4 RBFNN示意圖

2.2 拉丁超立方抽樣方法

經(jīng)典抽樣方法主要包含蒙特卡洛抽樣和拉丁超立方抽樣方法(Latin hypercube sampling,LHS)。蒙特卡洛抽樣是完全隨機(jī)的,在求解空間內(nèi)樣本可以落在任何位置,而且樣本常落于高分布概率區(qū)域,當(dāng)執(zhí)行迭代次數(shù)較少時(shí)會(huì)產(chǎn)生聚集問(wèn)題。與蒙特卡洛方法相比,拉丁超立方抽樣被設(shè)計(jì)成通過(guò)較少次數(shù)的抽樣來(lái)準(zhǔn)確地重建輸入分布,其關(guān)鍵是對(duì)輸入概率分布進(jìn)行分層,可以隨機(jī)而且系統(tǒng)地分布在求解空間內(nèi),所取得的樣本更具有代表性。

本文中使用拉丁超立方抽樣方法進(jìn)行采樣,并與隨機(jī)采樣(random sampling,RS)方法得到的數(shù)據(jù)分布進(jìn)行對(duì)比,結(jié)果如圖5所示。以二維空間、樣本點(diǎn)為10為例,相比隨機(jī)采樣,拉丁超立方抽樣方法采樣每行每列有且僅有一個(gè)樣本點(diǎn),保證了樣本散落均勻,能更廣泛地覆蓋變量空間。

圖5 LHS采樣與RS采樣對(duì)比

2.3 SHPSO算法流程

SHPSO是基于傳統(tǒng)粒子群(PSO)算法與社會(huì)學(xué)習(xí)粒子群算法SLPSO分層組合進(jìn)行優(yōu)化的。基本公式與傳統(tǒng)PSO類似,包含速度更新公式與位置更新公式。與傳統(tǒng)PSO不同的是其考慮了RBF代理模型輔助優(yōu)化而產(chǎn)生的最優(yōu)值sr,即RBFNN尋找到的最優(yōu)解。如果sr對(duì)應(yīng)的真實(shí)值優(yōu)于全局最優(yōu)值sg,則速度公式改為用sr代替sg求解。同理,在更新位置全局最優(yōu)sg和個(gè)體最優(yōu)sp時(shí)也采用對(duì)應(yīng)的sr,基本公式為

(4)

其中

式中,x為粒子位置;v為粒子速度;r為0到1的隨機(jī)數(shù);c1=c2=2.05。

內(nèi)部SLPSO用于求解RBFNN模型的最優(yōu)解,外部PSO基于更新的樣本庫(kù)迭代求解。其基礎(chǔ)框架如圖6所示。

圖6 SHPSO循環(huán)示意圖

樣本庫(kù)的更新來(lái)源為RBFNN求得的最優(yōu)解以及外層PSO中優(yōu)于sp的個(gè)體。算法的詳細(xì)流程如圖7所示。

圖7 SHPSO步驟示意圖

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

4口井分布于長(zhǎng)2 500 m、寬1 500 m、厚60 m的頁(yè)巖氣藏中,各井坐標(biāo)分別由井1的x、y坐標(biāo)以及橫向井間距L1、縱向井間距L2確定。裂縫等距分布在長(zhǎng)1 100 m的水平井內(nèi),間距是裂縫條數(shù)的函數(shù),模型如圖3所示。

模型基礎(chǔ)數(shù)據(jù)參照Ma等[2]歷史擬合過(guò)的Barnett頁(yè)巖,初始?xì)獠貕毫?0.6 MPa,井底流壓為3.5 MPa,生產(chǎn)年限為20 a,天然裂縫滲透率為0.000 1×10-3μm2,天然裂縫孔隙度為0.000 05,水力裂縫導(dǎo)流能力為0.3×10-3μm2·m,Langmuir壓力為4.48 MPa,Langmuir體積為0.002 7 m3/kg,SRV區(qū)域滲透率為0.000 5×10-3μm2;基質(zhì)滲透率為(0.000 1～0.000 3)×10-3μm2,基質(zhì)孔隙度為0.04～0.08,孔滲性質(zhì)示意圖如圖8所示。

圖8 基質(zhì)孔隙度與滲透率分布

受模型限制,變量取值范圍及經(jīng)濟(jì)參數(shù)取值[15]:裂縫條數(shù)為10～40,裂縫半長(zhǎng)為100～300 m,橫向井間距為1 100～1 500 m,縱向井間距為200～800 m,井1的橫坐標(biāo)為10～200 m,縱坐標(biāo)為100～900 m。氣體價(jià)格為0.791元/m3,折扣率為12.5%,豎直段鉆井成本為140萬(wàn)元,井日常操作費(fèi)用為420元。

裂縫造價(jià)參照馮其紅等[16]回歸分析得到的計(jì)算公式,

(5)

式中,Lxf為裂縫半長(zhǎng),m。

3.1 傳統(tǒng)單因素

針對(duì)裂縫半長(zhǎng)、裂縫條數(shù)、橫向和縱向井間距進(jìn)行傳統(tǒng)單因素分析。為簡(jiǎn)化計(jì)算,分析裂縫半長(zhǎng)與條數(shù)時(shí)僅取單井進(jìn)行模擬,平均基質(zhì)孔隙度為0.050 7,平均基質(zhì)滲透率為0.000 148×10-3μm2。圖9為不同裂縫條數(shù)下凈現(xiàn)值與裂縫半長(zhǎng)的關(guān)系曲線。由圖9可知,在不同裂縫條數(shù)下,凈現(xiàn)值均隨著裂縫半長(zhǎng)的增大而先增后減,最優(yōu)半長(zhǎng)為230 m。這主要是由于半長(zhǎng)值較小時(shí),增加半長(zhǎng)帶來(lái)的收益高于成本;隨著半長(zhǎng)增加對(duì)壓裂設(shè)備性能要求提高,成本上升。對(duì)比不同裂縫條數(shù)的曲線可知,在相同裂縫半長(zhǎng)的條件下,隨著裂縫條數(shù)的增加,凈現(xiàn)值先增加后減少,且最優(yōu)裂縫條數(shù)為26條。此時(shí)單條裂縫的控制區(qū)域可達(dá)最優(yōu)。

圖9 凈現(xiàn)值與裂縫半長(zhǎng)的關(guān)系

固定裂縫條數(shù)N=26、半長(zhǎng)Lxf=230 m和井1坐標(biāo)(x=50 m,y=400 m)條件下,分別求得均質(zhì)和非均質(zhì)時(shí)凈現(xiàn)值對(duì)橫、縱向井間距的關(guān)系曲線,如圖10所示。均質(zhì)條件下凈現(xiàn)值與橫、縱向井間距呈相同的變化趨勢(shì),均為先增大后減小。這是由于間距過(guò)小時(shí)井間干擾嚴(yán)重;間距過(guò)大時(shí),氣井生產(chǎn)受邊界影響變大。非均質(zhì)條件下,凈現(xiàn)值與縱向井間距變化規(guī)律不變,而隨橫向井間距的減小而增大,考慮其原因是受地層非均質(zhì)的影響,使布井傾向于豐度更高的位置。

圖10 凈現(xiàn)值與井間距的關(guān)系

3.2 基于SHPSO的參數(shù)優(yōu)化

采用數(shù)值模擬軟件ECLIPSE與MATLAB結(jié)合,基于SHPSO對(duì)裂縫半長(zhǎng)、裂縫條數(shù)、井坐標(biāo)和橫、縱向井距進(jìn)行優(yōu)化?；緟?shù)設(shè)置:外層粒子群種群數(shù)為50,允許最大調(diào)用ECLIPSE次數(shù)為250次;內(nèi)層粒子群種群數(shù)為100,最大迭代數(shù)為500,終止誤差為1×10-6,初始樣本庫(kù)數(shù)量設(shè)為80。圖11為樣本庫(kù)隨算法迭代而更新、排序及篩選的過(guò)程。

由圖11可知,首先基于LHS生成了80個(gè)遍布變量空間的隨機(jī)解,其凈現(xiàn)值均值為2.40×107元;隨迭代樣本庫(kù)數(shù)量增加至125、154、187、212、250,對(duì)應(yīng)凈現(xiàn)值均值分別為2.96×107、3.13×107、3.41×107、3.57×107和3.74×107元。樣本庫(kù)的更新來(lái)源只有2個(gè),即SLPSO尋得的最優(yōu)解sr和外層粒子大于個(gè)體最優(yōu)的解,均為篩選對(duì)比之后的變量組合。由散落的樣本數(shù)據(jù)可知,隨樣本庫(kù)數(shù)量的增加,新增的解逐步收斂至最優(yōu)值。

圖11 樣本選取示意圖

對(duì)于RBFNN模型訓(xùn)練集為每代更新后的前50個(gè)樣本,測(cè)試集為全體變量集合。一般來(lái)講,使用的訓(xùn)練樣本數(shù)越多,RBFNN模型越準(zhǔn)確;同時(shí)存在一閾值,超過(guò)該值會(huì)出現(xiàn)過(guò)擬合現(xiàn)象。樣本庫(kù)隨迭代實(shí)時(shí)更新,令訓(xùn)練集始終保持為50個(gè),既可以避免過(guò)擬合,同時(shí)也保證了訓(xùn)練精度和質(zhì)量。

裂縫分布如圖12所示。求解結(jié)果為井1坐標(biāo)x=40 m,y=700 m,4口井裂縫條數(shù)分別為26、24、27、25,裂縫半長(zhǎng)分別為230、220、230和230 m,橫向井間距為1 110 m,縱向井間距為530 m,最優(yōu)解為5.781×107元。分別比較單因素分析與整體分析、SHPSO與普通智能算法對(duì)比的結(jié)果。

圖12 裂縫分布

3.2.1 與傳統(tǒng)單因素分析對(duì)比

單因素分析是在固定其他參數(shù)時(shí)得到某一變量對(duì)整體凈現(xiàn)值的影響,在變量維度小時(shí)有效可靠;但當(dāng)多口井進(jìn)行組合優(yōu)化時(shí),變量維數(shù)增多,若再以非整體的思想去分析問(wèn)題,得到的結(jié)論不一定是合理的,尤其在非均質(zhì)條件下優(yōu)化布井問(wèn)題時(shí),井位的選擇需要結(jié)合每口井的裂縫參數(shù)去優(yōu)化考慮。根據(jù)最終結(jié)果可知,并非每口井參數(shù)都相同;當(dāng)4口井條數(shù)均為26,半長(zhǎng)均為230 m時(shí)的凈現(xiàn)值為5.541×107元,相比之下最優(yōu)解提高了4.3%(約240萬(wàn)元)?？梢?從成本的角度來(lái)講,整體優(yōu)化要比單因素分析結(jié)果更合理。

3.2.2 與普通智能算法對(duì)比

以調(diào)用ECLIPSE次數(shù)為橫坐標(biāo),縱坐標(biāo)為當(dāng)前全局最優(yōu)。分別采用標(biāo)準(zhǔn)PSO和SHPSO對(duì)模型進(jìn)行求解,結(jié)果如圖13所示。

圖13 SHPSO與PSO對(duì)比尋優(yōu)曲線

相比PSO算法,在調(diào)用ECLIPSE次數(shù)相同的前提下SHPSO算法能更快地捕捉到最優(yōu)值的鄰域,而且最終結(jié)果更優(yōu)。SHPSO平均最優(yōu)凈現(xiàn)值為5.74×107元,PSO平均最優(yōu)凈現(xiàn)值為5.37×107元,同比提高了6.89%。拉丁超立方抽樣數(shù)設(shè)為80,可以最大限度地遍歷變量空間;RBFNN輔助預(yù)測(cè)可以為算法提供潛在最優(yōu)方向。程序絕大部分時(shí)間用于ECLIPSE求解產(chǎn)能,相比之下代理模型的訓(xùn)練和尋優(yōu)可以忽略不計(jì),故將運(yùn)行時(shí)間看作與調(diào)用次數(shù)成正比。雖然智能算法可以有效避免求解梯度問(wèn)題,但對(duì)于多維數(shù)、單次運(yùn)行時(shí)間長(zhǎng)的情形難以有效解決,而代理模型的結(jié)合很好地彌補(bǔ)了這個(gè)缺陷。相比之下,SHPSO比普通智能算法擁有更快的收斂速度和結(jié)果準(zhǔn)確度。

4 結(jié) 論

(1)MATLAB與ECLIPSE的結(jié)合使用實(shí)現(xiàn)了自動(dòng)尋求最優(yōu)參數(shù)組合。優(yōu)化方法采用基于代理模型輔助優(yōu)化的分層粒子群算法SHPSO,相比普通智能算法,SHPSO擁有更快的收斂速度。

(2)相比傳統(tǒng)單因素分析,從整體的角度去優(yōu)化結(jié)果更加合理。井工廠模式同步壓裂并非每口井簡(jiǎn)單地復(fù)制疊加,單井的設(shè)計(jì)參數(shù)不能照搬到井工廠模式下,而是需要作為一個(gè)整體去優(yōu)化考慮,這樣才能更符合現(xiàn)場(chǎng)實(shí)際。