基于細胞自動機的斷控巖溶過程數(shù)值模擬方法

劉彥鋒, 張文彪, 段太忠, 趙華偉

(中國石化石油勘探開發(fā)研究院， 北京 100083)

隨著塔河油田勘探開發(fā)的不斷深入，發(fā)現(xiàn)在塔河油田托甫臺以及順北地區(qū)縫洞儲集體的發(fā)育與大型走滑斷裂之間具有較好的匹配關(guān)系，且溶蝕發(fā)育具有向“縱深”發(fā)展的趨勢，呈現(xiàn)出“斷控巖溶”特征。據(jù)此，魯新便等[1-2]提出“斷溶體”圈閉的概念，并通過多年的礦場實踐總結(jié)出斷溶體油藏高產(chǎn)井主要沿大型斷裂帶分布的特征。近年來，中國一些專家學者針對斷溶體的描述、分類、地震預測及三維表征等做了較多的工作。周文等[3]根據(jù)斷裂與溶洞發(fā)育關(guān)系，提出“控洞斷裂”概念，認為塔河地區(qū)大-中型斷裂對巖溶洞穴的控制作用明顯，奠定了斷控巖溶作用的基礎(chǔ)。商曉飛等[4]根據(jù)斷控巖溶的發(fā)育特點，按照溶蝕強度的變化劃分為不同的溶蝕相帶，指出溶蝕過程可以根據(jù)主控因素進行追蹤和描述。魯新便等[5]、胡向陽等[6-7]提出分級相控巖溶建模的思路，采用多尺度分類建模的思路，將不同發(fā)育規(guī)模的縫洞儲集體進行融合。已有關(guān)于斷溶體的描述及表征主要集中在對斷溶體圈閉形態(tài)的刻畫、連通性分析、成因模式總結(jié)等，對于斷溶體的地質(zhì)成因過程方面主要停留在定性描述階段，缺乏一種定量化的模擬方法來推動斷控巖溶過程的定量分析。已有巖溶縫洞型儲層地質(zhì)建模技術(shù)，目前仍主要以傳統(tǒng)的地質(zhì)統(tǒng)計學和地震雕刻法為主，尚缺乏一種基于巖溶形成過程的模擬技術(shù)來刻畫儲層的分布形態(tài)。

關(guān)于巖溶過程模擬，已有研究多數(shù)集中在小尺度孔隙的溶蝕過程物理模擬和數(shù)值計算，缺乏對大尺度特別是油藏尺度溶蝕孔洞發(fā)育規(guī)律的模擬研究。何治亮等[8]、彭軍等[9]根據(jù)深層高溫高壓的特點，通過物理模擬實驗的方式探索了白云巖與灰?guī)r溶蝕速率的變化，指出隨著溫度壓力增加，溶蝕速率均明顯增強的特征。范維等[10]通過室內(nèi)試驗研究差異溶蝕作用，表明除了溫度、壓力和水溶液介質(zhì)等因素的，碳酸鹽巖因為巖石成分和結(jié)構(gòu)的差異也會被表現(xiàn)出不同的溶蝕速率。Kang等[11]提出基于格子玻爾茲曼數(shù)值模擬方法進行孔隙級別的溶蝕-沉淀過程模擬，在溶蝕過程數(shù)值模擬方面邁出了重要一步，但該方法還主要停留在微觀孔隙的溶蝕模擬階段，對于大型溶蝕孔洞級別的模擬未曾涉及。

細胞自動機，作為一種自然界中復雜系統(tǒng)的完全離散化數(shù)學模型，非常適合描述在離散時間維上演化的動力學過程。Burgess[12]以碳酸鹽巖為例，分別考慮其平面演化及縱向生長的總體特征，并加入沉積物搬運及海平面變化等因素，采用細胞自動機的基本規(guī)則，再現(xiàn)了碳酸鹽巖沉積演化過程；Liang等[13-14]提出了基于細胞自動機的地貌動力學模型模擬河流-三角洲沉積過程，演化規(guī)則考慮了水體、沉積物、搬運、剝蝕等因素。說明細胞自動機方法在沉積物演化方面具有較好的適用性，表明溶蝕改造過程也可借鑒該方法進行數(shù)值模擬研究。

基于以上，現(xiàn)將細胞自動機方法引入到碳酸鹽巖溶蝕過程模擬研究，在斷控巖溶過程主控因素分析基礎(chǔ)上，設計溶蝕過程基本規(guī)則，形成基于細胞自動機的斷控巖溶數(shù)值模擬方法。

1 斷控巖溶過程地質(zhì)特征

斷溶體的形成過程是大氣水沿近垂直的斷裂系統(tǒng)下滲至碳酸鹽巖地層，巖溶水優(yōu)先沿寬度大的主干斷裂下滲，同時以較低的速度向側(cè)向滲流，表現(xiàn)為以深大斷裂帶為核心發(fā)生溶蝕擴大作用[15-16]?；诖罅恳巴饴额^、鉆井、巖心、測錄井及地震資料，分析斷溶體發(fā)育位置、形態(tài)、規(guī)模，與大型斷裂、裂縫、地層巖性之間的控制作用關(guān)系，確定斷控巖溶過程劃分為4個階段(圖1)。

圖1 斷控巖溶演化過程示意圖Fig.1 Schematic diagram of fault controlled karst evolution process

(1)斷溶初期：受走滑構(gòu)造運動影響，形成從地層深部到淺部的立體斷裂結(jié)構(gòu)，走滑斷裂縱向延伸長度較大，在主干斷裂附近會伴隨一系列小規(guī)模的斷裂，整體呈“V”字形，是斷控巖溶過程最重要的格架基礎(chǔ)。

(2)斷溶早期：大氣水沿近似垂直的斷裂系統(tǒng)下滲，優(yōu)先沿著多組裂縫交匯處、裂縫寬度大、密度大的區(qū)域下滲，裂縫的擴溶作用是斷控巖溶的最重要表現(xiàn)形式，斷溶體發(fā)育形成雛形。

(3)斷溶中期：隨著巖溶作用增強，孔洞規(guī)模逐漸加大，裂縫不斷溶蝕而拓寬，在裂縫密集的地方，孔洞開始連片、合并，形成縱向連通的溶蝕孔洞和小規(guī)模溶洞；除了斷裂系統(tǒng)，背景巖性也對溶蝕過程有一定影響。

(4)斷溶晚期：溶洞繼續(xù)擴大，周圍溶蝕孔洞也繼續(xù)增多增大，最后演化至現(xiàn)今大型洞穴，是巖溶最發(fā)育的區(qū)域，周圍發(fā)育大規(guī)模溶蝕孔和溶蝕擴大縫，向外逐漸過渡到裂縫儲集區(qū)，最后伴隨上覆地層沉積覆蓋，地表水引起的溶蝕作用逐漸停止。

2 細胞自動機巖溶模擬方法

2.1 細胞自動機基本原理

細胞自動機是一種離散迭代演化模型，可用于模擬復雜動態(tài)系統(tǒng)隨時間的變化，具有運算速度快、程序編碼簡單的優(yōu)勢[17-19]。該演化模型由無限(或有限)個相同的網(wǎng)格組成，每個網(wǎng)格處于一種有限狀態(tài)，每個網(wǎng)格在t時刻的狀態(tài)由t-1時的有限集合的狀態(tài)確定，且每次演化時，每個網(wǎng)格遵從同一規(guī)則一齊演進，大量細胞通過領(lǐng)域細胞的簡單相互作用而構(gòu)成動態(tài)系統(tǒng)的演化，它可以通過簡單的演化規(guī)則模擬宏觀的復雜系統(tǒng)。目前細胞自動機已成為一個應用廣泛的研究領(lǐng)域，在圖像處理、加密、壓縮、計算機仿真模擬等方面都有廣泛的應用[20-22]。

針對特定的物理過程需要在細胞自動機的基本假設下設計相應的演化模型。細胞自動機通過刻畫最小單元的局部作用及其多個最小單元的相互作用模擬實際物理過程。一個完整的細胞自動機模型由細胞、細胞空間、鄰域和細胞規(guī)則4部分組成。細胞是離散的組成細胞自動機的最小單元；全部的細胞組成了細胞空間，代表所描述的整個物理系統(tǒng)；鄰域表示相互作用的局部細胞團體；細胞演化規(guī)則刻畫鄰域內(nèi)細胞的作用機理，是細胞自動機的靈魂。對特定物理過程的細胞自動機模擬的好壞關(guān)鍵在于演化規(guī)則是否反映了客觀事物內(nèi)在的本質(zhì)特征。

2.2 巖溶過程細胞自動機建模

細胞自動機模擬斷溶體儲層的形成以走滑斷裂形成的裂縫網(wǎng)絡為初始條件，模擬地層在有一定酸濃度的大氣淡水在地層巖性差異、裂縫網(wǎng)絡和重力等因素的共同作用下的演化過程。細胞自動機采用自下而上的方法建模，對溶蝕過程建模的過程主要包括確定溶蝕過程的性質(zhì)、對斷控巖溶系統(tǒng)進行網(wǎng)格狀分割、確定細胞的初始狀態(tài)、設計溶蝕過程的演化規(guī)則。

(1)確定溶蝕過程的性質(zhì)。細胞自動機最早提出時用于對自組織系統(tǒng)的模擬，不受外界因素影響，只考慮局部因素，目前細胞自動機已經(jīng)可以用來描述含有外部因素的系統(tǒng)[23]。斷溶體形成過程中不斷有巖溶水從地表向下供給，因此，溶蝕過程是具有外部因素影響的三維動態(tài)系統(tǒng)。

(2)斷控巖溶系統(tǒng)網(wǎng)格狀分割。采用地質(zhì)模型網(wǎng)格的方法對系統(tǒng)進行三維結(jié)構(gòu)化網(wǎng)格劃分，整個網(wǎng)格系統(tǒng)組成細胞空間，每個網(wǎng)格代表一個待模擬的細胞單元，后續(xù)描述中的網(wǎng)格即代表細胞，演化規(guī)則考慮三維鄰域狀態(tài)。

(3)確定細胞初始狀態(tài)。初始狀態(tài)是走滑作用下形成的斷裂系統(tǒng)，該狀態(tài)由裂縫網(wǎng)絡和基質(zhì)網(wǎng)格組成，裂縫發(fā)育的網(wǎng)格f=1表示已溶蝕，基質(zhì)處溶蝕相f為0，表示未溶蝕。溶蝕狀態(tài)的取值只能為0或1，中間狀態(tài)通過下面敘述的其他參數(shù)表示。

(4)設計溶蝕演化規(guī)則。細胞自動機演化規(guī)則的確定要考慮細胞鄰域、響應、屬性及時間和空間上的相互影響。根據(jù)斷溶體地質(zhì)研究認識，斷控巖溶以斷裂破碎帶為核心發(fā)生溶蝕擴大，其溶蝕程度的差異受斷裂分布、裂縫密度、地層巖性、大氣淡水或底部熱液接觸時間等綜合因素影響。由于這些因素涉及地質(zhì)歷史時期的地質(zhì)過程，參數(shù)多，實際測量難度大，且量綱不同，因此，采用基于概率的這種無量綱方法設計溶蝕演化規(guī)則。

溶蝕演化規(guī)則的基本假設：以自身網(wǎng)格為基礎(chǔ)，三維空間相鄰的網(wǎng)格為鄰域，每個網(wǎng)格的屬性有酸濃度c、溶蝕狀態(tài)f、連通概率pd、不同溶蝕狀態(tài)接觸時連通概率pf、裂縫發(fā)育程度pp、背景巖石相的溶蝕容易程度pr，所有參數(shù)均為0～1的概率值。

演化規(guī)則：如果待模擬網(wǎng)格周圍的酸濃度比自己大，則該網(wǎng)格的濃度變大，周圍網(wǎng)格的濃度變小，同時滿足質(zhì)量守恒要求。酸的擴散考慮不同方向的連通概率、網(wǎng)格的溶蝕狀態(tài)、裂縫發(fā)育概率、背景巖石相溶蝕容易程度。通過不同類型、不同方向的系數(shù)體現(xiàn)不同類型地質(zhì)因素對溶蝕過程的影響，基于這些基本規(guī)則可以設計具體的溶蝕過程算法。

2.3 溶蝕過程算法設計

假設模擬的總時間步為N，前N1個時間步地層處于暴露狀態(tài)，通過大氣降水向地層表面提供含有一定濃度酸的巖溶水，從N1到N時間步為非暴露，不再有地表巖溶水供給，地層內(nèi)部的巖溶水繼續(xù)溶蝕和演化，一直達到平衡狀態(tài)，不再溶蝕。

為了描述具體的細胞自動機演化算法，需要進一步定義一些參數(shù)。假設某個網(wǎng)格酸濃度為c0(無量綱方式，最大值為1，最小值為0)，溶蝕狀態(tài)為f0，累積酸溶蝕為acu，其鄰域為前、后、左、右、上、下6個方向網(wǎng)格，上部網(wǎng)格酸濃度為c1，下部網(wǎng)格濃度為c2，左側(cè)網(wǎng)格酸濃度為c3，右側(cè)網(wǎng)格酸濃度為c4，前方網(wǎng)格酸濃度為c5，后方網(wǎng)格酸濃度為c6，上述幾個方向的連通概率分別為pd1、pd2、pd3、pd4、pd5、pd6，上述幾個方向不同溶蝕狀態(tài)接觸時的連通概率分別為pf1、pf2、pf3、pf4、pf5、pf6(圖2)。不同方向的連通系數(shù)和不通過溶蝕狀態(tài)接觸的連通概率需要預先給定，由于重力的影響，通常設置為pd1>pd3=pd4=pd5=pd6>pd2，背景巖石相影響的溶蝕概率為pr，已溶蝕與未溶蝕細胞接觸時連通概率遠小于已溶蝕與已溶蝕細胞接觸時的連通概率。每個時間步，暴露在地表的網(wǎng)格酸濃度重置為1，表示有充足的酸供給量，然后其他網(wǎng)格按演化規(guī)則同步更新。任意網(wǎng)格在時間維度上累積酸濃度表示為該網(wǎng)格的溶蝕程度，是模擬的主要目標。

圖2 不同方向酸濃度、連通概率及相接觸連通性示意圖Fig.2 Schematic diagram of acid concentration, connectivity coefficient and contact connectivity in different directions

演化的核心規(guī)則是在每個時間步每個網(wǎng)格計算酸濃度的增加量dc并更新累計酸溶蝕和溶蝕相f。酸濃度增加量的初始值為dc=0，依次考慮6個方向的濃度影響：

如果c1>c0，則dc+=pd1pp1pf1pr1c1，c1=(1-pd1pp1pf1pr1)c1；

如果c2>c0，則dc+=pd2pp2pf2pr2c2，c2=(1-pd2pp2pf2pr2)c2；

如果c3>c0，則dc+=pd3pp3pf3pr3c3，c3=(1-pd3pp3pf3pr3)c3；

如果c4>c0，則dc+=pd4pp4pf4pr4c4，c4=(1-pd4pp4pf4pr4)c4；

如果c5>c0，則dc+=pd5pp5pf5pr5c5，c5=(1-pd5pp5pf5pr5)c5；

如果c6>c0，則dc+=pd6pp6pf6pr6c6，c6=(1-pd6pp6pf6pr6)c6。

更新c0、acu、f：

c0+=dc，acu+=c0；

如果c0>1，則f+=1。

基于細胞自動機的斷溶體溶蝕模擬方法的步驟總結(jié)如下(圖3)。

圖3 基于細胞自動機的斷溶體模擬主要步驟Fig.3 Main steps of fault-controlled karst simulation based on cellular automata

步驟1設置初始溶蝕狀態(tài)、各個方向連通系數(shù)、各個方向不同溶蝕相接觸時的連通性、每個網(wǎng)格的裂縫發(fā)育密度，模擬總的時間步數(shù)N，以及暴露的時間步數(shù)N1，其中N1

步驟2在時間步t，設置表層網(wǎng)格的酸濃度，如果時間步t

步驟3生成隨機路徑。

步驟4沿隨機路徑訪問所有網(wǎng)格，如果是表層網(wǎng)格跳過，則否執(zhí)行下列步驟。

步驟5根據(jù)上述演化規(guī)則計算該網(wǎng)格的濃度c0，并更新累計酸溶蝕acu及溶蝕狀態(tài)f。

步驟6判斷是否訪問完所有節(jié)點，如果是則執(zhí)行步驟7，否則再次執(zhí)行步驟4～步驟6。

步驟7判斷時間步是否等于N，如果是則結(jié)束模擬，否則繼續(xù)執(zhí)行步驟2～步驟7。

3 應用實例

3.1 地質(zhì)背景

研究聚焦在塔河油田10區(qū)西南部與托甫臺區(qū)東北部交匯處的一個典型斷溶體單元，深大斷裂與局部構(gòu)造共同控制，儲集體整體發(fā)育，規(guī)模大，油氣富集程度高。該井組經(jīng)過多年的開采，表現(xiàn)出累產(chǎn)高、水體能力強、連通關(guān)系復雜的特征?；谇捌诳坍嫷膬w形態(tài)，雖對井間關(guān)聯(lián)性有了一定的認識，但仍存在許多目前地質(zhì)認識無法解釋的動態(tài)現(xiàn)象。為進一步挖掘剩余油，優(yōu)化開發(fā)生產(chǎn)對策，開展溶蝕過程模擬有助于明確縫洞儲集體空間配置關(guān)系和連通樣式，為方案調(diào)整提供地質(zhì)依據(jù)。

3.2 溶蝕過程模擬

地質(zhì)模型維數(shù)為170×251×410，網(wǎng)格大小為50 m×50 m×1 m。通過地震螞蟻體、斷裂系統(tǒng)分析等手段獲取離散裂縫網(wǎng)絡，并把離散裂縫粗化到網(wǎng)格系統(tǒng)作為初始狀態(tài)，設置裂縫所穿越網(wǎng)格為溶蝕相，其他位置為非溶蝕相，然后通過基于細胞自動機的方法模擬斷溶體的溶蝕過程演化。

設置總的模擬時間步為N=300，前期暴露時間步為N1=250，后期50個時間步處于埋藏狀態(tài)。設置不同方向的巖溶水溶蝕性連通概率，pd1=0.9、pd2=0.001、pd3=0.1、pd4=0.1、pd5=0.1、pd6=0.1，設置溶蝕相與周圍網(wǎng)格的連通概率pf1=0.9、pf2=0.01、pf3=0.5、pf4=0.5、pf5=0.5，pf6=0.5，設置非溶蝕相與周圍網(wǎng)格的連通性均為0.005。模擬結(jié)果基本可以再現(xiàn)巖溶水在斷裂系統(tǒng)等因素控制下對碳酸鹽巖地層的溶蝕過程，溶蝕程度演化過程復合地質(zhì)規(guī)律，可以清楚反映不同溶蝕程度儲層的發(fā)育位置，對認識優(yōu)勢儲層分布、儲集體連通關(guān)系具有重要意義(圖4，圖5)。

圖4 不同時間步剖面演化特征Fig.4 Profile characteristic of different evolution step

圖5 不同模擬時間步平面演化特征Fig.5 Plane characteristic of different evolution step

4 結(jié)論

(1)根據(jù)礦場實踐特征總結(jié)了斷控巖溶演化規(guī)則，并首次將細胞自動機的方法引入斷控巖溶過程模擬，設計規(guī)則均基于實際地層演化特征，能夠客觀反映形成過程本質(zhì)，具有較強的實用性。

(2)傳統(tǒng)地質(zhì)統(tǒng)計學方法相比，該方法考慮了成因過程，對于井資料少、地震資料品質(zhì)低的斷溶體發(fā)育區(qū)，具有預測準確度高的優(yōu)勢。

(3)初始條件及溶蝕參數(shù)對模擬結(jié)果有重要影響，為了使模擬結(jié)果與實際數(shù)據(jù)吻合，需要借助反演或優(yōu)化的方法，進行下一步深入研究。