鲜水河断裂带位于青藏高原东缘，是巴颜喀拉块体与川滇块体的边界断裂.断裂带形变观测资料丰富，历史强震频发，尤其是近年来断裂周缘陆续发生汶川、玉树、芦山等系列强震，是探讨形变变形机理与强震关系的理想场所.自2013年芦山M_S7.0地震发生到2014年1月，断裂带北西段跨断层测点短基线测项均出现不同程度的异常变化，尤其是侏倭AB测段、虚墟AC测段变化较为显著.通过对测点周边环境、点位稳定性实地考察及复测后，分析认为短基线场地附近无明显的环境、人为干扰，水准结果无明显变化，基线异常属实¹⁾.那么，该异常变化是否与芦山地震的发生存在关联性？值得深入研究.

1)四川省地震局测绘工程院.2014.2014年3月鲜水河断裂北段跨断层短基线异常核实报告.

从力学机理分析，强震是孕震体应力(能量)的积累、进入临界状态、最终失稳的过程.强震的过程可以清楚地划分为震间、同震和震后三过程(Scholz, 1998).近些年，随着大地测量、地球物理观测技术的发展和震例的积累，对于同震和震后地表变形的演化特征及物理机制研究获得系列认识.已有的研究表明，下地壳或上地幔的黏弹性松弛效应(Savage and Prescott, 1978; Reilinger, 1986; Pollitz and Sacks, 1992; Shen et al., 1994; Ueda et al., 2003; Paul et al., 2007; 邵志刚等，2008；朱守彪和蔡永恩，2009; Rollins et al., 2015; Shao et al., 2016)、震后余滑(Savage et al., 1994; Hsu et al., 2002)、孔隙弹性回跳(Peltzer, 1996; Pollitz et al., 2001; Fialko, 2004)是产生强震后地表变形的主要原因.Savage和Prescott (1978)构建包含岩石圈及软流圈的二维理想模型，分析软流圈黏弹性效应对地表变形的影响，结果认为当震源区深度与岩石圈厚度相当时影响较大.Shen等(1994)基于位错理论的分析认为来自地壳深部的震后余滑及黏弹性松弛效应能够较好地解释Landers地震后远场GPS观测资料的变形特征.针对昆仑山M_S8.1地震后东昆仑断裂两侧形变差异及同向运动特征，邵志刚等(2008)利用有限元方法分析介质的非均匀性、黏弹性松弛及孔隙流体作用对震后形变的影响.结果认为断裂两侧下地壳黏滞系数的差异造成震后形变的非对称性和同向运动.Savage等(1994)认为Loma Prieta地震震后余滑所产生的右旋逆冲运动可以较好地解释震后地表的变形特征.Fialko (2004)研究认为Landers地震与Hector Mine地震间近7年的地表变形特征主要来自于震后孔隙流体作用.相反，Barbot等(2008)则认为震后孔隙流体作用、震后余滑、黏弹性松弛效应都不是造成震后地表变形的主要原因，而余震的发生可以较好地解释震后的变形特征.另外，Masterlark等(2001)基于三维有限元模型，以墨西哥Colima-Jalisco M_W8.0地震为例，分析地壳介质各向同性或各向异性对同震及震后变形的影响.结果认为介质各向异性时模拟结果能够较好地解释实际GPS观测资料的变化.

本文基于区域地质构造、地球物理、大地测量、断层力学等多学科资料，构建鲜水河断裂带及邻区三维非线性黏弹性有限元模型.以GPS观测资料、芦山同震静态滑移量作为约束条件，通过开展多组数值模拟计算，探讨下地壳不同流变特征下，芦山地震发生前后不同阶段，不同跨断层基线的时序动态变化特征.通过与实际观测资料进行定性、半定量对比分析，尝试性地探讨基线异常变化机理.

以鲜水河断裂、小江断裂及龙门山断裂为主的“Y”字形构造体，将川西及邻区划分为川滇菱形块体、巴颜喀拉块体及华南块体三大活动块体(图 1).作为巴颜喀拉块体与川滇菱形块体的边界断裂，鲜水河断裂运动特征的研究对于理解青藏高原东向运移的物质分配具有重要意义.断裂包含结构特征简单的北西段和由若干分支断层组成的南东段，其中北西段由炉霍断裂、道孚断裂、乾宁断裂3条次级剪切断裂呈左阶羽列组合而成(闻学泽等，1989).侏倭、虚墟观测场地位于炉霍次级断裂，均有超过30年的观测历史，观测资料较多地用于探讨形变异常特征与地震活动的关系(吕弋培等，1997；李瑞莎等，2009；牛安福等，2013；苏琴等，2014).

图 1

Fig. 1

图 1 川西地区地质构造及形变测点分布图 Fig. 1 Map showing tectonics of western Sichuan and distribution of deformation measurement sites

虚墟、侏倭场地均为基线、水准综合观测场地，观测资料连续性好，观测精度高.侏倭场地自1980年4月正式观测，已有30多年的观测历史.2013年4月20日芦山M_S7.0地震发生后，该场地AB基线长度减小1.71 mm，从2013年5月份开始到2014年1月，基线长度快速增加，变化量达3.54 mm，2014年2月至今处于高值波动状态(图 2a).虚墟场地自1973年2月正式观测，已有40多年的观测历史.芦山地震发生后该场地AC基线长度迅速减小3.57 mm，从2013年5月份开始到2014年1月，基线累积增加量达3.93 mm，2014年2月至今处于高值波动且有增加的趋势(图 2b).针对两个测点的异常变化特征，通过开展实地场地环境、仪器工作状态核实及重复测量工作，没有发现明显的干扰现象，认为两场地观测资料异常确实存在¹⁾.

图 2

Fig. 2

图 2 基线测点异常变化特征 Fig. 2 Anomalous characteristics of baseline observation stations

由于巴颜喀拉块体受华南块体的阻挡作用，从而引起巴颜喀拉块体与川滇菱形块体南东向运动的差异，造成块体边界鲜水河断裂的左旋走滑运动.而作为巴颜喀拉块体与华南块体的边界断裂，龙门山断裂带上芦山地震的发生势必会对巴颜喀拉块体与川滇菱形块体的运动差异产生影响，从而造成鲜水河断裂运动特征的变化.因此，鲜水河断裂带的短基线观测表现出对应芦山M_S7.0地震的协同变化.

跨断层观测资料不仅包含断层活动的构造信息，也包含着非构造活动的信息(刘冠中等，2013).本文进一步利用多项式分段曲线拟合，获得去趋势后残差值分布图(图 3).从图 3可以看出，芦山地震发生前基线长度残差变化幅度较为平稳，而在地震发生后半年多的时间内基线出现快速的伸长变化，这种短时间内大幅度的变化可能蕴含着与芦山地震相关的信息.

图 3

Fig. 3

图 3 多项式分段曲线拟合残差值 Fig. 3 Residual values derived from polynomial subsection curve fitting

基于区域地质构造、地球物理场、大地测量等多学科资料，充分考虑对区域构造变形起主要影响的活动断裂，建立鲜水河断裂带及邻区三维地质模型(图 4a).在地质模型的基础上，利用单元离散化及断层区网格优化技术，建立目标区三维黏弹性有限元模型(图 4b).模型范围为99°E-105°E，26°N-33°N，深度30 km，包含上地壳和下地壳两层各15 km；鲜水河断裂、安宁河断裂、则木河断裂及龙门山断裂均处理为接触摩擦单元，其中龙门山断裂处理为倾滑断层，其他近似处理为垂直断层；模型各分区物性参数主要依据研究区已有的研究成果(Liu et al., 2014).上地壳、下地壳弹性模量变化范围分别为0.8×10¹¹~1×10¹¹Pa，0.75×10¹¹~1.3×10¹¹Pa，泊松比统一取0.28；断层杨氏模量为周围介质的1/10，泊松比0.3，摩擦系数0.4；有限元模型包含1060783节点、683034单元.模型边界约束参考2004-2007年GPS观测结果，模型表面自由，底部法向约束(图 4a).有限元分析中，应力函数以积分形式给出，在小应变理论下，各向同性材料的黏弹性本构关系表示为：

图 4

Fig. 4

图 4 鲜水河断裂带及邻区三维有限元模型(a)三维地质模型；(b)三维有限元模型. Fig. 4 Three-dimensional finite element model of the Xianshuihe fault zone and adjacent areas (a) Three-dimensional geological model; (b) Three-dimensional finite element model.

(1)

其中，σ为Cauchy应力，G(t)为剪切松弛核函数，K(t)为体积松弛核函数，e为偏应变，Δ为体应变，t为当前时间，τ为过去时间，I为单位张量.

黏弹性积分核函数G(t)和K(t)主要有广义Maxwell和Prony级数两种形式，本文采用Prony级数来表征介质的黏弹性特征，核方程如下：

(2)

其中，G_i和G_∞为剪切模量，K_i和K_∞为体积模量，τ_i^G和τ_i^K为各Prony级数分量的松弛时间，n_G和n_K为Prony单元的个数.

对于常载荷作用下的黏弹性问题，当变形时间远远小于松弛时间时，岩体没有屈服情况下以弹性变形为主，在岩体屈服时以塑性变形为主；当加载时间远远大于松弛时间后，岩体进入静岩压力状态，不产生剪切变形.本文研究岩体在不同阶段的变化特征，重点分析来自下地壳的流变性质对地表观测的影响.

为更好地理解芦山地震发生前后不同阶段，不同跨断层基线的异常变化特征，结合前人有关该区黏滞系数的研究成果，分别建立模型1(巴颜喀拉块体、川滇块体下地壳黏滞系数为10¹⁹Pa·s，华南块体为10²¹Pa·s)与模型2(巴颜喀拉块体、川滇块体下地壳黏滞系数为10¹⁸Pa·s，华南块体为10²¹Pa·s)，分析不同跨断层基线的时序动态变化特征.模型1中各点等效应力随时间的演化见图 5，其中AB、CD测线分别对应侏倭BA、虚墟AC基线.

图 5

Fig. 5

图 5 模型中测线分布及A~D点应力变化曲线 Fig. 5 Baseline distribution and stress change curves of points A~D

本文采用鲜水河断裂及邻区的三维有限元模型进行计算，计算的总时间步长为10.5年，包含3个时间步，21个时间子步.首先GPS边界位移约束随时间逐步加载到模型中，获得区域的相对初始应力场；在此基础上，以芦山地震的同震静态滑移量作为约束(王卫民等，2013)，通过将滑移量加载到模型中相应节点上，模拟芦山地震的发生；最后分析地震发生后半年时间尺度的黏弹性松弛效应.

Von Mises应力是基于剪切应变能的一种等效应力，单位体积形状改变的弹性位能(弹性形变能)达到某一常数时，材料开始屈服，可以表述为三个主应力的函数:

(3)

从模型1中A~D点的等效应力演化曲线可以看出(图 5)，介质经历地震发生前的应力积累、发生时瞬时的应力增加及地震后黏弹性松弛效应造成的应力减小的变化过程.

为分析构造作用下鲜水河断裂的运动学特征，利用跨越断层的AB、CD两条测线(测线分布见图 5)，通过分析断层面两侧不同方向的错动量变化，确定断层的活动方式.从图 6两条测线的变化特征可以看出，在以GPS观测数据为约束的块体间差异性运动的动力学边界作用下，鲜水河断裂表现为左旋走滑运动，进一步表明模型的合理性.

图 6

Fig. 6

图 6 基于模拟结果的断层运动性质 Fig. 6 Fault kinematics from numerical modeling

在构造运动等动力作用下，断层两侧的点位会发生变化.以侏倭测点为例(图 2)，假设在t时刻A、B两点的位移分别为(x_A，y_A，z_A)、(x_B，y_B，z_B).AB之间的长度为L_AB，伸长量变化为ΔL_AB.那么，在t时刻，AB间位移的变化量可近似表述为：

其中，伸长量变化增加时，表示断层两侧为拉张运动；减小时则表示为压缩运动.

从图 7可以看出，鲜水河断裂两侧的巴颜喀拉块体及川滇菱形块体在下地壳黏滞系数η=10¹⁹Pa·s时，基线AB间距离变化表现为地震发生前长度增加，发生时减小及地震后逐渐恢复的过程.地震发生后2个月内由于下地壳黏弹性松弛效应的影响，造成AB间距离快速增加到0.1916 m，此后进入缓慢增加的变化趋势.而在下地壳黏滞系数η=10¹⁸Pa·s时，基线AB间距离变化表现出相似的变化趋势，且黏弹性松弛效应的作用造成AB间距离变化快速增加到0.1951 m.

图 7

Fig. 7

图 7 下地壳不同黏滞系数时基线间距离变化 Fig. 7 Baseline length changes with different lower crust viscosities

而对于CD基线，下地壳黏滞系数η=10¹⁹Pa·s或者10¹⁸Pa·s时，基线间距离变化均表现为地震发生前长度增加，发生时减小及地震后逐渐恢复的过程.且黏弹性松弛效应的作用造成CD间距离变化分别增加到0.1935 m及0.1977 m.

可以看出，鲜水河断裂带两侧地块的下地壳黏滞系数在合理取值范围内，在构造应力作用下，基线AB、CD间长度存在增加的趋势；地震发生时基线长度减小；地震后2月短时间内，由于下地壳黏弹性松弛效应造成基线AB、CD间距离快速增加.

已有研究结果表明，黏弹性松弛效应、震后余滑、孔隙弹性作用是产生强震后地表变形的主要原因.震后余滑和孔隙弹性作用主要影响近场，而黏弹性松弛效应影响范围较广(谭凯等，2005).考虑到形变观测点与芦山地震震中较远，本文仅分析黏弹性松弛效应对地表变形的影响.研究结果表明，实际观测资料所反映的震后基线长度快速变化的特征可能是下地壳黏弹性效应引起的应力调整与构造作用共同作用的结果.

岩石圈流变结构影响介质中应力传递，进而影响地表变形.Flesch等(2001)利用偏应力和应变速率资料计算青藏高原及周边地区的平均等效黏滞系数，结果显示中下地壳等效黏滞系数在10¹⁸~10²¹Pa·s；石耀霖和曹建玲(2008)给出的川滇菱形块体25 km深处等效黏滞系数在10¹⁹~10²¹Pa·s.张晁军等(2008)基于PSGRN/PSCMP的多次模拟实验认为，鲜水河断裂带周缘30 km深处下地壳黏滞系数为10¹⁹Pa·s时，模拟结果与跨断层实测结果吻合较好.Wen等(2012)认为中下地壳黏滞系数为2×10¹⁹~5×10¹⁹Pa·s时，可以较好地解释InSAR反演的昆仑山口西地震震后变形特征.邵志刚等(2007)认为，用Maxwell体或Kelvin体模拟震后形变，黏滞系数较小时，黏弹性松弛效应很快平衡.而黏滞系数较大又无法解释震后应变率随时间衰减的现象.近年来，大地电磁测深及地球物理场波速结构的研究表明川西地区下地壳流的存在(Unsworth et al., 2005; Liu et al., 2014).为确定不同黏滞系数的选取对跨断层基线变化特征的影响，本文进一步分析巴颜喀拉块体、川滇块体及华南块体的下地壳黏滞系数分别为10¹⁹Pa·s、10¹⁸Pa·s和10²¹Pa·s时，不同跨断层基线的变化特征(图 8).

图 8

Fig. 8

图 8 下地壳不同黏滞系数时基线距离变化 Fig. 8 Baseline length changes with different lower crust viscosities

可以看出，巴颜喀拉块体、川滇地块下地壳黏滞系数分别为10¹⁹Pa·s及10¹⁸Pa·s时，地震发生后短时间内由于黏弹性松弛效应的作用造成AB、CD基线长度快速增加.且黏滞系数的选取不会造成基线长度变化趋势的改变.

基于有限元模型的模拟结果显示，地震后较短的时间尺度内，由于黏弹性松弛效应的影响，地表跨断层基线长度存在快速增加的特征.以模型中AB基线为例，从图 9地震发生前后实际观测资料与模拟结果的对比分析来看，两者所反映的基线长度的变化趋势是一致的.受模型精细程度影响，模拟给出的变化量很难与实际观测结果完全吻合.但从变化趋势上来看，地震后短时间内实际观测资料的快速变化可能是黏弹性松弛效应与构造作用共同作用的结果.考虑到地震后地表变形机制的复杂性，其他诸如震后余滑、孔隙流体作用等都不同程度地影响地表变形，因此在获取可靠的多学科观测资料后，需进一步开展震后地表变形机制的综合分析.

图 9

Fig. 9

图 9 地震发生前后实测资料与模拟结果的对比分析 Fig. 9 Comparison between observed and simulated results before and after the earthquake

芦山M_S7.0地震发生前，侏倭AB基线同样存在显著的观测异常，而处于同一次级断裂的虚墟基线长度变化则不明显.分析原因认为，首先这种变化可能与断裂的无震滑移有关.从活动构造来看，该区域主要发育有以左旋走滑为主的鲜水河-安宁河-则木河断裂系和右旋走滑兼逆冲性质的龙门山断裂带.如果是左旋走滑型活动断裂发生无震滑移，以鲜水河断裂带南东段康定次级断层为例，从图 10可以看出，该次级断层发生无震滑移时，AB基线变化所反映的断层运动性质与断层自身运动性质是一致的，即同样为左旋运动，这与实际观测不符；如果是右旋兼逆冲性质的龙门山断裂发生无震滑移事件，从本文针对芦山地震的研究结果来看，地震发生时两条基线的变化所反映的断层运动性质为右旋运动，可以较好地解释芦山地震前侏倭AB基线的变化特征，而虚墟基线长度变化不明显的现象可能存在其他原因，例如可能与断层不同位置断层自身的力学状态有关.如果虚墟测线所跨断层区摩擦系数较高、甚至处于闭锁状态，可以较好地解释虚墟基线长度变化不明显的问题.

图 10

Fig. 10

图 10 康定断裂无震滑移对AB基线长度变化的影响红色实线为破裂带. Fig. 10 Effect of aseismic slip on the Kangding fault on the change of AB baseline length Red solid line denotes rupture zone.

从断层不同摩擦系数对模拟结果影响的分析来看(图 11)，当断层摩擦系数较大时，基线间长度变化明显减小.由于目前缺乏关于断层力学性质的岩石物理实验结果，因此只能通过数值模拟的方法进行分析推测.

图 11

Fig. 11

图 11 断层不同摩擦系数时AB基线长度变化 Fig. 11 Length changes of baseline AB with different fault frictional coefficients

基于鲜水河断裂带及邻区三维黏弹性有限元模型，通过开展多组数值模拟实验，探讨下地壳不同流变特征下，芦山地震发生前后不同阶段，不同跨断层基线的时序动态变化特征.初步研究结果表明，侏倭、虚墟基线测点原始观测资料表现出地震发生时基线长度减小，发生后短时间内反向快速恢复及增加的协同变化.在以GPS观测数据为约束的块体间差异性运动的动力学边界条件下，鲜水河断裂表现为左旋走滑运动.芦山地震发生时，测点间距离迅速减小，表现为断层右旋的特征，这与实际观测资料显示的结果一致.鲜水河断裂带两侧地块的下地壳黏滞系数分别取10¹⁸Pa·s、10¹⁹Pa·s时，由于黏弹性松弛效应的影响，造成基线长度在地震后短时间尺度内增加.对比分析认为，地震发生后实测资料所反映的基线长度快速变化的特征可能是黏弹性效应与构造作用共同作用的结果.

感谢中国地震局地震预测研究所邵志刚研究员、中国地震局地壳应力研究所陆远忠研究员、刘冠中博士的指导与帮助，感谢匿名审稿专家提出的宝贵建议，部分图件使用GMT软件包绘制(Wessel and Smith, 1995)，在此一并表示感谢.