Geodynamical simulation of the effects of ridge subduction on the scale of the seismogenic zone south of Chile Triple Junction
-
摘要: 智利三联点以北地震较多,智利三联点以南地震很少且热异常显著。为探究洋脊俯冲对由温度定义的孕震区的影响,以智利三联点区域的地质背景为基础建立二维有限元数值模型,对洋脊俯冲过程进行数值模拟,并将俯冲角度和汇聚速率等因素对孕震区的影响进行了对比。结果表明,洋脊俯冲过程中孕震区宽度减小,导致发生在智利三联点以南的地震远少于智利三联点以北。剖面附近的观测数据与数值模拟结果的对比表明,数值模拟可以大致反映智利三联点区域板块间的孕震区宽度和地表热流特征。当俯冲汇聚量相同时,板块间的汇聚速率越大,洋脊俯冲过程中孕震区则越宽且其下边界越深,海沟附近的地表热流越高。与汇聚速率相比,俯冲角度等因素对地表热流的影响较小。俯冲角度越大,洋脊俯冲过程中孕震区越窄。当数值模型包含剪切生热时,洋脊俯冲过程中孕震区宽度大约可减小至15 km,且孕震区很浅。这可导致智利三联点以南部分区域难以发生地震,出现观测不到和达-贝尼奥夫带的现象。Abstract: Earthquakes are much more frequent to the north of Chile Triple Junction than to the south, where the thermal anomaly is also more significant. To study the effects of ridge subduction on the thermally defined seismogenic zone, two-dimensional finite element models were established based on the geology of Chile Triple Junction, the process of ridge subduction was simulated, and the effects of the initial slab dip and the convergence rate on the seismogenic zone were compared. The results show that the width of the seismogenic zone decreases during the ridge subduction, inducing earthquakes to occur much less to the south than to the north of Chile Triple Junction. By comparing the observed data in the vicinity of the profiles with the numerical simulation results, we find that the numerical simulation can roughly reflect the width of the interplate seismogenic zone and the surface heat flow in the area of Chile Triple Junction. At the same convergence distance, a larger convergence rate comes with a wider seismogenic zone and the deeper the downdip limit of the seismogenic zone, the higher the surface heat flow in the vicinity of the trench. Compared with the convergence rate, factors such as slab dip have little effect on the surface heat flow. In the process of ridge subduction, larger the slab dips leads to narrower seismogenic zones. When the effect of shear heating is included in the simulation, the width of the seismogenic zone in the process of ridge subduction can shrink to about 15 km and the depth of the seismogenic zone is small. Such a narrow and shallow seismogenic zone makes it hard for earthquakes to occur and for the Wadati-Benioff plane to be observed in some areas south of Chile Triple Junction.
-
Keywords:
- Chile Triple Junction /
- ridge subduction /
- seismogenic zone /
- numerical simulation
-
引言
中国大陆中部102°E与106°E之间突显一条纵贯南北的构造带,从滇西南经四川和甘肃,直至贺兰山,被称为南北构造带。张文佑(1959)指出,以南北构造带为界,中国东部的盆山构造体系受NNE和NNW向断裂控制,而西部受ENE和WNW向断裂控制。马杏垣(1989)在新构造图中明确标示了该大型构造带。南北构造带还是中国大陆浅源地震最活跃、地震分布带状明显的巨型地震带,被称为南北地震带。可见,对南北地震带的研究具有重要的科学意义。
Zhang (2013)指出南北地震带北起宁夏贺兰山,跨越西秦岭,穿过龙门山、小江、红河等断裂带向南延伸至缅甸境内,成为分隔中国大陆东部相对稳定的鄂尔多斯地块、四川盆地、华南地块与西部地震活动性强烈的青藏高原之间的重要边界带。长期以来,南北地震带区域范围的划定方法并没有严格的、统一的规定,按照不同的观点和划分依据,对其展布范围的认识也不尽相同。
深部地球物理探测表明,南北地震带为非常明显的重力梯度带和地壳厚度梯度带(马杏垣,1989),同时也是岩石圈速度结构、地壳密度、深部电性结构等横向显著变化的地带(王椿镛等,2015)。GPS测量得到的中国大陆速度场结果(Wang et al,2001)显示,我国东部和西部地壳运动速度场存在明显的差别,南北地震带位于西部与东部之间的过渡区,具有复杂的变化形态。在不同地球动力学过程的作用下,晚新生代和现今构造变形在南北地震带两侧的差异非常明显,因此南北地震带也是中国大陆地应力场一级分区的边界(谢富仁等,2004)。
地壳的应力状态是地壳最重要的性质之一。地壳表面和内部发生的各种构造现象(包括断层的形成与运动)及其伴生的各种物理化学现象均与地壳应力作用密切相关。区域构造应力场反映了地壳应力的分布状态,是揭示地球动力过程、研究断层相互作用的重要内容。目前,国内外针对构造应力场的研究方法种类众多,如地貌学及测量学方法、浅部地应力测量方法、地质构造方法和地球物理方法等。近几十年来,许多研究人员使用不同的方法对中国大陆构造应力场进行了大量研究(鄢家全等,1979;许忠淮等,1989;谢富仁等,2004;徐纪人等,2008;Wan,2010)。其中,利用震源机制资料研究区域构造应力场是最为重要的手段之一。针对南北地震带及其周缘的震源机制解及应力场研究已取得了不少成果,例如:崔效锋和谢富仁(1999)将滑动方向拟合法引入应力分区的研究中,使用震源机制解资料把中国西南及邻区分为5个应力区;Zhao等(2013)使用2008—2011年云南地区的174次中小地震的震源机制解及全球矩心矩张量(global centroid moment tensor,简写为GCMT)结果,运用阻尼线性反演方法,得到该地区的应力场结果;王晓山等(2015)利用震源机制解数据研究了南北地震带2°×2°网格大小的构造应力场特征;Xu等(2016)利用云南地区686个震源机制解资料,使用阻尼线性反演方法,给出了该区域的应力场特征;Luo等(2016)利用CAP (cut and paste)方法获取了川滇地区2007—2014年173次中小地震的震源机制解,并补充使用1976—2007年132次GCMT资料,给出了青藏高原东南缘P,T轴分布及构造应力场特征。但受数据资料的数量、空间分布等条件的限制,很多研究得到的区域构造应力场的空间尺度较大。
南北地震带频繁的地震活动和特殊的构造环境使其成为研究现今区域构造应力场的理想场所。本研究拟使用南北地震带及邻近区域2009年1月—2017年8月的466次ML≥3.5地震的震源机制解,并补充GCMT数据库公布的1976年1月—2017年8月发生在该区域的259次M≥4.5地震的震源机制解数据,基于区域应力张量阻尼反演方法(Hardebeck,Michael,2006),获得南北地震带的区域构造应力场,为该区域的动力学特征研究提供依据。
1. 研究方法
Michael (1984)基于地震断层面上的剪切应力方向与断层的滑动矢量方向一致且应力场的方向在一定区域内是一致的基本假设,构建单个区域震源机制解断层参数与应力张量参数的线性方程,提出一种利用一组震源机制解最佳拟合单个应力张量的测定方法。Hardebeck和Michael (2006)在该方法的基础上,通过引入一组可调整的阻尼参数抑制相邻子单元之间应力模式的突变,提出了区域应力张量阻尼反演方法。该方法对所有子区域的应力张量同时进行反演,很好地弥补了以往利用震源机制解反演区域构造应力场依赖于人为分区所带来的不足,从而能够更为准确地得到区域范围内的应力变化的特点。该方法的基本原理为:通过最小化每个断层面上滑动矢量与剪切应力矢量的差异,建立单个区域震源机制解断层面参数与应力张量参数的线性方程(Michael,1984)
$ {{G}}{{m}} {\text{=}} {{d}}{\text{,}} $
(1) 式中,矩阵m为应力张量分量的模型矢量,有6个独立分量。对于一组地震震源机制解,当应力张量不存在各向同性分量时,应力张量的迹为0,即σ33=-(σ11+σ22),此时模型矢量m具有5个独立参量,即
$ {{m}} {\text{=}}\left(\!\!{\begin{array}{*{20}{c}} {{\sigma _{11}}} \\ {{\sigma _{12}}} \\ {{\sigma _{13}}} \\ \begin{gathered} {\sigma _{22}} \\ {\sigma _{23}} \\ \end{gathered} \end{array}} \!\!\right){\text{,}} $
(2) 由模型矢量m可以确定应力张量的3个主应力方向以及描述主应力相对大小的量R。
式(1)中矩阵d为数据矢量,由K个地震震源机制解的3K个滑动矢量分量组成,可以表示为
$ {{d}} {\text{=}}\left(\!\!{\begin{array}{*{20}{c}} {{s_{11}}} \\ {{s_{12}}} \\ {{s_{13}}} \\ \begin{gathered} \; \vdots \\ {s_{K1}} \\ {s_{K2}} \\ {s_{K3}} \\ \end{gathered} \end{array}} \!\!\right){\text{,}} $
(3) This page contains the following errors:
error on line 1 at column 1: Start tag expected, '<' not foundBelow is a rendering of the page up to the first error.
数据核矩阵G由每个地震震源机制解的断层面法向矢量得到,其表达式为
$ {{G}} {\text{=}}\left(\!\!\!\!{\begin{array}{*{20}{c}} {{n_{11}} {\text{-}} n_{11}^3 {\text{+}} {n_{11}}n_{13}^2}\!\! & \!\!{{n_{12}} {\text{-}} 2n_{11}^2{n_{12}}}\!\! & \!\!{{n_{13}} {\text{-}} 2n_{11}^2{n_{13}}}\!\! & \!\!{{n_{11}}n_{13}^2 {\text{-}} {n_{11}}n_{12}^2}\!\! & \!\!{ {\text{-}} 2{n_{11}}{n_{12}}{n_{13}}} \\ {{n_{12}}n_{13}^2 {\text{-}} {n_{12}}n_{11}^2}\!\! & \!\!{{n_{11}} {\text{-}} 2n_{12}^2{n_{11}}}\!\! & \!\!{ {\text{-}} 2{n_{11}}{n_{12}}{n_{13}}}\!\! & \!\!{{n_{12}} {\text{-}} n_{12}^3 {\text{+}} {n_{12}}n_{13}^2}\!\! & \!\!{{n_{13}} {\text{-}} 2n_{12}^2{n_{13}}} \\ {n_{13}^3 {\text{-}}{n_{13}}n_{11}^2 {\text{-}} {n_{13}}}\!\! & \!\!{ {\text{-}} 2{n_{11}}{n_{12}}{n_{13}}}\!\! & \!\!{{n_{11}}{\text{-}}2n_{13}^2{n_{11}}}\!\! & \!\!{n_{13}^3 {\text{-}} {n_{13}}n_{12}^2 {\text{-}} {n_{13}}}\!\! & \!\!{{n_{12}} {\text{-}} 2n_{13}^2{n_{12}}} \\ \vdots \!\! & \!\! \vdots \!\! & \!\! \vdots \!\! & \!\!\vdots \!\! & \!\! \vdots \\ {{n_{K1}} {\text{-}} n_{K1}^3 {\text{+}} {n_{K1}}n_{K3}^2}\!\! & \!\!{{n_{K2}} {\text{-}} 2n_{K1}^2{n_{K2}}}\!\! & \!\!{{n_{K3}} {\text{-}} 2n_{K1}^2{n_{K3}}}\!\! & \!\!{{n_{K1}}n_{K3}^2 {\text{-}} {n_{K1}}n_{K2}^2}\!\! & \!\!{ {\text{-}} 2{n_{K1}}{n_{K2}}{n_{K3}}} \\ {{n_{K2}}n_{K3}^2 {\text{-}}{n_{K2}}n_{K1}^2}\!\! & \!\!{{n_{K1}} {\text{-}} 2n_{K2}^2{n_{K1}}}\!\! & \!\!{ {\text{-}} 2{n_{K1}}{n_{K2}}{n_{K3}}}\!\! & \!\!{{n_{K2}} {\text{-}}n_{K2}^3 {\text{+}} {n_{K2}}n_{K3}^2}\!\! & \!\!{{n_{K3}} {\text{-}} 2n_{K2}^2{n_{K3}}} \\ {n_{K3}^3 {\text{-}} {n_{K3}}n_{K1}^2{\text{-}}{n_{K3}}}\!\! & \!\!{ {\text{-}} 2{n_{K1}}{n_{K2}}{n_{K3}}}\!\! & \!\!{{n_{K1}} {\text{-}} 2n_{K3}^2{n_{K1}}}\!\! & \!\!{n_{K3}^3 {\text{-}} {n_{K3}}n_{K2}^2 {\text{-}} {n_{K3}}}\!\! & \!\!{{n_{K2}} {\text{-}} 2n_{K3}^2{n_{K2}}} \end{array}} \!\!\!\!\right) {\text{,}} \!\!\!\!\!\!\! $
(4) This page contains the following errors:
error on line 1 at column 1: Start tag expected, '<' not foundBelow is a rendering of the page up to the first error.
式(1)的解可用最小二乘法反演给出,即
$ {{{G}}^{\rm{T}}}{{Gm}}{\text{=}} {{{G}}^{\rm{T}}}{{d}}{\text{,}} $
(5) 通过高斯消元法来求解(Michael,1984),反演的不确定性可用Bootstrap重采样技术来估计。
Hardebeck和Michael (2006)将这种方法扩展到区域应力场反演,将研究区划分为I×J个二维空间网格,第i行第j列网格点上的模型矢量由mij表示,数据矢量由dij表示,核矩阵记为Gij。
构建长度为M=5I×J的新的模型矢量mall
$ {{{m}}_{{\rm{all}}}}{{ {\text{=}} }}\left(\!\!{\begin{array}{*{20}{c}} {{{{m}}_{11}}} \\ {{{{m}}_{12}}} \\ \vdots \\ {{{{m}}_{21}}} \\ {{{{m}}_{22}}} \\ \vdots \\ {{{{m}}_{IJ}}} \end{array}} \!\!\right){\text{,}} $
(6) 长度为N的新的数据矢量dall为
$ {{{d}}_{{\rm{all}}}}{{ {\text{=}} }}\left(\!\!{\begin{array}{*{20}{c}} {{{{d}}_{11}}} \\ {{{{d}}_{12}}} \\ \vdots \\ {{{{d}}_{21}}} \\ {{{{d}}_{22}}} \\ \vdots \\ {{{{d}}_{IJ}}} \end{array}} \!\!\right){\text{,}} $
(7) This page contains the following errors:
error on line 1 at column 1: Start tag expected, '<' not foundBelow is a rendering of the page up to the first error.
$ {{{G}}_{{\rm{all}}}}{{ {\text{=}} }}\left({\begin{array}{*{20}{c}} {{{{G}}_{{{11}}}}}&{}&{}&{}&{}&{}&{} \\ {}&{{{{G}}_{{{12}}}}}&{}&{}&{}&{}&{} \\ {}&{}& \ddots &{}&{}&{}&{} \\ {}&{}&{}&{{{{G}}_{{{21}}}}}&{}&{}&{} \\ {}&{}&{}&{}&{{{{G}}_{{{22}}}}}&{}&{} \\ {}&{}&{}&{}&{}& \ddots &{} \\ {}&{}&{}&{}&{}&{}&{{{{G}}_{IJ}}} \end{array}} \right){\text{,}} $
(8) 则在区域应力场反演中,所有网格点的数据矢量与模型矢量的关系可用线性方程描述为
$ {{{G}}_{{\rm{all}}}}{{{m}}_{{\rm{all}}}} {\text{=}} {{{d}}_{{\rm{all}}}}{\text{,}} $
(9) 用最小二乘法求解,可以得到每个网格点的应力张量。
为了最小化相邻(上下左右相邻)各网格点的应力张量的差异,通过引入阻尼系数矩阵实现相邻格网间的“平滑”约束,阻尼系数矩阵D可表示为
$ {{D}} {\text{=}} \left({\begin{array}{*{20}{c}} {{I}}&{ {\text{-}} {{I}}}& \cdots &{{0}} &{{0}}& \cdots \\ {{I}}&{{0}}& \cdots &{ {\text{-}} {{I}}}&{{0}}& \cdots \\ {{0}}&{{I}}& \cdots & {{0}} &{ {\text{-}} {I}}& \cdots \\ {\vdots}&{\vdots}& {} &{\vdots}&{\vdots}&{} \\ {{0}}& {{0}}& \cdots &{{I}}&{ {\text{-}} {{I}}}& \cdots\\ {\vdots}&{\vdots}& {} &{\vdots}&{\vdots}&{} \end{array}} \right){\text{,}} $
(10) 式中,I为5×5阶单位矩阵。模型长度可用二阶范数表示为
$ {\left. {\left\| {{{D}}{{{m}}_{{\rm{all}}}}} \right.} \right\|_2}{\text{,}} $
(11) 其可表示每对相邻格点应力张量之间的差异度。数据与模型理论值的残差用二阶范数表示为
$ {\left. {\left\| {{{{G}}_{{\rm{all}}}}{{{m}}_{{\rm{all}}}}{\text{-}} {{{d}}_{{\rm{all}}}}} \right.} \right\|_2}{\text{,}} $
(12) 则引入阻尼系数矩阵D后,式(9)的阻尼最小二乘问题的解,需要同时使得模型长度与数据残差的加权和达到最小,即
$ \left. {\left\{ {{{G}}_{\rm{all}}^ {\rm{T}}{{{G}}_{\rm{all}}} {\text{+}} {e^2}{{{D}}^ {\rm{T}}}{{D}}} \right.} \right\}{{{m}}_{{\rm{all}}}} {\text{=}} {{G}}_{\rm{all}}^ {\rm{T}}{{{d}}_{{\rm{all}}}}{\text{,}} $
(13) This page contains the following errors:
error on line 1 at column 1: Start tag expected, '<' not foundBelow is a rendering of the page up to the first error.
2. 数据
选择中小地震活动频繁的南北地震带及邻近区域(21°N—42°N,97°E—107°E)作为研究区域,搜集整理该区域2009年1月—2017年8月发生的ML≥3.5地震的波形数据(郑秀芬等,2009),利用CAP 方法(Zhao,Helmberger,1994;Zhu,Helmberger,1996),基于足够多的信噪比较高且分布均匀的台站数据,得到466个波形拟合情况较好的地震震源机制解。在此基础上,补充1976年1月—2017年8月研究区内259次M≥4.5地震的GCMT数据,则共计725个地震震源机制解作为本研究参与应力场反演的数据。参照Zoback (1992)的震源断层分类原则,将震源机制类型分为6种:走滑型(SS)、逆断型(TF)、逆走滑型(TS)、正断型(NF)、正走滑型(NS)和无法确定型(U)。对这725个震源机制解进行分类,绘制出本研究所用到的震源机制解空间分布图(图1)。鉴于震源机制解的主压应力(P)轴和主张应力(T)轴的方向与区域构造应力场的方向有关,将725个震源机制解的P,T轴分布显示在图2中,其中线段的长短与P,T轴仰角的大小有关,线段越长,仰角越小。
将研究区域(21°N—42°N,97°E—107°E)按照1.0°×1.0°大小划分网格,逐一整理每一网格内的地震震源机制解数据,形成参与本研究区域构造应力场反演的输入数据。
3. 结果与分析
在应用区域应力张量阻尼方法对式(13)进行反演计算时,阻尼系数e越小,数据方差在反演中占有的相对权重越大,但模型长度将变得复杂;阻尼系数e越大,模型长度在反演中占有的相对权重越大,各子区域的应力张量越平滑,但数据方差值将变大。通过数据方差和模型长度均衡曲线选择最佳阻尼系数e,使得二者之间达到平衡,该值即为均衡曲线拐点处的值。本研究最佳阻尼系数取为e=1.0,如图3所示。选取最佳阻尼系数后,使用时空应力反演(spatial and temporal stress inversion,简写为SATSI)算法程序进行区域应力张量反演,得到研究区内1.0°×1.0°网格大小的区域应力场空间分布图(图4)。
如图4所示,南北地震带及邻近区域最大主应力由北至南变化显著,大体可以分为3段:在29°N—34°N范围内,除龙门山断裂带NE段最大主应力取向为NW向外,其它地段最大主应力方向近EW向;34°N以北最大主应力方向逐渐由近EW向旋转至NE向;而29°N以南最大主应力方向逐渐由近EW向变化至NW或NE向,再至近NS向。整体而言,南北地震带及邻近区域最大主应力方向由北到南发生了顺时针转动。这种顺时针旋转可理解为:在印度板块向NNE方向楔入青藏高原和高原重力势能东向推挤的共同作用下,青藏高原东北缘受到阿拉善地块、蒙古地块及西伯利亚地块的阻挡,导致地壳向SE方向挤出,当东部受到华南地块的阻挡后,转变为向南运动,这使得川滇地区围绕喜马拉雅东构造结顺时针转动,经滇西南向印度洋方向扩展(王阎昭等,2008)。此外,各地块内部最大主应力方向基本一致,不同地块内最大主应力方向不同,相邻地块间应力方向转换带与块体边界对应较好。例如,川滇菱形地块内部的最大主应力为NNW向,应力方向转换带的东边界与川滇菱形地块的东边界基本一致,其东部最大主应力方向为NW向。南北地震带及邻近区域中南段东、西两部分的最大主应力方向特征差异显著,大体以金沙江断裂、红河断裂为界,西部为NNE向,东部为NNW向,转换带与川滇菱形块体西边界基本一致。徐纪人等(2008)通过对中国大陆应力场的研究认为,南北地震带南段西部和东部的应力场方向特征分界线是印度板块与菲律宾板块在大陆内部应力控制区的分界线。从最大主应力的仰角看,南北地震带及邻近区域绝大多数地区近水平,川滇菱形地块西部至西北部的仰角较陡,反映出该区域存在一定的拉伸机制。
4. 讨论与结论
由上所述,所谓的“南北地震带”,不同地段区域构造应力场最大主应力σ1的取向差别显著,与已有其它研究(鄢家全等,1979;许忠淮等,1987;谢富仁等,2004;徐纪人等,2008;Wan,2010;郝平等,2012;王晓山等,2015)所得到的应力场结果具有很好的一致性。由图4可见,南北地震带不同地段最大主应力方向的差别突出地表现在:秦岭以北,总体上呈NE取向;川滇菱形地块的东边界,即安宁河、则木河断裂和小江断裂,总体上呈现为NW取向;中间地段,即龙门山断裂带及周边,σ1的取向较为复杂,该复杂性已有不少研究人员从汶川地震余震序列震源机制、应力特征等方面给予了分析研究(Wang et al,2009;Zheng et al,2009,2010)。这里要特别指出的是,秦岭以南和以北地区,区域构造应力场的方向明显有别,其原因仍有待深入研究。长期以来,在地震预测等地震科学研究中许多人把我国大陆地质构造形态(断裂带、山脉展布等)、构造应力场和地震活动特征归于太平洋板块朝亚欧板块俯冲及印度板块与亚欧板块碰撞的共同作用。并认为我国大陆东部地区主要受太平洋板块俯冲作用的影响,而西部地区则主要受控于印度板块与亚欧板块的碰撞作用。但也有些研究(陈章立,李志雄,2013)强调仅依据印度板块与亚欧板块碰撞的远程效应难以对昆仑山—秦岭以北的西北地区的地质构造、区域构造应力场和地震活动特征作出合理的解释,认为由于北太平洋板块朝亚欧板块东北部的俯冲以及大西洋海底扩张对亚欧板块的推挤作用,可能使近刚性的西伯利亚块体向南运动,昆仑山—秦岭以北的西北地区受西伯利亚块体向南运动的影响较大。我们认为上述动力学问题值得进一步深入研究。这里要着重强调的是,仅就秦岭以南和以北的区域构造应力场的显著差异来说,长期以来俗称的所谓“南北地震带”,难以将其作为一个统一的地震带应用于中长期地震预测的研究与实践中。
综上所述,本研究通过使用震源机制解数据,利用区域应力张量阻尼方法反演获得了南北地震带及其邻近区域构造应力场特征,主要结论如下:
1) 南北地震带作为青藏高原的东边界,由于其动力环境和力学特征复杂,其内部最大主应力方向具有明显的空间差异性。南北地震带北段最大主应力方向为NE向;南北地震带中段及周边除龙门山断裂带NE段最大主应力为NW-NNW向外,其它地段最大主应力近EW向;南北地震带中南段最大主应力方向逐渐由近EW向变化至NW或NE向,再至近NS向。整体而言,南北地震带及邻近区域最大主应力方向由北到南发生了顺时针旋转。
2) 川滇菱形地块内部的最大主应力方向为NNW向,应力方向转换带与地块边界基本一致,其东边界以东最大主应力方向为NW向,西边界以西为NNE向。
3) 从区域构造应力场的角度分析,难以将“南北地震带”作为一个统一的地震带应用于中长期地震预测的研究与实践中。
本研究所用地震波形数据由中国地震局地球物理研究所“国家数字测震台网数据备份中心”提供,审稿人为本文的修改提出了宝贵的意见与建议,作者在此一并表示衷心的感谢。
-
图 1 智利三联点区域构造背景及观测数据
(a) 智利三联点区域1906—2021年地震分布,数据来自国际地震中心(Bondár,Storchak,2011;Storchak et al,2017,2020);(b) 智利三联点区域热流数据分布(Lucazeau,2019)
Figure 1. Geotectonic background of the Chile Triple Junction and the observed data
(a) Earthquake distribution in the area of CTJ from 1906 to 2021,the seismicity data are from the International SeismologicalCenter (Bondár,Storchak,2011;Storchak et al,2017,2020);(b) The heat flow data in the area of CTJ (Lucazeau,2019)
图 2 参考模型初始设置及边界条件
这些参数由干橄榄岩(Hirth,Kohlstedt,2004)的位错蠕变和扩散蠕变的流变实验获得
Figure 2. Initial setups and boundary conditions of the reference model
All these parameters come from the rheological experiments on dislocation creep and diffusion creep of dry olivine (Hirth,Kohlstedt,2004)
图 7 海沟附近地表热流分布
(a) 剖面1各模型热流结果;(b) 剖面1不同初始俯冲角度模型热流结果(汇聚速率=4.2 cm/a);(c) 剖面2不同初始俯冲角度模型热流结果(汇聚速率=1.8 cm/a)
Figure 7. Surface heat flow in the vicinity of the trench
(a) Heat flow results of different models of profile 1;(b) Heat flow results of models with different initial slab dips of profile 1 (convergence rate=4.2 cm/a);(c) Heat flow results of models with different initial slab dips of profile 2 (convergence rate=1.8 cm/a)
图 8 剖面1(左)和剖面2(右)的P波速度扰动(Simmons et al,2010)
Figure 8. Vertical cross sections of P wave velocity perturbations (Simmons et al,2010) along the profile 1 (left) and the profile 2 (right)
表 1 数值模型的相关参数
Table 1 Relevant parameters of numerical model
物质名称 厚度
/km热扩散率κ
/10−6 (m2·s−1)比热容Cp
/(J·kg−1· K−1)密度ρ
/(kg·m−3)热膨胀系数
α/10−5 K−1内摩擦角
φ/°黏聚力
C/106 Pa生热率A
/(W·m−3)流变
准则①晶粒大小
指数m地幔 660 0.989 750 3 300 2 20 20 干橄榄岩 3 大陆上地壳 20 1.21 750 2 800 2 20 20 10−6 湿石英 1 大陆下地壳 20 1.15 750 2 900 2 20 20 4×10−7 湿钙长石 1 岩石层地幔 0—87 9.87 750 3 300 2 20 20 干橄榄岩 3 沉积层 4 1.21 750 3 000 2 5 10 辉长岩 1 大洋地壳 8 1.15 750 3 100 2 20 10 辉长岩 1 薄弱带 1.21 750 3 300 2 0.03 1 辉长岩 1 物质名称 前因子② (Pa−n·m−p·s−1) 应力指数n 活化能E/103 (J·mol−1) 活化体积V/10−6 (m3·mol−1) 扩散蠕变 位错蠕变 扩散蠕变 位错蠕变 扩散蠕变 位错蠕变 扩散蠕变 位错蠕变 地幔 2.37×10−15 6.52×10−16 1 3.5 375 530 4 13 大陆上地壳 1×10−50 8.57×10−28 1 4.0 0 223 0 0 大陆下地壳 1×10−50 7.13×10−18 1 3.0 0 345 0 0 岩石层地幔 2.37×10−15 6.52×10−16 1 3.5 375 530 4 18 沉积层 1×10−50 1.12×10−10 1 3.4 0 497 0 0 大洋地壳 1×10−50 1.12×10−10 1 3.4 0 497 0 0 薄弱带 1×10−50 1.12×10−10 1 3.4 0 497 0 0 注:① 流变准则:干橄榄岩来自Hirth和Kohlstedt (2004);湿石英来自Gleason和Tullis (1995);湿钙长石来自Rybacki (2006);辉长岩来自Wilks和Carter (1990)。② 黏滞系数的前因子来自平面应变单轴试验结果(Naliboff,Buiter,2015;刘梦雪等,2019). 表 2 参考模型孕震区宽度演化
Table 2 Evolution of the width of the seismogenic zone of the reference model
时间/Ma 闭锁区宽度/km
x(350℃)-x(100℃)过渡带宽度/km
x(450℃)-x(350℃)孕震区宽度/km
x(350℃)-x(100℃) +0.5× [ x(450℃) -x(350℃) ]洋脊与海沟横
坐标差/km0 86 36 104 −164 3 70 31 86 −60 6 89 34 106 −13 10 68 44 90 32 15 52 39 72 78 19.5 49 34 66 112 表 3 不同初始俯冲角度模型孕震区宽度演化
Table 3 Evolution of the width of the seismogenic zone of models of different initial slab dips
时间/Ma 孕震区宽度/km 洋脊与海沟横坐标差/km 15° 30° 45° 15° 30° 45° 0 111 104 66 −145 −164 −158 3 96 86 76 −67 −60 −36 6 85 106 74 −33 −13 8 10 94 90 60 14 32 36 15 90 72 51 46 78 67 19.5 80 66 52 77 112 92 表 4 部分模型孕震区深度演化
Table 4 Evolution of the depth of the seismogenic zone for some models
时间/Ma 参考模型 汇聚速率=6.6 cm/a 俯冲角度=45° 剪切生热 孕震区上边界
深度/km孕震区下边界
深度/km孕震区上边界
深度/km孕震区下边界
深度/km孕震区上边界
深度/km孕震区下边界
深度/km孕震区上边界
深度/km孕震区下边界
深度/km0 7 47 7 47 8 39 8 43 3 23 70 23 70 12 74 24 67 6 14 74 18 94 7 71 11 73 10 17 71 17 91 10 60 11 66 15 17 69 − − 5 51 9 51 19.5 19 70 − − 1 47 19 26 -
贾鸿瑞,魏东平. 2021. 智利三联点相关的板块相对运动及其地球动力学意义[J]. 地球物理学报,64(10):3567–3575. doi: 10.6038/cjg2021O0470 Jia H R,Wei D P. 2021. Relativemotion of plates related to the Chile triple junction and geodynamic significance[J]. Chinese Journal of Geophysics,64(10):3567–3575 (in Chinese).
李忠海,石耀霖. 2016. 三维板块几何形态对大陆深俯冲动力学的制约[J]. 地球物理学报,59(8):2806–2817. doi: 10.6038/cjg20160808 Li Z H,Shi Y L. 2016. Constraints of 3-D plate geometry on the dynamics of continental deep subduction[J]. Chinese Journal of Geophysics,59(8):2806–2817 (in Chinese).
刘梦雪,魏东平,史亚男. 2019. 俯冲初始时板块分界面形状对俯冲过程的影响[J]. 地球物理学报,62(1):78–87. doi: 10.6038/cjg2019L0717 Liu M X,Wei D P,Shi Y N. 2019. Effect of plate interface geometry on the evolution of subduction[J]. Chinese Journal of Geophysics,62(1):78–87 (in Chinese).
史亚男,魏东平,皇甫鹏鹏,李忠海,刘梦雪. 2019. 海洋板块俯冲作用下上覆大陆岩石层减薄机制的动力学模拟[J]. 地球物理学报,62(1):63–77. doi: 10.6038/cjg2019L0785 Shi Y N,Wei D P,Huangfu P P,Li Z H,Liu M X. 2019. Dynamics of thinning of overriding continental lithosphere induced by oceanic plate subduction:Numerical modeling[J]. Chinese Journal of Geophysics,62(1):63–77 (in Chinese).
沈晓明,张海祥,马林. 2010. 洋脊俯冲及其在新疆阿尔泰地区存在的可能证据[J]. 大地构造与成矿学,34(2):181–195. doi: 10.3969/j.issn.1001-1552.2010.02.004 Shen X M,Zhang H X,Ma L. 2010. Ridge subduction and the possible evidences in Chinese Altay,Xinjiang[J]. Geotectonica et Metallogenia,34(2):181–195 (in Chinese).
王振山,魏东平. 2018. 全球板块运动三联点形成与演化规律的研究进展[J]. 地球物理学进展,33(5):1834–1843. doi: 10.6038/pg2018BB0286 Wang Z S,Wei D P. 2018. Research progress on the formation and evolution of triple junctions of global plate motions[J]. Progress in Geophysics,33(5):1834–1843 (in Chinese).
徐佳静,王振山,王少坡,魏东平. 2019. 智利三联点南部扩张洋脊俯冲区域岩石层热结构的数值模拟[J]. 地球物理学报,62(12):4729–4737. doi: 10.6038/cjg2019M0592 Xu J J,Wang Z S,Wang S P,Wei D P. 2019. Numerical simulation of the lithospheric thermal structure in the subduction zone of the South Chile triple junction[J]. Chinese Journal of Geophysics,62(12):4729–4737 (in Chinese).
张克亮,魏东平. 2011. 双地震带的影响因素探讨[J]. 地球物理学报,54(11):2838–2850. doi: 10.3969/j.issn.0001-5733.2011.11.014 Zhang K L,Wei D P. 2011. On the influence factors of double seismic zones[J]. Chinese Journal of Geophysics,54(11):2838–2850 (in Chinese).
Agurto-Detzel H,Andreas R,Klaus B,Miller M,Iwamori H,Priestley K. 2014. Seismicity distribution in the vicinity of the Chile Triple Junction,Aysén Region,southern Chile[J]. J South Am Earth Sci,51:1–11. doi: 10.1016/j.jsames.2013.12.011
Assumpção M,Mei F,Andrés T,Julià J. 2013. Models of crustal thickness for South America from seismic refraction,receiver functions and surface wave tomography[J]. Tectonophysics,609:82–96. doi: 10.1016/j.tecto.2012.11.014
Bagherbandi M,Bai Y,Sjöberg L E,Tenzer R,Abrehdary M,Miranda S,Alcacer Sanchez J M. 2017. Effect of the lithospheric thermal state on the Moho interface:A case study in South America[J]. J South Am Earth Sci,76:198–207. doi: 10.1016/j.jsames.2017.02.010
Bangerth W, Dannberg J, Gassmoeller R, Heister T. 2020. ASPECT: Advanced solver for problems in Earth’s convection, user manual[BE/OL]. [2021-10-21]. https://www.math.clemson.edu/~heister/manual.pdf.
Bohm M,Lüth S,Echtler H,Bataille K,Group I W. 2002. The Southern Andes between 36° and 40° S latitude:Seismicity and average seismic velocities[J]. Tectonophysics,356(4):275–289. doi: 10.1016/S0040-1951(02)00399-2
Bondár I,Storchak D A. 2011. Improved location procedures at the International Seismological Centre[J]. Geophys J Int,186:1220–1244. doi: 10.1111/j.1365-246X.2011.05107.x
Bourgois J,Michaud F. 2002. Comparison between the Chile and Mexico triple junction areas substantiates slab window development beneath northwestern Mexico during the past 12−10 Myr[J]. Earth Planet Sci Lett,201:35–44. doi: 10.1016/S0012-821X(02)00653-2
Breitsprecher K,Thorkelson D J. 2009. Neogene kinematic history of Nazca-Antarctic Phoenix slab windows beneath Patagonia and the Antarctic Peninsula[J]. Tectonophysics,464(1/2/3/4):10–20. doi: 10.1016/j.tecto.2008.02.013
Contreras-Reyes E,Flueh E R,Grevemeyer I. 2010. Tectonic control on sediment accretion and subduction off south central Chile:Implications for coseismic rupture processes of the 1960 and 2010 megathrust earthquakes[J]. Tectonics,29:1–27.
Eakin C M,Obrebski M,Allen R M,Boyarko D C,Brudzinski M R,Porritt R. 2010. Seismic anisotropy beneath Cascadia and the Mendocino triple junction:Interaction of the subducting slab with mantle flow[J]. Earth Planet Sci Lett,297:627–632. doi: 10.1016/j.jpgl.2010.07.015
Fraters M,Thieulot C,van den Berg A,Spakman W. 2019. The geodynamic world builder:A solution for complex initial conditions in numerical modeling[J]. Solid Earth,(10):1785–1807.
Gleason G C,Tullis J. 1995. A flow law for dislocation creep of quartz aggregates determined with the molten salt cell[J]. Tectonophysics,247(1/2/3/4):1–23. doi: 10.1016/0040-1951(95)00011-B
Goddard A,Fosdick J C. 2019. Multichronometer thermochronologic modeling of migrating spreading ridge subduction in southern patagonia[J]. Geology,47(6):555–558. doi: 10.1130/G46091.1
Gresho P M,Lee R L,Sani R L,Maslanik M K,Eaton B E. 1987. The consistent Galerkin FEM for computing derived boundary quantities in thermal and/or fluids problems[J]. Int J Numer Meth Fl,7(4):371–394.
Guo C,Sun P,Wei D. 2021. Numerical simulation of the effects of wedge subduction on the lithospheric thermal structure and the seismogenic zone south of Chile Triple Junction[J]. Front Earth Sci,9:782458. doi: 10.3389/feart.2021.782458
Hamza V M,Dias F,Gomes A,Terceros Z. 2005. Numerical and functional representations of regional heat flow in South America[J]. Phys Earth Planet Inter,152:223–56. doi: 10.1016/j.pepi.2005.04.009
Heister T,Dannberg J,Gassmöller R,Bangerth W. 2017. High accuracy mantle convection simulation through modern numerical methods,II:Realistic models and problems[J]. Geophys J Int,210(2):833–851. doi: 10.1093/gji/ggx195
Hirth G, Kohlstedt D. 2004. Rheology of the upper mantle and the mantle wedge: A view from the experimentalists[C]// Geophysical Monograph Series. Washington D. C.: American Geophysical Union: 83–105.
Kirby S, Engdahl E R, Denlinger R. 1996. Intermediate-depth intraslab earthquakes and arc volcanism as physical expressions of crustal and uppermost mantle metamorphism in subducting slabs[C]// Geophysical Monograph Series. Washington D C: American Geophysical Union: 195–214.
Klotz J, Abolghasem A, Khazaradze G, Heinze B, Vietor T, Hackney R, Bataille K, Maturana R, Viramonte J, Perdomo R. 2006. Long-Term Signals in the Present-Day Deformation Field of the Central and Southern Andes and Constraints on the Viscosity of the Earth’s Upper Mantle[M]. Berlin: Springer: 65–89.
Kronbichler M,Heister T,Bangerth W. 2012. High Accuracy mantle convection simulation through modern numerical methods[J]. Geophys J Int,191(1):12–29. doi: 10.1111/j.1365-246X.2012.05609.x
Lagabrielle Y,Christèle G,René C M,Bourgois J,Martin H. 2000. Magmatic–tectonic effects of high thermal regime at the site of active ridge subduction:The Chile Triple Junction model[J]. Tectonophysics,326:255–268. doi: 10.1016/S0040-1951(00)00124-4
Lange D,Rietbrock A,Haberland C,Bataille K,Dahm T,Tilmann F,Flüh E R. 2007. Seismicity and geometry of the south Chilean subduction zone (41.5°S—43.5°S):Implications for controlling parameters[J]. Geophys Res Lett,34:L06311.
Lucazeau F. 2019. Analysis and mapping of an updated terrestrial heat flow data set[J]. Geochem Geophys Geosyst,20(8):4001–4024. doi: 10.1029/2019GC008389
Maksymowicz A,Eduardo C,Ingo G,Flueh E. 2012. Structure and geodynamics of the post-collision zone between the Nazca–Antarctic spreading center and South America[J]. Earth Planet Sci Lett,345-348:27–37.
Murdie R E,Prioe D J,Styles P,Flint S S,Agar S M. 1993. Seismic response to ridge-transform subduction:Chile triple junction[J]. Geology,21:1095–1098.
Naliboff J,Buiter S J H. 2015. Rift reactivation and migration during multiphase extension[J]. Earth Planet Sci Lett,421:58–67. doi: 10.1016/j.jpgl.2015.03.050
Oleskevich D A,Hyndman R D,Wang K. 1999. The updip and downdip limits to great subduction earthquakes:Thermal and structural models of Cascadia,south Alaska,SW Japan,and Chile[J]. J Geophys Res:Solid Earth,104:14965–14991. doi: 10.1029/1999JB900060
Richards F D,Hoggard M J,Cowton L R,White N J. 2018. Reassessing the thermal structure of oceanic lithosphere with revised global inventories of basement depths and heat flow measurements[J]. J Geophys Res:Solid Earth,123:9136–9161. doi: 10.1029/2018JB015998
Rose I,Buffett B,Heister T. 2017. Stability and accuracy of free surface time integration in viscous flows[J]. Phys Earth Planet Inter,262:90–100. doi: 10.1016/j.pepi.2016.11.007
Rybacki E. 2006. Influence of water fugacity and activation volume on the flow properties of fine-grained anorthite aggregates[J]. J Geophys Res:Solid Earth,200-201:1–9.
Scherwath M,Flueh E,Grevemeyer I,Tillman F,Contreras-Reyes E,Weinrebe W. 2006. Investigating subduction zone processes in Chile[J]. Eos Trans Am Geophys Un,87(27):265–272.
Scherwath M,Contreras-Reyes E,Flueh E R,Grevemeyer I,Krabbenhoeft A,Papenberg C,Petersen C J,Weinrebe R W. 2009. Deep lithospheric structures along the southern central Chile margin from wide-angle P-wave modelling[J]. Geophys J Int,179:579–600. doi: 10.1111/j.1365-246X.2009.04298.x
Shi Y,Wei D,Li Z,Liu M Q,Liu M. 2018. Subduction mode selection during slab and mantle transition zone interaction:Numerical modeling[J]. Pure Appl Geophys,175(2):529–548. doi: 10.1007/s00024-017-1762-0,175(2):529—548
Simmons N A,Forte A M,Boschi L,Grand S P. 2010. GyPSuM:A joint tomographic model of mantle density and seismic wave speeds[J]. J Geophys Res,115:B12310. doi: 10.1029/2010JB007631
Sisson V B,Pavlis T L,Roeske S M,Thorkelson D J. 2003. Introduction:An overview of ridge-trench interactions in modern and ancient settings[J]. Geol Soc Am Spec Papers,371:1–18.
Stein C A,Stein S. 1992. A model for the global variation in oceanic depth and heat flow with lithospheric age[J]. Nature,359(6391):123–129. doi: 10.1038/359123a0
Storchak D A,Harris J,Brown L,Lieser K,Shumba B,Verney R,Giacomo D D,Korger E I M. 2017. Rebuild of the bulletin of the International Seismological Centre (ISC),part 1:1964–1979[J]. Geosci Lett,4:32. doi: 10.1186/s40562-017-0098-z
Storchak D A,Harris J,Brown L,Lieser K,Shumba B,Di Giacomo D. 2020. Rebuild of the bulletin of the International Seismological Centre (ISC),part 2:1980–2010[J]. Geosci Lett,7:18. doi: 10.1186/s40562-020-00164-6
Tebbens S F,Cande S C,Kovacs L,Parra,J C,LaBrecque J L,Vergara H. 1997. The Chile ridge:A tectonic framework[J]. J Geophys Res:Solid Earth,102:12035–12059. doi: 10.1029/96JB02581
Tichelaar B W,Ruff L J. 1991. Seismic coupling along the Chilean subduction zone[J]. J Geophys Res,96(B7):11997. doi: 10.1029/91JB00200
Tetreault J L,Buiter S J H. 2012. Geodynamic models of terrane accretion:Testing the fate of island arcs,oceanic plateaus,and continental fragments in subduction zones[J]. J Geophys Res:Solid Earth,117(B8):B08403.
van der Hilst R,Hoop M. 2005. Banana-doughnut kernels and mantle tomography[J]. Geophys J Int,163:956–961. doi: 10.1111/j.1365-246X.2005.02817.x
Völker D,Ingo G,Michael S,Wang K,He J. 2011. Thermal control of the seismogenic zone of southern central Chile[J]. J Geophys Res,116:1–20.
Wei D, Seno T, 1998. Determination of the Amurian plate motion[C]// Mantle Dynamics and Plate Interactions in East Asia: Geodynamics Series. San Francisco: AGU: 337–346.
Wilks K R,Carter N L. 1990. Rheology of some continental lower crustal rocks[J]. Tectonophysics,182(1/2):57–77. doi: 10.1016/0040-1951(90)90342-6
Zhang K,Wei D. 2012. Correlation between plate age and layer separation of double seismic zones[J]. Earthquake Science,25:95–101. doi: 10.1007/s11589-012-0835-5