论浮力对热液喷发型(SEDEX)矿床成矿的作用：以澳大利亚北部为例（杨建文,冯佐海,罗先熔,陈远荣《中国科学D辑》2009.5）.pdf

中国科学 D 辑 地球科学 2009 年 第 39 卷 第 5 期 594 601 中国科学杂志社 SCIENCE IN CHINA PRESS 引用格式引用格式 Yang J W, Feng Z H, Luo X R, et al. On the role of buoyancy force in the ore genesis of SEDEX deposits Example from Northern Australia. Sci China Ser D-Earth Sci, 2009, 524 452460, doi 10.1007/s11430-009-0043-x 论浮力对热液喷发型SEDEX矿床成矿的作用 以澳大利亚北部为例 杨建文 ①②*, 冯佐海①, 罗先熔①, 陈远荣① ① 桂林理工大学资源与环境工程系, 桂林 541004; ② Department of Earth and Environmental Sciences, University of Windsor, Ontario, N9B 3P4, Canada * E-mail jianwenyuwindsor.ca 收稿日期 2008-07-27; 接受日期 2008-10-08 广西人才小高地创新计划、加拿大国家自然科学基金编号 RGPIN261283和国家自然科学基金批准号 40772126资助 摘要 采用有限元法模拟了北澳大利亚 Mount Isa 盆地古流体迁移及温度场的分布, 特别评价了浮力对沉积喷发型 SEDEX 矿床成矿的作用. 数值模拟结果表明浮力驱动的 热液循环主要由断层深度及其与含水层的空间关系所控制. 海水从一条断层向下迁移, 然后沿着主要的含水层横向流动. 在此过程中地热导致温度逐渐增高, 最后在浮力作用 下加热后的富含金属元素的流体从另外一条断层向上运动并喷发到海底. 计算表明喷发 到海底的热液温度范围在 115160℃之间, 并在一百万年的时间内保持了 2.64.1 ma−1 的流速. 这些条件在适当的化学圈闭环境下有利于形成 Mount Isa 规模的 SEDEX 矿床. 因此, 浮力是盆地范围热液循环的一个重要驱动机制, 其强度之大足以形成类似于北澳 大利亚 Mount Isa 规模的超巨型热液喷发型矿床. 关键词 热液流体 有限元模拟 SEDEX 矿床 北澳 Mount Isa 盆地 热液喷发型SEDEX矿床是铅和锌的主要来源, 同时也是银的重要来源, 从全球来看, 它占到锌产量 的 40和铅产量的 60[1]. 地质及成矿研究表明此类 矿床可能是由上升到海底的富含金属卤水在热液喷 发过程中形成的. 尽管有此共识, 地质学家仍在继续 争论有关流体迁移、 热传导和物质运移的机理及对成 矿的作用. 通常认为地形、压实作用、浮力和构造形变可以 导致盆地范围的热液循环[2]. 其中, 浮力驱动具有非 常重要的意义, 因为它可以引起大范围的热液对流 和物质迁移. 其影响要比扩散作用大, 这是因为扩散 是由分子碰撞所引起的, 因而代表的是一种微观过 程. 在某些条件下, 浮力在促进高密度流体和低密度 流体的动态混合方面甚至比其他作用更为重要[3]. 流体的温度及盐度在空间和时间上的变化会导 致其密度发生变化, 从而产生浮力. 因为缺乏外部的 输入和输出, 所以浮力驱动的对流通常被称为自由 对流[4]. 此类流体的重要性在研究洋中脊的热液循环 系统[5,6]、 海水入侵[7,8]和溶质迁移[9,10]的过程中就已经 被认识到了. 近年来, 浮力驱动的自由对流模型正引 起越来越多的关注, 成为解释SEDEX矿床成因的最 具前景的水动力学机理之一[11]. 例如, Garven等[12]探 讨了自由对流在HYC矿床形成过程中可能的作用, 他们的数值结果揭示了一个由断层及含水层所控制 中国科学 D 辑 地球科学 2009 年 第 39 卷 第 5 期 的“补给-排出”流体模型. Yang等[13]模拟了随时间变 化的成矿流体在盆地中的迁移, 并提出了一系列圈 定沉积喷发型矿床的野外勘探准则. 近来, Yang等[14] 探讨了卤水盐度对成矿过程及热液迁移的影响. 本文的主要目的是评估浮力驱动的流体和热传 导是否是 SEDEX 矿床成矿的主要机制. 将针对北澳 大利亚的 Mount Isa 盆地建立一个古水动力系统理论 模型, 然后进行一些列数值模拟计算. 1 控制方程和有限元模拟 从数学上讲, 沉积盆地中的热液流体是由流体 连续性方程和热能守恒方程控制的, 这两个方程通 过Darcy定律耦合在一起. 此外, 热液流体也与状态 方程有关. 状态方程包括流体密度和黏度, 他们是温 度和压力的函数. 流体连续性方程在变密度系统中 可用等价纯水水头来表达[7]. 这些方程的详细推导可 参见文献[6, 15]. 如果多孔隙介质是不可形变的, 且 在液相和固相间存在热平衡, 那么流体连续性方程 和热能守恒方程可以表达如下 , rs hhh KKKKS h xxyyzzt ρ ⎛⎞∂∂∂∂∂∂∂ ⎛⎞⎛⎞ ⎜⎟⎜⎟⎜⎟ ∂∂∂∂∂∂∂ ⎝⎠⎝⎠ ⎝⎠ 1 和 mmmww TTT cq T xxyyzzx λλλρ ⎛⎞∂∂∂∂∂∂∂ ⎛⎞⎛⎞ − ⎜⎟⎜⎟⎜⎟ ∂∂∂∂∂∂∂ ⎝⎠⎝⎠ ⎝⎠ x wwywwz cq Tcq T yz ρρ ∂∂ −− ∂∂ [1 wwss T cc t ρ φρφ ∂ − ∂ ], 2 式中, qx , qy和 qz分别是 x, y, z方向上的 Darcy通量, K 是渗透率, φ是孔隙度, 1 mws φφ λλλ − λw和λs代表液 相和固相的热传导率, cw和 cs分别为液相和固相的 热容, ρw和ρs分别为液相和固相的密度, ρr为相对液 体密度, 定义为ρrρw/ρ0−1 ρ0是纯水参考密度, Ss为 储水率它依赖于流体密度和压缩率, 及岩石的压缩 率和孔隙度, t 为时间, h 为等价纯水水头, T 为温度. 这些方程构成了一个随时间变化的非线性耦合系统. 即使对于一个非常简单的一维地质问题, 它们的解 也非常棘手. 本文采用标准的Galerkin有限元法计算这些方程 的数值解[16,17]. 简言之, 从初始温度和初始流体条件 出发, 我们首先用前一步得到的最新温度数据来计 算新的流体密度和黏度; 然后解流体连续性方程求 得纯水水头分布; 再基于Darcy定律和新的流体参数 计算流体速度; 最后, 通过求解热能守恒方程来确定 温度分布. 上述这些步骤重复进行直到满足预设的 收敛准则. 我们采用非直角的四边形进行有限元网 格剖分, 这样可以更好地适应沉积盆地中不同岩石 单元和断层的几何形态. 每个网格单元都赋予包括 渗透率、孔隙度、热导率在内的各种控制流体和热传 导的物性参数. 具体计算方法及步骤可以参见文献[6, 18]. 2 理论模型 根据沉积盆地的基本特征, 特别是北澳大利亚 Mount Isa盆地古剖面的重建, 我们建立了如图 1 所示的二维古水动力学理论模型[19,20]. 它代表了在 Mount Isa盆地早期控制热液喷发型矿床成矿的高度 简化的地下层序和构造. 如图 1a所示, 该理论模型 涉及一个低渗透性的火山岩基底单元 4, 一个高渗 透性的砂岩含水层单元 3, 一个具有中等渗透性的 中间层单元 2, 以及一个由页岩和粉砂岩组成的顶 部容矿盖层单元 1. 两条渗透性更强的断层断层 1 和断层 2从顶部盖层一直深入基底. 根据地表露头 观测及地震剖面资料[21], 可判断断层向西陡倾, 断距 约为 1 km. 断层网格划分如图 1b所示. 类似的 Mount Isa-McArther盆地一般地质构造在先前的流体 模型中已经被采用过[1214,18,22,23]. 和先前的研究相似, 左边的断层1等价于Mount Isa断层体系, 右边的断层 2 模拟了一个远处北部的断裂带, 如在澳北世纪矿床 Century附近的Termite Range断层带. 理论模型尺寸 为 60 km20 km, 用 52 行和 19 列的二维非直角四边 形进行剖分如图 1b所示. 这两条向西陡倾的断层切断了砂岩含水层单元 3, 形成了一个有利于区域范围流体迁移的水动力框 架. 表 1 给出了地层单元和断层的渗透率和热传导率, 这些物性参数基于先前在Mount Isa-McArthur盆地所 进行的数值模拟研究[13,14,18,22,23]. 考虑到沉积岩的成 层特性, 我们假定其垂向渗透率比横向渗透率低两 个数量级. 另外, 我们假定cw4174 Jkg℃−1, 595 杨建文等 论浮力对热液喷发型SEDEX矿床成矿的作用 596 图图 1 a 基于北澳大利亚 Mount Isa 盆地建立的高度理想化的二维水动力学模型; b 用于有限元数值模拟的网格划分. 断层和岩层单元的物性参数见表 1 表表 1 断层和围岩的主要物性参数断层和围岩的主要物性参数 水文地质单元 横向渗透率 /m2 垂向渗透率 /m2 热传导率 /Jms℃−1 断层 410−14 410−14 2.0 单元 1 110−15 110−17 2.0 单元 2 210−15 210−17 2.5 单元 3 410−14 410−16 3.0 单元 4 210−16 210−18 3.0 cs 800 Jkg℃−1, λw 0.5 Wm℃−1, ρs2630 kgm−3, ρ01000 kgm−3, 且φ10. 根据先前的数值模拟研究, 断层之外的顶部边 界温度为 20℃, 处于海底, 流体可以通过. 在断层的 顶部, 垂向温度梯度固定为零, 这是因为流体不断从 其顶部高速通过, 因而靠近断层顶部流体基本等温. 底部边界温度为 450℃, 这由北澳大利亚的热流数据 所约束与 S. McLaren 及 M. Sandiford 私人通讯. 模 型的底部位于低渗透率的火山岩基底, 因而假定流 体不可渗透. 按照流体模拟惯例, 假定侧面边界是绝 热和非渗透的. 整个求解区域的流体初始速度设为 零, 初始温度随深度呈线形变化. 下面的数值模拟将 分别采用不同的断层深度和边界条件. 3 数值模拟结果 第一个数值模拟计算对应于上述边界条件和初 始条件, 模型几何形状和物性参数如图1和表1所示. 图 2 给出了温度和流体速度在不同时间段上的分布. 我们的数值结果表明冷的海水沿断层 2 向下流动, 然 后沿着透水性好的含水层横向迁移, 仅有少量的流 体渗透到透水性差的深部基底, 在此过程中流体被 中国科学 D 辑 地球科学 2009 年 第 39 卷 第 5 期 图图 2 数值模拟数值模拟 1对应于图对应于图 1 和表和表 1 所给出的参数所给出的参数 a 55000 年时温度等值线; b 100 万年时温度等值线, 此后热液系统达到稳定状态; c 稳态条件下的流体速度矢量 加热. 最后在浮力作用下, 温度抬高的流体沿着断层 1 上升并最终喷发到海底, 在那里形成类似于位于 Mount Isa 断层附近的 Mount Isa SEDEX 矿床. 另外, 在含水层内靠近边界的地方也形成了两个对流单元, 这可能是由假定的侧部边界不透水条件引起的. 对 比不同时间段的图2a和b可以发现, 当更多的海水 597 杨建文等 论浮力对热液喷发型SEDEX矿床成矿的作用 通过断层 2 渗入到盆地时, 整个系统冷却下来, 因而 从断层 1 喷发到海底的热液温度从 55000 年的 240℃ 降低至 100 万年的 160℃. 此后, 整个热液循环系 统达到稳定状态即, 温度和流体速度不再随时间变 化. 稳定状态下热液喷发到海底时的速度为 4.1 ma−1. 在第二个数值模拟计算中, 其他条件均与第一 个情况相同, 只是断层顶部被盖层封闭起来, 目的是 检验在低透水性岩层覆盖的情况下, 铅锌矿床能否 通过与黑色页岩的化学置换形成[24]. 如图 3 所示, 在 这种条件下, 没有明显的海水渗入盆地, 也没有热液 喷发到海底. 流体主要在断层及含水层内运移, 并形 成了 3 个 20 km规模的热液循环单元 一个在断层之 间, 两个在近边界处. 另外, 从图 3b还可以清楚地 看到流体在上覆的岩层下有明显的横向流动. 这说 明, 在适当的地球化学条件下比如, 当富含碳质的 粉砂岩或碳酸盐存在时, 此种流体机制有可能形成 层状的铅锌矿床. 但是, 现在的流体速度和温度要明 显低于第一种情况, 表明这种机制只能形成小规模 和低品位的矿体. 在第三个数值模拟计算中, 断层 1 的深度为 10 km仅到达含水层, 断层 2 的深度不变, 其他所有的 条件与第一个数值计算相同. 如图 4 所示, 数值计算 结果与第一种情况非常相似参见图 2b和c. 断层 1 作为流体排出通道, 断层 2 作为注入通道. 不同的 是, 热液喷发温度在稳态时为 115℃, 比第一种情况 低了45℃. 这显然是因为排出通道只有10 km深因而 不能抵达温度很高的火山岩基底. 在稳定状态下喷 图图 3 数值模拟数值模拟 2对应于断层被低渗透率顶部盖层封闭的情形对应于断层被低渗透率顶部盖层封闭的情形 a 稳态条件下温度等值线; b 稳态条件下流体速度矢量. 注意没有明显的热液喷发到海底 598 中国科学 D 辑 地球科学 2009 年 第 39 卷 第 5 期 图图 4 数值模拟数值模拟 3对应于断层对应于断层 1 仅穿透仅穿透 10 km 的情形的情形 a 稳态条件下温度等值线; b 稳态条件下流体速度矢量. 注意喷发热液的温度比数值模拟 1 中的低 发到海底的热液速度为 2.6 ma−1. 在最后一个数值模拟计算中, 我们假定两个断 层的深度都只有 10 km, 其他所有条件都与第一种情 况相同, 数值结果如图 5 所示. 可以看到, 与前面的 结果都不同, 现在形成了一个‘倒转’的流体充注-排 放模式. 断层 1 成为海水注入通道, 断层 2 成为热液 排出通道. 冷的海水沿断层 1 向下流动, 当到达含水 层时开始分散. 少量的海水在断层的左侧沿着含水 层向左迁移形成顺时针对流单元, 而大部分海水则 沿着透水性良好的含水层向右横向迁移. 同时从下 部边界获得的热能导致流体温度升高, 并驱动热液 沿断层 2 上升最后喷发到海底, 在那里有可能形成 SEDEX 型的矿床. 稳定状态时喷发到海底的热液温 度和速度分别为 142℃和 3.5 ma−1. 4 讨论和结论 本文利用有限元数值模拟方法探讨了浮力在沉 积喷发型 SEDEX 矿床成矿中的作用, 并以北澳大利 亚 Mount Isa 盆地为基础建立了高度理想化的理论模 型. 我们的数值模拟结果揭示断层和高渗透性碎屑 质的含水层之间的相互作用是控制盆地范围热液循 环的一个主要因素. 海水通过一条断层注入盆地, 然 后通过碎屑质的含水单元横向迁移, 同时由于地热 梯度而被加热并从围岩中滤出金属元素, 最后通过 另外一条断层向上流动并喷发到海底. 不同的断层 间隔和深度会形成不同的流体类型. 在前 3 个数值模 拟计算中, 断层 2 代表海水注入通道而断层 1 则是热 液喷发通道. 但当两个断层的深度都减少为 10 km 时, 则形成了‘倒转’的流体充注-排放模式. 导致这种变 599 杨建文等 论浮力对热液喷发型SEDEX矿床成矿的作用 图图 5 数值模拟数值模拟 4对应于两个断层均仅穿透对应于两个断层均仅穿透 10 km 的情形的情形 a 55000 年时温度等值线; b 稳态条件下流体速度矢量. 注意断层 1 为海水注入通道, 而断层 2 为热液喷发通道 化的原因还不完全清楚, 不过, 它们很可能受到盆地 结构及其与断层的空间关系所控制例如, Yang等[13]. 尽管已注意到一条断层为海水注入通道而另外一条 为热液喷发通道, 但本文并不打算深究是什么具体 原因导致了特定的流体迁移模式, 因为本文的重点 是探讨浮力在沉积喷发型SEDEX矿床成矿过程中的 潜在作用. 上述所有的数值模拟研究证实了浮力是盆地范 围热液迁移的一个有效的驱动机制. 我们需要回答 的关键问题是, 数值模拟所预测的流体速度及温度 是否有潜力形成一个类似于Mount Isa矿床规模的超 巨型矿床 我们在不同条件下的数值模拟结果表明 在一百万年的时间内, 热液喷发到海底时的速度在 2.64.1 ma−1之间. 因而, 单位长度断层的年平均流 体喷发通量约为 300 m3注意平均孔隙度为 10. 如 果Mount Isa矿床是由来自 2 km长的断层中的流体沉 积而成的, 且热液含有 100 ppm锌[25], 那么模拟的热 液喷发流体有能力在一百万年的时间内形成约 6 千 万吨规模的锌矿床假如热液流体中的锌能够全部沉 积成矿. 已知Mount Isa矿床的锌储量约为 1.57 千万 吨[1], 因而, 锌的有效沉积率只需 26, 即可形成该 矿床. 考虑到Mount Isa矿床的成矿温度[2123]及我们 数值模拟预测的结果, 上述第一、第三及第四种数值 模型似乎均能较好的解释Mount Isa类矿床的形成. 本文的数值模拟研究的最终结论如下 浮力是 一个盆地尺度范围流体迁移及热传导的重要的驱动 机制, 其强度之大足以形成类似于北澳大利亚 Mount Isa 规模的超巨型热液喷发型矿床. 600 中国科学 D 辑 地球科学 2009 年 第 39 卷 第 5 期 致谢 感谢两位匿名评审人提出的建设性意见. 参考文献 1 Goodfellow W D, Lydon J W, Turner R J W. Geology and genesis of strati sediment-hostedSedex zinc-lead-silver sulfide deposits. In Kirkham R V, Sinclair W D, Thorpe R I, et al, eds. Mineral Deposit Modelling. Geol Assoc Can Spec Publ, 1993, 40 201253 2 Garven G. Continental-scale groundwater flow and geologic processes. Annu Rev Earth Planet Sci, 1995, 24 89117[doi] 3 Raffensperger J P, Vlassopoulos D. The potential for free and mixed convection in sedimentary basins. Hydrol J, 1999, 7 505520 4 Sharp J M, Fenstemaker T R, Simmons C T, et al. Potential salinity-driven free convection in a shale-rich sedimentary basin example from the Gulf of Mexico basin in south Texas. Am Assoc Pet Geol Bull, 2001, 85 20892110 5 Fisher A T, Becker K. Correlation between seafloor heat flow and basement relief observational and numerical examples and implica- tions for upper crustal permeability. J Geophys Res, 1995, 100 1264112657[doi] 6 Yang J, Latychev K, Edwards R N. Numerical computation of hydrothermal fluid circulation in fractured earth structures. Geophys J Int, 1998, 135 627649 [doi] 7 Frind E O. Simulation of long-term transient density-dependent transport in groundwater. Adv Water Res, 1982, 5 7388 [doi] 8 Voss C I, Souza W R. Variable density flow and solute transport simulation of regional aquifers containing a narrow freshwater-salt- water transition zone. Water Resour Res, 1987, 23 18511866[doi] 9 Therrien R, Sudicky E A. Three-dimensional analysis of variably-saturated flow and solute transport in discretely-fractured porous media. J Contam Hydrol, 1996, 23 144[doi] 10 Shikaze S G, Sudicky E A, Schwartz F W. Density-dependent solute transport in discretely-fractured geological media is prediction pos- sible J Contam Hydrol, 1998, 34 273291 [doi] 11 Solomon M, Heinrich C A. Are high-heat-producing granites essential to the origin of giant lead-zinc deposits at Mount Isa and McArthur River, Australia Explor Min Geol, 1992, 1 8591 12 Garven G, Bull S W, Large R R. Hydrothermal fluid flow models of strati ore genesis in the McArthur Basin, Northern Territory, Australia. Geofluids, 2001, 1 289312[doi] 13 Yang J, Large R R, Bull S W. Factors controlling free thermal convection in faults in sedimentary basins implications for the ation of zinc-lead mineral deposits. Geofluids, 2004, 4 237247 [doi] 14 Yang J, Bull S W, Large R R. Numerical investigation of salinity in controlling ore-ing fluid transport in sedimentary basins exam- ple of the HYC deposit, Northern Australia. Miner Deposit, 2004, 39 622631[doi] 15 Bear J. Dynamics of Fluids in Porous Media. New York American Elsevier Publishing Company, 1972 16 Huyankorn P S, Pinder G F. Numerical s in Subsurface Flow. Toronto Academic Press, 1983 17 Kinzelbach W K H. Groundwater Modellingan Introduction with Sample Programs in BASIC. Amsterdam Elsevier Science Pub- lisher, 1986 18 Yang J. Finite element modeling of transient saline hydrothermal fluids in multi faulted sedimentary basins implications for ore-ing processes. Can J Earth Sci, 2006, 43 13311340[doi] 19 O’dea M G, Lister G S, MacCready T, et al. Geodynamic evolution of the Proterozoic Mount Isa terrain. Geol Soc Spec Pub, 1997, 121 99122[doi] 20 Betts P G, Giles D, Lister G S. Tectonic environment of shale-hosted massive sulfide Pb-Zn-Ag deposits of Proterozoic northeastern Australia. Econ Geol, 2003, 98 557576[doi] 21 Bierlein F P, Betts P G. The Proterozoic Mount Isa fault zone, northeastern Australia is it really a ca. 1.9 Ga terrane-bounding suture Earth Planet Sci Lett, 2004, 225 279294[doi] 22 Oliver N H S, McLellan J G, Hobbs B E, et al. Numerical models of extensional deation, heat transfer, and fluid flow across base- ment-cover interfaces during basin-related mineralization. Econ Geol, 2006, 101 131[doi] 23 McLellan J G, Oliver N H S, Hobbs B E. The relative effects of deation and thermal advection on fluid pathways in basin-related mineralization. J Geochem Explor, 2006, 89 271275 [doi] 24 Broadbent G C, Myers R E, Wright J V. Geology and origin of shale-hosted Zn-Pb-Ag mineralization at the Century deposit, northwest Queensland, Australia. Econ Geol, 1998, 93 12641294 25 Cooke D R, Bull S W, Large R R, et al. The importance of oxidized brines for the ation of Australian Proterozoic strati sedi- ment-hosted Pb-Zn Sedex deposits. Econ Geol, 2000, 95 118 601