2017, Vol. 60
2016年3月2日印尼苏门答腊西南海域发生了MS7.8地震，该地震震中坐标为(94.236°E，4.905°S)，震源深度为20 km(http://www.cea-igp.ac.cn/[2016-01-25]).根据该地区周边历史地震资料，USGS(美国地质调查局)向印度洋各国发布海啸预警.虽然此次走滑地震未造成海啸，但仍有多次强烈的M5级左右余震.
苏门答腊地区位于印度—澳大利亚板块向欧亚板块斜俯冲的边界，是法向和切向滑动分配(slip partitioning)的经典区域(见图 1a)，被称为地质财富宝库，是研究重要地质过程的天然实验室(McCaffrey，2009).印度—澳大利亚板块以约6 cm·a-1向欧亚板块斜向俯冲(Khan et al.，2005)，在缅甸微板块(Burma microplate)两侧解耦(Pubellier et al.，2004)为不同的断层运动形式: 西侧以俯冲为主(表现为沿巽他海沟发育的低角度逆冲断层)和东侧发生大规模右旋走滑(表现为苏门答腊走滑断层)，图 1b展示了这种解耦活动.古地磁研究表明印度—澳大利亚板块向北运动分量存在差异，在东印度洋海岭(也称东径90°海岭)东侧呈现出东高西低的现象(Socquet et al.，2006)，导致巽他海沟西侧探测者破裂区(Investigator Fracture Zone)发育多条N向和NNE向左旋走滑断裂(Royer and Gordon，1997)且历史上有多次大地震发生，见图 1a.例如，2012年4月11日发生的MW8.6和MW8.2地震及2012年1月10日发生MW7.2地震.
|
图 1 苏门答腊地区构造背景及历史(MW≥5)地震分布.红色震源球为该区域1976年至今历史地震(http://www.globalcmt.org/[2016-01-25])，黑色震源球为本次地震; 印度—澳大利亚板块相对巽他板块的运动速率瞿武林等(2016)和Subarya等(2006)，探测者破裂区断层数据引自Lin等(2014).(b)巽他海沟板块斜向俯冲的断裂解耦活动(McCaffrey，2009) Fig. 1 Tectonic background of Sumatra and distribution of historical earthquakes.Red-white focal spheres present historical earthquakes since 1976(http://www.globalcmt.org/[2016-01-25])and black-white focal sphere is the MS7.8 Sumatra earthquake.The movement velocity of the India-Australia Plate relative to the Sunda Plate is from the studies of Qu et al.(2016)and Subarya et al.(2006).The fault data is from Lin et al.(2014).(b)The decoupling of plates oblique subduction in the Sunda Trench(McCaffrey，2009) |
大地震的发生会引起区域位移场和应力场发生变化，进而改变区域内及临近断层的应力状态，可能增加断层滑动危险性或者使得断层趋于更加稳定，从而影响区域及邻区的地震活动性.因此，本文数值计算此次苏门答腊MS7.8大地震引起的同震变形和应力场及库仑应力变化，分析了库仑应力变化与余震分布的关系，讨论此次地震对周围断层的影响，为该区域的地震危险性评估提供一些参考.基于PREM分层模型，建立了全球等效体力分层球型模型(张贝等，2015); 地表给予自由边界条件，核幔边界施加弹簧边界，合理地回避了一般对有限区域进行计算时侧边界和地面边界的不确定性.通过采用断层附近自适应网格加密方法，有效地计算了2016年3月2日苏门答腊西南海域MS7.8地震引起的同震位移和应力变化，并计算周围区域的库仑应力变化，进而探讨了此次地震对周围区域的影响.
2 计算方法在应用弹性回跳理论解释地震发生的背景下，Steketee(1958)首先引入弹性位错理论来处理地震问题且迅速得到了应用.之后，众多学者发展了不同的位错模型，对大地震同震变形等问题进行了研究.Okada(1985，1992)在前人工作的基础上，给出了均匀弹性半无限空间中三种位错源引起的地表与地球内部的同震位移场及应力场的解析解，这种方法已成为计算弹性半无限空间同震效应的经典方法.Wang(2005a，2005b)考虑构造变形和物质迁移产生的重力变化，将Okada的方法拓展到分层弹性半无限空间.地球曲率对计算地震同震效应的影响不可忽略，Sun和Okubo(2002)建立了球形地球模型，将其与半无限空间介质进行对比，发现曲率与层状构造的影响不容忽视.特别指出的是，在从位移计算应力应变必须在球坐标下进行，而采用直角坐标常常会引起较大的误差(石耀霖和朱守彪，2006)，例如，2015年尼泊尔MS8.1大地震，处于北纬28.147°N，南北向位移达3 m多时，直角坐标计算误差可达数kPa，Cheng等(2016)通过数值模型和Okada模型对比分析得出采用直角坐标计算同震应变和应力时在30°和60°纬度最大相对误差分别可达到30%和60%.
Burridge和Knopoff(1964)提出弹性介质中的位错源与体力源完全等价.张贝等(2015)利用该理论将位错等效为体力，添加到平衡方程的右端项中，将地震位错问题转化为传统有限元问题，避免了常规有限元法处理不连续面时遇到的困难.本文基于该计算方法，在划分网格时使用自适应加密技术，采用PREM分层模型; 全球计算网格共3349696个节点和2661307个单元，发震断层附近网格密度约2.0 km(图 2a).计算中断层滑动模型(图 2b)采用中国地震局陈运泰院士研究小组提供的断层滑动模型(http://www.cea-igp.ac.cn/tpxw/273697.shtml[2016-01-25])，走向185°，倾角84°.
|
图 2 网格剖分和地震滑动断层模型 Fig. 2 Computational model with 3349696 grids and 2661307 elements，and fault slip model of the MS7.8 earthquake in this numerical computation |
图 3a给出了2016年3月2日苏门答腊MS7.8地震引起的近场同震位移，位移变化整体上呈现对称分布，1 m以上的水平滑移量集中在发震断层两侧.计算得到的最大水平位移量大约6.74 m，与断层滑动模型给出的相一致，位移量随着离开断层距离的增加而迅速减小.根据地表位移矢量图(图 3b)可知，断层上盘向南运动，下盘向北运动，呈现左旋走滑性质，位移变化量10 cm的范围约为154.25 km.图 3(c-e)给出了同震地表南北向和东西向及垂直方向位移量，地震造成南北方向水平位移较大，断层下盘向北运动最大位移约为6.14 m; 断层上盘向南最大位移为6.72 m，这与断层近南北向左旋走滑性质相关.在东西方向上，水平位移量相对较小，最大位移量约0.97 m，断层下盘由南向北先西向运动再东向运动，断层上盘则由北向南先东向运动再西向运动.在垂直方向上，在震中及周边海域，垂直位移较小，最大位移量约0.25 m，断层下盘由南向北先向下运动再向上运动，断层上盘则相反，因此发生海啸的可能性较低.
|
图 3 地表同震水平位移和位移矢量及各个方向同震位移分量 其中黑色震源球为主震，红色震源球为余震.(a)水平位移量;(b)水平位移矢量;(c)南北向位移分量;(d)东西向位移分量;(e)垂直位移分量. Fig. 3 Co-seismic surface horizontal displacements，vectors and co-seismic displacement components Black-white focal sphere represents the main shock and red-white focal spheres represent aftershocks.(a)The co-seismic horizontal displacements;(b)The co-seismic displacement vectors;(c)The co-seismic displacements in north-south direction;(d)The co-seismic displacements in east-west direction;(e)The co-seismic vertical displacements. |
根据上述可知，2016年3月2日苏门答腊地震是一个以左旋走滑为主的地震，与苏门答腊海域发育多条左旋走滑性质断层的构造背景相一致; 垂向同震位移较小，触发海啸几率较小.
3.2 同震应力场2016年3月2日苏门答腊MS7.8地震属于走滑型地震，其垂直方向应力变化很小，在此我们重点分析此次地震引起的12 km深度(CMT余震震源深度)同震水平方向近场应力，见图 4，红色表示应力值为正，蓝色表示应力值为负(引用弹性力学定义——张应力为正，应力增加即表示压应力减小).可以看出应力变化区域主要集中在地震断层两侧，同震应力随着离开断层距离的增加而迅速减小.地震造成断层两侧同震应力变化量在1.0 MPa以上，而在苏门答腊岛变化量级在10.0 kPa以下.东西向最大应力变化约-6.59 MPa，南北向最大应力变化约10.15 MPa，最大水平剪应力变化约31.52 MPa，南北向和东西向水平应力表现为正负相间，水平剪应力在断层面两侧主要表现为增加.
|
图 4 12 km深度同震水平应力变化 其中黑色震源球是本次地震，红色震源球为余震.(a)南北向应力分量;(b)东西向应力分量;(c)水平剪应力分量. Fig. 4 The co-seismic horizontal stress changes at 12 km depth Black-white focal sphere represents the MS7.8 earthquake and red-white focal spheres represent aftershocks.(a)The co-seismic horizontal stress changes in north-south direction;(b)The co-seismic horizontal stress changes in east-west direction;(c)The co-seismic horizontal shear stress changes. |
20世纪80年代，研究者发现大震后余震的发生与地震造成的剪应力变化或库仑应力变化有关(Stein and Lisowski，1983).众多学者据此利用库仑应力变化评估区域地震危险性，例如，Stein等(1997)研究1930—1992年土耳其北安纳托利亚断层10次M≥6.7地震后，指出1999年伊兹米特地区发生大地震的发震概率为12%.Nalbant等(1998)在研究土耳其西北部和爱琴海地区29次地震之间的应力触发问题后，指出伊兹米特未来可能发生地震，在1999年土耳其伊兹米特地区发生了M7.4大地震.万永革等(2007)在黏弹性介质中，考虑1920年海原地震以来大地震(M≥7.0)产生的应力变化和根据GPS速度场得出的长期构造加载共同作用的结果，发现20次M≥7.0地震中有17次大地震发生在库仑破裂应力变化为正的区域，触发率达85%.更多的研究表明历史地震序列与库仑应力变化有着相关性，且地震发生在库仑应力增加区域的可能性大，现在库仑应力变化已被广泛用于研究断层的危险性(单斌等，2012; Xiao and He，2015).
库仑应力变化(Coulomb Failure Stress changes(ΔCFS))表达式为
|
(1) |
其中Δτs为地震引起的剪应力变化(与滑动方向一致为正)，Δσn为正应力变化(以张应力为正)，ΔP为断裂带孔隙压力变化(已压缩为正)，μ为摩擦系数(范围为0.0~1.0)，在考虑孔隙压力时，常用有效摩擦系数μ′代替μ(Ma et al.，2005; Cattin et al.，2009)，即
|
(2) |
当ΔCFS为正时，有利于断层错动，使地震危险性增加; ΔCFS为负时，地震发生的危险性相应减小.
有效摩擦系数一般取值为0.0~0.75，在大多数计算中取值为0.4(Stein et al.，1992; King et al.，1994).有效摩擦系数受断层孔隙压力、累积滑动量和产状的影响.Parsons等(1999)认为高滑移积累量的断层有效摩擦系数较小(0.2~0.4)，而低滑移积累量的断层有效摩擦系数较高(0.6~0.8).Ali等(2008)和Parsons等(1999)认为对于走滑断层，剪应力变化对地震活动性具有重要影响，采用较小有效摩擦系数; 而对于逆断层，正应力变化的影响较大，需采用较大的有效摩擦系数.据此，本文中探测者破裂区(IFZ)左旋走滑断层有效摩擦系数取0.2.
单斌等(2012)认为左旋走滑地震ΔCFS随着深度的增加变化并不明显，在断层破裂面正下方及周围区域为ΔCFS增加区域.苏门答腊地震破裂过程以左旋走滑运动为主，且断层破裂面深度为45 km，因此本文利用上述公式(2)选取12 km深度计算苏门答腊地震及其周边区域ΔCFS，分析余震分布与ΔCFS关系.接受断层参数根据CMT提供的苏门答腊地震及余震震源机制解(表 1)选取，即走向6°，倾角84°，滑动角3°.
|
表 1 2016年3月2日苏门答腊地震及余震震源机制解(引自CMT) Table 1 The mechanism of this MS7.8 Sumatra earthquake on March 2，2016，and two big aftershocks(from CMT) |
图 5a给出了12 km深度库仑应力变化，ΔCFS增加量最大约0.38 MPa，减小量最大约32.26 MPa，库仑应力变化幅度较大处均集中在震中附近.其中震中断层两侧ΔCFS呈现负值为主且达MPa量级，ΔCFS增加呈4个扇区，其中两个扇区位于主震破裂面的两侧，另外两个扇区位于主震破裂面的两端，ΔCFS随距离增加迅速衰减，周边区域ΔCFS正负相间，在苏门答腊断层(Sumatra fault)上ΔCFS有正有负，但变化量均较小.
|
图 5 12 km深度库仑应力变化与余震分布 黑色震源球代表主震，红色震源球代表余震震源机制解，红色点表示余震，圆点大小表示震级.(a)震中及周边ΔCFS与USGS余震分布;(b)近场ΔCFS与CMT余震震源机制解;(c)近场ΔCFS与USGS余震分布. Fig. 5 The distribution of ΔCFSs at 12 km depth The black-white focal sphere represents earthquake mechanism of main shock and the red-white focal spheres represent mechanism of aftershocks. The red spots represent aftershock magnitude with different sizes.(a)The distribution of ΔCFSs of the epicenter and its surrounding region，as well as the distributions of aftershocks from USGS;(b)The distribution of ΔCFSs in near-field and the distribution of aftershocks from CMT;(c)The distribution of ΔCFSs in near-field and the distribution of aftershocks from USGS. |
库仑应力的增加使断层上受力增加，当ΔCFS超过应力触发值0.01 MPa时(Harris and Simpson，1998)，触发地震的可能性增大; 另外，也有一些地震在小于0.01 MPa区域触发(Ziv and Rubin，2000).现有研究中阀值0.01 MPa仅是地震应力触发的一种评估标准，地震触发主要还是与ΔCFS正负相关.库仑应力减弱区域(应力影区)，即断层上负荷部分卸载，地震危险性降低.图 5(b-c)和表 2给出了12 km库仑应力变化与余震震中的对应关系，结果显示仅约23.53%地震位于ΔCFS为正的区域，远场地区余震均分布在库仑应力增加区域且小于0. 01 MPa，但近场大部分余震处于库仑应力变化小于零的区域，仅3月4日发生的余震处于库仑应力增加区域.这似乎与库仑应力触发机制相违背，导致这种分布特征的原因值得分析.
|
|
表 2 余震震源处ΔCFS(余震数据来源于USGS(M≥4)) Table 2 The distribution of ΔCFSs(locations and magnitudes of aftershocks from USGS(M≥4)) |
在很多研究中，ΔCFS随摩擦系数增加而增加，即断层危险性随摩擦系数的增加而增大，这与常识相违背，朱守彪和缪淼(2016)提出其中的原因是在计算中未考虑构造应力的作用，即忽略了(3)式中的第一项，而仅考虑最后一项的影响.
|
(3) |
其中μ′0为震前有效摩擦系数，μ′1为震后有效摩擦系数.在计算ΔCFS时假设有效摩擦系数为常数，但μ′并非常数.根据速率-状态相关的断层摩擦本构关系(如Dieterich-Ruina定律)，速度阶跃会对摩擦系数产生直接效应和演化效应，使摩擦系数产生瞬时增加，然后在特征距离内下降直到达到稳态，导致速度强化或速度弱化.另外，Wang和Chen(2001)认为主震应力阶跃引起的孔隙流体压力变动导致μ′随时间变化，同震破裂作用造成断层的渗透性增高，但在海底受水头压力和孔隙压力扩散效应影响，孔隙压力沿断裂、裂隙等结构向四周扩散，可能导致余震震中孔隙压力不断增大.同震流体的热压作用、水-岩相互作用等物理化学反应也会造成矿物组成和化学成分发生变化，使断层弱化引起内摩察系数变化，Chen和Talwani(2001)认为断层摩擦系数的变化影响很大，值得讨论.上述因素均有可能造成有效摩擦系数不断变化.Johnston等(2011)根据Ghosh和Lithgow-Bertell and Guynn应力模型得到IFZ断层法向压应力均达10 MPa量级，若μ′减小0.05，则ΔCFS0会相应增加MPa数量级，使原来利用式(2)计算得到的ΔCFS小于零的部分区域变为正值，这或许能够解释主震断层周边大部分余震分布在ΔCFS为负的区域，但这种解释仅能适用于主震近场区域，因为上述影响有效摩擦系数的因素的影响范围也仅限主震断层周边区域.
计算ΔCFS时，一般不考虑区域构造应力场，仅计算地震引起的断层面上正应力和剪应力变化，但区域构造应力场的大小及方向也会影响ΔCFS(周宇明等，2008).投影断层的走向、倾角和滑动角也是计算库仑应力变化中重要参数，评估苏门答腊地区地震活动性需详细了解该地区断层构造.以下讨论在考虑区域构造应力场情况下，ΔCFS与余震的分布关系，以及苏门答腊周边断层ΔCFS.
4 讨论 4.1 假定区域应力已知时，靠近断层区域的库仑应力变化改进计算方法地震的发生对远场构造应力场影响较小，但在近场地震能改变区域构造应力方向.石耀霖和曹建玲(2010)考虑区域构造应力，改进了库仑应力变化计算方法(式(4)—(6)).计算中正应力变化仅与地震引起的应力变化有关，剪应力变化则考虑背景构造应力场，不限定断层面错动方向，计算断层面剪应力绝对值变化.
|
(4) |
|
(5) |
|
(6) |
探测者破裂区左旋走滑活动断层表明中间主应力轴σ2垂直，σ1和σ3水平，根据世界应力图(WSM)，最大水平应力由西北-东南向变为南北向，其方向与沃顿海盆运动方向近正交并平行于海沟(Johnston et al.，2011).Müller等(2015)认为沃顿海盆中最大水平应力近东-南向，大小约80 MPa; Johnston等(2011)根据Ghosh和Lithgow-Bertell and Guynn应力模型得到东印度洋最大压应力分别达到93 MPa和88 MPa，方向为近东-南向，南北向断层法向压应力分别为27 MPa和80 MPa; Sandiford等(2005)认为印度—澳大利亚板块中走滑断层的孕育表明最小水平应力与垂直应力近似相等，且最大水平应力约为100 MPa.
为探讨构造应力场影响，本文假定σ2垂直，根据5 km水深设为约-50 MPa; σ1与σ2近似相等，约-50 MPa，方向为东-北; σ3约-90 MPa，方向为东-南，接受断层参数不变，利用上述改进方法计算了12 km深度近场区域库仑应力变化与近场余震分布(图 6).与传统方法相比，改进方法改变了近场库仑应力变化的大小和分布，ΔCFS增加量最大约0.37 MPa，减小量最大约-22.09 MPa，其中2016年3月2日16时08分余震ΔCFS由-22.0158 MPa增加到-6.4754 MPa，变化突出，但其余余震变化并不明显，仅3月4日发生的余震分布在库仑应力增加区域.但若能详细了解区域构造应力场，改进方法或许能更好地反映近场余震与库仑应力变化的对应关系.
|
图 6 12 km深度库仑应力变化与余震分布 黑色震源球代表主震，红色震源球代表余震震源机制解，红色点表示余震，圆点大小表示震级. (a)近场库仑应力变化与CMT余震震源机制解.(b)近场库仑应力变化与USGS余震分布. Fig. 6 The distribution of ΔCFSs at 12 km depth The black-white focal sphere represents earthquake mechanism of MS7.8 main shock and the red-white focal spheres represent earthquake mechanism of aftershocks. The red spots represent aftershock magnitude with different sizes.(a)The distribution of ΔCFSs in near-field and the distribution of aftershocks from CMT;(b)The distribution of ΔCFSs in near-field and the distribution of aftershocks from USGS. |
在断层不发育或不清楚情况下，考虑区域构造应力场，库仑应力变化的计算也可以进行最优破裂方向投影.所谓最优破裂方向是指在计算库仑应力变化时，最优破裂方向上的库仑应力变化大于同一地点其他产状断层面上的库仑应力变化.King等(1994)提出了考虑构造应力场的最优破裂方向选取方式，具体方式如下，根据库仑破裂准则:
|
(8) |
其中σβ为断层面正应力，τs为断层面上剪应力.假设断层面与最大主应力σ1轴夹角为β，则
|
(9) |
|
(10) |
则方程(8)可变为
|
(11) |
方程(11)中，当tan2β=-1/μ时，σf取得最大值.在计算最优断层取向时，将已计算出的地震引起的应力σije与区域构造应力σijr叠加得到总应力σijt:
|
(12) |
在二维中，根据确定的最大主应力轴与轴的夹角为
|
(13) |
因此，最优破裂方向为θ±β.
利用4.1中假定的构造应力场，计算各点最优破裂方向上库仑应力变化.结果表明最优破裂方向上库仑应力变化显著(图 7)，震源附近原来ΔCFS减小区域大部分变为了增加区域，ΔCFS增加量最大约8.39 MPa，减小量最大约-32.74 MPa; 约76.47%的余震分布在库仑应力变化增加区域，而近场区域达到92.31%且大于0.01 MPa，与根据震源机制解求得的ΔCFS相比有很大提高.因此，断层参数对库仑应力变化的影响较大，在不清楚断层参数的情况下，考虑构造应力场利用最优破裂方向计算库仑应力变化能够较好地反映余震分布，验证了库仑应力变化触发作用.但由于余震震源机制解资料有限，无法判断余震破裂方向是否与最优破裂方向一致，且库仑应力变化对构造应力场较敏感，因此需要更详细的震源机制解和构造应力场资料才能获得更精确的结果.
|
图 7 12 km深度最优取向库仑应力变化与余震分布 黑色震源球代表主震，红色震源球代表余震震源机制解，蓝色点表示余震，圆点大小表示震级.(a)震中及周边ΔCFS与USGS余震分布.(b)近场ΔCFS与CMT余震震源机制解.(c)近场ΔCFS与USGS余震分布. Fig. 7 The distribution of ΔCFSs at 12 km depth using optimal orientation methods The black-white focal sphere represents the mechanism of main MS7.8 earthquake and the red-white focal spheres represent mechanism of aftershocks. The blue spots represent aftershock magnitude with different size.(a)The distribution of ΔCFSs of the epicenter and its surrounding region，as well as the distributions of aftershocks from USGS;(b)The distribution of ΔCFSs in near-field and the distribution of aftershocks from CMT;(c)The distribution of ΔCFSs in near-field and the distribution of aftershocks from USGS. |
库仑应力变化受断层几何参数和滑动角的影响，对苏门答腊地区断层做详细的调查研究有助于更准确地了解苏门答腊地区地震活动性.Cattin等(2009)给出苏门答腊及周边地区断层走向、倾角和滑动角数据，断层以高角度走滑为主.利用Cattin给出的断层数据，计算了断层上12 km深度库仑应力变化(图 8).
|
图 8 苏门答腊周边断层12 km深度库仑应力变化 Fig. 8 The distribution of ΔCFSs of Sumatra surrounding faults at 12 km depth |
计算结果表明，在苏门答腊周边断层上ΔCFS约在-594~241 Pa，量级较小，库仑应力增加区域主要集中在苏里曼断层(Seuliman fault)和苏门答腊断层系统(Sumatra fault system).根据0.01 MPa应力触发值，苏门答腊地震对周边断层的地震触发影响可能较小.但在长时间尺度上，研究苏门答腊地震对苏门答腊—安德曼断层地震活动性的影响还需要考虑脆性上地壳与黏弹性下地壳和上地幔耦合作用形成的黏弹性松弛效应.
4.4 介质结构复杂性与断层结构复杂性的影响相关构造物理实验和数值模拟结果表明，余震受控于断层面上的构造和介质的不均匀性(曲均浩和蒋海昆，2012).介质的不均匀性会产生局部应力调整，特别是震源区存在主震断层带、裂隙等，介质不均匀性较大，从而导致应力不均匀性，影响库仑应力变化的计算结果.另外，在库仑应力变化计算中，断层的几何参数会影响断层面上正应力和剪应力，因此，断层几何结构也是其影响因素，正如上述中根据主震区域周边断层的几何参数计算得到的库仑应力变化存在区别一样.
5 结论印度—澳大利亚板块整体向北运动，由于俯冲速度存在差异，在巽他海沟西侧洋壳发育多条北向和北北东向左旋走滑断裂，导致苏门答腊地区地震活动强烈且历史上发生多次大地震.本文建立全球等效体力法同震位错计算模型和地震断层滑动模型，计算苏门答腊地震同震位移、应力以及周边区域库仑应力变化，得出以下结论.
(1) 苏门答腊地震造成最大水平位移量约6.74 m，其中南北方向水平位移最大，东西方向水平位移量相对较小，垂直位移最小; 位移矢量表现为断层下盘向北运动，断层上盘向南运动.进一步说明苏门答腊地震为典型的左旋走滑为主的地震，引发海啸的可能性非常低.
(2) 有效摩擦系数变化对断层附近库仑应力变化影响可能较大.有效摩擦系数不是常数，滑动速率变化、孔隙流体压力和同震流体的物理化学作用导致它随时间变化.IFZ断层法向压应力达到10 MPa量级，若μ′减小0.05，则ΔCFS0会相应增加MPa数量级，使原来计算得到的库仑应力变化小于零的部分区域变为正值.
(3) 区域构造应力场对库仑应力变化影响较大.考虑区域构造应力场，改进方法会改变库仑应力变化的大小和分布，但总体并不明显.利用最优破裂方向计算得到的库仑应力变化与余震分布关系能很好地吻合库仑应力变化触发机制，验证了库仑应力触发作用.但最优破裂方向对区域构造应力场较敏感并且苏门答腊地区构造应力场数据有限，对库仑应力变化计算结果的精确性产生了一定的影响，另外余震震源机制解资料有限，较难判断余震破裂是否与最优方向一致.
(4) 苏门答腊地区周边断层库仑应力变化较小，对地震活动性影响可能较小.基于苏门答腊地区详细断层参数，计算得到的库仑应力变化在百帕以内，量级较小，对苏门答腊地区大地震触发作用有限.但在长时间尺度上，研究苏门答腊地震对苏门答腊周边断层地震活动性的影响需要考虑黏弹性松弛效应.
本文所用区域构造应力场资料有限，需进一步详细研究苏门答腊地区应力场，此外本文中采用的是线弹性本构关系，未考虑黏弹性效应，后续研究将弥补这一缺点.
致谢感谢中国地震局陈运泰院士研究小组提供的断层滑动模型和Cattin等(2009)提供的苏门答腊地区及其周边断层数据.感谢本文编辑和审稿专家针对本文提出的建设性意见和建议.
| Ali S T, Freed A M, Calais E, et al. 2008. Coulomb stress evolution in Northeastern Caribbean over the past 250 years due to coseismic, postseismic and interseismic deformation. Geophysical Journal International, 174(3): 904-918. DOI:10.1111/gji.2008.174.issue-3 | |
| Burridge R, Knopoff L. 1964. Body force equivalents for seismic dislocations. Bulletin of the Seismological Society of America, 54(6A): 1875-1888. | |
| Cattin R, Chamot-Rooke N, Pubellier M, et al. 2009. Stress change and effective friction coefficient along the Sumatra-Andaman-Sagaing fault system after the 26 December 2004 (Mw=9.2) and the 28 March 2005 (Mw=8.7) earthquakes. Geochemistry, Geophysics, Geosystems, 10(3): Q03011. | |
| Chen L, Talwani P. 2001. Mechanism of initial seismicity following impoundment of the Monticello Reservoir, South Carolina. Bulletin of the Seismological Society of America, 91(6): 1582-1594. DOI:10.1785/0120000293 | |
| Cheng H, Zhang B, Zhang H, et al. 2016. Necessity of using heterogeneous ellipsoidal Earth model with terrain to calculate co-seismic effect.//EGU General Assembly Conference. EGU General Assembly Conference Abstracts, 2512. | |
| Harris R A, Simpson R W. 1998. Suppression of large earthquakes by stress shadows:A comparison of Coulomb and rate-and-state failure. Journal of Geophysical Research:Solid Earth, 103(B10): 24439-24451. DOI:10.1029/98JB00793 | |
| Johnston M D, Long M D, Silver P G. 2011. State of stress and age offsets at oceanic fracture zones and implications for the initiation of subduction. Tectonophysics, 512(1-4): 47-59. DOI:10.1016/j.tecto.2011.09.017 | |
| Khan S A, Gudmundssonó. 2005. GPS analyses of the Sumatra-Andaman earthquake. Eos Transactions American Geophysical Union, 86(9): 89-94. | |
| King G C P, Stein R S, Lin J. 1994. Static stress changes and the triggering of earthquakes. Bulletin of the Seismological Society of America, 84(3): 935-953. | |
| Lin J Y, Sibuet J C, Hsu S K, et al. 2014. Could a Sumatra-like megathrust earthquake occur in the south Ryukyu subduction zone?. Earth, Planets and Space, 66: 49. DOI:10.1186/1880-5981-66-49 | |
| Ma K F, Chan C H, Stein R S. 2005. Response of seismicity to Coulomb stress triggers and shadows of the 1999 Mw=7.6 Chi-Chi, Taiwan, earthquake. Journal of Geophysical Research:Solid Earth, 110(B5): B05S19. | |
| McCaffrey R. 2009. The tectonic framework of the sumatran subduction zone. Annual Review of Earth & Planetary Sciences, 37: 345-366. | |
| Müller R D, Yatheesh V, Shuhail M. 2015. The tectonic stress field evolution of India since the Oligocene. Gondwana Research, 28(2): 612-624. DOI:10.1016/j.gr.2014.05.008 | |
| Nalbant S S, Hubert A, King G C P. 1998. Stress coupling between earthquakes in northwest Turkey and the north Aegean Sea. Journal of Geophysical Research Solid Earth, 103(B10): 24469-24486. DOI:10.1029/98JB01491 | |
| 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. | |
| Okada Y. 1992. Internal deformation due to shear and tensile faults in a half-space. Bulletin of the Seismological Society of America, 82(2): 1018-1040. | |
| Parsons T, Stein R S, Simpson R W, et al. 1999. Stress sensitivity of fault seismicity:A comparison between limited-offset oblique and major strike-slip faults. Journal of Geophysical Research:Solid Earth, 104(B9): 20183-20202. DOI:10.1029/1999JB900056 | |
| Pubellier M, Monnier C, Maury R, et al. 2004. Plate kinematics, origin and tectonic emplacement of supra-subduction ophiolites in SE Asia. Tectonophysics, 392(1-4): 9-36. DOI:10.1016/j.tecto.2004.04.028 | |
| Qu J H, Jiang H K. 2013. Review of aftershock mechanism research. Earthquake Research in China (in Chinese), 28(2): 109-120. | |
| Qu W L, Zhang B, Huang L Y, et al. 2016. Comparisons of global coseismic displacements from several fault slip models for the 2004 Sumatra earthquake. Chinese Journal of Geophysics (in Chinese), 59(8): 2843-2858. DOI:10.6038/cjg20160811 | |
| Royer J Y, Gordon R G. 1997. The motion and boundary between the Capricorn and Australian plates. Science, 277(5330): 1268-1274. DOI:10.1126/science.277.5330.1268 | |
| Sandiford M, Coblentz D, Schellart W P. 2005. Evaluating slab-plate coupling in the Indo-Australian plate. Geology, 33(2): 113-116. DOI:10.1130/G20898.1 | |
| Shan B, Li J H, Han L B, et al. 2012. Coseismic Coulomb stress change caused by 2010 MS=7. 1 Yushu earthquake and its influence to 2011 MS=5.2 Nangqên earthquake. Chinese Journal of Geophysics (in Chinese), 55(9): 3028-3042. DOI:10.6038/j.issn.0001-5733.2012.09.021 | |
| Shi Y L, Zhu S B. 2006. Discussion on method of calculating strain with GPS displacement data. Journal of Geodesy and Geodynamics (in Chinese), 26(1): 1-8. | |
| Shi Y L, Cao J L. 2010. Some aspects in static stress change calculation——case study on Wenchuan earthquake. Chinese Journal of Geophysics (in Chinese), 53(1): 102-110. DOI:10.3969/j.issn.0001-5733.2010.01.011 | |
| Socquet A, Vigny C, Chamot-Rooke N, et al. 2006. India and Sunda plates motion and deformation along their boundary in Myanmar determined by GPS. Journal of Geophysical Research:Solid Earth, 111(B5): B05406. | |
| Stein R S, Lisowski M. 1983. The 1979 Homestead Valley Earthquake Sequence, California:Control of aftershocks and postseismic deformation. Journal of Geophysical Research:Solid Earth, 88(B8): 6477-6490. DOI:10.1029/JB088iB08p06477 | |
| Stein R S, King G C P, Lin J. 1992. Change in failure stress on the Southern San Andreas fault system caused by the 1992 magnitude=7.4 Landers earthquake. Science, 258(5086): 1328-1332. | |
| Stein R S, Barka A A, Dieterich J H. 1997. Progressive failure on the North Anatolian fault since 1939 by earthquake stress triggering. Geophysical Journal International, 128(3): 594-604. DOI:10.1111/gji.1997.128.issue-3 | |
| Steketee J A. 1958. On Volterra's dislocations in a semi-infinite elastic medium. Canadian Journal of Physics, 36(2): 192-205. DOI:10.1139/p58-024 | |
| Subarya C, Chlieh M, Prawirodirdjo L, et al. 2006. Plate-boundary deformation associated with the great Sumatra-Andaman earthquake. Nature, 440(7080): 46-51. DOI:10.1038/nature04522 | |
| Sun W K, Okubo S. 2002. Effects of earth's spherical curvature and radial heterogeneity in dislocation studies-for a point dislocation. Geophysical Research Letters, 29(12): 46-1. | |
| Wan Y G, Shen Z K, Zeng Y H, et al. 2007. Evolution of cumulative Coulomb failure stress in northeastern Qinghai-Xizang (Tibetan) Plateau and its effect on large earthquake occurrence. Acta Seismologica Sinica (in Chinese), 29(2): 115-129. | |
| Wang R J. 2005a. The dislocation theory:a consistent way for including the gravity effect in (visco) elastic plane-earth models. Geophysical Journal International, 161(1): 191-196. DOI:10.1111/gji.2005.161.issue-1 | |
| Wang R J. 2005b. On the singularity problem of the elastic-gravitational dislocation theory applied to plane-Earth models. Geophysical Research Letters, 32(6): L06307. | |
| Wang W H, Chen C H. 2001. Static stress transferred by the 1999 Chi-Chi, Taiwan, earthquake:Effects on the stability of the surrounding fault systems and aftershock triggering with a 3D fault-slip model. Bulletin of the Seismological Society of America, 91(5): 1041-1052. | |
| Xiao J, He J K. 2015. 3D finite-element modeling of earthquake interaction and stress accumulation on main active faults around the Northeastern Tibetan Plateau edge in the past~100 years. Bulletin of the Seismological Society of America, 105(5): 2724-2735. DOI:10.1785/0120140342 | |
| Zhang B, Zhang H, Shi Y L. 2015. Equivalent-body force approach on modeling elastic dislocation problem using finite element method. Chinese Journal of Geophysics (in Chinese), 58(5): 1666-1674. DOI:10.6038/cjg20150518 | |
| Zhou Y, Shan B, Xiong X. 2008. Parameters sensitivity analysis of Coulomb stress change in static stress triggering. Journal of Geodesy and Geodynamics (in Chinese), 28(5): 21-26. | |
| Zhu S B, Miao M. 2016. On the study of earthquake triggering:Solution to paradox that Coulomb stresses increase with frictional coefficients. Chinese Journal of Geophysics (in Chinese), 59(1): 169-173. DOI:10.6038/cjg20160114 | |
| Ziv A, Rubin A M. 2000. Static stress transfer and earthquake triggering:No lower threshold in sight?. Journal of Geophysical Research:Solid Earth, 105(B6): 13631-13642. DOI:10.1029/2000JB900081 | |
| 曲均浩, 蒋海昆. 2012. 余震活动机理研究综述. 中国地震, 28(2): 109–120. | |
| 瞿武林, 张贝, 黄禄渊, 等. 2016. 2004年苏门答腊地震的几个断层滑动模型的全球同震位移对比. 地球物理学报, 59(8): 2843–2858. | |
| 单斌, 李佳航, 韩立波, 等. 2012. 2010年MS7.1级玉树地震同震库仑应力变化以及对2011年MS5.2级囊谦地震的影响. 地球物理学报, 55(9): 3028–3042. | |
| 石耀霖, 朱守彪. 2006. 用GPS位移资料计算应变方法的讨论. 大地测量与地球动力学, 26(1): 1–8. | |
| 石耀霖, 曹建玲. 2010. 库仑应力计算及应用过程中若干问题的讨论——以汶川地震为例. 地球物理学报, 53(1): 102–110. | |
| 万永革, 沈正康, 曾跃华, 等. 2007. 青藏高原东北部的库仑应力积累演化对大地震发生的影响. 地震学报, 29(2): 115–129. | |
| 张贝, 张怀, 石耀霖. 2015. 有限元模拟弹性位错的等效体力方法. 地球物理学报, 58(5): 1666–1674. | |
| 周宇明, 单斌, 熊熊. 2008. 静态应力触发中影响库仑应力变化的参数敏感性分析. 大地测量与地球动力学, 28(5): 21–26. | |
| 朱守彪, 缪淼. 2016. 地震触发研究中库仑应力随摩擦系数增加而增大的矛盾及其解决. 地球物理学报, 59(1): 169–173. | |