李彥澤 段 彬 商 琳，2，3 孫彥春 龔麗榮 辛春彥 邢建鵬

（1.中國石油冀東油田公司，河北 唐山 063000；2.中國石油冀東油田公司博士后科研工作站，河北 唐山 063000；3.中國石油勘探開發(fā)研究院博士后流動站，北京 100083；4.中國石油大學（北京）安全與海洋工程學院，北京 102200）

0 引 言

油藏數(shù)值模擬技術自20世紀50年代誕生，到21世紀的今天已發(fā)展為一項較為成熟的常用技術，作為一種定量預測技術，具有精度高、三維可視化、能夠直觀揭示油氣流在儲層內的滲流過程等特點。在理論研究上，學者將其應用于探索多孔介質儲層的復雜滲流問題；在油田開發(fā)中，數(shù)值模擬技術也廣泛應用于開發(fā)方案、開發(fā)調整方案、油氣藏提高采收率方案等，為油田開發(fā)技術政策的制定提供了重要支撐［1-6］。目前，Eclipse、CMG、VIP 等是目前應用比較廣泛的商業(yè)數(shù)值模擬軟件，都以達西滲流模型為核心，在中高滲透油藏中成熟應用。流體在特低滲透、致密砂巖等多孔介質儲層中的流動，由于邊界效益表現(xiàn)為存在最小啟動壓力梯度的低速非線性滲流，達西線性滲流難以表征此類儲層中流體的滲流過程［7-9］。因此，目前市面上常用的商業(yè)數(shù)值模擬軟件在特低滲透、致密砂巖油藏數(shù)值模擬中應用誤差相對增大。

自20世紀90年代以來，眾多學者開展了低滲透、致密砂巖儲層非線性滲流規(guī)律研究，也探索了相應的數(shù)值模擬方法［10-16］，主要認識有3 個方面：一是普遍認識到流體在特低滲透、致密砂巖等多孔介質儲層中的流動，存在最小啟動壓力梯度，且為低速非達西滲流，達西滲流模型在一定程度上夸大了儲層的滲流能力；二是部分學者提出的擬啟動壓力梯度模型，不能表征低速非線性滲流階段，表現(xiàn)出儲層滲流阻力增大，造成采用擬啟動壓力梯度模型的數(shù)值模擬軟件的預測結果誤差較大；三是部分學者提出變滲透率的數(shù)值模擬方法，該方法可以準確地描述特低滲透、致密儲層多孔介質的非達西滲流，但未能給出連續(xù)光滑的滲透率狀態(tài)方程，而且現(xiàn)有的變滲透率數(shù)值模擬軟件沒有能夠給出壓力梯度場和滲透率場的變化規(guī)律，因此變滲透率數(shù)值模擬方法未得到很好地應用。

為解決非線性滲流模擬問題，通過理論模型推導，建立非線性滲流模型，并求解非線性滲流系數(shù)。鑒于tNavigator 中的啟動壓力梯度限制為定值，無法體現(xiàn)非線性滲流規(guī)律，本文研究考慮對滲透率進行時變，即通過對滲透率進行修正來反映啟動壓力梯度的變化，從而反映非線性滲流?，F(xiàn)有軟件無法實現(xiàn)壓力梯度的網格表征，通過前期網格篩選，在定壓條件下進行模擬。

1 非線性滲流模型構建

原油成分極其復雜，并不是單純的牛頓流體，其存在屈服應力，而孔喉的微小使得原油的非牛頓性質更強，不可以忽略。以經典的毛管模型為例，流體勻速通過毛管時力學平衡關系式為

式中：Δp——毛管兩端的驅動壓差，MPa；

r——毛管半徑，μm；

τ0——流體屈服應力，MPa；

μ——流體黏度，mPa·s；

v——流體滲流速度，cm/s；

L——毛管長度，μm。

Δpπr2是驅動力，是內摩擦力，τ0體現(xiàn)出流體的非牛頓性對流體滲流的影響計算公式為

式中 ?p——驅動壓力梯度，MPa/m。

原油中的極性成分在滲流通道固相表面的吸附，導致邊界層的存在是必然的。從前述實驗可知，低滲透油藏孔喉細微，比表面大，邊界層的影響是不可以忽略的，考慮邊界層后流量的計算公式為

式中：q——單根毛細管流量，cm3；

δ——邊界層厚度，μm。

δ體現(xiàn)出固體與液體分子、液體與液體分子之間相互作用對流體滲流的影響；油藏巖心相當于一組半徑均為r的平行毛細管埋置于固體之中。如果與流動方向相垂直的每單位橫截面面積上有N根這樣的毛細管，則通過該多孔介質塊的流量計算公式為［17-18］

式中：Q——截面流量，cm3；

K——滲透率，10-3μm2；

A——截面積，cm2。

式（4）相對于達西定律更符合實際情況，但形式太過復雜，多位學者研究了邊界流體對低滲透、致密儲層滲流規(guī)律的影響，認為邊界層厚度隨驅動壓力梯度的增大而減小。龍騰［15］測定了去離子水通過微管的流量與驅替壓力梯度之間的關系，并計算出去離子水在微管中邊界層厚度隨驅替壓力梯度的變化曲線。閆慶來等［5］指出隨著驅替壓力梯度的增大，邊界層流體厚度呈指數(shù)減薄，說明邊界層厚度是驅替壓力梯度的函數(shù)。理論與實驗證實驅替壓力梯度越大，儲層的邊界層越薄。由前述邊界層實驗有邊界層比例與壓力梯度呈冪函數(shù)關系，與黏度呈線性關系，對于單根毛細管，有

式中a，b——擬合系數(shù)，常量。

式（5）非線性化較高，對于后面模型推廣與實用較差，因此，通過進一步結合實驗數(shù)據，對于壓力梯度采用反比例函數(shù)表示，即

由邊界層厚度定量表征可知，在壓力梯度過大時，邊界層比例減少較小，趨于穩(wěn)定，并不完全符合整體函數(shù)關系，去除掉壓力梯度過大后的點所得到的2 種表征方式系數(shù)與誤差如表1所示，轉換為反比例函數(shù)后誤差為15.17%，相對于冪函數(shù)關系12.77%的誤差，誤差增加不大，且可極大減小模型復雜程度。因此，在后續(xù)推導過程中使用簡化過后的邊界層厚度定量表征公式。

表1 邊界層公式簡化前后誤差對比Table 1 Error comparison of boundary layer formula before and after simplification

忽略高階無窮小項后，有

最終有

2 非線性滲流數(shù)值模擬方法

外部軟件計算模型擬合參數(shù)，非線性滲流模型里的a1、a2兩個參數(shù)，對于不同類型儲層，值均不同，通過現(xiàn)有的巖心數(shù)據，利用matlab 進行非線性擬合。通過對研究區(qū)塊進行全面、系統(tǒng)、綜合研究，應用油氣勘探、開發(fā)過程中獲取的地震、測井、測試、鉆井以及開發(fā)動態(tài)等資料，建立更加貼近于客觀沉積演化規(guī)律的精細化三維地質模型。對油藏儲層進行精細描述和預測，為非線性滲流模擬提供可靠的地質依據。對于采用的非線性滲流模擬方法，滲透率場和壓力場有直接影響，防止場數(shù)據的極值化，將滲透率場與壓力場初始化。對于井周網格進行平面方向的逐個篩選，并將其賦值，便于直接調用。將建立的非線性滲流模型通過滲透率時變的方式導入，輸出考慮非線性滲流的模擬結果，并與其他模型比較。導入歷史數(shù)據，通過不同算法對模擬結果進行智能歷史擬合，反饋數(shù)據，調整非線性滲流參數(shù)（圖1）。

圖1 致密砂巖油藏非線性滲流模擬流程Fig.1 Workflow of nonlinear flow simulation in tight sandstone reservoir

3 應用實例

3.1 模擬流程

由于現(xiàn)場壓力梯度測量較為困難，因此根據實驗室已開展的巖心驅替實驗數(shù)據進行非線性滲流參數(shù)的擬合，根據儲層分類結果，選取冀東油田南堡403X1 斷塊典型井組進行計算，并選取滲透率在1×10-3～3×10-3μm2內的巖心數(shù)據進行非線性擬合。南堡403X1 斷塊Ed2儲層為一套辮狀河三角洲前緣砂體沉積，儲層平均滲透率為2.8×10-3μm2，屬典型的特低滲油藏和層狀構造巖性油藏，采用壓裂注水開發(fā)技術進行開發(fā)，經過近10 a 滾動開發(fā)，已經進入高含水階段，注入水沿人工裂縫竄流嚴重，剩余油分布復雜，急需探索剩余油分布規(guī)律，為開發(fā)調整提供依據。結合鉆井資料、物性參數(shù)、生產動態(tài)等資料對M 區(qū)塊進行屬性建模，平均滲透率為1.29×10-3μm2，平均孔隙度為11.66%，平均深度為3 324.2 m，含油飽和度為46.2%，含水飽和度為53.8%。

3.1.1 初始化場圖

根據推導的非線性滲流模型，在采取滲透率時變的過程中，時變的滲透率是每一個時間步在上一個時間步的基礎上進行倍乘，不加以控制會導致滲透率無限倍乘出現(xiàn)極值，需對其進行初始化，防止極值情況出現(xiàn)。利用滲透率場初始化、壓力場初始化模型（圖2），使每步調用初始滲透率參與時變計算。

圖2 初始化模型Fig.2 Initialization model

3.1.2 井周網格篩選

tNavigator 無法做到壓力梯度的表征，因此根據泄油半徑距離，通過篩選網格的方式進行壓力梯度的模擬，通過Voronoi 圖指定半徑的方式進行網格篩選。可將井周范圍內的網格進行篩選，便于之后對壓力梯度進行表征。以網格長度為最小單位步長進行篩選（圖3）。

圖3 井周網格篩選Fig.3 Grid screening around the well

將井附近網格篩選過后，對于井周圍不同距離的網格可近似表示壓力梯度。

3.1.3 滲透率時變

tNavigator 中的啟動壓力梯度限制為定值，無法體現(xiàn)非線性滲流規(guī)律，采用對滲透率進行時變的方法，即通過對滲透率進行修正來反映啟動壓力梯度的變化，從而反映非線性滲流。對于井周不同距離的網格進行非線性滲流模型的導入，改變滲透率，在井底流壓大于0 MPa 的情況下，將推導的非線性滲流模型導入，結合實際生產情況，對井周不同范圍內的網格進行滲透率時變。

3.2 結果分析

歷史擬合情況及對比，對導入考慮非線性滲流影響的模型作智能歷史擬合，采用拉丁超立方實驗算法對日產油量、累計產油量進行擬合（圖4、圖5），選取擬合程度較高且參數(shù)均科學合理的擬合模型，并與現(xiàn)有的擬啟動壓力梯度模型進行比較。

圖4 不同模型日產油量對比Fig.4 Comparison of daily oil production of different models

圖5 不同模型累計產油量對比Fig.5 Comparison of cumulative oil production of different models

模擬結果表明，考慮非線性滲流的模型，模擬日產油量與實際產油量相關性、趨勢性較好。累計產油量由非線性模型擬合結果為12.085 9×104m3，與實際產油量11.424 8×104m3相比，誤差為5.47%，擬啟動模型產油量為9.826 8×104m3，誤差為13.99%。由此不難發(fā)現(xiàn)，非線性模型擬合日產油量和累計產油量與實際生產情況更接近，符合率更高，由此表明該方法適用于特低滲、致密砂巖油藏數(shù)值模擬。

4 結 論

（1）基于非線性滲流理論及啟動壓力實驗、邊界層實驗，推導了不同類型儲層的非線性滲流新模型，綜合考慮了啟動壓力梯度、邊界層厚度對非線性滲流的影響，新模型考慮了流體黏度的影響。

（2）基于非線性滲流模型，開展商業(yè)軟件二次開發(fā)，通過篩選井周網格及對滲透率場與壓力場初始化，結合實驗數(shù)據，對滲透率進行時變反映非線性滲流對儲層的影響，形成了非線性滲流數(shù)值模擬方法。

（3）利用本文方法對典型井組歷史擬合，通過對比來看擬啟動壓力梯度模擬方法夸大了地層阻力，數(shù)值模擬結果低于實際產量，非線性滲流模型更符合實際生產情況，誤差更小，由此證明該方法適用于特低滲致密砂巖油藏數(shù)值模擬。