基于多次數據吸收集合平滑算法的自動油藏歷史擬合研究

王澤龍，劉先貴，唐海發，呂志凱，劉群明

(1.中國科學院大學，北京 100049；2.中國科學院大學滲流流體力學研究所, 河北 廊坊 065007； 3.中國石油國際勘探開發有限公司，北京 100034；4.中國石油勘探開發研究院,北京 100083)

0 引 言

油藏數值模擬技術是室內油藏工程研究的重要手段。傳統油藏數值模擬技術是將油藏地質模型與流體模擬程序相結合，再輸入流體性質參數和生產制度，經過大量計算從而得到預測數據。地質模型靜態參數主要包括滲透率、孔隙度、初始溫度、壓力等，動態數據來自現場油井的實測生產數據。但研究過程中要求具備較高的地質學和油藏工程專業素養，才能較好地完成地質模型的調整工作。在實際工作中，如果油田處于開發初期，或者從業人員經驗尚淺，人為手動修正的地質模型往往不能夠真實反映出油藏地質特征[1]。此外，現代油藏管理需要進行不確定性分析，純粹的人工歷史擬合已難以適應發展要求。

自動油藏歷史擬合技術自20世紀90年代進入快速發展階段，大量學者開展了相關研究[2-5]。GeirEvensen通過改進卡爾曼濾波線性問題的局限，提出了擴展卡爾曼濾波，后為克服其敏感性矩陣的大量計算問題，提出集合卡爾曼濾波[6]，即EnKF，并開始在天氣預測領域應用。N?vda等[7]將EnKF引用到石油工程領域，用于更新油井附近的油藏模型，Lorentzen[8]、Oilver[9]、Chen[10]、Reynolds[11]、Jafarpour[12]、王玉斗[13]、薛亮[14]等學者均在EnKF引入到油藏模擬應用領域進行了大量研究，但均未徹底解決EnKF算法在油藏歷史擬合中計算效率低，且更新后的狀態參數與實際的狀態參數不一致的問題。Reynolds[15]發現對相同實測數據的多次數據吸收可以改善EnKF的效果，Emerick與Reynolds[16]結合集合平滑算法，提出了多次數據吸收集合平滑算法，即ES-MDA。在前人研究成果基礎上，進一步推導出ES-MDA算法的核心公式，開發了自動歷史擬合軟件，確定了多次數據吸收的選擇參數，并應用于海相砂巖油藏的歷史擬合，取得了較好的應用效果。

1 數據吸收算法的數學描述

EnKF和ES都是利用狀態參數(如孔隙度、滲透率、壓力、飽和度等)的集合，近似表達先驗概率的期望值和方差，通過求解后驗概率的期望值，實現后驗概率期望值與真值的誤差最小化、后驗概率方差最小化的計算過程。在此基礎上，通過求得的后驗概率，得到更新后的狀態參數，從而獲得更新后的油藏模型。2種算法均屬于蒙特卡洛模擬的范疇，基于貝葉斯概率論來計算后驗概率。

不同點在于：EnKF是在整個時間序列中一步一步地吸收實際測量數據，逐步更新后驗概率估計；ES則是卻在整個時間、空間域中，一次性吸收所有實測數據，直接更新整個油藏模型的后驗概率估計。對于線性問題，ES簡單高效。但對于非線性和混沌動態系統，EnKF相對于ES要更加優越，狀態集合的迭代更新過程符合實測數據所代表的真實情況。以下簡要對EnKF與ES算法進行數學描述。

1.1 EnKF與ES算法的數學描述

1.1.1 EnKF算法

EnKF是目前油藏歷史擬合中應用最廣泛的數據吸收算法。定義擴展狀態向量xk為時間序列k時刻估算的狀態參數，包括地質模型的靜態參數、油藏的主要動態參數以及實測數據。

(1)

式中：m為地質模型的靜態參數向量在k時刻的估值；p為動態參數向量(如油藏壓力、油氣水飽和度等)；d為實測數據向量(如油氣產量、產水量、井底壓力等)；下標k為時間步；x為擴展狀態向量。

EnKF的預測方程為：

(2)

(3)

EnKF的更新方程為：

(4)

H=[OI]

(5)

為方便計算，定義卡爾曼系數矩陣K為：

(6)

則，EnKF的更新方程可表達為：

(7)

同理，在k+1時刻及未來任意時刻，可通過連續循環以上過程，不斷更新地質模型集合。EnKF算法需要實時輸入狀態參數，同步進行正向的模擬計算，因此，EnKF更新方程中，模型參數和狀態參數同時更新，會造成更新后的狀態參數與實際的狀態參數不一致。如油藏的狀態壓力、溫度以及流體飽和度等參數的更新，與后驗地質模型的模擬結果相互矛盾。

1.1.2 ES算法的數學描述

ES可以全局性吸收所有時間序列的實測數據，并一步完成后驗概率的計算。對于一套靜態地質模型，只需運行一次完整時間序列的油藏流體模擬即可完成。因此，可直接由靜態地質參數矩陣mk代替擴展狀態參數矩陣xk，預測方程可表示為：

(8)

(9)

ES算法的更新方程為：

(10)

ES算法能夠一次性批量完成數據吸收并一次性完成所有模型集合的更新，避免了EnKF算法中需要頻繁調用正演模型文件獲得狀態參數的缺點，但由于油藏模擬是一個非常復雜的非線性計算過程，ES算法的更新不能給出準確預測結果，故不能應用于油藏歷史擬合。

1.2 ES-MDA算法的數學描述

借鑒EnKF算法多次循環迭代的思路，減小ES算法單次更新的步長，在小步長的情況下，對油藏模擬進行線性處理，通過多次小步長的更新完成計算，即為ES-MDA算法。

首先，人為選擇數據吸收的次數Na，并且針對每次數據的吸收過程定義系數αi，滿足如下關系：

(11)

在每次數據吸收過程中，進行以下步驟操作。

(1) 運行整套靜態地質模型集合的油藏模擬程序，獲得全時段的預測數據。

(2) 對于集合中任意一個模型子集，按照式(12)為實測數據添加誤差。

(12)

式中：dobs為實測數據。由于測量數據存在誤差，為修正誤差引入白色噪音向量zd，zd符合高斯分布zd～N(0,I)。

(3) 得到ES-MDA更新方程。

(13)

ES-MDA算法的流程框圖見圖1。

圖1 ES-MDA算法的流程框圖Fig.1 The flow chart of ES-MDA algorithm

相比EnKF，ES-MDA按照平滑算法的方式，一次性完成對所有歷史數據的吸收，進而減少了油藏狀態參數的代入需要，避免了油藏狀態參數更新異常問題。同時，ES-MDA是基于卡爾曼濾波在線性問題上構建貝葉斯概率理論似然函數的算法，其直接適用條件是線性問題或近似線性問題，并且參數的概率分布為高斯分布，初始模型的參數集合足夠大。大量文獻[10-13,17]證明EnKF作為集合算法的一種典型，當地質模型參數集合足夠大(如100組)時，即可取得較好的歷史擬合效果，故可推斷ES-MDA算法適應于油氣藏歷史擬合。按照上述方法，通過Matlab編程形成了一套自動歷史擬合程序，可用于油藏歷史擬合的實例研究。

2 實例應用

2.1 油藏地質及開發概況

北海布倫特油田總體為西東向帶狀半圓頂背斜構造，構造面積為31.5 km2，內部發育一條近東西走向斷層，北部邊界為閉合斷層。主力油藏埋深為1 956～2 210 m，油水界面在背斜頂部以下90 m處，油藏厚度約為15 m，初始含油飽和度約為80%，油藏類型為無氣頂邊底水油藏。沉積物源方向為近西東向，儲層巖相為砂巖，伴隨少量泥巖夾層，具有較強的非均質性。儲層孔隙度為5.1%～34.2%，平均為18.5%；儲層水平方向平均滲透率為325 mD，垂直方向平均滲透率為26 mD，為常規高孔高滲儲層。地層壓力系數為1.11，屬正常壓力系統。原始地層壓力為18 MPa，地層溫度為76 ℃。地面原油密度為0.71 g/cm3，76 ℃原油黏度為0.31 mPa·s。

本研究中,對照組給予常規用藥治療,觀察組則給予常規用藥聯合無創呼吸機治療。結果顯示,觀察組COPD并呼吸衰竭治療效果高于對照組,P<0.05;觀察組呼吸衰竭糾正的時間、住院時間優于對照組,P<0.05;治療前兩組PCO2、PO2兩項血氣指標和心率相近,P>0.05;治療后觀察組PCO2、PO2兩項血氣指標和心率優于對照組,P<0.05。

該油藏于2008年開始投入開發，以四點法部署開發井網，平均井距為1 000 m，部署油井20口、注水井10口，均為直井，所有井均配備數據傳感器，用于實時收集井口動態生產數據。截至目前，區塊累計產油1 165×104m3，采出程度約為32%。

2.2 ES-MDA自動油藏歷史擬合

2.2.1 模型建立

油藏構造模型由74×24×7個正交立方網格組成，網格大小為150.0 m×150.0 m×4.5 m。鑒于儲層巖相主要為海相砂巖，且孔隙度、滲透率呈高斯分布，采用序貫高斯模擬(簡稱SGSim)進行隨機儲層建模。以井口測井解釋孔隙度及鉆井取心測量滲透率作為隨機儲層建模硬數據約束，結合儲層沉積成藏認識，選取以下參數進行SGSim建模：水平方向最長相關長度為5 000 m，方位角為80°W，其正交方向相關長度為500 m。為定量描述儲層參數的置信度，通過SGSim生成100個初始地質模型集合。以滲透率為例表述靜態地質參數的更新調整，圖2a為100個初始地質模型的滲透率平均值分布，圖2b為其中4個初始地質模型子集的滲透率分布。由圖2可知，100個初始地質模型展現出100種油藏的非均質情況，但在滲透率平均值模型中，每個網格的數值都接近地質模型中滲透率的數學期望值，且變化趨勢光滑，因此，僅以井口采樣點數據作為約束條件的建模結果置信度較低。在初始的油藏認識情況下，一般用100個初始地質模型集合的平均值作為制訂開發方案的依據。按照閉合油藏管理理念，在獲得足夠油田開發動態數據后，進行油藏歷史擬合，對初始地質模型進行修正，以獲得更準確的更新地質模型。

圖2 初始地質模型集合的滲透率模型Fig.2 The permeability model of the initial geological model ensemble

2.2.2 歷史擬合

利用該地質模型及自主編寫的ES-MDA油藏自動歷史擬合程序對油田內20口油井、10口注水井進行歷史擬合，程序輸入文件為井口生產動態數據以及100個等概率油藏初始模型構成的集合。設置ES-MDA數據吸收次數為4，每次數據吸收的系數分別為α1=α2=α3=α4=4。每次數據吸收對單個油藏模型均調用一次油藏模擬器，并按照油田生產制度進行油藏全時段數值模擬，共進行400次油藏數值模擬運算，耗時2 d。自動歷史擬合程序輸出結果分別為更新后的地質模型集合以及動態數據擬合圖表，擬合結果(圖3～7)較好，各項預測指標與完測數據的誤差均控制在10%以內，較好地反映了油藏實際。

圖3 更新地質模型集合的滲透率模型Fig.3 The permeability model of the updated geological model ensemble

圖4～6為實測與模型預測的生產參數對比，由圖4～6可知：原始模型的預測數據具有很大的分散性，主要原因是地質參數的不確定性；而經過歷史擬合后的更新模型預測數據集中程度很高，且基本與實測數據吻合，再次證明更新地質模型的置信度很高。

圖4 生產井原油產量擬合結果

圖5 生產井產水量擬合結果Fig.5 The matching results of water yield of production well

圖6 注水井井底壓力擬合結果Fig.6 The matching results of downhole pressure of water injection well

為定量分析歷史擬合的置信度，定義油藏流體數值模擬的預測相對誤差為e，用于評判擬合后的模擬數據與實際測量數據的偏離程度。計算方法為：

(14)

式中：N為模擬數據的個數。

通過ES-MDA算法更新地質模型，地質模型集合中各項地質參數接近，模擬結果具有很高的集中度(表1～3)，油井的產油量以及注水井的井底流壓擬合結果較好，誤差在1%左右；油井的產水量擬合結果稍差，但誤差也低于10%。

表1 北海布倫特油田某區塊產油量歷史擬合統計Table 1 The history matching statistics of oil yield in a block of the Brent Oilfield in the North Sea

表2 北海布倫特油田某區塊產水量歷史擬合統計Table 2 The history matching statistics of water yield in a block of Brent Oilfield in the North Sea

表3 北海布倫特油田某區塊井底流壓歷史擬合統計Table 3 The history matching statistics of downhole flowing pressure of a block of the Brent Oilfield in the North Sea

在效率方面，自動歷史擬合優勢明顯，大尺度油藏模擬，若采用人工歷史擬合需要數月時間，而應用ES-MDA油藏自動歷史擬合，相同條件下效率可提高數百倍。與EnKF算法相比，ES-MDA算法效率也可提高5倍以上。

3 結 論

(1) 對于線性問題，ES算法能夠一次吸收所有時間序列數據，模型一次更新到位；但對于非線性問題，比如復雜油藏的模擬過程，由于ES算法計算步長過大，存在模型修正失真的問題。

(2) 通過有效結合ES算法和EnKF算法，推導出ES-MDA算法的核心公式，形成了自動油藏歷史擬合軟件。在實際應用中，地質模型集合中各個子項參數均需呈高斯分布，否則需要對模型參數進行高斯化處理。

(3) 應用ES-MDA算法對北海布倫特油田進行自動歷史擬合，結果表明更新后地質模型集合的子集模型之間滲透率分布具有更高的一致性，能夠有效表征真實油藏地質特征，且通過其開展的油藏數值模擬達到預期效果，預測數據與實測數據匹配程度達到90%以上。