銅官山礦田深部構(gòu)造的三維計算模擬及其找礦意義

劉亮明, 曹　偉(中南大學(xué) 計算地球科學(xué)研究中心, 教育部 有色金屬成礦預(yù)測重點(diǎn)實(shí)驗室, 湖南 長沙 410083)

銅官山礦田深部構(gòu)造的三維計算模擬及其找礦意義

劉亮明, 曹偉
(中南大學(xué) 計算地球科學(xué)研究中心, 教育部 有色金屬成礦預(yù)測重點(diǎn)實(shí)驗室, 湖南 長沙 410083)

計算模擬已成為研究礦田深部三維構(gòu)造重要可行手段, 對深部找礦具有重要意義。基于銅官山礦田地質(zhì)背景及所獲勘查數(shù)據(jù)特征, 我們采用知識驅(qū)動下多級約束和多源數(shù)據(jù)融合建模方法模擬礦田范圍內(nèi)–2000 m 以上主要地質(zhì)界面。以逆沖推覆導(dǎo)致斷展–斷彎褶皺理論模式解釋大地測深的電阻率數(shù)據(jù), 推斷地質(zhì)界面位置, 在 Micromine- GOCAD耦合平臺上在多級約束下融合多源的界面位置數(shù)據(jù), 以離散光滑插值算法和Delaunay三角網(wǎng)剖分法重構(gòu)主要地質(zhì)界面, 在計算機(jī)虛擬空間再現(xiàn)礦田三維構(gòu)造格架。結(jié)果顯示: 銅官山巖體和天鵝抱蛋巖體三維形態(tài)變化復(fù)雜, 總體上向西傾, 向深部巖體規(guī)模快速變小。位于銅官山巖體東接觸帶的銅官山銅礦田深部無找礦潛力, 原因是接觸帶產(chǎn)狀變陡, 并且是巖體直接與上泥盆統(tǒng)砂巖接觸; 礦田深部主要構(gòu)造為走向NE并沿傾向彎曲的逆沖推覆斷層; 銅官山背斜東翼含礦的中–上石炭統(tǒng)層位向深部倒轉(zhuǎn)后又在–1200 m至–1500 m深度返轉(zhuǎn)成向東南緩傾的正常層位, 其底部與逆沖推覆斷層重合, 是深部找礦的有利靶區(qū)。

深部構(gòu)造; 三維模型; 計算模擬; 銅官山礦田

0　引　言

隨著淺部礦產(chǎn)資源日益減少, 深部找礦不僅成為找礦勘查重要增儲方向, 也成為最具挑戰(zhàn)性方向(翟裕生等, 2004; 劉亮明等, 2005; Williams, 2005; 劉亮明和蔡愛良, 2009)。任何礦床空間定位都與構(gòu)造密切相關(guān)(陳國達(dá), 1985)。深部構(gòu)造的深度隱伏使得深部構(gòu)造的探測和分析存在相當(dāng)大難度, 深部構(gòu)造研究實(shí)際也就成了深部找礦所面對主要挑戰(zhàn)之一。

在大多數(shù)礦田中, 深部構(gòu)造組成要素在三維空間內(nèi)復(fù)雜多變, 傳統(tǒng)地質(zhì)方法顯然難以充分揭示勘查區(qū)域內(nèi)地質(zhì)要素的特征及其時空規(guī)律, 這已成為制約成礦理論預(yù)測能力主要原因(劉亮明, 2007)。要依據(jù)并不充分信息, 在三維空間內(nèi)盡可能準(zhǔn)確表征深部構(gòu)造和了解深部地質(zhì)過程, 計算模擬已成為一種必不可少的有效方法(趙義來和劉亮明, 2011; Liu et al., 2012)。隨著計算機(jī)技術(shù)、計算圖形學(xué)和信息科學(xué)發(fā)展, 復(fù)雜地質(zhì)體形態(tài)三維計算模擬在圖形生成算法、建模方法、不確定性分析、信息融合以及勘查應(yīng)用等方面都取得令人矚目的發(fā)展(Lorensen and Cline, 1987; Mallet, 1989, 2002; Victor, 1993; Houlding, 1994; Gong et al., 2004; Xue et al., 2004; Wu et al., 2005; Turner, 2006; Frank et al., 2007; Bistacchi et al., 2008; Calcagno et al., 2008; Kaufmann and Martin, 2008; Caumon et al., 2009; Zanchi et al., 2009; Ming et al., 2010; Wellmann et al., 2010; Wellmann and Regenauer-Lieb, 2012; Bouza-Rodríguez et al., 2014; Elsheikh A H and Elsheikh M, 2014; Saalmann and Laine, 2014), 這些成果為三維計算模擬復(fù)雜的深部構(gòu)造創(chuàng)造了有利條件。

銅官山礦田是安徽銅陵地區(qū)一個主要銅金多金屬礦田, 礦田內(nèi)主要礦床有銅官山銅礦、金口嶺銅礦、天馬山硫金礦和黃獅澇金礦。雖然這個礦田已經(jīng)歷長期開發(fā)(最早可追溯到唐代), 地表和淺部資源已消耗殆盡, 但近年來的深部找礦勘查工作顯示深部可能還存在相當(dāng)可觀的隱伏礦體。由于銅官山礦田礦床類型主要是與巖體相關(guān)的矽卡巖型和層控矽卡巖型, 長期以來, 地質(zhì)研究重點(diǎn)主要是礦床成因(郭文魁, 1957)、成礦巖體(周泰禧等, 1987; 楊學(xué)明和林文通, 1988; 杜楊松等, 2004, 2007; 吳才來等, 2010; 趙乘乘等, 2012)和礦床地球化學(xué)(田世洪等, 2005; 李紅陽等, 2008; 楊秋榮等, 2010)。對于礦田構(gòu)造研究并不深入, 主要集中于地表和淺部揭露構(gòu)造的幾何學(xué)及其與礦體空間關(guān)系方面(薛虎等, 1983;雷世和, 1985; 萬仁虎, 2004)。對深部構(gòu)造了解甚少,更沒有涉及到深部三維構(gòu)造格架, 這一局面無疑會嚴(yán)重制約深部礦體預(yù)測和勘查。本文綜合集成有關(guān)銅官山礦田深部地質(zhì)和地球物理信息, 以計算模擬方法再現(xiàn)礦田的深部三維構(gòu)造格架, 并對深部找礦有利靶區(qū)進(jìn)行預(yù)測分析。

1　地質(zhì)背景

銅官山礦田位于銅陵礦集區(qū)最西部, 西鄰長江。礦田最基本地質(zhì)特征是在晚古生代以碳酸鹽巖為主的沉積巖中侵入以閃長巖為主要成分的淺成侵入巖體, 區(qū)內(nèi)所有礦床都與這些長英質(zhì)淺成巖體密切相關(guān)。地表出露巖體主要有3個: ①銅官山巖體;②天鵝抱蛋巖體; ③金口嶺巖體(圖1)。銅官山巖體地表呈橢圓形, 出露面積為 1.5 km2, 巖性主要為石英閃長巖。鋯石 SHRIMP U-Pb同位素年齡為141.8±1 Ma(吳才來等, 2010), 為燕山早期。天鵝抱蛋巖體地表呈不規(guī)則圓形, 面積約 0.8 km2, 巖性主要為閃長巖, 鋯石 SHRIMP U-Pb同位素年齡為141.3±1.3 Ma(吳才來等, 2010)。金口嶺巖體地表呈一長軸近南北不規(guī)則橢圓形, 出露面積約 5 km2。巖體礦物成分和化學(xué)成分與銅官山巖體相似, 主要為石英閃長巖, 其黑云母 Ar-Ar同位素年齡為 137 Ma(周泰禧等, 1987)。考慮到銅官山巖體也獲得136 Ma的角閃石 Ar-Ar同位素年齡(周泰禧等, 1987),明顯低于鋯石 U-Pb同位素年齡。金口嶺巖體與銅官山巖體、天鵝抱蛋巖體應(yīng)同樣為燕山早期巖漿活動產(chǎn)物。

圍繞3個巖體出露地層為志留系?三疊系(圖1)。主要地層單元有: ①上志留統(tǒng)茅山組(S3m), 為砂巖、粉砂巖夾砂頁巖; ②上泥盆統(tǒng)五通組(D3w), 主要為砂巖, 有少量礫巖; ③中–上石炭統(tǒng)(C2+3), 主要為白云巖、白云質(zhì)灰?guī)r和灰?guī)r, 已變質(zhì)成大理巖;④下二疊統(tǒng)棲霞組(P1q), 主要為生物碎屑灰?guī)r, 已變質(zhì)成大理巖; ⑤下二疊統(tǒng)孤峰組(P1g), 主要為硅質(zhì)巖和硅質(zhì)灰?guī)r; ⑥上二疊統(tǒng)龍?zhí)督M(P2l), 為炭質(zhì)頁巖和粉砂頁巖, 夾煤層; ⑦上二疊統(tǒng)大隆組(P2d),為硅質(zhì)頁巖和硅質(zhì)灰?guī)r; ⑧下三疊統(tǒng)殷坑組(T1y),主要為夾粉砂質(zhì)頁巖的灰?guī)r; ⑨下三疊統(tǒng)和龍山組(T1h), 為含泥質(zhì)條帶灰?guī)r; ⑩下三疊統(tǒng)南陵湖組(T1n), 為中厚層灰?guī)r。

銅官山礦田地表主要構(gòu)造為軸向 NE銅官山背斜, 其次為走向NE和NW斷層(圖1)。銅官山背斜為向北東傾伏的不對稱倒轉(zhuǎn)背斜, 背斜軸線整體方向為42°左右。以銅官山制高點(diǎn)為分界, 整個背斜構(gòu)造分為兩段, 北東段背斜轉(zhuǎn)折端向北東方向傾伏,傾伏角近 10°, 形成較長傾伏端, 地表為第四系覆蓋,下伏地層從五通組到殷坑組都有。其北西翼地層向北西方向傾, 傾角為30°～50°, 南東翼地層在淺部傾向東南, 傾角為 25°左右, 向深變陡, 并逐漸轉(zhuǎn)為傾向北西。背斜西南段, 轉(zhuǎn)折端向上拱, 呈穹灣狀, 其核部地層為志留系砂頁巖, 兩翼地層從泥盆系到三疊系都有, 其北西翼地層向北西傾, 傾角為 40°～60°,其東南翼地層也向北西傾, 傾角大于60°。

銅官山礦田范圍內(nèi)已發(fā)現(xiàn)礦體主要為4種類型:①接觸帶矽卡巖型銅(鐵、硫)礦體, 產(chǎn)于巖體與碳酸鹽巖的接觸帶上, 成礦巖體主要為銅官山巖體, 其次為金口嶺巖體; ②層控矽卡巖型(鐵、硫)礦體, 礦體從巖體的接觸帶順巖層延伸; ③層控硫金礦體,含礦層位主要為中–上石炭統(tǒng)碳酸鹽巖; ④脈狀銅(鉬、金)礦體, 規(guī)模小, 主要產(chǎn)于靠近接觸帶的巖體和碳酸鹽巖中。

礦田內(nèi)已發(fā)現(xiàn)礦體空間分布及其與其他地質(zhì)要素之間關(guān)系顯示: 控制礦體空間位置與規(guī)模的兩個主要地質(zhì)因素是成礦巖體和中–上石炭統(tǒng)碳酸鹽巖,所以通過深部構(gòu)造研究揭示巖體在深部形態(tài)變化(包括是否存在未發(fā)現(xiàn)的隱伏巖體)以及中–上石炭統(tǒng)在深部的展布對預(yù)測和勘查深部礦體具有十分重要的意義。

圖1　銅官山礦田地質(zhì)簡圖(據(jù)銅陵有色金屬集團(tuán)股份有限公司礦產(chǎn)資源中心, 2011修改)Fig.1　Geological sketch map of the Tongguanshan ore field

2　地質(zhì)構(gòu)造三維形態(tài)的計算模擬

2.1三維形態(tài)計算模擬方法

依據(jù)低維地質(zhì)數(shù)據(jù)構(gòu)建三維地質(zhì)體模型計算模擬方法很多, 概括起來可分為 3類: ①面元模擬(surface-based modeling); ②體元模擬(volume-based modeling)和③混合模擬(hybrids modeling)(Shi, 1996; Gong et al., 2004; Bouza-Rodríguez et al., 2014)。面元模擬是根據(jù)地質(zhì)體邊界點(diǎn)和線構(gòu)建地質(zhì)體邊界面,依生成面的方法不同模擬方法不同, 主要有: 不規(guī)則三角網(wǎng)(triangulated irregular network, 簡稱 TIN)模擬(De Floriani, 1989)、邊界再現(xiàn)(boundary representation, 簡稱B-Rep)模擬(Stroud, 2006)、非均勻有理 B樣條曲線(non-uniform rational B-splines,簡稱NURBS)模擬(Fisher and Wales, 1991)等。體元模擬是用各種形狀體元來擬合地質(zhì)體, 主要有體素(volumetric pixel, 簡稱 Voxel)模擬(Marschallinger, 1996)、八叉樹(octree)模擬(Dunstan and Mill, 1989)、四面體格網(wǎng)(tetrahedron network, 簡稱 TEN)模擬(Pilout et al., 1994)、三棱柱(tri-prism, 簡稱TP)模擬(Houlding, 1994)等。混合模擬則將兩種方法混合使用, 如Octree-TEN(李德仁和李清泉, 1997)、TIN-Octree (Shi, 1996)等。

由于地質(zhì)要素中有不少是二維延展, 如地層界面、斷層、不整合面等, 面元模擬在地質(zhì)體三維模擬中有著比體元模擬更加明顯優(yōu)勢(Houlding, 1994; Caumon et al., 2009)。TIN已成為構(gòu)建三維地質(zhì)模型最為廣泛采用方法, 其基本原理是將空間上無重復(fù)的散亂數(shù)據(jù)點(diǎn)集按剖分規(guī)則進(jìn)行三角剖分, 使這些散亂點(diǎn)形成連續(xù)但不重疊的不規(guī)則三角面片網(wǎng), 并以此來描述三維物體的表面。TIN模擬中最廣泛使用的剖分規(guī)則是Delaunay三角網(wǎng)規(guī)則, 這最早是由俄羅斯數(shù)學(xué)家 Delaunay(1934)提出并最終以他的名字命 名 的 剖 分 規(guī) 則 (Victor, 1993; George and Borouchaki, 1998), 它可以使剖分出的TIN最小內(nèi)角最大化, 這就能最大限度地避免“瘦長”三角形, 從而自動向等邊三角形靠近, 因此, Delaunay剖分被認(rèn)為是最優(yōu)的三角網(wǎng)剖分法之一。 Delaunay三角網(wǎng)最重要性質(zhì)是空圓特性, 即三角網(wǎng)中任意一個三角形外接圓范圍內(nèi)必然不包含其他三角形頂點(diǎn), 這條特性成為構(gòu)建Delaunay三角網(wǎng)的首要準(zhǔn)則和自動剖分算法的基礎(chǔ)。基于Delaunay法則創(chuàng)造的一些算法已成三維形態(tài)模擬的精典算法, 如 Lawson算法、Bowyer-Watson算法(Rice and Usowt, 1968; Rebay, 1993; George and Borouchaki, 1998)。

由于空間數(shù)據(jù)不充分和不均勻分布, 直接由原始數(shù)據(jù)點(diǎn)經(jīng)Delaunay三角剖分生成三維曲面偏離真實(shí)且極不美觀(如圖2a), 還必須進(jìn)一步優(yōu)化。優(yōu)化途徑有兩條, 一是由一些插值的算法加密數(shù)據(jù)點(diǎn); 二是人工插入一些關(guān)鍵控制點(diǎn)。通過插值加密數(shù)據(jù)點(diǎn)算法中以離散光插值(discrete smooth interpolation,簡稱DSI)算法應(yīng)用最廣, 同樣數(shù)據(jù)經(jīng)過DSI處理后的效果顯然好很多(如圖2b)。

DSI算法是Mallet (1989, 1992, 2002)發(fā)明的, 其基本內(nèi)容是: 對一個離散化自然體構(gòu)建其各個部分相互聯(lián)系的網(wǎng)格(曲面)模型, 根據(jù)約束條件已知某些點(diǎn)位置或?qū)傩灾? 則未知節(jié)點(diǎn)上的相應(yīng)值可以通過解線性方程得到。設(shè)網(wǎng)格所有節(jié)點(diǎn) k構(gòu)成的集合為?, φ(k)為定義在?上的函數(shù); 設(shè)對于?的子集L中節(jié)點(diǎn)l, φ(l)已知, 對于?的子集I中節(jié)點(diǎn)i, φ(i)未知; L+I=?。在?上, 對已知點(diǎn){φ(l): l∈L}可以構(gòu)建無窮多個插值函數(shù) φ(k), 而離散光滑插值的目的即是從這些插值函數(shù)中找到一個使給定準(zhǔn)則 R*(φ)有最小值的函數(shù)。R*(φ)的表達(dá)式為:

R*(φ)= R(φ)+ρ(φ)

其中, R(φ)為函數(shù) φ(k)的全局粗糙度, ρ(φ)表征函數(shù)φ(k)與它應(yīng)滿足的約束條件之間的偏差值(約束違反度), 約束條件即為從原始數(shù)據(jù)提取的各種點(diǎn)狀、線狀約束。可見DSI插值生成的曲面需要滿足擬合約束及光順兩個要求。即通過DSI插值能夠在滿足所有約束條件下使生成曲面盡可能地光順美觀,光順美觀曲面是否代表真實(shí)地質(zhì)界面則取決于約束條件多少。由于礦田深部有工程控制的約束條件稀少, 所以人工插入控制點(diǎn)就顯得尤為重要。人工插入控制點(diǎn)可靠性和效果取決于建模者對模擬對象真實(shí)形狀了解。

基于銅官山礦田實(shí)際情況, 我們對礦田深部構(gòu)造計算模擬采用基于Delaunay三角剖分的TIN面模擬法模擬礦田構(gòu)造中主要地質(zhì)界面, 對缺失數(shù)據(jù)一些關(guān)鍵點(diǎn)通過地質(zhì)構(gòu)造理論對電阻率數(shù)據(jù)的解釋來決定人工插值, 最終由 DSI算法進(jìn)行優(yōu)化處理。采用計算模擬軟件為 Micromine和 GOCAD。Micromine是由澳大利亞 Micromine公司開發(fā)專門進(jìn)行礦床勘探和開發(fā)設(shè)計的商業(yè)軟件(見 http: // www.micromine.com/), 在勘探數(shù)據(jù)錄入和管理方面以及圖形操作方面具獨(dú)特優(yōu)勢。GOCAD是由法國Nancy大學(xué)開發(fā)用于地質(zhì)對象三維建模的軟件(見http: //www.gocad.org/), 由于軟件采用DSI算法, 在生成面模型時具有獨(dú)特優(yōu)勢。所以我們在建模時,采用 Micromine進(jìn)行勘探數(shù)據(jù)錄入和處理, 而用GOCAD來生成最終的面模型。

圖2　DSI優(yōu)化效果對比Fig.2　Efficiency of DSI optimization

2.2模擬數(shù)據(jù)及其所存在不足

銅官山礦田長期勘探和開發(fā)積累了大量地質(zhì)勘探資料, 但由于幾個礦山都經(jīng)歷過關(guān)破、轉(zhuǎn)產(chǎn)或準(zhǔn)備轉(zhuǎn)產(chǎn)等, 不少資料已丟失。早期勘探鉆孔和坑道測量、編錄和分析數(shù)據(jù)基本都已丟失, 只能找到部分勘探線剖面圖。目前可利用數(shù)據(jù)包括如下4種類型: ①早期勘探剖面和中段平面圖, 不完整, 共計 77幅剖面圖、46幅中段及地表圖; ②近年來施工的地質(zhì)探礦鉆孔資料, 孔口坐標(biāo)、測斜數(shù)據(jù)、編錄數(shù)據(jù)和分析數(shù)據(jù)齊全, 共計447個鉆孔, 99812 m巖性編錄; ③包含部分老鉆孔的新近完成鉆孔勘探線剖面圖, 共計26幅; ④我們新近完成的大地電磁測深電阻率數(shù)據(jù), 共計441個測點(diǎn)的9475個電阻率數(shù)據(jù), 探測深度到–2000 m。

上述看似龐大數(shù)據(jù)其實(shí)對建立如此大面積的礦田深部構(gòu)造模型來說還是非常不充足的, 表現(xiàn)在:①由工程控制的地質(zhì)數(shù)據(jù)分布范圍太局限, 平面上限于礦體分布區(qū), 深度上大都在–500 m以上, 超過1000 m鉆孔只有5個, 超過500 m鉆孔也只有50個;②深部的地質(zhì)構(gòu)造主要靠大地電磁測深獲得的電阻率數(shù)據(jù)來推斷。電阻率數(shù)據(jù)至少存在3個方面問題:第一、由于地表建筑和工農(nóng)業(yè)設(shè)施, 很多地方無法進(jìn)行大地電磁測深; 第二、地表和地下電力設(shè)施的干擾使得獲得的電阻率數(shù)據(jù)不可避免地含有假異常;第三、依電阻率數(shù)據(jù)解釋深部地質(zhì)構(gòu)造本身毫無疑問具有很強(qiáng)多解性。

正是由于模擬數(shù)據(jù)存在上述不足, 可以肯定將上述數(shù)據(jù)輸入計算模擬軟件平臺不可能由計算機(jī)能自動生成一個合理的礦田構(gòu)造格架的三維模型。計算模擬遠(yuǎn)非僅僅依靠計算機(jī)就能勝任工作, 還需要地質(zhì)人員應(yīng)用地質(zhì)知識從理論上進(jìn)行推斷和選擇,以彌補(bǔ)模擬數(shù)據(jù)不足。

2.3知識驅(qū)動多級約束和多源數(shù)據(jù)融合模擬途徑

我們以面元模型模擬銅官山礦田的深部地質(zhì)構(gòu)造, 需要重構(gòu)礦田范圍內(nèi)–2000 m深度以上的主要地質(zhì)界面, 就是將工程控制的地質(zhì)邊界和推斷的地質(zhì)界面數(shù)據(jù)通過三維面模型形式在計算機(jī)虛擬空間中表達(dá)出來。由模擬流程來分析(圖3), 模擬的最終結(jié)果可靠與否主要取決于輸入的界面數(shù)據(jù)。由于工程控制范圍十分有限, 大部分空間地質(zhì)界面位置需要推斷。推斷的證據(jù)主要來自大地電磁測深所獲得的電阻率數(shù)據(jù), 所以由電阻率數(shù)據(jù)推斷的深部地質(zhì)界面對本區(qū)深部構(gòu)造三維模擬至關(guān)重要。

圖3　三維地質(zhì)構(gòu)造形態(tài)模擬流程圖Fig.3　Flow chart for the 3D geological structural modeling

對于大地電磁測深獲得電阻率數(shù)據(jù), 通過Kriging插值獲得了電阻率在三維空間分布的塊體模型(圖4)。這個模型非常形像地展示電阻率在三維空間分布, 總體上與銅官山背斜構(gòu)造形態(tài)是吻合的,背斜核部未蝕變的志留系砂巖表現(xiàn)為高阻, 背斜翼部未蝕變的晚古生代碳酸鹽巖也表現(xiàn)為高阻(圖 4),但低阻卻具有很強(qiáng)的多解性, 既可能是 P1q含炭質(zhì)灰?guī)r和P1l含煤頁巖, 也能是蝕變的碳酸鹽巖和砂巖或蝕變巖體, 還能是含水破碎帶, 甚至還可能是因為干擾而引起。所以, 僅僅依靠這個電阻率三維空間分布模型無法正確推斷無工程控制部位的地質(zhì)邊界。這就需要知識驅(qū)動, 即應(yīng)用適合本地地質(zhì)構(gòu)造理論(概念)模式, 充分考慮地表及地下工程揭露地質(zhì)體產(chǎn)狀和空間形態(tài)變化趨勢來解釋電阻率分布三維塊體模型, 推斷地質(zhì)界面的空間位置(圖3)。地質(zhì)構(gòu)造理論模式是推斷深部地質(zhì)構(gòu)造非常重要理論依據(jù)(劉亮明和蔡愛良, 2009), 而銅官山礦田需要推斷地質(zhì)界面空間范圍要遠(yuǎn)遠(yuǎn)大于已有工程控制范圍, 所以這種依構(gòu)造理論模式推斷深部地質(zhì)構(gòu)造的知識驅(qū)動作用對本區(qū)深部構(gòu)造的計算模擬具有十分重要的意義。

圖4　銅官山礦田電阻率三維分布塊體模型Fig.4　3D block model of resistivity in the Tongguanshan ore field

以知識驅(qū)動推斷地質(zhì)(體)邊界位置正確與否主要取決于所應(yīng)用地質(zhì)構(gòu)造概念模式是否真正反映研究區(qū)地質(zhì)構(gòu)造形成的動力機(jī)制(劉亮明和蔡愛良, 2009)。在銅官山礦田, 地質(zhì)構(gòu)造概念模式必須從區(qū)域地殼演化的動力學(xué)解釋形成銅官山背斜獨(dú)特特征的機(jī)制和原因。綜合分析銅官山背斜地質(zhì)特征以及區(qū)域地殼動力學(xué)背景, 我們認(rèn)為銅官山背斜是由逆沖推覆而派生的斷展–斷彎褶皺(圖 5)。地層在沖斷演化過程中, 推覆斷層的滑動完全被褶皺作用所吸收, 沒有滑動傳遞到褶皺以外, 形成斷展褶皺(fault-propagation fold)。斷層繼續(xù)滑動可以追蹤層間薄層帶形成彎曲的推覆斷層, 形成斷彎褶皺(fault bend fold)。Suppe (1983)、Suppe and Medwedeff (1990)和Mercier et al. (1997)等創(chuàng)立的這種斷展–斷彎褶皺理論模式適用于銅官山礦田, 其證據(jù)有如下 3個方面: ①銅官山背斜是東南翼局部倒轉(zhuǎn)的不對稱背斜,而且這種不對稱褶皺是銅陵–安慶地區(qū)褶皺的普遍特征, 這也是逆沖推覆形成褶皺的典型特征; ②銅陵–安慶地區(qū)褶皺的古生代地層中存在多個層間滑脫帶, 如: D3w與C2+3之間、P2l和T1h之間; ③逆沖推覆是整個銅陵–安慶地區(qū)印支期造山的主要構(gòu)造型式(董樹文和邱瑞龍, 1993)。基于上述考慮, 我們將應(yīng)用逆沖推覆引起斷展–斷彎褶皺理論模式來解釋銅官山礦田大地電磁測深所獲得的電阻率數(shù)據(jù),推斷礦田深部的主要地質(zhì)界面。

圖5　逆沖推覆斷層派生斷展褶皺和斷彎褶皺概念模式(據(jù)Suppe, 1983; Mercier et al., 1997修改)Fig.5　Conceptual model for thrust induced fault-propagation fold and fault-bend fold

將推斷地質(zhì)界面數(shù)據(jù)和控制地質(zhì)界面數(shù)據(jù)輸入到一個由Micromine和GOCAD相互耦合的計算平臺, 采用知識驅(qū)動下多級約束和多源數(shù)據(jù)融合方法進(jìn)行三維地質(zhì)界面模擬, 基本流程見圖3。整個模擬過程關(guān)鍵操作是人工插入控制點(diǎn), 多級約束和多源數(shù)據(jù)融合是為了保證人工插入控制點(diǎn)的合理性和可靠性, 主要內(nèi)容和操作包括如下4個方面。

(1) 用來保證插入點(diǎn)合理性和可靠性約束條件是多級的, 包括: ①來自工程點(diǎn)位置嚴(yán)格約束, 即插入點(diǎn)必須嚴(yán)格地與工程控制點(diǎn)相協(xié)調(diào); ②來自地質(zhì)現(xiàn)象的產(chǎn)狀和形態(tài)變化趨勢的約束, 即插入點(diǎn)必須符合這種產(chǎn)狀和形態(tài)變化趨勢; ③來自概念模式輪廓約束, 即插入點(diǎn)與已有控制點(diǎn)的空間展布應(yīng)該與概念模式所反映的輪廓相吻合; ④來自三維電阻率場的塊段約束, 即插入點(diǎn)與已有點(diǎn)所包圍塊體屬性應(yīng)該與電阻率場所反映屬性一致。所有插入邊界點(diǎn)要與這4級約束相符, 且同一屬性插入點(diǎn)和控制點(diǎn)分布要大致均勻。

(2) 用于確定邊界點(diǎn)位置和屬性的數(shù)據(jù)是多源的, 包括: ①鉆孔和坑道編錄確定的邊界數(shù)據(jù); ②地表地質(zhì)圖確定的邊界數(shù)據(jù); ③勘探線剖面圖和中段平面地質(zhì)圖圈定的邊界數(shù)據(jù); ④電阻率數(shù)據(jù)以及根據(jù)電阻率推斷的邊界數(shù)據(jù)。這些不同來源數(shù)據(jù)形式和可靠程度都不一樣。將上述各類數(shù)據(jù)按照Micromine的點(diǎn)(位置, 屬性)格式統(tǒng)一于同一坐標(biāo)空間內(nèi)進(jìn)行編輯和操作。根據(jù)原始數(shù)據(jù)約束條件差異, 可將整個模擬區(qū)域分成多個子區(qū)域, 在各個子區(qū)域中通過控制點(diǎn)、解譯線與輔助面混合使用, 插入盡可能合理邊界點(diǎn), 然后再將各子區(qū)域融合成一個統(tǒng)一的點(diǎn)數(shù)據(jù)庫。

(3) 將所有點(diǎn)數(shù)據(jù)按 GOCAD格式要求導(dǎo)入GOCAD進(jìn)行DSI處理, 然后再由Delaunay算法生成相應(yīng)的曲面。

(4) 在整個模擬過程中, 知識驅(qū)動主要體現(xiàn)在兩個方面: ①根據(jù)電阻率變化推斷地質(zhì)解釋的選擇,由于電阻率數(shù)據(jù)的地質(zhì)解釋具有多解性, 根據(jù)概念模式來選擇一種合理解釋; ②松散約束條件下多種可能模型選擇。由于礦田深部嚴(yán)格工程約束稀少, 主要是來自概念模式和三維電阻率場的松散約束, 可能構(gòu)建多種可選模型, 最后只能根據(jù)構(gòu)造動力學(xué)背景知識和已清楚的礦田構(gòu)造特征選擇一個合理模型。

3　礦田構(gòu)造三維模擬結(jié)果及找礦意義分析

應(yīng)用上述計算模擬方法, 我們最終模擬銅官山礦田地表至–2000 m范圍的下列地質(zhì)界面: ①銅官山巖體的邊界面; ②天鵝抱蛋巖體邊界面; ③金口嶺巖體的東部邊界面; ④D3w與C2+3分界面; ⑤P1q與 P1g分界面; ⑥天馬山礦床主要礦體的邊界面;⑦銅官山礦床主要礦體的邊界面; ⑧深部逆沖推覆斷層, 模擬結(jié)果見圖6。

綜合分析地質(zhì)構(gòu)造三維計算模擬結(jié)果(圖 6)和電阻率分布三維計算模擬結(jié)果(圖4), 揭示銅官山礦田構(gòu)造在深部空間展布具有如下規(guī)律, 這些規(guī)律對指示深部找礦具有重要意義。

(1) 銅官山巖體和天鵝抱蛋巖體作為銅官山礦田兩個最主要控礦因素, 其空間形態(tài)變化復(fù)雜, 總體特征是位于銅官山背斜西北翼靠核部, 向西傾,并且向深部快速變小, 在–2000 m深度以上沒有與這兩個小巖株相連的巖漿房。

(2) 銅官山礦田最主要已有礦床是與銅官山巖體相關(guān)的銅官山銅礦床以及與天鵝抱蛋巖體相關(guān)的天馬山硫金礦床。兩個礦床的控礦構(gòu)造及空間變化差異很大, 銅官山銅礦床礦體連續(xù)性較好, 主要受控于圍巖巖性和接觸帶形態(tài), 淺部接觸帶平緩, 圍巖為中–上石炭統(tǒng)碳酸鹽巖, 成礦好。向深部, 接觸帶變陡, 而且?guī)r體直接與泥盆系砂巖接觸, 成礦條件差, 所以銅官山巖體東接觸帶深部應(yīng)該沒有找礦潛力。天馬山礦床礦體主要位于與天鵝抱蛋巖體相接觸的中–上石炭統(tǒng)中, 礦體規(guī)模小, 連續(xù)性差, 礦體主要受控于層位, 而不是巖體的接觸帶, 礦體定位空間與D3w-C2+3分界面形態(tài)密切相關(guān), D3w-C2+3分界面呈現(xiàn)出一系列沿傾向延伸的凹槽和凸起, 礦體趨向定位于凹槽靠凸起的部位, 向深部至–1200 m到–1500 m仍具成礦條件, 具有找礦潛力。

(3) 銅官山背斜東南翼地層在深部向北西倒轉(zhuǎn),但在–1200 m到–1500 m深度, 倒轉(zhuǎn)地層又會返轉(zhuǎn)成向東南緩傾的正常產(chǎn)狀, 返轉(zhuǎn)部位是北高南低, 意味著本區(qū)最有利成礦的中–上石炭統(tǒng)和下二疊統(tǒng)棲霞組在天馬山礦床以東呈現(xiàn)向東南緩傾的產(chǎn)狀, 位于–1200 m到–1500 m的深度。

(4) 深部最主要構(gòu)造是彎曲的逆沖推覆斷層,這個斷層向北西陡傾部分切穿地層, 而東南緩傾部分則平行于地層, 而且是追蹤D3w與C2+3的分界面。沿此斷層陡傾部分是巖體侵位的有利場所, 電阻率三維分布模型也顯示這一構(gòu)造部位可能存在隱伏巖體。所以, 在礦田東部深部緩傾的 C2+3和P1q應(yīng)該是找礦有利靶區(qū)。

4　結(jié)　論

基于銅官山礦田勘探資料及大地電磁測深電阻率數(shù)據(jù), 采用知識驅(qū)動下多級約束和多源數(shù)據(jù)融合方法模擬礦田范圍在–2000 m 以上主要地質(zhì)界面,形象地展示礦田深部三維構(gòu)造格架。結(jié)果顯示: 銅官山巖體和天鵝抱蛋巖體形態(tài)變化復(fù)雜, 總體上向西傾, 向深部巖體規(guī)模快速變小。位于銅官山巖體東接觸帶的銅官山銅礦深部找礦潛力小, 原因是接觸帶產(chǎn)狀變陡, 并且圍巖也變成上泥盆統(tǒng)砂巖; 礦田深部最主要構(gòu)造為走向 NE并沿傾向彎曲的逆沖推覆斷層, 斷層?xùn)|南段向東南緩傾; 銅官山背斜東翼含礦層位中–上石炭統(tǒng)向深部倒轉(zhuǎn)后又在–1200 m至–1500 m深度返轉(zhuǎn)成向東南緩傾的正常層位, 其底部與深部逆沖推覆斷層重合, 是深部找礦有利靶區(qū)。

圖6　礦田構(gòu)造三維模型Fig.6　3D model of the geological structures in the ore field

致謝: 研究工作得到了國家自然科學(xué)基金和銅陵有色金屬集團(tuán)股份有限公司的資助, 銅陵有色金屬集股份有限公司還對現(xiàn)場研究工作給予了大力支持,席振銖和范福來等承擔(dān)了大地電磁測深工作, 在此深表感謝！感謝審稿人合肥工業(yè)大學(xué)袁峰教授和中國地質(zhì)大學(xué)(武漢)左仁廣教授所提的修改意見！

陳國達(dá). 1985. 成礦構(gòu)造研究法(第二版). 北京: 地質(zhì)出版社: 1–4.

董樹文, 邱瑞龍. 1993. 安慶–月山地區(qū)構(gòu)造作用與巖漿活動. 北京: 地質(zhì)出版社: 1–58.

杜楊松, 李順庭, 曹毅, 秦新龍, 樓亞兒. 2007. 安徽銅陵銅官山礦區(qū)中生代侵入巖的形成過程: 巖漿底侵、同化混染和分離結(jié)晶. 現(xiàn)代地質(zhì), 21(1): 71–77.

杜楊松, 秦新龍, 田世洪. 2004. 安徽銅陵銅官山礦區(qū)中生代巖漿–熱液過程: 來自巖石包體及其寄主巖的計據(jù). 巖石學(xué)報, 20(2): 339–350.

郭文魁. 1957. 論安徽銅官山銅礦成因. 地質(zhì)學(xué)報, 37(3): 317–323.

雷世和. 1985. 銅官山銅礦主控礦斷裂構(gòu)造的特征及找礦意義. 河北地質(zhì)學(xué)院學(xué)報, 30(2): 10–21.

李德仁, 李清泉. 1997. 一種三維 GIS 混合數(shù)據(jù)結(jié)構(gòu)研究. 測繪學(xué)報, 26(2): 128–133.

李紅陽, 李英杰, 康桂玲, 曹積飛, 楊竹森, 侯增謙, 蒙義峰. 2008. 安徽銅官山銅–鐵–金–硫礦床的地球化學(xué)特征. 地質(zhì)科學(xué), 43(2): 370–376.

劉亮明. 2007. 成礦理論的預(yù)測能力及其改善途徑. 地學(xué)前緣, 14(5): 82–91.

劉亮明, 蔡愛良. 2009. 深部找礦: 地質(zhì)理論、勘查技術(shù)、投資決策的難題及對策. 地質(zhì)通報, 28(7): 923–932.

劉亮明, 呂俊武, 彭省臨, 王國平, 方維萱, 邵擁軍. 2005.成熟勘探區(qū)新一輪找礦: 勘查戰(zhàn)略創(chuàng)新及銅陵礦集區(qū)找礦實(shí)例. 地質(zhì)論評, 51(3): 325–333.

田世洪, 丁悌平, 侯增謙, 楊竹森, 謝玉玲, 王彥斌, 王訓(xùn)誠. 2005. 安徽銅陵小銅官山銅礦床稀土元素和穩(wěn)定同位素地球化學(xué)研究. 中國地質(zhì), 32(4): 604–613.

銅陵有色金屬集團(tuán)股份有限公司礦產(chǎn)資源中心. 2011. 安徽省銅陵市銅官山礦區(qū)地形地質(zhì)圖(1: 5000).

萬仁虎. 2004. 銅官山礦田金礦床類型和構(gòu)造控礦作用分析. 地質(zhì)力學(xué)學(xué)報, 10(2): 122–128.

吳才來, 高前明, 國和平, 郭祥炎, 劉良根, 郜源紅, 雷敏, 秦海鵬, 陳其. 2010. 銅陵地區(qū)銅官山礦田侵入巖鋯石SHRIMP定年. 地質(zhì)學(xué)報, 84(12): 1746–1758.

薛虎, 車樹政, 董茂興, 張富陶. 1983. 安徽銅官山–馬山金屬礦田的構(gòu)造研究. 中國地質(zhì)科學(xué)院南京地質(zhì)礦產(chǎn)研究所所刊, 4(3): 1–24.

楊秋榮, 王金芳, 馮景志, 侯增謙, 楊竹森, 蒙義峰, 王聰穎, 陳文韜. 2010. 安徽銅陵天馬山金硫礦床地質(zhì)地球化學(xué)特征. 地球?qū)W報, 31(2): 203–208.

楊學(xué)明, 林文通. 1988. 銅官山火成雜巖體成巖機(jī)理研究.地質(zhì)論評, 34(1): 25–35.

翟裕生, 鄧軍, 王建平, 彭潤民, 劉家軍, 楊立強(qiáng). 2004.深部找礦研究問題. 礦床地質(zhì), 23(2): 142–149.

趙乘乘, 楊曉勇, 馮敏, 孫衛(wèi)東, 洪長春, 孫健, 胡俊杰,唐樹平. 2012. 銅官山巖體礦物學(xué)–礦物化學(xué)特征:巖漿結(jié)晶動力學(xué)意義. 地質(zhì)學(xué)報, 86(11): 1748–1760.

趙義來, 劉亮明. 2011. 復(fù)雜形態(tài)巖體接觸帶成礦耦合動力學(xué)三維數(shù)值模擬: 以安慶銅礦為例. 大地構(gòu)造與成礦學(xué), 35(1): 128–136.

周泰禧, 李學(xué)明, 趙俊深, 張富陶. 1987. 安徽銅陵銅官山礦田火成巖的同位素地質(zhì)年齡. 中國科學(xué)技術(shù)大學(xué)學(xué)報, 17(3): 403–407.

Bistacchi A, Massironi M, Dal Piaz G, Dal Piaz G, Monopoli B, Schiavo A and Toffolon G. 2008. 3D fold and fault reconstruction with an uncertainty model: An example from an Alpine tunnel case study. Computers and Geosciences, 34(4): 351–372.

Bouza-Rodríguez B, Comesa?a-Campos A, Menéndez-Díaz A and Garcia-Cortes S. 2014. A novel geometric approach for 3-D geological modelling. Bulletin of Engineering Geology and the Environment, 73(2): 551–567.

Calcagno P, Chilès J P, Courrioux G and Guillen A. 2008. Geological modelling from field data and geological knowledge Part I: Modelling method coupling 3D potential-field interpolation and geological rules. Physics of the Earth and Planetary Interiors, 171(1): 147–157.

Caumon G, Collon-Drouaillet P, Le Carlier, de Veslud C, Viseur S and Sausse J. 2009. Surface-based 3D modeling of geological structures. Mathematical Geosciences, 41(8): 927–945.

De Floriani L. 1989. A pyramidal data structure for trianglebased surface description. IEEE Computer Graphics and Application, 8(2): 67–78.

Delaunay B. 1934. Sur la sphere vide. Izv. Akad. Nauk SSSR, Otdelenie Matematicheskii Estestvennyka Nauk, 7: 793–800.

Dunstan S P and Mill A J B. 1989. Spatial indexing of geological models using linear Octrees. Computers and Geosciences, 15(8): 1291–1301.

Elsheikh A H and Elsheikh M. 2014. A reliable triangular mesh intersection algorithm and its application in geological modelling. Engineering with Computers, 30(1): 143–157.

Fisher T R and Wales R Q. 1991. Three-dimensional solid modeling of geological objects using non-uniformrational B-splines (NURBS) // Turner A K. Three dimensional modelling with geoscientific information systems. Dordrecht: Kluwer: 85–105.

Frank T, Tertois A L and Mallet J L. 2007. 3D-reconstruction of complex geological interfaces from irregularly distributed and noisy point data. Computers and Geosciences, 33(7): 932–943.

George P L and Borouchaki H. 1998. Delaunay Triangulation and Meshing: Application to Finite Elements. Paris: Hermes: 1–410.

Gong J, Cheng P and Wang Y. 2004. Three-dimensional modeling and application in geological exploration engineering. Computers and Geosciences, 30(4): 391–404.

Houlding S W. 1994. 3D Geoscience Modeling: Computer Techniques for Geological Characterization. Berlin: Springer: 1–301.

Kaufmann O and Martin T. 2008. 3D geological modelling from boreholes, cross-sections and geological maps, application over former natural gas storages in coal mines. Computers and Geosciences, 34(3): 278–290.

Liu L M, Zhao Y L and Sun T. 2012. 3D computational shape-and cooling process-modeling of magmatic intrusion and its implication for genesis and exploration of intrusion-related ore deposits: An example from the Yueshan intrusion in Anqing, China. Tectonophysics, 526–529: 110–123.

Lorensen W E and Cline H E. 1987. Marching cube: A high resolution 3d surface construction algorithm. Computer Graphics, 21(4): 163–169.

Mallet J L. 1989. Discrete smooth interpolation. ACM Transactions on Graphics, 8(2): 121–144.

Mallet J L. 1992. Discrete smooth interpolation in geometric modelling. Computer-aided design, 24(4): 178–191.

Mallet J L. 2002. Geomodeling. Oxford: Oxford Press: 1–600.

Marschallinger R. 1996. A voxel visualization and analysis system based on Autocad. Computer and Geosciences, 22(4): 379–386.

Mercier E, Outtani F and Lamotte D F D. 1997. Later-stage evolution of fault-propagation folds: Principles and example. Journal of Structural Geology, 19(2): 185–193.

Ming J, Pan M, Qu H and Ge Z. 2010. GSIS: A 3D geological multi-body modeling system from netty cross-sections with topology. Computers and Geosciences, 36(6): 756–767.

Pilout M, Tempfli K and Molenaar M. 1994. A tetrahedron-based on 3D vector data model for geoinformation // Molenaar M, de Hoop S. Advanced Geographic Data Modeling, Publications on Geodesy Vol. 40. Delft: Netherlands Geodetic Commission: 129–140.

Rebay S. 1993. Efficient unstructured mesh generation by means of Delaunay triangulation and Bowyer-Waston algorithm. Journal of Computational Physics, 106(1): 125–138.

Rice J R and Usowt K H. 1968. The Lawson algorithm and extensions. Mathematics of Computation, 22(101): 118–127.

Saalmann K and Laine E L. 2014. Structure of the Outokumpu ore district and ophiolite-hosted Cu-Co-Zn-Ni-Ag-Au sulfide deposits revealed from 3D modeling and 2D high-resolution seismic reflection data. Ore Geology Reviews, 62: 156–180.

Shi W Z. 1996. A hybrid model for 3D GIS. Geoinformatics, 1: 400–409.

Stroud I. 2006. Boundary representation modelling techniques. Berlin: Springer: 1–808.

Suppe J. 1983. Geometry and kinematics of fault-bend folding. American Journal of Science, 283(7): 684–721.

Suppe J and Medwedeff D A. 1990. Geometry and kinematics of fault-propagation folding. Eclogae Geologicae Helvetiate, 83(3): 409–459.

Turner A K. 2006. Challenges and trends for geological modelling and visualisation. Bulletin of Engineering Geology and the Environment, 65(2): 109–127.

Victor J D. 1993. Delaunay triangulation in TIN creation: An overview and a linear-time algorithm. International Journal of Geographical Information Systems, 7(6): 501–524.

Wellmann J F, Horowitz F G, Schill E and Regenauer-Lieb K. 2010. Towards incorporating uncertainty of structural data in 3D geological inversion. Tectonophysics, 490(3): 141–151.

Wellmann J F and Regenauer-Lieb K. 2012. Uncertainties have a meaning: Information entropy as a quality measure for 3-D geological models. Tectonophysics, 526–529: 207–216.

Williams N. 2005. Australian mineral exploration in the 21st century: Are we up to the challenge? Perth: Presentation at AMEC National Mining Congress.

Wu Q, Xu H and Zou X. 2005. An effective method for 3Dgeological modeling with multi-source data integration. Computers and Geosciences, 31(1): 35–43.

Xue Y, Sun M and Ma A. 2004. On the reconstruction of three-dimensional complex geological objects using Delaunay triangulation. Future Generation Computer System, 20(7): 1227–1234.

Zanchi A, Francesca S, Stefano Z, Simone S and Graziano G. 2009. 3D reconstruction of complex geological bodies: Examples from the Alps. Computers and Geosciences, 35(1): 49–69.

Computational 3D Modeling on Deep Structure Architecture and Implication for Ore Exploration in the Tongguanshan Ore Field, Tongling, China

LIU Liangming and CAO Wei
(Computational Geosciences Research Center, MOE Key Laboratory of Nonferrous Metallogenic Prediction, Central South University, Changsha 410083, Hunan, China)

Computational modeling has become an important virtual technique for revealing the deep structures in 3D, which is critical for deep ore exploration. In view of the geological setting and the exploration data available for the Tongguanshan ore field, we simulate the major geological interfaces within the ore field deep to ?2000 m by applying the knowledge-driven modeling method, during which multi-sourced data are integrated under multi-scaled constraints. The concept model for thrust induced fault-propagation fold and fault-bend fold is used to interpret the MT resistivity data, and the inferred positions of geological interfaces are used where no drills or tunnels are available. The multi-sourced data about the geological interfaces are integrated under multi-scaled constraints in the soft platform of coupled Micromine and GOCAD. The algorithms of DSI and Delaunay TIN are applied to reconstruct the major geological interfaces for virtually revealing the structure architecture of the ore field. The final modeling results demonstrate that the Tongguanshan and Tian’ebaodan intrusions, with very complex spatial variation in 3D geometric shape, are generally west-dipping and diminished sharply at depth. The majority of the ore bodies of the Tongguanshan deposit occur in the shallow part of the east contact zone of the intrusion, however, the economic potential of ore mineralization at depth is permissive because the steeper and directly contact with the Upper Devonian sandstone. The principal structure deep in the ore field is the NE-trending bend thrust fault. The ore-hosting strata, the Mid- to Upper Carboniferous sequences, in the southeastern limb of the Tongguanshan upturn anticline will be appear as gently dipping strata and cut by the bend thrust fault at the depth of ?1200 m to ?1500 m, and thus is likely favorable target for ore prospecting.

deep structure; 3D model; computational modeling; Tongguanshan ore field

P542

A

1001-1552(2016)05-0928-011

10.16539/j.ddgzyckx.2016.05.003

2015-03-05; 改回日期: 2015-04-01

項目資助: 國家自然科學(xué)基金項目(41240017和41372338)資助。

劉亮明(1964–), 男, 教授, 主要從事礦床地質(zhì)和構(gòu)造地質(zhì)研究。Email: lmliu@csu.edu.cn