山东地区地壳三维纵横波速度和泊松比结构及其对地震活动性的影响＊

范建柯丁志峰郭慧丽苏道磊张斌，

1） 中国北京 100081 中国地震局地球物理研究所

2） 中国山东青岛 266071 中国科学院海洋研究所

3） 中国山东青岛 266071 中国科学院海洋地质与环境重点实验室

4） 中国济南 250001 济南市地震监测中心

引言

山东地区位于中国大陆东部，自中生代以来经历了复杂的地质构造运动，如中生代早期华北地块与扬子地块的碰撞、中生代中晚期的华北克拉通破坏以及新生代以来的拉张变形等.陆地区域内发育沂沭断裂带（郯庐断裂带山东段）和聊考断裂带两大断裂系统，北部海域发育张家口——蓬莱断裂带，导致地震活动频繁，并伴有大地震的发生，如1668年郯城M8.5地震，造成了重大的人员伤亡和经济损失.

2003年6月5日和2020年2月18日，山东省最大的两个城市青岛和济南分别发生了一次M4.1地震（图1），虽然未导致较大的人员伤亡和经济损失，但由于其发生于人口均近千万的两大城市，且有感范围较大，余震次数较多，而上述两个地区通常被认为属于“弱震”区，因此，这两次地震仍然引发了非常大的关注.前人综合余震重定位结果、震源机制解等数据，推测这两次地震可能属于相对完整岩体条件下的一次新破裂活动（潘元生等，2004）或区域构造应力作用下附近断裂或次级派生断裂活动的结果（张斌等，2020）.另外，自2000年以来，山东半岛地区相继发生了崂山震群、乳山震群和长岛震群，但目前上述发震区尚未有高精度的层析成像数据，速度结构与发震机理仍不清楚.因此，有必要查明上述地区的浅部及深部速度结构形态，探讨速度结构与地震发震机理之间的关系，为后期的防震减灾提供理论支持.

图1 山东及周边地区构造分区图（修改自苏道磊等，2016）Fig.1 Tectonic settings in and around Shandong area（modified from Su et al，2016）

1 数据与方法

本文所用数据来源于山东及邻区96个台站（图2a）记录的区域内地震事件的P波和S波到时数据，主要分为两部分：第一部分为1975年1月至2014年1月发生的天然地震事件（苏道磊等，2016）；第二部分为2016年1月至2019年12月发生的天然地震事件.原始数据共包含7 271 个地震事件.

图2 台站（a）及地震（b）分布图Fig.2 Distribution of seismic stations （a） and local earthquakes （b） used in this study

为保证反演结果的准确性，需要对地震事件进行严格挑选，设定标准如下：① 每个地震事件至少被4个台站接收；② 震相走时残差绝对值小于3.0 s；③ 重定位前后，发震时刻偏差小于2.0 s，水平向偏差小于6 km，震源深度偏差小于8 km.最终筛选出4 652个地震事件，包括3万6 482个P波震相到时和3万2 600个S波震相到时用于成像反演（图2b）.

本文采用TOMOG3D三维层析成像反演方法（Zhaoet al，1992）.该方法采用三维网格节点表示空间的速度分布，空间内每一点的速度值由周围八个节点的速度值进行线性插值获得，允许三维空间内存在间断面，通过伪弯曲算法快速准确地计算射线路径和走时，并能够同时处理近远震及后续震相到时数据.

合适的一维初始速度模型对最终反演结果的准确性至关重要.苏道磊等（2016）分别测试了三个初始速度模型：① 鲁西地区地壳速度模型，同时考虑莫霍面起伏（嘉世旭，张先康，2005）；② 山东地震台网定位用地壳速度模型；③ 华北地区地壳速度模型（陈立华，宋仲和，1990）.测试结果表明，鲁西地区地壳速度模型的走时残差均方根最小，也更符合真实地层情况，而华北地区地壳速度模型走时残差远大于其它两个速度模型.因此，本文仅对比前两个初始速度模型的走时残差，并同时考虑莫霍面起伏的影响.其中，莫霍面埋深数据来自于CRUST1.0模型（Laskeet al，2013）和郑宏等（2021）利用接收函数反演得到的山东地区莫霍面埋深等（图3）.初始S波速度模型由P波速度除以1.732得到.

图3 基于不同模型反演所得的莫霍面埋深（a，b）及本研究所用莫霍面埋深（c）Fig.3 Moho depths inverted from different models （a，b） and used in this study （c）

经过计算，上述两个初始速度模型的残差均方根分别为0.841 s和0.860 s，走时残差分布如图4所示.对比发现，鲁西地区速度模型要优于山东地震台网定位用地壳速度模型，与苏道磊等（2016）的结论相似.因此，后续的层析成像反演采用鲁西地区速度模型，并考虑莫霍面起伏的影响.

图4 使用鲁西地区地壳速度模型（a）和山东地震台网定位用地壳速度模型（b）在地震重定位前后所得的走时残差分布图Fig.4 Travel time residual distribution based on crustal velocity model in Luxi area （a） and Shandong seismic network （b） before and after relocation

本研究中，初始模型三维网格节点横向间隔为0.2°×0.2°，纵向上分别在1，10，20和30 km深度处设置节点层.反演前，基于初始速度模型和原始到时数据对所有地震事件进行了重新定位.重定位前后的走时残差分布如图4所示，重定位后的总体走时残差均方根由0.841 s降低到 0.649 s，其中 P 波走时残差均方根由 0.792 s降低到 0.620 s，S 波走时残差均方根由0.893 s降低到0.680 s，说明地震定位精度有了较大幅度的提高.走时残差绝大多数位于±2 s以内，因此选择走时残差在±2 s以内的震相参与最终的成像反演.反演采用带阻尼和平滑因子的最小二乘（least squares QR-factorization，缩写为 LSQR）方法（Paige，Saunders，1982）得到最终的三维速度结构.通过大量的测试，P波和S波速度反演采用的最优阻尼和平滑因子均为 5.0 和 50.0 （图5）.

图5 P 波（a，b）和 S 波（c，d）成像中不同阻尼和平滑因子对应的三维速度模型标准差与走时残差均方根关系曲线Fig.5 Trade-off curves between the standard deviation of the 3-D velocity model and the root-mean-square of travel time residual according to the damping （a and c） and smoothing （b and d） parameters for P- （a，b） and S-wave （c，d） tomographies

2 分辨率测试

在分析成像结果之前，需要对成像结果的可靠性进行评估，一般采用检测板（checkerboard）方法（Zhaoet al，1992）.在该方法中，首先将三维空间内相邻网格节点分别设置±3%的速度扰动，然后利用相同的地震和台站分布计算理论走时，并在计算理论走时过程中加入标准差为0.1 s的随机误差来检测计算稳定性.随后对得到的理论走时在初始一维速度模型的基础上进行反演，通过对比三维网格节点处扰动值的恢复情况对成像结果的可靠性进行评估.如果反演后扰动值与反演前扰动值的分布相似，则说明分辨率较好.

图6显示了横向间隔分别为 0.33°×0.33°，0.4°×0.4°和 0.5°×0.5°的分辨率测试结果.结果显示，对于P波和S波成像结果，研究区内大部分地区的分辨率能达到0.4°×0.4°，部分地区如鲁西南地区和山东半岛地区可达到0.33°×0.33°的分辨率.其中：网格间距为0.33°×0.33°时，P波和S波检测板扰动值恢复准确性分别达到74%和81%，振幅恢复达到70%的比例约为47%和59%，振幅恢复达到100%的比例约为32%和41%；网格间距为0.4°×0.4°时，P波和S波检测板扰动值恢复准确性分别达到76%和82%，振幅恢复达到70%的比例约为52%和62%，振幅恢复达到100%的比例约为36%和41%；网格间距为0.5°×0.5°时，P波和S波检测板扰动值恢复准确性分别达到77%和82%，振幅恢复达到70%的比例约为55%和64%，振幅恢复达到100%的比例约为37%和42%.

图6 检测板分辨率测试结果Fig.6 Results of checkerboard resolution tests

3 结果与讨论

经过反演，P波和S波走时残差均方根由反演前的0.620 s和0.680 s分别降低到0.378 s和0.417 s，减少近40%.利用P波、S波反演结果计算得到了泊松比成像结果，如图7所示.苏道磊等（2016）的P波层析成像结果揭示出：在1——10 km深度切片上，沂沭断裂带沿构造走向表现出强烈的横向不均一性，高低速异常交替出现，胶东半岛北部海域、胶莱盆地和济阳坳陷主要表现为低速异常；20——30 km深度，鲁西地区存在较大范围的低速异常（苏道磊等，2016）.上述结果与本文的P波成像结果非常一致（图7a）.另外，本文的S波成像结果与Li等（2018）基本一致，如在10 km深度上，沂沭断裂带北部以低速异常为主（图7b），20——30 km山东半岛表现为大范围的低速异常（图7a）.但由于本文应用了更多的P波和S波到时数据，因此本文的成像结果具有更高的分辨率，如在1 km和10 km深度，沂沭断裂带南部西侧的两个低速异常区相互分离，断裂带东侧表现为高速异常（图7a），这在苏道磊等（2016）的研究中揭示得并不明显.同时，本文的检测板测试结果也要优于苏道磊等（2016），鲁西南和山东半岛地区分辨率能够达到0.33°×0.33°.山东半岛地区中下地壳存在大范围低速异常，说明该地区可能存在较强烈的地幔上涌（李志伟等，2006）.本文同时利用Liu和Zhao （2018）提出的方法计算了1——30 km的P波、S波和泊松比的平均值.平均泊松比异常与郑宏等（2021）利用接收函数反演得到的山东地区的泊松比分布基本一致，进一步说明了本文成像结果的可靠性.

图7 山东地区不同深度 P 波速度扰动（a）、S 波速度扰动（b）及泊松比相对扰动（c）结果底部3幅图分别为1——30 km的P波、S波速度扰动及泊松比相对扰动的平均值Fig.7 Map views of P-wave （a），S-wave （b） and Poisson’s ratio （c） perturbations in Shandong area The averages of P-wave，S-wave and Poisson’s ratio anomalies from 1 to 30 km are shown at the bottom

泰山作为山东地区的最高点，自新生代以来经历了多期快速抬升（李理，钟大赉，2006）.接收函数研究结果显示，泰山地区具有较薄的地壳厚度（32 km左右）和较大的泊松比（0.27），且没有明显的方位各向异性（郑宏等，2021）.本文的成像结果显示，泰山地区下方存在明显的低速异常（图8左下四幅图），这说明该地区仍存在较强的地幔上涌，导致泰山新生代以来出现显著抬升，并且现在仍处于抬升阶段（郑宏等，2021）.另外，在浅部1——10 km深度（图8左上四幅图），泰山北部为高速异常，南部为低速异常，这与泰山的岩性和地质构造是一致的：以泰山山前断裂为界，断裂北侧为泰山主体，主要出露前寒武纪泰山变质杂岩，以隆升为主；山前断裂南侧以沉降为主，主要为泰安——莱芜盆地巨厚的第三系碎屑岩及第四系沉积物（李理，钟大赉，2006）.

图8 长清地震及泰山周边不同深度处 P 波速度扰动（a）、S 波速度扰动（b）及泊松比相对扰动（c）结果Fig.8 Map views of P-wave （a），S-wave （b） and Poisson’s ratio （c） perturbations around Changqing earthquake and Mount Tai

2020年2月18日，济南长清发生了一次MS4.1地震，震源深度约为2.7 km，并引发近40次的余震活动（张斌等，2020）.本文的成像结果显示，该地震震中位于P波、S波高低速异常和高低泊松比异常过渡带，可能与该地震有关的长清断裂也位于P波、S波和泊松比高低异常过渡带（图8上半部和图9）.震中东部的高速异常与该地区的背景噪声成像结果一致，可能与济南侵入岩体有关（雷霆，2020）.震源机制解显示，该地震具有正断兼走滑性质，其余震多为走滑性质（张斌等，2020）.另外，余震重定位结果显示，主震两侧的破裂呈明显不对称分布，以西北侧破裂为主（张斌等，2020）.因此，济南长清MS4.1地震可能是区域构造应力下长清断裂发生左旋走滑运动的结果（张斌等，2020），这与本文的层析成像结果相一致.另外，主震的西北向不对称破裂的形成可能是由于震中东侧高速异常体的存在阻碍了地震的东向破裂.

图9 过长清地震的剖面 AA′（a）的 P 波速度扰动（b）、S 波速度扰动（c）和泊松比扰动（d）Fig.9 Vertical cross section of P-wave （b），S-wave （c） and Poisson’s ratio （d） perturbations along a profile AA′ （a） across the Changqing earthquake shown in the inset map

2003年6月5日青岛崂山地区发生了ML4.1地震（图10a左侧），2004年11月1日附近地区又发生一次ML3.6地震，这两次地震都形成了震群序列.根据定位结果，这两次震群序列大致呈北西走向，与附近的主要断裂近似垂直、与次级断裂展布方向基本一致但并不重合，可能属于相对完整岩体条件下的一次新的破裂（潘元生等，2003，2005）.根据本文的成像结果，这两次地震及震群序列发生在P波高低速异常过渡带、S波低速异常区和高低泊松比过渡区且偏向于高泊松比区域（图10a——c 1km），同时深部存在明显的P波和S波低速异常（图10a，b 10——30 km），这与1995年日本神户地震（ML7.2）震源区的速度和泊松比结构非常相似（Zhaoet al，1996）.另外，震中附近水资源丰富，同时为花岗岩侵入形成地区，且构造裂隙发育（赵广涛等，1996）.因此本文认为，可能是深部的流体填充相对完整岩体内的裂隙并引发破裂，从而导致地震的发生.需要注意的是，虽然这两次地震及震群序列与本地区北东走向的主断层没有直接关系，但本文的成像结果显示，该地区北东走向的主断层基本位于低速区内以及高低泊松比异常的过渡带（图10），因此仍需警惕该地区发生震级更大的中强地震的可能性.

乳山地区位于大别——苏鲁超高压变质带的东部，可能存在大量的中生代岩浆侵入体（郭敬辉等，2005）.自2013年10月开始，记录到的乳山震群地震数已超过1万余次，最大震级为M5.0，目前震群活动可能已结束.震群周边的断裂大多以NE-NNE向为主，而近期的重定位结果和震源机制解显示，乳山地震序列主要呈NW向展布，绝大多数地震发生于中上地壳，发震断层可能为倾角近直立的左旋走滑断层，与区域内的主要断裂并不一致（张斌等，2017）.本文的成像结果显示，乳山震群周边1 km深度以高速异常和高泊松比异常为主，10 km深度主要为高低速异常过渡区，20——30 km深度以低速异常和低泊松比异常为主（图10）.另外，震群附近地热资源较为丰富（田粟，2012）.因此，本文认为，乳山地区深部存在热地幔物质上涌，内部含有的流体注入浅部相对较完整的侵入体内或侵入体之间，在区域应力场的作用下导致侵入体的破裂或侵入体之间隐伏断裂的活动，从而引发乳山震群.这与前人对乳山震群进行的震源谱参数反演、震中空间分布、震中随时间的演化规律以及震源区应力状态的研究结果是一致的（郑建常等，2016；王鹏，2019）.

图10 山东半岛地区不同深度 P 波速度扰动（a）、S 波速度扰动（b）及泊松比相对扰动（c）结果Fig.10 Map views of P-wave （a），S-wave （b） and Poisson’s ratio （c） perturbations in and around Shandong Peninsula

长岛震群的发生可能主要受控于NWW向的张家口——蓬莱断裂，这是一条深大断裂，可能已切穿莫霍面甚至岩石圈，成为地幔热物质或基性物质上涌的通道（张岭等，2007），断裂主要以正断兼走滑运动为主（索艳慧等，2013）.历史上，附近海域曾发生过M6.0和M7.0左右的大地震（王志才等，2006）.本文的成像结果显示（图10），张家口——蓬莱断裂带的地壳速度结构特征与沂沭断裂带具有非常大的相似性（苏道磊等，2016），地壳速度和泊松比结构在1——10 km深度处具有非常强烈的横向不均一性，断裂一侧高低速异常交替分布，断裂带整体位于高低速异常与高低泊松比异常的过渡带，20 km以下以低速和低泊松比异常为主，但部分地区仍有显著的高泊松比异常，可能反映了地幔热物质的上涌.速度和泊松比异常在10 km左右发生明显变化，这与前人得到的b值在8.5 km左右发生转折是一致的（申金超等，2019）.本文认为，深部地幔热物质沿断裂带上涌所产生的构造应力，和/或地幔热物质内部含有的流体沿断裂上涌或侵入裂隙导致了长岛震群、甚至周边地区的强震活动.

4 结论

利用山东及邻区的地震台站记录的P波和S波到时数据反演获得了研究区内高精度的纵横波速度结构和泊松比异常分布形态.研究结果揭示了研究区内的地壳结构具有强烈的横向不均一性.2020年济南长清MS4.1地震可能是区域构造应力下长清断裂发生左旋走滑运动的结果，震中东侧的高速异常体可能阻碍了地震的东向破裂.2003年青岛崂山ML4.1地震可能是由深部的流体填充相对完整岩体内的裂隙并引发破裂所致.崂山震群、乳山震群和长岛震群可能与深部流体有非常强的相关性.