目前采用地震波数据确定地壳应力场大体可分为两种方法:(1)先确定地震的震源机制，然后采用多个地震的震源机制确定应力场的方向和相对大小^{[1， 2]}.(2)采用一个地区多个地震的大量P 波初动确定应力场^[3~7].采用第一种方法必须首先确定研究区域的多个震源机制，假定这些震源机制的空间取向是随机的，然后根据震源机制确定该区域的地壳应力场.目前确定震源机制的常规方法包括利用P波初动信息^[8]、P/S 振幅比信息^[9~13]和波形信息^[14~17]，然而震源机制解目前仍然存在很大的不确定性^[18]，而且能够确定震源机制的地震数目有限，这就使得地壳应力场的求解存在较大的不确定性.采用第二种确定地壳应力场的方法虽然只利用P波初动信息，但它可以利用大量地震(同样假定震源机制的空间取向是随机的)，且包括只有少量地震观测资料而无法确定其震源机制的小震观测信息，数据量大大增加使得应力场求解结果比较稳定.

虽然目前利用大量地震的P 波初动信息已经求得了某些地区的地壳应力场^{[6， 7]}，但是目前采用P波初动信息求解应力场的方法仅考虑了应力方向导致的大量地震产生的P 波初动的辐射花样^[6](以下简称综合P 波辐射花样)，并没有考虑到产生地震的应力场主应力值的相对大小(即假定应力场主应力值等差排列，仅考虑某一应力状态下的最大剪切应力的平面为破裂面)及与地震错动有关的摩擦系数(假定摩擦系数为零)的影响.

应力相对大小可以根据震源机制或断层滑动资料求得^{[1， 2， 19]}，是解释地球动力学模型的重要参数^[20].摩擦系数也是影响地震滑动的一项重要参数^[21~24].这两个参数会影响地壳中位错破裂的空间形态，进而影响大量地震产生的综合P 波辐射花样.欲研究根据大量P 波初动资料确定地壳应力方向、应力相对大小及摩擦系数问题，必须先研究应力方向及其相对大小和摩擦系数对综合P 波辐射花样的影响.

本研究拟通过模拟方式给出不同应力场状态和摩擦系数导致大量均匀取向断层破裂产生的综合P波辐射花样的空间模式，以期在P 波初动资料足够多的情况下，研究是否能同时确定应力场主应力的相对大小和摩擦系数.

为了探讨应力场对综合P 波辐射花样的影响，首先研究给定应力场在均匀取向断层面上的剪应力、摩擦力和库仑破裂应力的分布.由于库仑破裂应力(摩擦系数为零时退化为剪应力，这时又称为粘聚力)表示了断层滑动的可能性及与研究区域的地震活动性相联系^{[25， 26]}同时被地质断层活动状态所验证^[27]，因此可将按照库仑破裂准则破裂的地震产生的P 波辐射花样(辐射花样大小与库仑破裂应力大小成正比)叠加得到研究区域综合P波辐射花样.

假定地壳应力场的最大主应力、中间主应力和最小主应力分别为σ₁，σ₂，σ₃，拉张为正.在由地震震源机制解求取应力场方向时还可以得到相对应力大小R，Gephart和Forsyth^[1]给出的R定义式为

(1)

通常很难确定地下的绝对应力^[28]，在研究地壳应力场和地震破裂的关系时通常假定σ₃ =-1，这里考虑的是应力各分量的相对值，忽略单位.则根据(1)式可得^[29]

(2)

(3)

在该应力场作用下，假定各种空间取向断层均可能发生破裂.断层面法向与σ₃轴夹角呈δ，与σ₁ 轴夹角呈Φ 的方向(图 1)可表示为

(4)

图 1

Fig. 1

图 1 地震断层和应力方向示意图 Fig. 1 Skepch map of seismic fault and stress direction

应力在此断层面上产生的剪切应力方向为

(5)

与n和t垂直的方向为

(6)

应力张量在法向为n的地震断层面上的应力矢量(s1，s2，s3)计算公式为^[30]

(7)

法向应力pv的计算公式为

(8)

在t方向上产生的剪切应力p_t为

(9)

在b方向上产生的剪切应力p_b为

(10)

则应力场在此断层面上产生的合成剪切应力为

(11)

应力场在此断层面上产生的摩擦应力为

(12)

此处，μ 为摩擦系数.则此断层面上产生的库仑破裂应力为

(13)

按照上面的计算步骤，计算了R=0.0，0.25，0.5，0.75，1.0五种情况下的剪应力、摩擦应力和库仑破裂应力的空间分布.将剪应力、摩擦应力和库仑破裂应力的大小表示在空间取向断层面的法向方向上.注意这里的应力量为相对大小，摩擦系数取0.6.可以看到，在R=0.0(即σ₂=σ₃)的情况下，剪应力在法向与σ₁ 轴夹角相同的断层面上大小相同，因此具有相同的破裂可能性(图 2a)；同时，由于σ₂ =σ₃=-1(压应力)，此时法向在σ₂σ₃平面上断层面的摩擦应力均为负，阻碍了断层面破裂(图 2b).计算的库仑破裂应力关于σ₁ 轴对称，法向在σ₂σ₃平面及其邻域内的断层面上，由于摩擦应力的存在，破裂不可能发生，而在法向与σ₁ 轴(张轴)夹角约30°的断层面上(断层面与σ₂σ₃平面夹角约60°)库仑破裂应力最大，这可以根据摩擦系数为0.6时与压轴夹角的计算公式$\frac{1}{2}$tan^-1(1/μ)=29.5°^[31]所验证.随着R的增大(σ₂ 值逐渐远离σ₃值)，剪应力逐渐呈现四象限分布，摩擦应力逐渐表现为在σ₁ 轴上达到最大值，在σ₃轴上达到最小值.到R=0.5(σ₁，σ₂，σ₃等差排列)时剪应力达到完全四象限分布，此时剪应力最大的断层面与σ₁、σ₃轴夹角为45°(图 2g)，库仑破裂应力最大的断层面与σ₃轴的夹角为60°(法向与σ₃轴的夹角约30°，图 2i).随着R的进一步增大(σ₂ 逐渐靠近σ₁)，剪应力、摩擦应力和库仑破裂应力逐渐变为以σ₃为对称轴的轴对称分布(图 2j，k, l)，到R=1.0(σ₁=σ₂)时，剪应力(图 2m)、摩擦应力(图 2n)和库仑破裂应力(图 2o)在断层面法向上的分布完全以σ₃轴对称.

图 2

Fig. 2

图 2 典型应力状态(R=0.0，0.25，0.5，0.75，1.0)下的剪应力、摩擦系数为0.6的摩擦应力和库仑破裂应力的空间分布 应力相对大小表示在所取断层面的法向上.(a)、(d)、(g)、(j)、(m)分别为R=0.0，0.25，0.5，0.75，1.0情况下的剪应力空间分布；(b)、(e)、(h)、(k)、(n)分别为R=0.0，0.25，0.5，0.75，1.0的摩擦应力空间分布；(c)、(f)、(i)、(l)、(o)分别为R=0.0，0.25，0.5，0.75，1.0的库仑破裂应力空间分布.蓝色到绿色表示应力为负，绿色到红色表示应力为正. Fig. 2 The spatial distribution of the shear stress, friction stress (The friction coefficient is choose as 0.6) and theCoulomb failure stress of the fault varies as a function of its orientation with respect to different stress ratios CR=0.0，0.25，0.5，0.75，1.0) The relative stress magnitude is expressed on the normal direction of the fault plane.Figures (a)，(d)，(g)，()and (m) show shear stress distributions correspond to R=0.0，0.25，0.5，0.75，1.0 respectively； Figures (b)，(e)，(h)，(k) and (n) show the spatial friction stress distributions correspond to R=0.0，0.25，0.5，0.75，1.0 respectively； Figures (c)， (!)， (i)，(l) and (o) are the spatialcoulomb failure stress distributions correspond to R=0.0，0.25 » 0.5，0.75，1.0 respectively.Blue to green express the negative stress, green to red express the positive stress.

由此分析可知，在三轴应力状态下，断层面破裂的可能情况不仅与三轴的空间方位有关，而且与三轴应力的相对大小有关，通常情况下所说的最大剪应力分布在与σ₁ 和σ₃轴夹角为45°的断层面上，只有在理想情况下(σ₁<σ₂<σ3)才是完全正确的.在σ₁=σ₂或σ₂=σ₃的情况下最大剪应力或库仑破裂应力分布在一个环形圆周上(图 2a，c和m, o).在分析地球内部处于三轴应力状态时要特别注意.

根据以上应力状态可以给出每个断层的破裂可能性，库仑破裂应力(剪应力可看成摩擦系数为零的库仑破裂应力)越大的断层破裂的可能性越大，单个地震产生的P波辐射花样越大.将地震断层面在空间均匀取样，求取每个断层面上的库仑破裂应力大小从而确定单个地震断层面滑动的P 波辐射花样大小，将单个地震断层面产生的辐射花样叠加得到综合P 波辐射花样的空间分布，研究其在空间的变化.

球坐标系Φ 在0~360°范围内均匀取样，δ 在0~180°范围内均匀取样可以得到覆盖整个断层面取样空间的应力状态分布.然而，这样取值导致了δ趋近0或180°时取值较密集，而δ趋近90°取值较稀疏^[32]，断层面空间取向取值不均匀，从而使综合P波辐射花样出现偏差.为了改善这种情况，我们采用在不同的δ 取不同的Φ 采样间隔:dΦ/sinδ，dΦ 为取样间隔.这样可使得断层面空间取向的均匀.

将第2节求得的p_b、p_t及摩擦系数μ，求得剪应力或库仑破裂应力的大小，其方向均按照(14)式计算:

(14)

给出每个地震断层面错动产生的相对于主轴的立体辐射花样^[33]，按其大小与库仑破裂应力数值成比例的方式在某一方向进行叠加，可得到综合P 波辐射花样(图 3).

图 3

Fig. 3

图 3 均匀空间取向的单个断层面错动产生的P波辐射花样叠加得到的综合P波辐射花样. (a)，(d)，(g)，(j)，(m)分别为R=0.0，0.25，0.5，0.75，1.0、摩擦系数为零(即剪应力)产生的综合P波辐射花样；(b)，(e)，(h)，(k)，(n)分别为R=0.0，0.25，0.5，0.75，1.0、摩擦系数为0.6产生的综合P波辐射花样；(c)，(f)，(i)，(l)，(o)分别为R=0.0，0.25，0.5，0.75，1.0、摩擦系数为1.2产生的综合P波辐射花样.蓝色到绿色表示P波初动为负，绿色到红色表示P波初动为正. Fig. 3 The composite P wave radiation patterns by superposition of the individual P wave radiation patterns with a uniform orientation distribution of fault planes. Figure (a)，(d)，(g)，(j) and (m) are the composite P wave radiation patterns correspond to ^ = 0.0 » 0.25 » 0.5 » 0.75，1.0 respectively, with the friction coefficient is 0 ； Figure (b)，(e)，(h)，(k) and (n) are the composite P wave radiation patterns correspond to R = 0.0，0.25，0.5，0.75，1.0 respectively, with the friction coefficient is0.6； Figure (c)，(f)，()，(l) and (o) are the composite P wave radiation patterns correspond to R=0.0，0.25，0.5，0.75，1.0 respectively, with the friction coefficient is 1.2； Blue to green express the negative P wave first motion, green to red express the positive P wave first motion.

由图 3(a, b, c)可知，在R=0(σ₂=σ3)时库仑破裂应力产生的P波辐射花样在σ₁ 轴及其周围为正，且随着远离σ₁ 轴逐渐减少，在与σ₁ 轴夹角为45°的圆周上减少至零，以后随着逐渐靠近σ₂σ₃平面，P波振幅逐渐增大，但极性为负，在σ₂σ₃平面上极性为负的振幅达到最大.R=0的综合P波辐射花样很难识别出摩擦应力所起的作用，即μ=0.0(a)，0.6(b)，1.2(c)产生的综合P波辐射花样看不出有明显的模式变化.随着R值的增大(σ₂ 值逐渐远离σ₃值)，库仑破裂应力产生的综合P 波辐射花样有所区别，表现为P波初动振幅在σ₂ 上的分布不同和远离σ₃轴和σ₂ 轴的振幅衰减情况不同.R=0.5(σ₁，σ₂，σ₃等差排列)、摩擦系数为零产生的综合P 波辐射花样呈现完全的四象限分布(图 3g)，但μ =0.6(图 3h)和1.2(图 3i)产生的综合P波辐射花样在σ₂ 上具有负振幅初动(振幅大小随摩擦系数的增大而增大).随着R的进一步增大(σ₂ 值逐渐靠近σ₁ 值)，摩擦系数为零的库仑破裂应力产生的综合P 波辐射花样在σ₂ 轴上逐渐呈现正振幅初动，振幅大小随R的增大而增大.在R>0.5的情况下，摩擦应力使得σ₂ 轴上正初动振幅减小(图 3k，l).当R=1.0(σ₂ =σ₁)时，库仑破裂应力产生的综合P波辐射花样在σ₃轴及其周围为负，且随着远离σ₃轴逐渐减少，在与σ3轴夹角为45°的圆周上减少至零，以后随着逐渐靠近σ₁σ₂ 平面，正极性的P波振幅逐渐增大，在σ₁σ₂ 平面上极性为正的振幅达到最大，呈现以σ₃轴为对称轴的轴对称分布，并且不同摩擦系数产生的综合P波辐射花样没有明显变化，因此在这种情况下，也无法看出摩擦系数的影响.

综合考虑图 3的各种情况可以看到，R值和摩擦系数均会对综合P 波辐射花样产生影响，在R=0.0和R=1.0 的情况下看不出摩擦系数导致的综合P波辐射花样的差别，即看不出摩擦系数产生的影响.当R值处于两者之间时，摩擦所起的作用不可忽略.但在这里需要注意的是，虽然R值和摩擦系数均对P波综合辐射花样造成影响，但是同时精确确定这两个参数却是困难的，例如，图 3g和图 3l分别为R=0.5、μ=0.0和R=0.75、μ=0.6时产生的综合P波辐射花样没有明显变化，即使稍微有变化也可以通过调整摩擦系数或R值而达到完全一致，图 3c和图 3f、图 3e和图 3i也具有类似情况.如果将R=0.5、μ=1.2 情况下的摩擦系数增大(图3i)，则库仑破裂应力产生的辐射花样会逐渐趋近于R=0.25，μ=0.6 的综合P 波辐射花样(图 3e).这说明，在一定的允许范围内，R值的减小和摩擦系数的增大具有相同的效应，仅仅根据P 波辐射花样完全确定R 值和摩擦系数是困难的.

为研究不同应力状态下综合P 波辐射花样的空间形态以及摩擦系数对综合P 波辐射花样的影响，我们确定了空间任意取向断层上的剪应力、摩擦应力和库仑破裂应力相对大小，并绘制了三维图形.结果表明在三轴应力下，剪应力、摩擦应力和库仑破裂应力的空间分布随R 值的不同有很大不同，并且摩擦系数在库仑破裂应力分布中也起了一定作用.

根据不同应力状态下的剪应力和库仑破裂应力大小，假定地震在各个断层取向上均匀发生，我们给出了考虑摩擦系数的不同应力状态下的综合P 波辐射花样.研究发现综合P 波辐射花样随R值和摩擦系数不同而有较大的差别.综合P 波辐射花样只有在特定情况下才呈现规则的四象限分布(如R=0.5、μ=0.0).在R=0.0和1.0情况下产生的综合P波辐射花样分别以σ₁ 轴和σ₃轴为对称轴的轴对称分布，并且不能识别出摩擦系数对P 波综合辐射花样的影响.在其他情况下摩擦系数和R值同时影响P波辐射花样，并且R值的增大和摩擦系数的增大具有相同的效应.这说明仅根据P 波初动符号完全确定R值和摩擦系数是困难的.这也说明了前人根据综合P 波初动确定应力场假定R=0.5，μ=0.0^[3~6]的相对合理性.要完全确定应力相对大小和摩擦系数，还需要其他资料，如岩石力学实验资料的约束.虽然如此，摩擦是岩石力学实验所验证的一种因素^[31]，如果相对合理的假定摩擦系数大小(如μ=0.4^[34])，则可以根据综合P波初动资料识别出应力相对大小R值，这是我们下一步的工作.

本研究假定地震在各个方位取向断层面上均匀发生，如果地震沿均一断层面破裂，本研究摩擦系数对P 波综合辐射花样影响的结论将不成立.然而，目前的观测表明，即使大震破裂后的余震区，小地震震源机制也呈现非常大的多样性^{[35， 36]}，如果对于没有大地震破裂的地区，地震更可能在各个取向断层面上均匀破裂.这可能是由于地球内部具有很大压力，致使地震在各个取向断层面均可发生的缘故.

本研究将空间取向断层面上的单个地震破裂产生的P波辐射花样按照库仑破裂应力大小进行加权，叠加得到综合P 波辐射花样，这只是考虑到一个统计平均.其实，也许地震并不是库仑破裂应力越大，产生的地震越大，但这里考虑的综合效应是断层面上库仑破裂应力越大，发生同等地震的概率越大，同样会得到对该断层面上地震的P 波辐射花样加权大的效果，本模拟结果不会因此而改变.

Xu等^[6]在综合震源机制确定平均P、T轴的模拟研究中假定断层上的剪应力超过最大剪应力的一半时才会使断层破裂，进而根据大量地震的P 波初动信息反推平均P、T轴.本研究没有考虑此假定.如果考虑此假定，需要研究剪应力究竟达到多大才会发生地震，目前这个问题还没有定论.这里考虑较小的库仑破裂应力的地震破裂会产生较小的P 波辐射花样，其叠加的效果也不会由于考虑了Xu等^[6]的假定而明显改变本研究的模拟结果.

本文的模拟结果为采用大量地震产生的P 波初动信息反演应力场提供了考虑摩擦系数和应力相对大小的思路.随着地震资料的积累和地震台站的加密，越来越多的P 波初动资料可用于应力场反演，为本文模拟结果的应用提供了广阔的应用前景.基于此模拟结果编制考虑摩擦系数的应力场(包括相对应力大小)反演程序，进而用于精确应力场的反演是我们下一步的工作目标.