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

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

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

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

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

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

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

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

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

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

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

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

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

2 細(xì)胞自動(dòng)機(jī)巖溶模擬方法

2.1 細(xì)胞自動(dòng)機(jī)基本原理

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

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

2.2 巖溶過(guò)程細(xì)胞自動(dòng)機(jī)建模

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

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

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

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

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

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

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

2.3 溶蝕過(guò)程算法設(shè)計(jì)

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

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

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

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

如果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。

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

圖3 基于細(xì)胞自動(dòng)機(jī)的斷溶體模擬主要步驟Fig.3 Main steps of fault-controlled karst simulation based on cellular automata

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

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

步驟3生成隨機(jī)路徑。

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

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

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

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

3 應(yīng)用實(shí)例

3.1 地質(zhì)背景

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

3.2 溶蝕過(guò)程模擬

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

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

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

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

4 結(jié)論

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

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

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