杜 川, 喻思羽, 李少華, 方 紅
(1.長江大學 地球科學學院,武漢 430100;2.中國石油 天然氣遼河油田分公司 勘探開發研究院,盤錦 124010)
隨著多點地質統計學(multi-point statistics, MPS)受到廣泛的關注,涌現出許多相對成熟的MPS建模算法。相對于傳統的兩點地質統計學,多點地質統計學利用訓練圖像代替變差函數表達地質變量的空間結構性[1],從中獲取更高階的樣本統計量[2],MPS建模優點在于既能夠模擬地質體復雜的幾何結構,也能忠實于條件數據。Guardiano等[2]提出的首個MPS算法每模擬一個網格節點都需要遍歷整個訓練圖像,計算成本非常高,使得該算法難以實際應用。之后,更好的算法被開發出來改善這些缺點。多點統計技術包括基于象元的建模算法SNESIM[3]、基于樣式的建模算法FILTERSIM[4]、SIMPAT[5]、DISPAT[6]、基于直接采樣的算法DS[7]、IMPACA[8]、HOSIM[9]、CCSIM[10]、基于圖像縫合技術的Image Quilting[11]。這些MPS算法使用不同的方法來讀取訓練圖像、存儲信息和復制模式。Mariethoz[12]對這些多點地質建模方法進行了定性比較。多點地質統計學建模已用于許多應用之中,例如常見的在儲層地質建模研究領域,利用Simpat模擬河流相儲層分布[13];多點地質統計學在河流相儲層建模中的應用[14];在三角洲相儲層建模領域的應用,如多點地質統計學建模方法在復雜疊置樣式砂體表征中的應用[15]。
在兩點地質統計學建模中,通過變差函數考慮統計數據之間的空間相關性[16]。類似的在多點地質統計學建模過程,需要通過給定的數據事件對訓練圖像進行掃描來獲取一個統計先驗模型,相關的統計數據取決于給定數據事件的參數,通過人為的參數調整,找到合適的參數設置是獲取良好建模效果的前提,不同的參數不僅影響建模的質量,也涉及CPU及內存開銷問題[17]。為此前人也進行了大量的探究,常見的參數敏感性分析,如Liu[18]提供了關于SNESIM參數以及每個參數對模型再現質量和CPU成本的完整指南;Meerschman[19]也為DS直接采樣算法做了類似的工作;Honarkhah[6]提出了DISPAT,它對用戶提供的參數不太敏感。在該方法中,通過計算不同樣板尺寸的圖案的平均熵,然后使用所謂的肘形圖的拐點,自動選擇最佳樣板尺寸。Strebelle[20]提出了一種通過優化數據樣板的大小來最小化SNESIM的計算成本,同時保持模式再現質量的方法。在該方法中,最小可接受的數據樣板大小是通過閾值來選擇的,超過該閾值,額外的條件數據不會提高估計的條件概率。Kolbj?rnsen[21]開發了一種自動方法來定量確定馬爾可夫網格建模中使用的多個網格的最佳數量,該方法使用從訓練圖像獲得的方向相關函數;Bai[22]還提出了一種量化方法來估計SNESIM算法的多個網格的最小數量。它使用由改進的聯合計數統計量化的空間關聯度來尋找最粗略的尺度,并因此基于模擬尺度應不小于目標尺度的標準來估計多個網格的最小數量。
傳統建模參數的評價主要依靠人工識別,即首先給定一組(有序)參數集,使用該參數集里的每個參數模擬一組隨機模型,然后通過人工判別優選與訓練圖像相似度較為合適的模型,以該模型對應的參數作為最優參數。人工識別的精確性取決于建模工作者的經驗,具有較強主觀性,同時人工識別效率低,不適應于現代自動化生產的需要。為了解決這一問題,我們引入了一種基于連通性函數[23-24]的多點地質統計學建模參數敏感性分析方法,通過計算條件數據下生成的模擬圖像的連通性概率曲線,并把訓練圖像與模擬圖像連通性概率曲線之間的相似性用余弦相似度來量化,得到基于連通性函數的空間相關性評價指標與建模參數的關系曲線,通過對量化數據進行分析,找到相似度曲線趨于平緩的拐點值,以此拐點值作為最優參數,大于最優值參數將變得不敏感,模型效果不會隨著參數值增加而得到明顯改善。
連通性函數被定義為網格上的兩個點屬于同一個連接區域的概率[24],利用連通性函數計算區域化變量指示圖的連通性,可以反映區域化變量在某一方向某一距離上連通性概率的變化程度。理論上隨機模型與訓練圖像的連通性差異越小,則表明二者沉積相的空間相關性及位置分布特征越相近,也說明了用于生成該隨機模型的建模參數越可靠。Eulogio[25]詳細地描述了連通性概率計算原理,其計算原理如下:
圖1定義工區網格為矩形G,G的網格節點(或象元)屬于子集S(砂巖相)或補集Sc(泥巖相)
G=SUSc
(1)
G的象元取值函數定義為指示函數I(u)
(2)
如果u和u′=u+h屬于相同的S,用u?u′表示,連通性函數τ(h)定義為屬于S的兩點連通的概率:
τ(h)=P(u?u+h│u,u+h∈S)
(3)
其近似計算公式為式(4)。
(4)

以砂巖相連通性為例,如圖1所示u和u′有三種情況,①u和u+h有一個或者兩個都屬于泥巖相,例如u3與u3+h或者u4與u4+h,則不參與砂巖相連通計算;②u和u+h都是砂巖,但不屬于相同砂巖體,例如u2與u2+h,則增加N(u,u+h∈S)數值;③u和u+h都是砂巖,屬于相同砂巖體,例如u1與u1+h則增加#N(u?u+h│u,u+h∈S)和#N(u,u+h∈S)值。

圖1 連通性示意圖Fig.1 Schematic diagram of connectivity
基于連通性函數的多點地質統計學建模參數優選方法是建立在一個前提條件下:以經典多點算法Snesim的一個重要建模參數--“搜索節點數”為例,基于連通性函數優選搜索節點數;隨著搜索節點的增加,隨機模型與訓練圖像的空間相關性、幾何形態越來越相似。
在此前提條件下采用有序排列的建模參數集建立一組隨機模型集,然后計算基于連通性函數和余弦相似度CosSim函數,評價建模參數集的多點地質統計隨機模型與訓練圖像的空間相關性及結構特征相似性,進而建立基于連通性函數的空間相關性評價指標與建模參數的關系曲線,選取評價指標開始平穩、進入平臺區域時的拐點所對應的參數值作為最優參數。
以建模參數“搜索節點數”為例,基于連通性函數的建模參數優選方法具體流程及流程圖(圖2)如下:

圖2 優選方法流程圖Fig.2 Flow chart of optimization method
1)輸入訓練圖像TI,給定MPS的關鍵建模參數——搜索節點的K個連續單調遞增數值(例如:1,2,…,K)組成的集合S={1,2,…,K}。
2)采用MPS,使用TI和S生成隨機模型集M,其中第k個參數Sk相應的隨機模型集定義為Mk,每個Mk包含n個隨機模型,則M中共計有N=K×n個模型。


(5)

(6)
計算Mk所有模型與訓練圖像TI關于連通性函數參數組的余弦相似度,然后取其數學期望作為建模參數Sk對應的空間相關性評價指標ξk
ξk=E[cosθk(Ai,Bi)]
(7)
5)k增加1,如果k小于K-1,則進入步驟4),否則進入步驟6)。
6)建立參數k與相似度的關系曲線,分析曲線變化趨勢,把曲線變化趨于平緩時的拐點建模參數k(i)作為最優參數。隨著建模參數值的遞增,基于建模參數值的隨機模型與訓練圖像之間相似度逐步增加,并最終達到難以區分的程度,即距離也來越小,相似度越來越高,同時對應的參數越不敏感。
訓練圖像是地質先驗模型,其樣式、大小等都會對模擬結果產生重大影響,以一個網狀河訓練圖像channel為例,圖3訓練圖像大小為101×101的網格,黑色部分代表河道砂,白色部分代表背景泥,訓練圖像數字化的結果是一個101×101的矩陣,“1”代表河道,“0”代表背景泥。其中背景泥占比0.7,河道占比0.3。

圖3 訓練圖像channelFig.3 Training image channel
2.2.1 搜索節點數
利用Snesim基于像元的多點地質建模對訓練圖像進行非條件模擬,為了找出最優建模參數值,這個過程中其他參數保持不變,只對Snesim重要建模參數“搜索節點”數進行更改,設置搜索節點數5到75,步長10。其他參數分別設置為:分類變量比例0.7、0.3,伺服器系統參數0.5;搜索半徑hx=50,hy=50,搜索角度angel_1=0, angel_2=0,多重網格層數為3。對應每個搜素節點數下進行模擬生成300個隨機模型。圖4是不同搜索節點數所模擬出來的一個模型(300個實現模型中的一個),隨著搜索節點數量的增加,所模擬出來的河道形態越來越好,與訓練圖像越來越相似,當增加到25和35以后,設置不同搜索節點數所模擬出來的模型的空間特征相似用肉眼就無法直接區分開來,同時隨著搜索節點增加,該參數敏感性降低,所生成的隨機模型效果改善不明顯甚至幾乎難以區分,往往在考慮計算效率以及實現效果上,通過人為肉眼觀察給定一個最優參數值往往具有誤差,會給建模中帶來不確定性。

圖4 不同搜索節點數對應的實現模型Fig.4 Implementation models corresponding to different node numbers(a)Node=5;(b)Node=15;(c)Node=25;(d)Node=35;(e)Node=45;(f)Node=55;(g)Node=65;(h)Node=75
根據圖2提出的優選方法流程,計算隨機模型與訓練圖像砂巖泥巖相在x、y方向上不同距離下的連通概率,基于連通性概率函數計算公式(4)計算訓練圖像channel砂巖相泥巖相在x、y方向上的連通概率,圖5是訓練圖像在x和y方向上砂巖泥巖在不同滯后距下的連通概率,隨著滯后距的增加砂巖泥巖相在x、y方向上的連通性總體上不斷趨于減小,其中砂巖的變化有周期化的上下波動,并且其變化規律能大致反映砂巖的空間分布特征。

圖5 訓練圖像channel在x、y方向上的連通性概率曲線Fig.5 Connectivity probability curve of training image channel in x and y directions
不同搜索節點所實現的隨機模型與訓練圖像在空間型態結構特征上存在差異,隨著搜索節點數量的增加,隨機模型的形態結構等特征與訓練圖像越來越相似,差異會越來越小,這種差異也同時體現在了連通性概率曲線上如圖6所示,其中mod1代表搜索節點數等于5的隨機模型(綠色),Mod2是搜索節點數等于75的隨機模型(藍色),實線與虛線分別代表泥巖相與砂巖相在x或y方向上的連通概率,由圖6可知,模型與訓練圖像的形態特征和空間結構越相似,二者之間的連通性概率曲線就越相近相似。(圖6曲線紅色和藍色相近,紅色與綠色差異大)

圖6 訓練圖像及其隨機模型的連通性概率曲線對比圖Fig.6 Comparison of the connectivity probability curves of the training image and its random model(a)x方向;(a)y方向
利用余弦相似度cosSimilarity函數計算訓練圖像與隨機模型連通性概率曲線的相似度,得到隨機模型(搜索節點數從5~75)與訓練圖像關于連通性概率的余弦相似度,然后取其數學期望作為建模參數“搜索節點”對應的空間相關性評價指標ξk建立相似度與搜索節點數的關系曲線(圖7),隨著搜索節點數從5開始增加,300個隨機模型與訓練圖像的相似度期望不斷增加,表明模型與訓練圖像越來越相似(根據余弦相似度定義,余弦值越接近1,表示越相似),當搜索節點數達到35時,空間相關性評價指標進入一個平臺區域,參數變得不敏感,模型質量將不再有明顯提升,是本例中的優選參數。

圖7 建模參數(搜索節點數)基于連通性函數的余弦相似度變化曲線Fig.7 Cosine similarity change curve based on connectivity function for modeling parameters (node numbers)
2.2.2 多重網格層數
多重網格層數對模擬的效果影響很大,但是隨之而來的是消耗更大的計算性能,利用上述分析方法,可以找到一個最優值,既能達到良好建模效果,同時也節約計算性能。對“多重網格重數”這一參數利用相同方法測試,設置多重網格重數從1到8。
保持其他參數不變設置搜索節點數35、分類變量比例0.7、0.3、伺服系統參數0.5,搜索半徑hx=50,hy=50、搜索角度angel_1=0, angel_2=0。利用snesim多點地質建模對訓練圖像進行非條件模擬,在每個多重網格重數下進行模擬生成300個隨機模型(圖8顯示對應網格層數下所模擬的300個實現模型中的一個)。

圖8 不同網格重數對應的實現模型Fig.8 Implementation models corresponding to different multigrids(a)Multigrids=1;(b)Multigrids=2;(c)Multigrids=3;(d)Multigrids=4;(e)Multigrids=5;(f)Multigrids=6;(g)Multigrids=7;(h)Multigrids=8
利用余弦相似度cosSimilarity函數計算訓練圖像與隨機模型連通性概率曲線的相似度,得到隨機模型(多重網格重數1-8)與訓練圖像關于連通性概率的余弦相似度,然后取其相似度均值建立其對應多重網格重數變化的曲線(圖 9)。

圖9 建模參數(網格重數)基于連通性函數的余弦相似度變化曲線Fig.9 Cosine similarity change curve based on connectivity function for modeling parameters (multigrids)
隨著多重網格數從1開始增加,300個隨機模型與訓練圖像的相似度期望不斷增加,表明模型與訓練圖像越來越相似(根據余弦相似度定義,余弦值越接近1,表示越相似),當網格重數到3時,空間相關性評價指標進入一個平臺區域,模型質量將不再有明顯提升,參數敏感性降低。這是因為當網格很粗時,相應的用來掃描訓練圖像的樣板規模會大于訓練圖像的規模,導致大的數據事件很難被找到。此時數據事件就會相應縮小,而這就等同于使用較小的數據樣板進行掃描,一般來說,網格重數設置多少取決于最大尺度結構信息的大小,在本例中101×101的網格訓練圖像中,利用連通性方法定量化地給出最佳網格重數為3,超過3以后,模擬效果不會有明顯改善,該實驗結果也符合前人所提出的認識。
采用基于連通性函數的多點地質統計學建模參數敏感性分析方法,分別對兩個關鍵建模參數(搜索節點數、多重網格層數)進行敏感性測試,得出以下認識:
1)隨著搜索節點數和多重網格層數增加,模型與訓練圖像的形態視覺特征越來越相似,利用連通性函數來量化這一空間相關性特征,結合余弦相似度度量方法對其進行定量的分析,找出價指標開始平穩、進入平臺區域時的拐點所對應的參數值作為最佳參數。
2)針對本例101×101的河道訓練圖像,在snesim模擬過程中搜索節點數不能過少,過少捕捉不到河道的結構特征,過多并不能明顯改善建模效果,參數在到達一定值后變得不敏感,反而會增加計算量,降低效率,多重網格層數也是如此。利用連通性的方法定量化的分析參數敏感性,找到搜索節點數設置在35左右,多重網格層數設置為3層,能獲取良好模擬效果,也能減少不必要的運算,提高效率。