劉偉,趙輝,曹琳,張凱

(1.長江大學(xué)石油工程學(xué)院,湖北 武漢 430100;2.中國石油大學(xué)(華東)石油工程學(xué)院,山東 青島 266580)

0 引言

油藏?cái)?shù)值模擬歷史擬合是油藏?cái)?shù)值模擬中的關(guān)鍵環(huán)節(jié),為油藏未來開發(fā)動態(tài)預(yù)測和開發(fā)方案制定提供正確模型[1-3].傳統(tǒng)人工歷史擬合速度慢,準(zhǔn)確性低.若依賴于油藏工作者對油藏的分析認(rèn)識和擬合經(jīng)驗(yàn),就無法滿足油田開發(fā)高效化、實(shí)時(shí)化的要求.油藏模擬輔助歷史擬合是基于最優(yōu)化理論,并借助計(jì)算機(jī)實(shí)現(xiàn)對油藏模型參數(shù)的自動調(diào)整,相比于人工歷史擬合,可以極大地提高擬合效率和模型的準(zhǔn)確度[4-6].

歷史擬合通常是在一定網(wǎng)格區(qū)域內(nèi)調(diào)整模型參數(shù),得到準(zhǔn)確的單井敏感性區(qū)域,能提升歷史擬合質(zhì)量.單井敏感性區(qū)域是劃定一個(gè)臨界區(qū)域,認(rèn)為單井的生產(chǎn)動態(tài)數(shù)據(jù)只與臨界區(qū)域內(nèi)網(wǎng)格的模型參數(shù)存在相關(guān)性,臨界區(qū)域外的網(wǎng)格對生產(chǎn)動態(tài)數(shù)據(jù)沒有影響[6-8].人工歷史擬合是基于已有的油藏認(rèn)知或以往經(jīng)驗(yàn),劃定一個(gè)不精確的區(qū)域.在輔助歷史擬合研究中,研究者通?；诘刭|(zhì)信息(如沉積相態(tài)等),結(jié)合距離截?cái)喾椒╗6-8]得到單井敏感性區(qū)域,但距離截?cái)啻嬖谝欢ǖ娜藶橐蛩?不夠準(zhǔn)確.另外,歷史擬合需要反演的油藏模型參數(shù)通常數(shù)以十萬計(jì),應(yīng)用主流無梯度優(yōu)化算法[9-12]的計(jì)算代價(jià)大.參數(shù)化[13-15]是一種提高計(jì)算效率的途徑,其依據(jù)特定方式對油藏模型重新參數(shù)化,以減少模型參數(shù)的數(shù)量.

Fast-Marching Method(FMM)是一種可以快速高效追蹤邊界運(yùn)移情況的方法,Sharifi等[16]應(yīng)用FMM來標(biāo)定油藏中單井控制儲量.本文將FMM引入歷史擬合單井敏感性區(qū)域的確定中,充分考慮油藏模型靜態(tài)參數(shù)信息,準(zhǔn)確獲得單井敏感性區(qū)域,且不需要進(jìn)行油藏模擬計(jì)算,網(wǎng)格計(jì)算效率高.獲得單井敏感性區(qū)域后,提出一種簡單的參數(shù)化方法,將單井敏感性區(qū)域內(nèi)模型參數(shù)的平均值,作為歷史擬合目標(biāo)函數(shù)中的參數(shù)變量,極大降低了變量的維數(shù),提高了計(jì)算效率.

1 單井敏感性區(qū)域定量表征方法

FMM是通過求解波動傳播問題的程函方程獲知邊界的運(yùn)移,程函方程的一般形式為

式中:F為速度,cm/s0.5;τ為飛行時(shí)間,s0.5;x為飛行空間位置.

在油藏多孔介質(zhì)中,井點(diǎn)處油藏壓力傳播方程也可以表示為程函方程的形式[17],在靜態(tài)地質(zhì)參數(shù)場下,可利用FMM追蹤計(jì)算壓力波從井點(diǎn)穿過每個(gè)網(wǎng)格的飛行時(shí)間,利用追蹤的飛行時(shí)間劃定油藏三維單井敏感性區(qū)域.

對于油藏網(wǎng)格系統(tǒng),假設(shè)油藏綜合壓縮系數(shù)和流體黏度是常數(shù).壓力波穿過每個(gè)網(wǎng)格的速度定義為

式中:K為滲透率,μm2;μ為流體黏度,mPa.s;Ct為油藏綜合壓縮系數(shù),MPa-1;φ為孔隙度.

考慮油藏的非均質(zhì)性,基于逆風(fēng)差分近似求解程函方程,計(jì)算飛行時(shí)間的方程為

其中

式中:Fx,Fy,Fz分別為網(wǎng)格x,y,z方向上的速度,cm/s0.5;Δx,Δy,Δz分別為 x,y,z方向上的網(wǎng)格尺寸,cm;i,j,k分別為油藏模型網(wǎng)格體系中x,y,z方向上的坐標(biāo).

以二維非均質(zhì)油藏模型為例,油藏網(wǎng)格飛行時(shí)間的追蹤過程如圖1所示.

1)將某一確定網(wǎng)格作為初始起點(diǎn)(τ=0),標(biāo)記為凍結(jié)點(diǎn)(紅色網(wǎng)格)(見圖1a的網(wǎng)格13).

2)找到與凍結(jié)點(diǎn)相鄰的網(wǎng)格,標(biāo)記為相鄰點(diǎn)(綠色網(wǎng)格)(見圖1b的網(wǎng)格8,12,14,18).

3)利用式(3)計(jì)算壓力波從凍結(jié)點(diǎn)到達(dá)相鄰點(diǎn)的飛行時(shí)間,選出飛行時(shí)間值最小的相鄰點(diǎn),標(biāo)記為凍結(jié)點(diǎn)(見圖1c的網(wǎng)格12).

4)將網(wǎng)格12作為新的初始起點(diǎn),并找到網(wǎng)格12的相鄰點(diǎn)(不包括已經(jīng)計(jì)算過飛行時(shí)間的網(wǎng)格)(見圖1d的網(wǎng)格7,11,17).

5)計(jì)算壓力波從網(wǎng)格12到達(dá)網(wǎng)格7,11,17的飛行時(shí)間,并找到所有相鄰點(diǎn)中飛行時(shí)間值最小的點(diǎn),標(biāo)記為凍結(jié)點(diǎn)(見圖1e的網(wǎng)格18).

6)重復(fù)步驟4)-5),直到所有的網(wǎng)格被標(biāo)記為凍結(jié)點(diǎn)(見圖1f).

圖1 飛行時(shí)間追蹤過程示意

在計(jì)算整個(gè)油藏網(wǎng)格的飛行時(shí)間過程中,通常將井點(diǎn)處網(wǎng)格作為初始起點(diǎn).離井點(diǎn)越近,飛行時(shí)間值越小;離井點(diǎn)越遠(yuǎn),飛行時(shí)間值越大.確定單井敏感性區(qū)域的具體處理方法為:當(dāng)油藏存在多口井時(shí),每一個(gè)網(wǎng)格對每口井會得到不同的飛行時(shí)間,選出最小的飛行時(shí)間值,將該網(wǎng)格劃歸于最小飛行時(shí)間值所對應(yīng)井的敏感性區(qū)域;基于得到的單井敏感性區(qū)域,將區(qū)域內(nèi)油藏模型參數(shù)的平均值,作為歷史擬合目標(biāo)函數(shù)中的變量元素,不以單個(gè)網(wǎng)格的模型參數(shù)為變量反演對象.

2 歷史擬合數(shù)學(xué)模型及優(yōu)化求解

2.1 數(shù)學(xué)模型的建立及油藏模型的參數(shù)化

以擬合歷史動態(tài)生產(chǎn)數(shù)據(jù)為目標(biāo),歷史擬合目標(biāo)函數(shù)定義為[11]

變量u的元素包含孔隙度、滲透率、飽和度等油藏模型參數(shù),向量dobs的元素包含壓力、含水率、產(chǎn)油量等歷史動態(tài)生產(chǎn)數(shù)據(jù).對于油藏模擬歷史擬合問題,本質(zhì)是利用最優(yōu)化方法求解變量u,使得目標(biāo)函數(shù)Z()u取得最小值.

對于實(shí)際油藏歷史擬合,以單個(gè)網(wǎng)格的模型參數(shù)為反演對象,變量u的維數(shù)巨大.本文提出一種近似的參數(shù)化方法,采用單井敏感性區(qū)域內(nèi)網(wǎng)格模型參數(shù)的平均值,作為變量u的元素.由于油藏模型參數(shù)在每一層存在差異,因此模型參數(shù)的平均值是對每一層的敏感性區(qū)域進(jìn)行計(jì)算.變量u的維數(shù)僅與油藏井?dāng)?shù)、油藏模型的縱向?qū)訑?shù)和反演參數(shù)的種類數(shù)有關(guān),與油藏模型的網(wǎng)格數(shù)量沒有關(guān)系.參數(shù)化后,變量維數(shù)極大降低,計(jì)算效率提高.

依據(jù)擬合前后單井敏感性區(qū)域內(nèi)網(wǎng)格模型參數(shù)平均值的倍數(shù)變化關(guān)系,更新油藏模型參數(shù),以孔隙度為例,參數(shù)化公式為

式中:φI,J,K,φI,J,K′分別為擬合前、后第I井敏感性區(qū)域內(nèi)第J層中第K個(gè)網(wǎng)格的孔隙度;φave,I,J,φave,I,J′分別為擬合前、后第I井敏感性區(qū)域內(nèi)第J層中網(wǎng)格孔隙度的平均值;Nw為總井?dāng)?shù);Nz為油藏模型縱向?qū)訑?shù);nI,J為第I井敏感性區(qū)域內(nèi)第J層的網(wǎng)格總數(shù).

2.2 近似擾動梯度升級算法

本文采用趙輝等[18]提出的近似擾動梯度升級算法,對歷史擬合目標(biāo)函數(shù)進(jìn)行優(yōu)化求解,該方法是目前常用歷史擬合算法SPSA的一種改進(jìn)算法.對于目標(biāo)函數(shù)Z(u),考慮在第l個(gè)迭代步的最優(yōu)模型參數(shù)變量為ul,在其周圍進(jìn)行擾動生成N個(gè)模型參數(shù)變量:

式中:ul,i為第l步第i個(gè)模型參數(shù)變量;γ為擾動步長;Ll,i為服從某分布的擾動向量.

式中:L為N個(gè)擾動向量組成的矩陣;ΔZ為N維向量;a為常數(shù);M為一個(gè)N維下三角矩陣.

向量ΔZ的元素為模型參數(shù)變量ul,i對應(yīng)的目標(biāo)函數(shù)值與當(dāng)前最優(yōu)目標(biāo)函數(shù)值的差值.另外,當(dāng)Ll,i服從伯努利分布,M為單位矩陣時(shí),式(7)即是SPSA算法的梯度計(jì)算公式.

由于式(10)中除M外都為確定量,不牽涉油藏模擬器運(yùn)算,可以使用傳統(tǒng)的優(yōu)化算法[19]對其進(jìn)行優(yōu)化迭代,如擬牛頓法.每個(gè)迭代步都對M進(jìn)行優(yōu)化,然后采用優(yōu)化后的M構(gòu)造梯度,這樣得到的梯度更加接近真實(shí)梯度.獲得近似擾動梯度后,更新模型參數(shù)變量:

式中:ul+1為第l+1迭代步參數(shù)向量;ul為第l迭代步參數(shù)向量;αl為迭代步長;為的無窮范數(shù).

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

運(yùn)用本文提出的方法對某碳酸鹽巖油藏進(jìn)行了歷史擬合測試.該油藏非均質(zhì)性極強(qiáng),油藏模型網(wǎng)格劃分為76X115X20,總網(wǎng)格數(shù)為174 800,共有油水井29口.擬合的歷史動態(tài)生產(chǎn)數(shù)據(jù)是油藏生產(chǎn)4 564 d的區(qū)塊累計(jì)產(chǎn)油量和單井日產(chǎn)油量,反演的油藏模型參數(shù)為滲透率和孔隙度,參數(shù)化后變量的維數(shù)為1 160.油藏模型第1層網(wǎng)格飛行時(shí)間的計(jì)算結(jié)果見圖2,根據(jù)飛行時(shí)間確定第1層單井敏感性區(qū)域(見圖3).

圖2 第1層網(wǎng)格飛行時(shí)間

圖3 第1層網(wǎng)格單井敏感性區(qū)域

歷史擬合目標(biāo)函數(shù)優(yōu)化結(jié)果如圖4所示,在經(jīng)過約35次迭代步后,目標(biāo)函數(shù)逐漸收斂.油藏的區(qū)塊累計(jì)產(chǎn)油量和單井日產(chǎn)油量擬合結(jié)果如圖5、圖6所示,擬合后模型的生產(chǎn)數(shù)據(jù)較擬合前更接近實(shí)際值,擬合效果好.從歷史擬合前后油藏模型第1層網(wǎng)格滲透率場分布可以看出,擬合前后滲透率變化較大(見圖7、圖8).圖中滲透率取對數(shù)值.

圖4 目標(biāo)函數(shù)優(yōu)化結(jié)果

圖5 區(qū)塊累計(jì)產(chǎn)油量擬合結(jié)果

圖6 單井日產(chǎn)油量擬合結(jié)果

圖7 第1層網(wǎng)格擬合前平面滲透率

圖8 第1層網(wǎng)格擬合后平面滲透率

4 結(jié)論

1)針對歷史擬合中單井敏感性區(qū)域的確定,基于油藏中壓力傳播方程,采用FMM求解其程函方程,提出了定量表征油藏中單井敏感性區(qū)域的計(jì)算方法.該方法可以快速準(zhǔn)確地確定油藏三維單井敏感性區(qū)域,具有穩(wěn)定性好、計(jì)算效率高的優(yōu)點(diǎn),且適用于大型實(shí)際油藏.

2)針對實(shí)際油藏應(yīng)用中歷史擬合目標(biāo)函數(shù)的參數(shù)變量維數(shù)大的問題,提出以油藏各層中單井敏感性區(qū)域內(nèi)網(wǎng)格模型參數(shù)的平均值作為歷史擬合目標(biāo)函數(shù)中參數(shù)變量的元素,大幅度降低了參數(shù)變量的維數(shù),提高了計(jì)算效率.

3)油藏測試結(jié)果表明,基于FMM確定的單井敏感性區(qū)域,采用近似擾動梯度升級算法,對油藏區(qū)塊和單井動態(tài)生產(chǎn)指標(biāo)的擬合取得了較好的效果,驗(yàn)證了本文提出方法的可靠性.