2024, Vol. 44
2. 武汉大学地球空间环境与大地测量教育部重点实验室，武汉市珞喻路129号，430079;
3. 自然资源部地球物理大地测量重点实验室，武汉市珞喻路129号，430079;
4. 中国地震局地震研究所，武汉市洪山侧路40号，430071
研究大地震的震后时变机制，对认识地震破裂与发震环境的关系，探索地震断层的震后运动行为、摩擦性质及岩石圈流变结构等具有重要意义[1]。地震断层发生同震错动后，断层面的震后余滑会导致地壳产生形变，许多地震的早期震后形变可以用余滑模型来解释。目前用于解释震后余滑形变的模型包括运动学余滑模型和基于速率-状态摩擦准则的应力驱动余滑模型。2021-05-22青海玛多发生MW7.4地震，震中为98.34°E、34.59°N，震源深度约17 km。余震精定位结果显示，地震发生在巴颜喀拉块体北部的NW-W向昆仑山口-江错断裂带[2-3]。布设在玛多地区的GPS连续运行参考站记录了地震的震后高时空分辨率形变信息[4-5]。学者们利用大地测量观测数据，通过运动学余滑模型反演震后余滑[3, 6-7]，同时利用应力驱动余滑模型解释震后形变时空特征[7]。本文采用玛多地震震后160 d的GPS形变时间序列数据，利用基于速率-状态摩擦准则的应力驱动余滑模型[8]，正演此次玛多地震的震后余滑分布，以分析该地震的震后余滑特征及断层摩擦属性。
1 数据和方法 1.1 震后位移和断层模型玛多地震后，学者们使用不同的观测数据给出各自的断层滑动分布模型[5-6]，所得滑移大小和地震矩基本相同。震后早期余滑(约160 d)形变由位于发震断层约50 km内的16个连续GPS (cGPS)站记录得到，其中13个站点在主震后10 d内建立，另外3个在2021年之前建立(图 1)。本文使用的数据为Xiong等[5]剔除16个cGPS站点数据的震间形变场和季节性变形后得到的震后早期形变时间序列，形变表现出对称特征，符合高倾角断层及以走滑为主的滑动机制。值得注意的是，4454站点仅记录到微弱的震后位移，表明该站点附近几乎没有发生余滑，这一结论同样体现在InSAR观测和运动学反演结果中[6-7]。
|
图 1 GPS参考站分布 Fig. 1 Distribution of GPS reference stations |
本文采用Jin等[9]利用InSAR观测数据反演的玛多地震同震滑动模型作为同震输入模型，通过应力驱动余滑模型拟合Xiong等[5]处理得到的震后早期形变时间序列数据。
1.2 震后余滑模型使用Relax程序[10]模拟玛多地震震后余滑引起的地表形变。应力驱动余滑模型描述的是地震断层在受同震应力加载和断层面摩擦参数控制下的震后滑动行为。根据速率-状态摩擦准则，以同震破裂在断层面上产生的剪应力变化为应力初始条件，当给定断层的摩擦参数(a－b)σ和参考速率V0时，可计算断层面的震后余滑和地表形变随时间的演化。与运动学模型相比，应力驱动余滑模型顾及了同震滑动产生的应力变化和断层面的摩擦参数对震后余滑行为的控制作用。在应力驱动余滑模型中，同震库仑应力的增加促进了断层继续滑动[11]。断层滑动速率可以表示为：
| V=2V0sinh[Δτ(a−b)σ] | (1) |
式中，V为断层的震后余滑滑动速率，Δτ为同震剪应力变化，V0为参考滑动速率，(a－b)σ为0.1~10 MPa范围内的本构参数，σ为断层面上的有效正应力。
本文取剪切模量为30 GPa，泊松比为0.25，摩擦系数为0.6。采用Barbot等[12]的方法计算应力变化的同震滑动模型，即假设同震滑动量较大的区域是速度弱化区域，不参与震后蠕滑。为使同震滑动和余滑之间的重叠区域允许发生余滑，去除同震滑动峰值30%以下的滑动量，以集中同震滑动分布；然后将其余的滑动量放大，以保持地震矩不变[12-13]。
实际操作中，由于V0和(a－b)σ之间通常存在很强的相关性，不可能独立确定。尽管改变(a－b)σ的值会导致搜索结果不同，但并不能通过数据拟合来区分不同的(a－b)σ的值[14]，这种相关性正是Barbot等[12]认为不能赋予V0任何特定物理意义的部分原因。因此，本文遵循以往使用Relax软件的方法，将(a－b)σ的值固定为5 MPa，以降低搜索的复杂程度[11]。
2 结果分析震后形变机制主要包括同震应力驱动的震后余滑、下地壳和上地幔粘弹性松弛和地壳浅部的孔隙弹性回弹。为评估除震后余滑外的震后形变机制是否参与了震后早期形变，参考Huang等[15]研究巴颜喀拉块体内断层相互作用时使用的下地壳和上地幔粘滞系数(下地壳1018~1019 Pa·s，上地幔1020 Pa·s)，使用Relax软件计算玛多地震震后160 d粘弹性松弛引起的累积地表形变，预测水平位移小于6 mm，远小于GPS观测值。采用Jin等[9]的同震滑动模型，利用地震前后泊松比的变化计算地表形变场，并对其进行差分，得到平衡状态的孔隙弹性响应形变，预测水平位移不足1 mm，表明孔隙回弹引起的形变可忽略不计。
2.1 余滑模型1假定震后断层余滑面与同震断层滑动面完全重合。由于玛多地震以走滑分量为主，震后GPS形变观测数据的EW分量整体远大于NS分量，故拟合GPS观测数据时仅使用EW分量计算残差，而用NS分量作定性比较。通过改变V0的值对地表震后EW向形变时间序列进行拟合，得到最优拟合参数为V0=1.5 m/a，RMS=0.76 cm，从而得到GPS站点震后形变时间序列(图 2和3)和余滑分布(图 4(a))。可以看出，震后余滑主要分布在同震断层显著滑动的下倾和两侧区域，最大滑移量为0.67 m。余滑模型1所得的地表形变在距离断层迹线较近的GPS站点处均大于观测数据，而较远站点处则明显小于观测值，且大部分站点具有比模型1更明显的衰减趋势。根据前人研究结果[9, 16]，玛多地震发震断层的摩擦属性在深度上可能存在差异，因此进一步构建在深部和浅部具有不同摩擦属性的应力驱动余滑模型(模型2)。
|
图 2 应力驱动余滑模型模拟的震后地表EW向形变与GPS观测数据对比 Fig. 2 Comparison of post-seismic surface deformation at EW direction simulated by stress-driven afterslip model and GPS observation data |
|
图 3 应力驱动余滑模型模拟的震后地表NS向形变与GPS观测数据对比 Fig. 3 Comparison of post-seismic surface deformation at NS direction simulated by stress-driven afterslip model and GPS observation data |
|
图 4 最佳拟合余滑分布 Fig. 4 Best fitted afterslip distribution |
基于余滑模型1的结果，通过大量正演实验测试余滑界面位置和摩擦属性分层深度的影响，最终将震后余滑对应的不同摩擦属性深度分界设定为13 km(图 5)。余滑模型2以该深度为界，对不同深度的余滑界面设置不同的参考滑动速率V0[17]，并参考运动学余滑结果[5-7]约束断层两端7 km深度以上的浅部余滑。
|
图 5 不同摩擦属性分界深度对应的最小残差 Fig. 5 The minimum residual corresponding to the depth of boundary of different friction properties |
通过对地表形变时间序列进行拟合，对0~13 km和13 km以下的
|
图 6 余滑模型2对应的拟合残差 Fig. 6 Fitting residuals corresponding to afterslip model 2 |
为测试不同同震模型对计算结果的影响，在其余参数不变的情况下，分别使用Wang等[18]和Zhao等[7]的同震滑动模型计算余滑模型2的震后地表LOS向形变，并与本文结果进行对比。由图 7可知，不同模型结果的主要特征比较接近，均沿走向分布在断层迹线两侧，但地表形变分布的集中程度不同，差异主要在断层附近。本文得到的模拟结果在形变分布和量级上介于Wang等[18]和Zhao等[7]的结果之间。
|
图 7 不同同震模型的LOS向形变模拟结果 Fig. 7 LOS direction deformation simulation results of different coseismic models |
为评估均匀弹性模型对搜索结果的影响，在最优拟合模型的基础上分别对断裂南北两侧GPS站点的观测数据进行拟合。结果显示，南侧的最优拟合参数
Xiong等[5]和He等[6]利用运动学余滑模型反演得到的滑动模型均能较好地解释玛多地震的震后短期大地测量观测数据。Xiong等[5]认为，余滑主要发生在同震破裂下侧，同震破裂区域基本没有发生滑移，最大滑动量约为0.3 m，约位于14 km深度处。He等[6]认为，余滑主要发生在同震破裂下倾区域，余滑峰值为0.3 m，约位于20 km深度处。运动学震后余滑特征的差异可能与采用的数据不同有关。Zhao等[7]采用InSAR形变约束的应力驱动余滑结果表明，余滑主要分布在发震断层的深部，浅部几乎没有发生余滑，滑动峰值为0.8 m，位于约10 km深度处，与本文余滑集中区域一致，量级接近余滑模型1的结果。
在断层西侧，He等[6]及Zhao等[7]给出的余滑量较小，而本文和Xiong等[5]的余滑量较大。这是由于本文使用的的同震滑动模型在断层西侧引起的应力变化较大。本文给出的余滑集中分布在5~15 km深度处，而Xiong等[5]及He等[6]给出的余滑分布较广，在远离同震破裂区域(深度>20 km)仍大量存在。这是由于应力驱动的余滑往往集中分布在同震滑动峰值区附近，远离同震显著滑动区无法产生余滑[19]。本文余滑模型2计算得到的震后余滑量最大值为1.24 m，远大于Xiong等[5]及He等[6]通过运动学余滑模型反演的结果；Zhao等[7]通过应力驱动余滑模型得到的震后余滑量最大值介于本文余滑模型1和余滑模型2之间。本文余滑模型2的结果中参考滑移速率
本文采用应力驱动余滑模型对2021年青海玛多MW7.4地震的震后余滑进行模拟，较好地解释了玛多地震震后160 d的GPS形变观测数据，获得了地震的震后余滑特征。玛多地震的震后余滑主要发生在同震断层显著滑动的下倾和两侧区域，集中分布在5~15 km深度(占总量80%)，最大滑移量达1.24 m，累积释放地震矩1.08×1019 Nm，相当于MW6.62地震。以13 km为界，浅部和深部在最优拟合下对应的V0分别为6.50 m/a和0.15 m/a，表明断层摩擦属性存在深度差异。
| [1] |
Bürgmann R, Dresen G. Rheology of the Lower Crust and Upper Mantle: Evidence from Rock Mechanics, Geodesy, and Field Observations[J]. Annual Review of Earth and Planetary Sciences, 2008, 36(1): 531-567 DOI:10.1146/annurev.earth.36.031207.124326
(  1) |
| [2] |
王未来, 房立华, 吴建平, 等. 2021年青海玛多MS7.4地震序列精定位研究[J]. 中国科学: 地球科学, 2021, 51(7): 1 193-1 202 (Wang Weilai, Fang Lihua, Wu Jianping, et al. Study on Precise Location of Maduo MS7.4 Earthquake Sequence in Qinghai in 2021[J]. Scientia Sinica Terrae, 2021, 51(7): 1 193-1 202)
( 1) |
| [3] |
程佳, 刘杰, 甘卫军, 等. 1997年以来巴颜喀拉块体周缘强震之间的粘弹性触发研究[J]. 地球物理学报, 2011, 54(8): 1 997-2 010 (Cheng Jia, Liu Jie, Gan Weijun, et al. Coulomb Stress Interaction among Strong Earthquakes around the Bayan Har Block since the Manyi Earthquake in 1997[J]. Chinese Journal of Geophysics, 2011, 54(8): 1 997-2 010)
( 2) |
| [4] |
Masuti S. Inferring Upper Mantle Rheology from Numerical Modeling of Experimental and Geodetic Observations[D]. Singapore: Nanyang Technological University, 2020
( 1) |
| [5] |
Xiong W, Chen W, Wang D Z, et al. Coseismic Slip and Early Afterslip of the 2021 MW7.4 Maduo, China Earthquake Constrained by GPS and InSAR Data[J]. Tectonophysics, 2022, 840
( 11) |
| [6] |
He L J, Feng G C, Wu X X, et al. Coseismic and Early Postseismic Slip Models of the 2021 MW7.4 Maduo Earthquake(Western China) Estimated by Space-Based Geodetic Data[J]. Geophysical Research Letters, 2021, 48(24)
( 8) |
| [7] |
Zhao D Z, Qu C Y, Chen H, et al. Tectonic and Geometric Control on Fault Kinematics of the 2021 MW7.3 Maduo(China) Earthquake Inferred from Interseismic, Coseismic, and Postseismic InSAR Observations[J]. Geophysical Research Letters, 2021, 48(18)
( 10) |
| [8] |
Marone C. The Effect of Loading Rate on Static Friction and the Rate of Fault Healing during the Earthquake Cycle[J]. Nature, 1998, 391(6 662): 69-72
( 1) |
| [9] |
Jin Z Y, Fialko Y. Coseismic and Early Postseismic Deformation Due to the 2021 M7.4 Maduo(China) Earthquake[J]. Geophysical Research Letters, 2021, 48(21)
( 3) |
| [10] |
Barbot S. Relax v1.0.7[Z]. Computational Infrastructure for Geodynamics, 2014
( 1) |
| [11] |
Wang K, Fialko Y. Observations and Modeling of Coseismic and Postseismic Deformation Due to the 2015 MW7.8 Gorkha(Nepal) Earthquake[J]. Journal of Geophysical Research: Solid Earth, 2018, 123(1): 761-779 DOI:10.1002/2017JB014620
( 2) |
| [12] |
Barbot S, Fialko Y, Bock Y. Postseismic Deformation Due to the MW6.0 2004 Parkfield Earthquake: Stress-Driven Creep on a Fault with Spatially Variable Rate-and-State Friction Parameters[J]. Journal of Geophysical Research, 2009, 114(B7)
( 3) |
| [13] |
Jin Z Y, Fialko Y, Zubovich A, et al. Lithospheric Deformation Due to the 2015 M7.2 Sarez(Pamir) Earthquake Constrained by 5 Years of Space Geodetic Observations[J]. Journal of Geophysical Research: Solid Earth, 2022, 127(4)
( 1) |
| [14] |
Tian Z, Freymueller J T, Yang Z Q. Postseismic Deformation Due to the 2012 MW7.5 Craig Earthquakes and Its Implications for Regional Rheological Structure[J]. Journal of Geophysical Research: Solid Earth, 2021, 126(2)
( 1) |
| [15] |
Huang L Y, Zhang B, Shi Y L. Stress Transfer at the Northeastern End of the Bayan Har Block and Its Implications for Seismic Hazards: Insights from Numerical Simulations[J]. Earth and Space Science, 2021, 8(12)
( 1) |
| [16] |
Fang J, Ou Q, Wright T J, et al. Earthquake Cycle Deformation Associated with the 2021 MW7.4 Maduo(Eastern Tibet) Earthquake: An Intrablock Rupture Event on a Slow-Slipping Fault from Sentinel-1 InSAR and Teleseismic Data[J]. Journal of Geophysical Research: Solid Earth, 2022, 127(11)
( 1) |
| [17] |
Wang K, Bürgmann R. Probing Fault Frictional Properties during Afterslip Updip and Downdip of the 2017 MW7.3 Sarpol-e Zahab Earthquake with Space Geodesy[J]. Journal of Geophysical Research: Solid Earth, 2020, 125(11)
( 2) |
| [18] |
Wang M, Wang F, Jiang X, et al. GPS Determined Coseismic Slip of the 2021 MW7.4 Maduo, China, Earthquake and Its Tectonic Implication[J]. Geophysical Journal International, 2021, 228(3): 2 048-2 055 DOI:10.1093/gji/ggab460
( 2) |
| [19] |
Huang M H, Bürgmann R, Freed A M. Probing the Lithospheric Rheology across the Eastern Margin of the Tibetan Plateau[J]. Earth and Planetary Science Letters, 2014, 396: 88-96 DOI:10.1016/j.epsl.2014.04.003
( 1) |
2. Key Laboratory of Geospace Environment and Geodesy, Ministry of Education, Wuhan University, 129 Luoyu Road, Wuhan 430079, China;
3. Key Laboratory of Geophysical Geodesy, MNR, 129 Luoyu Road, Wuhan 430079, China;
4. Institute of Seismology, CEA, 40 Hongshance Road, Wuhan 430071, China