地球物理学进展  2017, Vol. 32 Issue (6): 2394-2400   PDF    
引用本文
杨九元, 张永志, 祝意青, 等. 2017. 尼泊尔MW 7.9地震同震及震后地表形变及重力变化模拟. 地球物理学进展, 2017, 32(6): 2394-2400, doi: 10.6038/pg20170614  
YANG Jiu-yuan, ZHANG Yong-zhi, ZHU Yi-qing, et al. 2017. Simulation of co-seismic and post-seismic surface deformation and gravity changes of Nepal MW 7.9 earthquake. Progress in Geophysics (in Chinese), 32(6): 2394-2400, doi: 10.6038/pg20170614   
尼泊尔MW 7.9地震同震及震后地表形变及重力变化模拟
杨九元1, 张永志1, 祝意青2, 姜永涛3, 杨珍1, 徐云马2     
1. 长安大学地质工程与测绘学院, 西安 710054
2. 中国地震局第二监测中心, 西安 710054
3. 南阳师范学校环境科学与旅游学院, 河南南阳 473061
收稿日期：2017-04-22; 修回日期：2017-10-12.
基金项目：国家自然科学基金项目（41374028，61627824，41274083，41374026）和南阳师范学院博士专项基金（ZX2016007）联合资助.
作者简介：杨九元, 男, 1990年生, 硕士生, 主要研究方向为地壳形变与地球动力学.(E-mail:291813306@qq.com)
摘要：基于USGS公布滑动分布模型，本文利用用地壳分层模型，考虑到自重及黏弹特性，采用数值模拟的方法，对尼泊尔MW 7.9大地震同震及黏弹松弛效应引起的震后形变场及重力变化进行模拟计算，并与发震断层地表投影面附近区域实测同震GPS水平形变数据对比，结果显示：（1）同震形变场及重力变化均显示该次地震主要表现为逆冲型，且主要形变及重力变化主要发生在发震断层地表投影区域，离断层越远区域，形变量及重力变化值越小；（2）震后形变量及重力变化均增大，且影响范围逐渐增大，其中，垂直形变变化趋势与重力变化相反，表明地表高程变化与重力变化具密切的联系；（3）模拟同震结果与实测同震水平形变结果比较符合，少数台站点差异较大，但运动趋势与实测结果基本一致，一定程度上验证了模拟同震及震后形变量及重力变化的可靠性.
关键词：USGS    地壳分层模型    尼泊尔地震    黏弹松弛效应    GPS数据    
Simulation of co-seismic and post-seismic surface deformation and gravity changes of Nepal MW 7.9 earthquake
YANG Jiu-yuan1 , ZHANG Yong-zhi1 , ZHU Yi-qing2 , JIANG Yong-tao3 , YANG Zhen1 , XU Yun-ma2     
1. School of Geology Engineering and Geomatics, Chang'an University, Xi'an 710054, China
2. Second Crust Monitoring and Application Center, China Earthquake Administration, Xi'an 710054, China
3. School of Environment Science and Tourism, Nanyang Normal University, Henan Nanyang 473061, China
Abstract: Based on the sliding distribution model provided by the USGS, this paperuses the earth's crust layered model and the method of numerical simulation, considering the gravity and sticky elastic characteristics, to calculatethe co-seismic deformation field and the gravity change and the post-seismic deformation field and the gravity change, then, this paper compares the co-seismic deformation with theco-seismic GPS date, Result show:(1) Co-seismic deformation field and gravity change shows that the earthquake belongs to thrust type, and the main deformation and gravity change occurs mainly in the seismogenic fault surface projection region, the farther from the fault zone, the smaller form variables and gravity change value; (2) The post-seismic deformation and gravity change all increase, and the influence range is increased gradually; Among them, the vertical deformation trend is the opposite of gravity change trend, it show that the surface height and gravity changes have close connection; (3)The simulation co-seismic results are in good agreement with the measured co-seismic deformation results, a few sites are different, but the trend of movement and the measured results are basically identical, it partly verifies the reliability of the co-seismic and post-seismic deformation field and the gravity change by the simulation calculation.
Key words: USGS     earth's crust layered model     Nepal earthquake     sticky relaxation effect     GPS data    
0 引言

尼泊尔地处欧亚块体与印度块体的主要撞裂带上，紧邻喜马拉雅山，地形比较险峻，沿该撞裂带历史上发生了多次强震，如：MW 8.7察隅地震、MW 7.6 Kashimir地震及MW 6.9 sikkim地震等，而2015年4月25日尼泊尔发生近81年来最强地震(MW 7.9)，震中位于加德满都西北约80 km处，地理坐标为(28.15°N, 84.71°E), 震源深度为15 km(USGS公布)，属于危害较大的浅源地震.由于该次地震破坏性较大，引起了国内许多研究者关注，张勇等(2015)基于远场地震波数据在震后2 h 9 min快速反演出该次地震的断层破裂模型；王卫民等(2015)利用实时地震波数据与GPS同震数据联合反演出该地震震源机制及破裂模型，为研究地质构造、预估地震危害提供依据；谭凯等(2016)等考虑到地壳垂直分层及横向差异性，以GPS和INSAR同震形变数据为约束，反演尼泊尔MW 7.9主震及MW 7.3余震滑动破裂分布模型，为研究青藏高原动力学特性奠定基础；在已有断层滑动分布模型研究的基础上，熊维等(2015)考虑到岩石圈分层及震后黏弹松弛效应，对尼泊尔地震同震及震后应力场效应进行了深入的研究；张贝等(2015)采用PREM地球分层模型，以美国地质调查局(USGS)和中国地震局(CEA)公布的断层滑动模型为基础，计算了尼泊尔地震同震形变及应力场变化，并探讨该次地震震后较短时间内地震活动特性，在研究中考虑到地球曲率的影响，但其未考虑黏弹特性在震后引起形变效应；为了进一步了解尼泊尔MW 7.9大地震破裂过程，基于USGS已公布不同深度的断层滑动分布，本文利用地壳分层模型，考虑到自重及区域流变系数的影响，采用数值模拟的方法，对该次地震同震及黏弹松弛效应引起的震后长时间内形变及重力变化效应进行计算模拟，并与收集到的实测同震GPS水平形变数据进行对比，详细探讨尼泊尔地震同震及震后形变及重力变化趋势，为以后预防该区域孕震及震后形变与重力测量提供一定的参考依据.

1 计算方法及模型

地震的发生往往会引发物理场的变化，其中包括同震引起的位移场、重力场、应力变化等弹性响应及震后黏弹性松弛效应引起的形变场及应力变化等；而位错作为描述地震断层模型与地表形变场及重力变化之间关系的重要理论，自应用于地震学计算形变场以来，已从简单的均匀弹性半空间模型逐渐发展到考虑分层、不均匀性、横向差异性、黏弹性、自重等影响因素使位错理论模型趋于更加完善.继Wang等(2003)推算出基于弹性半空间分层模型计算同震及震后形变场及重力变化的软件EDGRN/EDCMP后，又开发出基于黏弹性半空间分层模型计算同震及震后形变及重力变化的软件PSGRN/PSCMP(Wang et al., 2006), 且考虑到自重的影响；其中，PSGRN程序主要用来计算描述黏弹分层模型中不同深度处四个基本位错源的时变格林函数，而PSCMP程序基于PSGRN的计算结果，将断层破裂面离散成点源位错，通过线性叠加的方法计算同震及震后不同时间段内的形变场、重力变化及大地水准面变化等.

为了进一步了解孕震、地震触发、断层破裂过程及滑动分布等问题，本文基于Wang等(2006)开发的地震形变场计算程序(PSGRN/PSCMP)对尼泊尔地震同震及震后形变效应进行模拟计算，且计算中考虑介质自重影响及震后黏弹性松弛效应，分层半空间中的黏滞部分一般采用MAXWELL流体来近视表示.其中，地壳分层结构模型参考已有研究成果(熊维等，2015)，如表 1.

表 1 地壳分层模型 Table 1 The earth's crust layered model

地震发生后美国地质调查局(USGS)快速公布的断层滑动分布模型主要根据地震台站的远场地震波数据反演确定，与前人(刘刚等，2015；王卫民等，2015；谭凯等，2016)研究结果基本一致，均显示出该次地震破裂主要以逆冲滑动为主.本文研究采用USGS已公布的断层滑动模型，该模型断层面长220 km，宽165 km，分别沿走向和倾向方向以步长20 km、15 km均匀离散为11×11个有限断层面，其中，走向角为295°，倾角为10°，最大滑动量为3.1 m.如图 1所示.

图 1 断层滑动分布模型(USGS提供) Figure 1 Fault slip distribution model(provided by USGS)
2 尼泊尔地震同震及震后形变及重力变化模拟 2.1 同震形变及重力变化模拟

同震变形是地壳介质对地震破裂的弹性响应，为更好分析地震同震引起区域形变及重力效应，本文主要对26N°~30N°，82E°~88E°区域形变场进行模拟计算，采用上文中地壳分层结构及断层滑动模型，运用PSGRN/PSCMP软件包进行同震计算，结果如图 2.

图 2 尼泊尔同震水平经向位移(a)、纬向位移(b)、垂直位移(c)及重力变化(d) Figure 2 The coseismic level meridional displacement(a)、zonal displacement(b)、vertical displacement (c) and gravity change of Nepal(d)

图 2a为尼泊尔震中附近区域同震水平经向位移场变化图，以震中破裂带为界，断层面附近西北区域表现为负值，最大负值达到-357.1 mm，东南区域表现为正值，最大正值达到362.5 mm，整体上看，等值线具有一定的对称特性，且附近位移量变化较大值主要位于震中断层投影面区域，断层面近场区域位移量变化量大，断层面远场区域，位移量变化较慢，影响范围较大.

图 2b为震中附近区域同震纬向位移场变化图，整体上，纬向位移变化值主要表现为负值，最大负值达到-1255 mm，断层面震中附近纬向位移变化迅速，离断层面越远，位错变化速度越慢，整个模拟区域纬向位移变化范围为-1255~20.62 mm.

图 2c为震中附近同震垂直位移图，可知，断层面投影区域以外垂直位移变化表现为较明显的对称现象，断层面以南地区主要表现为地壳上升，最大上升值为520.4 mm，断层面以北喜马拉雅山区域主要表现为地壳沉降，最大沉降值为-204.4 mm.

图 2d为震中附近同震重力场变化图，对比2b可知，重力场变化图与垂直位移场变化正负值刚好相反，且除断层面区域之外也表现出一定的对称现象，证明重力场变化与区域地壳台上或沉降有着一定的联系.模拟区域重力变化范围为-111.5~36.38，断层面以北喜马拉雅山区域重力变化值为正，最大值为36.38，断层面以南区域重力变化主要表现为负，最大值为-111.5.

图 3为尼泊尔地震震中附近区域地形图及模拟同震水平位移场矢量图，可知，发震断层上盘向上运动，下盘向下运动，说明此次为逆冲型地震.其中，水平形变值较大主要集中在27.6°N~29.2°N，84.2°E~86.4°E区域及周边.

图 3 尼泊尔同震水平位移矢量图(紫色方框为断层地表投影面) Figure 3 The vector map of Nepal coseismic horizontal displacement (Purple box in picture for illustrating the fault surface plane)
2.2 黏弹性松弛效应引起震后地表形变及重力场模拟

震后形变过程主要由震后余滑、上/下地壳黏弹松弛、空系流体等多种复杂构造机制组成，其中，震后余滑和空系流体作用在震后较短时间内快速消失，而黏弹松弛效应在震后数十年对地壳形变场仍有很大影响.选取黏滞系数的大小不同虽然对震后水平形变及垂直形变场变化趋势影响不大，但对于震后模拟形变场值大小影响较大(邓明莉等，2014)，因此，黏滞系数的选取尤为重要.根据徐晶等(2013)研究可知，青藏高原中下地壳黏滞系数在1019~1020 Pa·s之间，而青藏高原上地幔黏滞系数可以等效为1020 Pa·s(万永革等，2008)，本文参考前人对地壳速度结构模型的研究成果，设定地壳分层模型上地壳(0~20 km)为弹性体，中下地壳(20~51 km)及上地幔(51 km以下)用MAXWELL流体来近似表示，其中，中下地壳黏滞系数取1.0×1019 Pa·s，上地幔黏滞系数取1.0×1020 Pa·s，运用PSGRN/PSCMP软件对尼泊尔地震震中附近区域震后形变及重力场变化进行模拟，进一步了解震后形变变化趋势，对该区域今后震后防震防灾工作具有一定的参考价值.

图 4a和b分别为模拟区域震后10年、70年水平形变值，与图 3同震水平形变相比，在方向上整体趋势不变，断层投影面区域形变值随着黏弹松弛时间的增加，形变值逐渐减小，而断层面远场区域西北处及东南处形变值逐年增大.其中，震后10年水平经向形变及纬向形变值范围分别为-313.1~296.5 mm、-1106~46.98 mm，震后70年水平经向形变及纬向形变值范围分别为-197~221 mm、-741~179 mm，相对于同震形变值变化范围均较小.

图 4 (a)-(b)、(c)-(d)、(e)-(f)分别为尼泊尔震后水平形变，垂直形变、重力变化图(包括同震形变效应) 图(a)、(c)、(e)为震后10年变化，图(b)、(d)、(f)为震后70年变化. Figure 4 (a)-(b)、(c)-(d)、(e)-(f)are respectively horizontal deformation, vertical deformation, gravity variation of Nepal earthquake Fig(a)、(c)、(e)are postseismic10years changes, Fig(b)、(d)、(f)are postseismic 70years changes.

图 4c和d分别为模拟区域震后10年、70年垂直形变，与图 2c同震垂直形变相比可知，震后10年垂直形变整体分布趋势基本没变，但断层投影面附近以北区域垂直形变值沉降值减下，震后70年垂直形变分布趋势有较大改变，断层投影面以北附近地区沉降范围逐渐减小，并逐渐表现为抬升，且地壳抬升范围逐渐增大，而断层面以南由垂直抬升逐渐表现为沉降，地壳抬升区域范围变大.

图 4e和f分别为模拟区域震后10年、70年重力场变化，与震后垂直形变场对比可知，震后重力变化区域与垂直分布趋势相反，且同震重力变化与同震垂直变化分布趋势也相反，因此，地壳高程的变化对重力变化有较大影响，具有一定的关联性.

为研究尼泊尔地震长期时变效应，本文在断层投影面两侧任取4个点，通过计算发现该次地震在震后300年后地壳形变和重力变化基本达到稳定，因此，为了深入了解该次尼泊尔大震震后地壳总形变量及总重力变化，对震后400年(去除同震效应)形变场及重力变化进行了模拟计算，对震后危害性评估具有一定的意义，结果如图 5所示.

图 5 (a)-(c)分别为震后400年水平形变、垂直形变及重力场变化 (注去除同震效应) Figure 5 (a)-(c) are respectivelyhorizontal deformation、vertical deformation and gravity changesin 400years after the earthquake (note: Remove the coseismic effects)

黏弹松弛效应引起的震后400年(去除同震效应)断层投影面区域水平形变形变值(图 5a)变化最大，主要表现向北扩展趋势，与同震引起断层地表投影面区域的变化趋势相反，而断层面以北水平形变值大小逐年增大，整体趋势变化与同震引起的形变趋势基本一致，但影响范围增大，断层面以南亦如此；黏弹松弛效应引起的震后400年垂直形变变化(图 5b)在模拟区大部分范围表现沉降，只有在断层面以北区域表现为抬升，变化幅度为-369.09~622.9 mm，震后400年重力变化趋势(图 5c)与垂直变化趋势仍然相反，变化幅度为-124.85~78.007 μGal.

2.3 讨论

限于资料，本文仅收集到了尼泊尔MW 7.9地震同震水平形变数据(李瑜等，2015)，在一定程度上可以验证模拟区域同震及震后形变效应的可靠性，由于PSGRN/PSCMP软件包是基于平面位错的基础之上考虑到自重及黏弹性分层，而未考虑到地球曲率的影响，对远场台站的计算误差相对较大，因此，本文主要选取离震中断层面地表投影区域较近的GPS连续观测站台进行模拟计算及分析，结果表 2所示.

表 2 GPS实测与模拟结果 Table 2 GPS measurement and simulation results

由表中模拟结果可知，除CHLM、KKN4、NAST及SNDL站台与实测值差异较大外，其他GPS站台实测值与模拟计算值虽有一定的差异，但基本在误差范围内；CHLM、KKN4及NAST站台模拟计算结果相对于实测值均较小，而SNDL站台模拟计算值相对于实测值较大，但该4个站台点模拟结果运动趋势均与实测值保持一致.总的来说，模拟计算值相较于实测值得到了比较满意的结果.

3 结论

本文基于尼泊尔地震滑动分布模型(USGS提供)，考虑到区域流变系数及自重的影响，采用数值模拟的方法对尼泊尔地震同震及震后形变及重力变化效应进行模拟计算，并与实测同震数据进行对比分析，结果表明：

(1) 尼泊尔地震为逆冲型地震，该次地震同震形变量及重力变化主要表现在发震断层地壳投影面区域，离断层越远区域，形变量及重力变化值越小；其中，水平经向位移、纬向位移幅度范围分别达到-357.1~362.5 mm、-1255~20.62 mm，垂直形变变化在断层面以南地区主要表现为地壳上升，断层面以北喜马拉雅山区域主要表现为地壳沉降，幅度范围为-204.4~520.4 mm，重力变化与垂直形变变化趋势相反.

(2) 震后10年、70年水平形变量在断层地表投影区域逐渐减小，而在断层远场区域逐渐增大，通过震后400年(去除同震效应)总的形变及重力变化可知，这种情况是由于断层面投影区域震后总的水平形变变化趋势与同震水平形变变化趋势相反，而断层远场震后总的水平形变变化趋势与同震水平形变变化趋势一致导致的；垂直形变变化在断层投影面以北附近地区沉降范围逐渐减小，并逐渐表现为抬升，且地壳抬升范围逐渐增大，而断层面以南由垂直抬升逐渐表现为沉降，地壳抬升区域范围变大；重力变化与垂直变化趋势相反，表明地壳垂直变化与重力变化具有密切的联系.

(3) 模拟结果与实测同震GPS数据对比可知，大部分台站点与实测数据比较符合，均在一定的误差范围内，CHLM、KKN4及NAST站台模拟计算结果相对于实测值均较小，而SNDL站台模拟计算值相对较大，但4个站台点模拟结果运动趋势均与实测值保持一致，在一定程度上验证了模拟同震及震后形变及重力变化模拟的可靠性.

致谢 感谢审稿专家提出的修改意见和编辑部的大力支持！
参考文献
[] Deng M L, Sun H P, Xu J Q, et al. 2014. Theoretical simulation of co-seismic and post-seismic deformations and gravity changes of Lushan earthquake[J]. Earth Science——Journal of China University of Geosciences, 39(9): 1373–1382.
[] Li Y, Liu J, Liang H, et al. 2015. Co-seismic displacement field associated with the 25 April, 2015 MW 7.8 Nepal earthquake recorded by Global Positioning System[J]. Chin. Sci. Bull., 60(36): 3606–3616. DOI:10.1360/N972015-00602
[] Liu G, Wang Q, Qiao X J, et al. 2015. The 25April 2015 Nepal MS 8.1 earthquake slip distribution from joint inversion of teleseismic, static and high-rate GPS data[J]. Chinese Journal of Geophysics, 58(11): 4287–4297. DOI:10.6038/cjg20151133
[] Tan K, Zhao B, Zhang C H, et al. 2016. Rupture models of the Nepal MW 7.9 earthquake and MW 7.3 aftershock constrained by GPS and InSAR coseismic deformations[J]. Chinese Journal of Geophysics, 59(6): 2080–2093. DOI:10.6038/cjg20160614
[] Wan Y G, Shen Z K, Sheng SZ, et al. 2009. The influence of 2008 Wenchuan earthquake on surrounding faults[J]. Acta Seismologica Sinica, 31(2): 128–139.
[] Wang R., Lorenzo Martı'n F., Roth F.. 2003. Computation of deformation induced by earthquakes in a multi-layered elastic crust-FORTRAN programs EDGRN/EDCMP[J]. Computers & Geosciences, 29(2): 195–207.
[] Wang R., Roth F., Lorenzo-Martin 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[J]. Computers and Geosciences, 32(4): 527–541. DOI:10.1016/j.cageo.2005.08.006
[] Wang W M, Hao J L, He J K, et al. 2015. Rupture process of the MW 7.9 Nepal earthquake April 25, 2015[J]. Science China:Earth Sciences, 58(10): 1895–1900. DOI:10.1007/s11430-015-5170-y
[] Xiong W, Tan K, Liu G, et al. 2015. Coseismic and postseismic Coulomb stress changes on surrounding major faults caused by the 2015 Nepal MW 7.9 earthquake[J]. Chinese Journal of Geophysics, 58(11): 4305–4316. DOI:10.6038/cjg20151135
[] Xu J, Shao Z G, Ma H S, et al. 2013. Evolution of Coulomb stress and stress interaction among strong earthquakes along the Xianshuihe fault zone[J]. Chinese Journal of Geophysics(in chinese), 56(4): 1146–1158. DOI:10.6038/cjg20130410
[] Zhang B, Cheng H H, Shi Y L. 2015. Calculation of the co-seismic effect MS 8.1earthquake, Apirl 25, 2015, Nepal[J]. Chinese Journal of Geophysics, 58(5): 1794–1803. DOI:10.6038/cjg20150529
[] Zhang Y, Xu L S, Chen Y T. 2015. Rupture process of the 2015 Nepal MW 7.9 earthquake:Fast inversion and preliminary joint inversion[J]. Chinese Journal of Geophysics (inChinese), 58(5): 1804–1811. DOI:10.6038/cjg20150530
[] 邓明莉, 孙和平, 徐建桥, 等. 2014. 芦山地震同震和震后地表形变及重力变化的理论模拟[J]. 地球科学——中国地质大学学报, 39(9): 1373–1382.
[] 李瑜, 刘静, 梁宏, 等. 2015. 全球定位系统测定的尼泊尔MW 7.8级地震同震位移[J]. 科学通报, 60(36): 3606–3616. DOI:10.1360/N972015-00602
[] 刘刚, 王琪, 乔学军, 等. 2015. 用连续GPS与远震体波联合反演2015年尼泊尔中部MS 8.1地震破裂过程[J]. 地球物理学报, 58(11): 4287–4297. DOI:10.6038/cjg20151133
[] 谭凯, 赵斌, 张彩红, 等. 2016. GPS和InSAR同震形变约束的尼泊尔MW 7.9和MW 7.3地震破裂滑动分布[J]. 地球物理学报, 59(6): 2080–2093. DOI:10.6038/cjg20160614
[] 万永革, 沈正康, 盛书中, 等. 2009. 2008年汶川大地震对周围断层的影响[J]. 地震学报, 31(2): 128–139.
[] 王卫民, 郝金来, 何建坤, 等. 2015. 2015年4月25日尼泊尔MW 7.9级地震震源过程[J]. 中国科学:地球科学, 45(9): 1421–1426.
[] 熊维, 谭凯, 刘刚, 等. 2015. 2015年尼泊尔MW 7.9地震对青藏高原活动断裂同震、震后应力影响[J]. 地球物理学报, 58(11): 4305–4316. DOI:10.6038/cjg20151135
[] 徐晶, 邵志刚, 马宏生, 等. 2013. 鲜水河断裂带库仑应力演化与强震间关系[J]. 地球物理学报, 56(4): 1146–1158. DOI:10.6038/cjg20130410
[] 张贝, 程惠红, 石耀霖. 2015. 2015年4月25日尼泊尔MS 8.1大地震的同震效应[J]. 地球物理学报, 58(5): 1794–1803. DOI:10.6038/cjg20150529
[] 张勇, 许力生, 陈运泰. 2015. 2015年尼泊尔MW 7.9地震破裂过程:快速反演与初步联合反演[J]. 地球物理学报, 58(5): 1804–1811. DOI:10.6038/cjg20150530