NUMERICAL MODELING FOR MECHANISM OF YUSHU MS7.1 EARTHQUAKE
-
摘要: 根据玉树地区的地应力场、速度场和断层展布, 对青海玉树2010年4月14日MS7.1级地震发震机理进行了数值模拟。将围岩看成弹性体, 断层看成具有应变软化的弹塑性体, 断层和围岩组成统一的地质介质系统。在地应力、孔隙压力及边界位移的作用下, 应力逐渐积累, 当达到断层摩擦破坏强度时, 断层产生应变软化, 断层突然滑动, 能量突然释放, 应力突然下降, 形成地震。模拟结果表明:玉树7.1级地震是在印度板块向北推挤, 青藏高原向东南侧向挤压, 在玉树地区形成主压应力为北东80°方向的水平应力场, 使甘孜-玉树断裂带产生左旋走滑错动形成的。计算结果给出了应力降、能量释放量、断层走滑错动量、地震复发周期、应力积累速度等重要参数, 模拟结果与野外调查资料具有较好的一致性。Abstract: In this paper, strain softening model is used to establish the modeling for Yushu MS7.1 earthquake. The results shows that Yushu MS7.1 earthquake occurred as a combination of northward compressing of the Indian Plate and southeastward compressing of the Qinghai-Tibet Plateau. The stress drop, energy release amount, the magnitude of dislocation of fault, and earthquake recurrence intervals are given by the modeling. The modeling results show a good match with the geologic survey in the field.
-
0. 前言
本世纪以来, 继2001年11月14日昆仑山口西8.1级强震和2008年5月12日汶川8级[1~8]强震之后, 于2010年4月14日7时49分40.7秒在我国青海省玉树县(33.23°N, 96.6°E)发生MS 7.1级地震[9~10], 震源深度14km, 该地震发生在甘孜—玉树断裂带上, 震源机制得出的断层面走向为119°/倾角8 3°/滑动角-2°。这次地震破裂过程持续了约23s, 包括两次主要的子事件, 第一次是在初始破裂后0~5.5s, 第二次在5.5~23s左右[11]。地表破裂主要集中在距震中东南约10~30km处, 破裂长度延续23km, 破裂错动以反扭走滑为主, 错动量可达1.75m[12]。地震的破裂主要由北西向东南方向扩展, 破裂滑移量最大的区域位于玉树城区附近。
玉树地震位于巴颜喀拉地块南边界, 汶川地震属于巴颜喀拉地块东南边界。他们属于同一地块不同的边界上, 具有相同的构造动力学背景。
地震造成严重破坏及人员伤亡和财产损失。震后第一时间中国地质科学院地质力学研究所专家赶赴现场进行考查, 为灾后重建, 减灾救灾和地震科学研究提供了重要资料。在此基础上, 我们对玉树地震发生的机理进行了数值分析, 为地震机理的研究和预测预报提供依据。
我们将围岩看成弹性介质, 断层看成具有软化特性的弹塑性介质, 采用非稳定断层滑动准则, 利用FLAC3D软件进行地震发生机理的数值模拟。
在地应力场的作用下, 应力逐渐积累, 当应力达到断层的摩擦破坏强度时, 断层发生软化, 产生突然滑动, 储存在围岩中的能量突然释放, 形成地震[13~15]。
1. 玉树地区的应力场和速度场
玉树地震区及其附近地区地应力场最大水平主应力方向为北东到北东东向如图 1[16]。玉树MS 7.1级地震的震源机制解得出的p轴方位为74°[17]。
震源机制解得出的p轴方位与主应力方位一般情况下是不同的。他们之间的区别与岩石或断层的内摩擦有关。内摩擦角越大, 差别越大。p轴方位是在不考虑内摩擦的情况下确定的, 它与断层面的夹角为45°。但是, 当存在内摩擦时, 主应力方向与断层面的夹角小于45°。考虑岩石和断层内摩擦角的存在, 将最大水平主应力方向调整为北东80°, 与断层走向的夹角为40°。
据本区域地应力测量结果统计最大水平主应力σH与最小水平主应力σh之比在1.27~1.98之间, 多数比值为1.45左右[18]。地应力随深度可视为线性增加。并可用下式表示:
(1) (2) 式中z为深度(m)。
GPS测量得到的速度场如图 2[19].可以看出, 玉树附近的速度方向为北东70°到北东90°。西部速度快, 东部速度慢。由东向西有速度梯度。正因为存在速度梯度, 才会产生应力和应变。如果没有速度梯度, 速度处处相同, 则只有物体的整体移动, 而无变形, 因而无应变, 也无应力。玉树地区的应力场和速度场的特点, 正好有利于甘孜—玉树断层产生左旋走向滑动, 形成走滑地震。
2. 计算模型、边界条件和模型网格
模型的位置取在玉树附近。断层产状为直立断层, 走向120°。模型长度80km, 方向北东80°。模型宽度15km, 深度27km如图 3。模型中的断层长度23km, 与地表破裂长度相当。模型受到地应力和孔隙压力的作用。最大水平主应力σH的方向与模型长轴平行, 与断层走向夹角为40°。最小水平主应力的方向与长轴垂直。因为断层是直立的, 故不考虑重力作用。最大水平主应力和最小水平主应力皆随深度变化, 其表达式如式(1)及(2)。孔隙压力为静水压力, 其表达式为: p=0.01z(MPa), z为深度(m)。
两个水平主应力和孔隙压力均作为初始应力和初始孔隙压力施加在模型上。同时, 在模型的左边界施加速度, 右边界横向约束。模型的底边界z方向约束。
图 4是模型网格。对于围岩使用六面体网格, 共划分1500个单元。断层使用不连续界面单元。
3. 本构方程
3.1 断层的本构方程
(3) (4) (5) 式中: kt, kn为断层的切向和法向刚度; τzx, τzy, τ, dτzx, dτzy, dσz为沿断层面的剪应力和法向应力分量及增量;
; ux, uy, uz, 为断层面在x, y, z方向的位移间断量; μ为摩擦系数; c′和μ′为内聚力和摩擦系数随塑性内变量的变化, 表征介质的强化或软化特性; b为断层厚度。当考虑孔隙压力时, 上式中的应力皆为有效应力。3.2 围岩的本构方程
围岩为弹性体, 他的本构方程为:
(6) 式中: De为弹性矩阵, ε为应变。
4. 断层非稳定滑动的能量准则
断层和围岩组成一个地质介质系统, 在断层失稳滑动前的时刻, 断层处在非稳定平衡状态, 一个小的扰动, 如位移、应力、孔隙压力等的扰动, 系统的弹性应变能将发生改变, 变化的能量, 除了驱动断层滑动外, 还能产生地震波的动能。下式表示断层非稳定滑动的能量准则[22, 24]:
(7) 当上式满足时, 断层发生非稳定滑动。式中ΔU是系统弹性应变能的变化, ΔE是系统的耗散能; 即断层非稳定滑动时, 系统弹性应变能的改变大于断层破坏滑动过程所需的能量, 其中的差值用于产生地震波。
5. 介质参数的选择
我们参考文献[22], [25], [26]模拟地震时的摩擦系数和内聚力以及他们的软化情况进行了参数选择如表 1。
表 1 断层和围岩参数Table 1. Parameters of fault and surroundings6. 地震应力降和能量释放
在水平地应力σH、σh和孔隙应力p的作用下, 随着远场位移(边界位移)的增加, 应力逐渐加大。当平均应力达到396.35MPa时, 断层突然错动(图 5), 能量突然释放(图 6), 应力突然下降(图 7), 此时发生地震。地震造成的应力降为1.1MPa, 释放的能量密度为9×104尔格/cm3。应力降和能量释放的大小与断层位错量有关, 断层错动量越大, 能量释放越大, 应力降也大。
图 5 地震时断层的突然走向滑动(监测点为图 8 a点)Figure 5. Sudden sliding of faultwhen the earthquake happened平均应力是指模型深度范围内的水平应力平均值, 由下式得出:
(8) 式中: H0为深度, H0=27km。
7. 断层的错动
图 8为地震时, 断层在地表的相对错动。可以看出, 断层为逆时针走向滑动。断层两盘的相对滑动量最大为1.75m。与现场调查结果一致。地表破裂长度约23km。
图 9为断层位移随离开断层距离的变化。可以看出, 断层线附近断层相对位移最大, 随着离开断层线距离的加大, 位移量变小。在远处出现少量反向位移。
8. 地震复发周期
在当前的应力状态和地层参数不变的情况下, 应力增加1.1MPa后, 地应力又达到临界状态, 可再发生7.1级地震。地震的复发周期与地应力的积累速度有关。应力的积累速度越快, 复发周期越短。地应力的积累速度可由图 7求出。即由图 7的应力-位移曲线求出远场每单位位移引起的应力变化, 而单位位移所需要的时间可由速度求出, 进而求出地应力的积累速度。这里所指的速度是模型东西两侧的相对速度, 即两侧的速度差。由图 2可求得模型两侧的速度差为1.36~3.36mm/a, 可求得应力积累速度为1.23~3×10-3MPa/a。同样, 可求得应变速率, 其结果为1.7~4.2×10-8/a。可见, 应力和应变积累速度是非常缓慢的。
在这样的应力积累速度下, 应力增加1.1MPa的时间为360~900年。即玉树7.1级地震复发周期为360~900年。
9. 结论
1.玉树地震是印度板块持续向北挤压, 青藏高原物质向东南方向流动, 在玉树地区形成主压应力为北东80°方向的水平应力场和速度场。在这样的水平地应力和速度场作用下, 应力逐渐积累, 再加上孔隙流体压力的作用, 使甘孜—玉树断裂的应力达到其破裂强度, 断层突然滑动, 发生地震。
2.随着水平地应力的缓慢增加, 当达到396.35MPa时(平均值), 断层失稳滑动, 形成地震。地震的应力降为1.1MPa, 能量释放密度为9×104尔格/cm3。断层产生反扭走向滑动, 相对滑动量可达1.75m。
3.玉树地震区应力的积累速度为1.23~3×10-3 MPa/a。应变积累速度为1.7~4.2×10-8/a。应力和应变积累速度是缓慢的。
4.在当前玉树地震区速度梯度的情况下。地震复发周期为360~900年。当速度梯度改变或断层参数改变时, 地震复发周期将改变。
-
图 5 地震时断层的突然走向滑动(监测点为图 8 a点)
Figure 5. Sudden sliding of faultwhen the earthquake happened
表 1 断层和围岩参数
Table 1. Parameters of fault and surroundings
-
[1] 谭成轩, 孙叶, 吴树仁, 等."5.12"汶川MS 8.0大地震后关于我国区域地壳稳定性评价的思考[J].地质力学学报, 2009, 15(2):142~150. doi: 10.3969/j.issn.1006-6616.2009.02.004TAN Cheng-xuan, SUN Ye, WU Shu-ren, et al.A consideration on regional crustal stability assessment after MS 8.0 Wenchuan strong earthquake in China[J].Journal of Geomechanics, 2009, 15(2):142~150. doi: 10.3969/j.issn.1006-6616.2009.02.004 [2] 姚鑫, 张永双.基于差分干涉雷达的汶川地震同震形变特点[J].地质力学学报, 2009, 15(2):151~161. doi: 10.3969/j.issn.1006-6616.2009.02.005YAO Xin, ZHANG Yong-shuang.Co-seismic deformation of 5.12 Wenchuan earthquake based on D-INSAR[J].Journal of Geomechanics, 2009, 15(2):151~161. doi: 10.3969/j.issn.1006-6616.2009.02.005 [3] 张永双, 石菊松, 孙萍, 等.汶川地震内外动力耦合及灾害实测[J].地质力学学报, 2009, 15(2):131~141. doi: 10.3969/j.issn.1006-6616.2009.02.003ZHANG Yong-shuang, SHI Ju-song, SUN Ping, et al.Coupling between endogenic and exogenic geological processes in the Wenchuan earthquake and example analysis of geo-hazards[J].Journal of Geomechanics, 2009, 15(2):131~141. doi: 10.3969/j.issn.1006-6616.2009.02.003 [4] 张永双, 雷伟志, 石菊松, 等.四川5.12地震次生地质灾害的基本特征初析[J].地质力学学报, 2008, 14(2):109~116. doi: 10.3969/j.issn.1006-6616.2008.02.002ZHANG Yong-shuang, LEI Wei-zhi, SHI Ju-song, et al.General characteristics of 5.12 earthquake-induced geo-hazards in Sichuan[J].Journal of Geomechanics, 2008, 14(2):109~116. doi: 10.3969/j.issn.1006-6616.2008.02.002 [5] 王连捷, 崔军文, 周春景, 等.汶川5.12地震发震机理的数值模拟[J].地质力学学报, 2009, 15(2):105~113. doi: 10.3969/j.issn.1006-6616.2009.02.001WANG Lian-jie, CUI Jun-wen, ZHOU Chun-jing, et al.Numerical modeling for Wenchuan earthquake mechanism[J].Journal of Geomechanics, 2009, 15(2):105~113. doi: 10.3969/j.issn.1006-6616.2009.02.001 [6] 王连捷, 周春景, 孙东升, 等.汶川5.12地震引起的库仑变化及其对周边地震活动的影响[J].地质力学学报, 2008, 14(3):193~200. doi: 10.3969/j.issn.1006-6616.2008.03.001WANG Lian-jie, ZHOU Chun-jing, SUN Dong-sheng, et al.Coulomb stress changes caused by Wenchuan earthquake and its influence on seismic activity in the adjacent area[J].Journal of Geomechanics, 2008, 14(3):193~200. doi: 10.3969/j.issn.1006-6616.2008.03.001 [7] 彭华, 马秀敏, 姜景捷, 等.龙门山北端青川断裂附近应力测量与断层稳定性[J].地质力学学报, 2009, 15(2):114~130. doi: 10.3969/j.issn.1006-6616.2009.02.002PENG Hua, MA Xiu-min, JIANG Jing-jie.Stability and stress measurement near the Qingchuan fault in the northern Longmen Mountains[J].Journal of Geomechanics, 2009, 15(2):114~130. doi: 10.3969/j.issn.1006-6616.2009.02.002 [8] 彭华, 马秀敏, 姜景捷, 等.山丹地应力监测站应变仪的地震效应[J].地质力学学报, 2008, 14(2):97~108. doi: 10.3969/j.issn.1006-6616.2008.02.001PENG Hua, MA Xiu-min, JIANG Jing-jie.Analysis of the volume strain data from the Shandan in-situ stress Monitoring Station[J].Journal of Geomechanics, 2008, 14(2):97~108. doi: 10.3969/j.issn.1006-6616.2008.02.001 [9] 马寅生, 张永双, 胡道功, 等.玉树地震地表破裂与宏观震中[J].地质力学学报, 2010, 16(2):115~124. doi: 10.3969/j.issn.1006-6616.2010.02.002MA Yin-sheng, ZHANG Yong-shuang, HU Dao-gong.et al.The Surface Ruptures And The Macroscopical Epicenter Of Yushu MS 7.1 Earthquake[J].Journal of Geomechanics, 2010, 16(2):115~124. doi: 10.3969/j.issn.1006-6616.2010.02.002 [10] 姚鑫, 张永双, 杨农, 等.玉树地震地表变形观测及初步分析[J].地质力学学报, 2010, 16(2):129~136. doi: 10.3969/j.issn.1006-6616.2010.02.003YAO Xin, ZHANG Yong-shuang, YANG Nong, et al.D-InSAR Deformation Observation and Preliminary Analysis of the MS 7.1 Yushu Earthquake[J].Journal of Geomechanics, 2010, 16(2):129~136. doi: 10.3969/j.issn.1006-6616.2010.02.003 [11] 张永, 许力生, 陈运泰. 2010年4月14日青海玉树MS 7. 1级地震[EB/OL]. http://www.cea-igp.ac.cn/近期显著地震研究/标题. htm, 2010-04-17.ZHANG Yong, XU Li-sheng, CHEN Yun-tai. MS 7. 1 earthquake in Yushu of Qinghai Province on April 12 th 2010[EB/OL]. http://www.cea-igp.ac.cn/, 2010-04-17 [12] 马寅生, 张永双, 胡道功. 地质力学研究所专家组青海玉树地震断裂调查[EB/OL]. http://www.geomech.ac.cn/xinwen/2010/yushudizhen/yushu0417.htm, 2010-04-17.MA Yin-sheng, ZHANG Yong-shuang, HU Dao-gong. Experts from Institute of Geomechanics investigated earthquake faults in Yushu, Qinghai Province[EB/OL]. http://www.geomech.ac.cn/xinwen/2010/yushudizhen/yushu0417.htm, 2010-04-17. [13] 陈颙, 地壳岩石的力学性能[M].北京:地震出版社, 1988.CHEN Yong.Mechanical properties of Crustal rocks[M].Beijing:Seismological Publishing House, 1988. [14] J C耶格, N G W库克. 岩石力学基础[M]. 中国科学院工程力学研究所(译), 北京: 地震出版, 1981.Jaeger J C, Cook N G W. Basis of rock mechanics[M]. Institute of engineering mechanics, Chinese Academy of Science, (Trans. ). Beijing: Seismological Publishing House, 1981. [15] Rice J R.Theory of precursory processes in the inception of earthquake rupture[J].Gerlands Beitraege Fuer Geophysik, 1979, 88, 91~127. https://www.researchgate.net/publication/220010732_Theory_of_precursory_process_in_the_inception_of_earthquake_rupture [16] 谭成轩, 马秀敏, 杨农, 等. 地质力学研究所地应力监测台站获得青海玉树MS 7. 1级地震信息[EB/OL]. http://www.geomech.ac.cn/xinwen/2010/yushudizhen/yushujianbao.htm, 2010-04-15TAN Cheng-xuan, MA Xiu-min, YANG Nong, et al. Information of MS 7. 1 earthquake in Yushu of Qinghai Province was captured by the monitor station of crustal stress from the Institute of Geomechanics[EB/OL]. http://www.geomech.ac.cn/xinwen/2010/yushudizhen/yushujianbao.htm, 2010-04-15 [17] 美国地调局. Magnitude 6. 9-Southern Qinghai, China[EB/OL]. http://earthquake.usgs.gov/earthquakes/recenteqsww/Quakes/us2010vacp.php, 2010-04-13.U S Geological Survey. Magnitude 6. 9-Southern Qinghai, China[EB/OL]. http://earthquake.usgs.gov/earthquakes/recenteqsww/Quakes/us2010vacp.php, 2010-04-13. [18] 彭华, 马秀敏, 李金锁, 等.南水北调西线一期工程地壳稳定性研究[M].北京:地震出版社, 2009.PENG Hua, MA Xiu-min, LI Jin-suo, et al.Study on crustal stability of West Route of South-to-North Water Transfer Project[M].Beijing:Seismological Publishing House, 2009. [19] 赖锡安, 黄立文, 徐菊生, 等.中国大陆现今地壳运动[M].北京:地震出版社, 2004.LAI Xi-an, HUANG Li-wen, XU Ju-sheng, et al.Present crustal movement in China[M].Beijing:Seismological Publishing House, 2004. [20] Goodman R E, Taylor R L, Brekke T L.A model for the mechanics of jointed rock[J].Journal of Soil Mech.Found.Div., ASCE, 1968, 94(SM3):637~659. http://d.old.wanfangdata.com.cn/Periodical/zgkydxxb201106009 [21] Zienkiewicz O C, Taylor R L.The finite element method (Fifth Edition), Volume 1:The basis, Volume 2:Solid mechanics, Volume 3:Fluid dynamics, Oxford:Butterworth-Heineman, 2000. [22] 殷有泉, 张宏.模拟地震的应变软化的数学模型[J].地球物理学报, 1982, 25(5):414~423. doi: 10.3321/j.issn:0001-5733.1982.05.004YIN You-quan, ZHANG Hong.A mathematical model of strain softening in simulating earthquake[J].Chinese Journal of Sinica, 1982, 25(5):414~423. doi: 10.3321/j.issn:0001-5733.1982.05.004 [23] 殷有泉, 张宏.断裂带内介质的软化特性和地震的非稳定模型[J].地震学报, 1984, 6(2):135~145. http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=QK000000375734YIN You-quan, ZHANG Hong.The softening behavior of fault zone medium and an instability model of earthquakes[J].Acta Seismologica Sinica, 1984, 6(2):135~145. http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=QK000000375734 [24] 殷有泉.非线性有限元基础[M].北京:北京大学出版社, 2007.YIN You-quan.Nonlinear finite element basis[M].Beijing:Peking University Press, 2007. [25] Stuart W D.Quasi-static earthquake mechanics[J].Reviews of Geophysics and Space Physics, 1979, (17):115~1120. http://d.old.wanfangdata.com.cn/Periodical/ytlx201012023 -