盧異，胡浩 ，成亞斌，夏國朝，任光文

1.中國石油大港油田勘探開發(fā)研究院，天津濱海新區(qū)300280

2.中海福陸重工有限公司，廣東珠海519090

3.中國石油大港油田公司資源評價處，天津濱海新區(qū)300280

引言

隨著數(shù)值模擬理論的不斷完善以及計算機技術的快速發(fā)展，油氣藏數(shù)值模擬逐漸由人工歷史擬合向自動歷史擬合發(fā)展[1-7]，由傳統(tǒng)的單模型預測向多模型預測改進[8-11]。多模型自動歷史擬合主要分為兩個部分：多初始模型的建立與挑選以及多模型的自動歷史擬合。由于所建立的地質模型即使經過歷史擬合，也僅能說明模型計算結果和已知的歷史動態(tài)基本一致，無法保證能反映地下儲層的真實分布，同時也無法保證對未來的預測準確。而采用多個特征各異卻又符合歷史動態(tài)的模型同時進行預測，提供一個未來生產情況的參考區(qū)間而非參考值，將更能夠降低開發(fā)風險。

1 方法原理

1.1 主成分分析法降維

實際油藏建模中，即使經過粗化，大型油藏模型的網格數(shù)量仍可達數(shù)十萬以上，數(shù)據(jù)量龐大，導致模型后期處理以及自動歷史擬合效率低下，因此，需要對模型數(shù)據(jù)進行降維[12-14]。主成分分析（Principal Component Analysis，PCA）[15]是由Hotelling[16]提出的數(shù)據(jù)降維方法，通過對原始變量的相關矩陣或協(xié)方差矩陣結構的研究，將多個隨機變量轉換為少數(shù)幾個新的隨機變量（同時保留原始變量絕大部分特征信息），從而達到降維目的。

設隨機實現(xiàn)樣本為M={m1,m2,···,mNr}，其中，mi表示第i個實現(xiàn)的參數(shù)向量，Nr則是實現(xiàn)個數(shù)。

根據(jù)最優(yōu)重建準則，PCA 目標函數(shù)為

式中：

降維矩陣和原始參數(shù)向量相乘可將模型的參數(shù)向量由Nm維轉化為n維（通常Nm?n），樣本矩陣M和單個實現(xiàn)數(shù)據(jù)向量m的降維可由式（4）計算

式中：

1.2 改進的K 中心點聚類

多個地質模型經過聚類，將屬性相近的地質模型歸并成一類，可以減少用于歷史擬合的模型數(shù)量，本次采用K 中心點聚類方法進行聚類[17]，主要原理如下：

設數(shù)據(jù)集合X中有f條a維數(shù)據(jù)，每條數(shù)據(jù)xq為一個對象。當e個初始中心點對象隨機或按照一定依據(jù)選定后，將剩余的f?e個非中心點劃分為e個組，分組的規(guī)則是將非中心點對象化歸到離它最近的中心點。準確地說：如果Oq是一個非代表對象，Ow是一個代表對象（中心點），并且Oq與Ow的距離d(Oq,Ow)是Oq與所有中心點距離最近的，則稱Oq屬于Ow所代表的簇。兩個對象的距離為

式（5）表示的是曼哈頓距離，根據(jù)不同情況可以選擇歐氏距離等其他形式的距離。

初始中心點確定并劃分好簇后，隨后每一步隨機使用一個非中心點對象Or替換現(xiàn)有中心點對象Ow，提高聚類質量。為了確保中心點Ow能夠被Or很好地替代，在每一次替換時需要考慮4種情況。

對于所有f?e個非代表對象Oq，一次替代所造成的總代價為

式中：

若Ct小于0，則認為本次替換可以提高聚類質量，認可本次替換，反之則拒絕替換。

常規(guī)K 中心點聚類時間復雜度很高，最多需要e(f?e)2次計算判斷，尤其對初始中心點的選擇異常敏感?；陬I域的K 中心點算法[18]通過對中心點設定鄰域半徑，選取相互距離較遠的e個處于樣本分布密集區(qū)域的數(shù)據(jù)作為K 中心點算法的初始聚類中心，即可以在初期找到佳中心點或到達最佳中心點附近，減少中心點替換次數(shù)，又能避免初始中心點處于同一簇造成重復計算。

數(shù)據(jù)對象的鄰域半徑R定義為

對于任意數(shù)據(jù)對象xw，以xw為中心，半徑為R的圓形區(qū)域內的數(shù)據(jù)對象稱為數(shù)據(jù)對象xw的鄰域，用δw表示。

在鄰域的限制下，K 中心點聚類效率得到顯著提高。

1.3 自動歷史擬合模型參數(shù)

在自動歷史擬合研究中，需要將參數(shù)表示為向量的形式。如同時對孔隙度、滲透率、飽和度、相對滲透率和水體參數(shù)進行調整，則油藏參數(shù)組成的控制變量h可表示為

為了得到統(tǒng)一的歷史數(shù)據(jù)格式，即每項數(shù)據(jù)在每個時間點都有數(shù)值，對應時間步測量則為測量值，沒有測量也必須填寫任意值，通常為0。

1.4 目標函數(shù)

近年來，基于貝葉斯理論目標函數(shù)[11]被廣泛用于自動歷史擬合領域，不僅可以考慮動態(tài)歷史和模型響應之間的偏差，同時可以充分利用隨機實現(xiàn)中的先驗信息，使得擬合后的模型更為符合實際地質統(tǒng)計規(guī)律。

目標函數(shù)的表達式為

將式（5）代入式（10），通過PCA 進行變換可得模型降維后目標函數(shù)式（11），可有效減少求解的計算量。

1.5 SPSA 優(yōu)化算法及改進

SPSA 算法[19]可以一次同時擾動所有控制變量，每次迭代僅需擾動兩次即可求取隨機梯度，隨機梯度的期望為真實梯度，且恒為上山方向，在解決多變量的歷史擬合問題時，效率高且收斂速度快。

本次在SPSA 算法的基礎上，為了使隨機梯度更接近真實梯度，提高算法穩(wěn)定性，采用在每個迭代步中求取多個隨機梯度，以其均值作為搜索方式

目前普遍認可實際油藏相鄰網格的參數(shù)存在一定的相關性，參數(shù)場中相鄰網格之間為漸變關系，引入控制變量協(xié)方差矩陣來指導生成擾動向量[20-21]

2 多模型自動歷史擬合流程

多模型歷史擬合完整流程如圖1所示，主要包括5個步驟。

圖1 多模型歷史擬合完整流程Fig.1 Process of multi-model history matching

（1）整理靜態(tài)地質資料和動態(tài)觀測數(shù)據(jù)，依據(jù)靜態(tài)資料使用隨機模擬建立多個實現(xiàn)；

（2）使用PCA 算法對模型數(shù)據(jù)進行降維；

（3）使用改進的K 中心點聚類算法挑選指定數(shù)量的隨機實現(xiàn)作為初始模型；

（4）設定SPSA 算法參數(shù)，對所有初始模型進行自動歷史擬合；

（5）使用擬合后的模型進行動態(tài)預測，評價開發(fā)不確定性。

3 實例驗證

3.1 數(shù)據(jù)準備及隨機實現(xiàn)的生成

以反五點井網模型為例，工區(qū)面積為1000m×1000m，網格劃分為25×25×1=625個，網格尺寸為40 m×40m×20m，包含1口注水井（I）和4口采油井（A、B、C 和D 井）。研究區(qū)儲層參數(shù)如表1 所示，5口井的物性參數(shù)見表2。

表1 儲層參數(shù)Tab.1Properties of reservoir

表2 井的物性參數(shù)Tab.2 Physical properties of five wells

依據(jù)5口井的地質資料，使用序貫高斯模擬隨機生成701組孔隙度、滲透率和凈毛比的模型。選擇其中一組作為參考模型，將其進行數(shù)值模擬運算后的生產數(shù)據(jù)作為油藏歷史數(shù)據(jù)。

圖2為參考模型的孔隙度、滲透率及凈毛比的空間分布。該模型共模擬7200d，前6000d 作為歷史數(shù)據(jù)，用于歷史擬合，將后1200d 假設為未來真實動態(tài)，用于和經歷史擬合確定后的模型參數(shù)預測結果做對比，分析預測的不確定性。

模型中的4口采油井采用定液量的方式生產，產油量變化如圖3所示。從圖3中可以看出，4口井的穩(wěn)產時間不同，D 井的穩(wěn)產時間最長、B 井的穩(wěn)產時間最短。

圖3 參考模型生產井產油量曲線Fig.3 Oil productionrate of four production wellssimulated by reference geologic model

3.2 模型降維及初始模型的聚類挑選

使用PCA 算法對700個隨機實現(xiàn)進行降維，每個實現(xiàn)包含孔隙度、滲透率和凈毛比3類共625×3=1875個數(shù)據(jù)，經降維計算，前221個主成分即可表達原始數(shù)據(jù)90%以上的特征信息。由此每個實現(xiàn)的參數(shù)數(shù)據(jù)可由1875維降至221維。主成分數(shù)據(jù)對總樣本數(shù)據(jù)的特征貢獻率如圖4所示。降維后的模型數(shù)據(jù)量大大減少，可以明顯減少聚類計算量以及擬合階段需要調整的參數(shù)數(shù)量。

圖4 主成分特征貢獻率Fig.4 Contribution rateof principal component

使用改進的K 中心點聚類算法，在降維后的數(shù)據(jù)中挑選出5個實現(xiàn)作為初始模型，初始模型的滲透率分布如圖5所示，觀察可知，5個初始模型的滲透率分布具有明顯差別，各初始模型能很好地代表一類實現(xiàn)，能夠更廣泛地包含儲層物性的分布。

圖55 個初始模型的滲透率參數(shù)場Fig.5 Permeability distributionof five initial models

3.3 自動歷史擬合及預測

本例中僅考慮孔隙度、滲透率和凈毛比3類變量對歷史擬合的影響；經敏感性分析，發(fā)現(xiàn)滲透率對擬合效果產生較大影響，孔隙度和凈毛比影響可以忽略，且自動歷史擬合迭代超過30次后，目標函數(shù)已趨于收斂。

將滲透率作為擬合參數(shù)，設定最大擾動次數(shù)為30次。其余參數(shù)設置為：隨機梯度個數(shù)為3，增益系數(shù)為10（增益系數(shù)是期望迭代次數(shù)的1/10或者更少），初始擾動幅度為0.6，初始搜索步長為0.3。歷史擬合后的模型滲透率分布如圖6所示，生產井A歷史擬合前后的產油量對比見圖7。

圖6 擬合后的滲透率模型Fig.6 Permeability distribution after history matching

圖7 不同初始模型下A 井產油量對比Fig.7 Contrast of oil productionrate of Well A simulated by five initial geologic models

由圖7可以看出，A 井經過生產歷史擬合后，5個模型的產油量與歷史數(shù)據(jù)的符合率顯著提高，表明擬合效果較好。同時，從圖6可以看出，盡管各模型的滲透率分布存在差異，但都能獲得相似的生產歷史擬合效果，這也反映出生產歷史擬合存在多解性，即多個地質模型的分布可以達到較為相似的生產效果。

另外，從A 井6 000～7 200 d 的生產預測結果可以看出，使用任何一個經歷史擬合后的模型預測A井的未來產量，都與A 井真實產量存在一定差異，但A 井的真實產量包含在5個模型預測的產量范圍內（圖7）。這表明，使用單一模型預測的生產效果不佳，而采用多個模型預測的生產效果接近真實情況，因此，在今后的油藏數(shù)值模擬中，有必要采用多個模型進行生產效果預測。

4 結論

（1）地質模型經過PCA 降維，可以保留地質模型中大部分特征的信息，并大大降低參數(shù)維度，減少聚類計算量，提高擬合效率。

（2）提出了利用改進K 中心點聚類方法為基礎的初始模型挑選方法，使挑選出的模型具有很好的代表性，更廣泛地反映儲層的可能分布狀況。

（3）在自動歷史擬合中將基于貝葉斯理論的目標函數(shù)和PCA 算法相結合，提高了目標函數(shù)的計算效率。

（4）以理想反五點井網模型為研究對象，進行模型降維、聚類挑選及自動歷史擬合的完整流程，驗證了多模型自動歷史擬合的可行性，實例中的整體動態(tài)數(shù)據(jù)的匹配度有明顯提高，效果顯著。