2019, Vol. 62
2. 中国科学院海洋研究所海洋地质与环境重点实验室, 青岛 266071;
3. 中国科学院大学, 北京 100049
2. Key Laboratory of Marine Geology and Environment, Institute of Oceanology, Chinese Academy of Sciences, Qingdao 266071, China;
3. University of Chinese Academy of Sciences, Beijing 100049, China
在地质时间尺度上(>106 yr), 岩石圈板块通常被当作“漂浮”于软流圈之上的弹性板块(Turcotte and Schubert, 2002).其厚度可以通过弹性板在受到载荷时发生的弹性挠曲变形来确定(Watts, 2001)，对应的厚度称为岩石圈有效弹性厚度(effective elastic thickness) (Te)，是表征岩石圈强度的参数，对认识岩石圈的力学性质具有重要意义(Forsyth, 1985; 付永涛等, 2002; Hunter and Watts, 2006; 杨安等, 2016).自20世纪70年代以来，国内外学者采用不同的方法计算了不同区域的岩石圈有效弹性厚度，取得了巨大的成果.目前，针对表面载荷(海山或山脉)引起板块挠曲，主要根据重力与地形的均衡响应函数进行研究.如Kalnins和Watts(2009)通过滑动窗口导纳技术(MWAT)估算了西太平洋海区的岩石圈有效弹性厚度，并探讨了海洋岩石圈有效弹性厚度与大洋年龄之间的关系.赵俐红等(2010)根据三维频率域中地形和重力的关系，计算了中西太平洋海山区的有效弹性厚度，并推测它们形成于白垩纪期间(约130~90 Ma).苏达权(2012)通过三维导纳法计算出我国南海中央海盆和南沙海域的岩石圈有效弹性厚度分别为6~7 km和10 km.利用相似的方法，胡敏章等(2015)对西北太平洋的岩石圈有效弹性厚度进行了研究，估算出其平均厚度为13.2 km；杨安等(2016)对西太平洋卡罗琳板块及其附近区域的岩石圈有效弹性厚度进行了研究，认为其变化范围为1~34 km.
但是，上述的重力导纳法只适用于由表面载荷(q)引起的板块挠曲的计算，对于海沟处俯冲板片的挠曲却不适用(Kalnins and Watts, 2009).因为重力导纳法总是假设岩石圈的挠曲变形是来源于上覆的载荷，而海沟处的板片挠曲主要来源于俯冲板片的拖曳，在模拟时通常将其作为一种边界载荷考虑(如海沟轴部的弯矩-M0和垂向剪力-V0，取向上为正).Caldwell等(1976)通过拟合垂直于海沟的深度剖面，分别计算了马里亚纳海沟、伊豆-小笠原海沟、阿留申海沟和千岛海沟的有效弹性厚度，并给出了适用于俯冲带的二维板片挠曲的解析解.利用相似的方法，McNutt(1984)计算了秘鲁和智利海沟的挠曲参数；Levitt和Sandwell(1995)模拟了全球大多数海沟的挠曲形态，并探讨了大洋岩石圈有效弹性厚度与洋壳年龄之间的关系，认为大洋岩石圈有效弹性厚度与洋壳年龄呈正相关.上述工作均将二维剖面的有效弹性厚度作为定值来进行反演，而最近的研究表明当俯冲板片在靠近海沟时，由于挠曲产生的应力大于岩石圈的屈服极限，会在海沟的隆起前缘形成一系列平行于海沟走向的挠曲相关正断层(bending-related normal faults)(Ranero et al., 2003).地震数据揭示挠曲相关正断层可以切至海底6 km的深度，到达上地幔的底部，并可以将水带入引起地幔的蛇纹石化(Ranero et al., 2005).这个过程可能导致大洋岩石圈的刚度降低，因此在模拟海沟俯冲板片的挠曲时用固定的有效弹性厚度并不准确(Contreras-Reyes and Osses, 2010).Contreras-Reyes和Osses(2010)通过对智利海沟大洋岩石圈的挠曲模拟，提出了模拟二维变有效弹性厚度的板片挠曲的有限差分法，并指出俯冲板块有效弹性厚度的减小与大洋年龄之间的关系不明显. Zhang等(2014)利用该方法计算了马里亚纳海沟大洋岩石圈的挠曲参数，并认为挑战者深渊之所如此深，是因为其具有最大的轴部垂向剪力和有效弹性厚度的减小量.但是这些研究都是以二维薄板模型为基础，忽略了沿海沟走向的影响.关于海沟俯冲板片三维挠曲模拟的研究相对较少.Wessel(1996)给出了适用于俯冲带的岩石圈三维挠曲变形的解析解，但是他并没有考虑有效弹性厚度在平面空间的可变性.Manríquez等(2014)利用有限元方法求解了三维厚板挠曲方程，并模拟了智利海沟大洋岩石圈的三维挠曲变形，为在三维尺度研究俯冲带岩石圈的挠曲变形提供了新思路.然而该研究只考虑了沿着垂直于海沟方向有效弹性厚度的变化，并没有考虑沿海沟走向的变化.
为了能够考虑任意方向挠曲参数的变化(Te, q, -M0和-V0)，本文利用有限差分方法求解三维薄板挠曲方程，建立了适用于俯冲-碰撞带的岩石圈挠曲模型.在反演挠曲参数时，因参数量大而无法运用常规的反演方法，我们提出了基于粒子群优化(Particle Swarm Optimization，PSO)算法的俯冲带三维有效弹性厚度反演.并模拟了南海马尼拉海沟的岩石圈挠曲形态，反演了其挠曲参数.
1 方法原理 1.1 二维挠曲模型海沟处俯冲板片挠曲的基本模型通常如图 1所示.根据Kirchhoff-Love薄板理论，在弹性域内，当板块在海沟处受到俯冲板块的拖曳作用，会产生挠曲变形，其变形形态ω可在二维情况下通过下列方程表示(Turcotte and Schubert, 2002)：
|
(1) |
|
图 1 俯冲板块挠曲模型.M0和V0为海沟轴部的弯矩和剪力.ρm，ρc和ρw分别表示地幔、地壳和海水的密度 Fig. 1 Schematic model of flexural bending of the subducting plates. M0 and V0 represent the bending moment and vertical force applied at trench axis respectively. ρm, ρc and ρw are densities of mantle, crust and sea water, respectively |
其中，Δρ=ρm-ρw，表示地幔(ρm)和海水(ρw)的密度差.q(x)为由海山或者沉积物造成的表面载荷，由于沉积载荷相对于边界剪力可以忽略，通常可以令q(x)为0.M为弯矩，可由挠曲量w表示：
|
(2) |
其中α为挠曲波长，
|
(3) |
二维挠曲模型忽略了沿海沟走向挠曲参数的变化的影响，其假设在沿海沟走向上挠曲参数是不变的.但是由于岩石圈的不均匀性及构造运动的差异性，挠曲参数通常是变化的，二维挠曲模型无法真实反映俯冲带岩石圈的挠曲变形.因此，我们将二维挠曲模型扩展为三维模型，以更好的模拟岩石圈在俯冲带的挠曲变形.
1.2 三维挠曲模型考虑薄板的x和y两个维度时，板块的挠曲刚度D(x, y)即成为空间的函数(Garcia et al., 2015).微分方程(1)则可展开为如下形式(van Wees and Cloetingh, 1994; Braun et al., 2013)：
|
(4) |
同样，俯冲板片的挠曲不考虑水平力和上覆载荷，即令F和q(x, y)等于零.此时，弯矩M和剪力V与挠曲量w之间的关系分别为：
|
(5) |
|
(6) |
(6) 式可以看出，在三维情况下，Vx和Vy均受到了来自横向的扭矩
|
图 2 本文采用的三维挠曲模型.Γ1—Γ4表示不同的边界条件(见公式(7)) Fig. 2 Schematic 3-D model used in this study. Γ1—Γ4 represent different boundary conditions (See formula (7)) |
|
(7) |
此时轴部弯矩(-M0)和剪力(-V0)都不再是单一值，而变为向量.根据有限差分方法原理，将矩形板划分为间隔分别为dx和dy的网格(本文为200×200个网格)，并将偏微分方程改写为差分方程：AW=Q，其中A为系数矩阵，W为挠曲向量，Q为载荷向量.求解线性方程组后即可得到板块的三维挠曲形态.为了验证该方法的正确性，我们将三维模型沿海沟走向的挠曲参数设置为不变，并与前人的二维模拟结果进行对比：
(1) 固定M0，V0和Te(图 3a)：分别设置轴部剪力和弯矩为V0=0.19×1012 N·m-1，M0=1.37×1016 N和Te=25.88 km(与Contreras-Reyes和Osses(2010)文中剖面06相同)，在该条件下，二维挠曲形态可直接由方程(2)得出.从三维结果中抽取一条剖面与二维结果对比，两者结果相差在0.2%以内(图 3a).
|
图 3 本文三维模型与前人工作的对比 (a)有效弹性厚度为定值时本文与前人工作的对比；(b)可变有效弹性厚度时，本文与Contreras-Reyes和Osses(2010)工作的对比；(c)同时加载边界载荷和上覆载荷时，本文与Contreras-Reyes和Osses(2010)工作的对比. Fig. 3 Comparisons between 3-D model in this paper with the previous work (a) Comparisons between our model with the 2-D analytic solution with uniform Te; (b) Comparisons between our model with the 2-D model of Contreras-Reyes and Osses (2010) with variable Te; (c) Comparisons between our model with the 2-D model of Contreras-Reyes and Osses (2010) with boundary loading and surface loading. |
(2) 固定M0和V0，Te可变(图 3b)：与实验1一样，分别设置轴部剪力和弯矩为V0=0.19×1012 N·m-1，M0=1.37×1016 N，Te设置为在边界附近由25.88 km线性减小至14.2 km.从三维结果中抽取一条剖面与Contreras-Reyes和Osses(2010)的二维模拟结果相比较，相差在1%以内(图 3b).
(3) 固定M0，V0和Te，同时考虑上覆载荷q(图 3c)：M0, V0和Te的值与实验1设置为一样，同时令上覆载荷q= -5.23×106 N·m-1.将三维模拟结果与Contreras-Reyes和Osses(2010)的二维模拟结果相比较，其相差在0.3%以内(图 3c).
上述正演数值实验可以看出，我们的三维模型具有较好的稳定性和精度.
1.3 粒子群优化算法的基本原理粒子群优化算法(PSO)是由Eberhart和Kennedy在1995年提出的一种基于群体的随机优化技术.它将在空间搜索最优解的过程类比为鸟类在飞行空间寻找食物的过程：将每个候选解抽象为无质量无体积，具有一定运动速度的微粒，所有粒子都通过一个适应度函数来决定粒子位置的优劣.粒子之间并不相互独立，可以通过学习群体和自身的经验来不断地调整自己的速度和位置.
设在一个M维目标搜索空间中，有n个粒子组成的一个群体.在t时刻，第i个粒子的空间位置为Xi=(xi1, xi2, xi3, , …, xiM)，i=1, 2, 3, …, n，速度为Vi=(vi1, vi2, vi3, …, viM).每个粒子的位置都是一个可行解，将其代入一个目标函数计算出其适应值，并根据适应值的大小来判断解的优劣.记第i个粒子到目前为止经历的最优位置为Pi=(Pi1, Pi2, Pi3, …, PiM)，整个粒子群获得的最优位置为Pgi=(Pg1, Pg2, Pg3, …, PgM).粒子通过如下公式来更新自身的速度和位置(Shi and Eberhart, 1998)：
|
(8) |
|
(9) |
其中，c1和c2称为学习因子，为非负的常数，r1和r2为相互独立的，在0到1之间均匀分布的随机数.vim为粒子的飞行速度，通常限制在一定的范围里面，即vim∈[-vmax, vmax]，vmax为最大速度.如果搜索空间为[-xmax, xmax]，那么可以设vmax=kxmax，其中0.1≤ k≤ 1.0.ω为非负数，称为惯性权重，代表了前一时刻速度对当前速度的影响.当ω较大时，前一时刻的速度影响较大，粒子的全局搜索能力较强；当ω较小时，前一时刻的速度影响小，局部搜索能力较强.众多学者讨论了ω对寻找最优解的影响(Shi and Eberhart, 1998; 岳碧波等, 2009)，其中较常用的是线性递减的权重PSO算法.即在搜索的过程中，让惯性权重从最大的ωmax线性减小至ωmin，其公式为：
|
(10) |
式中，ωmax，ωmin分别代表惯性权重的最大值和最小值，t代表当前迭代步数，tmax表示最大迭代步数.这种设置满足了在搜索初期，具有较大的惯性权重，便于全局搜索；而在搜索后期惯性权重变小，利于局部搜索.
2 利用PSO反演挠曲参数如前所述，在三维挠曲模拟过程中，轴部弯矩(-M0)和剪力(-V0)均为向量，有效弹性厚度(Te)为与节点大小相同的矩阵，待反演的参数量非常巨大.为了简化反演过程，我们主要考虑挠曲参数在沿海沟走向上变化，假设垂直于海沟走向上Te不变.此时待反演的Te也简化为向量.将-M0，-V0和Te合并为一个向量：X= (-M0，-V0，Te)，X即为PSO反演中的粒子.
2.1 适应度计算在粒子群反演中，重要的一步是根据目标函数计算每个粒子的适应度值，然后根据适应度值的大小来衡量解的优劣.本文将实测挠曲形态数据与模拟挠曲形态数据之间的最小二乘误差作为适应度函数，具体如下：
|
(11) |
其中N为节点数，ωobs为观测的挠曲形态数据，ωcalc为根据公式(5)，(6)和(7)计算的挠曲形态数据.从公式(11)可以看出，粒子搜索到的模拟挠曲形态越接近实测挠曲形态，该粒子的适应度值就越小.
2.2 PSO有效弹性厚度反演测试为了验证反演方法的可行性，我们先将-M0，-V0设置为定值，以Te可变为例，描述粒子群算法的有效弹性厚度反演过程.具体的测试步骤为：
(1) 首先构造出一个三维挠曲形态，作为“观测值”：给定-M0，-V0的值，并设置可变的Te，利用公式(5)，(6)和(7)构造出一个挠曲形态，作为公式(11)中的ωobs.为简化测试模型，本文设置600 km×600 km的方形板，网格大小为200×200.测试的目的就是利用PSO方法反演出变化的Te.
(2) 设置PSO算法的参数，本文设置：学习因子c1=c2=2，速度系数k=0.15，惯性权重ωmin=0.4，ωmax=0.9.
(3) 在约束条件下随机设置粒子群中每个粒子的初始速度和位置.本文为了加快搜索速度，先利用二维挠曲模型反演挠曲参数，然后以该结果上下波动30%作为三维反演的搜索空间，以约束求解空间的大小.
(4) 利用公式(8)，(9)对粒子群进行处理，不断更新其速度和位置.以公式(11)作为判断运算结束的条件，最终输出最优解.
为了对比二维和三维反演的结果，我们设置不同的Te变化模式和M0，V0变化模式，并分别利用传统的二维反演方法和本文的三维反演方法进行挠曲参数的反演(图 4，5和6).结果显示，Te沿边界走向呈突变时，基于PSO算法的三维有效弹性厚度反演可以识别有效弹性厚度发生突变的位置，反演结果与初始设置条件的误差在5%以内，局域较好的反演效果.而二维反演则无法识别有效弹性厚度突变的位置；当Te沿边界走向呈线性变换时，二维和三维反演效果均较好，与初始设置条件的误差在5%以内(图 4).图 5分别将Te设置为抛物线型和周期型变化，其中抛物线型变化时，二维和三维的反演结果误差在5%以内，而周期型变化时，二维反演的结果误差较大，可达近30%，我们认为这和周期变化的波长(L)与板块的特征波长(α)的比值有关：即L/α的值越大，反演的误差越小.图 6将Te设置为定值(25 km)，边界弯矩M0和剪力V0设置为周期型变化，二维和三维反演的Te结果误差在2%以内.
|
图 4 有效弹性厚度沿走向呈突变和线性变化时，二维与三维方法反演有效弹性厚度的结果差别 (a)三维模型参数；(b)三维正演结果；(c)二维和三维反演方法结果差别；(d)三维反演方法的RMS值. Fig. 4 The difference between Te estimated by 2-D and 3-D model with lateral sharp and linear variable Te (a) 3-D models parameters; (b) 3-D forward results; (c) The difference between Te estimated by 2-D, 3-D inversion results; (d) The RMS of 3-D inversion results. |
|
图 5 有效弹性厚度沿走向呈抛物线型和周期型变化时，二维与三维方法反演有效弹性厚度的结果差别 (a)三维模型参数；(b)三维正演结果；(c)二维和三维反演方法结果差别；(d)三维反演方法的RMS值. Fig. 5 The difference between Te estimated by 2-D and 3-D model with lateral parabola and periodic variable Te (a) 3-D models parameters; (b) 3-D forward results; (c) The difference between Te estimated by 2-D, 3-D inversion results; (d) The RMS of 3-D inversion results. |
|
图 6 边界沿走向呈周期性变化时，二维与三维方法反演有效弹性厚度的结果差别 (a)三维模型参数；(b)三维正演结果；(c)二维和三维反演方法结果差别；(d)三维反演方法的RMS值. Fig. 6 The difference between Te estimated by 2-D and 3-D model with periodic variable M0 and V0 (a) 3-D models parameters; (b) 3-D forward results; (c) The difference between Te estimated by 2-D, 3-D inversion results; (d) The RMS of 3-D inversion results. |
此外，由图 4, 5和6的迭代曲线看出，RMS经常是一个突然下降，然后进入一段迭代收敛的平缓期，然后又突然下降，这是由粒子群算法自身的搜寻方式导致的.因为在真正寻得最优解之前，所有粒子的位置会随着迭代次数的增加而不断地调整，RMS值也会因粒子越来越接近最优解而一直下降，直到搜寻到最优解.所以即便是在迭代的后期，当某次迭代出现更接近最优解的粒子时，也有可能出现RMS值的突然下降.
3 马尼拉海沟有效弹性厚度反演马尼拉海沟位于南海海盆东缘，由南海海盆俯冲至菲律宾海板块之下形成，其向北与台湾碰撞造山带相连，向南与民都洛地震构造带相连，是一条正在活动的板块汇聚边界(图 7)(李家彪等, 2004).研究认为马尼拉海沟的形成时间为中中新世(Wolfe, 1988)或晚渐新世(Hayes and Lewis, 1984; Wang and Li, 2009).南海岩石圈向东俯冲形成贝尼奥夫带(Benioff zone)，其俯冲倾角从北部吕宋岛的约40°向南逐渐变为近垂直(Hayes and Lewis, 1984).近年研究表明，在15°N—16°N附近存在古扩张脊的俯冲，俯冲板片在100 km左右的深度可能存在撕裂(詹文欢等, 2017).这种深部的板片撕裂是否会对浅部的岩石圈挠曲变形产生影响？还需进一步探讨.
|
图 7 南海马尼拉海沟海底地形图 海底地形参考30 arc seconds水深数据(Becker et al., 2009)，白色虚线为洋陆边界(COB)(参考Briais et al., 1993)，黑色方框为本文研究区域. Fig. 7 Regional bathymetry map of the Manila trench and the South China Sea Bathymetry data are from Becker et al. (2009). The white dashed line represents the continent oceanic boundary (COB) (from Briais et al., 1993) and the black box is the area of this study. |
通过对马尼拉海沟外缘隆起形态的研究，Chang等(2012)利用二维挠曲方法反演了马尼拉海沟的有效弹性厚度为13~16 km，并认为南海海盆岩石圈的Te与海盆年龄的算术平方根呈线性关系.Zhang等(2018)通过求解二维薄板挠曲方程，在可变有效弹性厚度情况下，计算马尼拉海沟的有效弹性厚度为4~28 km.但是上述研究均局限于二维的挠曲参数反演，并没有考虑来自横向的参数变化和影响.
本文选取经度为116°E—119°E，纬度为14°N—17°N的矩形作为我们的研究区域(见图 7黑色方框)，利用PSO方法对该区域的挠曲参数进行了反演：为了有效地提取南海俯冲板片的挠曲形态，用海底地形数据减去了海洋沉积物厚度.同时，因为研究区的海山(洋中脊)较多，而常规的数字滤波很难完全将海底的局部地形干扰消除(Zhang et al., 2018)，因此我们根据地形数据解释了长波长的构造挠曲形态(图 8b)，其结果作为观测数据ωobs，然后利用公式(8)—(11)进行反演.
|
图 8 马尼拉海沟三维挠曲模拟及反演结果 (a)海底地形(Becker et al., 2009)；(b)解释的挠曲形态；(c)三维模拟挠曲形态；(d)三维反演的挠曲参数，其中黄色区域为洋中脊的位置. Fig. 8 3-D flexural modelling and inverted results of the Manila trench (a) Bathymetry map of study area (Becker et al., 2009); (b) Interpreted flexural shape; (c) 3-D flexural modelling result; (d) 3-D inverted flexural parameters, and the yellow area represents the position of the middle ocean ridge. |
结果显示：南海中央海盆岩石圈有效弹性厚度(Te)在洋中脊以南的区域为10~15 km，在洋中脊以北的区域为5~10 km，Te随距离洋中脊的距离增大而增加，且对于南海而言，距洋中脊向南Te增加更快.这与前人研究认为大洋岩石圈有效弹性厚度与洋壳年龄成正相关的结论一致(Chang et al., 2012; Hunter and Watts, 2016).而海沟轴部弯矩(-M0)则在洋中脊两侧存在分段现象，即在洋中脊以南的区域，随着距洋中脊的距离的增加，M0的值由-0.2×1016 N增加至2.3×1016 N；而在洋中脊北部，M0的值在-0.2×1016 N到0.1×1016 N之间.我们认为这可能与南海板块在洋中脊附近发生撕裂有关，板片在深部撕裂后，由深部传递至浅部的拖曳力发生改变，从而对浅部岩石圈的挠曲形态产生影响.海沟轴部剪力(-V0)的变化形态与海沟深度变化形态具有较好的耦合，我们认为V0主要反映了海沟深度的变化.
4 结论本文通过求解三维薄板挠曲方程，提出了适用于俯冲-碰撞的岩石圈三维挠曲正演模型.同时发展了基于粒子群算法的俯冲带三维有效弹性厚度反演方法，并用于对南海马尼拉海沟岩石圈挠曲的模拟，主要得到以下结论：
(1) 南海中央海盆岩石圈有效弹性厚度随着距洋中脊距离的增加而增大，且向南增加更快；
(2) 马尼拉海沟轴部弯矩变化在洋中脊两侧存在分段性，即洋中脊以南的值大于以北的值.这种分段性可能与南海板片在深部沿洋中脊发生撕裂有关.
Becker J J, Sandwell D T, Smith W H F, et al. 2009. Global bathymetry and elevation data at 30 arc seconds resolution:SRTM30_PLUS. Marine Geodesy, 32(4): 355-371. DOI:10.1080/01490410903297766 |
Braun J, Deschamps F, Rouby D, et al. 2013. Flexure of the lithosphere and the geodynamical evolution of non-cylindrical rifted passive margins:results from a numerical model incorporating variable elastic thickness, surface processes and 3D thermal subsidence. Tectonophysics, 604: 72-82. DOI:10.1016/j.tecto.2012.09.033 |
Briais A, Patriat P, Tapponnier P. 1993. Updated interpretation of magnetic anomalies and seafloor spreading stages in the South China Sea:implications for the Tertiary tectonics of Southeast Asia. J. Geophys. Res., 98(B4): 6299-6328. DOI:10.1029/92JB02280 |
Caldwell J G, Haxby W F, Karig D E, et al. 1976. On the applicability of a universal elastic trench profile. Earth Planet. Sci. Lett., 31(2): 239-246. DOI:10.1016/0012-821X(76)90215-6 |
Chang J H, Yu H S, Lee T Y, et al. 2012. Characteristics of the outer rise seaward of the Manila Trench and implications in Taiwan-Luzon convergent belt, South China Sea. Mar. Geophys. Res., 33(4): 351-367. DOI:10.1007/s11001-013-9168-6 |
Contreras-Reyes E, Osses A. 2010. Lithospheric flexure modelling seaward of the Chile trench:implications for oceanic plate weakening in the Trench Outer Rise region. Geophys. J. Int., 182(1): 97-112. |
Forsyth D W. 1985. Subsurface loading and estimates of the flexural rigidity of continental lithosphere. J. Geophys. Res., 90(B14): 12623-12632. DOI:10.1029/JB090iB14p12623 |
Fu Y T, Li A C, Qin Y S. 2002. Effective elastic thickness of the oceanic and continental marginal lithospheres. Marine Geology & Quaternary Geology (in Chinese), 22(3): 69-75. |
Garcia E S, Sandwell D T, Luttrell K M. 2015. An iterative spectral solution method for thin elastic plate flexure with variable rigidity. Geophys. J. Int., 200(2): 1012-1028. DOI:10.1093/gji/ggu449 |
Hayes D E, Lewis S D. 1984. A geophysical study of the Manila Trench, Luzon, Philippines:1. Crustal structure, gravity, and regional tectonic evolution. J. Geophys. Res., 89(B11): 9171-9195. |
Hu M Z, Li J C, Li H, et al. 2015. The lithosphere effective elastic thickness and its tectonic implications in the Northwestern Pacific. Chinese J. Geophys. (in Chinese), 58(2): 542-555. DOI:10.6038/cjg20150217 |
Hunter J, Watts A B. 2016. Gravity anomalies, flexure and mantle rheology seaward of circum-Pacific trenches. Geophys. J. Int., 207(1): 288-316. DOI:10.1093/gji/ggw275 |
Kalnins L M, Watts A B. 2009. Spatial variations in effective elastic thickness in the Western Pacific Ocean and their implications for Mesozoic volcanism. Earth Planet. Sci. Lett., 286(1-2): 89-100. DOI:10.1016/j.epsl.2009.06.018 |
Levitt D A, Sandwell D T. 1995. Lithospheric bending at subduction zones based on depth soundings and satellite gravity. J. Geophys. Res., 100(B1): 379-400. DOI:10.1029/94JB02468 |
Li J B, Jin X L, Ruan A G, et al. 2004. Indentation tectonics in the accretionary wedge of middle Manila trench. Chinese Science Bulletin, 49(12): 1279-1288. DOI:10.1360/03wd0412 |
Manríquez P, Contreras-Reyes E, Osses A. 2014. Lithospheric 3-D flexure modelling of the oceanic plate seaward of the trench using variable elastic thickness. Geophys. J. Int., 196(2): 681-693. DOI:10.1093/gji/ggt464 |
McNutt M K. 1984. Lithospheric flexure and thermal anomalies. J. Geophys. Res., 89(B13): 11180-11194. DOI:10.1029/JB089iB13p11180 |
Ranero C R, Morgan J P, McIntosh K, et al. 2003. Bending-related faulting and mantle serpentinization at the Middle America trench. Nature, 425(6956): 367-373. DOI:10.1038/nature01961 |
Ranero C R, Villaseñor A, Morgan J P, et al. 2005. Relationship between bend-faulting at trenches and intermediate-depth seismicity. Geochem. Geophys. Geosyst., 6(12): Q12002. DOI:10.1029/2005GC000997 |
Shi Y, Eberhart R. 1998. A modified particle swarm optimizer.//Proceeding of 1998 IEEE International Conference on Evolutionary Computation. Anchorage, AK, USA: IEEE, 69-73.
|
Su D Q. 2012. A study of the effective elastic thickness of the oceanic lithosphere. Chinese J. Geophys. (in Chinese), 55(10): 3259-3265. |
Turcotte D L, Schubert G. 2002. Geodynamics. 2nd ed. Cambridge: Cambridge University Press.
|
van Wees J D, Cloetingh S. 1994. A finite-difference technique to incorporate spatial variations in rigidity and planar faults into 3-D models for lithospheric flexure. Geophys. J. Int., 117(1): 179-195. DOI:10.1111/j.1365-246X.1994.tb03311.x |
Wang P X, Li Q Y. 2009. Oceanographical and geological background.//Wang P X, Li Q Y eds. The South China Sea. Dordrecht: Springer, 25-73.
|
Watts A B. 2001. Isostasy and Flexure of the Lithosphere. Cambridge: Cambridge University Press.
|
Wessel P. 1996. Analytical solutions for 3-D flexural deformation of semi-infinite elastic plates. Geophys. J. Int., 124(3): 907-918. DOI:10.1111/j.1365-246X.1996.tb05644.x |
Wolfe J A. 1988. Arc magmatism and mineralization in North Luzon and its relationship to subduction at the East Luzon and North Manila Trenches. J. Southeast Asian Earth Sci., 2(2): 79-93. DOI:10.1016/0743-9547(88)90011-6 |
Yang A, Fu Y T, Li A C. 2016. The effective elastic thickness of the Caroline plate and its adjacent areas. Chinese J. Geophys. (in Chinese), 59(9): 3280-3290. DOI:10.6038/cjg20160913 |
Yue B B, Peng Z M, Hong Y G, et al. 2009. Wavelet inversion of pre-stack seismic angle-gather based on particle swarm optimization algorithm. Chinese J. Geophys. (in Chinese), 52(12): 3116-3123. |
Zhan W H, Li J, Tang Q Q. 2017. Subduction of the paleo-spreading-ridge in eastern South China Sea. Marine Geology & Quaternary Geology (in Chinese), 37(6): 1-11. |
Zhang F, Lin J, Zhan W H. 2014. Variations in oceanic plate bending along the Mariana Trench. Earth Planet. Sci. Lett., 401: 206-214. DOI:10.1016/j.epsl.2014.05.032 |
Zhang F, Lin J, Zhou Z Y, et al. 2018. Intra-and intertrench variations in flexural bending of the Manila, Mariana and global trenches:implications on plate weakening in controlling trench dynamics. Geophys. J. Int., 212(2): 1429-1449. DOI:10.1093/gji/ggx488 |
Zhao L H, Jin X L, Gao J Y, et al. 2010. The effective elastic thickness of lithosphere in the mid-west Pacific and its geological significance. Earth Science-Journal of China University of Geosciences (in Chinese), 35(4): 637-644. DOI:10.3799/dqkx.2010.078 |
付永涛, 李安春, 秦蕴珊. 2002. 大洋和大陆边缘岩石圈有效弹性厚度的研究意义. 海洋地质与第四纪地质, 22(3): 69-75. |
胡敏章, 李建成, 李辉, 等. 2015. 西北太平洋岩石圈有效弹性厚度及其构造意义. 地球物理学报, 58(2): 542-555. DOI:10.6038/cjg20150217 |
李家彪, 金翔龙, 阮爱国, 等. 2004. 马尼拉海沟增生楔中段的挤入构造. 科学通报, 49(10): 1000-1008. DOI:10.3321/j.issn:0023-074X.2004.10.015 |
苏达权. 2012. 海洋岩石圈板块有效弹性厚度研究. 地球物理学报, 55(10): 3259-3265. DOI:10.6038/j.issn.0001-5733.2012.10.008 |
杨安, 付永涛, 李安春. 2016. 卡罗琳板块及其附近地区的岩石圈有效弹性厚度. 地球物理学报, 59(9): 3280-3290. DOI:10.6038/cjg20160913 |
岳碧波, 彭真明, 洪余刚, 等. 2009. 基于粒子群优化算法的叠前角道集子波反演. 地球物理学报, 52(12): 3116-3123. DOI:10.3969/j.issn.0001-5733.2009.12.021 |
詹文欢, 李健, 唐琴琴. 2017. 南海东部古扩张脊的俯冲机制. 海洋地质与第四纪地质, 37(6): 1-11. |
赵俐红, 金翔龙, 高金耀, 等. 2010. 中西太平洋海山区的岩石圈有效弹性厚度及其地质意义. 地球科学——中国地质大学学报, 35(4): 637-644. |