銅陵獅子山礦田矽卡巖銅(金)成礦系統三維動力學計算模擬及其意義

周飛虎, 劉亮明*, 劉洪升, 盛 歡, 伍紅志, 閔令帥, 畢晨曦

銅陵獅子山礦田矽卡巖銅(金)成礦系統三維動力學計算模擬及其意義

周飛虎1, 劉亮明1*, 劉洪升1, 盛 歡2, 伍紅志1, 閔令帥1, 畢晨曦1

(1. 中南大學 計算地球科學研究中心, 教育部有色金屬成礦預測與地質環境監測重點實驗室,湖南 長沙 410083; 2. 銅陵有色金屬集團股份有限公司冬瓜山銅礦, 安徽 銅陵 244000)

獅子山礦田產有銅陵地區乃至整個長江中下游成礦帶最大的矽卡巖型銅(金)礦床, 由于成礦特征的復雜性, 礦田內巨大的層狀矽卡巖礦床的定位機制一直存在爭議。為了探討本區成礦的動力學機制及其控制因素, 本次研究基于礦田內全部的勘探資料建立了礦田三維地質模型, 在此基礎上建立動力學模型, 在FLAC3D平臺進行巖體冷卻過程力?熱?流耦合動力學計算模擬。模擬結果顯示在拉張構造應力環境中, 在淺部東、西獅子山礦床產出部位只在巖體接觸帶上形成了規模不大的高擴容區(體積應變值≥1.2%), 而在深部層狀礦體產出的地層中形成了大規模的層狀高擴容變形區, 同時流體也向高擴容區匯聚。這些流體匯聚的高擴容區與已知礦體的容礦空間高度吻合, 由此推斷礦體的定位受控于巖體冷卻過程中產生的擴容和流體流動。石炭系中流體匯聚的高擴容空間的規模遠大于已知的成礦規模, 表明礦田深部的冬瓜山礦床仍存在較大的增儲潛力。改變模型的邊界條件和模型組成單元的特征參數進行假設模型計算實驗, 結果表明拉伸的構造應力環境為大規模匯流容礦空間的形成提供了有利條件, 并且賦礦地層中更弱的抗剪切強度和上覆地層低滲透的孔隙特征也是控制礦體定位空間和成礦規模的重要因素。

成礦動力學; 計算模擬; 矽卡巖礦床; 獅子山礦田; 銅陵

0 引 言

獅子山礦田位于長江中下游成礦帶銅陵礦集區, 是整個長江中下游成礦帶中規模最大的矽卡巖銅(金)礦田(常印佛等, 1991; Pan and Dong, 1991)。礦田內礦體的形態和空間分布規律復雜多變, 眾多不同規模的不規則狀、囊狀和層狀礦體圍繞冬瓜山閃長巖體不均勻分布。其中最引人注目的特征是冬瓜山深部巨大的層狀矽卡巖型礦體, 僅I號礦體的銅儲量就超過了100萬噸。由于礦床的巨大規模以及復雜的、多層分布的大規模層狀矽卡巖銅礦體, 吸引了眾多學者的關注, 已開展研究包括礦床的特征、控礦因素與成礦模式(常印佛和劉學圭, 1983; 儲國正, 1992; 黃許成, 1993; 陸建軍等, 2008)、礦床與成礦巖體的地球化學特征和成巖成礦演化過程(肖建新等, 2002; 邱士東等, 2007; 徐兆文等, 2007; 徐曉春等, 2008; 侯增謙等, 2011; Xu et al., 2011; Yang et al., 2011; 郭維民等, 2013; Wang et al., 2015; Cao et al., 2017; Liu et al., 2018; Xie et al., 2020; Xiao et al., 2021)、成礦機制(於崇文等, 1995; 凌其聰等, 2002; 徐文藝等, 2004; Liu et al., 2014)以及成礦的區域動力學背景等(吳淦國, 2003; 毛景文等, 2009; 董樹文等, 2011; Deng et al., 2011), 這些研究促進了對獅子山礦田和冬瓜山礦床的成礦規律和成礦機制的認識。但關于巨大的層狀矽卡巖礦床的成因機制仍存在爭議, 有噴流沉積成因說(顧連興和徐克勤, 1986; Gu et al., 2007)、噴流沉積+巖漿熱液(矽卡巖)疊加改造成因說(肖新建等, 2002; 徐兆文等, 2007; 陸建軍等, 2008; 侯增謙等, 2011)和巖漿熱液(矽卡巖)成因說(毛景文等, 2009; Liu et al., 2014; Cao et al., 2017; Xie et al., 2021)。盡管越來越多的研究支持巖漿熱液(矽卡巖)成因說, 但并沒有對成礦的動力學過程和關鍵控制因素進行詳細的論證, 對于形成如此大規模的層狀矽卡巖型礦體的關鍵控制因素仍不清楚。其主要原因在于控制礦體空間定位和規模的動力學機制復雜, 傳統的地質研究方法是依據礦床的地質地球化學特征來反推礦床形成的動力學過程, 其結論帶有極大的不確定性, 并難以驗證。

計算機和計算科學的發展使計算模擬成為研究復雜地質系統動力學的有效手段(Hobbs et al., 2000; Price and Stoke, 2002; Zhao et al., 2009; Liu et al., 2010)。利用計算機強大的計算能力和虛擬時空表達能力再現復雜地質系統的時空演化, 通過不斷試錯使模擬的復雜動力學系統更接近真實。因此, 與傳統的地質研究方法相比, 動力學計算模擬在結論的確定性和可靠性方面具有明顯的優勢, 因此在復雜的巖漿熱液成礦系統的動力學研究方面得到了廣泛的應用, 深化了對巖漿熱液成礦動力學機制的認識, 特別是巖漿熱液在不同的位置形成不同規模、不同類型礦體的動力學過程及關鍵控制因素(Norton and Taylor, 1979; Hanson, 1996; Eldursi et al., 2008; Ingesbritsen et al., 2010; Liu et al., 2011, 2012, 2014; Weis et al., 2012; Sch?pa et al., 2017; Yapparova et al., 2019; Korges et al., 2020; Lamy-Chappuis et al., 2020; Vigneresse and Truche, 2020; Annen and Burgisser, 2021)。近年, 動力學計算模擬技術得到了快速發展, 由早期二維簡單模型的單一過程模擬發展成三維復雜模型的多過程耦合動力學模擬(Ingebritsen et al., 2010; Liu et al., 2012; Vigneresse and Truche, 2020)。能模擬復雜的動力學系統的商業和科研軟件有: Itasca公司的有限差分軟件FLAC與FLAC3D軟件(http://www.itasca.com/)、Ansys公司的ANSYS和FLUENT系列軟件(http://www.ansys.com/)、COMSOL公司的COMSOL系列軟件(http://www.comsol.com/)、蘇世理工的CSMP++(Yapparova et al., 2019)、美國Lawrence Berkeley國家實驗室的TOUGHREACT軟件等(Xu et al., 2011)。利用這些專業的軟件, 學者們對成礦過程中的構造變形、流體演化及成礦預測等開展了卓有成效的研究(Zhang et al., 2012; 賈蔡等, 2014; Li et al., 2017; 袁峰等, 2019; Li et al., 2019; Vigneresse and Truche, 2020; Liu et al., 2021; 肖凡和王愷其, 2021)。

本文以FLAC3D軟件為工具, 模擬獅子山礦田成礦巖體冷卻過程中的力?熱?流耦合動力學, 并根據模擬結果分析礦體空間定位的動力學機制及關鍵控制因素。

1 礦田地質背景與礦床特征

獅子山礦田位于銅陵礦集區中部, 距銅陵市中心東南方向約7 km。從地表往深部, 礦田內地層依次為中三疊統月山組, 下三疊統南陵湖組、和龍山組和殷坑組, 上二疊統大隆組和龍潭組, 下二疊統孤峰組和棲霞組, 中?上石炭統和上泥盆統五通組(圖1)。這些地層已發生褶皺變形, 形成軸向NE的舒緩背斜(圖2)。從區域地殼構造特征分析, 褶皺變形是中三疊世?中侏羅世印支運動的產物(董樹文等, 2011), 早于本區的礦床和成礦巖體的年齡。

圖1 獅子山礦田地層層序及賦存的礦床規模(據安徽321地質隊, 2014)

圖2 銅陵獅子山礦田礦體的水平投影、動力學模型范圍(a; 據銅陵冬瓜山礦床內部資料)和典型礦床剖面圖(b)

與銅陵礦集區的其他礦田一樣, 獅子山礦田最重要的地質要素也是巖體。雖然地表看似有青山腳、獅子山、大團山、胡村、白芒山、曹山、龍虎山等多個巖體(圖2), 但實際上, 這些巖體在三維空間上是相連成一體的。巖體在地表的出露范圍并不大, 其三維空間形態極為復雜, 呈現為NE和NW兩個方向巖墻交織的復雜網狀形態(圖2a)。巖性主要為閃長巖和石英閃長巖, 部分地段發育隱爆角礫巖(圖2a)。巖石地球化學特征表明巖體具有高鋁和富鉀的特征, 巖漿來源單一, 主要為巖石圈地幔和地殼物質混染的產物(Wang et al., 2015; Liu et al., 2018b),并且具有較高的氧逸度, 有利于形成矽卡巖型礦床(Xie et al., 2012; Zhang et al., 2021)。鋯石U-Pb測年結果顯示成巖年齡集中在136～148 Ma(徐曉春等, 2008; 郭維民等, 2013; Liu et al., 2018b), 巖漿作用發生于燕山早期, 此時區域地殼由擠壓轉為伸展(董樹文等, 2011)。礦石的Re-Os同位素年齡為138～140.8 Ma(毛景文等, 2005), 表明成礦時間稍晚于成巖時間。與所有的矽卡巖礦床一樣, 獅子山礦田的成礦作用也發生在巖漿結晶冷卻的退化蝕變階段。

獅子山礦田成礦巖體的規模并不大, 但成礦規模卻很廣, 礦體的投影范圍遠大于巖體出露的面積(圖2a)。成礦作用圍繞巖體不均勻發育, 形成了東獅子山、西獅子山、大團山、老鴉嶺、花樹坡和冬瓜山礦床(圖2)。但是礦床不是全部集中在巖體接觸帶上, 而是分布于不同的地層: 大團山的礦體產于下三疊系殷坑組; 老鴉嶺的礦體產于上二疊統大隆組和龍潭組; 花樹坡的礦體產于下二疊統棲霞組, 部分跨至孤峰組底部; 冬瓜山的礦體產于中?上石炭統(圖1、2c)。隨著賦礦層位深度的增加, 礦床的規模逐漸增大, 其中冬瓜山I號礦體的規模最大, 銅儲量在100萬噸以上, 約占整個礦田儲量的60%(圖1)。從淺部到深部, 礦石的結構構造發生了規律性的變化, 淺部的東獅子山和西獅子山礦床的礦石主要呈角礫狀構造, 角礫有矽卡巖、閃長巖和大理巖, Cu主要賦存于黃銅礦中, 與其他熱液礦物充當角礫的膠結物, 或產于裂隙之中。大團山和老鴉嶺礦床也發育角礫狀礦石, 但不是主要類型。再往深部, 花樹坡和冬瓜山礦床, 已沒有角礫狀礦石, 但在塊狀硫化物礦石中發育圍巖角礫, 受稠密裂隙控制的浸染狀礦石在深部也成為一種重要類型。這種礦石構造的變化規律反映了形成容礦空間的巖石破裂機制和方式的變化, 與應力狀態及孔隙流體壓力的變化密切相關。除了不同礦床間具有垂向差異特征外, 礦床內部的礦石類型也表現出明顯的分帶特征。冬瓜山礦床主礦體的組成由頂板到底板分別為: ①含銅矽卡巖→②含銅黃鐵礦礦石+少量含銅磁黃鐵礦礦石→③含銅矽卡巖+少量含銅硬石膏巖→④含銅磁黃鐵礦礦石+少量含銅黃鐵礦礦石→⑤含銅矽卡巖+含銅硬石膏巖+含銅磁鐵礦礦石→⑥含銅磁黃鐵礦礦石+少量含銅黃鐵礦礦石→⑦含銅蛇紋巖→⑧含銅粉砂巖(唐永成等, 1998; 侯增謙等, 2011; 徐曉春等, 2014)。這種垂直于層理方向的分帶規律性, 不僅表明礦體上下盤的圍巖參與了成礦過程的水?巖反應, 也表明在成礦過程中流體不僅從巖體往圍巖流動, 也存在垂直于圍巖層理方向的流動。

2 成礦系統三維模型及特征

2.1 三維建模方法及數據

采用知識驅動多級約束多源數據融合的顯式建模方法(liu et al., 2012; 劉亮明和曹偉, 2016), 利用MicroMine和Gocad三維建模軟件, 融合礦田內所有勘探數據, 建立了成礦巖體的接觸面、主要地層界面(D3-C2+3界面、C2+3-P1界面、P1-P2界面、P2-P2界面、P2-T界面)及礦體三維面元模型。地質界面以不規則三角網(TIN)表示。所使用的勘探數據包括: 2159個鉆孔總進尺347567 m的巖心編錄; 939832 m巷道編錄; 58523個樣品的基本數據; 125張較老的沒有編錄與分析數據庫的剖面和中段地質圖(圖3)。建模所用數據基本覆蓋了整個三維建模研究區, 此次三維建模過程中, 地質界線不需要做過多的外推。

圖3 建模所依據的鉆孔、坑道及地質圖的分布(未顯示全部地質圖)

2.2 成礦系統的三維空間結構特征

獅子山礦田的三維模型(圖4)反映成礦系統的三維空間結構具有如下特征和規律:

(a) 獅子山礦田現代地表下三維地質模型; (b) 巖體與礦體模型; (c) 東、西獅子山礦床; (d) 大團山礦床; (e) 老鴉嶺、花樹坡礦床; (f) 冬瓜山礦床。

(1) 巖體的形態極為復雜, 總體上表現為由南東向北西侵入的巖墻狀。巖體的北西側接觸帶的形態也比其東南側接觸帶更復雜(圖4b)。

(2) 所有已知礦床均位于一個平緩的背斜褶皺中, 深部的冬瓜山礦床的礦體在背斜兩翼和轉折端都有, 并連成一個呈背斜形態的大礦體(圖4f), 其余礦床的礦體均產出于背斜的南東翼, 而且規模也比冬瓜山礦體小得多(圖4b、e)。

(3) 只有淺部的東獅子山和西獅子山礦床的礦體是產在巖體接觸帶上的(圖4c), 往深部主要是受層位控制的似層狀和層狀礦體, 而且礦體規?？傮w具有明顯往深部變大的趨勢(圖4d、f)。

(4) 不同產出層位的礦體與巖體的空間關系也不相同, 深部的冬瓜山礦床主要產出于巖體的北西盤(下盤)圍巖中, 其他礦床則主要產出于巖體的南東盤(上盤)圍巖中(圖4b)。

(5) 冬瓜山礦床的層狀礦體產于青山背斜的軸部及近軸部兩翼, 且其總體形態也呈背斜狀(圖4f)。正是因為這一特征, 許多研究者總是試圖證明礦床是沉積成因(顧連興和徐克勤, 1986; Gu et al., 2007),或者是在沉積的基礎上疊加巖漿熱液改造成因(肖新建等, 2002; 徐兆文等, 2007; 陸建軍等, 2008; 侯增謙等, 2011)。實際上, 層狀礦體的邊界并非像層理一樣平整, 而是極不平整的鋸齒狀, 不僅在塊狀硫化物礦體中發育有圍巖角礫, 也有不規則細脈狀礦化從主礦體一直延伸到圍巖中, 與層狀主礦體連為一體, 礦體的后生充填特征十分明顯(Liu et al., 2014)。

上述特征表明, 在獅子山礦田內, 不但礦體定位的空間規律非常復雜, 礦床形成的動力學過程及控制因素亦十分復雜, 這也正是本文開展動力學計算模擬研究的原因之所在。

3 動力學計算模擬的原理、方法及條件

3.1 動力學系統的數學描述

本次研究采用Itasca公司的FLAC3D軟件(Itasca, 2012)模擬巖體冷卻過程中的力?熱?流耦合動力學過程。模擬的地質體被看成是流體飽和的多孔彈塑性介質, 其本構關系符合Mohr-Coulomb定律。對于飽和流體多孔的Mohr-Coulomb材料, 其耦合的變形、流體流動和傳熱的動力學過程, 計算公式如下:

3.2 模擬的方法及離散化模型的建立

動力學計算模擬就是用數值的方法求解動力學方程組。FLAC3D軟件是有限差分軟件, 即利用差分的方法求解偏微分方程的數值解, 先通過模型的離散化將連續變量變成離散變量, 然后利用Taylor展開式求取偏微分方程的近似解(Quarteroni and Valli, 1994)。所以進行成礦動力學計算模擬的第一步是建立成礦系統的地質模型, 然后將模型離散化, 再對模型的特征參數、初始條件和邊界條件賦值, 最后用數值方法求解。

對區域地殼構造演化(常印佛等, 1991; 董樹文等, 2011)、成礦巖體的成巖年代(郭維民等, 2013)與成礦年代(毛景文等, 2009)的研究, 表明成礦作用發生至今被剝蝕掉的地層主要為侏羅系, 而冬瓜山礦床的成礦深度約為2.09～5.19 km(杜楊松等, 2004; 陸三明, 2007)。因此, 我們確定了動力學模型的總深度為2.5 km, 并在地質體三維面模型的基礎上恢復了被剝蝕的部分侏羅系和三疊系, 將礦田的三維模型還原為成礦作用發生時成礦系統模型的初始狀態。然后利用我們自主開發的剖分工具(Liu et al., 2012)將2000 m′3300 m′2500 m的三維模型剖分成478168個四面體單元(圖5)。

圖5 動力學數值模擬用的三維實體模型及力學邊界條件

3.3 模型的特征參數及模擬的初始條件與邊界條件

根據獅子山礦田各地質體單元的組成及其與成礦的時空關系, 我們將成礦系統的動力學模型劃分成巖體和7個地層單元(圖5, 表1)。各單元的密度、體積模量、孔隙度、滲透率、導熱率等特征參數是根據各地質體中巖石的組分、結構而選定。巖石物理實驗手冊(Sch?n, 1998)和FLAC3D應用手冊(Itasca Consulting Group, 2012)給出的各種巖石的物理特征參數一般有一個較大的可選范圍, 我們通過控制模型的參數變量, 對參數范圍內的值進行試算, 淘汰了使模型計算結果不穩定的參數值, 最后優選出最合理的特征參數組合(表1)。

表1 動力學模型的巖石物理及水文模型參數

模型初始孔隙壓力的選擇依據已普遍接受的Fournier模型(Fournier, 1999), 將巖體中孔隙壓力設為靜巖壓力, 圍巖中的孔隙壓力設置為靜水壓力。四周及底部的流體邊界條件設置為不透水邊界, 頂部設置為透水邊界以模擬含礦流體與大氣降水的共同作用。將地表溫度固定為25 ℃, 地層中按20 ℃/km的地溫梯度增溫。巖體底部的初始溫度設為650 ℃, 并且向上按10 ℃/km的梯度遞減。

已有的研究成果表明整個長江中下游成礦帶燕山期的斑巖?矽卡巖成礦作用發生在區域地殼由擠壓轉換成伸展的構造環境下(毛景文等, 2009; 董樹文等, 2011), 因此我們將模型的力學邊界條件設置為NE和SW方向的對稱拉伸(圖5), 拉伸速率設為1.65′10?8m/s, 這個速率高于真實的構造應變速率,主要是為了縮短計算的時間。

4 模擬結果及其意義

4.1 真實模型的模擬結果及其意義

對于上述根據獅子山礦田的地質特征及演化歷史而建立的真實的動力學模型(圖5), 進行同伸展巖體冷卻演化的動力學計算模擬, 即以顯式時間控制求解力?熱?流耦合動力學方程不同時間的數值解, 在計算機的虛擬時空中再現巖體冷過程中的力?熱?流耦合動力學過程的演化。模擬的結果展現出的力?熱?流空間變化規律, 能為了解本區的構造?巖漿控礦的動力學機制提供有力的支撐。

模擬所產生的最有意義的變化是體積應變, 巖體在同伸展的冷卻過程中能產生具有特征性的擴容應變分布, 最大最強的擴容應變區分布于中?上石炭統, 其次是下二疊統棲霞組和孤峰組及上二疊統大隆組, 而且這些地層中的擴容應變的最強部分在巖體接觸帶上(圖6a、d)。礦田中礦體的空間分布規律與擴容空間的分布規律相關(圖6、7), 提取體積應變值為1.2%的等值面與礦體空間位置對比, 可以發現工程揭露的礦體全部被等值面包裹, 無礦的空間則被排除在外(圖7)。體積應變＞1.2%的等值面則無法完全包裹已知礦體, ＜1.2%的等值面又會把許多工程揭露的無礦空間包裹在內, 我們認為體積應變值≥1.2%的高擴容空間對確定礦體空間定位和預測找礦靶區最具意義。這些高擴容空間也是流體匯流空間, 特別是中?上石炭統的高擴容帶, 進一步證實不同來源流體于容礦空間中混合反應是礦石沉淀的重要機制(劉亮明等, 2008)。

箭頭表示流體向量; (a) 模型體積應變場; (b) 巖體?P2d-C2+3體積應變場; (c) 三疊系中體積應變場; (d) 石炭系中體積應變場; (e) 巖體?P2d-C2+3孔隙壓力場; (f) 巖體?P2d-C2+3溫度場。

圖7 獅子山礦田內高擴容應變區(≥1.2%)分布規律及與礦體的空間關系

模擬結果也顯示, 在巖體冷卻過程中, 超壓的高溫巖體在降壓降溫的同時, 導致近巖體的圍巖升溫(圖6f), 但沒有導致近巖體的圍巖明顯的升壓, 而在靠近巖體的擴容空間反而形成明顯的低壓帶(圖6e), 高壓的巖漿流體向低壓擴容空間流動必將導致巨大的流體泄壓, 巖漿流降壓沸騰, 這也是促使礦石快速沉淀的一種重要機制。前人對本區的不同礦石類型中流體包裹體的組合、溫度、鹽度和壓力的研究也證實存在這一機制(肖新建等, 2002)。

動力學計算模擬結果能為礦體的預測提供有益的支撐。體積應變值≥1.2%的擴容帶是礦體定位的有利部位, 還沒開展過鉆探的高擴容空間可能是有潛力的容礦空間, 具有發現新礦體或擴大已發現礦體的可能。模擬結果顯示, 中?上石炭統的擴容空間大于已發現的冬瓜山礦床的規模(圖7), 特別是在青山背斜的南東翼還有部分區域沒有鉆探, 這些區域具有擴大礦體的可能。

4.2 假設模型的模擬結果及其意義

為了探索制約擴容應變和流體流動的因素和條件, 找出本區礦體定位規律的原因, 在保持模型幾何結構不變的情況下, 通過改變邊界條件和模型組成單元的特征參數進行假設模型的動力學計算模擬實驗。共設定了4個假設模型, 具體參數見表2。

表2 假設模型的參數及模擬條件

假設模型1施加NE-SW向擠壓: 結果以負體積應變(壓縮變形區)為主, 主要分布在上泥盆統五通組砂巖及頂部侏羅系蓋層中。石炭系靠近巖體接觸帶的位置出現正體積應變, 但很小(圖8a)。相比于拉伸模型中形成的大規模層狀匯流擴容區, 擠壓模型的擴容強度及規模則顯得比較局限, 并且在其他碳酸鹽巖地層中沒有形成明顯擴容區, 也沒有明顯的流體匯流規律。

箭頭表示流體向量; (a) NE-SW向擠壓; (b) 無構造邊界條件; (c) 增大了石炭系(C2+3)模型的剪切模量(4.4′1010 Pa); (d) 增大了上二疊統龍潭組(P2l)和三疊系模型的滲透率(25′10?12 m2)。

假設模型2無構造外力加載: 結果亦未形成明顯的匯流擴容區, 中?上石炭統的擴容很弱, 匯流不明顯(圖8b), 沒能產生有利于冬瓜山礦體定位的條件。

這兩個假設模型的模擬結果表明拉張的構造應力環境對獅子山礦田, 特別是對冬瓜山礦床的成礦作用的重要性。

假設模型3將中?上石炭統的剪切強度調大:整個獅子山礦田模型的正體積應變值變小, 并且高擴容應變在巖體和上二疊統大隆組(P2)產生, 而中?上石炭統中沒有形成層狀的高擴容區(圖8c), 亦沒能產生有利于冬瓜山礦體定位的條件。以上結果表明碳酸鹽巖地層中更弱的抗剪切強度是造成地層中體積應變不均勻分布的關鍵因素, 也說明容礦地層的力學強度對成礦空間的形成至關重要, 弱的剪切強度更容易使地層巖石發生剪切破裂, 形成層狀容礦空間, 從而有利于形成冬瓜山式的層狀礦體。

假設模型4將上二疊統龍潭組和三疊系的滲透率調大:巖體主要成了強擴容區, 所有地層擴容均小, 沒有層狀的高擴容區形成。流體運動更劇烈, 表現為由巖體流向圍巖, 但沒有明顯的匯流規律(圖8d)。以上結果表明地層的水文條件對形成層狀矽卡巖礦床的重要性, 碳酸鹽巖地層之上存在低滲透率地層更有利于形成層狀矽卡巖礦床。

對比假設模型與真實模型模擬結果的差異(圖8、6), 反差最大的是擠壓構造環境的模擬(圖8a)和無低滲透性地層的模擬(圖8d)。這也表明, 除了成礦巖體外, 拉伸的構造應力環境以及賦礦地層之上低滲透擋水層對形成層狀矽卡巖礦床也有至關重要的控制作用。

5 結 論

根據獅子山礦田成礦系統在三維形態模擬基礎之上的動力學計算模擬研究, 獲得了如下認識:

(1) 獅子山礦田的矽卡巖礦床表現出非常特殊的成礦規律, 即規模并不大、形態極其復雜的成礦巖體形成規模超大的礦體群, 雖然礦體都是圍繞巖體分布的, 但主要礦體都是受層位控制的層狀礦體, 并且位于最深部的層狀礦體的規模遠遠大于其上的各個礦體。造成這一規律的主要動力學原因是在伸展構造背景下, 巖體冷卻過程中因力?熱?流耦合而在接觸帶及圍巖中現在礦體定位的部位產生的匯流的高擴容空間。

(2) 中?上石炭統中匯流高擴容空間的規模要大于現在已發現的冬瓜山礦體, 表明冬瓜山礦床還具備擴大儲量的空間。

(3) 除了成礦巖體外, 區域地殼所處的拉張的環境, 賦礦巖層的低剪切強度以及賦礦巖層之上的低滲透性擋水層也是形成獅子礦田特殊成礦規律的重要因素。

致謝:感謝銅陵有色金屬集團股份有限公司冬瓜山銅礦的資助和支持, 特別是提供有關的勘探資料。感謝中山大學肖凡副教授和另一位匿名審稿人對本文提出的寶貴修改意見, 向兩位審稿專家致以誠摯的謝意。

安徽321地質隊. 2014. 安徽省銅陵市冬瓜山銅礦床資源儲量復核地質報告: 1–72.

常印佛, 劉湘培, 吳言昌. 1991. 長江中下游銅鐵成礦帶. 北京: 地質出版社: 1–379.

常印佛, 劉學圭. 1983. 關于層控式矽卡巖型礦床——以安徽省內下揚子坳陷中一些礦床為例. 礦床地質, 2(1): 11–20.

儲國正. 1992. 銅陵獅子山礦田構造及其控巖控礦作用研究. 安徽地質, 2(2): 1–14.

董樹文, 馬立成, 劉剛, 薛懷民, 施煒, 李建華. 2011. 論長江中下游成礦動力學. 地質學報, 85(5): 612–625.

杜楊松, 秦新龍, 李鉉具. 2004. 安徽銅陵地區中生代幔源巖漿底侵作用——來自礦物巨晶和巖石包體的證據. 巖石礦物學雜志, 23(2): 109–116.

顧連興, 徐克勤. 1986. 論長江中、下游中石炭世海底塊狀第一線化物礦床. 地質學報, 60(2): 176–188.

郭維民, 陸建軍, 蔣少涌, 章榮清, 招湛杰. 2013. 安徽銅陵獅子山礦田巖漿巖年代學、Hf同位素、地球化學及巖石成因. 中國科學(D輯), 43(8): 1268–1286.

侯增謙, 楊竹森, 呂慶田, 曾譜勝, 謝玉玲, 蒙義峰, 田世洪, 徐文藝, 李紅陽, 姜章平, 王訓成, 姚孝德. 2011.安徽銅陵冬瓜山大型銅礦: 海底噴流?沉積與矽卡巖化疊加復合成礦過程. 地質學報, 85(5): 659–686.

黃許成. 1993. 銅陵獅子山礦田多位一體(多層樓)模式. 礦床地質, 12(3): 221–230.

賈蔡, 袁峰, 張明明, 李曉暉, 周濤發, 邵尉, 鄭通科, 高道明. 2014. 寧蕪盆地白象山鐵礦床成礦作用過程數值模擬. 巖石學報, 30(4): 1031–1040.

凌其聰, 劉叢強, 鮑征宇. 2002. 地質流體的輸運?化學反應與成礦作用微觀機制研究——以安徽冬瓜山層控夕卡巖型銅礦床為例. 大地構造與成礦學, 26(1): 55–61.

劉亮明. 2011. 淺成巖體引發的流體超壓與巖石破裂及其對成礦的制約. 地學前緣, 18(5): 78–89.

劉亮明, 曹偉. 2016. 銅官山礦田深部構造的三維計算模擬及其找礦意義. 大地構造與成礦學, 40(5): 928–938.

劉亮明, 疏志明, 趙崇斌, 萬昌林, 蔡愛良, 趙義來. 2008. 矽卡巖礦床的匯流擴容空間控礦機制及其對深部找礦的意義: 以銅陵?安慶地區為例. 巖石學報, 24(8): 1848–1856.

陸建軍, 郭維民, 陳衛鋒, 蔣少涌, 李娟, 顏曉蓉, 徐兆文. 2008. 安徽銅陵冬瓜山銅(金)礦床成礦模式. 巖石學報, 24(8): 1857–1864.

陸三明. 2007. 安徽銅陵獅子山銅金礦田巖漿作用與流體成礦. 合肥: 合肥工業大學博士學位論文: 1–158.

毛景文, 邵擁軍, 謝桂青, 張建東, 陳毓川. 2009. 長江中下游成礦帶銅陵礦集區銅多金屬礦床模型. 礦床學報, 28(2): 109–119.

邱士東, 謝玉玲, 徐九華, 王葆華, 楊竹森, 蒙義峰. 2007. 安徽銅陵冬瓜山銅礦床成礦流體特征及演化. 礦床地質, 26(2): 204–212.

唐永成, 吳言昌, 儲國正, 邢鳳鳴, 王永敏, 曹奮揚, 常印佛. 1998. 安徽沿江地區銅金多金屬礦床地質. 北京: 地質出版社: 1–351.

吳淦國, 張達, 臧文拴. 2003. 銅陵礦集區構造滑脫與分層成礦特征研究. 中國科學(D輯), 33(4): 300–308.

肖凡, 王愷其. 2021. 德興斑巖銅礦床斷裂與侵入體產狀對成礦的控制作用: 從力?熱?流三場耦合數值模擬結果分析. 地學前緣, 28(3): 190–207.

肖新建, 顧連興, 倪培. 2002. 安徽銅陵獅子山銅?金礦床流體多次沸騰及其與成礦的關系. 中國科學(D輯), 32(3): 199–206.

徐文藝, 楊竹森, 蒙義峰, 曾普勝, 史大年, 田世洪, 李紅陽. 2004. 安徽銅陵礦集區塊狀硫化物礦床成因模型與成礦流體動力學遷移. 礦床地質, 23(3): 353–364.

徐曉春, 范子良, 何俊, 劉雪, 劉曉燕, 謝巧勤. 陸三明, 樓金偉. 2014. 安徽銅陵獅子山礦田銅金多金屬礦床的成礦模式. 巖石學報, 30(4): 1054–1074.

徐曉春, 陸三明, 謝巧勤, 柏林, 儲國正. 2008. 安徽銅陵獅子山礦田巖漿巖鋯石SHRIMP定年及其成因意義. 地質學報, 82(4): 500–509.

徐兆文, 陸現彩, 高庚, 方長泉, 王云健, 楊小男, 蔣少涌, 陳幫國. 2007. 銅陵冬瓜山層狀銅礦同位素地球化學及成礦機制研究. 地質論評, 53(1): 44–51.

於崇文, 蔣耀松, 肖正域. 1995. 安徽銅陵層控夕卡巖型銅礦床的成礦作用動力學. 地質學報, 69(3): 243–254.

袁峰, 張明明, 李曉暉, 葛粲, 陸三明, 李建設, 周宇章, 蘭學毅. 2019. 成礦預測: 從二維到三維. 巖石學報, 35(12): 3863–3874.

Annen C, Burgisser A. 2021. Modeling water exsolution from a growing and solidifying felsic magma body., 402, 105799.

Cao W, Liu L M, Liu H S, Lai F. 2020. Investigating the irregular localization of skarn orebodies by computational modeling in the Fenghuangshan ore field, Tongling district, Anhui Province, China., 29(5): 2967–2988.

Cao Y, Gao F P, Du Y S, Du Y L, Pang Z S. 2017. Genesis of the Datuanshan stratabound skarn Cu(-Mo) deposit, Middle- Lower Yangtze Valley, Eastern China—Constraints from geology, Re-Os geochronology, mineralogy, and sulfur isotopes., 52(3): 443–462.

Deng J, Wang Q F, Xiao C H, Yang L Q, Liu H, Gong Q J, Zhang J. 2011. Tectonic-magmatic-metallogenic system, Tonglingore cluster region, Anhui Province, China., 53(5–6): 449–476.

Eldursi K, Branquet Y, Guillou-Frottier L, Marcoux E. 2008. Numerical investigation of transient hydrothermal processes around intrusions: Heat-transfer and fluid circulation controlled mineralization patterns., 288(1–2): 70–83.

Fournier R O. 1999. Hydrothermal processes related to movement of fluid from plastic into brittle rock in the magmatic-epithermal environment., 94(8): 1193–1211.

Gu L X, Zaw K, Hu W X, Zhang K J, Ni P, He J X, Xu Y T, Lu J J, Lin C M. 2007. Distinctive features of Late Palaeozoic massive sulphide deposits in South China., 31(1–4): 107–138.

Hanson R B. 1996. Hydrodynamics of magmatic and meteoric fluids in the vicinity of granitic intrusions.:, 87(1–2): 251–259.

Hobbs B E, Zhang Y, Ord A, Zhao C. 2000. Application of coupled deformation, fluid flow, thermal and chemical modeling to predictive mineral exploration.,69–70: 505–509.

Ingebritsen S E, Geiger S, Hurwitz S, Driesner T. 2010. Numerical simulation of magmatic hydrothermal systems., 48(1), RG1002.

Itasca Consulting Group, Inc. 2012. FLAC3DUser’s Guide, version 5.0. Minneapolis.

Korges M, Weis P, Andersen C. 2020. The role of incremental magma chamber growth on ore formation in porphyry copper systems., 552, 116584.

Lamy-Chappuis B, Heinrich C A, Driesner T, Weis P. 2020. Mechanisms and patterns of magmatic fluid transport in cooling hydrous intrusions., 535, 116111.

Li X H, Yuan F, Zhang M M, Jowittd S M, Ord A, Zhou T F, Dai W Q. 2019. 3D computational simulation-based mineral prospectivity modeling for exploration for concealed Fe-Cu skarn-type ineralization within the Yueshan orefield, Anqing district, Anhui Province, China., 105: 1–17.

Li Z H, Chi G X, Bethune K M, Thomas D, Zaluski G. 2017. Structural controls on fluid flow during compressional reactivation of basement faults: Insights from numerical modeling for the formation of unconformity-related uranium deposits in the Athabasca Basin, Canada., 112(2): 451–466.

Liu L M, Sun T, Zhou R C. 2014. Epigenetic genesis and magmatic intrusion’s control on the Dongguashan strata bound Cu-Au deposit, Tongling, China: Evidence from field geology and numerical modeling., 144: 97–114.

Liu L M, Wan C L, Zhao C B, Zhao Y L. 2011. Geodynamic constraints on orebody localization in the Anqing orefield, China: Computational modeling and facilitating predictive exploration of deep deposits., 43: 249–263.

Liu L M, Zhao Y L, 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., 526–529: 110–123.

Liu L M, Zhao Y L, Zhao C B. 2010. Coupled geodynamics in the formation of Cu skarn deposits in the Tongling- Anqing district, China: Computational modeling and implications for exploration., 106(1–3): 146–155.

Liu X C, Xiao C H, Zhang S H, Chen B L. 2021. Numerical modeling of deformation at the Baiyun Gold Deposit, northeastern China: Insights into the structural controls onmineralization., 32(1): 174– 184.

Liu Z F, Shao Y J, Wang C, Liu Q Q. 2018a. Genesis of the Dongguashan skarn Cu-(Au) deposit in Tongling, EasternChina: Evidence from fluid inclusions and H-O-S-Pb isotopes., 104: 462–476.

Liu Z F, Shao Y J, Zhang Y, Wang C. 2018b. Geochemistry and geochronology of the Qingshanjiao granites: Implications for the genesis of the Dongguashan copper (gold) ore deposit in the Tongling ore district, Eastern China., 99: 42–57.

Norton D, Taylor H P. 1979. Quantitative simulation of the hydrothermal systems of crystallizing magmas on the basis of transport theory and oxygen isotope data: An analysis of the Skaergaard intrusion., 20(3): 421–486.

Pan Y M, Dong P. 1999. The Lower Changjiang (Yangzi/ Yangtze River) metallogenic belt, east central China: Intrusion- and wall rock-hosted Cu-Fe-Au, Mo, Zn, Pb, Ag deposits., 15(4): 177–242.

Price G P, Stoker P. 2002. Australian geodynamics cooperative research centre’s integrated research program delivers a new minerals exploration strategy for industry., 49(4): 595–600.

Quarteroni A, Valli A. 1994. Numerical Approximation of Partial Differential Equations. Berlin, Heidelberg: Springer: 1–537.

Sch?n J H. 1998. Physical properties of rocks. Fundamentals and principles of petrophysics. Oxford: Elsevier: 1–360.

Sch?pa A, Annen C, Dilles J H, Sparks R S J, Blundy J D. 2017. Magma emplacement rates and porphyry copper deposits: Thermal modeling of the Yerington batholith, Nevada., 112(7): 1653–1672.

Vigneresse J L, Truche L. 2020. Modeling ore generation in a magmatic context., 116, 103223.

Wang S W, Zhou T F, Yuan F, Fan Y, Zhang L J, Song Y L. 2015. Petrogenesis of Dongguashan skarn-porphyry Cu-Au deposit related intrusion in the Tongling district, eastern China: Geochronological, mineralogical, geochemical and Hf isotopic evidence., 64: 53–70.

Weis P, Driesner T, Heinrich C A. 2012. Porphyry-copper ore shells form at stable pressure-temperature fronts with dynamic fluid plumes., 338(6114): 1613–1616.

Xiao X, Zhou T F, White N C, Zhang L J, Fan Y, Chen X F. 2021. Multiple generations of titanites and their geochemicalcharacteristics record the magmatic-hydrothermal processes and timing of the Dongguashan porphyry-skarn Cu-Au system, Tongling district, Eastern China., 56(1): 363–380.

Xie J C, Ge L K, Qian L, Li Q Z, Sun W D. 2020. Trace element characteristics of pyrite in Dongguashan Cu(Au) deposit, Tongling region, China., 5(4): 233–246.

Xie J C, Yang X Y, Sun W D, Du J G. 2012. Early cretaceous dioritic rocks in the Tongling region, eastern China: Implications for the tectonic settings., 150: 49–61.

Xu T F, Spycher N, Sonnenthal E, Zhang G X, Zheng L G, Pruess K. 2011. TOUGHREACT Version 2.0: A simulator for subsurface reactive transport under non-isothermal multiphase flow conditions., 37(6): 763–774.

Xu X C, Zhang Z Z, Liu Q N, Lou J W, Xie Q Q, Chu P L, Frost R. 2011. Thermodynamic study of the association and separation of copper and gold in the Shizishan ore field, Tongling, Anhui Province, China., 43(1): 347–358.

Yang X N, Xu Z W, Lu X C, Jiang S Y, Ling H F. 2011. Porphyry and skarn Au-Cu deposits in the Shizishan ore field, Tongling, East China: U-Pb dating andHf isotope analysis of zircon and petrogenesis of associated granitoids., 43(1): 182–193.

Yapparova A, Miron G D, Kulik D A, Kosakowski G, Driesner T. 2019. An advanced reactive transport simulationscheme for hydrothermal systems modelling., 78: 138–153

Zhang J D, Liu L, Yu Z Q, Yang X Y, Liu Z F, Li H. 2021. Petrogenesis and metallogeny of the Dongguashan Cu-Au deposit in the Tongling ore-cluster region, the Lower Yangtze River Metallogenic Belt: Constraints from geochemistry and geochronology., 81(3), 125822.

Zhang Y, Karrech A, Schaubs P M, Regenauer-Lieb K, Poulet T, Cleverley J S. 2012. Modelling of deformation around magmatic intrusions with application to gold-related structures in the Yilgarn Craton, Western Australia., 526–529: 133–146.

Zhao C B, Hobbs B E, Alison A. 2009. Fundamentals of Computational Geoscience. Berlin: Springer: 1–241.

3D Computational Geodynamic Modeling of Skarn Copper (Gold) Mineralization System and its Implication for Ore Formation in the Shizishan Orefield, Tongling District, Eastern China

ZHOU Feihu1, LIU Liangming1*, LIU Hongsheng1, SHENG Huan2, WU Hongzhi1, MIN Lingshuai1, BI Chenxi1

(1.MOEKey Laboratory of Metallogenic Prediction of Nonferrous Metals and Geological Environment Monitoring; Computational Geosciences Research Centre, Central South University; Changsha 410083, Hunan, China; 2. Dongguashan Copper Mine, Tongling Nonferrous Metals Group Holding Co., Ltd., Tongling 244000, Anhui, China)

The Shizishan orefield hosts the largest Cu (Au) skarn deposit, the Dongguashan deposit, in the Tongling district, even in the whole Lower Yangtze Metallogenic Belt. The major orebody of the Dongguashan deposit is strataboud and stratiform Cu (Au) skarn which has more than 1 Mt Cu reserves, accounting about 60% of the total reserves of the orefield. The localization mechanism of this significantly large stratiform orebody is still actively debated because of the complexity of the metallogenic characteristics. For better understanding the dynamic mechanism and the controlling factors for ore formation, we constructed a 3D geological model of the orefield by integrating all available exploration data, and then carried 3D computational geodynamic modelling in the FLAC3Dplatform, to simulate the coupled mechano-thermo-hydrological (MTH) processes during the cooling of the intrusion and mineralization in the Shizishan orefield. The results of the computational geodynamic modelling show that, in the extensional tectonic regime, the shallow high-dilation zones (volumetric strain ≥1.2%) are not significantly large and located in the contact zone that hosts orebodies of the East and West Shizishan deposit; while the deep high-dilation zones, where the fluids from the intrusion and sedimentary rocks are focusing, are the largest and most intensive zones located in the carbonate rocks of the mid to upper Carboniferous that hosts the Dongguashan ore deposit. Because these fluid focusing high-dilation zones are spatially consistent with the known orebodies, we conclude that the fluid focusing high-dilation zones were produced during the intrusion cooling processes and controlled the ore localization in the Shizishan orefield. Furthermore, the scale of the fluid focusing high-dilation zones in the mid to upper Carboniferous formation is larger than that of the known Dongguashan deposit, indicating that the mid to upper Carboniferous formation are favorable for the Cu (Au) mineralization. For investigating ore formation-controlling factors, we conducted dynamic simulation experiments on 4 different scenarios with different boundary conditions and material parameters: scenario 1, under NE-SW compression; scenario 2, without tectonic-load on boundaries; scenario 3, increasing the shear modulus of the C2+3to 4.4×1010Pa; scenario 4, increasing the permeability of the P2and T to 25×10?12m2. The simulation experiments of these 4 scenarios show that the patterns of deformation and fluid flow are very different from that occurred in the former dynamic simulation experiment with the credible boundary conditions and material parameters supported by the geological investigation in the orefield. The most distinct differences occurred in the scenarios under compression and with increased permeability of P2and T. The results of these 4 scenarios, especially for the scenarios 1 and 4, are unfavorable for formation and localization of the stratiform skarn orebodies like the Dongguashan deposit. Based on overall dynamic modelling results, we conclude that the fluid focusing dilatant deformation spaces resulted from the coupled MTH processes during cooling of the intrusion controlled the location and scale of the skarn orebodies in the Shizhishan orefield, whereas the tectonic regime, mechanic, and hydrological parameters are important constraints on ore formation, the extentional tectonic load and the lower permeability of the overlying strata, in particular, are more favorable for ore formation.

metallogenic dynamics; computational modeling; skarn ore deposit; Shizishan orefield; Tongling

2021-11-14;

2022-03-16

國家自然科學基金項目(41772351)資助。

周飛虎(1998–), 男, 碩士研究生, 地質資源與地質工程專業。E-mail: fhzhou_csu@163.com

劉亮明(1964–), 男, 教授, 主要從事礦床地質和構造地質及其計算模擬研究。E-mail: lmliu@csu.edu.cn

P628; P612

A

1001-1552(2023)06-1242-014

10.16539/j.ddgzyckx.2023.06.004