断裂两盘岩性差异对汶川地震的影响

姚 琪，邢会林，徐锡伟，张 微

1 活动构造与火山重点实验室，中国地震局地质研究所，北京 100029

2 中国地震台网中心，北京 100045

3 Earth Systems Science Computational Centre（ESSCC），The University of Queensland，St.Lucia，QLD 4072.Australia

4 国土资源厅航遥中心，北京 100083

1 引 言

造山楔与前陆盆地之间往往存在巨大的岩性差异，而随着前陆盆地的演化，早期沉积层逐渐被卷入逆冲褶皱带中，进一步增大了造山楔与前陆之间的岩性差异.板块构造背景和多套滑脱层固然是控制造山带断层活动的决定性因素，然而岩性影响也不可忽视.岩性差异导致层间变形差异是常见的构造变形现象，能够造成有限应变状态差异、劈理折射、香肠构造、褶皱形态差异等.在较短的地质时期内，譬如一个地震复发周期，岩性也能影响断层位移，以及位移与断裂带长度的线性关系[1］，而在强地震动的传播中，如果发震断层位于岩性差异的界限，那么物性软弱的一盘则会产生更大的变形[2-4］.

位于青藏高原东缘的龙门山断裂带为巴颜喀拉块体向南东滑动与扬子地块发生碰撞而形成，代表了青藏亚板块与华南板块的碰撞边界.该断裂带中段北川—映秀断裂虹口至清平地区，断裂上盘为坚硬的前震旦系褶皱基底，下盘则为软弱的前陆盆地沉积物.2008年汶川8.0级地震震后地质考察研究表明，映秀—北川地表破裂带龙门山镇—清平段，沿着北川—映秀断裂分布，即沿着彭灌杂岩与紧邻的上三叠统须家河组的边界分布（图），而伴随地震断层出露地表的滑动面大多沿炭质泥岩和煤层发育[5］，汶川地震特有的两条平行破裂同时参与同震破裂的段落[6］基本上与彭灌杂岩长轴位置大致相当.近地表三维同震滑移矢量显示2008年MS7.9汶川地震沿着北川—映秀断裂分布的地表破裂带和沿着灌县—江油断裂分布的地表破裂带是深部斜滑断裂在上地壳脆性域发生应变分解的结果[7］.北川—映秀断裂上下盘岩性的差异对汶川地震的孕震、发震以及传播过程是否有影响，汶川地震多重地表破裂带与彭灌杂岩的空间对应关系是否属于偶然，岩性差异在断裂活动的应变分解中是否有作用，作用多大，均是值得我们关注的问题.

基于龙门山地区的人工地震反射剖面，震源机制解以及汶川地震地表破裂带，建立上地壳范围内北川—映秀断裂2D构造模型，以及北川—映秀断裂与灌县—江油断裂共存的双断坡模型，用基于R-minimum的有限元算法对脆性地层中的弹性断层活动进行非线性摩擦接触模拟，分别用杨氏模量、泊松比的差异来代表岩性差异，尝试了24对上下盘材料参数共存的情况，获得了在一个地震复发间隔内单条断层的破裂时间、破裂过程、下盘变形强度随着断层两盘岩性差异的变化，并对比单条断层和双断坡构造在上下盘岩性一致或是同时具有强硬上盘与软弱下盘情况下的断层行为与块体变形特征，来分析上下盘之间岩性差异对汶川地震双断坡构造的影响.

2 区域地质概况

龙门山断裂带由后山断裂、中央断裂、前山断裂和山前隐伏断裂这4条近于平行的断裂组成（图1，图2），其中后山断裂以西为松潘一甘孜褶皱带，后山断裂带与中央断裂带之间为韧性变形带，中央断裂带与前山断裂带之间为基底卷入冲断带，前山断裂带与山前隐伏断裂系之间为前缘—褶皱冲断带，山前隐伏断裂系以东为前陆坳陷带.

自南而北，分别有宝兴杂岩、彭灌杂岩和轿子顶杂岩出露于龙门山断裂带之中.这些杂岩为混合岩化的变质岩和混合岩，变质岩主要为斜长角闪岩和变粒岩，混合岩则主要为黑云斜长混合岩、角闪斜长混合岩等，质地较古生代和中生代基岩更为坚硬.其中，彭灌杂岩作为龙门山地区出露面积最大的杂岩，其南边界与北川—映秀断裂一致，大部分北边界与汶川—茂汶断裂一致，而以含碳泥岩为主的上三叠统须家河组地层分布于映秀至汉旺一带，长轴范围与紧邻的彭灌杂岩体大致相当.彭灌杂岩内部普遍而广泛地分布着一系列近叠瓦状的走向大致北东、倾向以北西向为主的断裂系统，且在愈靠近茂汶断裂，断裂发育密度愈大，杂岩体变形的韧性程度也愈高，而近映秀断裂一侧则多具脆性断裂[10］，反映出该杂岩体与龙门山断裂活动具有密切的联系.

在彭灌杂岩的前缘，则是四川盆地广泛的沉积地层，包括以泥质岩和碳酸盐岩建造为主的志留系，碳酸盐岩与碎屑岩不等厚互层的泥盆系，以及以石英砂岩与页岩互层为主，并夹煤层的上三叠统须家河组.这些沉积地层属于“非能干层”，在同样的应力作用下，成为易于滑动，有利于断层活动的软弱层，也是龙门山断裂带山前隐伏薄皮构造发育的滑脱层[9，11］.

图1 龙门山地区区域构造简图（a）及地层柱状图[6，8-9］（b）GJF：灌县—江油断裂；BYF：北川—映秀断裂；MWF：茂汶—汶川断裂.Fig.1 Simplified geologic maps of the Longmenshan thrust belt（a）and adjacent regions stratigraphic column（b）GJF：Guanxian-Jiangyou Fault；BYF：Beichuan-Yingxiu Fault；WMF：Wenchuan-Maowen Fault.

图2 （a）汶川地震发震构造模型；（b）人工地震剖面；（c）横跨龙门山地区的深地震反射剖面[6，11-12］（图a位置见图c的黑框范围，图b位置见图a黑框范围）Fig.2 （a）Seismogenic structure of the Wenchuan earthquake，（b）petrol seismic reflection profile，（c）and deep seismic reflection profile（black frame in figure（c）refers to the location of figure（a），and the black frame in figure（a）refers to the location of figure（b））

石油地震剖面（图2）及上三叠统须家河组地层等厚图（图1）揭示，彭灌杂岩的前渊是川西晚三叠纪时期的周缘前陆盆地的沉降中心，且在此基础上又叠加了晚白垩纪—早第三纪再生前陆盆地的相关沉积，印—藏碰撞的持续挤压作用又使得晚新生代构造变形不断向东扩展进入川西盆地南部[9］，导致该处沉积地层厚度较厚，卷入层次较深.深地震反射剖面[12］显示，北川—映秀断裂两侧的岩性差异直至深度20km的范围内都是存在的，而在龙门山断裂带的西侧，地层仍旧保持了较高的P波速度，表明在深达20km的范围内，不仅仅是出露的彭灌杂岩这一局部范围，龙门山断裂带西侧的地层较之东侧的地层更为强硬.

由此可见，虽然汶川地震是较软的巴颜喀喇块体在印度板块向亚洲板块俯冲的作用下的向东挤压受到冷硬的扬地块体阻挡而长期积累的应力释放而成，而在上地壳，则是较硬的受强烈的褶皱变形的结晶基底冲断席在块体运动的作用下沿着巴颜喀拉块体深部约20km处的低阻低速层向东滑动，作用于较软的四川盆地前寒武纪基底和上覆的古生代和中生代沉积盖层组成，而龙门山断裂带中央断裂则是这两个地层岩性差异巨大块体的分界.

3 摩擦接触有限单元数值模拟方法

大部分的地震发生在先存的断裂构造上，地震的孕育、发生过程中也就是断层闭锁、滑动、解锁过程，其中断层面上的摩擦关系起着决定性的作用.然而在断层的活动中，断层面上的摩擦系数数值较大（0.2～0.7），对环境条件也更为敏感（如速度、压力等），常用于计算断层非线性摩擦问题的静态隐式算法可能会导致不收敛问题，建立复杂地质模型所需的光滑技术可以减少结果的多样化，但仍不能保证收敛，并可能黏滑失稳时停止计算.

Xing and Makinouchi（2002，2003）[13-14］应用统一的数学公式表述了速率相关的摩擦接触中黏着（sticking）和滑移（sliding）这两种不同的运动状态；有限元计算中采用静力显示的时间积分方法，基于R最小策略，控制时间步长以保持力学状态变化稳定，从而保证有限元计算过程平稳、收敛.该程序能够计算在外边界的均匀载荷下，模型内部的接触面上摩擦状态黏着和滑移两种过程的变化，即断层上闭锁和解锁两种不同运动状态的变化特征，也就是地震的孕育、发生的过程.

有限元公式则是采用更新拉格朗日列式算法来描述非线性接触问题，平衡方程可以表述如下：

公式中上标的（～）表示变形块体之间的相对值，且l，m=1，2；i，j，k=1，2，3.其中V 和S 分别表示t时刻变形体B的体积和表面积；ρ为密度，SF是表面S 上外载的作用边界；虚拟速度场，在速度边界上满足是 柯 西 应 力Jaumann率；σij为Kirchhoff应力张量；L为速度梯度张量，且Lij=∂vi／∂xj，δLij是虚拟速度梯度张量；D 和W是张量L的对称和反对称部分；是接触表面Sc上接触力的速率，是接触对（点与点接触）之间的相对滑动速度.

在法向力的计算中，则运用罚函数方法处理法向约束，使之在接触发生时法向不产生穿透，即：

f为法向接触力，En为罚因子，n为从接触体表面的外法线方向，gn为法向穿透距离.

采用Mohr-Coulomb摩擦模型来描述接触面黏着（sticking）和滑移（slipping）的摩擦行为，类似于理想弹塑性材料的屈服函数，这里将黏着状态类比为“弹性”，滑移过程比作“塑性”.为了避免状态变量计算导致的计算困难，采用速率相关的Coulomb摩擦准则，使用将摩擦中黏着与滑移两部分分开的办法处理，而不考虑状态量变化的效应，若也忽略温度影响，则摩擦力可表示为[13-15］：

在式（3）（4）中，μ为摩擦系数，它可能与多种变量有关，如法向接触力fn，等效切向速度状态变量φ，即是沿接触面切向的一个常数；F是临界摩擦力始值；是的滑动部分（不可逆），则是的闭锁部分（可逆）.

文中采用显式的时间积分方法，通过时间步长增量的比例控制系数Rmin（＜1）来控制载荷增量的大小.通过选取适当的时间步长来保证每个增量步内单元的力学状态和界面的接触状态平稳变化.

用式（5）中的增量替代方程（1）中的所有率变化量，则公式（1）可表示为式（6）的形式

式中K为整个物体的标准刚度矩阵，ΔF为力边界的外荷载增量，ΔFf为接触力增量，Kf为所有接触单元的接触刚度矩阵，Δu为节点位移增量.

4 模型设置

根据USGS地震记录，汶川地震起始破裂处距北川—映秀地表破裂带垂向距离约17km，而地震震源 深 度 约 19km（USGS，http：／／earthquake.usgs，gov／earthquakes／eqinthenews／2008／us2008ryan），若简化断面为平直的面，则断裂倾角为48°.汶川震源机制解则显示最佳双力偶解的节面I倾角可能为33°、59°、32°、39°，与之对应的节面II倾角为70°、47°、63°、57°[19-20］.而在静水孔隙压力条件下，低倾角逆冲断裂是指倾角30°及其以下的断裂，高倾角逆冲断裂是指45°以上的断裂.由此，本文以45°倾角的平直断面做为北川—映秀断裂的主要几何形态，并用与之相切的曲面将主体断面与深部的滑脱构造连接起来（模型I，图3），连接位置为19km深处，即汶川地震震源深度.

Xu Xiwei[6］和Jia Dong等[11］（2010）通过石油人工地震勘探剖面解释和断层相关褶皱模型提出北川—映秀断裂与灌县—江油断裂在地下可能相连，但连接处深度不超过10km（图2a，b）.但是由于人工地震勘探的深度与范围的限制，这两条断裂在深部的连接方式和共同破裂机制尚存在争议.王卫民（2009）[21］用深部不相连的双断层模型模拟表明地震波能够从发震断层传递到不相连的其他断层上.考虑到龙门山构造带的具有沿倾向层次渐浅、强度递减、卷入层位变新的趋势，呈现前展式扩展的特征[22］，灌县—江油断裂的倾角应小于或等于北川—映秀断裂的倾角，为了避免不同断层倾角对结果造成混淆，本文将灌县—江油断裂简化为倾角45°，向下延伸至10km，与北川—映秀断裂相距14km的次级断层，建立双断坡模型（模型II，图3）.

龙门山断裂带距其西北侧NE走向的龙日坝断裂带约185km，距其东南侧NE走向的龙泉山背斜约100km.右旋走滑兼逆冲性质的龙日坝断裂带可能分担了一部分巴颜喀啦块体向四川盆地的挤压，而川滇地区的地形局部高程差显示主要的高程变化出现在龙日坝断裂带与龙门山断裂带之间[23］.龙泉山背斜是川西前陆盆地沉积边缘，代表了迄今造山带对前陆的最大影响范围.由此，为了避开龙日坝断裂带活动的影响，并包含龙门山断裂活动对四川盆地的作用，模型中断层上盘长约60～88km，下盘长度约112～140km，模型总长度200km（图3）.

图3 几何模型及边界条件设置Uxy为所有方向均固定，Ux为X方向固定，Uy为Y方向固定，DUx为Y方向固定，而在X方向以恒定的速率3×10-4 Vref施加位移边界条件，Vref为相对速度.Fig.3 Models and boundary conditions Here Uxy-fixed at all the coordinates；Ux-fixed at Xcoordinates；Uy-fixed at Ycoordinates；DUx-displacement applied along the X coordinate by a constant velocity of 3×10-4 Vref，but fixed at the Ycoordinate；Vrefis the reference velocity for the entire process.

鉴于巴颜喀喇块体之下的低阻低速层对龙门山断裂活动的重要作用，在断层上盘之下设置一2km厚的软弱层，在下盘则设置2km厚的强硬层，而上下盘下部的边界条件均设置为纵向固定，横向自由.这样不仅可以使得上盘沿着该软弱薄层发生滑动并止步于强硬薄层，使得应力能够在断层上积累，也能够限制上下盘的运动方向，避免畸变或是块体扭曲发生，由此，模型总高度26km，如图3所示.

有限元建模及计算结果后处理利用MSC.PAtran商业软件进行，采用Paver网格对模型进行划分，网格长度1km，其中模型I中有8节点六面体单元21512个，节点数为28250，模型II中有8节点六面体单元21648个，节点数为28530.模型南东边界设置为横向固定，纵向自由，模型南东末端节点固定，模型上表面为自由边界.震前的GPS观测表明，横跨整个龙门山断裂带的滑动速率不超过～2mm／a，单条断裂的活动速率不超过～1mm／a[24］，由此将3×10-4Vref的持续位移直接加载在模型北西边界（图3）.

虽然青藏高原内部具有弱负均衡重力异常，而龙门山地区为正均衡重力异常，四川盆地则为负均衡重力异常区[25-26］，这种重力不均衡的分布表现了龙门山地区的高地形，以及对应的地壳厚度自龙门山地区向四川盆地减薄，可能导致断层复发间隔增大且破裂过程变长.然而本文主要模拟对象为深度仅20余km的上地壳的弹性变形，并不包括厚达46～60km的整个地壳[27］，中上地壳之间的低阻低速层分隔了上地壳变形与中下地壳变形，上地壳的重力不均衡更多的来自于地形的变化，此外，千余年尺度的大地震复发间隔[28］对百万年尺度的重力均衡所需时间[29］而言极小.本文计算中初始地形为水平面，上下盘选择一致的密度，且主要模拟目标是一个地震复发周期内侧向挤压条件下的断层行为，重力作用的影响是极小，因此本文不考虑重力作用.

接触面上的摩擦依速度而定：μ=0.60+（0.01-0.025）ln（V／Vref），μ为摩擦系数，V 和Vref分别为滑动速度及参考速度.

上地壳材料假定为线弹性，密度取为2.60g／cm3，考虑到彭灌杂岩和四川盆地沉积地层的岩性，参考四川地区水电站建设所作的岩石力学特性试验，在不考虑岩性差异的模型中（模型I-23；模型II-1），断层两侧杨氏模量取为44.80GPa（花岗闪长岩或石英细砂岩[30］），泊松比取为0.22（花岗闪长岩[30］），略大于地震台阵观测所的龙门山断裂附近地壳平均泊松比（0.20）[27］，稍小于粉砂岩的（0.25）[30］.

据GPS观测资料，汶川地震的成因可能是较为柔软的巴颜喀喇块体向东运动，受到冷硬的四川盆地的阻挡而积累大量应力，导致块体边界的逆断层破裂引发大地震[24，31］，然而深部地球物理勘探显示在上地壳范围内，断裂上盘巴颜喀拉地区较为强硬，下盘四川盆地较为软弱（图3），低阻低速层分隔了中下地壳与上地壳，中下地壳则是巴颜喀拉地区较软弱而四川盆地较强硬[12，27］.因此，本文在计算断层两侧岩性差异对上地壳断层行为的影响时，采用强硬上盘和软弱下盘的材料参数（模型I-1—22，模型II-2，表1）.

弹性材料属性参数常用的有杨氏模量（Young′s modulus）与泊松比（Poisson ratio），其中杨氏模量是描述固体材料抵抗形变能力的物理量，是物体弹性变形难易程度的表征，只与材料的化学成分有关，与其组织变化无关，而泊松比则是反映材料横向变形的弹性常数，对普通材料而言，这两者没有直接的关系.为了弄清这两个参数的变化分别对断层行为的影响，本文首先通过改变上下盘杨氏模量比值，并维持上下盘一致的泊松比（如模型II的无岩性差异模型，保持0.22）来计算杨氏模量的变化对断层行为的影响（模型I-1—16），然后通过改变上下盘的泊松比，并维持上下盘一致的杨氏模量（如模型II的无岩性差异模型，保持44.80GPa），来计算泊松比的变化对断层行为的影响（模型I-17—22），从而分别得到两个最重要的弹性材料参数的影响，并于无岩性差异的模型（模型I-23）进行对比，以确定龙门山断裂两盘材料参数（模型I-24），并用于双断坡模型中（模型II-2）.其中，根据龙门山断裂带上下盘的岩性（图1），杨氏模量的上限取为73.98GPa，即混合岩或是闪长岩、花岗岩，下限则取为7.69GPa，即泥质灰岩或泥质粉砂岩[30］.泊松比的下限略小于地震台阵观测所的龙门山地区地壳平均泊松比（0.20），取为0.19，泊松比上限则略大于断层下盘出露的奥陶系粉砂岩和志留系页岩泊松比（0.27），取为0.30[30］，如表1所示.

5 计算结果

5.1 上下盘岩性差异的影响

软弱的下盘能够通过自身的变形来释放积累的应力，使得断层不易滑动，延长地震复发间隔，而较为强硬的上盘难以变形，更容易积累应力，使得断层更易滑动，缩短地震复发间隔.而在这两者的共同作用下，岩性的差异对断层行为的影响是否有变化，如何变化，是值得探讨的问题.模型I的计算结果给出了22个强硬上盘与软弱下盘情况下的断层摩擦行为与块体变形特征，其中模型I-1—16为改变上下盘杨氏模量，保持泊松比，模型I-17—22则改变泊松比，保持杨氏模量，模型I-23则是无岩性差异的情况.

表1 材料参数列表Table 1 Material parameters

（1）上下盘存在杨氏模量差异

断层面上接触点的滑动时间表明，在强硬上盘与软弱下盘共存的情况下，上下盘杨氏模量之比越大，断层滑动时间越晚，滑动过程花费时间越长（图4），且无论是对整条断层还是某一接触点来说（图5），断层滑动时间与上下盘杨氏模量比值均为正相关关系.当比值小于1.58时，最大破裂时间与18.65km处接触点的破裂时间一致，即最晚破裂点位于18.65km深处，最大破裂时间与杨氏模量比值呈线性关系，而当比值大于1.58时，最大破裂时间与18.65km处接触点的破裂时间发生了一定的偏差，且断层破裂时间与杨氏模量比值的关系转而呈现出近于指数增长的关系.

上下盘杨氏模量之比的变化并没有改变断层的破裂过程，最大滑动时间仍发生在15～19km处，该处也是应力积累最大的地方.唯一的例外就是杨氏模量比值为9.622时，具有最大滑动速度的接触点深度小于15km，即上下盘岩性差异过大的情况下，应力积累的位置能够发生改变.当杨氏模量比值小于或等于4.97时，破裂过程可在1×3／Vref时间内完成，且破裂过程呈现阶段性，即有若干接触点同时破裂，表明能够产生较大的地震，且地震发生前有若干次较大的小震.而比值小于等于1.58时，破裂过程曲线整体差别不大，仅最晚破裂时间仍有少许的差异.

图6 下盘距断层面4km处接触点的累积垂向位移.图中右侧数字前面为上下盘杨氏模量之比，4km表示距断层面距离Fig.6 Accumulated vertical displacements of the nodes 4km away from the fault plane in the footwall using different ratio of the Young′s modulus

靠近断层面的地方变形较大，而远离断层面的地方变形较小，但由于断层面是倾斜的，横坐标相同的节点的垂向位移并不能用来相互对比，因此本文采用距离断层面垂向距离为4km、12km节点的、在破裂过程结束后积累的垂向位移进行对比（如图3所示）.

在距断层4km处（图6），累积垂向位移在深度5km以上为最大值，深度5～10km处位移量迅速减小，深度10km以下则缓慢减小，近地表处（深度＜5km）的垂向位移曲线斜率较大，而在10km以下垂向位移斜率减小，表明浅部的变形较大，且从地表向地下迅速减小，这与前山断裂带的向下延伸范围一致.累积垂向位移不仅随着深度的增大而减小，也随着上下盘杨氏模量比值的减小而减小，两者具有正相关关系，且杨氏模量比值越大，浅部地层与深部地层的变形差异越大.在距断层12km处（图7），累计垂向位移与上下盘杨氏模量比值仍然保持了正相关关系，垂向位移在近地表处较小，以近乎线性的方式随着深度的增加逐渐降低，在深度10km仍然是垂向位移曲线斜率变化的转折点.

图7 下盘距断层面12km处接触点的累积垂向位移，右侧数字含义同图6Fig.7 Accumulated vertical displacements of the nodes 12km away from the fault plane in the footwall using different ratio of the Young′s modulus

当杨氏模量比值为9.62时，其垂向位移曲线在深度10～18km处保持了较高的值，这表明下盘在深度10～18km发生了较大的变形，这与断层面上接触点的破裂过程（图5）一致，即下盘在深度10～15km的通过地层的较大变形吸收了应力，导致应力累积中心由深部约19km处向浅部转移，并积累更多的应力，花费更多的时间，最后在15km附近的位置产生大破裂.与断层破裂时间一致的是，当杨氏模量比值小于3.62时，断层下盘的地层变形开始减弱.而当杨氏模量比值小于1.58时，下盘的变形曲线越发趋于和缓.从图6和图7可以看出，杨氏模量比值为1.58和1.06的两条变形曲线的形态基本类似，表明此时断层下盘的变形也是类似的.

（2）上下盘存在泊松比差异

在上下盘杨氏模量相同，泊松比有差异的情况下，断层破裂时间与上下盘泊松比差异并没有明显的线性关系，泊松比差值在0.01～0.11的区间时，断层破裂时间仅差距约9.82%，而在泊松比差值为0.01～0.09区间，断层破裂时间仅相差3.86%.在泊松比差值为0.05和0.11时，最大破裂时间与18.65km深处接触点的破裂时间发生了偏差（图8）.断层面上各接触点的相对破裂时间（图9）表明上下盘的泊松比差值变化并没有对断层的破裂过程产生较大影响，各计算结果显示断层破裂过程所花费时间不超过0.67×3／Vref，接近于断裂两盘岩性无差异的情况下断层破裂过程花费的时间（0.47×3／Vref），远小于杨氏模量有差异的情况下的断层破裂过程花费时间（（0.54～4480）×3／Vref），表明泊松比对断层行为的影响远远小于杨氏模量的影响.下盘地层的垂直累积位移（图10）同样显示上下盘泊松比差异为断层的垂向变形影响不大.

与其他计算结果不同的是，破裂过程大部分呈现阶段性，这种阶段性的破裂主要发生在深度3～19km处的接触点上，而在19km以下有部分接触点破裂时间较晚，但最晚破裂的接触点深度仍为12～19km.在断层下盘距断层4km处节点的水平累计位移分布（图11）可以看出，水平累积位移与泊松比差值没有线性关系，但在3km和19km处，可见水平累积位移均发生了较小的不连续.这表明泊松比之间的差异导致的横向变形的差异主要体现在3km处和19km处，3km处体现了无重力作用下地表自由面的影响范围，而19km则为弧形断面和平直断面的连接处，这表明在泊松比影响下平直断面与弧形断面（图3）的相连处成为了较快滑动的地点.

（3）岩性差异的影响

图12 单条断层（模型I）与双断坡（模型II）滑动时间对比（加载时间为相对滑动时间×3／Vref）Fig.12 Relative slip time of contact nodes of the single fault（Model I）and the two thrust faults（Model II）（Loading time=relative slip time×3／Vref）

模型I-1—22的计算结果分别显示了单条断层行为与杨氏模量比值、泊松比差值之间的关系，从断层破裂时间、破裂过程以及下盘的变形来看，杨氏模量比值1.58是断层行为特征的分界点，当比值小于或等于1.58时，断层破裂过程迅速，且主震释放能量大，而下盘的变形量较小，而当比值大于1.58时，则会造成破裂过程长，能量缓慢均匀地释放，且下盘变形量过大的结果，这与汶川地震下盘变形量小，震级大的特征不符.上下盘泊松比差值对断层行为的影响不是很大，但当差值为0.07时，断层破裂的阶段性更强，12～19km深度的接触点能够同时破裂，可以造成更大的地震.

由此，本文将选取上下盘杨氏模量比值2.77，泊松比差值0.07这两个关键参数的模型（表1，模型I-24）与无岩性差异的模型（表1，模型I-23）进行对比.其中，岩性有差异的模型I-24是以岩性无差异模型I-23为基础，分别增加上盘的刚度，减弱下盘的刚度，达到所需的材料参数（表1）.计算结果显示（图12），在单条断层的情况下，上下盘的岩性差异推迟了断层破裂时间，延长了断层破裂过程，但仍表现出鲜明的阶段式破裂的特征，且最晚破裂点在15～17km范围内.

张勇[20］等利用全球地震台网（Global Seismographic Network，简写为GSN）的长周期数字地震资料反演的汶川大地震的震源机制和动态破裂过程显示破裂开始于汶川县的映秀镇地面下方约15km处，王卫民（2008）[21］反演的汶川地震破裂过程显示断层在深度15.5km处达到最大错动，USGS则指出汶川地震震源深度为14km或19km，这与我们模拟结果一致.由此认为模型I-24的材料参数是可行的.

5.2 岩性差异对双断坡构造的影响

基于上述计算分析结果，我们将模型I-24的材料参数用于双断坡构造的计算中，对比双断坡模型无岩性差异（模型II-1）和有岩性差异（模型II-2）的断层破裂过程，以及垂直于地表方向的累积位移滑动时间.

北川—映秀断裂的破裂过程（图12）结果显示，在上下盘岩性一致的情况下，双断坡构造推迟了主断层的滑动时间，延长了破裂过程，破裂过程曲线类似于单断层在上下盘杨氏模量比值较高的情况下的断层行为特征，推测这是由于主断层前缘次级断层能够造成应力分解，通过断层破裂来释放主断层下盘累积的能量，在上下盘岩性一致的情况下，其存在相当于减小了下盘的刚度.在上下盘存在岩性差异的情况下，双断坡构造中主断层12～20km深处的接触点破裂时间也呈现出阶段性破裂特征，最晚破裂时间与无岩性差异的情况下的双断坡最晚破裂时间，以及有岩性差异情况下的单条断层最晚破裂时间相近，为无岩性差异条件下单条断层最晚破裂时间的1.18～1.27倍，也就是说在主断层具有强硬上盘的情况下，下盘变形量更大，位于下盘的次级断裂的破裂与变形对主断层的断层行为造成的影响是有限的，强硬上盘对主断层的破裂起了决定性作用.但与单断层模型不同的是，在深度小于12km的接触点，其破裂过程较为连续，而非阶段性破裂，表明这些接触点逐次破裂，释放的能量较小，可能是受到了主断层前缘次级断层的影响.

双断坡构造对主断层下盘的变形具有显著影响，距主断层面4km处和12km处节点的累积垂向位移（图13）显示，无论有无岩性差异，双断坡构造主断层下盘的变形在近断层处不再具有深浅部变形差异，而表现为自地表向深处逐渐减小的特征，在距主断层面12km处，即次级断层面的下盘，仍然具有较小的深浅部变形差异，尤其是在有岩性差异的情况下，在近地表的差异较大，可能是近地表的自由面影响.

主断层破裂时累积垂向位移分布及地表垂向位移曲线（图14），即地表在这一个地震复发间隔内的隆升总量（不计震后弹性回跳）也显示有较大的变形发生在主断层和前缘的次级断层之间.虽然在有岩性差异和没有岩性差异的情况下主断层的最晚破裂时间近似，但是在有岩性差异的条件下，造成的地表累积位移更大，主破裂与次级破裂的位移量差距更大，下盘变形更集中于断层周边.

图13 单条断层（模型I）与双断坡（模型II）下盘距断层4km和12km处的垂向变形量Fig.13 Accumulated vertical displacements respectively 4km and 12km away from the fault plane in the footwall of the single fault（Model I）and the two thrust faults（Model II）

图14 双断坡构造发生滑动时的累计垂向位移量（单位为mm）Fig.14 Accumulated vertical displacement（mm）of two thrust faults（Model II）when the fault slipped

由此可见，双断坡构造在下盘浅部能够产生滑移分解，10km以上的较大变形转化为次级断层的破裂，而北川—映秀断裂上下盘巨大的岩性差异是促进双断坡构造发育的一个重要因素，也是汶川地震多条地表破裂带发育的原因之一.这与陈桂华等（2009）[7］通过利用多种地质地貌标志测绘分析得到了汶川地震发震断裂的近地表三维同震滑移矢量结果一致，也与InSAR图像显示的汶川地表破裂带前缘的同震位移量分布特征吻合[32］.

6 结论和讨论

通过基于R-minimum的有限元算法对上地壳弹性断层活动进行非线性摩擦接触模拟，我们取得了对一个地震复发间隔内单条断层的破裂过程、下盘变形强度随着岩性差异增大的变化，以及双断坡构造在强硬上盘与软弱下盘共存情况下的断层行为与块体变形的计算结果，取得了以下结论：

（1）上下盘杨氏模量之比的变化不仅能够影响单条断层的破裂时间，也能影响断层下盘的变形范围与变形量.其中，当上下盘杨氏模量之比大于1.58时，断层破裂时间随着比值的增大呈近于指数增长，破裂过程相应的延长，能量释放较为缓慢，地震震级较小，而当杨氏模量之比小于1.58时，断层破裂时间则随着杨氏模量比值的增大则呈线性增长，断层破裂过程呈现出鲜明的阶段性，即能量能够集中释放，导致较大的地震的产生.在大部分情况下，最晚破裂的接触点位于15～19km处，这与地震破裂过程反演结果近于一致[20-21］，但当杨氏模量比值为9.622时，具有最大滑动速度的接触点深度小于15km，即上下盘岩性差异过大的情况下，应力积累的位置能够发生改变.因此，北川—映秀断裂上盘的彭灌杂岩与下盘含煤层的沉积地层的共同作用可延迟断层的破裂时间，延长破裂过程，但是上下盘的杨氏模量比值应该不会太大.

（2）上下盘泊松比的差异并没有对断层行为和块体变形产生较大影响，断层破裂时间与泊松比差值之间没有明显的相关性，且破裂过程仍然保持了阶段性破裂的特征，最晚破裂的接触点位于12～19km处.

（3）单条断层的下盘的变形量随着上下盘杨氏模量比值的增大而增大，下盘地层的变形在深度5km以上较大，深度5～10km则迅速减小，深度10km以下则缓慢减小，表明单条断层能够造成浅部10km以上地层的强烈变形，杨氏模量比值越大，浅部地层与深部地层的变形差异越大.地球物理勘探显示，龙门山断裂带为前展式的逆冲断裂带，位于北川—映秀断裂前缘的灌县—江油断裂带向下延伸的深度不超过10km[9］.这与我们的模拟结果一致，即北川—映秀断裂的活动造成了断层下盘10km以上和10km以下地层变形的差异，而断裂两盘不同的岩性增大了这种深浅部变形的差异，有利于在断层前缘产生新断层.

（4）双断坡构造能够削弱单条断层下盘深度10 km以上的大变形，将近主断层的变形转化为次级断层的应力积累和断层破裂，使深部斜滑断裂在上地壳脆性域发生应变分解.

（5）在上下盘岩性一致的情况下，双断坡构造相当于减小了主断层下盘的刚度，延迟了断层的破裂时间，延长了破裂过程，而在强硬上盘和软弱下盘共存的情况下，也就是说在主断层具有强硬上盘的情况下，下盘变形量更大，位于下盘的次级断裂的破裂与变形对主断层的断层行为造成的影响是有限的，强硬上盘对主断层的破裂起了决定性作用.但是北川—映秀断裂上下盘巨大的岩性差异是促进双断坡构造发育的一个重要因素，也是汶川地震多条地表破裂带发育的原因之一.

由于单条断层深部地层和浅部地层的变形差异随着上下盘杨氏模量比值的增大而增大，双断坡构造中次级断裂相应地也受到更大变形量的推动，但由于块体越软，发育其中的断层越不容易破裂，软弱下盘中的次级破裂对主断层面的破裂时间和破裂过程的影响需要更多的工作来研究.

本文的研究着重于上地壳脆性地层中断层的弹性变形行为，且为了简化模型起见，断层两侧采用了相同的密度，不考虑重力的作用，仅以侧向挤压力推动块体造成变形.龙门山地区的实际情况可能更为复杂，希望本文关于该地区材料参数对断层行为的研究能够有助于更深入研究的开展.

致 谢 本研究是在澳大利亚昆士兰大学地球系统科学计算中心的超级计算机上利用该中心邢会林领导研发的ESyS＿Crustal软件进行的，计算中得到工作人员的大力协助，在此表示感谢.

（References）

[1]Gudmundsson A.Effects of Young′s modulus on fault displacement.Comptes Rendus Géosciences，2004，361（1）：85-92.

[2]Ben-Zion Y，Andrews D J.Properties and implications of dynamic rupture along a material interface.Bulletin of the Seismological Society of America，1998，88（4），1085-1094.

[3]Ben-Zion Y.Dynamic ruptures in recent models of earthquake faults.Journal of the Mechanics and Physics of Solids，2001，49（9）：2209-2244.

[4]Ben-Zion Y，Shi Z Q.Dynamic rupture on a material interface with spontaneous generation of plastic strain in the bulk.Earth and Planetary Science Letters，2005，236（1-2）：486-496.

[5]王萍，付碧宏，张斌等.汶川8.0级地震地表破裂带与岩性关系.地球物理学报，2009，52（1）：131-139.Wang P，Fu B H，Zhang B，et al.Relationships between surface ruptures and lithologic characteristics of the Wenchuan Ms8.0earthquake.Chinese J.Geophys.（in Chinese），2009，52（1）：131-139.

[6]Xu X W，Wen X Z，Yu G H，et al.Co-seismic reverse and oblique-slip surface faulting generated by the 2008 Mw7.9 Wenchuan earthquake，China.Geology，2009，37（6），515-518.

[7]陈桂华，徐锡伟，于贵华等.2008年汶川MS8.0地震多断裂破裂的近地表同震滑移及滑移分解.地球物理学报，2009，52（5）：1384-1394.Chen G H，Xu X W，Yu G H，et al.Co-seismic slip and slip partitioning of multi-faults during the Ms8.0 2008Wenchuan earthquake.Chinese J.Geophys.（in Chinese），2009，52（5）：1384-1394.

[8]Li Y，Allen P A，Densmore A L，et al.Evolution of the Longmen Shan foreland basin（Western Sichuan，China）during the late Triassic Indosinian Orogeny.Basin Research，2003，15：117-138.

[9]贾东，陈竹新，贾承造等.龙门山前陆褶皱冲断带构造解析与川西前陆盆地的发育.高校地质学报，2003，9（3）：402-410.Jia D，Chen Z X，Jia C Z，et al.Structural features of the Longmen Shan Fold and thrust belt and development of the western Sichuan foreland basin，central China.Geological Journal of China Universities （in Chinese），2003，9（3）：402-410.

[10]林茂炳，马永旺.论龙门山彭灌杂岩体的构造属性.成都理工学院学报，1995，22（1）：42-46.Lin M B，Ma Y W.A discussion on the tectonic characteristics of Peng-Guan Complex in Longmen Mountains.Journal of Chengdu Institute of Technology （in Chinese），1995，22（1）：42-46.

[11]Jia D，Li Y Q，Lin A M，et al.Structural model of 2008 Mw7.9Wenchuan earthquake in the rejuvenated Longmen Shan thrust belt，China.Tectonophysics，2010，491（1-4）：174-184.

[12]Li Q S，Gao R，Wang H Y，et al.Deep background of Wenchuan Earthquake and the upper crust structure beneath the Longmen Shan and adjacent areas.Acta Geologica Sinica，2009，83（4）：733-739.

[13]Xing H L， Makinouchi A.Finite-element modeling of multibody contact and its application to active faults.Concurrency and Computation：Practice and Experience，2002，14（6-7）：431-450.

[14]Xing H L， Makinouchi A.Finite element modelling of frictional instability between deformable rocks.Int.J.Numer.Anal.Meth.Geomech，2003，27（12）：1005-1025.

[15]Xing H L， Makinouchi A， MORA P.Finite element modeling of interacting fault systems.Physics of the Earth and Planetary Interiors，2007，163（1-4）：106-121.

[16]Dieterich J H.Modeling of rock friction.1.Experimental results and constitutive equations.J.Geophys.Res.，1979，84（B5）：2161-2168.

[17]Ruina A.Slip instability and state variable friction laws.J.Grophys.Res.，1983，88（B12）：10359-10370.

[18]Scholz C H.Earthquakes and friction laws.Nature，1998，391（6662）：37-42.

[19]刘超，张勇，许力生等.一种矩张量反演新方法及其对2008年汶川Ms8.0地震序列的应用.地震学报，2008，30（4）：329-339.Liu C，Zhang Y，Xu L S，et al.A new technique for moment tensor inversion with applications to the 2008Wenchuan MS8.0earthquake sequence.Acta Seismologica Sinica（in Chinese），2008，30（4）：329-339.

[20]张勇，许力生，陈运泰.2008年汶川大地震震源机制的时空变化.地球物理学报，2009，52（2）：379-389.Zhang Y，Xu L S，Chen Y T.Spatio-temporal variation of the source mechanism of the 2008great Wenchuan earthquake.Chinese J.Gephys.（in Chinese），2009，52（2）：379-389.

[21]王卫民，赵连锋，李娟等.四川汶川8.0级地震震源过程.地球物理学报，2008，51（5）：1403-1410.Wang W M，Zhao L F，Li J，et al.Rupture process of the Ms8.0Wenchuan earthquake of Sichuan，China.Chinese J.Geophys.（in Chinese），2008，51（5）：1403-1410.

[22]刘树根，李智武，曹俊兴等.龙门山陆内复合造山带的四维结构构造特征.地质科学，2009，44（4）：1151-1180.Liu S G，Li Z W，Cao J X.et al.4-D textural and structural characteristics of Longmen intracontinental composite orogenic belt，southwest China.Chinese Journal of Geology （in Chinese），2009，44（4）：1151-1180.

[23]刘静，曾令森，丁林等.青藏高原东南缘构造地貌、活动构造和下地壳流动假说.地质科学，2009，44（4）：1227-1255.Liu J，Zeng L S，Ding L，et al.Tectonic geomorphology，active tectonics and lower crustal channel flow hypothesis of the southeastern Tibetan Plateau.Chinese Journal of Geology （in Chinese），2009，44（4）：1227-1255.

[24]张培震，徐锡伟，闻学泽等.2008年汶川8.0级地震发震断裂的滑动速率、复发周期和构造成因.地球物理学报，2008，51（4）：1066-1073.Zhang P Z，Xu X W，Wen X Z，et al.Slip rates and recurrence intervals of the Longmenshan active fault zone，and tectonic implications for the mechanism of the May 12 2008，Wenchuan earthquake，Sichuan，China.Chinese J.Geophys.（in Chinese），2008，51（4）：1066-1073.

[25]李勇，徐公达，周荣军等.龙门山均衡重力异常及其对青藏高原东缘山脉地壳隆升的约束.地质通报，2005，24（12）：1162-1168.Li Y，Xu G D，Zhou R J，et al.Isostatic gravity anomalies in the Longmen Mountains and their constraints on the crustal uplift below the mountains on the eastern margin of the Qinghai—Tibet Plateau.Geological Bulletin of China （in Chinese），2005，24（12）：1162-1168.

[26]王谦身，滕吉文，张永谦等.四川中西部地区地壳结构与重力均衡.地球物理学报，2009，52（2）：579-583.Wang Q S，Teng J W，Zhang Y Q，et al.The crustal structure and gravity isostasy in the middle western Sichuan area.Chinese J.Geophys.（in Chinese），2009，52（2）：579-583.

[27]刘启元，李昱，陈九辉等.汶川Ms8.0地震：地壳上地幔S波速度结构的初步研究.地球物理学报，2009，52（2）：309-319.Liu Q Y，Li Y，Chen J H，et al.Wenchuan Ms8.0 earthquake：preliminary study of the S-wave velocity structure of the crust and upper mantle.Chinese J.Geophys.（in Chinese），2009，52（2）：309-319.

[28]Ran Y K，Chen L C，Chen J，et al.Paleoseismic evidence and repeat time of large earthquakes at three sites along the Longmenshan fault zone.Tectonophysics，2010，491（1-4）：141-153.

[29]朱守彪，张培震.2008年汶川MS8.0地震发生过程的动力学机制研究.地球物理学报，2009，52（2）：418-427.Zhu S B，Zhang P Z.A study on the dynamical mechanisms of the Wenchuan MS8.0earthquake，2008.Chinese J.Geophys.（in Chinese），2009，52（2）：418-427.

[30]叶金汉.岩石力学参数手册.北京：水利水电出版社，1991.Ye J H.Manual of Rock Mass Parameters （in Chinese）.Beijing：China Water Power Press，1991.

[31]张培震，闻学泽，徐锡伟等.2008年汶川8.0级特大地震孕育和发生的多单元组合模式.科学通报，2009，54（7）：944-953.Zhang P Z，Wen X Z，Xu X W，et al.Tectonic model of the great Wenchuan earthquake of May 12，2008，Sichuan，China.Chinese Science Bulletin （in Chinese），2009，54（7）：944-953.

[32]Hao K X，Si H J，Fujiwara H，et al.Coseismic surfaceruptures and crustal deformations of the 2008Wenchuan earthquake Mw7.9，China.Geophys.Res.Lett.，2009，36：L11303，doi：10.1029／2009GL037971.