非均質(zhì)低滲透油藏非線性流動數(shù)值模擬

姚 軍， 黃 濤， 黃朝琴

（中國石油大學（華東）石油工程學院，山東青島 266580）

非均質(zhì)低滲透油藏非線性流動數(shù)值模擬

姚 軍， 黃 濤， 黃朝琴

（中國石油大學（華東）石油工程學院，山東青島 266580）

提出一種適用于非均質(zhì)低滲透油藏非線性流動數(shù)值模擬的方法.通過引入單元邊界壓力和速度兩組物理量，在單個網(wǎng)格單元中構造數(shù)值計算格式，可精確地獲取每個網(wǎng)格單元的壓力梯度.以單相非線性滲流問題為出發(fā)點，詳細闡述該方法的基本原理；通過數(shù)值算例驗證方法的正確性和適用性.

非均質(zhì)；低滲透油藏；非線性滲流；油藏數(shù)值模擬

0 引言

隨著我國油氣田勘探開發(fā)的不斷深入，低滲透油藏已成為我國增儲上產(chǎn)的主要領域［1］.大量低滲透滲流實驗和生產(chǎn)實踐表明［2－4］：油氣流動存在最小啟動壓力梯度，其滲流曲線偏離了經(jīng)典的達西定律，具有顯著的非線性特征.許多專家學者對其流動機理進行了大量研究，提出了多種非線性滲流模式，包括擬啟動壓力梯度模型、三階段模型、非線性連續(xù)性模型等［5－7］，并進行了相應的數(shù)值模擬研究［8－11］.現(xiàn)有的低滲透油藏非線性流動數(shù)值模擬，大都是采用傳統(tǒng)的有限差分法或控制體積有限差分.此類方法均基于相鄰單元中心點的物理量值來構造數(shù)值計算格式，求取驅替壓力梯度，進而結合滲流曲線來進行判斷和求解.然而，對于非均質(zhì)低滲透油藏，不同區(qū)域甚至不同網(wǎng)格單元的滲流模式不盡相同，數(shù)值模擬中必須獲取各網(wǎng)格單元的壓力梯度，并選擇相應的滲流模式，此時上述傳統(tǒng)的數(shù)值計算方法已無法滿足要求.

對此，在控制體積有限差分法的基礎上，通過引入界面中心點壓力，提出了一種新的數(shù)值算法.該算法基于單個網(wǎng)格單元來構造數(shù)值計算格式，能夠精確地獲取每一網(wǎng)格單元的壓力梯度，單元間則通過流量連續(xù)條件聯(lián)系起來.并以此為基礎，建立了一套非均質(zhì)低滲透油藏非線性流動數(shù)值模擬方法.

1 滲流數(shù)學模型

由于受啟動壓力梯度的影響，低滲透油藏流體流動存在明顯的動邊界［12－13］.實驗結果表明［14］：低滲透油藏滲流模式一般為如圖1所示的滲流曲線（實線OA－AE），A點表示孔隙介質(zhì)中最大孔道中的流體流動時的啟動壓力梯度，B點表示最小孔道中流體流動時所需的啟動壓力梯度，曲線AD上一點Di的平均啟動壓力梯度為Ci，C為直線DE與橫坐標軸的交點，即擬啟動壓力梯度，D為曲線轉直線的折點.

之前很多研究采用簡化的擬啟動壓力梯度滲流模式，即用直線OC、CE來描述滲流過程.擬啟動壓力梯度模型在計算時對于在低壓力梯度時阻力較小的大孔道中的流動估計偏低，因而綜合經(jīng)濟技術指標會偏低.考慮某一封閉低滲透油藏區(qū)域Ω中的不可壓縮單相液體的等溫滲流過程，為簡單起見，本文僅研究各向同性問題，并假設巖石微可壓縮，對于擬啟動壓力梯度其數(shù)學模型如下：

圖1 典型非達西滲流曲線Fig.1 Non-Darcy flow curve

其中，υ為滲流速度，K為滲透率張量，μ為流體粘度，p為油藏壓力（或流動勢），λ為啟動壓力梯度，φ為巖石孔隙度，q表示源匯項.

對于實際非線性滲流曲線（圖1實線部分），文獻［15－17］建議采用動態(tài)滲透率的概念，即滲透率和啟動壓力梯度是隨著壓力梯度的增大而增大的，當壓力梯度大于某一臨界值時滲透率和啟動壓力梯度不再發(fā)生變化，如圖2所示.此時，可求出每個驅替壓力梯度對應的滲透率及啟動壓力梯度，具體如下：①由巖心做出低滲介質(zhì)滲流曲線圖；②由驅替壓力梯度，找到對應于滲流曲線上的滲流點；③對該滲流點沿滲流曲線方向做切線，該切線的斜率即為滲透率，而切線與驅替壓力梯度坐標軸的交點即為該壓力梯度時孔道的啟動壓力梯度.顯然，對于動態(tài)滲透率模型，每一個壓力梯度值對應一個擬啟動壓力梯度模型，因此，本文將重點闡述擬啟動壓力梯度滲流數(shù)學模型及其數(shù)值求解方法.對于圖2實線所示的非線性滲流模式，其滲流數(shù)學模型僅需改變方程（1），其中，在曲線AD段滲透率是驅替壓力梯度函數(shù)，具體如下

圖2 滲透介質(zhì)的滲透率隨壓力梯度的變化Fig.2 Mutative permeability vs.pressuregradient

2 數(shù)值計算格式構造

假設油藏區(qū)域Ω∈R2被一組互不重疊的多邊形網(wǎng)格Ω＝｛Ωi｝所剖分，如圖3所示.

圖3 網(wǎng)格單元分析示意圖Fig.3 Schematic of grid analysis

由方程（1）易知邊界面上的法向滲流速度vf與壓力梯度Δp、啟動壓力梯度λ之間有

然，方程（15）中的第一行對應于方程（1）Darcy定律，第二行對應于方程（2）質(zhì)量守恒定律，第三行則表征了單元邊界面上的法向速度連續(xù)性條件.上述方程的系數(shù)矩陣表達式：

其中，Ne為網(wǎng)格單元總數(shù)；Ii＝Em．從方程（16）可知：方程（15）中的各系數(shù)矩陣僅與網(wǎng)格單元的幾何信息和油藏參數(shù)有關，而對網(wǎng)格的幾何形狀沒有特殊要求，求解方便，原則上適用于任何復雜網(wǎng)格.并且方程（14）中啟動壓力項只與每個網(wǎng)格單元有關，對于非均質(zhì)油藏只需改變每個單元啟動壓力梯度，所以具有很好的適應性.

3 數(shù)值算例

通過簡單的二維單相流算例對低滲透油藏數(shù)值模擬進行分析.考慮如圖4（a）所示二維1／4五點井網(wǎng)模型，油藏長100 m，寬100 m，在（0，0）點有一注入井，在（100，100）點有一采出井，注入和采出的速度均為17.28 m3·d－1.采用三角網(wǎng)格剖分，節(jié)點數(shù)為190，單元數(shù)為330，如圖4（b）所示.

圖4 油藏幾何信息圖Fig.4 Reservoir geometry information

3.1 算例1

擬啟動壓力梯度模型：設滲透率Kx＝Ky＝0.1 md（10－3μm2），Kxy＝Kyx＝0，初始油藏壓力Pe＝Pf＝10 MPa，擬啟動壓力梯度λ＝0.1 MPa·m－1，巖石壓縮系數(shù)Cf＝8.64×10－4MPa－1，流體粘度μ＝1 mPa·s，t＝200 d（天）.

從圖5可以看出，在擬啟動壓力梯度模型滲流過程中，存在明顯的動邊界，當驅替壓力小于啟動壓力時，油藏中的流體是不流動的，油藏壓力保持初始值不變.當壓力傳播到整個油藏時，會形成優(yōu)勢通道，流體流速較快.

圖5 算例1壓力及速度分布Fig.5 Distributions of pressure and velocity in Example 1

3.2 算例2

非線性滲流模型：該模型與上述擬啟動壓力梯度模型基本相同，唯一區(qū)別是，滲透率和啟動壓力梯度隨驅替壓力梯度增大而增大，當達到某一臨界值時保持不變，其對應關系見表1.圖6為非線性滲流模型不同時刻壓力及速度分布圖.

圖6 算例2壓力及速度分布Fig.6 Distributions of pressure and velocity in Example 2

如圖7所示，對比兩個模型不同時間步t＝80、160、200 d注采井點連線之間壓力分布，可看出擬啟動壓力梯度模型在同一時間步中，壓力傳播范圍較非線性滲流模型小，且靠近源點處相對于初始地層壓力壓差大；而在非線性滲流模型中，壓力傳播較快，油藏較早達到穩(wěn)態(tài).這是因為擬啟動壓力梯度模型是一種簡化的滲流模式，模型假設了一個啟動壓力梯度（圖1中的C點），當驅替壓力梯度小于該啟動壓力梯度時，地層中的流體是不流動的.而實驗結果表明實際滲流模式應為圖1實線所示，當驅替壓力梯度大于某值（圖1中的A點）時，流體已經(jīng)開始流動.所以說，對比2類模型，由于流體流動所需要的啟動壓力大，擬啟動壓力梯度模型在計算時對于在低壓力梯度時阻力較小的大孔道中的流動估計偏低，因而綜合經(jīng)濟技術指標會偏低，非線性滲流模型近似圖1實線部分所示滲流模式，計算時更符合實際.

圖7 不同時刻兩種模型井間連線壓力值Fig.7 Pressures between two model's wells

3.3 算例3

非線性滲流模式算例中均采用了相同的滲透率、啟動壓力梯度隨驅替壓力梯度變化曲線，對于非均質(zhì)油藏，由于計算過程中針對每個網(wǎng)格單元均求出了驅替壓力梯度，所以只要給定區(qū)塊甚至網(wǎng)格單元的變化曲線，都可以求出網(wǎng)格單元相應的滲透率以及啟動壓力梯度，這也是該方法的優(yōu)勢所在，為了簡單起見，均采用同一變化曲線，見表1.

考慮與算例2中非線性滲流相同的模型，唯一區(qū)別是其滲透率分布不同，在此算例中滲透率Kx、Ky通過MatLab的隨機函數(shù)randn生成，Kxy＝Kyx＝0，如圖8（a）所示，圖8（b）、8（c）、8（d）為不同時刻壓力及速度分布.

低滲透油藏非均質(zhì)性極強，網(wǎng)格之間滲透率差別較大，如圖8（a）所示.從數(shù)值模擬結果看出，流體主要沿著滲透率較大的網(wǎng)格流動、流速較大，壓力傳播快，符合實際，并且本文提出的方法適用于非均質(zhì)性油藏.

圖8 算例3壓力及速度分布Fig.8 Distributions of pressure and velocity in Example 3

4 結論

1）提出一種適用于非均質(zhì)低滲透油藏流動數(shù)值模擬方法，以單相滲流問題為出發(fā)點，詳細闡述了方法的基本原理，并通過數(shù)值算例驗證了方法的正確性.

2）新方法僅基于單個網(wǎng)格單元的信息來構造數(shù)值計算格式，理論上適用于任何復雜網(wǎng)格，可精確地獲取每個網(wǎng)格單元的壓力梯度，適用于非均質(zhì)低滲透油藏數(shù)值模擬，尤其是不同區(qū)域或網(wǎng)格單元具有完全不同滲流模式的情形.

3）本文僅對非均質(zhì)低滲透油藏的單相滲流問題予以研究，下一步將該方法推廣至兩相、多相滲流以及三維問題.

［1］Wang Guangfu，Liao Rongfeng，Li Jianglong，et al.The development situation and future of low permeability oil reservoirs of SINOPEC［J］.Petroleum Geology and Recovery Efficiency，2007，14（3）：84－89.

［2］Miller R J，Low P F.Threshold gradient for water flow in clay systems［J］.Soil Science Society of America Proceedings，963，27（6）：605－609.

［3］Lv Chengyuan，Wang Jian，Sun Zhigang，et al.An experimental study on starting pressure gradient of fluids flow in lowpermeability sandstone porous media［J］.Petroleum Exploration and Development，2002，29（2）：86－89.

［4］Wan Quanlin，Tan Hai，Lv Dongliang，et al.An experimental study on threshold pressure gradient in low permeability reservoir［J］.Petroleum Geology and Recovery Efficiency，2011，18（1）：97－100.

［5］Huang Yanzhang.Nonlinear percolation feature in low permeability reservoir［J］.Special Oil＆Gas Reservoirs，1997，4（1）：9－14.

［6］Yao Yuedong，Ge Jiali.Study on non Darcy flow pattern in low permeability oil reservoir［J］.Xinjiang Petroleum Geology，2000，21（3）：213－215.

［7］Xiao Zengli，Pu Chunsheng，Qin Wenlong，et al.Nonlinear percolation feature of low permeability reservoirs and its effect ［J］.Oil Drilling＆Production Technology，2007，29（3）：105－107.

［8］Guo Yongcun，Lu Detan，Ma Linxiao，et al.Numerical simulation of fluid flow in low permeability reservoir using finite difference method［J］.Journal of Hydrodynamics，2004，19（3）：288－293.

［9］Han Hongbao，Cheng Linsong，Zhang Minglu，et al.Physical simulation and numerical simulation of ultra-low permeability reservoir in consideration of starting pressure gradient［J］.Journal of University of Petroleum，China（Natural Science Edition），2004，28（6）：49－53.

［10］Xue Guoqing，Li Ming，Wang Daocheng，et al.Numerical simulation of low permeability and heterogeneous reservoir［J］. Fault-Block Oil and Gas Field，2007，14（4）：43－44.

［11］Yang Zhengmin，Yu Rongcheng，Su Zhixin，et al.Numerical simulation of the nonlinear flow in ultra-low permeability reservoirs［J］.Petroleum Exploration and Development，2010，7（1）：94－98.

［12］Guo Yongcun，Zeng Qinghong，Wang Zhongxun.Numerical simulation of low permeability flow with moving-boundary using meshless methods［J］.Engineering Mechanics，2006，23（11）：188－192.

［13］Liu Wenchao，Yao Jun，Wang Jianzhong.Front tracking for moving boundary of non-Darcy seepage flows in low-permeability porous media［J］.Chinese Journal of Computational Physics，2012，29（6）：823－827.

［14］Yan Qinglai，He Qiuxuan，Yu Ligang.A laboratory study on percolation characteristics of single phase flow in low-permeability reserviors［J］.Journal of Xi'an Petroleum Institute，1990，5（2）：1－6.

［15］Lei Qun，Xiong Wei，Yuan Jiangru，et al.Behavior of flow through low-permeability reservoirs［A］.SPE 113144，2008.

［16］Yao Jun，Liu Shun.Well test interpretation model based on mutative permeability effects for low-permeability reservoir［J］. Acta Petrolei Sinica，2009，30（3）：430－433.

［17］Yin Zhilin，Sun Wenjing，Yao Jun.Numerical simulation of the 3D oil-water phase dynamic permeability for low-permeability reservoirs［J］.Acta Petrolei Sinica，2011，32（1）：117－121.

Numerical Simulation of Nonlinear Flow in Heterogeneous and Low-permeability Reservoirs

YAO Jun，HUANG Tao，HUANG Zhaoqin
（School of Petroleum Engineering，China University of Petroleum，Qingdao 266580，China）

A numerical method is developed to model nonlinear flow in heterogeneous and low-permeability reservoirs.In order to obtain pressure gradient of each gridcells precisely，we introduced two sets of physical quantities，cell face pressures and surface efflux for each gridcell.Principle of the method is elaborated based on single-phase fluid flow problem.Several numerical examples are presented to demonstrate efficiency and applicability of the scheme.

heterogeneity；low-permeability reservoirs；nonlinear porous flow；numerical simulation

date：2013－10－15；Revised date：2014－02－19

TE348；TE319；O241.82

A

2013－10－15；

2014－02－19

國家重點基礎研究發(fā)展計劃（2011CB201004），國家自然科學基金（51234007），國家重大科技專項（2011ZX05014-005-003HZ）和中央高?；究蒲袠I(yè)務費專項資金（13CX05007A，14CX02042A，13CX05017A，14CX06091A）資助項目

姚軍（1964－），男，山東臨沂，博士生導師，主要從事油氣田開發(fā)工程教研工作，E-mail：RCOGFR_UPC＠126.com

1001-246X（2014）05-0552-07