2008年5月12日汶川8.0级地震发生在巴颜喀拉与华南地块交界的北东向龙门山断裂带上.地震形成了长200多公里､南西-北东向延伸的地表破裂带.相对四川盆地，龙门山沿该地表破裂带既有逆冲运动，又有北东向走滑运动^[1].多位学者利用多种资料初步判定汶川地震的发震构造是龙门山中央断裂中北段和前山断裂中段，同时基底滑脱带在地震孕育过程中起到重要作用^[2-3].汶川地震是有历史记录以来首例发生在大陆内部低滑动速率断裂上的高角度逆冲型8级强震^[4]，其强度远远超出了龙门山断裂带上的最大历史地震，因此仅根据数百年至一两千年的历史地震记载远不足以正确评估一条较低速率的大型活动断裂带的潜在地震危险性^[4-5].

根据现代大地测量资料，反演断层深浅部的闭锁程度和滑动亏损分布特征是判断断裂带未来一段时期内地震危险性的重要手段.发生汶川地震的龙门山断裂带及其附近地区积累了丰富的GPS 观测结果，为研究该断裂带震前变形提供了可能.本文利用1999-2007 期GPS 水平速度场结果，采用Defnode负位错反演程序估算了龙门山断裂在汶川地震前的断层闭锁程度､滑动亏损速率分布等，结合龙门山断裂带附近地表水平应变率场结果，分析了震前断层深浅部地壳变形特征.另外，根据不同学者由数字地震波形资料和同震GPS 位移数据反演给出的最大同震滑动量的可能范围，讨论了汶川地震可能的复发间隔.

McCaffrey^[6]给出了Defnode负位错反演程序，该程序假定块体内部点的运动为块体旋转与块体边界由于断层闭锁产生的滑动亏损而引起的地表弹性变形之和.程序可以综合利用GPS 数据､地质数据等结果，基于Okada模型^[7-8]，采用网格搜索和模拟退火方法反演出块体旋转运动的欧拉极和块体边界断层的闭锁程度(或滑动亏损).McCaffrey等^{[6, 9-11]}运用此程序计算了苏门答腊岛､哥斯达黎加和美国西海岸等地区的断层闭锁和块体旋转与变形；Wallace等^[12-14]计算了新西兰南岛､北岛和巴布亚新几内亚等地区的俯冲带断层闭锁与构造块体旋转；还有其他一些学者对海洋块体俯冲带的断层闭锁与块体旋转进行了计算，但是基于此程序计算大陆块体内部高角度逆冲区域的断层闭锁与滑动亏损分布还未见相关文章发表.

为了约束块体旋转和断层闭锁程度，利用GPS水平速度场数据进行反演，所用公式^[6]如下:

(1)

其中X为GPS测站的位置，V_i(X)为测站X的速度，B为块体数目，Δ_b为b块体模型区域的子集(如果点X在块体b范围内H=1，否则H=0)，i为i方向速度单位矢量(x或者y)，_RΩ_b为块体b相对于参考框架的欧拉极，_hΩ_f=_hΩ_R-_fΩ_R为断层下盘f相对于上盘h的欧拉极，F为断层数量，N_k为定义断层k的节点数，X_nk为断层k上节点n的位置，Φ_nk为断层k上节点n的闭锁程度，j为断层面上j方向的单位矢量(垂直等深线方向或者沿走向方向)，G_ij(X，X_nk)为响应函数，它表示断层面上节点X_nk处j方向上的单位滑动速率引起的地表位置X处i方向的速度分量.公式右端前半部分代表块体旋转，用块体旋转极求出块体上所有点的刚体运动速率；后半部分代表断层闭锁产生的影响，用欧拉极(_hΩ_f)求出下盘与上盘的相对滑动速率，该值与断层闭锁程度乘积可计算得到滑动亏损速率.

根据地质结果和现今大地测量结果，以龙门山断裂带和鲜水河断裂带为界将本文研究区域划分为3个地块-巴东地块､华南地块和川西北地块(图 1).在反演时，共用了1999-2007 期的98 个GPS水平速度场数据，其中巴东块体32个､华南块体44个､川西北块体22个，速度场采用GAMIT/GLOBK^[15-16]和QOCA^[17]软件解算得到(图 1).该速度场结果共利用了4 期观测数据(1999 期､2001 期､2004 期和2007期)，平均误差为1.7mm/a.

图 1

Fig. 1

图 1 研究区GPS站水平运动速度场(相对于华南块体)与活动断裂分布 Fig. 1 The GPS velocity field of horizontal movement(related to South China Block) and active faults in the study region

反演过程中，沿龙门山断裂走向方向节点之间的距离为50~55km，每条等深线上有9个节点；节点深度方向从地表起依次为0.1km､6km､12km､16km､21km 和24km，垂直等深线方向有6排节点；根据文献的研究结果^{[1-3, 18-21]}，每两排节点之间的断层倾角从地表至深部依次为55°､55°､20°､7°和7°，断裂总长度为420km，断层面总宽度为92km(图 2).利用Defnode程序可计算得到每个节点处断层闭锁程度，进而通过双线性插值方法计算相邻节点之间断层网格(沿走向方向长4km，深度方向宽4km)区域的闭锁程度.对断层节点Φ 值所加约束如下:Φ值从地表沿垂直等深线向下单调递减(0.1km与6km 处的Φ 值相等)，龙门山断裂北段､中段和南段各自的节点在同一等深线上的Φ 值相等，并且在24 km 处整条断裂的Φ 值均为0(完全蠕滑).

图 2

Fig. 2

图 2 断层模型与节点设置 (a)断层面结构与节点位置；(b)断层剖面结构. Fig. 2 Fault model and node setting (a) Fault surface structure and node setting；(b) Fault profile structure.

在逆冲断层中，Φ 值可以随深度发生变化，变化规律可以是单调递减函数，也可以是其它的一些深度的函数.McCaffrey^[6]曾做实验证明，如果随深度增加断层滑动亏损速率变大，那么位错模型就可能得到错误的地表变形，另外Wang等^[22]根据逆冲断层的热性能也得到过类似结论.因此本文假定Φ 值随深度增大而减小，在达到一定深度后为0.

反演所求参数拟合的好坏可由公式(2)表征:

(2)

其中n为所有观测数据的数量，dof为自由度(所有观测数据数量-自由参数数量)，ri为观测数据残差，σi为数据标准差，f为数据误差权重因子，它可以平衡各种数据(如GPS 水平速度场和GPS 垂直速度场)的影响，一般GPS 水平速度场数据误差的f取值范围为1~5^[23].反演时通过不断改变其大小进行逐步择优，寻求最恰当的f值，进而保证x_n²≈1，此时模型能够准确地解释观测数据.反演过程中，通过使用网格搜索和模拟退火方法同时反演块体旋转运动的欧拉极和块体边界断层的闭锁程度，并使x_n²最小.

通过大量试算，我们得到了最优模型，其中GPS水平速度场数据误差权重因子f取为1.25，2=1.008 (观测值个数为196，自由度为166)，根据前面对f和2 的分析表明模型是有效的.输入数据的残差分布可以反映模型拟合程度的好坏，我们在图 3中给出了GPS 速度残差分布.结果显示，只有巴东块体西端远离断层的几个测站速度残差值稍大；而其它测站特别是龙门山断裂附近测站的速度残差值很小，基本处在误差范围之内，表明模型拟合较好.

图 3

Fig. 3

图 3 模型速度残差 Fig. 3 Residual velocities of the model

表 1给出了计算得到的巴东块体和川西北块体相对于华南块体的欧拉极.结果表明巴东和川西北块体均相对于华南块体做顺时针旋转，其中川西北块体的旋转速率明显高于巴东块体；川西北块体相对于华南块体顺时针旋转运动的欧拉极距离川西北块体很近，而巴东块体相对于华南块体顺时针旋转运动的欧拉极距离巴东块体较远.

表 1(Table 1)

表 1 块体欧拉矢量 Table 1 Residual velocities of the model

表 1 块体欧拉矢量 Table 1 Residual velocities of the model

因为基准参数不包含地壳相对运动与变形信息，所以基准转换过程保持了地壳相对运动的不变性.这表明尽管同一地区相同资料来源给出的不同基准速度场显示出不同的图像，但其中实际包含的地壳变形信息是完全一致的^[24].因此参考基准的不同并不影响表 1结果.

根据前面的最优模型结果，我们得到了汶川地震前龙门山断裂的闭锁程度，如图 4 所示，其中起始破裂深度取19km^[25].结果表明，震前龙门山断裂中段和北段处于完全闭锁状态，而且闭锁很深，达到21km(闭锁比例0.99)左右；南段只有地表以下12km 闭锁程度较高(闭锁比例0.99)，在12~16km处闭锁比例约为0.83，在16~21km 处约为0.75，较中北段同一深度的闭锁程度稍弱.在21~24km 处整条断裂均逐步由闭锁转变为完全蠕滑.

图 4

Fig. 4

图 4 震前龙门山断裂闭锁程度 (a)断层闭锁地面投影；(b)断层闭锁分布. Fig. 4 Fault locking of the Longmenshan fault before the Wenchuan earthquake (a) The projection of fault locking；(b) Fault locking distribution.

断层间的相对运动由断层两侧块体的欧拉极决定，断层滑动亏损速率为闭锁程度与两个块体间断层上的相对运动矢量的乘积.图 5 给出了龙门山断裂的滑动亏损分布，结果表明震前龙门山断裂中北段垂直断层方向的挤压滑动亏损速率与平行断层方向的右旋滑动亏损速率均比南段大.图 5a表明龙门山断裂中北段从地表至21km处垂直断层方向滑动亏损速率约为2.2mm/a(沿走向方向变化范围为1.8~2.5mm/a，最北端最大，往南西逐渐减小)；南段从地表至12km 处滑动亏损速率约为1.4mm/a，12~16km 处滑动亏损速率约为1.2mm/a，16~21km处滑动亏损速率约为1.1mm/a.图 5b表明龙门山断裂中北段从地表至21km 处平行断层方向滑动亏损速率约为4.6 mm/a；南段从地表至12km 处滑动亏损速率约为4.6mm/a，12~16km 处滑动亏损速率约为3.8mm/a，16~21km 处滑动亏损速率约为3.5 mm/a.在21~24km 处整条断裂均逐步转变为完全蠕滑.

图 5

Fig. 5

图 5 震前龙门山断裂滑动亏损分布 Fig. 5 Fault slip deficit rate of the Longmenshan fault before the Wenchuan earthquake

GPS揭示的汶川地震大区域水平同震位移场最显著的特征是以龙门山断裂带的映秀-北川断裂为核心的相向运动和地壳缩短，兼有右旋走滑^[26].逆冲的位移量与右旋走滑位移量之比在破裂带西北侧的龙门山地区约为2∶1.大地震所释放的弹性应变是在震前长期积累形成的，因此震前龙门山断裂西北侧的垂直断层的应变积累与平行断层的应变积累也应该约为2∶1.但是GPS 资料揭示，震前龙门山断裂中北段西侧巴颜喀拉东部500km 尺度变形是以右旋剪切变形为主，垂直断层的地壳缩短量与平行断层的右旋扭动量之比约为1∶3，因此震前8年间弹性应变积累与同震能量释放是不对应的^[27].本文根据GPS资料反演所得震前龙门山断裂中北段垂直断层方向滑动亏损速率约为2.2 mm/a，平行断层方向滑动亏损速率约为4.6 mm/a，二者比例约为1∶2，相对来说更接近震前GPS资料反映的垂直断层的地壳缩短量与平行断层的右旋扭动量之比1∶3，这个比例与同震能量释放比例之间的不对应关系是一个值得深思的问题.文献^[27]提出了一种可能的分析解释，即大地震的弹性应变积累是一个极其漫长的过程，很可能在此次地震孕育的最后阶段龙门山断裂带附近孕震区域的地壳弹性变形已经趋于极限，因此处于难以再发生显著变形的相持阶段.这个解释与本文所得结果---震前龙门山断裂中北段处于完全闭锁状态，而且闭锁很深(达到21km)是一致的.

图 6为龙门山断裂中段和南段GPS 水平观测值与模型值的剖面拟合结果(剖面位置见图 1)，其中剖面长700km，中段宽180km，南段宽120km.图中灰色线为平行断层运动的模型值，斜率为正表示左旋，反之为右旋；黑色线为垂直断层运动的模型值，斜率为正表示拉张，反之为挤压.图 6结果表明，只有巴东块体西端几个GPS 测站的速度与模型值差别稍大(3~4 mm/a)，可能源于这些测站距离断层较远，含有较大的模型误差或块体内其它断层上存在变形，影响了拟合效果；其它GPS 测站的速度与模型值差值基本保持在2mm/a以内，拟合较好.另外，各测站垂直断层方向的速度拟合结果好于平行断层方向上的.

图 6

Fig. 6

图 6 横跨龙门山断裂中段和南段的GPS站速度剖面和拟合结果(相对于华南地块) (a)龙门山断裂中段结果；(b)龙门山断裂南段结果. Fig. 6 GPS velocity profiles across the middle segment and the southern segment of the Longmenshan fault and the fitting results (related to South China Block) (a) Results of the middle segment of the Longmenshan fault；(b) Results of the southern segment of the Longmenshan fault.

龙门山断裂中段，GPS测站平行与垂直断裂方向上的速度分量均逐渐下降并接近四川盆地的速度分量；而龙门山断裂南段，靠近断裂两侧的GPS测站速度显示出一定的杂乱性，没有明显的规律.龙门山断裂中段的模型速度受断层闭锁的影响比南段稍大一些，不仅表现为断层附近速度分量受闭锁的影响明显，而且其受影响范围也大一些.上述现象表明龙门山断裂中段应变积累程度高，存在显著的闭锁作用，且闭锁较深；而南段闭锁程度较北段稍弱，依然可以通过断裂近场区域的变形来积累应变和调节､吸收能量.

由于GPS观测主要反映了地壳上部的变形特征，因此汶川地震前龙门山断裂带及其附近地区的地壳变形特征也应该可以反映该地震的孕育特征.选取图 1中(99°E-107°E；25°N-34°N)GPS 水平速度场数据(共253 个测站)，采用最小二乘配置方法计算得到了该区域的应变率分布(图 7).关于最小二乘配置应变率计算方法可参考文献^{[24, 28-29]}.图 7a结果表明发生汶川地震的龙门山断裂带附近的变形具有明显的分区特征，汶川地震震源区南西侧断裂的应变积累速率明显大于北东侧，整条断裂带表现为挤压兼右旋剪切变形特征，其中南段 的挤压变形幅度大于中北段.震源北东侧断层的北西侧变形幅度表现为从远离断裂带较大到靠近断裂带逐渐减小的特征，表明其在震前较长时期已经积累了较高的应变能，致使有GPS观测以来无法观测到断裂带显著变形，该特征更利于发生破裂滑动.汶川地震的破裂过程和余震分布特征(沿主震发生地点向北东向破裂､余震主要分布于主震的北东方向断层上)也证明了这一点^[30].图 7b结果表明受青藏高原东向运动影响，在汶川地震震源区的南西方向形成了显著的东西向挤压应变集中区；以震源位置以北的龙门山断裂带为界，断裂带的北西侧挤压变形明显，而断裂带的南东侧则属于东西向拉张变形区，该结果表明龙门山断裂带作为东西向变形的分界带，具备长期挤压应变积累背景，并且闭锁程度很高.

图 7

Fig. 7

图 7 1999-2007期GPS高分辨率应变率场结果 (a)主应变率分布；(b)东西向应变率分布. Fig. 7 Results of GPS strain rate in high resolution for1999-2007 (a) Distribution of principle strain rate；(b) Distribution of strain rate in the east-west direction.

(1) 断层模型的不同参数化讨论

反演过程中，断层模型的不同参数化对结果有一定影响，因此对于断层模型参数的设置需要尽可能接近断层实际情况.Wang等^[18]已经利用同震位错反演对断层结构做了大量研究，提出的断层模型也对资料做出了很好的解释，因此采用其断层模型参数进行震前负位错反演应当是一个合理的选择.

(2) 反演结果不确定性讨论

一般GPS水平速度场数据误差的权重因子f取值范围为1~5^[23].反演时通过不断改变其大小进行逐步择优，寻求最恰当的f值.本文通过大量试算，获取了f值的最优解(f=1.25)，进而得到了最优模型.我们对f取值1.15和1.35进行了反演计算，结果表明断层的闭锁程度与最优模型结果基本一致，说明了反演结果具有一定的稳定性.另外，最优模型中GPS数据的残差分布基本处在误差范围之内，也表明模型拟合较好.

(3)汶川地震复发间隔讨论

龙门山断裂带在有历史记录以来没有发生过7级以上地震，因此，对类似于汶川8.0级这样的地震在龙门山断裂带的复发间隔的认识是我们探索地震发生规律必须思考的一个问题.地震复发间隔时间T可以表示为T=S/V，其中S为断层同震滑动，V为断层滑动亏损(或者断层长期滑动)速率^[31].同震破裂释放的弹性应变是在震前长期积累形成的，汶川大地震以逆冲运动为主，因此震前长时间尺度的应变积累也应该以逆冲积累为主.虽然震前8 年间GPS观测到的平行断层方向右旋剪切变形速率高，而垂直断层方向挤压缩短变形速率低，但这可能并不代表长期(整个震间期)应变积累的情况，而计算地震复发间隔更应当考虑同震破裂的结果.另外，如果单考虑平行断层方向的右旋滑动亏损速率，根据同震位移推算的地震复发间隔偏短(600~1200a)，与具有1700年记录以来龙门山断裂带无强震发生的观测不符合.因此，我们主要参考垂直断层方向的挤压滑动亏损速率计算了地震复发间隔.Wang等^[18]､王卫民等^[21]､Shen 等^[31]､赵翠萍等^[32]､张勇等^[33]､Nishimuru等^[34]､王敏^[35]和Wang等^[36]曾对汶川地震的最大同震滑动量做过相关研究，将他们给出的最大同震滑动量分别除以断层的滑动亏损速率2.5mm/a(同震逆冲的位移量与右旋走滑位移量之比为2∶1，因此由垂直断层方向滑动亏损速率2.2mm/a可以求得整个断层的滑动亏损速率)，我们推算了龙门山地区8.0级地震可能的复发间隔为2500~5000a(表 2)，这与张培震等^[1]､李勇等^[37]､Burchfiel等^[38]所得结果基本一致.

表 2(Table 2)

表 2 龙门山断裂带8.0级地震复发间隔估计 Table 2 8.0 earthquake recurrence estimates of the Longmenshan fault zone

表 2 龙门山断裂带8.0级地震复发间隔估计 Table 2 8.0 earthquake recurrence estimates of the Longmenshan fault zone

(1) 巴东块体的运动包括它相对于华南块体的块体旋转与块体边界由于断层闭锁产生的滑动亏损而引起的地表弹性变形之和.震前龙门山断裂中北段处于完全闭锁状态，闭锁深度达到21km(闭锁比例0.99)左右；南段只有地表以下12km 闭锁程度较高(闭锁比例0.99)，在12~16km 处闭锁比例约为0.83，在16~21km 处约为0.75，较中北段同一深度的闭锁程度稍弱.在21~24km 处整条断裂均逐步转变为完全蠕滑.上述特征表明龙门山断裂中北段的地震危险程度高于南段，该结果可以较好地解释汶川地震破裂方式及强余震发生位置，从另一角度说明本文反演结果较为可靠.

(2) 模型值剖面与GPS 速度剖面拟合较好，各测站垂直断层方向的速度拟合结果好于平行断层方向上的.龙门山断裂中段，GPS 测站平行与垂直断裂方向上的速度分量均逐渐下降并接近四川盆地的速度分量；而龙门山断裂南段两侧区域的GPS测站速度显示出一定的杂乱性，没有明显的规律.龙门山断裂中段的模型速度受断层闭锁的影响比南段稍大一些，表现为断层附近速度分量受闭锁的影响明显，且受影响范围大一些.这些现象表明龙门山断裂中段存在显著的闭锁作用，且闭锁较深，而南段闭锁程度较北段稍弱，依然可以积累应变.

(3) 最小二乘配置应变率结果显示，龙门山断裂带中北段的应变积累速率明显小于南段，该结果表明中北段已经积累了较高的应变能且闭锁程度很高，而南段应变积累速率较大，闭锁程度较北段弱.该分析结果与Defnode负位错反演结果相吻合.

(4) 龙门山断裂中北段应变积累程度很高､闭锁程度很高､闭锁深度很深，已经趋于应变积累的极限阶段，当震源区发生破裂时断裂近场区域无法再通过变形来调节和吸收能量，所以发生了破裂；龙门山断裂南段闭锁程度较北段稍弱､应变积累程度也稍弱，还没达到应变积累的极限，震源区发生破裂时依然可以通过断裂近场区域的变形来调节和吸收能量，因此地震过程中并没有发生破裂.

感谢中国地壳运动监测网络中心为本研究提供数据；感谢中国地震局地质研究所王敏研究员､中国地震局地震预测研究所邵志刚副研究员为本研究提供的帮助；感谢两位审稿专家提出的中肯意见.