铜陵狮子山矿田矽卡岩铜(金)成矿系统三维动力学计算模拟及其意义

周飞虎, 刘亮明*, 刘洪升, 盛 欢, 伍红志, 闵令帅, 毕晨曦

铜陵狮子山矿田矽卡岩铜(金)成矿系统三维动力学计算模拟及其意义

周飞虎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), 往深部主要是受层位控制的似层状和层状矿体, 而且矿体规模总体具有明显往深部变大的趋势(图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