引言

地震通常都伴随断层的滑移，长期构造运动产生地壳变形，当变形超过断层的承受临界值时，断层就会发生滑动破坏。断层滑移是理解震源机制和发震机理的重要参数（柳畅等，2014）。目前对断层滑移的研究大多集中于同震滑移，由于地表断层破裂数据不能反映断层的真实滑移，李志才等（2009）和王阅兵等（2012）根据GPS观测的地表同震变形数据，通过对断层同震滑移的反演计算，分析地震的发震类型和破坏特征，并据此理解震源机制。震前断层滑动是与地震的孕育和发生直接关联的物理现象。所以，震前断层滑移与同震参数的比较分析是理解和进行地震预测的一个具有明显物理背景的分析方法。

基于此，本文以汶川地震为背景，利用中国地壳运动观测中心提供的1999年至2001年、2001年至2004年、2004年至2007年共3期的153个固定站点GPS测量数据（水平精度2mm，垂直精度6mm）分析震前地表变形的变化趋势，采用遗传算法反演震前3期断层滑移参数，并将震前断层滑移参数与同震测量的38个GPS站点数据反演得到的同震断层参数进行比较，分析震前断层滑移与同震断层滑移的关联。

1.   断层位置与GPS测站布置

汶川地震震后地质考查表明（焦青等，2008；刘静等，2008；李勇等，2008），此次地震由中央断裂带、山前江油-灌县断裂带和山后茂汶-汶川断裂带等多条断层同时参与地震断层破裂过程，但中央断裂带为主要断裂带，其走滑和倾滑位移较大，山前和山后断裂带倾滑运动很小，可以忽略，仅考虑其走滑运动，因此本文将中央断裂带作为模型研究对象。断层位置与GPS测点分布如图 1所示，根据断层位置，将其拟合为一条直线，将断层视为平面，且断层面上的位移是均匀分布的。资料显示（陈运泰，2008），断层破裂长度在300km左右，通过地质考察等方式（李海兵等，2008；张培震等，2008）得到汶川地震主断裂带断层近地面倾角为60°—70°，断层在20km处断层随深度向下逐渐变缓。根据这些背景资料，本文选取断层长度L=300km，倾角θ=60°，断层深度D=20km作为模型形状参数。

图 1  断层位置与GPS测点分布图

Figure 1.  Location offault and GPS measurements

下载: 全尺寸图片 幻灯片

2.   震前3期断层滑移参数

2.1   震前断层参数计算方法

断层滑移参数的计算通常采用反演算法（徐果明，2003）。采用图 2所示的断层位错模型，图中竖直方向坐标轴的原点z=0定义在地表平面，d为震源深度。断层长、宽分别为L和W，倾角为θ。x轴平行于断层走向。U₁、U₂、U₃分别代表断层上盘相对于下盘的走滑、倾滑和张引位错分量。

图 2  断层位错模型

Figure 2.  Fault dislocation model

下载: 全尺寸图片 幻灯片

为了提高断层位移反演的精度，本文采用所有GPS测量点的3个方向测量值和模型值差的绝对值作为目标函数值，即：

式中，M为目标函数值，n为GPS测点数量，u_xⁱ、u_yⁱ、u_zⁱ为GPS测量的实际地表位移数据。u_x0、u_y0、u_z0是基于图 2的断层位错模型，由断层初始参数U₁、U₂、U₃计算得到的各GPS测点对应的地表位移值。所以，目标函数值M越接近于0，模型计算出的地表位移值u_x0、u_y0、u_z0越接近GPS实测位移值u_xⁱ、u_yⁱ、u_zⁱ，反演得到的断层滑移参数U₁、U₂、U₃就越精确。本文采用遗传算法（周明，1999；雷英杰，2014）求解断层滑移参数U₁、U₂、U₃。遗传算法是一种仿生算法，它把待求参数编制成一组二进制代码，将这组编码看作一个个体，并随机产生包含一定数量个体的种群，编制目标函数并将目标函数值经过变换得到适应度值，即评价种群中个体值与真实值的差异程度的标准。将种群中适应度值高的个体保留下来，对这些个体进行选择、交叉、变异等运算形成新的种群，并用种群中适应度高的个体代替原种群适应度低的个体，经过多次重复迭代使种群中个体值接近真实值。

在整个算法中的关键步骤是基于位错模型（Steketee，1958；Okada，1985；1992）计算断层滑动参数U₁、U₂、U₃对应的地表变形u_x0、u_y0、u_z0。根据点源位错模型，断层内一点位移对地表产生的影响可以表示为：

式中，u_i为地表某点位移，下标i=x、y、z，Δu_j为位错点位移，下标j=1、2、3，k=3，n=1、2，j、k、n上标与下标含义相同，1，2，3分别对应于U₁、U₂、U₃三个方向（如图 2所示）。δ_jk是Kronecker张量（当k=j时δ_jk=1，当k≠j时，δ_jk=0），λ、μ为拉梅参数，v_k为断层面法向分量方向余弦，$ u_{i}^{j}\left(u_{i}^{n}, u_{i}^{k} \right) $，表示在点（ξ₁，ξ₂，ξ₃）处振幅为F的点震源的j（n，k）分量在点（x₁，y₁，z₁）处产生的i分量位移。

由点源位错计算公式（2）沿断层面积分，可得有限矩形位错公式。

2.2   模型验证

下面，通过一个模拟断层来校验本文算法。已知断层滑移等参数列于表 1，其中L为断层长度，W为断层宽度，d为断层深度，θ为断层倾角，U₁、U₂、U₃为断层3个方向的滑移。据此计算地表位移分量，选取其中5组数据（表 2）进行分析。然后，基于这5组数据，利用遗传算法反演断层滑动参数，结果列于表 3。计算结果与给出的断层参数一致。为更进一步说明计算效果，图 3示出了反演计算中目标函数演化过程。可以看出，目标函数收敛迅速，较好地趋近0且并未有波动趋势，表明本文算法的稳定性和较好的收敛效果。

表 1  模拟断层参数

Table 1.  Simulation of fault parameters

断层参数	L/km	W/km	d/km	θ/°	U1/m	U2/m	U3/m
数值	9	6	7	85	0.5	0.2	0.1

下载: 导出CSV 

| 显示表格

表 2  断层正演结果

Table 2.  Forwardsolution of fault parameters

编号	x/m	y/m	ux/m	uy/m	uz/m
1	1000	2000	-0.0951	-0.0304	-0.0664
2	2000	2000	-0.0948	-0.0333	-0.0666
3	5000	3000	-0.0905	-0.0465	-0.0589
4	-2000	3000	-0.0949	-0.0192	-0.0575
5	-1000	1000	-0.0426	-0.0069	-0.0409

下载: 导出CSV 

| 显示表格

表 3  断层反演结果

Table 3.  Inversion result of fault parameters

编号	目标函数值	U1'/m	U2'/m	U3'/m
1	6.35E-04	0.5001	0.2001	0.1001
2	1.70E-05	0.5000	0.2002	0.0995
3	5.92E-05	0.5000	0.2002	0.0994
4	3.50E-07	0.5000	0.2001	0.0996
5	1.00E-03	0.5001	0.1997	0.1007
平均值	-	0.5000	0.2001	0.0999

下载: 导出CSV 

| 显示表格

图 3  目标函数的演化过程

Figure 3.  Evolution of objective function

下载: 全尺寸图片 幻灯片

2.3   震前三期断层滑移特征

将GPS测量的震前3期地表位移数据进行高斯换算和坐标旋转，把大地坐标系转换到以断层下缘线为水平轴的平面坐标系中。利用变换后的地表位移以及上述计算模型反演龙门山断层1999年到2007年震前3期滑动位移的累积量和断层同震滑移量。震前3期的时间（以年为单位）为：第一期1999.771—2001.678，时间间隔1.907年；第二期为2001.678—2004.583，时间间隔2.905年，第三期为2004.583—2007.579，时间间隔2.996年。

从震前3期的目标函数演化过程（图 4）可以看出，3期结果均收敛较好，其中3期的目标函数值分别为0.0231mm、0.0261mm、0.0283mm，平均每个测点每个方向的误差小于1×10^-4mm，同震目标函数值为3.93×10^-2m，每个测点每个方向的误差小于10^-3m，震前3期和同震反演误差均小于GPS测量数据的误差，故反演结果可靠。

图 4  震前三期和同震目标函数的演化过程

Figure 4.  Evolution of the pre-earthquake and coseismic objective function

下载: 全尺寸图片 幻灯片

同震断层附近地表位移如图 5所示，同震反演断层走滑位移为-3.029m，倾向位移为4.725m。震后地质考察结果显示（吴珍汉等，2008；刘健等，2012），中央断裂带中，映秀至北川断裂带平均逆冲位移为6.4m，北川至青川断裂带逆冲为5m，映秀至北川破裂带垂直位移量从0.2m到11m不等，伴随的右旋平均水平位移量为2—3m。所以，同震断层反演结果与地质观测结果相符合。陈运泰等（2008）通过地震波反演得出，断层平均逆冲位错为5m，断层平均走滑位错为2m，与本文的结果基本一致。

图 5  同震断层附近GPS测点位移矢量图

Figure 5.  Coseismic deformation displacement of GPS measurements near the fault

下载: 全尺寸图片 幻灯片

另一方面，根据断层滑移参数可以计算地震矩：

式中，μ为断层的刚性系数（本文取为3×10¹⁰Pa），S为断层面积（根据断层模型求出断层面积为$ S=L\times d/\sin \theta =300\times {{10}^{3}}\times 20\times {{10}^{3}}/\sin {{60}^{\circ }}=6.93\times {{10}^{9}} $m²），D为断层标量滑移。基于本文反演结果，$ D=\sqrt{{{3.029}^{2}}+{{4.725}^{2}}}=5.613 $m，于是求得$ {{M}_{0}}=1.17\times {{10}^{21}} $J，根据矩震级公式：

计算得矩震级M_W=7.985。所以，本文反演的结果与根据地震波计算的矩震级为M_W=7.9（张勇等，2008）基本相当。

为了与同震滑移数据进行对比，了解震前滑移与地震的关联，表 4给出了震前3期的断层滑移的反演结果。可以看出，震前3期总的平均走滑速率为-5.39mm/a，倾向速率为2.66mm/a，断层有逆冲和右旋走滑的趋势。由于走滑位移会受到断层两端的约束，而倾向位移不受地表约束更易发生滑动，所以从震前断层滑移参数分析可知，断层的倾向滑移以逆冲为主。这与震后观测结果及上述同震反演结果（逆冲和右旋运动）一致。这个结果可能预示着，可以通过震前断层的运动趋势来推测地震发生时断层的破坏形式。

表 4  龙门山断层震前3期反演结果

Table 4.  Inversion result of pre-earthquake fault parameters

量	U1'/m	U2'/m
震前第一期滑移	-10.98	4.37
震前第二期滑移	-15.80	7.73
震前第三期滑移	-14.89	9.06
震前第一期年均滑移	-5.76	2.29
震前第二期年均滑移	-5.44	2.66
震前第三期年均滑移	-4.97	3.02
震前3期平均滑移	-5.39	2.66

下载: 导出CSV 

| 显示表格

进一步地，断层震前各期的反演结果（表 4）表明，震前第一、二、三期的年均走向滑移分别为5.76mm/a、5.44mm/a、4.97mm/a，逐渐减小。而倾向的年均滑移分别为2.29mm/a、2.66mm/a和3.02mm/a，显示上升趋势。这表明断层的矢量运动方向从以走滑为主的右旋运动逐渐偏转到以倾向滑移为主的逆冲运动。该趋势预示着断层在未来地震时可能发生以逆冲为主的破坏，这与汶川地震时断层发生的巨大逆冲滑移相一致。

龙门山断层震前和同震的滑移参数可以反映地震的能量累计和释放量，进而能够估算汶川地震的复发周期。龙门山断裂带晚第四纪以来地貌学的研究结果显示（马保起等，2005），茂汶-汶川断裂、北川-映秀断裂和江油-灌县断裂晚第四纪逆冲滑动速率的总和约为1mm/a左右。龙门山断裂带历史上没有超过7级的地震，6级以上地震只有3次，由于地震震级相差一级释放的能量相差很大，所以忽略复发周期内地震释放的能量，由于临近地震，后2期的逆冲滑动有加速趋势，故将第一期逆冲滑动速率2.29mm/a作为长期积累量，每年平均释放量为1mm/a，用积累量与释放量之差计算净积累量为1.29mm/a，利用断层平均滑动位移4.725m与断层的滑动速率的净积累量之比计算得到地震复发周期为3700年，地质学估计的地震复发周期为2000—10000年（谢富仁等，2008；李玉江等，2012）。

3.   结论

（1）本文计算了震前3期断层滑移累积量，断层平均走滑速率为-5.39mm/a，倾滑向速率为2.66mm/a，与同震相比，断层运动形式一致，都为逆冲兼右旋的滑动。

（2）比较震前3期的滑移速率，从第一、二期到第二、三期逆冲滑移速度增加，说明逆冲方向有加速的趋势，断层滑移加速方向和地震发生时断层的滑动形式一致。

（3）根据震前、同震断层滑移和地质资料，利用断层同震逆冲滑移与震前断层滑动速率的净积累量推算了汶川地震的复发周期，结果表明汶川地震复发周期为3700年。