余 燕，周琳瑯，甘笑非，胡 燕，淦文杰，鄧 莊

（中國石油西南油氣田公司川中油氣礦，四川遂寧 629000）

0 引言

碳酸鹽巖儲集空間類型多樣、孔隙結(jié)構(gòu)復(fù)雜，儲層裂縫、孔洞發(fā)育，具有很強(qiáng)的非均質(zhì)性［1］。鄭松青等［2］針對碳酸鹽巖縫洞型儲層呈離散介質(zhì)的特征，提出縫洞網(wǎng)絡(luò)模型，并在此基礎(chǔ)上建立了油水兩相流數(shù)學(xué)模型，簡化了儲層描述參數(shù)，計(jì)算成本低，可以進(jìn)行大規(guī)模的工程應(yīng)用；賈永祿等［3］基于基質(zhì)、微裂縫、溶蝕孔洞和大裂縫組成的四重介質(zhì)，建立了縫洞型四重介質(zhì)的碳酸鹽巖試井模型；熊鈺等［4］針對大尺度的碳酸鹽巖縫洞型油氣藏，考慮縫洞的多尺度特征，建立了井打在大尺度裂縫上的裂縫線性流試井解釋模型；程飛［5］針對縫洞型碳酸鹽巖儲層的非均質(zhì)性和各向異性強(qiáng)的特征，建立了碳酸鹽巖油藏儲層類型識別模板，并劃分了3 種儲層類型；姜瑞忠等［6］針對低滲油藏存在的壓敏效應(yīng)和啟動(dòng)壓力梯度現(xiàn)象，建立了雙重介質(zhì)低滲油藏的斜井?dāng)?shù)學(xué)模型，并繪制了斜井試井曲線；杜鑫等［7］針對縫洞型油氣藏，建立了井筒與溶洞相連通情況下的流動(dòng)方程并給出了井筒和溶洞之間的壓力差表達(dá)式；劉建峰等［8］建立了孔洞和裂縫組成滲流空間的典型雙重介質(zhì)模型；王蓓等［9］針對磨溪龍王廟組縫洞型碳酸鹽巖儲層，建立了多尺度離散裂縫模型，明確了該氣藏高、低滲分布區(qū)域；姜瑞忠等［10］建立了縫洞型碳酸鹽巖油藏水平井分形非線性滲流模型，并結(jié)合測井?dāng)?shù)據(jù)將新模型應(yīng)用于礦場滲流參數(shù)的解釋，但這些有關(guān)碳酸鹽巖儲層的滲流模型均未考慮二次壓力梯度非線性項(xiàng)的影響。

聶仁仕等［11-13］建立了考慮井儲和表皮效應(yīng)的不穩(wěn)定滲流數(shù)學(xué)模型，同時(shí)，該模型考慮了二次壓力梯度非線性項(xiàng)的影響；Nie 等［14］針對底水油藏和厚油層油藏，建立了不同邊界條件下考慮二次壓力梯度非線性項(xiàng)的球面滲流方程；張林等［15］建立了不同邊界條件下考慮二次壓力梯度非線性項(xiàng)的三區(qū)復(fù)合滲流方程；張強(qiáng)等［16］建立了考慮啟動(dòng)壓力梯度和二次梯度非線性項(xiàng)影響的滲流模型；Wang 等［17］針對多層復(fù)合儲層，建立了考慮二次壓力梯度非線性項(xiàng)的多區(qū)復(fù)合非線性模型；Guo 等［18］針對地下地層中氣井生產(chǎn)，建立了考慮二次壓力梯度非線性項(xiàng)的直井模型；王美楠等［19］建立了考慮啟動(dòng)壓力梯度、動(dòng)邊界、二次梯度非線性項(xiàng)的低滲透變形介質(zhì)油藏滲流模型；Lu 等［20］建立了考慮二次壓力梯度非線性項(xiàng)的多孔介質(zhì)水平井滲流模型。目前，針對碳酸鹽巖縫洞型三孔介質(zhì)儲層，尚無有關(guān)二次壓力梯度三孔滲流模型及非線性滲流特征的研究。

流體在地下儲層中的滲流具有二次壓力梯度的非線性滲流特性［18-20］，縫洞型碳酸鹽巖儲層也不例外，采用常規(guī)的三孔介質(zhì)線性滲流模型無法精準(zhǔn)描述流體在縫洞型碳酸鹽巖儲層中的滲流特征。筆者首次建立二次壓力梯度三孔滲流模型，研究縫洞型碳酸鹽巖儲層受二次壓力梯度控制的非線性滲流模型的缺失問題，以期實(shí)現(xiàn)對縫洞型碳酸鹽巖儲層滲流特征的精準(zhǔn)描述。

1 滲流模型

1.1 物理模型

平面徑向流物理模型假設(shè)如下：①裂縫孔洞型三孔介質(zhì)水平等厚地層中心一口井以定產(chǎn)量q生產(chǎn)；②地層厚度為h，地層滲透率為常數(shù)且與壓力無關(guān)，原始地層壓力為pi；③井半徑為rw，井底壓力為pwf；④井到地層外邊界的距離為re，外邊界可以是無窮大地層、封閉地層或恒壓邊界；⑤在彈性驅(qū)動(dòng)方式下，黏度為μ的單相液體滲流符合等溫達(dá)西定律；⑥考慮天然裂縫與井筒溝通，溶洞系統(tǒng)和基質(zhì)系統(tǒng)中流體均向天然裂縫系統(tǒng)發(fā)生擬穩(wěn)態(tài)竄流，因溶洞滲透率大于基質(zhì)滲透率，故溶洞系統(tǒng)中流體優(yōu)先向天然裂縫系統(tǒng)竄流；⑦考慮井筒儲集效應(yīng)和表皮效應(yīng)的影響。

1.2 數(shù)學(xué)模型及其求解

1.2.1 數(shù)學(xué)模型

一個(gè)完整的數(shù)學(xué)模型由滲流控制微分方程、初始條件、內(nèi)外邊界條件構(gòu)成。對縫洞型碳酸鹽巖儲層，其滲流控制微分方程包括裂縫系統(tǒng)方程、基質(zhì)系統(tǒng)方程和溶洞系統(tǒng)方程。

（1）滲流控制微分方程

裂縫系統(tǒng)方程：

式中：pfD為無因次裂縫系統(tǒng)壓力；pmD為無因次基質(zhì)系統(tǒng)壓力；pvD為無因次溶洞系統(tǒng)壓力；rD為無因次徑向距離；β為無因次二次壓力梯度項(xiàng)系數(shù)；λmf為基質(zhì)系統(tǒng)向裂縫系統(tǒng)的竄流系數(shù)；λvf為溶洞系統(tǒng)向裂縫系統(tǒng)的竄流系數(shù)；ωf為裂縫系統(tǒng)的彈性儲容比；S為井的表皮系數(shù)；tD為無因次時(shí)間。下標(biāo)f，m，v，D 分別代表裂縫（fracture）、基質(zhì)（matrix）、溶洞（vug）和無因次（Dimensionless）。

基質(zhì)系統(tǒng)方程：

式中：ωm為基質(zhì)系統(tǒng)的彈性儲容比。

溶洞系統(tǒng)方程：

式中：ωv為溶洞系統(tǒng)的彈性儲容比。

（2）初始條件

（3）內(nèi)邊界條件

式中：CD為無因次井筒儲集系數(shù)。

（4）外邊界條件

無限大外邊界

恒壓邊界

式中：reD為無因次外邊界距離。

封閉邊界

上述數(shù)學(xué)模型的無因次參數(shù)定義如下：無因次壓力

式中：pf為裂縫系統(tǒng)壓力，MPa；pm為基質(zhì)系統(tǒng)壓力，MPa；pv為溶洞系統(tǒng)壓力，MPa；pi為初始油藏壓力，MPa；q為流量，m3/d；μ為黏度，mPa·s；B為地層流體的體積系數(shù)，無因次；kf為裂縫滲透率，D；h為儲層厚度，m。

基于有效井徑定義的無因次距離

式中：r為徑向距離，m；rw為井半徑，m。

無因次時(shí)間

式中：t為時(shí)間，h；φf為裂縫孔隙度；φv為溶洞孔隙度；φm為基質(zhì)孔隙度；Cft為裂縫彈性綜合壓縮系數(shù)，MPa-1；Cvt為溶洞彈性綜合壓縮系數(shù)，MPa-1；Cmt為基質(zhì)彈性綜合壓縮系數(shù)，MPa-1。

無因次井筒儲集系數(shù)

式中：Cs為井的彈性儲容系數(shù)，m3/MPa。

彈性儲容比

無因次二次壓力梯度項(xiàng)系數(shù)

式中：Cρ為地層流體的壓縮系數(shù)，MPa-1。

1.2.2 積分變換求解模型

引入基于tD的拉氏變換

式中：L［］為Laplace 算子為變量pD對應(yīng)的La‐place 空間變量；u為Laplace 變量。

對式（1）—（3）經(jīng)基于tD的Laplace 變換后，得

當(dāng)ωf=1 時(shí)，非線性模型則簡化為常規(guī)單一介模型，式（22）可以簡化為質(zhì)

換句話講，在相同的油井生產(chǎn)和地層邊界條件下，Laplace 空間下的三孔介質(zhì)模型和均質(zhì)模型的差異僅僅表現(xiàn)在f（u）的函數(shù)表達(dá)式不同。因此，如果首先在Laplace 空間中得到均質(zhì)模型的解，就可以直接通過改變f（u）來得到Laplace 空間下三孔介質(zhì)模型的解。

先令ωf=1，則考慮單一介質(zhì)情形，滲流控制方程中就不存在竄流項(xiàng)，則式（1）變?yōu)?/p>

作如下變量代換，將式（24）線性化

式中：x為無因次裂縫壓力的中間代換變量。

再經(jīng)基于tD的Laplace 變換后，可得到Laplace空間中的數(shù)學(xué)模型

1.2.3 數(shù)學(xué)模型的解

用式（22）替換式（27），就可得到經(jīng)變量代換線性化后的裂縫孔洞型三孔介質(zhì)滲流偏微分控制方程

（1）無窮大外邊界模型的解

故無窮大外邊界條件模型在Laplace 空間的通解為

井底rD=1，當(dāng)p=pwf時(shí)，pD=pwD，則x=xw，故無窮大地層模型井底動(dòng)態(tài)壓力在Laplace 空間為

（2）恒壓外邊界模型的解

將式（34）代入恒壓外邊界條件式（30），有

聯(lián)解式（35），（36），（39），將其寫為矩陣形式，利用克萊姆法則或高斯迭代消元法，即可解出系數(shù)A和B，從而得出外邊界條件為恒壓邊界時(shí)井底壓力動(dòng)態(tài)在Laplace 空間的解。

（3）封閉外邊界模型的解

將式（34）代入封閉外邊界條件式（31），有

聯(lián)解式（35），（36），（40），將其寫為矩陣形式，利用克萊姆法則或高斯迭代消元法，即可解出系數(shù)A和B，從而得出外邊界條件為封閉邊界時(shí)井底壓力動(dòng)態(tài)在Laplace 空間的解。

對上述三類外邊界條件下Laplace 空間的解，利用Stehfest 數(shù)值算法進(jìn)行數(shù)值反演，得到實(shí)空間的解xw，再根據(jù)壓力pD與x的變換關(guān)系式，可求出裂縫孔洞型三孔介質(zhì)非線性滲流模型實(shí)空間井底壓力動(dòng)態(tài)響應(yīng)的數(shù)值解，進(jìn)而編程繪制pwD～tD/CD與p'wDtD/CD～tD/CD的無因次雙對數(shù)滲流特征曲線。

2 非線性滲流特征分析

2.1 非線性滲流過程分析

圖1 為受非線性項(xiàng)系數(shù)影響的無窮大裂縫孔洞型三孔介質(zhì)地層擬穩(wěn)態(tài)竄流非線性滲流模型的典型曲線。該曲線主要受無因次二次梯度非線性項(xiàng)系數(shù)β控制。在其他參數(shù)一定時(shí)（CD=10，S=10，ωf=0.001，ωm=0.989，ωv=0.01，λmf=1×10-10，λvf=1×10-8），通過分別設(shè)定β值為0，0.01 和0.05，來分析非線性滲流特征曲線與線性滲流特征曲線的差異。當(dāng)β=0 時(shí)，該曲線則蛻化為線性模型的滲流典型曲線（圖1 曲線①）。從圖1 可以看出，非線性滲流特征曲線向下偏離線性滲流特征曲線，隨β增大，偏離量越大。圖1 主要有如下7 個(gè)流動(dòng)階段：

第Ⅰ階段：純井筒儲集階段。此流動(dòng)階段的曲線特征不受非線性二次壓力梯度項(xiàng)的影響。

第Ⅱ階段：井筒儲集和表皮效應(yīng)階段。在此階段，非線性二次壓力梯度項(xiàng)的影響開始出現(xiàn)，非線性滲流特征曲線也開始偏離線性滲流特征曲線。

第Ⅲ階段：早期裂縫徑向流階段。裂縫系統(tǒng)中的流體徑向流入井筒，而溶洞系統(tǒng)和基質(zhì)系統(tǒng)中的流體還未開始流動(dòng)。與第Ⅱ階段相比，線性模型和非線性模型典型曲線的差異更為明顯。

圖1 受非線性項(xiàng)系數(shù)（β）影響的三孔介質(zhì)擬穩(wěn)態(tài)竄流模型典型曲線Fig.1 Type curves of pressure transients controlled by varying nonlinear coefficient（β）

第Ⅳ階段：溶洞系統(tǒng)向裂縫系統(tǒng)竄流階段。壓力導(dǎo)數(shù)曲線呈“V 形”，反映了溶洞系統(tǒng)向裂縫系統(tǒng)的竄流。線性模型和非線性模型典型曲線的差異進(jìn)一步增加。

第Ⅴ階段：裂縫系統(tǒng)和溶洞系統(tǒng)整體徑向流階段。非線性二次壓力梯度項(xiàng)的影響越來越明顯，其壓力導(dǎo)數(shù)曲線不符合“0.5 線”定律，而是位于“0.5線”下方。

第Ⅵ階段：基質(zhì)系統(tǒng)向裂縫系統(tǒng)竄流階段。壓力導(dǎo)數(shù)曲線呈“V 形”，表明基質(zhì)系統(tǒng)開始向裂縫系統(tǒng)竄流。線性模型和非線性模型特征曲線的差異越來越大。

第Ⅶ階段：裂縫系統(tǒng)、溶洞系統(tǒng)和基質(zhì)系統(tǒng)整體徑向流階段。非線性二次壓力梯度項(xiàng)的影響越來越明顯，其壓力導(dǎo)數(shù)曲線不符合“0.5 線”定律，而是位于“0.5 線”下方，隨著時(shí)間的推移，壓力導(dǎo)數(shù)曲線逐漸偏離“0.5 線”，如圖1 中曲線②和曲線③所示。隨著時(shí)間的推移，非線性模型的壓力導(dǎo)數(shù)曲線逐漸偏離線性模型的壓力導(dǎo)數(shù)曲線。

非線性三孔介質(zhì)曲線的外邊界反映特征與常規(guī)線性模型曲線的外邊界反映特征類似，此處不再贅述。

該曲線除了受β控制，同時(shí)還受別的因素影響。圖2 為受溶洞向裂縫竄流的竄流系數(shù)（λvf）的非線性滲流模型典型曲線。在其他參數(shù)一定時(shí)（CD=10，S=1，ωf=0.001，ωm=0.989，ωv=0.01，β=0.05），分別設(shè)定λvf值為1×10-7，1×10-8和1×10-9。隨λvf增加，溶洞向裂縫竄流時(shí)間推遲，在圖2 中表現(xiàn)為壓力導(dǎo)數(shù)曲線的“V 形”向左移動(dòng)。

圖2 受溶洞向裂縫竄流的竄流系數(shù)（λvf）影響的典型曲線Fig.2 Type curves of pressure transients controlled by varying inter-porosity flow factor（λvf）

圖3 為受基質(zhì)向裂縫竄流的竄流系數(shù)（λmf）的非線性滲流模型典型曲線。在其他參數(shù)一定時(shí)（CD=10，S=1，ωf=0.001，ωm=0.989，ωv=0.01，β=0.05），分別設(shè)定λmf值為1×10-9，1×10-10和1×10-11。溶洞向裂縫竄流時(shí)間推遲，在圖3 中表現(xiàn)為壓力導(dǎo)數(shù)曲線的“V 形”向右移動(dòng)。

圖3 受基質(zhì)向裂縫竄流的竄流系數(shù)（λmf）影響的典型曲線Fig.3 Type curves of pressure transients controlled by varying inter-porosity flow factor（λmf）

圖4 為受裂縫系統(tǒng)彈性儲容比（ωf）、基質(zhì)系統(tǒng)彈性儲容比（ωm）影響的非線性滲流模型典型曲線。其他參數(shù)一定時(shí)（CD=10，S=1，ωv=0.1，λvf=1×10-10，λmf=1×10-8，β=0.05），分別設(shè)定ωf值為0.001，0.01 和0.08，相應(yīng)的ωm值為0.889，0.89 和0.82。隨ωf的增加，在圖4 中表現(xiàn)為壓力導(dǎo)數(shù)曲線的“V 形”變淺、變窄。

圖5 為受溶洞系統(tǒng)彈性儲容比（ωv）、基質(zhì)系統(tǒng)彈性儲容比（ωm）影響的非線性滲流模型典型曲線。其他參數(shù)一定時(shí)（CD=10，S=1，ωf=0.001，λvf=1×10-10，λmf=1×10-8，β=0.05），分別設(shè)定ωv值為0.01，0.01 和0.5，相應(yīng)的ωm值為0.989，0.889 和0.499。隨ωv的增加，在圖5 中表現(xiàn)為壓力導(dǎo)數(shù)曲線的第1個(gè)“V 形”變得更深、更寬，第2 個(gè)“V 形”變得更窄、更淺。

圖4 受裂縫系統(tǒng)彈性儲容比（ωf）影響的典型曲線Fig.4 Type curves of pressure transients controlled by fluid capacitance coefficient（ωf）

圖5 受溶洞系統(tǒng)彈性儲容比（ωv）影響的典型曲線Fig.5 Type curves of pressure transients controlled by fluid capacitance coefficient（ωv）

2.2 非線性滲流模型典型曲線的偏差分析

為了定量描述非線性項(xiàng)系數(shù)的影響，定義絕對偏差（DV）與相對偏差（RDV）：

其他參數(shù)一定時(shí)（CD=10，S=1，ωf=0.001，ωm=0.989，ωv=0.01，λmf=1×10-10，λvf=1×10-8），分別設(shè)定β值為0.01 和0.05，分析非線性項(xiàng)系數(shù)對滲流典型曲線的影響（表1—2）。

當(dāng)β＝0.01 時(shí)，由表1 中數(shù)據(jù)知：當(dāng)tD/CD＝103時(shí)，無因次壓力曲線的絕對偏差和相對偏差分別為0.579 8 和8.83%，無因次壓力導(dǎo)數(shù)曲線的絕對偏差和相對偏差分別為0.085 4 和17.37%；當(dāng)tD/CD＝107時(shí)，無因次壓力曲線的絕對偏差和相對偏差分別為1.4 和13.43%，無因次壓力導(dǎo)數(shù)曲線的絕對偏差和相對偏差分別為0.089 7 和25.40%。

表1 非線性與線性滲流典型曲線的定量偏差數(shù)據(jù)（β＝0.01）Table 1 Theoretical offset of type curves between the linear and nonlinear models（β＝0.01）

表2 非線性與線性滲流典型曲線的定量偏差數(shù)據(jù)（β＝0.05）Table 2 Theoretical offset of type curves between the linear and nonlinear models（β＝0.05）

當(dāng)β＝0.05，由表2 中數(shù)據(jù)知：當(dāng)tD/CD＝103時(shí)，無因次壓力曲線的絕對偏差和相對偏差分別為2.135 0 和32.52%，無因次壓力導(dǎo)數(shù)曲線的絕對偏差和相對偏差分別為0.272 4 和55.40%；當(dāng)tD/CD＝107時(shí)，無因次壓力曲線的絕對偏差和相對偏差分別為4.527 2 和43.43%，無因次壓力導(dǎo)數(shù)曲線的絕對偏差和相對偏差分別為0.241 1 和68.28%。

絕對偏差和相對偏差均隨著無因次時(shí)間的增加而增加。在同一時(shí)刻，壓力導(dǎo)數(shù)的相對偏差值大于壓力的相對偏差值，例如，當(dāng)tD/CD＝103，β＝0.01時(shí)，壓力導(dǎo)數(shù)的相對偏差值就（17.37%）大于壓力的相對偏差（8.83%）。絕對偏差和相對偏差值均隨β的增加而增加，例如，當(dāng)tD/CD＝103時(shí)，β＝0.01的壓力相對偏差值為8.83%，而β＝0.05 的壓力相對偏差值為32.52%。

3 結(jié)論

（1）針對裂縫孔洞型三孔介質(zhì)儲層建立的二次壓力梯度非線性滲流模型比不考慮二次壓力梯度的常規(guī)線性滲流模型能更精準(zhǔn)地描述儲層的滲流特征。

（2）隨油井生產(chǎn)時(shí)間的增加，裂縫孔洞型三孔介質(zhì)儲層非線性滲流特征曲線逐漸偏離常規(guī)線性滲流特征曲線，且低于常規(guī)線性滲流特征曲線。

（3）非線性與線性滲流特征曲線之間的偏離量可用“絕對偏差”和“相對偏差”等2 個(gè)參數(shù)來定量描述，其偏差隨著二次壓力梯度非線性系數(shù)的增加而增加。當(dāng)非線性系數(shù)為0 時(shí)，非線性三孔滲流模型可蛻化為常規(guī)線性三孔滲流模型。