基于三维粘弹性有限元研究汶川地震对川滇地区的震后影响
蒋锋云 朱良玉 李玉江
摘要：在考虑川滇地区地壳介质的横向分块、纵向分层特征，及汶川地震同震断层面静态位错上下盘非对称性的基础上，构建川滇地区地壳三维粘弹性有限元模型，研究了汶川地震对川滇地区的震后影响。结果表明：（1）震后l0 a由于介质粘弹性松弛效应产生的远场水平变形在川滇菱形地块内部约0～20mm，对其北部影响较大，南部较刁、，其产生的远场垂直变形较小，在川滇菱形地块大部分区域表现为0～4mm的上升：（2）离发震断层较近的鲜水河断裂、东昆仑断裂，震后水平和垂直形变呈现高梯度带：（3）为了和汶川地震之后川滇地区强震活动进行对比，计算了主要块体边界活动断裂水平应力状态和库仑破裂应力变化。从同震及震后位移场空间分布特征、应力状态反映的断层活动特征及库仑破裂应力变化反映的断层面应力加卸载情况来看，与该区域强震活动分布在空间上存在一定的相关性。
关键词：汶川地震；有限元；粘弹性松弛；川滇地区
中图分类号：P315.725 文献标识码：A 文章编号：1000-0666（2018）02-0233-11
0 引言
自前，普遍的认识是一次大地震震后短期发震断层近场的地壳变形主要由发生在断层面及其延伸部分的滑移（震后余滑）和孔隙介质反弹控制（Jouanne et al，2011；Jonsson et al，2003），而震后长期的远场变形则主要受脆性上地壳下部的粘性下地壳和上地幔内部同震应力变化导致的粘弹性松弛效应控制（Rice，Gu，1983）。而有关研究表明，由于介质粘弹性松弛效应引起的断层面应力变化（库仑破裂应力变化）甚至比同震影响更为显著，可以导致周围断层失稳，从而触发强震。如沈正康等（2003）模拟了东昆仑断裂带1937年以来的5次7级以上强震在黏弹性成层介质中地震断层错动产生的应力演化过程，并计算了在后续地震破裂面上产生的库仑破裂应力变化，发现前4次地震均使2001年昆仑山口西8.1级大震发震断层面库仑破裂应力增加，且中地壳和下地壳的黏弹性松弛效应使得库仑破裂应力场随着时间的推移而逐渐加强，前3个地震震后松弛效应相比同震影响更为显著。万永革等（2007）采用多层弹粘性介质中的位错模型，以M≥7.0大地震和GPS数据得出的长期构造加载为形变源，计算得到了地震产生的应力变化，给出了青藏高原东北部1920年以来积累库仑破裂应力演化，结果表明：17次大地震均发生在库仑破裂应力变化为正的区域，触发率达85%。
2008年5月12日汶川8.0级特大地震发生至今已接近10年，有关该地震引起的震后远场变形的研究不多，少数的研究（谈洪波等，2009；李强等，2013）局限于利用弹性位错模型，仅考虑介质的纵向分层而没有考虑研究区复杂的介质横向非均匀性，而川滇地区相比华南地块普遍存在流变性更强的中下地壳，对震后形变场产生重要的影响（季灵运等，2015）。相比位错模型，有限元技术由于能充分考虑地球介质的非均匀性、调整断层面上下盘同震静态位移而在数值模拟中得到广泛的应用。因此，本文借助有限元技术，综合考虑活动地块、地球物理、大地测量、地震活动等多学科资料，构建考虑地壳介质纵向分层、横向分块的三维粘弹性有限元模型，以USGS给出的汶川地震同震破裂（Ji，Hayes，2008）作为加载，研究汶川地震对川滇地区由介质粘弹性松弛效应引起的震后形变场、应力场的时空分配关系。
汶川地震之后，川滇地区中强地震活动相比地震前出现了一些显著的变化，如鲜水河断裂（3次MS≥5.0地震）、大凉山地块附近（除芦山地震外，发生7次MS≥5.0地震）强震活动较为活跃。本文利用有限元模拟结果计算了汶川地震同震及震后不同时段研究区域主要块体边界活动断裂应力状态和库仑应力变化，并结合震后形变应变场时空演化规律，进一步分析了汶川地震对川滇地区的震后影响和地震活动之间的相关性。
川滇地区三维粘弹性有限元模型
根据川滇地区全新世活动断裂、地壳上地幔三维波速结构等资料，充分考虑对地质构造运动和地震活动起决定作用的活动断裂和活动断块，建立川滇地区三维地质构造模型。
1.1 介质分区及物性参数
川滇地区的活动断裂控制了主要构造活动的发生和发展，同时成为活动地块的边界（张培震等，2003），依据川滇地区的构造特征模型划分出10个块体，具体见图1。
构建有限元几何模型时，模型横向分块主要依照前述10个块体划分结果（图1）设定，纵向分层及各块体介质参数主要根据川滇地区地震波速资料揭示的地壳上地幔结构和P波、S波速度结果设定（王椿墉等，2015；吴建平等，2006，2009），具体设定为：（1）吴建平等（2006，2009）研究结果表明川滇地区中下地壳普遍存在低速层，而华南地块并不明显。因此，在介质参数设定时，考虑了中下地壳低速体的存在及其横向差异性，将巴颜喀拉块体、川西北块体、滇中块体、大凉山块体、保山块体、景谷块体6个块体中下地壳粘滞系数设定为1.0×1018Pa·s，而稳定的华南地块设定为1.0×1021Pa·s，其他區域为1.0×1020Pa·s；（2）已有研究（王达镛等，2015；吴建平等，2006，2009）表明，南北地震带中南段上地壳统一为强度较高的介质，横向和构造存在显著的相关性。因此，上地壳统一设定为较硬且粘滞系数较高的脆性层，但在弹性模量的设置上体现了华南地块相比其他区域介质相对较硬的特性：（3）整个区域上地幔具有相同的介质参数，王辉等（2007）根据地震波速进行统计计算得，粘滞系数约为1.0×1020Pa·s：（4）由于需要在龙门山断裂上加载汶川地震的同震位移，而同震滑动可以看作无摩擦的瞬间滑动，因此龙门山发震断层采用无摩擦的接触面来设计，具体的断层深度、倾角的设定和滑动分布的加载依据Ji和Hayes（2008）给出的有限断层面解来进行：（5）为简化模型，研究区其他断裂均采用倾角为90°的直立断层，参考王辉等（2007）的做法，断层的设定用粘滞系数较小（1.0×1020Pa·s）的软弱体代替，具体模型分块及分层参数见表1。为了尽可能地减少模型边界效应对计算结果的影响，我们将研究区域向外围扩展了200km。在上述地质模型的基础上，借助ANSYS有限元软件建立川滇地区三维粘弹性有限元模型（图2）。模型研究范围为（94°～110°E，21°～34°N），垂直区域为地表至上地鳗顶部100km深处，采用SOLID185和SOLID186这2种实体单元进行划分，整个有限元模型共划分555330个单元，737989个节点。
1.2 汶川地震同震发震断层模型
有关汶川地震的同震断层破裂模型已发表的研究成果有很多，Ji和Hayes（2008）的模型仅考虑中远场地震波，断层仅为一个面，比较简单。陈运泰等（2008）的研究结果虽然也只有一个断层面，但在资料的使用上，增加了近场的地震波资料约束。而之后的地震破裂模型由于考虑了大地测量、地震波、InSAR等多种观测（毕研磊等，2017；Shen et al，2009；Wang et al，2011；Field-ing et al，2013；Wan et al，2017），断层分了很多段，且考虑了断层深浅部倾角的差异，相对比较复杂，也更为精确。谈洪波等（2009）对陈运泰等（2008）以及ii和Hayes（2008）研究中的不同模型進行比较，发现不同模型同震及震后影响差异主要集中在断层附近，中远场差异较小。考虑到本文模型空间尺度范围较大，且主要考查中远场同震及震后形变，为了使有限元模型建立过程中同震断层面处理更方便，本文采用Ji和Hayes（2008）的简单模型（图3），这种模型在断层附近和实际结果存在一定的差异，但对于中远场的分析而言，差异可以忽略。图3中箭头表示下盘相对上盘的运动方向，最优选取的发震断层参数为走向229°，倾角33°，分别沿走向和倾向均匀的分成21块和8块子断层，每个子断层的尺度为15km×5km，并给出了每个子断层面的滑动角和滑动量。
1.3 模型边界约束与同震加载方式
为了得到汶川地震对川滇地区真实的震后形变场，必须先得到和实际观测接近的同震位移场。Wang等（2011）为我们提供了可靠的、高密度的、大范围的汶川地震的同震观测水平位移场（图4a），其结果表明：（1）汶川地震同震水平形变存在明显的非对称性，即上盘（青藏高原一侧）形变强度和范围明显大于下盘（华南地块一侧），这种现象在INSAR同震近场形变观测中也可以看到（Shen et al，2009）：（2）汶川地震同震影响在川滇菱形地块东边界以西从南至北，具有和背景构造水平运动场接近的顺时针旋转特征，在滇西北地块表现为近南北向运动，而在滇中地块表现出明显的南西向整体运动，远至滇西北和滇中块体南部仍然有将近10mm的同震形变量。
相关研究（Wang et al，2007；Funning et al，2007；邵志刚等，2008；李宁等，2017；季灵运等，2017，焦佳爽等，2017）表明，同震形变场的非对称性成因机理主要有以下几种：（1）发震断层两盘介质性质存在显著差异：（2）震后粘弹性松弛引起的有效粘度发生横向变化：（3）断层存在多条近平行的剪切带：（4）断层产状，即断层并非直立断层，存在倾角。对于汶川地震的发震断层龙门山断裂而言，断层两侧青藏高原相对华南地块明显的介质性质差异和断层产状都可造成同震形变场的非对称性。而第2种情况，相对同震短时间而言，影响可忽略不计。由于龙门山断裂相对比较集中，第3种情况也可忽略不计。我们在有限元模型构建时就考虑发震断层330的北西倾向，也考虑了发震断层两盘介质性质的差异。在此基础上，在边界条件的设定上，参照李玉江等（2013）将模型四周和底部做法向约束，其他方向可自由移动。目前大多数同震静态位移断层面加载方式为上下两盘断层面对称加载（李玉江等，2013），即上下两盘同震静态位移各占一半。通过试算，我们发现这种加载方式对近场的拟合效果较好，而对中远场拟合效果较差，无法完全体现出上下两盘形变的非对称性。通过分析，我们认为除了上述4种情况，同震形变的非对称性可能与断层两盘本身的运动状态有关，即断层的主动盘与被动盘。青藏高原受到印度板块北东向推挤，在外围构造环境约束下，龙门山以西表现为垂直龙门山断裂近东西向的运动特征，而华南地块，特别是稳定的四川盆地，相对青藏高原而言是稳定不动的。这表明龙门山断裂上盘青藏高原一侧是主动盘，而华南地块一侧是相对稳定的被动盘。换句话说，华南地块和青藏高原本身的构造边界约束条件是不一致的，这种不一致性有可能造成上盘断层面同震滑移要比下盘断层面大，从而使得同震形变在更大范围内表现出上盘比下盘显著的特征。因此，我们拟通过调整断层上下两盘同震静态滑移比例来实现模拟同震水平位移（主要是中远场）和实际观测的最佳匹配，在实际给定比例时，要保证上盘一直大于下盘，这样才符合物理机理。试算结果表明上下两盘同震滑移比例为7.5：2.5时获得的同震水平移场（图4b）和实际观测结果（图4a）较为接近。
需要强调的是虽然受到有限元模型边界条件、地壳介质的简化以及发震断层产状与同震滑动分布简略等影响，模拟结果和观测结果仍然存在显著的差异，从模拟值与观测值矢量差可以更为清晰地看出（图4c），差异主要表现为：（1）在断层附近模拟值比实际观测值要小，这可能是由于使用的同震破裂模型是用远场地震波资料反演得到，而没有近场的约束引起的：（2）整个川滇地区和四川盆地远离发震断裂地区，模拟得到速度场比观测速度场量值大，且存在系统性偏差。其原因一方面可能是有限元模拟中边界条件仅做法向约束，和实际情况可能存在差异，另一方面可能是上下盘断层面同震静态位移比例调整，仅用一个比例系数，过于简单：（3）整个区域除了发震断层附近上下两盘一定范围内，大部分地区模拟得到的速度场和实际观测存在10mm以下的残差。
虽然模拟结果与观测结果存在上述明显的差异，但其反映的整个同震形变场的整体运动形态，特别是中远场形态和实际观测结果（图4a）具有较好的相似性。因此，利用上述得到的上下盘同震加载比例及边界约束方式进行汶川地震同震及震后形变场的研究具有一定合理性。
2 汶川地震对川滇地区震后影响
采用上述边界条件及断层面静态同震破裂加载，计算了汶川地震同震及震后1a、5a、10a这3个时间点由于介质粘弹性松弛效应引起的川滇地区位移场以及主要块体边界活动断裂水平主应力状态和库仑应力的变化，并结合川滇地区近年来强震活动特征，对这些物理量的时空演化特征进行分析验证。
从模拟同震水平位移等值线图（图5a）可以看出，和实际观测结果类似，汶川地震的同震影响主要沿垂直发震断裂两侧方向延伸，上盘形变强度和范围明显大于下盘，上下两盘存在明显的非对称性。同时平行断裂方向也存在一定分量，表现出同震形变主要以断裂逆冲为主，同时兼具右旋走滑的特征。需要说明的是模拟得到的汶川地震对川滇交界东部地区同震影响比实际观测结果量值偏大。
从震后水平形变场来看，震后1a（图5b）形变主要集中在发震断裂附近，且上盘变形强度明显大于下盘。但就影响范围来看，下盘影响范围比上盘要大。震后5a、10a形变（图5c、d）则在短期形变的基础上随时间推移向外围扩展。其中从川滇交界东部地区向滇西南、滇西北方向扩展特征更为明显，这和汶川地震之后川滇地区中强地震活动主要呈北東向条带分布一致。另外，在形变扩展过程中，沿着巴颜喀拉地块南北两侧边界带出现形变高梯度带，且随时间变化逐渐加强。在这一过程中，先后在其南边界鲜水河断裂上发生了2010年道孚5.4级、2011年炉霍5.7级、2014康定6.3级强震，近期又在其北边界东北角发生了2017年8月8日九寨沟7.0级强震，可能与汶川地震对该区域的震后影响有关。大凉山地块震后形变随着时间变化增强显著，水平变形量震后10a达到20mm，该区域也是汶川震后强震活跃显著的区域。
从垂直形变特征（图6）来看，和水平同震形变类似，同震垂直形变主要沿着龙门山断裂东西两侧垂直断层方向延伸，靠近断层附近下盘表现为下沉，上盘表现为上升。距断层一定范围内（上盘约300km，下盘约150km）表现为上升，之后向两侧延伸又表现为下沉。随着离开断层距离的增加，同震垂直形变衰减较快。断裂西南部表现为显著的上升，而断裂的东北部表现为上升-下沉交替变化特征。震后1a尺度垂直形变（扣除同震）则主要集中在断层附近较小的区域，靠近断层附近上盘和同震形变一致继续表现为上升，下盘仍然表现为下沉。相对远离断层区域和同震形变相反，表现为下沉区。随着时间的推移，震后垂直形变不断向外扩展，和水平形变类似，扩展方向除了沿着断裂两侧垂直断裂方向外，还具有从川滇交界东部向滇西北、滇西南方向扩展的特征。
从同震及震后主要块体边界水平向主应力变化特征（图7）来看：（1）甘孜—玉树断裂南段、鲜水河断裂炉霍—道孚—乾宁段以及金沙江断裂北段同震影响以垂直断裂方向张应力为主，震后随时间推移逐渐加强。尤其是鲜水河断裂，震后10a的张应力累积变化量已超过同震变化量。而该断裂上相继发生了2010年道孚5.4级、2011年炉霍5.7级、2014康定6.3级强震，随时间的推移震级呈增强特征，和汶川地震震后应力场变化具有显著相关性；（2）安宁河—则木河—小江断裂、马边断裂、红河断裂、金沙江断裂南段同震及震后影响均表现为与断层走向呈一定夹角的庄应力变化特征。安宁河—则木河—小江断裂主压应力方向（北北东）和背景构造主压应力方向（北北西向）相反，使断层表现为右旋，具有负影响，汶川震后强震活动较弱。红河断裂应力所反映的断层走滑活动特征和背景特征一致，表现为右旋走滑，但量值相对较小。金沙江断裂中南段以垂直断裂的压应力特征为主：（3）总体上，汶川地震同震及震后对鲜水河断裂影响最大，其次为金沙江断裂北段、安宁河—则木河断裂、马边断裂、鲁甸—昭通断裂、金沙江断裂中南段、红河断裂、丽江—小金河断裂，对小江断裂影响最小。从影响特征和影响强度2方面考虑，汶川震后川滇地区强震活动的主体区域和其具有较好的一致性。
为了进一步分析汶川地震对川滇地区主要块体边界活动断裂的影响，将汶川地震产生的应力变化投影到断裂滑动方向上，用库仑破裂应力变化来表示该地震对周围接收断裂产生的影响。具体计算公式参照Jaeger等（1969）的研究，在计算中有效摩擦系数参照Stein等（1992）和King等（1994）的研究，取值为0.4。有关接收断裂活动性质几何参数参考邓起东等（2002）结果给出。由于研究区域的地震震源深度绝大多数在0～25km范围内，因此我们研究汶川地震在10km深度上产生的库仑破裂应力变化，我们也计算了其他深度上产生的库仑破裂应力变化，发现模式和量值变化不大。
从库仑破裂应力变化结果（图8）来看：（1）同震产生的库仑破裂应力显著增加断裂主要包括龙门山断裂南端和北端、甘孜—玉树断裂、鲜水河断裂、白玉断裂、东昆仑断裂东段部分段落，这些部位库仑破裂应力增加值均超过0.01MPa，达到强震触发的阈值（King et al，1994；Harris，Simgson，1998）。震后随时间的变化，这些断裂库仑破裂应力大部分在同震的基础上有所增加，但增加量值不太显著。汶川地震之后，这些断裂相继有5.0级以上强地震发生，其中包括在甘孜一玉树断裂发生的2010年玉树7.1级、龙门山断裂南段发生的2013年芦山7.0级和鲜水河断裂发生的2014年康定6.3级地震。（2）同震库仑应力有所增加，但增加不是很显著的断裂，主要包括大凉山地块东边界马边断裂、大凉山地块南边界鲁甸一昭通断裂、金沙江断裂、丽江小金河断裂的南段、南汀河断裂，库仑应力增加值小于0.01MPa。鲁甸一昭通断裂库仑破裂应力震后随时间变化逐渐增强，震后10a达到0.01MPa，该断裂附近先后发生了几次5.0级以上强震，其中包括2014年8月3日鲁甸6.5级地震。丽江一小金河断裂震后随时间推移库仑破裂应力也逐渐增强，其南段震后10a库仑应力接近0.01MPa，在其附近发生了2012年宁蒗5.7级地震。金沙江断裂震后随时间变化库仑破裂应力有所增强，南汀河断裂震后随时间变化库仑应力有所增加，但并不显著。（3）同震库仑破裂应力减少的断裂包括龙门山断裂中段发震断裂、安宁河—则木河—小江断裂、丽江—小金河断裂北段、红河断裂及东昆仑断裂东段部分段落。除了丽江—小金河断裂北段在震后随时间推移，库仑破裂应力由同震的卸载转为震后的加载外，其他断裂震后随时间变化情况和同震效应一致，均存在不同程度的卸载增强的特征。总之，同震及震后断层库仑破裂应力反映的加卸载情况和强震活动空间分布具有一定的相关性。需要说明的是由于使用的模型及方法的不同，本文计算获得的同震库仑破裂应力变化和万永革等（2009）采用位错方法得到的结果存在明显的差异。
3 討论与结论
本文借助三维粘弹性有限元技术，以Wang等（2011）给出的大范围汶川地震同震水平位移场作为约束，模拟了汶川地震对川滇地区的同震及震后形变场、应力的时空演化特征，在模拟过程中，不但考虑了川滇地区地球介质的横向和纵向非均性，而且对同震上下盘断层面静态位移加载比例进行了调整，尽可能使我们的模拟结果在中远场和实际观测结果在形态上具有一致性。
模拟得到的川滇地区同震和震后形变场、应力状态反映的断层活动特征及库仑破裂应力变化反映的断层面应力加卸载情况，在空间上和强震活动具有较好的相关性。对于理解和解释汶川地震之后川滇地区出现的一些观测异常具有一定的参考意义。但本文模拟结果相比实际观测结果还存在许多不一致的地方，在模拟过程中仍然存在许多有待商榷和解决的问题。主要包括以下几个方面：（1）有限元边界约束条件对结果有着重要的影响，位错模型认为在离开断层无限远处和无限深处同震位移为0，即固定不动。在有限元模型中，由于模型范围总是有限的，我们无法做到这一点。本文采用大多数学者的做法，即模型边界法向约束，切向自由移动，这和实际情况存在明显的差异，也可能是本文模拟结果相比观测结果中远场同震形变相对较大的主要原因。因此探索更为有效、合理的边界条件约束是在模拟过程需要进一步解决的问题：（2）本文简单地通过采用上下盘断层面同震错动位移比例，来拟合测同震形变场，仅仅采用一个比例系数，过于简单，也可能是拟合效果存在差异的原因之一：（3）本文模型没有考虑重力、地壳背景构造应力场作用，重力、背景构造应力对同震及震后变形也存在一定的影响。另外，地壳介质本身具有复杂性（地形、断层的非连续性）等，这对模拟结果都会产生一定影响。在后续的数值模拟工作中需要进一步探索和研究。
参考文献：
毕研磊，张永志，曹海坤，等.2017.利用断层自剖分技术反演汶川MW7.9地震滑动分布[J].40（2）：211-215.
陈运泰，许力生，张勇，等.2008.2008年5月12日汶川特大地震震源特性分析报告[EB/OL].（2017-09-20）[2008-05].ht-tp：//www.csi.ac.cn/sichuan/.
邓起东，张培震，冉勇康，等.2002.中国活动构造基本特征[J].中国科学：地球科学，32（12）：1020-1030.
季灵运，刘传金，徐晶，等.2017.九寨沟MS7.0级地震的InSAR观测及发震构造分析[J].地球物理学报，60（10）：4069-4082.
季灵运，刘立炜，郝明.2015.利用InSAR技术研究滇西南镇康永德地区现今地壳形变特征[J].地震研究，38（1）：84-89.
焦佳爽，张永志，张凯南，等.2017.基于同震GPS数据的汶川大震同震应变特征研究[J].地震工程学报，39（3）：227-233.
李宁，赵强，李金.2017.基于InSAR形变观测反演2015年皮山MS6.5地震滑动分布[J].地震，37（2）：67-76.
李强，江在森，武艳强，等.2013.汶川地震震后短期青藏高原东北缘地壳变形特征成因讨论[J].地震，33（3）：124-132.
李玉江，陈连旺，杨树新.2013.基于应变能变化的芦山强震同震效应的数值模拟[J].震灾防御技术，8（4）：361-369.
邵志刚，傅容珊，薛霆虓，等.2008.昆仑山MS8.1级地震震后变形场数值模拟与成因机理探讨[J].地球物理学报，51（3）：805-816.
沈正康，万永革，甘卫军，等.2003.东昆仑活动断裂带大地震之间的黏弹性应力触发研究[J].地球物理学报，46（6）：786-795.
谈洪波，申重阳，李辉，等.2009.汶川大地震震后重力变化和形变的黏弹分层模拟[J].地震学报，31（5）：491-505.
万永革，沈正康，盛书中，等.2009.2008年汶川大地震对周围断层的影响[J]地震学报，31（2）：129-139.
万永革，沈正康，曾跃华，等.2007.青藏高原东北部的库仑应力积累演化对大地震发生的影响[J].地震学报，29（2）115-129.
主椿镛，杨文采，吴建平，等.2015.南北构造带岩石圈结构与地震的研究[J].地球物理学报，58（11）：3867-3901.
王辉，曹建玲，张怀，等，2007.川滇地区下地壳流动对上地壳运动变形影响的数值模拟[J].地震学报，29（6）：581-591.
吴建平，黄媛，张天中，等，2009.汶川MS8.0级地震余震分布及周边区域P波三维速度结构研究[J].地球物理学报，52（2）：320-328.
吴建平，明跃红，王椿镛.2006.川滇地区速度结构的区域地震波形反演研究[J].地球物理学报，49（5）：1369-1376.
张培震，邓起东，张国民，等.2003.中国大陆的强震活动与活动地块[J].中国科学：地球科学，33（增刊1）.12-20.
Fielding E J，Sladen A，Li Z，et aL 2013.Kipemalic fault slip evolutionsource models of he 2008 M7.9 Wenchuan.earthquake in Chinafrom SAIL interferometry，GPS and teleseismic analysis and implies-Lionsfor Lougmen Shan tectonics[J].Geophysical journal interns-tional，194（2）：1138-1166.
Funning G J，Parsons 13，Wright T J.2007.Fault slip in the 1997Manyi，Tibet earthquake from linear elastic modelling of InSAR displace-ments[J].Geophysical Journal International，169（3）：988-1008.
Ham's R A，Simpson R W.1998.Supp cession of large earquakes by stressshadow：A Conpaiison of Coulomb and rate-and-state failure[J].Journal of Geophysical Research solid Earth，103（B10）：24439-24451.
Jaeger J C，Cook N G W，Zimmerman R.1969.Fundamentals of rock me-chanics[M].Chapman and Hall'Landon.
Ji C，Hayes G.2008.NEIC，Preliminary result of theMay 12，2008MW7.9 eastern Sichuan，China earthquake[EB/OL].（2017-09-20）[2008-05].http；//earthquake.nags.gov/eqcenter/eqinthe-news/us2008ryan/.
Jonssons，Segall P，Pedersen R，et al.2003.Post-earthquake groundmovements correlated to pore-pressure transients[J].Nature，424（4945）：179-183：
Jouanne F，Awan A，Madji A，et al.2011.Postseismic deformation in Pa-kistan after the 8 October 2005 earthquake：evidence of afterslip a-long a flat north of the Balakot-Bagh thrust[J].Journal of Geo-physical Research：Solid Earth，116（B07401）：1-22.
King G C P，Stein Rs，Lin J.1994.Static stress changes and the trigge-ring of earthquakes[J].Bulletin of the Seismological Society of A-merica，84（3）：935-953.
Rice J R，Go J.1983.Earthquake after effects and triggered seismic phe-nomena[J].Pure and Applied Geophysics，121（2）：187-219.
Shen Z K，Sun J，Zhang P，et al.2009.Slip maxima at fault junctions andrupturing of barriers dining the 2008 Wenchuan earthquake[J].Na-ture geosoence，2009，2（10）：718-724.
Stein R S，King G C P，Lin J.1992.Change in failure stress on the south-ern San Andreas fault system caused by the 1992 magnitude=7.4Lenders earthquake[J].Science，258（5086）：1328-1332.
Wan Y，Shen Z K，Burgmann R，et at.2017.Fault geometry and slip dis-tribution of the 2008 MW7.9 Wenchuan，China earthquake，inferredfrom GPS and InSAR measurements[J].Geophysical Journal tnter-nationel，208（2）：748-766.
Wang，Xu C，Ge L.2007.Coseismic deformation and slip distribution ofthe 1997 MW7.5 Manyi，Tibet，earthquake from InSAR measure-ments[J].Journal of Geodynamics，44（3）：200-212.
Wang Q，Qiao X J，Len Q G，et al.2011.Rupture of deep faults in the2008 Wenchuan earthquake and uplift of the Longmen Shan[J].Na-Lure Ceascience，4（9）：634-640.
摘要：在考虑川滇地区地壳介质的横向分块、纵向分层特征，及汶川地震同震断层面静态位错上下盘非对称性的基础上，构建川滇地区地壳三维粘弹性有限元模型，研究了汶川地震对川滇地区的震后影响。结果表明：（1）震后l0 a由于介质粘弹性松弛效应产生的远场水平变形在川滇菱形地块内部约0～20mm，对其北部影响较大，南部较刁、，其产生的远场垂直变形较小，在川滇菱形地块大部分区域表现为0～4mm的上升：（2）离发震断层较近的鲜水河断裂、东昆仑断裂，震后水平和垂直形变呈现高梯度带：（3）为了和汶川地震之后川滇地区强震活动进行对比，计算了主要块体边界活动断裂水平应力状态和库仑破裂应力变化。从同震及震后位移场空间分布特征、应力状态反映的断层活动特征及库仑破裂应力变化反映的断层面应力加卸载情况来看，与该区域强震活动分布在空间上存在一定的相关性。
关键词：汶川地震；有限元；粘弹性松弛；川滇地区
中图分类号：P315.725 文献标识码：A 文章编号：1000-0666（2018）02-0233-11
0 引言
自前，普遍的认识是一次大地震震后短期发震断层近场的地壳变形主要由发生在断层面及其延伸部分的滑移（震后余滑）和孔隙介质反弹控制（Jouanne et al，2011；Jonsson et al，2003），而震后长期的远场变形则主要受脆性上地壳下部的粘性下地壳和上地幔内部同震应力变化导致的粘弹性松弛效应控制（Rice，Gu，1983）。而有关研究表明，由于介质粘弹性松弛效应引起的断层面应力变化（库仑破裂应力变化）甚至比同震影响更为显著，可以导致周围断层失稳，从而触发强震。如沈正康等（2003）模拟了东昆仑断裂带1937年以来的5次7级以上强震在黏弹性成层介质中地震断层错动产生的应力演化过程，并计算了在后续地震破裂面上产生的库仑破裂应力变化，发现前4次地震均使2001年昆仑山口西8.1级大震发震断层面库仑破裂应力增加，且中地壳和下地壳的黏弹性松弛效应使得库仑破裂应力场随着时间的推移而逐渐加强，前3个地震震后松弛效应相比同震影响更为显著。万永革等（2007）采用多层弹粘性介质中的位错模型，以M≥7.0大地震和GPS数据得出的长期构造加载为形变源，计算得到了地震产生的应力变化，给出了青藏高原东北部1920年以来积累库仑破裂应力演化，结果表明：17次大地震均发生在库仑破裂应力变化为正的区域，触发率达85%。
2008年5月12日汶川8.0级特大地震发生至今已接近10年，有关该地震引起的震后远场变形的研究不多，少数的研究（谈洪波等，2009；李强等，2013）局限于利用弹性位错模型，仅考虑介质的纵向分层而没有考虑研究区复杂的介质横向非均匀性，而川滇地区相比华南地块普遍存在流变性更强的中下地壳，对震后形变场产生重要的影响（季灵运等，2015）。相比位错模型，有限元技术由于能充分考虑地球介质的非均匀性、调整断层面上下盘同震静态位移而在数值模拟中得到广泛的应用。因此，本文借助有限元技术，综合考虑活动地块、地球物理、大地测量、地震活动等多学科资料，构建考虑地壳介质纵向分层、横向分块的三维粘弹性有限元模型，以USGS给出的汶川地震同震破裂（Ji，Hayes，2008）作为加载，研究汶川地震对川滇地区由介质粘弹性松弛效应引起的震后形变场、应力场的时空分配关系。
汶川地震之后，川滇地区中强地震活动相比地震前出现了一些显著的变化，如鲜水河断裂（3次MS≥5.0地震）、大凉山地块附近（除芦山地震外，发生7次MS≥5.0地震）强震活动较为活跃。本文利用有限元模拟结果计算了汶川地震同震及震后不同时段研究区域主要块体边界活动断裂应力状态和库仑应力变化，并结合震后形变应变场时空演化规律，进一步分析了汶川地震对川滇地区的震后影响和地震活动之间的相关性。
川滇地区三维粘弹性有限元模型
根据川滇地区全新世活动断裂、地壳上地幔三维波速结构等资料，充分考虑对地质构造运动和地震活动起决定作用的活动断裂和活动断块，建立川滇地区三维地质构造模型。
1.1 介质分区及物性参数
川滇地区的活动断裂控制了主要构造活动的发生和发展，同时成为活动地块的边界（张培震等，2003），依据川滇地区的构造特征模型划分出10个块体，具体见图1。
构建有限元几何模型时，模型横向分块主要依照前述10个块体划分结果（图1）设定，纵向分层及各块体介质参数主要根据川滇地区地震波速资料揭示的地壳上地幔结构和P波、S波速度结果设定（王椿墉等，2015；吴建平等，2006，2009），具体设定为：（1）吴建平等（2006，2009）研究结果表明川滇地区中下地壳普遍存在低速层，而华南地块并不明显。因此，在介质参数设定时，考虑了中下地壳低速体的存在及其横向差异性，将巴颜喀拉块体、川西北块体、滇中块体、大凉山块体、保山块体、景谷块体6个块体中下地壳粘滞系数设定为1.0×1018Pa·s，而稳定的华南地块设定为1.0×1021Pa·s，其他區域为1.0×1020Pa·s；（2）已有研究（王达镛等，2015；吴建平等，2006，2009）表明，南北地震带中南段上地壳统一为强度较高的介质，横向和构造存在显著的相关性。因此，上地壳统一设定为较硬且粘滞系数较高的脆性层，但在弹性模量的设置上体现了华南地块相比其他区域介质相对较硬的特性：（3）整个区域上地幔具有相同的介质参数，王辉等（2007）根据地震波速进行统计计算得，粘滞系数约为1.0×1020Pa·s：（4）由于需要在龙门山断裂上加载汶川地震的同震位移，而同震滑动可以看作无摩擦的瞬间滑动，因此龙门山发震断层采用无摩擦的接触面来设计，具体的断层深度、倾角的设定和滑动分布的加载依据Ji和Hayes（2008）给出的有限断层面解来进行：（5）为简化模型，研究区其他断裂均采用倾角为90°的直立断层，参考王辉等（2007）的做法，断层的设定用粘滞系数较小（1.0×1020Pa·s）的软弱体代替，具体模型分块及分层参数见表1。为了尽可能地减少模型边界效应对计算结果的影响，我们将研究区域向外围扩展了200km。在上述地质模型的基础上，借助ANSYS有限元软件建立川滇地区三维粘弹性有限元模型（图2）。模型研究范围为（94°～110°E，21°～34°N），垂直区域为地表至上地鳗顶部100km深处，采用SOLID185和SOLID186这2种实体单元进行划分，整个有限元模型共划分555330个单元，737989个节点。
1.2 汶川地震同震发震断层模型
有关汶川地震的同震断层破裂模型已发表的研究成果有很多，Ji和Hayes（2008）的模型仅考虑中远场地震波，断层仅为一个面，比较简单。陈运泰等（2008）的研究结果虽然也只有一个断层面，但在资料的使用上，增加了近场的地震波资料约束。而之后的地震破裂模型由于考虑了大地测量、地震波、InSAR等多种观测（毕研磊等，2017；Shen et al，2009；Wang et al，2011；Field-ing et al，2013；Wan et al，2017），断层分了很多段，且考虑了断层深浅部倾角的差异，相对比较复杂，也更为精确。谈洪波等（2009）对陈运泰等（2008）以及ii和Hayes（2008）研究中的不同模型進行比较，发现不同模型同震及震后影响差异主要集中在断层附近，中远场差异较小。考虑到本文模型空间尺度范围较大，且主要考查中远场同震及震后形变，为了使有限元模型建立过程中同震断层面处理更方便，本文采用Ji和Hayes（2008）的简单模型（图3），这种模型在断层附近和实际结果存在一定的差异，但对于中远场的分析而言，差异可以忽略。图3中箭头表示下盘相对上盘的运动方向，最优选取的发震断层参数为走向229°，倾角33°，分别沿走向和倾向均匀的分成21块和8块子断层，每个子断层的尺度为15km×5km，并给出了每个子断层面的滑动角和滑动量。
1.3 模型边界约束与同震加载方式
为了得到汶川地震对川滇地区真实的震后形变场，必须先得到和实际观测接近的同震位移场。Wang等（2011）为我们提供了可靠的、高密度的、大范围的汶川地震的同震观测水平位移场（图4a），其结果表明：（1）汶川地震同震水平形变存在明显的非对称性，即上盘（青藏高原一侧）形变强度和范围明显大于下盘（华南地块一侧），这种现象在INSAR同震近场形变观测中也可以看到（Shen et al，2009）：（2）汶川地震同震影响在川滇菱形地块东边界以西从南至北，具有和背景构造水平运动场接近的顺时针旋转特征，在滇西北地块表现为近南北向运动，而在滇中地块表现出明显的南西向整体运动，远至滇西北和滇中块体南部仍然有将近10mm的同震形变量。
相关研究（Wang et al，2007；Funning et al，2007；邵志刚等，2008；李宁等，2017；季灵运等，2017，焦佳爽等，2017）表明，同震形变场的非对称性成因机理主要有以下几种：（1）发震断层两盘介质性质存在显著差异：（2）震后粘弹性松弛引起的有效粘度发生横向变化：（3）断层存在多条近平行的剪切带：（4）断层产状，即断层并非直立断层，存在倾角。对于汶川地震的发震断层龙门山断裂而言，断层两侧青藏高原相对华南地块明显的介质性质差异和断层产状都可造成同震形变场的非对称性。而第2种情况，相对同震短时间而言，影响可忽略不计。由于龙门山断裂相对比较集中，第3种情况也可忽略不计。我们在有限元模型构建时就考虑发震断层330的北西倾向，也考虑了发震断层两盘介质性质的差异。在此基础上，在边界条件的设定上，参照李玉江等（2013）将模型四周和底部做法向约束，其他方向可自由移动。目前大多数同震静态位移断层面加载方式为上下两盘断层面对称加载（李玉江等，2013），即上下两盘同震静态位移各占一半。通过试算，我们发现这种加载方式对近场的拟合效果较好，而对中远场拟合效果较差，无法完全体现出上下两盘形变的非对称性。通过分析，我们认为除了上述4种情况，同震形变的非对称性可能与断层两盘本身的运动状态有关，即断层的主动盘与被动盘。青藏高原受到印度板块北东向推挤，在外围构造环境约束下，龙门山以西表现为垂直龙门山断裂近东西向的运动特征，而华南地块，特别是稳定的四川盆地，相对青藏高原而言是稳定不动的。这表明龙门山断裂上盘青藏高原一侧是主动盘，而华南地块一侧是相对稳定的被动盘。换句话说，华南地块和青藏高原本身的构造边界约束条件是不一致的，这种不一致性有可能造成上盘断层面同震滑移要比下盘断层面大，从而使得同震形变在更大范围内表现出上盘比下盘显著的特征。因此，我们拟通过调整断层上下两盘同震静态滑移比例来实现模拟同震水平位移（主要是中远场）和实际观测的最佳匹配，在实际给定比例时，要保证上盘一直大于下盘，这样才符合物理机理。试算结果表明上下两盘同震滑移比例为7.5：2.5时获得的同震水平移场（图4b）和实际观测结果（图4a）较为接近。
需要强调的是虽然受到有限元模型边界条件、地壳介质的简化以及发震断层产状与同震滑动分布简略等影响，模拟结果和观测结果仍然存在显著的差异，从模拟值与观测值矢量差可以更为清晰地看出（图4c），差异主要表现为：（1）在断层附近模拟值比实际观测值要小，这可能是由于使用的同震破裂模型是用远场地震波资料反演得到，而没有近场的约束引起的：（2）整个川滇地区和四川盆地远离发震断裂地区，模拟得到速度场比观测速度场量值大，且存在系统性偏差。其原因一方面可能是有限元模拟中边界条件仅做法向约束，和实际情况可能存在差异，另一方面可能是上下盘断层面同震静态位移比例调整，仅用一个比例系数，过于简单：（3）整个区域除了发震断层附近上下两盘一定范围内，大部分地区模拟得到的速度场和实际观测存在10mm以下的残差。
虽然模拟结果与观测结果存在上述明显的差异，但其反映的整个同震形变场的整体运动形态，特别是中远场形态和实际观测结果（图4a）具有较好的相似性。因此，利用上述得到的上下盘同震加载比例及边界约束方式进行汶川地震同震及震后形变场的研究具有一定合理性。
2 汶川地震对川滇地区震后影响
采用上述边界条件及断层面静态同震破裂加载，计算了汶川地震同震及震后1a、5a、10a这3个时间点由于介质粘弹性松弛效应引起的川滇地区位移场以及主要块体边界活动断裂水平主应力状态和库仑应力的变化，并结合川滇地区近年来强震活动特征，对这些物理量的时空演化特征进行分析验证。
从模拟同震水平位移等值线图（图5a）可以看出，和实际观测结果类似，汶川地震的同震影响主要沿垂直发震断裂两侧方向延伸，上盘形变强度和范围明显大于下盘，上下两盘存在明显的非对称性。同时平行断裂方向也存在一定分量，表现出同震形变主要以断裂逆冲为主，同时兼具右旋走滑的特征。需要说明的是模拟得到的汶川地震对川滇交界东部地区同震影响比实际观测结果量值偏大。
从震后水平形变场来看，震后1a（图5b）形变主要集中在发震断裂附近，且上盘变形强度明显大于下盘。但就影响范围来看，下盘影响范围比上盘要大。震后5a、10a形变（图5c、d）则在短期形变的基础上随时间推移向外围扩展。其中从川滇交界东部地区向滇西南、滇西北方向扩展特征更为明显，这和汶川地震之后川滇地区中强地震活动主要呈北東向条带分布一致。另外，在形变扩展过程中，沿着巴颜喀拉地块南北两侧边界带出现形变高梯度带，且随时间变化逐渐加强。在这一过程中，先后在其南边界鲜水河断裂上发生了2010年道孚5.4级、2011年炉霍5.7级、2014康定6.3级强震，近期又在其北边界东北角发生了2017年8月8日九寨沟7.0级强震，可能与汶川地震对该区域的震后影响有关。大凉山地块震后形变随着时间变化增强显著，水平变形量震后10a达到20mm，该区域也是汶川震后强震活跃显著的区域。
从垂直形变特征（图6）来看，和水平同震形变类似，同震垂直形变主要沿着龙门山断裂东西两侧垂直断层方向延伸，靠近断层附近下盘表现为下沉，上盘表现为上升。距断层一定范围内（上盘约300km，下盘约150km）表现为上升，之后向两侧延伸又表现为下沉。随着离开断层距离的增加，同震垂直形变衰减较快。断裂西南部表现为显著的上升，而断裂的东北部表现为上升-下沉交替变化特征。震后1a尺度垂直形变（扣除同震）则主要集中在断层附近较小的区域，靠近断层附近上盘和同震形变一致继续表现为上升，下盘仍然表现为下沉。相对远离断层区域和同震形变相反，表现为下沉区。随着时间的推移，震后垂直形变不断向外扩展，和水平形变类似，扩展方向除了沿着断裂两侧垂直断裂方向外，还具有从川滇交界东部向滇西北、滇西南方向扩展的特征。
从同震及震后主要块体边界水平向主应力变化特征（图7）来看：（1）甘孜—玉树断裂南段、鲜水河断裂炉霍—道孚—乾宁段以及金沙江断裂北段同震影响以垂直断裂方向张应力为主，震后随时间推移逐渐加强。尤其是鲜水河断裂，震后10a的张应力累积变化量已超过同震变化量。而该断裂上相继发生了2010年道孚5.4级、2011年炉霍5.7级、2014康定6.3级强震，随时间的推移震级呈增强特征，和汶川地震震后应力场变化具有显著相关性；（2）安宁河—则木河—小江断裂、马边断裂、红河断裂、金沙江断裂南段同震及震后影响均表现为与断层走向呈一定夹角的庄应力变化特征。安宁河—则木河—小江断裂主压应力方向（北北东）和背景构造主压应力方向（北北西向）相反，使断层表现为右旋，具有负影响，汶川震后强震活动较弱。红河断裂应力所反映的断层走滑活动特征和背景特征一致，表现为右旋走滑，但量值相对较小。金沙江断裂中南段以垂直断裂的压应力特征为主：（3）总体上，汶川地震同震及震后对鲜水河断裂影响最大，其次为金沙江断裂北段、安宁河—则木河断裂、马边断裂、鲁甸—昭通断裂、金沙江断裂中南段、红河断裂、丽江—小金河断裂，对小江断裂影响最小。从影响特征和影响强度2方面考虑，汶川震后川滇地区强震活动的主体区域和其具有较好的一致性。
为了进一步分析汶川地震对川滇地区主要块体边界活动断裂的影响，将汶川地震产生的应力变化投影到断裂滑动方向上，用库仑破裂应力变化来表示该地震对周围接收断裂产生的影响。具体计算公式参照Jaeger等（1969）的研究，在计算中有效摩擦系数参照Stein等（1992）和King等（1994）的研究，取值为0.4。有关接收断裂活动性质几何参数参考邓起东等（2002）结果给出。由于研究区域的地震震源深度绝大多数在0～25km范围内，因此我们研究汶川地震在10km深度上产生的库仑破裂应力变化，我们也计算了其他深度上产生的库仑破裂应力变化，发现模式和量值变化不大。
从库仑破裂应力变化结果（图8）来看：（1）同震产生的库仑破裂应力显著增加断裂主要包括龙门山断裂南端和北端、甘孜—玉树断裂、鲜水河断裂、白玉断裂、东昆仑断裂东段部分段落，这些部位库仑破裂应力增加值均超过0.01MPa，达到强震触发的阈值（King et al，1994；Harris，Simgson，1998）。震后随时间的变化，这些断裂库仑破裂应力大部分在同震的基础上有所增加，但增加量值不太显著。汶川地震之后，这些断裂相继有5.0级以上强地震发生，其中包括在甘孜一玉树断裂发生的2010年玉树7.1级、龙门山断裂南段发生的2013年芦山7.0级和鲜水河断裂发生的2014年康定6.3级地震。（2）同震库仑应力有所增加，但增加不是很显著的断裂，主要包括大凉山地块东边界马边断裂、大凉山地块南边界鲁甸一昭通断裂、金沙江断裂、丽江小金河断裂的南段、南汀河断裂，库仑应力增加值小于0.01MPa。鲁甸一昭通断裂库仑破裂应力震后随时间变化逐渐增强，震后10a达到0.01MPa，该断裂附近先后发生了几次5.0级以上强震，其中包括2014年8月3日鲁甸6.5级地震。丽江一小金河断裂震后随时间推移库仑破裂应力也逐渐增强，其南段震后10a库仑应力接近0.01MPa，在其附近发生了2012年宁蒗5.7级地震。金沙江断裂震后随时间变化库仑破裂应力有所增强，南汀河断裂震后随时间变化库仑应力有所增加，但并不显著。（3）同震库仑破裂应力减少的断裂包括龙门山断裂中段发震断裂、安宁河—则木河—小江断裂、丽江—小金河断裂北段、红河断裂及东昆仑断裂东段部分段落。除了丽江—小金河断裂北段在震后随时间推移，库仑破裂应力由同震的卸载转为震后的加载外，其他断裂震后随时间变化情况和同震效应一致，均存在不同程度的卸载增强的特征。总之，同震及震后断层库仑破裂应力反映的加卸载情况和强震活动空间分布具有一定的相关性。需要说明的是由于使用的模型及方法的不同，本文计算获得的同震库仑破裂应力变化和万永革等（2009）采用位错方法得到的结果存在明显的差异。
3 討论与结论
本文借助三维粘弹性有限元技术，以Wang等（2011）给出的大范围汶川地震同震水平位移场作为约束，模拟了汶川地震对川滇地区的同震及震后形变场、应力的时空演化特征，在模拟过程中，不但考虑了川滇地区地球介质的横向和纵向非均性，而且对同震上下盘断层面静态位移加载比例进行了调整，尽可能使我们的模拟结果在中远场和实际观测结果在形态上具有一致性。
模拟得到的川滇地区同震和震后形变场、应力状态反映的断层活动特征及库仑破裂应力变化反映的断层面应力加卸载情况，在空间上和强震活动具有较好的相关性。对于理解和解释汶川地震之后川滇地区出现的一些观测异常具有一定的参考意义。但本文模拟结果相比实际观测结果还存在许多不一致的地方，在模拟过程中仍然存在许多有待商榷和解决的问题。主要包括以下几个方面：（1）有限元边界约束条件对结果有着重要的影响，位错模型认为在离开断层无限远处和无限深处同震位移为0，即固定不动。在有限元模型中，由于模型范围总是有限的，我们无法做到这一点。本文采用大多数学者的做法，即模型边界法向约束，切向自由移动，这和实际情况存在明显的差异，也可能是本文模拟结果相比观测结果中远场同震形变相对较大的主要原因。因此探索更为有效、合理的边界条件约束是在模拟过程需要进一步解决的问题：（2）本文简单地通过采用上下盘断层面同震错动位移比例，来拟合测同震形变场，仅仅采用一个比例系数，过于简单，也可能是拟合效果存在差异的原因之一：（3）本文模型没有考虑重力、地壳背景构造应力场作用，重力、背景构造应力对同震及震后变形也存在一定的影响。另外，地壳介质本身具有复杂性（地形、断层的非连续性）等，这对模拟结果都会产生一定影响。在后续的数值模拟工作中需要进一步探索和研究。
参考文献：
毕研磊，张永志，曹海坤，等.2017.利用断层自剖分技术反演汶川MW7.9地震滑动分布[J].40（2）：211-215.
陈运泰，许力生，张勇，等.2008.2008年5月12日汶川特大地震震源特性分析报告[EB/OL].（2017-09-20）[2008-05].ht-tp：//www.csi.ac.cn/sichuan/.
邓起东，张培震，冉勇康，等.2002.中国活动构造基本特征[J].中国科学：地球科学，32（12）：1020-1030.
季灵运，刘传金，徐晶，等.2017.九寨沟MS7.0级地震的InSAR观测及发震构造分析[J].地球物理学报，60（10）：4069-4082.
季灵运，刘立炜，郝明.2015.利用InSAR技术研究滇西南镇康永德地区现今地壳形变特征[J].地震研究，38（1）：84-89.
焦佳爽，张永志，张凯南，等.2017.基于同震GPS数据的汶川大震同震应变特征研究[J].地震工程学报，39（3）：227-233.
李宁，赵强，李金.2017.基于InSAR形变观测反演2015年皮山MS6.5地震滑动分布[J].地震，37（2）：67-76.
李强，江在森，武艳强，等.2013.汶川地震震后短期青藏高原东北缘地壳变形特征成因讨论[J].地震，33（3）：124-132.
李玉江，陈连旺，杨树新.2013.基于应变能变化的芦山强震同震效应的数值模拟[J].震灾防御技术，8（4）：361-369.
邵志刚，傅容珊，薛霆虓，等.2008.昆仑山MS8.1级地震震后变形场数值模拟与成因机理探讨[J].地球物理学报，51（3）：805-816.
沈正康，万永革，甘卫军，等.2003.东昆仑活动断裂带大地震之间的黏弹性应力触发研究[J].地球物理学报，46（6）：786-795.
谈洪波，申重阳，李辉，等.2009.汶川大地震震后重力变化和形变的黏弹分层模拟[J].地震学报，31（5）：491-505.
万永革，沈正康，盛书中，等.2009.2008年汶川大地震对周围断层的影响[J]地震学报，31（2）：129-139.
万永革，沈正康，曾跃华，等.2007.青藏高原东北部的库仑应力积累演化对大地震发生的影响[J].地震学报，29（2）115-129.
主椿镛，杨文采，吴建平，等.2015.南北构造带岩石圈结构与地震的研究[J].地球物理学报，58（11）：3867-3901.
王辉，曹建玲，张怀，等，2007.川滇地区下地壳流动对上地壳运动变形影响的数值模拟[J].地震学报，29（6）：581-591.
吴建平，黄媛，张天中，等，2009.汶川MS8.0级地震余震分布及周边区域P波三维速度结构研究[J].地球物理学报，52（2）：320-328.
吴建平，明跃红，王椿镛.2006.川滇地区速度结构的区域地震波形反演研究[J].地球物理学报，49（5）：1369-1376.
张培震，邓起东，张国民，等.2003.中国大陆的强震活动与活动地块[J].中国科学：地球科学，33（增刊1）.12-20.
Fielding E J，Sladen A，Li Z，et aL 2013.Kipemalic fault slip evolutionsource models of he 2008 M7.9 Wenchuan.earthquake in Chinafrom SAIL interferometry，GPS and teleseismic analysis and implies-Lionsfor Lougmen Shan tectonics[J].Geophysical journal interns-tional，194（2）：1138-1166.
Funning G J，Parsons 13，Wright T J.2007.Fault slip in the 1997Manyi，Tibet earthquake from linear elastic modelling of InSAR displace-ments[J].Geophysical Journal International，169（3）：988-1008.
Ham's R A，Simpson R W.1998.Supp cession of large earquakes by stressshadow：A Conpaiison of Coulomb and rate-and-state failure[J].Journal of Geophysical Research solid Earth，103（B10）：24439-24451.
Jaeger J C，Cook N G W，Zimmerman R.1969.Fundamentals of rock me-chanics[M].Chapman and Hall'Landon.
Ji C，Hayes G.2008.NEIC，Preliminary result of theMay 12，2008MW7.9 eastern Sichuan，China earthquake[EB/OL].（2017-09-20）[2008-05].http；//earthquake.nags.gov/eqcenter/eqinthe-news/us2008ryan/.
Jonssons，Segall P，Pedersen R，et al.2003.Post-earthquake groundmovements correlated to pore-pressure transients[J].Nature，424（4945）：179-183：
Jouanne F，Awan A，Madji A，et al.2011.Postseismic deformation in Pa-kistan after the 8 October 2005 earthquake：evidence of afterslip a-long a flat north of the Balakot-Bagh thrust[J].Journal of Geo-physical Research：Solid Earth，116（B07401）：1-22.
King G C P，Stein Rs，Lin J.1994.Static stress changes and the trigge-ring of earthquakes[J].Bulletin of the Seismological Society of A-merica，84（3）：935-953.
Rice J R，Go J.1983.Earthquake after effects and triggered seismic phe-nomena[J].Pure and Applied Geophysics，121（2）：187-219.
Shen Z K，Sun J，Zhang P，et al.2009.Slip maxima at fault junctions andrupturing of barriers dining the 2008 Wenchuan earthquake[J].Na-ture geosoence，2009，2（10）：718-724.
Stein R S，King G C P，Lin J.1992.Change in failure stress on the south-ern San Andreas fault system caused by the 1992 magnitude=7.4Lenders earthquake[J].Science，258（5086）：1328-1332.
Wan Y，Shen Z K，Burgmann R，et at.2017.Fault geometry and slip dis-tribution of the 2008 MW7.9 Wenchuan，China earthquake，inferredfrom GPS and InSAR measurements[J].Geophysical Journal tnter-nationel，208（2）：748-766.
Wang，Xu C，Ge L.2007.Coseismic deformation and slip distribution ofthe 1997 MW7.5 Manyi，Tibet，earthquake from InSAR measure-ments[J].Journal of Geodynamics，44（3）：200-212.
Wang Q，Qiao X J，Len Q G，et al.2011.Rupture of deep faults in the2008 Wenchuan earthquake and uplift of the Longmen Shan[J].Na-Lure Ceascience，4（9）：634-640.
