地球物理学报  2020, Vol. 63 Issue (6): 2210-2220   PDF    
震后GPS观测数据揭示的日本MW9.0地震周边地区地幔黏滞性结构垂向变化
陈飞1, 刘泰2, 付广裕1, 佘雅文1,3     
1. 中国地震局地震预测研究所, 北京 100036;
2. 中国地震局地球物理研究所, 北京 100081;
3. 河北地震局, 石家庄 050021
摘要:本文利用大范围的震后GPS数据和黏弹性球形地球位错理论,定量研究了日本MW9.0地震周边地区地幔黏滞性结构的垂向变化.首先结合陆地和海底的GPS观测数据,以及基于球形地球位错理论格林函数和贝叶斯反演方法,反演了该地震的同震滑动分布,发现其最大错动量高达59 m.然后在均一地幔黏滞性结构的假设前提下,确定了震源周边地区地幔黏滞因子的最优解,发现依据该地幔黏滞因子获得的理论远场震后位移和GPS观测结果之间的均方根误差高达0.81 cm,不能解释远场观测结果.为解决上述问题,本文对震中周边地区地幔黏滞性结构沿垂向方向进行分层,建立了一个随深度变化的地幔黏滞性构造模型,然后综合利用远近场的GPS数据对该地区地幔黏滞因子进行反演研究,结果表明,震源周边地区岩石圈弹性层厚度最优解为40 km,40~220 km深度的地幔黏滞因子最优解为6×1018Pa·s,220~670 km深度之间的地幔黏滞因子最优解为1.5×1019Pa·s.上述地幔黏滞性构造使远场的均方根误差降为0.12 cm,仅为利用均一地幔黏滞性构造所得均方根误差值的15%,大大提高了远场模拟结果的准确性.最后,观测值和模拟值之间的均方根误差分析表明,近场震后形变数据主要约束浅层的地幔黏滞性结构,而远场震后形变数据主要约束深部的地幔黏滞性结构.
关键词: 日本MW9.0地震      球形地球位错理论      地幔黏滞因子(黏滞因子)      震后形变      岩石圈弹性层厚度(弹性厚度)     
Variation of the mantle viscosity around the Tohoku-Oki MW9.0 earthquake revealed by post-seismic GPS data
CHEN Fei1, LIU Tai2, FU GuangYu1, SHE YaWen1,3     
1. Institute of Earthquake Forecasting, China Earthquake Administration, Beijing 100036, China;
2. Institute of Geophysics, China Earthquake Administration, Beijing 100081, China;
3. Hebei Earthquake Agency, Shijiazhuang 050021, China
Abstract: In this paper, we research quantitatively the variation of the mantle viscosity around the Tohoku-Oki MW9.0 earthquake by using the post-seismic GPS data that is widely distributed, as well as the viscoelastic spherical earth dislocation theory. At first, combined land and seafloor GPS data, we inverse the co-seismic slip distribution of the Tohoku-Oki MW9.0 earthquake by using Bayesian inversion method as well as the Green's functions of the spherical earth dislocation theory. The maximum of the slip is about 59 m. Then, under the assumption of a uniform model, we inverse the viscosity of mantle around the Tohoku-Oki MW9.0 earthquake. The root mean square error between the calculated post-seismic displacements and the observed ones reaches up to 0.81 cm in the far field. Next, we establish a model of the mantle viscous structure varying with the change of the depth, and inverse the viscosity of mantle in this region by using far-field and near-field GPS data. As a result, the optimal solution of the lithosphere thickness around the epicenter is 40 km. The optimal viscosity of the mantle between 40 km to 220 km is about 6×1018 Pa·s, and changes to be 1.5×1019 Pa·s between 220 km to 670 km. Such a viscous structure of mantle makes the root mean square error between the observations and calculations in the far field to be as small as 0.12 cm, about 15% of the one under the assumption of uniform mantle viscosity. The accuracy of the theoretical result is improved. Finally, the analysis of root mean square error between the observed value and the simulated value indicates that the near-field data of post-seismic deformation mainly constrain the viscous structure of the shallow mantle, while the far-field data of post-seismic deformation mainly constrain the viscous structure of the deep mantle.
Keywords: The Tohoku-Oki MW9.0 earthquake    Spherical earth dislocation theory    Viscosity of the mantle    Post-seismic deformation    Elastic thickness of the lithosphere    
0 引言

2011年3月11日,日本本岛东北部海域发生了9.0级地震,震源周边地区板块汇聚速率高达8 cm·a-1(Taira, 2001).在距离海沟50 km的海底GPS站监测到31 m的同震水平位移(Kido et al., 2011),中国大陆及其周边的GPS测站也监测到了此次地震引起的远场同震形变(杨少敏等, 2011; Zhao et al., 2018).同时学者们基于地震波形数据、大地测量数据等各种观测数据分别给出了此次地震的同震滑动分布模型(Hayes, 2011; Shao et al., 2011; 邵志刚等, 2011; Wei et al., 2012; Pollitz et al., 2011; Diao et al., 2012),为研究该地区震后形变机理奠定了基础.

地震不仅可以引起显著的同震形变,还会产生随时间变化的震后效应(Pollitz, 1997; Wang et al., 2006).迄今为止,大量学者对日本MW9.0地震引起的震后形变进行了深入研究,包括该地震的震后余滑和黏滞性松弛效应(Wang et al., 2012; 刘泰等, 2017; Zhao et al., 2018; Wang et al., 2018).震后初期的研究(Ozawa et al., 2011; Silverii et al., 2014)忽略地幔黏滞性松弛效应的影响,发现震后余滑发生在同震破裂以下的断层深部区域.随着震后时间的推移,Diao等(2014)发现忽略黏滞性松弛的影响会过高的估计震后余滑.依据海底GPS观测数据,Sun等(2014)发现黏滞性松弛效应在震源区占主导地位,同时发现海沟地区的陆向移动主要源自应力松弛.Yamagiwa等(2015)研究了日本MW9.0地震震后2.5年的地表形变,基于弹性-黏弹性层状地球模型计算出震后2.5年内黏滞性松弛效应产生的累积位移量.刘泰等(2017)研究了该强震震后4.5年的地壳形变,分离了震后余滑和黏滞性松弛效应,发现震后余滑在震后初期占主导地位,随后快速减弱,并在2.5年后接近为零.在刘泰等(2017)研究的基础上,Zhao等(2018)研究了该强震在中国东北及其周边区域引起的震后形变,发现在远场的震后形变中,黏滞性松弛效应的影响尤为显著.

刘泰等(2017)基于近场GPS观测数据,反演获得的该区域地幔黏滞因子(黏滞因子)为6.0×1018 Pa·s;Zhao等(2018)基于远场GPS观测数据,得到黏滞因子为1.0×1019 Pa·s.刘泰等(2017)Zhao等(2018)采用完全相同的数据处理策略,却得到了明显不同的黏滞因子,表明均一地幔黏滞性构造无法同时兼顾近场和远场的震后变形.鉴于此,本文对日本MW9.0地震震中周边地区地幔黏滞性结构沿垂向方向进行分层,建立了一个随深度变化而变化的地幔黏滞性构造模型,然后基于大范围的震后GPS数据,利用Tanaka等(2006, 2007)的黏弹性球形地球位错理论,对该地区黏滞因子进行反演研究,定量给出该地震周边地区地幔黏滞性结构垂向变化,该结果对于研究震后变形机理、大地震循环周期等科学问题提供了深部构造模型依据,从而具有重要的科学价值.

1 日本MW9.0地震同震滑动分布模型

本节对日本MW9.0地震同震断层的滑动分布反演所采用的观测数据为日本地理空间信息局发布的日本MW9.0地震的同震数据,其中包含1232个遍布日本列岛的GPS同震观测结果以及5个海底GPS观测站的结果.首先进行棋盘检测,验证反演方法的可行性与GPS观测数据对日本MW9.0地震同震断层滑动分布模型的约束程度.验证过程中,采用如图 1a所示的由450个子断层组成的初始模型, 各子断层长度和宽度分别为25 m和20 m,每相邻的30个子断层的初始滑动量在0和5 m之间交替分布,滑动角统一设为90°,断层几何模型参考刘泰等(2017)断层滑动分布模型建立.利用图 1a所示的断层模型进行正演计算,获取上述1237个测站的同震理论位移,在理论位移的基础上加入2 mm的随机误差,然后对滑动分布进行计算,反演方法均采用贝叶斯反演方法(Zhou et al., 2014).图 1b是利用Okada(1985)的平面半空间位错理论得到位错格林函数,再通过贝叶斯反演方法获得的同震滑动分布模型,其对断层深部区域的滑动量能较为准确的再现,而靠近海沟的浅部断层再现程度不高,其再现的最大滑动量达7.8 m.图 1c是利用球形地球位错理论(Sun et al., 2009)格林函数和贝叶斯反演方法获得的滑动分布模型,它得到的模型滑动量比使用平面半空间位错理论格林函数反演出来的小,最大滑动量为6.7 m,更加接近于初始值,但断层浅部滑动分布再现程度略低于平面半空间位错理论格林函数反演的模型.总之,据图 1可知,靠近日本岛一侧的深部滑动,利用两种格林函数均可以获得较为可靠的结果,但海沟附近的反演结果可靠性相对较差.该棋盘检测表明,本文所使用的GPS数据可以较好地约束断层的同震滑动分布.

图 1 反演方法的棋盘检验 (a)初始滑动分布;(b)使用平面半空间位错理论的格林函数反演的滑动分布;(c)使用球形地球位错理论格林函数反演的滑动分布. Fig. 1 Checkerboard test of inversion method (a) Initial slip distribution; (b) Inversed slips by using Green′s functions of half-space dislocation theory; (c) Inversed slips by using Green′s functions of spherical earth dislocation theory.

接着,考虑到远场的震后形变,采用球形地球位错理论格林函数和贝叶斯反演方法,利用日本地理空间信息局发布1232个陆地GPS和5个海底GPS观测站的同震数据,反演了日本MW9.0地震的同震滑动分布,发现最大错动量为59 m.而远场GPS数据可视为上述模型是否可靠的外部检核条件.基于图 2给出的滑动分布模型,利用球形地球位错理论(Sun et al., 2009),正演出该地震远场的同震水平位移,其结果如图 3所示.正演数据包含中国大陆的东北部及韩国和俄罗斯部分的GPS观测站同震结果.远场GPS站点的同震水平位移计算值与观测值高度匹配,证明利用上述GPS数据获得的同震滑动分布模型在精度上满足本次研究的要求.

图 2 陆地GPS和海底GPS数据联合反演出的日本MW9.0地震同震破裂滑动分布模型 Fig. 2 Co-seismic slip model of the 2011 Tohoku-Oki MW9.0 earthquake inversed by land and seafloor GPS data
图 3 日本MW9.0地震远场GPS站点同震水平位移观测和模拟对比 蓝色箭头代表观测值,红色箭头代表模拟值. Fig. 3 Comparison of the GPS co-seismic horizontal displacements of the 2011 Tohoku-Oki MW9.0 earthquake and the calculated ones in far field The blue and red arrows represent the observed and calculated values, respectively.
2 震后GPS观测数据及其处理过程

本研究的目标是利用GPS观测数据约束日本MW9.0地震周边地区的地幔黏滞性构造,因此有必要选择信噪比较高的观测数据展开研究.总体来说,日本MW9.0地震是一个逆冲型地震,因此垂直于断层面分布的GPS观测数据信号较大(刘泰等,2017),信噪比较高.鉴于此,近场GPS数据的选取范围限定为35°N—41°N/138°E—144°E,共计430个GPS连续站点,数据来自美国内华达大学大地测量实验室网站.首先剔除28个震前2年内没有数据的点.接着剔除受其他地震影响较大的数据,如图 4所示,2011年6月22日在39.241°N,142.463°E处发生MW6.7的地震(蓝色五角星所示),2015年5月12日在38.906°N,142.032°E处发生MW6.8的地震(红色五角星所示),上述两地震周边数据予以剔除.以GPS站点J546为例(图 4中蓝色三角形),该观测站东西向的形变数据在2015年5月12日出现1.5 cm左右的阶跃,表明该观测站受到2015年MW6.8地震的影响较大,为了减少其他地震对本研究的影响,我们剔除距离这两次地震较近的GPS数据.但因为数据信号的问题仍需保留垂直于断层面的一些GPS站点,这些GPS站点受其他地震的影响相对较小.然后,因为海底GPS测站不是连续站点,不满足本文对震后形变研究的要求,予以剔除.其次,剔除日本MW9.0地震断层南北两端的观测数据,因为断层两端的GPS站点对日本MW9.0地震断层模型以及周边地区黏滞性结构的约束较低.最后,选取的震后形变观测站尽可能均匀分布,确保反演结果的稳定性.

图 4 用于黏滞因子反演研究的GPS观测站分布图 蓝色和红色五角星代表两个6级以上地震震中位置;黑色五角星为日本MW 9.0地震震中;蓝色框显示同震滑动分布模型,红色和蓝色三角形为GPS观测站. Fig. 4 Distribution of GPS stations that are used to inverse the viscosity of the mantle Blue and red stars represent the epicenter of two earthquakes (M>6). Black star represents the epicenter of the Tohoku-Oki MW9.0 earthquake. Blue frame shows the co-seismic slip distribution model of the great earthquake. Red and blue triangles represent GPS stations.

经过上述原则的筛选,本研究最终选取分布于日本本岛的23个信噪比较高的GPS观测站的震后观测数据,作为研究震源及周边地区地幔黏滞性结构的基础数据,此外,本研究还选取8个远场GPS观测站的观测结果(Zhao et al., 2018)进行研究,其详细分布见图 4.所使用的数据为震后2.5~4.5年期间观测到的水平位移.由于震后2.5年震后余滑已经基本停止,黏滞性松弛开始占据主导作用(刘泰等,2017),上述GPS观测数据不包含震后余滑的影响,适宜于进行黏滞因子反演研究.

在GPS数据处理过程中,首先扣除震前趋势项对震后形变的影响,之后对震后形变曲线进行拟合,拟合公式(1)来源于Heinkelmann等(2008),公式中包含震后余滑项和黏滞性松弛项.利用公式(1)进行拟合提取震后2.5年和4.5年的数据,得到震后水平位移.

(1)

式中,a1a2分别为对数和指数部分的振幅,δt为震后时间,τln为震后余滑的持续时间,τe为地幔黏滞性松弛的时间.

3 日本MW9.0地震周边地区均一地幔黏滞性结构

随着地幔深度的变化,黏滞因子也在不断变化.但是,依据GPS等大地测量观测数据,通过位错理论进行大地震周边地区黏滞因子的反演研究时,一般要求震后GPS观测数据的空间分布范围较大,信噪比较高,否则就无法约束较为深部的黏滞因子,因此,定量反演黏滞因子随深度的变化是一个较为困难的科学问题.鉴于上述原因,研究者一般假设地幔黏滞性结构是均一的,即不考虑地幔黏滞性随深度的变化,然后在此基础上进行定量的反演研究(Yamagiwa et al., 2015; Zhao et al., 2018; 刘泰等, 2017).本节将依照Zhao等(2018)刘泰等(2017)的研究思路,利用上文挑选出的GPS震后数据,研究日本MW9.0地震周边地区均一的地幔黏滞性构造.

本节采用黏弹性球形地球位错理论(Tanaka et al., 2006, 2007)模拟日本MW9.0地震引起的震后形变.模拟过程中,依次选取不同的岩石圈弹性层厚度(弹性厚度)和黏滞因子,计算相应参数组合对应的位错Love数,进而得到对应的位移格林函数;然后利用位移积分计算程序(付广裕和刘泰,2017)计算各个时间段对应的震后水平位移;接着利用均方根误差公式(2)计算不同模型组合对应的模拟值与GPS观测值之间的均方根误差;最后进行网格搜索,从统计学的角度定量估算模拟值和观测值之间的不符合程度,从而获取震源及周边地区黏滞因子的最优解.

(2)

其中,Uobsi为GPS观测值,Ucali为理论计算值,N为GPS观测站点数.

本节以近场(日本列岛)的GPS震后观测数据为基础展开反演研究,使用远场(中国东北、韩国和俄罗斯)的观测数据对反演结果进行检验,即利用远场的震后观测结果对仅利用近场数据反演出的均一的黏滞因子模型进行检验,探讨使用均一地幔黏滞性构造的有效性.在对震后水平位移进行模拟的过程中,弹性厚度的搜索区域设定为20~50 km,步长设定为10 km,黏滞因子的搜索区域设定为(3~10)×1018Pa·s,步长设定为1×1018Pa·s.由于正演计算量较大的缘故,所以选取了上述较大的步长.图 5给出了弹性厚度与黏滞因子不同组合对应的均方根误差曲线.反复拟合后发现,当弹性厚度取为40 km,黏滞因子取为7×1018Pa·s时,对应的均方根误差值为3.1 cm,是上述所有组合中的最优解,模拟值和GPS观测值符合程度最高,此时的反演结果最优.其次,当黏滞因子为4×1018Pa·s,弹性厚度为30 km时,均方根误差值为3.6 cm(图 5)是所有组合中第二小的值,也是一个相对比较理想的结果.从图 5中可以看出当弹性厚度为30 km,黏滞因子为4×1018Pa·s和5×1018Pa·s时,两点的均方根误差几乎相同,因此可能存在更小值,所以我们将步长缩小一半,计算出黏滞因子为4.5×1018 Pa·s时,均方根误差大小为3.4 cm,没有对最优解的选取造成影响.

图 5 日本MW9.0地震周边地区不同黏滞因子和弹性厚度组合对应的均方根误差 Fig. 5 The root mean square (RMS) error, based on different viscosity of the mantle and different elastic thickness around the region of the Tohoku-Oki MW9.0 earthquake

最后使用GPS远场观测数据,对上述网格搜索结果进行检核.即设定日本MW9.0地震周边地区的弹性厚度和黏滞因子分别为40 km和7×1018Pa·s,然后进行正演计算,获取远场的8个GPS观测站的理论震后水平位移,最后与观测结果进行比较.图 6a给出了日本列岛的震后观测值和模拟值,二者吻合度非常高,图 6b显示的残差较小,图 6c6d分别给出远场GPS观测与模拟结果,以及它们的残差分布.图 6清晰表明,二者之间的残差均方根高达0.81 cm,表明远场模拟值无法解释观测结果.该结果表明,远场和近场的GPS观测数据,可能由不同深度的地幔黏滞性分别控制.为了能够同时解释远场观测结果,下一节将对地幔黏滞性进行分层处理.至于分层构造下的反演结果检核,可先只拿一部分的近场数据进行反演,再利用未参与反演的另一部分近场数据进行验证,检验反演结果的正确性.远场亦是如此.

图 6 日本MW9.0地震周边地区均一地幔黏滞性结构假设前提下,GPS观测值(蓝色箭头)与最优模拟值(红色箭头),以及它们之间的残差分布 Fig. 6 The post-seismic horizontal displacements observed by GPS (blue arrows), the optimal calculated values (red arrows) under the assumption of uniform mantle viscous structure around the region of the 2011MW9.0 Tohoku-Oki earthquake, and the residuals
4 日本MW9.0地震周边地区地幔黏滞性结构的垂向变化

由于日本MW9.0地震震级足够大,不仅在日本列岛引起了显著的震后位移(刘泰等,2017),而且在中国东北、韩国和俄罗斯引起了可被GPS观测到的、信噪比较高的震后水平位移(Zhao et al., 2018),为研究日本MW9.0地震周边地区地幔黏滞性结构的垂向变化提供了可能.而刘泰等(2017)Zhao等(2018)给出的日本MW9.0地震周边地区地幔黏滞性结构的显著差异,也表明了研究该地震周边地区地幔黏滞性构造随深度变化的必要性.

为了同时解释GPS观测到的远场和近场震后位移数据,我们拟对日本MW9.0地震周边地区的地幔黏滞性结构进行垂向分层.考虑到GPS观测站位置的限制(图 4),只有日本列岛和欧亚大陆分布有GPS观测站,中间的日本海并没有GPS观测站,所以在进行垂向黏滞性构造反演时,无法反演出一个随地幔深度变化而变化的连续曲线,鉴于此,我们构造出一个阶梯型黏滞性结构模型.Pollitz等(2006)在研究苏门答腊MW9.3地震时,依据PREM地球模型(Dziewonski and Anderson, 1975)的垂向密度结构变化将整个地幔分成三层,220 km处为软流层的底部,670 km处为上下地幔交界处,同时这两处上下的密度也不连续.本研究将该模型作为地幔黏滞性结构模型划分的依据.此外,Pollitz等(2006)将670 km至核幔边界间即下地幔的黏滞因子设为1×1021Pa·s,因为图 5中的GPS观测数据对下地幔的黏滞因子基本没有影响,至于更远处的GPS数据又因为震后信号相对较小,信噪比低,不宜于对该层的黏滞因子展开反演研究,所以下地幔的黏滞因子采用Pollitz等(2006)的研究结果.至于上部地幔的黏滞因子,因为研究的区域不同而不同,可依据观测数据进行反演获得.参考第3节反演获得的均一地幔黏滞性构造,我们将日本MW9.0地震周边地区弹性厚度设为40 km,40~220 km和220~670 km两层的黏滞因子的初始值均设定为7×1018Pa·s,670 km至核幔边界间的黏滞因子根据Pollitz等(2006)设为1×1021Pa·s,结合上述选取结果作为地幔黏滞性结构垂向变化的初始模型.

利用上述模型对近场(日本列岛)的震后位移进行研究,模拟震后2.5~4.5年内的水平位移,将其和GPS观测结果进行比较,寻求最佳拟合解.如上文所述,均一地幔模型下黏滞因子的最优解为7×1018Pa·s,考虑到黏滞因子的分层结构后,本节将地幔深度为670 km以下介质的黏滞因子增大至1×1021 Pa·s,同时在保证对近场震后形变的模拟改变最小的情况下,需将670 km深度以上介质的黏滞因子微调至6×1018Pa·s, 并以此作为网格搜索初值,调整220~670 km深度范围内的黏滞因子.本文利用近场GPS数据主要控制浅层黏滞构造,利用远场GPS数据微调浅层黏滞构造和主要控制深层黏滞构造的研究思路,将处于40~220 km深度范围内的浅层黏滞因子减小至6×1018Pa·s,然后利用远场的震后形变对深部黏滞因子即深度为220~670 km处的黏滞因子进行调整,将步长设定为1×1018Pa·s,调整220~670 km深度范围内的黏滞因子使其从6×1018Pa·s依次增加,检验近场和远场震后2.5~4.5年震后位移场的观测值和模拟值的吻合程度.此处有必要说明的是,即使是正演计算,工作量与计算机耗时都比较大,不能展开全区域小步长搜索反演,而只能采用上述较大的步长进行拟合研究.

随着深度为220~670 km处的黏滞因子的不断增加,远场震后2.5~4.5年水平位移场的模拟值在不断减小,当其增加至1.5×1019Pa·s时,观测和模拟之间的均方根误差值为0.12 cm,是远场所有组合中的最优解,随后均方根误差又逐步变大(图 7b).这一现象表明远场的GPS数据对深部的地幔黏滞性结构有较好的约束效果.而近场的模拟值随着220~670 km处深部黏滞因子的不断增加改变很小(图 7a),总体呈现先减小后增加的趋势.黏滞因子在(8~15)×1018Pa·s之间变化时,近场GPS观测值与拟合值之间的均方根误差曲线呈现波动状态,这是因为220~670 km处的地幔是较为深层的地幔,深层地幔的黏滞性结构对近场震后形变控制能力有限,不能获取更高的分辨率,在观测误差的影响下,无法单调约束近场的均方根误差至最优解.图 7a中绿色曲线为均方根误差拟合值,可以看出如果没有GPS观测误差影响,近场的均方根误差曲线可以约束到最优解,证明了近场的震后形变受深层的地幔黏滞性结构影响较低.上述均方根误差分析表明,近场的GPS观测数据主要约束浅层地幔黏滞性结构,远场的GPS观测数据主要约束深部的地幔黏滞性结构.

图 7 日本MW9.0地震近场和远场GPS观测与拟合值之间的均方根误差曲线图 红色直线代表均方根误差值为3.6 cm,绿色虚线代表拟合结果.η40-220η670-2887分别代表 40~220 km深度和670~2887 km深度的黏滞因子. Fig. 7 Far-field and near-field Root Mean Square (RMS) error between the observed and calculated results of the Tohoku-Oki MW9.0 earthquake The red line represents the RMS of 3.6 cm. The dashed line in green represents the fitting result. η40-220 and η670-2887 represent respectively the viscosity of the mantle in depth between 40 km to 220 km, and between 670 km and 2887 km.

均一地幔黏滞性结构对应的均方根误差(图 5)表明,当近场的均方根误差值小于3.6 cm时,其和最优解3.1 cm最为接近,对应的模拟值可以较好的解释近场观测结果.该结果表明,图 7a中的均方根误差值低于红线即小于3.6 cm时,对应的黏滞因子都是可以接受的结果.换句话说,因为深部的地幔黏滞性结构不能将近场的观测值单调约束到最优,所以在此范围内任何一个结果都是可能的.鉴于此,我们依据远场震后形变的均方根误差进行分析,发现黏滞因子为1.5×1019Pa·s时均方根误差为最优解为0.12 cm,于是确定1.5×1019Pa·s为220~670 km处的黏滞因子最优解.最终,日本MW9.0地震周边地区地幔黏滞性结构垂向变化模型如下(图 8):0~40 km深度的地层设定为岩石圈弹性层,40~220 km深度的黏滞因子设定为6×1018Pa·s,其主要由近场GPS数据确定;220~670 km深度的黏滞因子设定为1.5×1019Pa·s,其主要由远场GPS数据确定;670 km深度至核幔边界之间地层的黏滞因子为1×1021Pa·s,该地层的黏滞因子依据Pollitz等(2006)的研究结果确定.

图 8 (a) 日本MW9.0地震周边地区地幔黏滞性结构垂向变化图,(b)为(a)中蓝色矩形框部分的地幔黏滞性结构(放大图) Fig. 8 The variation of the mantle viscous structure around the Tohoku-Oki MW9.0 earthquake. (b) A larger image shows the variation of the viscous structure of the mantle in the blue rectangular of (a)

图 9为在上述黏滞性结构模型条件下,远场和近场震后水平位移的模拟值和观测值对比结果,以及它们之间的残差分布.对比图 6可清晰发现,新模型下,近场GPS站点的观测值与模拟值之间的残差变化不大,但远场的观测值与模拟值之间的残差明显变小,均方根误差降至0.12 cm.该结果表明,黏滞因子垂向变化模型可以同时解释日本MW9.0地震引起的远场和近场震后位移.

图 9 地幔黏滞性结构垂向变化模型条件下,远近场模拟值最优解和GPS观测值对比图,以及它们之间的残差分布图 Fig. 9 The near-field and far-field GPS observations, the optimal calculated values under the condition of vertical variation model of the mantle viscous structure, and the residuals

Diao等(2014)参考Pollitz等(2006)Han等(2008)中的地幔黏滞性垂向结构,将日本MW9.0地震震源周边区域深度为220~670 km的地幔黏滞因子设定为1.0×1020Pa·s,利用近场的数据对浅层的地幔黏滞性结构进行约束,得到的弹性厚度为50 km,50~220 km深度的黏滞因子为8×1018Pa·s.梁明等(2018)利用平面半空间位错理论和近场GPS数据,反演了浅层的地幔黏滞性结构,利用GRACE数据修正深部的黏滞因子,得到深度40~120 km之间地层的黏滞因子大陆侧地幔顶层为1.5×1019Pa·s,海洋侧为6×1018Pa·s;120~220 km和220 km以下的黏滞因子分别为5×1019Pa·s和1.0×1021Pa·s.Hu等(2016)利用三维黏弹性有限元模型和日本岛的GPS观测数据对该地区的黏滞因子进行了研究,因为陆地上水平向海运动主要受地幔楔黏滞因子的控制,他们将地幔层分为海洋地幔和地幔楔,使用的三维模型中包含弹性层、黏弹性地幔层、弹性俯冲板和黏弹性剪切带等.Diao等(2014)梁明等(2018)Hu等(2016)和本文的研究表明,地幔黏滞性构造确实存在垂向变化,Zhao等(2018)刘泰等(2017)给出的均一地幔黏滞性构造过于简单.测震、形变与模拟研究结果表明,实际地球的黏滞因子是一个随深度变化而不断变化的物理量(Tanaka et al., 2007Dziewonski and Anderson, 1975石耀霖与曹建玲, 2008).本文依据日本MW9.0地震引起的大范围震后变形数据,反演了震中周边地区的地幔黏滞性垂向变化.虽然由于数据分布的局限(日本海缺乏GPS观测数据),反演结果依然简单,却利用大地测量数据证明了黏滞因子随深度变化这一事实.

最后,根据本文结果表明,均一地幔黏滞性结构不能同时解释远场和近场的震后形变.换句话说,GPS观测到的大范围震后位移数据,可用来约束日本MW9.0地震周边地区的地幔黏滞性随深度的变化.此外,位错理论的完备程度也会影响黏滞因子的反演精度.本文利用球形地球模型位错理论,依据大范围的GPS震后数据展开研究,结果应相对可靠.

5 结论

首先,本文基于同震GPS观测数据,利用球形地球位错理论格林函数和贝叶斯反演方法(Zhou et al., 2014),然后,利用第1节反演的断层滑动分布模型,通过黏弹性球形地球位错理论(Tanaka et al., 2007),正演计算出分布于日本列岛的23个GPS观测点、以及分布于中国东北、韩国和俄罗斯的8个GPS观测站的震后2.5~4.5年的水平位移,通过震后GPS观测值和模拟值的反复拟合,获取日本MW9.0地震周边地区最优的地幔黏滞因子(黏滞因子).在均一地幔黏滞性结构的假设前提下,当岩石圈弹性层厚度(弹性厚度)和黏滞因子分别为40 km和7×1018Pa·s时,近场的观测和模拟吻合程度最高,二者之间的均方根误差约为3.1 cm;但远场拟合值系统性高于GPS观测结果,均方根误差高达0.81 cm.随后,采用地幔黏滞性结构的垂向分层方法,对整个地幔细分为三层:深度为0~40 km之间的地层为弹性层;深度为40~220 km之间的地层为浅层地幔,其黏滞因子主要由近场GPS数据约束,最优解为6×1018Pa·s;深度为220~670 km之间的地层为深层地幔,其黏滞因子主要由远场GPS数据约束,最优解为1.5×1019Pa·s;670 km深度至核幔边界之间黏滞因子为1×1021Pa·s.利用上述分层模型得到的近场均方根误差为3.5 cm,接近均一模型下的最优解;远场的均方根误差大约为0.12 cm,明显优于均一模型的模拟结果.最后,根据观测值和模拟值之间的均方根误差分析表明,近场GPS观测数据主要约束浅层地幔黏滞性结构,远场GPS数据主要约束深层地幔黏滞性结构.

本研究给出的地幔黏滞性分层结构不仅可以解释近场震后GPS观测数据,也可以解释远场观测数据,对研究日本MW9.0地震周边地区大范围的震后形变分布形态及其形成机理具有重要作用.

致谢  本文采用了日本地理空间信息局发布的同震GPS数据以及美国内华达州大地测量实验室发布的震后GPS数据,使用了东京大学Tanaka博士提供的计算软件,同时两位匿名审稿老师的建议对提高论文质量有很大帮助,在此一并致谢.
References
Diao F Q, Xiong X, Yong Z. 2012. Static slip model of the MW9.0 Tohoku (Japan) earthquake:Results from joint inversion of terrestrial GPS data and seafloor GPS/acoustic data. Chinese Science Bulletin, 57(16): 1990-1997. DOI:10.1007/s11434-012-5014-5
Diao F Q, Xiong X, Wang R J, et al. 2014. Overlapping post-seismic deformation processes:Afterslip and viscoelastic relaxation following the 2011 MW9.0 Tohoku (Japan) earthquake. Geophysical Journal International, 196(1): 218-229. DOI:10.1093/gji/ggt376
Dziewonski A M, Anderson D L. 1981. Preliminary reference earth model. Physics of the Earth and Planetary Interiors, 25(4): 297-356. DOI:10.1016/0031-9201(81)90046-7
Fu G Y, Liu T. 2017. A user-friendly code for calculating post-seismic displacements and gravity changes on a symmetric viscoelastic spherical earth. Journal of Geodesy and Geodynamics (in Chinese), 37(7): 661-667.
Han S C, Sauber J, Luthcke S B, et al. 2008. Implications of postseismic gravity change following the great 2004 Sumatra-Andaman earthquake from the regional harmonic analysis of GRACE intersatellite tracking data. Journal of Geophysical Research, 113(B11).
Hayes G P. 2011. Rapid source characterization of the 2011 MW9.0 off the Pacific coast of Tohoku Earthquake. Earth, Planets and Space, 63(7): 529-534. DOI:10.5047/eps.2011.05.012
Heinkelmann R, Freymueller J, Schuh H. 2008. A postseismic relaxation model for the 2002 Denali earthquake from GPS deformation analysis applied to VLBI data. Ivs General Meeting Proceeding: 340-355.
Hu Y, Bürgmann R, Uchida N, et al. 2016. Stress-driven relaxation of heterogeneous upper mantle and time-dependent afterslip following the 2011 Tohoku earthquake. Journal of Geophysical Research:Solid Earth, 121(1): 385-411. DOI:10.1002/2015JB012508
Kido M, Osada Y, Fujimoto H, et al. 2011. Trench-normal variation in observed seafloor displacements associated with the 2011 Tohoku-Oki earthquake. Geophysical Research Letters, 38(24): 389-395.
Liang M, Wang W X, Zhang J. 2018. Post-seismic deformation mechanism of the MW9.0 Tohoku-Oki earthquake detected by GPS and GRACE observations. Chinese Journal of Geophysics (in Chinese), 61(7): 2691-2704. DOI:10.6038/cjg2018L0356
Liu T, Fu G Y, Zhou X, et al. 2017. Mechanism of post-seismic deformations following the 2011 Tohoku-Oki MW9.0 earthquake and general structure of lithosphere around the source. Chinese Journal of Geophysics (in Chinese), 60(9): 3406-3417. DOI:10.6038/cjg20170911
Okada Y. 1985. Surface deformation due to shear and tensile faults in a half-space. Bulletin of the Seismological Society of America, 75(4): 1135-1154.
Ozawa S, Nishimura T, Suito H, et al. 2011. Coseismic and postseismic slip of the 2011 magnitude-9 Tohoku-Oki earthquake. Nature, 475(7356): 373-376. DOI:10.1038/nature10227
Pollitz F F. 1997. Gravitational viscoelastic postseismic relaxation on a layered spherical Earth. Journal of Geophysical Research:Solid Earth, 102(B8): 17921-17941. DOI:10.1029/97JB01277
Pollitz F F, Bürgmann R, Banerjee P. 2006. Post-seismic relaxation following the great 2004 Sumatra-Andaman earthquake on a compressible self-gravitating earth. Geophysical Journal International, 167(1): 397-420. DOI:10.1111/j.1365-246X.2006.03018.x
Pollitz F F, Bürgmann R, Banerjee P. 2011. Geodetic slip model of the 2011 M9.0 Tohoku earthquake. Geophysical Research Letters, 38(7): L00G08. DOI:10.1029/2011GL048632
Shao G F, Li X Y, Ji C, et al. 2011. Focal mechanism and slip history of the 2011 MW9.1 off the Pacific coast of Tohoku Earthquake, constrained with teleseismic body and surface waves. Earth, Planets and Space, 63(7): 559-564. DOI:10.5047/eps.2011.06.028
Shao Z G, Wu Y Q, Jiang Z S, et al. 2011. The analysis of coseismic slip and near-field deformation about Japanese 9.0 earthquake based on the GPS observation. Chinese Journal of Geophysics (in Chinese), 54(9): 2243-2249. DOI:10.3969/j.issn.0001-5733.2011.09.006
Shi Y L, Cao J L. 2008. Effective viscosity of China continental lithosphere. Earth Science Frontiers (in Chinese), 15(3): 82-95. DOI:10.1016/S1872-5791(08)60064-0
Silverii F, Cheloni D, D'Agostino N, et al. 2014. Post-seismic slip of the 2011 Tohoku-Oki earthquake from GPS observations:implications for depth-dependent properties of subduction megathrusts. Geophysical Journal International, 198(1): 580-596. DOI:10.1093/gji/ggu149
Sun T, Wang K L, Iinuma T, et al. 2014. Prevalence of viscoelastic relaxation after the 2011 Tohoku-Oki earthquake. Nature, 514(7520): 84-87. DOI:10.1038/nature13778
Sun W K, Okubo S, Fu G Y, et al. 2009. General formulations of global co-seismic deformations caused by an arbitrary dislocation in a spherically symmetric earth model-Applicable to deformed earth surface and space-fixed point. Geophysical Journal International, 177(3): 817-833. DOI:10.1111/j.1365-246X.2009.04113.x
Taira K. 2001. Late Quaternary paleoceanography in eastern Asia. Journal of Coastal Research, 17(1): 114-117.
Tanaka T, Okuno J, Okubo S. 2006. A new method for the computation of global viscoelastic post-seismic deformation in a realistic earth model(Ⅰ)-vertical displacement and gravity variation. Geophysical Journal International, 164(2): 273-289. DOI:10.1111/j.1365-246X.2005.02821.x
Tanaka T, Okuno J, Okubo S. 2007. A new method for the computation of global viscoelastic post-seismic deformation in a realistic earth model (Ⅱ)-Horizontal displacement. Geophysical Journal International, 170(3): 1031-1052. DOI:10.1111/j.1365-246X.2007.03486.x
Wang K L, Hu Y, He J H. 2012. Deformation cycles of subduction earthquakes in a viscoelastic Earth. Nature, 484(7394): 327-332. DOI:10.1038/nature11032
Wang K L, Sun T, Brown L, et al. 2018. Learning from crustal deformation associated with the M9 2011 Tohoku-Oki earthquake. Geosphere, 14(2): 552-571. DOI:10.1130/GES01531.1
Wang R J, Lorenzo-Martín F, Roth F. 2006. PSGRN/PSCMP-a new code for calculating co-and post-seismic deformation, geoid and gravity changes based on the viscoelastic-gravitational dislocation theory. Computers & Geosciences, 32(4): 527-541.
Wei S J, Graves R, Helmberger D, et al. 2012. Sources of shaking and flooding during the Tohoku-Oki earthquake:A mixture of rupture styles. Earth and Planetary Science Letters, 333-334: 91-100. DOI:10.1016/j.epsl.2012.04.006
Yamagiwa S, Miyazaki S, Hirahara K, et al. 2015. Afterslip and viscoelastic relaxation following the 2011 Tohoku-Oki earthquake (MW9.0) inferred from inland GPS and seafloor GPS/Acoustic data. Geophysical Research Letters, 42: 66-73. DOI:10.1002/2014GL061735
Yang S M, Nie Z S, Jia Z G, et al. 2011. Far-field coseismic surface displacement caused by the MW9.0 Tohoku earthquake. Geomatics and Information Science of Wuhan University (in Chinese), 36(11): 1336-1339.
Zhao Q, Fu G Y, Wu W W, et al. 2018. Spatial-temporal evolution and corresponding mechanism of the far-field post-seismic displacements following the 2011 MW9.0 Tohoku earthquake. Geophysical Journal International, 214(3): 1774-1782. DOI:10.1093/gji/ggy226
Zhou X, Cambiotti G, Sun W K, et al. 2014. The coseismic slip distribution of a shallow subduction fault constrained by prior information:the example of 2011 Tohoku (MW9.0) megathrust earthquake. Geophysical Journal International, 199(2): 981-995. DOI:10.1093/gji/ggu310
付广裕, 刘泰. 2017. 基于粘弹性球体地球模型的震后位移与重力变化计算软件. 大地测量与地球动力学, 37(7): 661-667.
梁明, 王武星, 张晶, 等. 2018. 联合GPS和GRACE观测研究日本MW9.0地震震后变形机制. 地球物理学报, 61(7): 2691-2704. DOI:10.6038/cjg2018L0356
刘泰, 付广裕, 周新, 等. 2017. 2011年日本MW9.0地震震后形变机制与震源区总体构造特征. 地球物理学报, 60(9): 3406-3417. DOI:10.6038/cjg20170911
邵志刚, 武艳强, 江在森, 等. 2011. 基于GPS观测分析日本9.0级地震同震位错与近场形变特征. 地球物理学报, 54(9): 2243-2249. DOI:10.3969/j.issn.0001-5733.2011.09.006
石耀霖, 曹建玲. 2008. 中国大陆岩石圈等效粘滞系数的计算和讨论. 地学前缘, 15(3): 82-95. DOI:10.3321/j.issn:1005-2321.2008.03.006
杨少敏, 聂兆生, 贾志革, 等. 2011. GPS解算的日本MW9.0级地震的远场同震地表位移. 武汉大学学报, 36(11): 1336-1339.