儲層裂縫隨機建模方法研究進展

董少群 曾聯波* Xu Chaoshui曹 菡 王圣嬌 呂文雅(①中國石油大學(北京)地球科學學院，北京 102249； ②油氣資源與探測國家重點實驗室，北京102249； 阿德萊德大學土木、環境和采礦工程學院,澳大利亞阿德萊德 5005； 中國科學院對地觀測與數字地球科學中心，北京 100094)

1 引言

裂縫是指巖石受到構造變形發生破裂作用或物理成巖作用形成的沒有明顯位移的面狀不連續體，常呈網狀分布[1]。裂縫在很多領域均有涉及，其空間展布研究對于油氣勘探開發、地熱開發、采礦、巖土工程、水文地質等領域均有重要意義。在油氣田勘探開發領域常稱其為“裂縫”，在巖土工程、采礦、水文地質等領域則習慣稱其為“裂隙”，兩者中文名稱雖不同，但其所指一般等同[2]。為便于敘述，本文統一將其稱為裂縫。裂縫是巖體中流體移動的通道，如裂縫網絡是儲層中油氣、地熱系統中熱液、水文地質工程中地下水等流體運移的通道。

裂縫建模是一種常用的研究裂縫空間展布的手段，可分為確定性建模和隨機建模。確定性建模是根據已知信息建立確定的裂縫模型，如通過地震資料解釋出規模較大的裂縫，該方法一般不適用于規模較小的裂縫，且不能較好地綜合利用多種資料。裂縫隨機建模則是利用裂縫的先驗信息，通過隨機模擬方式生成可選的相同概率裂縫模型。這類方法不僅滿足已知點的裂縫統計學特征，而且承認未知區域裂縫發育的隨機性，較好地尊重了裂縫模擬不確定性的客觀事實[3]。近年來裂縫隨機建模方法應用較為廣泛，且已在油氣開采、地熱開發、采礦等領域均取得較好的應用效果。

目前常用的隨機裂縫建模方法大致可分為五類： ①基于空間剖分的裂縫建模； ②離散裂縫網絡建模； ③基于變差函數的裂縫建模； ④基于多點地質統計學的裂縫建模； ⑤基于分形特征迭代的裂縫建模。本文結合實例介紹五類方法的原理、優缺點及適用性等，同時討論這些方法在應用中存在的問題，并在此基礎上總結了裂縫隨機建模方法研究的四個發展趨勢。

2 基于空間剖分的裂縫建模

假設巖石中單個裂縫為平面多邊形，所有裂縫多邊形延伸后的平面便可將巖石剖分成多個巖塊，即進行空間剖分。裂縫建模時，若空間剖分合理，則這些平面上的某些多邊形即為裂縫?；诳臻g剖分的裂縫建模方法便是這種建模思路的一種實現，常用模型包括正交模型、馬賽克鑲嵌模型和Veneziano模型等[4]。

2.1 正交模型

正交模型是針對巖石裂縫系統進行幾何模擬的一種建模方法，其早期應用于水文地質等領域[5]。實際應用中的正交模型多為三維模型，初期的正交模型為無邊界正交模型，它將巖石裂縫系統中的裂縫簡化為三個方向上的無限正交面(圖1a的灰色平面)。不同方向裂縫面間距分別為S1、S2、S3，可通過隨機模擬獲得。圖1b為澳大利亞塔斯馬尼亞的棋盤石[6]，通過野外露頭測量統計不同方向的裂縫間距，分別獲得S1、S2、S3的分布函數，隨后便可通過蒙特卡洛隨機模擬得到隨機正交裂縫模型。該方法適用于正交的高角度裂縫系統建模，且裂縫規模相對于研究區范圍足夠大。該方法的優點是簡單易行，但實際裂縫系統并非都為如此規整正交裂縫網絡(圖1d中的延河剖面露頭、圖1f中的加利福尼亞的碳酸鹽巖露頭剖面[7])。針對此類問題，有學者[4，8，9]將無邊界正交模型進行了改進，提出有邊界正交模型(圖1c、圖1e)。將圖1a中剖分出的小網格進行編號，然后依據一定的規則選擇某些小網格或小網格的一部分作為裂縫。

圖1 正交模型及與其相似露頭剖面(a)無邊界正交模型； (b)塔斯馬尼亞較規整的棋盤石(正交裂縫)[6]； (c)有邊界正交模型；(d)延河剖面露頭； (e)有邊界正交模型； (f)加利福尼亞圣巴巴拉市碳酸鹽巖露頭剖面[10]

2.2 馬賽克鑲嵌模型

二維馬賽克鑲嵌模型，先在平面上生成一系列散點，然后通過相鄰兩點的中垂線將平面分割成彼此鑲嵌的凸多邊形“馬賽克”[11]，圖2a是其俯視圖。當給予每個多邊形一定高度后便可形成三維馬賽克鑲嵌模型(圖2a)，多邊形柱狀體相鄰面即為裂縫，這種方法可較好地模擬柱狀節理(圖2b)。另一種三維馬賽克實現方式如圖2c所示，它直接在空間中生成一系列散點，然后通過相鄰兩點的中垂面將三維空間剖分維彼此鑲嵌的凸多面體“馬賽克”。此時巖石被剖分成斷塊或基質巖塊，按一定規則選定凸多面體的部分面為裂縫，所形成的裂縫網絡即為所得的裂縫模型。野外露頭及巖心常見此類裂縫，如圖2d所示。在利用這種方法建立裂縫模型時，生成平面和空間中的散點十分關鍵，因為它直接決定了模型是否可用。

2.3 Veneziano模型

Veneziano模型也稱泊松平面模型(Poisson flat model)，是由Veneziano于1978年提出，其建立過程如圖3所示，可分為三步[4，12]： ①泊松面過程(圖3a)，在空間中生成隨機平面Ai(i=1,2,…,n)，n為預定義的參數； ②泊松線過程(圖3b)，通過泊松過程生成平面Ai上的隨機直線，將平面分割成多邊形集Pi(i=1,2,…,m)，m為第i個平面上的多邊形的個數； ③依概率選取多邊形作為裂縫面(圖3c)。整個模型生成過程由裂縫密度決定[12]，步驟①反映了主應力走向及研究區的整體裂縫密度[13]，由P30、P31或P33(表1)控制，步驟②、步驟③通過裂縫多邊形的大小和位置反映裂縫密度的變化，步驟②由P20、P21或P22(表1)控制。前面提及的正交模型、馬賽克鑲嵌模型均可視為Veneziano模型的一種特例。

圖2 馬賽克鑲嵌模型及與之相似的露頭剖面(a)二維馬賽克鑲嵌擴展模型； (b)濟州島柱狀節理； (c)三維馬賽克鑲嵌模型； (d)塔河油田一間房組的巖心

圖3 Veneziano模型(a)泊松面過程確定隨機面； (b)泊松線過程確定隨機線； (c)確定裂縫面； (d)三維裂縫模型

表1 裂縫密度、強度、孔隙度參數

該模型提出后，一些學者陸續對該模型做了改進。Dershowitz等[8，9]將該模型的后兩步進行了改進，提出了Dershowitz模型。Einstein、Ivanova等[4]在Veneziano模型(圖3d)的基礎上進行了改進，通過將多邊形進行旋轉及篩選等操作將Veneziano模型轉換為離散裂縫網絡模型。Ivanova等也嘗試對步驟②進行改進，使用1.2節中介紹的Voronoi圖法將平面進行剖分[13]，提高了生成符合條件的裂縫多邊形的效率。Mosser等[14]也對該模型進行了研究，認為直接利用泊松線過程建立的裂縫模型過度估計了裂縫的連通性，并且針對此問題提出了STIT方法。STIT方法在泊松線過程的基礎上，對其剖分出的網格進行進一步剖分，Mosser等[14]認為該方法可以緩解過度估計裂縫連通性的問題。

3 離散裂縫網絡建模

離散裂縫網絡(Discrete fracture networks，DFN)最早由Baecher等[12，15]引入。如圖4a所示，Baecher模型假設裂縫為一圓盤(圓形或橢圓形)，各圓盤之間可以交叉，是一種基于示性點過程的隨機建模方法[16，17]。其實現過程分為兩步：點過程、示性過程。點過程確定裂縫的中心位置； 示性過程確定點的屬性，如裂縫大小、傾角、傾向、開度、滲透率等屬性[17]。

自Baecher模型被引入以來，便得到了廣泛應用，并且出現了不同形式的改進，主要體現在以下三個方面。

(1)裂縫面形狀表示方法的改進。三維裂縫面最初假設為圓盤，因為裂縫產生時最初近似為圓形或橢圓形，隨著破裂過程的繼續，受層面等多方面限制，其形狀變得相對不規則，故后續有學者相繼將其改進為多邊形等，并對裂縫的期次(如裂縫多邊形之間的的切割、限制關系)做了相關改進(圖4b～圖4e)。圖4b為增強Baecher模型，它將橢圓換成近似橢圓的多邊形，而且裂縫可被其他裂縫截斷； 圖4c為BART(Baecher Algorithm Revised Terminations)模型，其對裂縫終止過程做了進一步改進[9]； 圖4d為泊松矩形模型，是一種簡化的增強Baecher模型，其裂縫通過矩形表示，而不再是原來近似圓形的多邊形[4]； 圖4e中模型為隨機多邊形模型，它將裂縫面定義為邊數固定的隨機多邊形[17，20]，目前很多裂縫建模商業軟件中已采用該方案，可預設多邊形邊數，如FracMan、Petrel等。

(2)裂縫中心位置的產生規則的改進。裂縫位置的確定一般通過穩態泊松過程來實現， 為擴展其適用范圍，其他方式的裂縫位置產生規則相繼被引入，如基于非穩態泊松過程、馬爾科夫鏈蒙特卡洛(Markov Chain MonteCarlo，MCMC)[21]、分形[18，22]、Cluster過程[17]、Cox過程[17]等(圖4f～圖4j)。圖4f中模型為最近鄰模型，是一種半概率的、基于模式的模型。裂縫易分布于斷層或大規模裂縫附近，該模型模擬該趨勢，優先在斷層或大規模裂縫附近生成裂縫[4]。圖4g中模型為戰區模型，也是一種半概率的、基于模式的模型。它通過優先在近似平行的裂縫之間生成裂縫來模擬剪切帶的幾何形態[23]。圖4h～圖4j分別對應Levy-Lee分形模型[24]、母女模型(Parent-daughter model)[25，26]、二元分形裂縫網絡(Binary fractal fracture network，BFFN)模型[18，22]，它們是考慮裂縫位置的分形特征的改進DFN模型。

圖4 離散裂縫網絡模型[4，9，11，13，17-19](a)Baecher模型； (b)增強Baecher模型； (c)BART模型； (d)泊松矩形模型； (e)隨機多邊形模型； (f)最近鄰模型； (g)戰區模型；(h)Levy-Lee分形模型； (i)母女模型； (j)二元分形裂縫網絡模型； (k)非平面帶模型； (l)GEOFRAC模型

(3)裂縫屬性分布的改進。很多學者對裂縫的走向、傾向、傾角、開度、孔隙度、滲透率、間距等參數的分布做了大量研究工作。如描述裂縫長度常使用指數分布[27]、對數正態分布[15，28]及伽馬分布[29]等； 描述裂縫走向的有von-Mises分布[17]； 描述裂縫傾向、傾角的有Fisher分布； 描述開度的有指數分布、對數正態分布[30]、與裂縫長度呈線性或非線性關系[31]等； 描述裂縫間距的有指數分布[12]、對數正態分布等； 描述裂縫跡長的有對數正態分布[12]、 指數分布等。

對DFN模型進行改進時，常兼顧上述多個方面，如圖4h中的非平面帶模型，是一種半隨機模型，其裂縫位置和走向隨設定曲面變化[4]，對裂縫面的形狀及裂縫中心位置的產生規則均有所改進。

建模方法的主要原理相同，即示性點過程。除了基于示性點過程的離散裂縫網絡建模外，Ivanova等[4，13，32]在Veneziano模型的基礎上進行了改進，也可以生成離散裂縫網絡模型，稱為GEOFRAC模型(圖4l)。限于篇幅，本文僅選取幾組有代表性的DFN模型進行介紹。

3.1 Baecher模型

Baecher模型也稱為隨機圓盤模型(Random-disk model)[12]，它通過點過程和示性過程共同確定裂縫網絡中的各條裂縫。原始Baecher模型有五個假設條件[4，12]： ①單個裂縫為平面圓盤； ②裂縫中心位置符合隨機分布，且各裂縫中心位置彼此獨立； ③圓盤半徑符合對數正態分布； ④裂縫半徑與傾角不相關(統計角度上獨立)； ⑤裂縫半徑與裂縫中心位置不相關(統計角度上獨立)。由于構造應力作用，實際儲層的裂縫常具有一個或多個優勢產狀(傾向、傾角)?；诖?，Baecher模型建立過程如圖5所示，分為四步。

(1)利用極點圖劃分裂縫組系。如圖5a所示，半徑對應裂縫法線與鉛垂線夾角(即裂縫傾角)的正弦，方位角對應裂縫法線在水平面的投影與正北方的夾角(即傾向)，顏色代表裂縫的數量，暖色調表示裂縫數量多，可以看出存在兩組高角度裂縫，走向分別為北北東和北北西。

(2)點過程確定裂縫位置。以第①組裂縫為例，通過下式中均勻分布[33]確定隨機數(x,y,z)，即為裂縫中心位置坐標

(1)

式中V為研究區體積。該組裂縫圓盤的中心位置如圖5b所示。

(3)示性過程確定每條裂縫的屬性(如傾向、傾角、走向、大小、開度等)。通過蒙特卡洛隨機模擬產生的隨機數作為分位數，從累計分布函數上獲得分位數對應的裂縫屬性值。裂縫產狀(傾向、傾角)模擬可通過模擬裂縫面法線來實現[17，33-36]，確定裂縫面法線后再將其轉換為裂縫面的傾向和傾角。除此之外，產狀也可通過直接模擬裂縫面的傾向、傾角來實現[17，19，33-35]，即首先擬合裂縫傾向、傾角所符合的分布，然后直接對兩個分布進行隨機模擬，獲取相應的傾向和傾角數值。

下面將以通過裂縫面法線確定裂縫產狀的方式進行介紹。利用von-Mises分布(式(2))擬合裂縫傾向，利用均勻分布(式(3))擬合裂縫傾角，利用對數正態分布(式(4))擬合裂縫圓盤半徑和裂縫開度。

(2)

(3)

式中a、b分別表示裂縫傾角的最大值和最小值。

(4)

式中：σ體現圓盤半徑或開度的變化程度，σ越大，分布函數越平緩；μ決定半徑或開度分布的峰值位置。

通過圖5a中第①組裂縫的參數統計結果進行概率密度估計可獲得式(2)～式(4)中的參數。隨后便可得到對應的累積概率密度曲線，如圖5c～圖5f所示。

通過蒙特卡洛隨機模擬產生四個隨機數作為四個分布的分位數，然后找到對應累積概率密度曲線的橫坐標即為裂縫屬性，如圖11中虛線所示，傾向、傾角、半徑、開度分別為61°、84°、1.2m、0.0192mm。重復上述操作3次后，生成的裂縫如圖5g所示，顏色表示開度的大小。繼續重復上述過程直至滿足終止條件即可生成離散裂縫網絡模型。終止條件多與裂縫密度有關，如表1所示。常用的終止條件有： ①單位體積中裂縫的條數P30達到預設值； ②單位體積中裂縫的面積P32達到預設值； ③單位體積中裂縫的體積P33達到預設值。這些預設值常通過巖心、測井、露頭等測得的裂縫密度給出。完成第一組裂縫模擬后，用同樣的方式可以實現另一組裂縫的建模。

圖5 Baecher模型建模流程(a)極點圖確定裂縫組系； (b)點過程確定裂縫中心位置； (c)示性過程確定裂縫傾向、傾角；(d)示性過程確定裂縫半徑長度； (e)示性過程確定裂縫開度； (f)第①組裂縫模擬結果

(4)質量控制(QC)。裂縫網絡模型建完后，很重要的一步就是質量控制，而裂縫密度則是常用的一類衡量參數，如統計井軌跡通過P10(單位長度上裂縫的條數)或P11(單位長度上裂縫的長度)是否與井數據相對吻合，如果不吻合，則需要對模型進行調整。

3.2 最近鄰模型

按裂縫形成時構造應力場作用范圍，可將裂縫分為區域構造裂縫和局部構造裂縫。局部構造裂縫又可分為與斷層相關的裂縫和與褶皺相關的裂縫。與斷層相關的裂縫存在于斷層附近，與之伴生，其發育的應力場與斷層一致，且隨著與斷層距離的增大，裂縫的線密度明顯降低[1，37](圖6[38])。最近鄰模型是一種改進的Baecher模型，它對裂縫位置做了改進，假設裂縫密度與裂縫和主裂縫的距離有關，它隨著距離的增大呈指數減小，從而使裂縫更易于分布于斷層或大規模裂縫附近[4]。

針對上述現象，高志勇等[37]對庫車前陸盆地裂縫與斷層關系進行了定量分析。圖7顯示裂縫線密度與距斷面距離成較好的指數關系，與最近鄰模型假設相近，故可用于此類裂縫建模及最后的裂縫模型效果檢驗。

圖6 裂縫分布與距斷層距離相關的剖面[38]

圖7 庫車前陸盆地裂縫與斷層關系[38，40]

使用最近鄰法建立與斷層相關的裂縫模型時，裂縫的性質可以參考與斷層相關的裂縫模式[39]。通過分析斷層的應力特征及發育裂縫的力學特征可以獲得裂縫的走向、傾向、傾角等相關信息，也可以輔助檢驗輸入模型參數合理性及建立的裂縫網絡模型的有效性。

在建立與斷層相關的裂縫模型時，也可以使用母女模型[25，26]，同樣能達到與最近鄰模型類似的效果，母女模型除了考慮斷層附近裂縫密度高的特性外，還考慮了裂縫的分形特征，本文將在3.3節中進行介紹。

3.3 分形模型

利用了裂縫分形特征的DFN模型均可看作是分形模型，認為裂縫長度符合冪律分布[33，41]。除考慮裂縫長度的分形特征外，裂縫位置的分形特征也常被引入到DFN模型中，如Levy-Lee分形模型[24](圖4h)、母女模型[25，26](圖4i)、二元分形裂縫網絡模型[18，22](圖4j)等。

實際應用中，裂縫的分形特征也被學者廣泛認同[42-45]。很多學者認為晶質巖中的裂縫系統常呈現分形特征，具體體現在兩個方面：裂縫空間分布符合冪律分布； 裂縫長度符合冪律分布。BFFN模型正是基于上述兩點提出的。模型中裂縫中心位置按照“分形破碎”的方法產生，裂縫長度利用隨機過程結合冪律定理產生[18，22]。其具體步驟如下。

(1)利用“分形破碎”的方法，根據D2(式(5))確定裂縫中心位置

(5)

(2)示性過程確定裂縫屬性，如圖8h所示。裂縫長度根據式(6)確定，傾向和傾角依據Fisher分布(式(7))獲得

lgN=alglmax-algl

(6)

(7)

式中：N為裂縫長度不小于l的裂縫數量；lmax為最大裂縫的長度；φ、θ分別為裂縫面的傾向、傾角。

3.4 Veneziano改進模型

Ivanova等在離散裂縫建模方法研究和改進方面未沿用Baecher模型中示性點過程的方法，而是在Veneziano模型[4，12]基礎上進行了改進，將圖3d中生成的裂縫模型進行了平移、旋轉，從而反映裂縫位置和走向的局部變化，并將其稱為GEOFRAC模型[13]，如圖4l所示。原始GEOFRAC模型在裂縫生成過程中存在一個問題，泊松線過程產生的多邊形中有大量因其形狀不能較好地反映裂縫形態而被刪除，即裂縫多邊形生成效率低，如圖3b、圖3c所示。針對這一問題，Ivanova等用Delauney三角剖分代替泊松線過程對原始模型進行了改進，有效地改善了裂縫多邊形生成效率低的問題，他們將改進后的模型稱為新GEOFRAC模型[13]。新GEOFRAC模型的生成過程如圖9所示。

圖8 三維BFFN模型構建過程[18](a)～(g)裂縫中心位置產生過程； (h)基于示性點過程生成裂縫模型

圖9 新GEOFRAC模型中裂縫生成過程[13](a)通過泊松過程生成平面； (b)平面網格剖分； (c)多邊形隨機平移； (d)多邊形隨機旋轉

4 基于變差函數的裂縫建模

基于變差函數的建模方法在沉積相、儲層物性等建模中廣泛應用[46]，但在裂縫建模中應用相對較少。通過變差函數可以表征空間中不同位置裂縫的相關性[16，47-49]，因此不少學者也嘗試將其應用于裂縫建模。大部分工作是在2.1節中Baecher模型的框架內展開的，即利用變差函數模擬裂縫密度或確定裂縫位置，即代替點過程確定裂縫位置[16]。Dowd等[16]、Xu等[47]利用截斷高斯指示模擬(Truncated Gaussian indicator simulation)和截斷多高斯場截斷模擬(Truncated Pluri Gaussian simulation)[50]對裂縫位置進行模擬，并結合示性過程進行裂縫建模，實驗結果表明截斷高斯指示模擬適用于裂縫系統相對簡單的模擬，截斷多高斯場截斷模擬可更好地對復雜裂縫系統中裂縫位置進行模擬。

變差函數除了可以用于代替點過程外，也可以用于代替示性過程。Koike等[48，49，51，52]利用通過指示變換和變差函數來模擬裂縫走向、傾向，其文章中將其方法也稱為GEOFRAC[48，49]，但根據發表先后，3.4節中方法要早，因此本文中GEOFRAC如果沒有特別注釋，則代指2.3節中的方法。Koike等提出的方法的操作步驟如圖10所示，可分為以下三步。

(1)確定小平面中心，如圖10a所示。該方法不直接模擬裂縫，而是先生成尺度較小的小平面，然后再將小平面連接成最終裂縫面。在確定裂縫中心時，先將儲層進行網格化，然后利用序貫高斯模擬(Sequential Gaussian Simulation，SGS)模擬每個網格的裂縫密度，結合裂縫密度通過隨機過程在各網格中確定是否有小平面中心，若有，則確定其中心位置。

(2)確定小平面走向和方位角(圖10b～圖10e)。將走向和傾角劃分為16個區域，并進行指示變換，I=(g1,g2,…,gL)，此時L為16，當走向和方位角位于第i個區域時，gi=1，其他元素為0。此處Koike等定義的走向符合右手定則，范圍為[0,2π]。通過普通克里金(Ordinary Kriging，OK)計算某個小平面中心位置處的I,如I=(0.2,0.9,…,0.05)，則其屬于第2個區域。隨后再根據該區域內裂縫的走向和方位的分布函數隨機生成裂縫走向和方位。在普通克里金插值時，為了減少計算量，該方法引入了主成分分析。

(3)將小平面連接形成裂縫面(圖10c～圖10f)。首先判斷兩個小平面是否連接，它可以通過兩個小平面的走向差、傾向差和小平面中心距離判斷是否連接，也可以通過投影后的走向差和中心距離來判斷兩個小平面是否連接。隨后將孤立的小平面轉換成圓盤，作為孤立裂縫，將相互連接的小平面連接成裂縫面，如圖10g，可以直接將小平面中心連接成不規則三角形格網(Triangulated Irregular Network，TIN)，也可以根據小平面的走向、傾角的均值計算裂縫面，并用小平面的包絡作為裂縫面的邊界。最終生成的裂縫模型，如圖10h所示。

圖10 基于變差函數的裂縫建模(據文獻[48， 49]，有改動)(a)確定小平面中心； (b)～(e)確定小平面走向和方位角； (f)連接小平面形成裂縫面； (g)、(h)離散裂縫網絡模型

5 基于多點地質統計學的裂縫建模

多點地質統計學較以變差函數為核心的兩點地質統計學而言，易于表征復雜的空間結構、再現復雜目標的幾何形態[46]。多點地質統計學包括迭代的和非迭代的兩大類方法。迭代方法有基于模擬退火的方法[53]、基于吉布斯取樣的后處理迭代方法等，但這類方法普遍存在迭代收斂問題[46]。文獻中提及的多點地質統計學方法多指非迭代的方法，該類方法最早由Guardiano等提出，其通過訓練圖像掃描獲取待估點的概率分布，進而模擬沉積相[54]，但在掃描圖像時運算量及內存消耗都較大，不同學者對其嘗試進行改進。目前常用的改進算法有SNESIM(Single Normal Equation Simulation)算法[55， 56]、SIMPAT算法[57-59]、FILTERSIM 算法[60]、DisPAT算法[61,62]等。SNESIM使用一個平衡方程進行概率估計，一次性將訓練圖像的條件概率存儲于“搜索樹”中，極大地減少了機時。SIMPAT利用訓練圖像獲取地下結構的模式，采用相似性方法對地下儲層進行圖像恢復和再現。FILTERSIM為基于模式濾波的多點地質統計模擬方法。DisPAT算法構建數據樣式的距離矩陣，應用多維尺度分析進行降維，再用K-means聚類建立樣式聚類，后續模擬與SIMPAT相同[62]。這些方法從不同角度對原始方法改進，以提高效率和拓展應用范圍。

多點地質統計學早期主要用于沉積相模擬，近年來正逐漸被應用于裂縫建模，對于裂縫建模，野外露頭剖面、地震解釋的斷層等都可以作為裂縫幾何形態的先驗認識，即訓練圖像，同時測井解釋或巖心觀察的裂縫可以為硬數據。Dowd等[16]嘗試使用SNESIM算法對尤卡山二維裂縫剖面進行模擬，并與其他方法進行對比，認為多點地質統計學模擬可以體現裂縫的空間相關性，并且能模擬復雜的裂縫網絡。謝青[63]嘗試使用多線程并行運算提高SNESIM算法效率，并對頁巖儲層天然裂縫進行了三維空間建模。Mohammadmoradi[64]、張烈輝等[65]、Jia等[66]針對FILTERSIM算法在裂縫建模中的問題分別進行了改進，并進行了二維裂縫建模。Liu等[67]使用SIMPAT算法對二維裂縫剖面進行建模并數模，結果表明所建立模型可以較好地描述裂縫系統。利用多點地質統計學進行裂縫建模，普遍存在穩態的問題，針對這一問題，Chugunova等[68]提出一種非穩態多點地質統計學的模擬方法，該方法除了使用裂縫圖像作為訓練圖像外，還加入了裂縫密度圖作為輔訓練圖像進行裂縫建模。隨后，Chugunova等[69]利用多點地質統計學建立裂縫網絡模型，并結合動態數值模擬對模型進行了驗證。其建模具體流程如圖11所示，具體分為三步：①訓練圖像構建(圖11a～圖11c),通過遙感圖像獲取巖性及裂縫，并選擇某一相同巖性的區域作為訓練圖像; ②確定硬數據； 利用遙感數據解釋明顯的斷層作為硬數據； ③結合訓練圖像和硬數據對待估點進行預測； ④通過生產數據驗證模型的有效性。

圖11 基于多點地質統計裂縫建模(據文獻[69]，有改動)(a)遙感圖像； (b)提取的不同巖性區域的裂縫； (c)訓練圖像砂巖4對應的裂縫； (d)硬數據(以(c)中

6 基于分形特征迭代的裂縫建模

基于分形特征[70，71]迭代的裂縫建模的常用方法有迭代函數系統(Iterated Function System， IFS)[72-74]、L-System、謝爾賓斯基三角形等，本章重點介紹前兩種。

IFS將待生成的圖像看作由許多與整體相似的(自相似)或經過一定變換與整體相似的(自仿射)小塊拼貼而成。而裂縫網絡常表現出自相似性或自仿射性，所以IFS也被用于裂縫網絡模擬。Acuna等[72，75]利用IFS對二維和三維裂縫系統進行建模，如圖12a、圖12b所示，該方法首先對初始點群進行一系列迭代數值變換，隨后每次迭代中，應用該系統中的函數對點群進行轉換、映射、旋轉、收縮、扭曲。當圖像中的點群符合分形目標時便終止迭代，此時即可獲得最終裂縫分布圖像[46，75]。周德華等[74]將IFS應用于油氣藏天然裂縫網絡建模，方法實施包括三部分： ①裂縫網絡的迭代生成，根據生成模型的分形維數Df與實際地層裂縫參數空間Dr的分形維數的差異決定迭代發生器的選擇和迭代的終止； ②裂縫網絡中流體流動的有限元模擬； ③通過與試井資料的匹配關系，利用模擬退火算法優化IFS參數，進而優化天然裂縫網絡的幾何結構。

L-System是描述植物生長的數學模型，其基本思想可解釋為理想化的樹木生長過程：從一條樹枝(種子)開始，發出新的芽枝，而發過芽枝的枝干又全部發新芽枝，直至最后長出葉子。該方法與裂縫生長過程較為相近，故可嘗試應用于裂縫建模。根據L-System 的基本思想，先建立二維樹的數學模型：從樹干開始，然后沿著樹干逐漸擴展到連接的樹枝，再以遞歸的方式進行同樣的過程，該過程持續到最終分枝。Zhou等[76，77]結合微地震數據、生產數據對非常規儲層中的裂縫系統利用L-system進行模擬，該方法將微地震數據解釋為與之匹配的分形網絡，進而解釋復雜裂縫網絡的分叉及多尺度特征。注水誘導裂縫的形成過程近似“從樹干逐漸擴展到連接的樹枝”，因此也有學者嘗試將其應用于注水誘導縫的模擬。Fan等[78]基于分形特征對分形誘導裂縫網絡進行模擬，并結合數值模擬進行了驗證。

圖12 基于分形特征迭代的裂縫模型(據文獻[72， 75， 76， 78]，有改動)(a)二維IFS模型； (b)三維IFS模型； (c)二維L-system模型； (d)三維分形誘導裂縫網絡示意圖

7 結論和討論

儲層裂縫分布的復雜性使得儲層裂縫建模成為一大難點，目前隨機裂縫建模方法是較常應用的一類建模方法。本文系統地調研了現有的裂縫隨機建模方法，將其分為五大類，他們的優缺點、適用性、發展趨勢概括如下。

(1)基于空間剖分的裂縫建模

基于空間剖分的裂縫建模原理簡單，便于操作，適用于高角度正交裂縫、柱狀節理、巖石破裂等，但不適用于復雜裂縫建模。無邊界正交模型和有邊界正交模型適合模式較為簡單的正交裂縫系統，如區域構造裂縫。馬賽克鑲嵌模型簡單易行，可較好地應用于玄武巖柱狀節理建模及巖石破碎模擬，但此類方法初始隨機點的確定準則是一個難點。Veneziano模型包含大部分的空間剖分的裂縫模型，但其生成符合裂縫形態的多邊形的效率較低。基于空間剖分的裂縫建模近年來發展緩慢，但仍有應用潛力。

(2)離散裂縫網絡建模

離散裂縫網絡模型中的Baecher模型是目前最常用的裂縫模型之一，簡單易行，便于后續數模。它可綜合多種信息將已有裂縫模式在研究區進行裂縫空間位置拓展，所建立的模型可包含較全的裂縫屬性(傾向、傾角、長度、開度等)，這也是其他模型所不具備的。該方法易于與地質模式相結合，進而形成更加符合地質認識的裂縫模型。例如在斷層附近的裂縫系統建模時可以通過改進點過程，使其符合遠離斷層裂縫密度低的規律；晶質巖中的裂縫系統常呈現分形特征，將分形特性融入到點過程和示性過程，即可較好地體現此類裂縫系統的特性；裂縫密度約束點過程可以使生成裂縫網絡更加逼近此類裂縫系統；此外也可將研究區應力狀態導致的裂縫模式融入。實際裂縫系統中裂縫展布常受層控，因為巖性變化對裂縫終止有較強影響，如沉積巖儲層中裂縫，DFN模型可以通過分層位、不同組系、不同規模裂縫的疊加來實現多尺度復雜裂縫系統的建模。

然而DFN模型在實際應用中也存在三個主要問題： ①需要預設裂縫的形狀、大小、傾向、傾角等的分布及其分布中的參數。這類參數可以通過露頭分析、巖心觀察、測井解釋等獲取，但常存在采樣不均一、不充分、數據需要校正等一系列問題； ②裂縫系統中裂縫在空間上常具有相關性，DFN模型未考慮這一點，這是該模型的不足，而這一點正是基于變差函數和多點地質統計學的裂縫建模的一大優點。③不同組系裂縫數量比例如何分配，這一問題也是前三類方法都面臨的問題。

近年來，DFN建模研究趨勢主要有以下幾個方面：①裂縫中心位置確定方法的研究，如點過程的改進(非穩態泊松過程、MCMC、Cluster過程、Cox過程等)、裂縫密度約束下的裂縫中心位置產生方法研究等；②示性過程中裂縫屬性的研究，如裂縫形狀的選擇、模型中裂縫參數(傾向、傾角、大小、開度等)的分布函數研究；③DFN模型中參數獲取問題的研究。獲取裂縫參數的手段越來越多，如通過遙感衛星數據獲取大尺度裂縫及斷裂帶信息[79，80]；基于三維激光掃描(如Light Detection and Ranging， LiDAR)技術對露頭進行掃描獲取三維數據[81-83]，并通過手動或自動的方式提取裂縫參數[80，84]；基于掃描線法的露頭裂縫信息獲取[85]；利用無人駕駛飛機(Unmanned Aerial Vehicle，UAV)獲取多尺度裂縫模式、長度等信息[86]；利用成像測井獲取單井裂縫信息[87]，如位置、傾角、開度等；通過巖心觀察統計單井裂縫信息；利用CT掃描獲取小尺度裂縫信息[88，89]，如裂縫幾何形態等。這些裂縫規模數千米到幾納米，不可避免的是裂縫多尺度耦合的問題，如何合理地將其轉換成研究問題模型中的參數，即參數校正成為一項重點工作[90]；④地質信息的融合，針對具體地質模式對DFN的改進，如與褶皺伴生的裂縫系統建模[91，92]、與斷層伴生的裂縫系統建?？紤]斷距的影響因素[93，94]等；⑤構建裂縫密度約束體的研究。建模時，裂縫參數的統計分布、露頭統計數據等可以反映裂縫的統計信息，但不能有效地指出裂縫更為發育的區域，而裂縫密度約束體正是起到這一作用，對裂縫中心位置的生成過程或裂縫方位確定過程進行約束。目前構建約束體有：通過疊前三維地震數據預測裂縫密度及方位[95，96]、弱度[97，98]、扁率、通過疊后地震數據邊緣檢測獲取裂縫帶的信息[99-101]、通過應力場數值模擬構建裂縫密度約束體[94，102]，井震結合地質建模等實現[93，103]；⑥裂縫模型有效性檢驗，即質量控制(QC)。如對生成的裂縫模型進行抽樣，統計一維密度、裂縫長度和二維密度，并與實測裂縫進行對比[104]，將裂縫連通性分析[20，105，106]結果與示蹤劑數據進行對比，并與實際生產數據進行對比[93]。根據對比結果對模型參數等進行調整，校正模型。重復上述過程直至滿足質量控制條件，否則模型不可用。

除上述基于示性點過程的離散裂縫網絡模型外，利用其他數學方法生成離散裂縫網絡的模型也逐漸得到重視，如Veneziano改進模型。

(3)基于變差函數的裂縫建模

基于變差函數的裂縫建模有兩類方法：一類是用沉積相建模流程直接對井上測井解釋的裂縫用地質統計學的方法進行，但由于裂縫非均質性遠超過沉積相的變化，因此不能較好地模擬裂縫形態，同時對裂縫屬性也不能較好地估計； 另一類研究多在DFN模型的框架下開展，因此它繼承了DFN模型的大部分特點。它用指示變換和變差函數來模擬裂縫屬性，其功能近似于Baecher模型中示性過程，但其考慮了裂縫的空間相關性，是對DFN模型的改進，而兩點地質統計學對復雜裂縫系統中裂縫的空間相關性表征不足是其面臨的一大問題。從發表成果看，近年來該方法在實用性方面有了較大改進。其研究趨勢除了DFN模型中提到的六點外，復雜裂縫系統空間相關性的描述也是亟待解決的問題。

(4)基于多點地質統計學的裂縫建模

基于多點地質統計學的裂縫建模，易于結合先驗裂縫模式信息，且忠實于硬數據，并在后續數值模擬方面具有較強優勢，適用于靜態裂縫建模。但訓練圖像獲取及算法效率提高是其發展過程中必須要解決的兩個問題。目前訓練圖像多為二維圖像，常通過衛星、露頭等獲得，如何獲取適合研究區裂縫建模的三維訓練圖像是近年來逐漸被關注的一個研究方向，如通過三維地震解釋斷層、通過應力場數值模擬三維斷層、通過三維CT獲取裂縫。獲取的訓練圖像普遍存在尺度耦合的問題，這也是該類方法需要解決的問題。同時訓練圖像穩定性問題是應用多點地質統計學進行模擬時要考慮的一個問題。算法效率提高也是多點地質統計學近年來的研究熱點，數據結構、降維、并行計算等手段相繼被引入，雖有所緩解，但仍待改進。

(5)基于分形特征迭代的裂縫建模

基于分形特征迭代的裂縫建模簡單易行，既可用于靜態裂縫模擬，也可用于動態模擬(如巖石破裂的模擬、注水誘導縫模擬等)。該方法僅適用于分形特征明顯的裂縫系統模擬。如何高效地生成符合約束條件(如微地震數據、裂縫密度數據)的分形模型是目前該方法的一個研究重點，如通過約束條件和優化算法改進迭代規則生成符合地質認識的裂縫模型。

綜上所述，近年來裂縫隨機建模方法的研究趨勢可歸納為四點： ①原有數學算法的改進及新算法的引入，如空間剖分方法中Delaunay三角剖分算法的改進、離散裂縫網絡模型中隨機過程穩定性的考慮、多點地質統計學中掃描效率的提高、分形特征迭代模型中迭代算法的改進等； ②將地質認識融合到數學模型中，建立針對某一特定地質體的裂縫系統，如與斷層伴生的裂縫系統建模、具有明顯分形特征的晶質巖裂縫等； ③確定模型參數的相關研究。首先是資料獲取技術的改進，如疊前三維地震、成像測井、CT掃描、三維激光掃描、無人駕駛飛機等手段的引入。其次是所獲取的不同尺度裂縫信息的耦合，將所獲取的裂縫信息轉換成模型中的參數； ④裂縫模型有效性驗證，即質量控制。

感謝毛哲、王小垚、曹東升等在本文成文過程中提供了幫助！

參考文獻

[1] 曾聯波．低滲透砂巖儲層裂縫的形成與分布．北京：科學出版社，2008．

Zeng Lianbo．Formation and Distribution of Fracture in Low-permeability Reservoir．Science Press，Beijing，2008.

[2] 邱殿明．斷裂、斷層、節理、劈理、裂隙、裂縫之間的關系小結．吉林大學學報(地球科學版)，2013，43(5)：1392.

[3] 潘欣．構造裂縫識別與建模研究——以鄂爾多斯盆地中西部定邊—華池地區延長組為例[學位論文]．陜西西安：西北大學，2010．

Pan Xin．Identifying and Modeling Research of Structural Fracture — A Case Study of the Yanchang Formation in Dingbian-Huachi Region of the Central-West Ordos Basin [D]．Northwest University，Xi’an，Shaanxi，2010.

[4] Ivanova V M．Three-dimensional Stochastic Modeling of Rock Fracture Systems [D]．Massachusetts Institute of Technology， 1995.

[5] Childs E C．The anisotropic hydraulic conductivity of soil．European Journal of Soil Science，1957，8(1)：42-47.

[6] Robinson N I，Sharp Jr J M，Kreisel I．Contaminant transport in sets of parallel finite fractures with fracture skins．Journal of Contaminant Hydrology，1998，31(1-2)：83-109.

[7] Dunne W M，North C P．Orthogonal fracture systems at the limits of thrusting：an example from southwestern Wales．Journal of Structural Geology，1990，12(2)：207-215.

[8] Dershowitz W S．Rock Joint Systems [D]．Massachusetts Institute of Technology，1984.

[9] Dershowitz W S，Einstein H H．Characterizing rock joint geometry with joint system models．Rock Mechanics and Rock Engineering，1988，21(1)：21-51.

[10] Bai T，Maerten L，Gross M R et al．Orthogonal cross joints：do they imply a regional stress rotation?Journal of Structural Geology，2002，24(1)：77-88.

[11] Alghalandis Y F．Stochastic Modelling of Fractures in Rock Masses [D]．University of Adelaide，2014.

[12] Baecher G B，Lanney N A，Einstein H H．Statistical description of rock properties and sampling．Proceedings of the 18th US Symposium on Rock Mechanics,Golden,CO,1977.

[13] Ivanova V M，Sousa R，Murrihy B et al．Mathematical algorithm development and parametric studies with the GEOFRAC three-dimensional stochastic model of natural rock fracture systems．Computers & Geosciences，2014，67：100-109.

[14] Mosser L J，Matth?i S K．Tessellations stable under iteration．International Discrete Fracture Network Engineering Conference，2014.

[15] Baecher G B．Statistical analysis of rock mass fracturing．Journal of the International Association for Mathematical Geology，1983，15(2)：329-348.

[16] Dowd P A，Xu C，Mardia K V et al．A comparison of methods for the stochastic simulation of rock fractures．Mathematical Geology，2007，39(7)：697-714.

[17] Xu C，Dowd P．A new computer code for discrete fracture network modelling．Computers & Geosciences，2010，36(3)：292-301.

[18] Nakaya S，Nakamura K．Percolation conditions in fractured hard rocks：A numerical approach using the three-dimensional binary fractal fracture network (3D-BFFN) model．Journal of Geophysical Research，2007，112(B12)：1-15.

[19] Alghalandis F Y．ADFNE：Open source software for discrete fracture network engineering，two and three dimensional applications．Computers & Geosciences，2017，102：1-11.

[20] Dong S，Zeng L，Dowd P et al．A fast method for fracture intersection detection in discrete fracture networks．Computers and Geotechnics，2018，98：205-216.

[21] Mardia K V，Nyirongo V B，Walder A N et al．Markov chain Monte Carlo implementation of rock fracture modelling．Mathematical Geology，2007，39(4)：355-381.

[22] Nakaya S，Yoshida T，Shioiri N．Percolation conditions in binary fractal fracture networks：Applications to rock fractures and active and seismogenic faults．Journal of Geophysical Research：Solid Earth，2003，108(B7)：1-13.

[23] Black J H．Hydrogeology of fractured rocks-a question of uncertainty about geometry．Applied Hydrogeology，1994，2(3)：56-70.

[24] Geier J，Lee K，Dershowitz W S．Field validation of conceptual models for fracture geometry．Transactions of American Geophysical Union，1988，69(44)：1177.

[25] Chiles J P．Fractal and geostatistical methods for modeling of a fracture network．Mathematical Geology，1988，20(6)：631-654.

[26] Billaux D，Fuller P．An algorithm for mesh simplification applied to fracture hydrology．Mathematical Geology，1989，21(2)：221-232.

[27] Kulatilake P H S W，Um J，Wang M et al．Stochastic fracture geometry modeling in 3-D including validations for a part of Arrowhead East Tunnel，California，USA．Engineering Geology，2003，70(1)：131-155.

[28] Zadhesh J，Jalali S-M E，Ramezanzadeh A．Estimation of joint trace length probability distribution function in igneous， sedimentary，and metamorphic rocks．Arabian Journal of Geosciences，2014，7(6)：2353-2361.

[29] Noroozi M，Kakaie R，Jalali S E．3D geometrical-stochastical modeling of rock mass joint networks：case study of the right bank of rudbar Lorestan dam plant．Journal of Geology and Mining Research，2015，7(1)：1-10.

[30] Watanabe N，Ishibashi T，Hirano N et al．Precise 3D numerical modeling of fracture flow coupled with X-ray computed tomography for reservoir core samples．SPE Journal，2011，16(3)：683-691.

[31] Gong J，Rossen W R．Modeling flow in naturally fractured reservoirs：Effect of fracture aperture distribution on dominant sub-network for flow．Petroleum Science，2017，14(1)：138-154.

[32] Ivanova V M．Geologic and Stochastic Modeling of Fracture Systems in Rocks [D]．Massachusetts Institute of Technology，1998.

[33] 鄭松青，姚志良．離散裂縫網絡隨機建模方法．石油天然氣學報，2009，31(4)：106-110．

Zheng Qingsong，Yao Zhiliang．Discrete fracture network stochastic modeling．Journal of Oil and Gas Technology，2009，31(4)：106-110.

[34] 于青春，大西有三．巖體三維不連續裂隙網絡及其逆建模方法．地球科學，2003，28(5)：522-526．

Yu Qingchun，Ohnishi Yuzo．Three-dimensional discrete fracture network model and its inverse method．Earth Science—Journal of China University of Geosciences，2003，28(5)：522-526.

[35] 宋曉晨，徐衛亞．裂隙巖體滲流模擬的三維離散裂隙網絡數值模型(Ⅰ)：裂隙網絡的隨機生成．巖石力學與工程學報，2004，23(12)：2015-2020．

Song Xiaochen，Xu Weiya．Numerical model of three-dimensional discrete fracture network for seepage in fractured rocks (Ⅰ)：generation of fracture network．Chinese Journal of Rock Mechanics and Engineering，2004，23(12)：2015-2020.

[36] Wettstein S J，Wittel F K，Araújo N A M et al．From invasion percolation to flow in rock fracture networks．Physica A：Statistical Mechanics and its Applications，2012，391(1)：264-277.

[37] 高志勇，朱如凱，馮佳睿等．中國前陸盆地構造—沉積填充響應與深層儲層特征．北京：地質出版社，2016．

Gao Zhiyong，Zhu Rukai，Feng Jiarui et al．Tectonic-sedimentary Filling Response and Deep Reservoir Characteristics of Foreland Basins in China．Geological Publishing House，Beijing，2016.

[38] Surrette M J，Allen D M．Quantifying heterogeneity in variably fractured sedimentary rock using a hydrostructural domain．Geological Society of America Bulletin，2008，120(1-2)：225-237.

[39] Lyu W，Zeng L，Liao Z et al．Fault damage zone cha-racterization in tight-oil sandstones of the Upper Triassic Yanchang Formation in the southwest Ordos Basin， China： Integrating cores， image logs， and conventional logs．Interpretation，2017，5(4)：SP27-SP39.

[40] Lyu W，Zeng L，Liu Z et al．Fracture responses of conventional logs in tight-oil sandstones：A case study of the Upper Triassic Yanchang Formation in southwest Ordos Basin，China．AAPG Bulletin，2016，100(9)：1399-1417.

[41] Tóth T M．Determination of geometric parameters of fracture networks using 1D data．Journal of Structural Geology，2010，32(7)：878-885.

[42] Davy P，Sornette A，Sornette D．Some consequences of a proposed fractal nature of continental faulting．Nature，1990，348(6296)：56-58.

[43] Bour O，Davy P．Connectivity of random fault networks following a power law fault length distribution．Water Resources Research，1997，33(7)：1567-1583.

[44] Bour O，Davy P．Clustering and size distributions of fault patterns：Theory and measurements．Geophysical Research Letters，1999，26(3)：2001-2004.

[45] Magnusdottir L，Horne R N．Inversion of time-lapse electric potential data to estimate fracture connectivity in geothermal reservoirs．Mathematical Geosciences，2015，47(1)：85-104.

[46] 吳勝和．儲層表征與建模．北京：石油工業出版社，2010．

Wu Shenghe．Reservoir Description and Modeling．Petroleum Industry Press，Beijing，2010.

[47] Xu C，Dowd P A，Mardia K V et al．A flexible true plurigaussian code for spatial facies simulations．Computers & Geosciences，2006，32(10)：1629-1645.

[48] Koike K，Kubo T，Liu C et al．3D geostatistical modeling of fracture system in a granitic massif to characterize hydraulic properties and fracture distribution．Tectonophysics，2015，660：1-16.

[49] Koike K，Liu C，Sanga T．Incorporation of fracture directions into 3D geostatistical methods for a rock fracture system．Environmental Earth Sciences，2012，66(5)：1403-1414.

[50] Sebacher B，Hanea R，Stordal A S．An adaptive pluri-Gaussian simulation model for geological uncertainty quantification．Journal of Petroleum Science and Engineering，2017，158：494-508.

[51] Koike K，Sakamoto H，Ohmi M．Detection and hydrologic modeling of aquifers in unconsolidated alluvial plains through combination of borehole data sets：A case study of the Arao area，Southwest Japan．Engineering Geology，2001，62(4)：301-317.

[52] Masoud A，Koike K．Tectonic architecture through Landsat-7 ETM+/SRTM DEM-derived lineaments and relationship to the hydrogeologic setting in Siwa region，NW Egypt．Journal of African Earth Sciences，2006，45(4-5)：467-477.

[53] Deutsch C V，Journel A G．Gslb：Geostatistical Software Library and User’s Guide．Oxford University Press，1992.

[54] Guardiano F B，Srivastava R M．Multivariate Geostatistics：Beyond Bivariate Moment．Springer，1993.

[55] Strebelle S．Conditional simulation of complex geolo-gical structures using multiple-point statistics．Mathematical Geology，2002，34(1)：1-21.

[56] 張文彪，段太忠，劉彥鋒等．綜合沉積正演與多點地質統計模擬碳酸鹽巖臺地——以巴西Jupiter油田為例．石油學報，2017，38(8)：925-934．

Zhang Wenbiao，Duan Taizhong，Liu Yanfeng et al．Integrated sedimentary forward modeling and multipoint geostatistics in carbonate platform simulation：A case study of Jupiter oilfield in Brazil．Acta Petrolei Sinica，2017，38(8)：925-934.

[57] Arpat G B．Sequential Simulation with Patterns [D]．Stanford University，2005.

[58] Arpat G B，Caers J．Conditional simulation with patterns．Mathematical Geology，2007，39(2)：177-203.

[59] Yin Y，Qin Z，Zuo R．Simulation of reservoir geometry using simpat method．Acta Geologica Sinica，2013，87(Z1)：607-609.

[60] Zhang T，Switzer P，Journel A．Filter-based classification of training image patterns for spatial simulation．Mathematical Geology，2006，38(1)：63-80.

[61] Honarkhah M，Caers J．Stochastic simulation of pat-terns using distance-based pattern modeling．Mathematical Geosciences，2010，42(5)：487-517.

[62] 喻思羽，李少華，何幼斌等．基于樣式降維聚類的多點地質統計建模算法．石油學報，2016，37(11)：1403-1409．

Yu Siyu，Li Shaohua，He Youbin et al．Multiple-point geostatistics algorithm based on pattern scale-down cluster．Acta Petrolei Sinica，2016，37(11)：1403-1409.

[63] 謝青．基于GPU的頁巖儲層裂縫建模及壓力模擬[學位論文]．安徽合肥：中國科學技術大學，2014．

Xie Qing．GPU based Fracture Modeling and Pressure Distribution Simulation of Shale Reservoirs [D]．University of Science and Technology of China，Hefei，Anhui，2014.

[64] Mohammadmoradi P．Facies and fracture network modeling by a novel image processing based method．Geomaterials，2013，3(4)：156-164.

[65] 張烈輝，賈鳴，張芮菡等．裂縫性油藏離散裂縫網絡模型與數值模擬．西南石油大學學報(自然科學版)，2017，39(3)：121-127．

Zhang Liehui，Jia Ming，Zhang Ruihan et al．Discrete fracture network modeling and numerical simulation of fractured reservoirs．Journal of Southwest Petroleum University(Science & Technology Edition)，2017，39(3)：121-127.

[66] Jia M，Zhang L，Guo J．Combining a connected-component labeling algorithm with FILTERSIM to simulate continuous discrete fracture networks．Environmental Earth Sciences，2017，76(8)：327.

[67] Liu X，Zhang C，Liu Q et al．Multiple-point statistical prediction on fracture networks at Yucca Mountain．Environmental Geology，2009，57(6)：1361-1370.

[68] Chugunova T L，Hu L Y．Multiple-point simulations constrained by continuous auxiliary data．Mathematical Geosciences， 2008，40(2)：133-146.

[69] Chugunova T，Corpel V，Gomez J P．Explicit fracture network modelling：From multiple point statistics to dynamic simulation．Mathematical Geosciences，2017，49(4)：541-553.

[70] Qiao J．Julia sets and complex singularities of free energies．Memoirs of the American Mathematical Society，2015，5：1-3.

[71] Qiao J．Julia sets and complex singularities in diamond-like hierarchical Potts models．Science China Mathematics，2005，48(3)：388-412.

[72] Acuna J A，Yortsos Y C．Application of fractal geome-try to the study of networks of fractures and their pressure transient．Water Resources Research，1995，31(3)：527-540.

[73] 袁杰，李霞麗，白水．迭代函數系統及其應用．中央民族大學學報(自然科學版)，2004，13(4)：317-321．

Yuan Jie，Li Xiali，Bai Shui．Iterated function system and its application．Journal of The Central University For Nationalities (Natural Science Edition)，2004，13(4)：317-321.

[74] 周德華，焦方正，葛家理．一種油氣藏裂縫網絡模擬優化生成技術．斷塊油氣田，2002，9(4)：30-32．

Zhou Dehua，Jiao Fangzheng，Ge Jiali．A new approach to optimistically model fracture networks of hydrocarbon reservoirs．Fault-block Oil & Gas Field，2002，9(4)：30-32.

[75] Acuna J A，Yortsos Y C．Numerical construction and flow simulation in networks of fractures using fractal geometry．Society of Petroleum Engineers，1991:1-7.

[76] Zhou Z，Su Y，Wang W et al．Application of the fractal geometry theory on fracture network simulation．Journal of Petroleum Exploration & Production Technology，2017，7(2)：487-496.

[77] Zhou Z，Su Y，Wang W et al．Integration of microseismic and well production data for fracture network calibration with an L-system and rate transient analysis．Journal of Unconventional Oil & Gas Resources，2016，15：113-121.

[78] Fan D，Ettehadtavakkol A．Semi-analytical modeling of shale gas flow through fractal induced fracture networks with microseismic data．Fuel，2017，193：444-459.

[79] Dasho O A，Ariyibi E A，Akinluyi F O et al． Application of satellite remote sensing to groundwater potential modeling in Ejigbo area，Southwestern Nigeria．Modeling Earth Systems and Environment，2017，3(2)：615-633.

[80] Cao T，Xiao A，Wu L et al．Automatic fracture detection based on Terrestrial Laser Scanning data：A new method and case study．Computers & Geosciences，2017，106：209-216.

[81] 曾慶魯，張榮虎，盧文忠等．基于三維激光掃描技術的裂縫發育規律和控制因素研究——以塔里木盆地庫車前陸區索罕村露頭剖面為例．天然氣地球科學，2017，28(3)：397-409．

Zeng Qinglu，Zhang Ronghu，Lu Wenzhong et al．Fracture development characteristics and controlling factors based on 3D laser scanning technology：An outcrop case study of Suohan village，Kuqa foreland area，Tarim Basin．Natural Gas Geoscience，2017，28(3)：397-409.

[82] Vazaios I，Vlachopoulos N，Diederichs M S．Integra-tion of lidar-based structural input and discrete fracture network generation for underground applications．Geotechnical and Geological Engineering，2017，35(5)：2227-2251.

[83] Casini G，Hunt D W，Monsen E et al．Fracture characterization and modeling from virtual outcrops．AAPG Bulletin，2016，100(1)：41-61.

[84] Wang X，Zou L，Shen X et al．A region-growing approach for automatic outcrop fracture extraction from a three-dimensional point cloud．Computers & Geosciences，2017，99：100-106.

[85] Santos R F V C，Miranda T S，Barbosa J A et al．Chara-cterization of natural fracture systems：Analysis of uncertainty effects in linear scanline results．AAPG Bulletin，2015，99(12)：2203-2219.

[86] Bisdom K，Nick H M，Bertotti G．An integrated workflow for stress and flow modelling using outcrop-derived discrete fracture networks．Computers & Geosciences，2017，103：21-35.

[87] 肖小玲，靳秀菊，張翔等．基于常規測井與電成像測井多信息融合的裂縫識別．石油地球物理勘探，2015，50(3)：542-547．

Xiao Xiaoling，Jin Xiuju，Zhang Xiang et al．Fracture identification based on information fusion of conventional logging and electrical imaging logging．OGP，2015，50(3)：542-547.

[88] Lai J，Wang G，Fan Z et al．Three-dimensional quantitative fracture analysis of tight gas sandstones using industrial computed tomography．Scientific Reports，2017，7(1)：1825.

[89] Karimpouli S，Tahmasebi P，Ramandi H L et al．Stochastic modeling of coal fracture network by direct use of micro-computed tomography images．International Journal of Coal Geology，2017，179：153-163.

[90] Huang L，Tang H，Tan Q et al．A novel method for correcting scanline-observational bias of discontinuity orientation．Scientific Reports，2016，6(1)：22942.

[91] Sanai L，Chenini I，Ben M A et al．Fracture related-fold patterns analysis and hydrogeological implications：Insight from fault-propagation fold in Northwestern of Tunisia．Journal of African Earth Sciences，2015，101：375-382.

[92] 祖克威，曾聯波，鞏磊．斷層相關褶皺概念模型中的裂縫域．地質科學，2013，48(4)：1140-1147．

Zu Kewei,Zeng Lianbo,Gong Lei.Fractures domains in conceptual models of fault-related folds．Chinese Journal of Geology，2013，48(4)：1140-1147.

[93] 郎曉玲，郭召杰．基于DFN離散裂縫網絡模型的裂縫性儲層建模方法．北京大學學報(自然科學版)，2013，？(6)：964-972．

Lang Xiaoling，Guo Zhaojie．Fractured reservoir modeling method based on discrete fracture network model．Acta Scientiarum Naturalium Universitatis Pekinensis，2013，？(6)：964-972.

[94] 張雨晴，王暉，范廷恩等．花崗巖潛山儲層裂縫建模表征方法——以渤海花崗巖潛山A油田為例．吉林大學學報(地)，2016，46(5)：1312-1320．

Zhang Yuqing，Wang Hui，Fan Ting’en et al．Granite buried hill reservoir characterization and modeling：Taking offshore an oilfield in Bohai Bay as an example．Journal of Jilin University (Earth Science Edition)，2016,46(5)：1312-1320.

[95] 熊曉軍，簡世凱，李翔等．基于標準差統計的窄方位角疊前裂縫預測方法．石油地球物理勘探，2017，52(1)：114-120．

Xiong Xiaojun，Jian Shikai，Li Xiang et al．Fracture prediction on prestack narrow azimuth data with standard deviation statistic.OGP,2017,52(1):114-120.

[96] 李慧瓊，張盟勃，蒲仁海等．黃257井區疊前縱波方位各向異性裂縫分布預測．石油地球物理勘探，2017，52(2)：350-359．

Li Huiqiong，Zhang Mengbo，Pu Renhai et al．Late Triassic fracture detection with seismic azimuth anisotropics in Huang257 survey，OrdosBasin．OGP，2017，52(2)：350-359.

[97] Xue J，Gu H，Cai C．Model-based amplitude versus offset and azimuth inversion for estimating fracture parameters and fluid content．Geophysics，2017，82(2)：M1-M17.

[98] 薛姣，顧漢明，蔡成國．準靜態孔縫介質廣義裂隙弱度研究．石油地球物理勘探，2015，50(6)：1146-1153．

Xue Jiao，Gu Hanming，Cai Chengguo．General fracture weaknesses for quasi-static porous fractured media．OGP，2015，50(6)：1146-1153.

[99] 袁曉宇，李映濤，葉寧等．基于頻譜分解的螞蟻追蹤裂縫檢測技術在玉北地區的應用研究．石油地球物理勘探，2015，50(4)：665-671．

Yuan Xiaoyu，Li Yingtao，Ye Ning et al．The application of ant-tracking fracture detection based on spectrum decomposition in Yubei area．OGP，2015，50(4)：665-671.

[100] Wilson T H，Smith V，Brown A．Developing a model discrete fracture network，drilling，and enhanced oil recovery strategy in an unconventional naturally fractured reservoir using integrated field，image log，and three-dimensional seismic data．AAPG Bulletin，2015，99(4)：735-762.

[101] 汪忠德，薛詩桂，李紅敬等．油氣地球物理解釋技術研究新進展．地球物理學進展，2016，31(3)：1187-1201．

Wang Zhongde，Xue Shigui，Li Hongjing et al．New progress in the technology of oil and gas geophysical interpretation．Progress in Geophysics，2016，31(3)：1187-1201.

[102] 劉建偉，張云銀，曾聯波等．非常規油藏地應力和應力甜點地球物理預測——渤南地區沙三下亞段頁巖油藏勘探實例．石油地球物理勘探，2016，51(4)：792-800．

Liu Jianwei，Zhang Yunyin，Zeng Lianbo et al．Geophysical prediction of stress and stress desserts in unconventional reservoirs：An example in Bonan area．OGP，2016，51(4)：792-800.

[103] 薛艷梅，夏東領，蘇宗富等．多信息融合分級裂縫建模．西南石油大學學報(自然科學版)，2014，36(2)：57-63．

Xue Yanmei，Xia Dongling，Su Zongfu et al．Fracture modeling at different scales based on convergent multi-source information．Journal of Southwest Petroleum University (Science & Technology Edition)，2014，36(2)：57-63.

[104] 張成，夏露，李英杰等．柴達木盆地東部石炭系頁巖三維裂縫網絡模型研究．地學前緣，2016，23(5)：184-192．

Zhang Cheng，Xia Lu，Li Yingjie et al．Study on three-dimensional fracture network models of Carboniferous shale in Eastern Qaidam Basin．Earth Science Frontiers，2016，23(5)：184-192.

[105] Xu C，Dowd P A，Mardia K V et al．A connectivity index for discrete fracture networks．Mathematical Geology，2007，38(5)：611-634.

[106] Alghalandis Y F，Dowd P A，Xu C．Connectivity field：a measure for characterising fracture networks．Mathematical Geosciences，2015，47(1)：63-83．