引言

在震源物理研究中，作为二阶对称张量的地震矩张量已被广泛接受并且得到了成功的应用 (Gilbert, Dziewonski, 1975; McCowan, 1976; Fitch et al, 1980; Dziewonski et al, 1981; Ekström, 1989).基于天然地震是发生于地球内部的震源 (内源) 这一前提，从角动量守恒可得出地震矩张量具有对称性的结论 (Backus, Mulcahy, 1976a, b; Aki, Richards, 1980; Kennett, 1983).然而，地震矩张量的对称性不过是在一定条件下引入的简化和限制，并不是绝对必要的 (Takei, Kumazawa, 1995).在平面断层的假设下，对称矩张量是对地震震源的很好近似.但是，若研究更为接近实际的震源模型，例如考虑震源区的质量流动 (Takei, Kumazawa, 1994)，或断层具有厚度和强度弱化区 (Knopoff, Chen, 2009)，则需要引入非对称矩张量以表示相应的地震效应.

若考虑到地震是发生在具有厚度的强度弱化区的剪切断层上，那么与其远场地震波对应的震源项，要比不考虑具有厚度的强度弱化区的剪切断层 (Burridge, Knopoff, 1964) 多出一个单力偶项 (Knopoff, Chen, 2009).由于与该单力偶项对应的地震矩张量是非对称张量，因而与计及断层厚度和强度弱化区域的剪切位错震源所对应的地震矩张量也是非对称张量.与对称地震矩张量不同的是，非对称地震矩张量中的反对称部分使得断层面解中方向相反的两对正交单力偶的矩不再相等.因此，根据具有厚度的断层模型，可判定与力偶矩较大的单力偶相联系的节面为断层面，另一节面则为辅助面 (刘超, 陈运泰，2014).

本文将在上述工作的基础上，进一步研究时间域内的非对称地震矩张量反演的理论与方法，特别是研究从对称地震矩张量反演到非对称地震矩张量反演，是否存在过度拟合，以及影响非对称地震矩张量反演的因素等问题；并将通过数值试验，检验非对称地震矩张量反演方法的可行性.

1.   理论与方法

1.1   理论

非对称地震矩张量反演的理论与方法和对称地震矩张量反演类似，只需在其基础上略加改动即可实现.参考对称地震矩张量反演的思路，可将作为观测点位置x和时间t函数的地面运动位移u_i(x, t) 表示为表征震源特性的、含时非对称地震矩张量M_jk(t) 与其相对应的格林函数G_{ij, k}(x, t; 0, 0) 卷积的线性组合 (刘超等, 2008; 刘超, 2011; 刘超, 陈运泰, 2014):

式中“*”表示卷积.在具有厚度的地震断层模型中 (Knopoff, Chen, 2009; 刘超, 2011; 刘超, 陈运泰, 2014)

式中，M_jk^u(t) 和M_jk^T(t) 分别为与断层面上的位错u和应力错T相联系的、含时地震矩张量,

式中, m_jk^u(ξ, t) 和m_jk^T(ξ, t) 分别为与断层面上的位错u和应力错T相联系的、作为位置ξ与时间t函数的地震矩张量密度 (又称“地震矩密度张量”).由式 (2) 可知，由于M_jk^u(t) 为二阶对称张量，M_jk^T(t) 为二阶非对称张量，所以M_jk(t) 为二阶非对称张量，即M_jk(t)≠M_kj(t).

由式 (2)，可将式 (1) 改写为

若以S^u(t) 和S^T(t) 分别表示与M_jk^u(t) 和M_jk^T(t) 对应的归一化震源时间函数 (以下简称震源时间函数), M_jk^u和M_jk^T分别表示与M_jk^u(t) 和M_jk^T(t) 对应的不含时非对称地震矩张量：

则

对于断层面上的某一点ξ，地震矩张量密度m_jk^u(ξ, t) 与m_jk^T(ξ, t) 的时间历史不同，即m_jk^u(ξ, t) 对应滑动开始至滑动停止的过程，持续时间为T_s；m_jk^T(ξ, t) 则对应强度弱化过程，持续时间为T′_s，且T_s>T′_s.对于某一次地震事件，M_jk^u(t) 的持续时间为整个断层面的破裂时间T_s^t，M_jk^T(t) 的持续时间为T_s^t－(T_s－T′_s).由于断层面上某一点的滑动持续时间T_s远小于整个断层面的破裂时间T_s^t(Heaton, 1990; Nielsen, Madariaga, 2003)，则有T_s≤T_s^t，所以T_s－T′_s≤T_s^t，也即M_jk^T(t) 的持续时间近似等于T_s^t.因此，M_jk^u(t) 与M_jk^T(t) 具有近似的、由断层面的破裂过程决定的震源时间函数：S^u(t)≈S^T(t)≈S(t)，故由式 (6) 可得

这样，利用上式即可将对含时非对称地震矩张量M_jk(t) 的反演转换为分别对不含时非对称地震矩张量M_jk和震源时间函数S(t) 进行的线性反演.若反演M_jk，则有9个待解参数；若反演S(t)，则当用N个离散点表示S(t) 时，便有N个待解参数.理论上, 震源时间函数S(t) 是由断层面的破裂过程决定的, 与波传播速度无关.实际上, 由于破裂传播效应, 在不同台站 (因而不同方位)、由不同波型的波 (如P波、S波) 反演得到的S(t) 是视震源时间函数，与波型 (波传播速度) 有关, 因此，我们分别称由P波和S波反演得到的震源时间函数为P波震源时间函数S^P(t) 和S波震源时间函数S^S(t).

1.2   方法

与对称地震矩张量反演对比可见，从对称地震矩张量反演到非对称地震矩张量反演，只是增加了3个待解参数.因此, 对称地震矩张量反演的方法和步骤可以应用于非对称地震矩张量反演，只是当地震矩张量对称时，M_jk= M_kj，地震矩张量的待解参数减为6个；当地震矩张量非对称时，地震矩张量的待解参数为9个.因此，只要将对称地震矩张量反演方法和程序作适当改动, 便可用于非对称地震矩张量反演.

在反演非对称地震矩张量时，首先利用对称地震矩张量反演，分别得到垂直向P波和水平向S波的震源时间函数，然后再联合反演非对称地震矩张量，具体步骤如下:

1) 选择9个地震矩张量作为基本地震矩张量M′_jk进行一次正演，得到与M′_jk相应的格林函数G_{ij, k}.这9个基本地震矩张量M′_jk为

2) 合成对称地震矩张量反演所需的格林函数G_{ij, k}(x, t; 0, 0)；

3) 利用P波进行对称地震矩张量反演，得到P波震源时间函数S^P(t)；

4) 利用S波进行对称地震矩张量反演，得到S波震源时间函数S^S(t)；

5) 用u_i^P和G_{ij, k}^P分别表示P波观测数据和相应的格林函数，用u_i^S和G_{ij, k}^S分别表示S波观测数据和相应的格林函数，令

将其代入式 (7), 可以看出, 此时在式 (7) 中, 仅有地震矩张量M_jk未知；

6) 将式 (7) 中的S(t) 与G_{ij, k}(x, t; 0, 0) 卷积得到

即可将式 (7) 进一步简化为

据此即可通过线性反演得出地震矩张量M_jk.

由式 (7) 可知，在反演非对称地震矩张量M_jk时，需要分别计算出M_jk中9个分量对应的格林函数G_{ij, k}用于构建反演问题的系数矩阵.非对称地震矩张量的各个分量能否被精确地反演，与位移数据是否包含有相应的震源信息密切相关.在对称地震矩张量反演时，利用垂直向P波位移数据即可反演其6个分量.然而，在非对称地震矩张量反演时，若仅用垂直向P波位移数据将会无法区分M_xy与M_yx这两个分量，需加入水平向信号进行联合反演 (刘超, 2014).这是非对称地震矩张量反演与对称地震矩张量反演的最大区别.

考虑到S波的水平向地震信号有较高的信噪比，因此可运用S波水平向数据与P波垂直向数据联合反演非对称地震矩张量.对于S波，可将计算的格林函数旋转至N--S分量和E--W分量后，再用于反演.不用SH波而分别选用SH波的N--S分量和E--W分量，是因为在数据处理实践中，直接用N--S分量和E--W分量可获得更多可靠的S波信号. SH波需要通过旋转N--S向和E--W向的水平记录得出，这要求两个水平分向的仪器响应参数、标定参数、极性和授时等均准确无误.然而在处理S波波形数据时，两个水平分量记录的不一致性很难识别 (如极性反转，标定不准确，仪器响应误差或授时不一致等)，使得部分错误的SH波会被用于反演，增大了反演的不确定性.另一方面，即使识别出不一致的水平向地震数据，也只能将两条数据都舍弃，从而减少了可用于反演的SH波数据量.但是，如果我们单独选择可用的N--S分量或E--W分量数据用于反演，则只需旋转理论计算的水平向格林函数，从而避免对N--S分量或E--W分量数据一致性的检查，并且使得可用的S波波形数据增多.在S波与P波的联合反演中，数据的反演权重会自动调整，使得每一条单独的P波或S波波形数据，均具有相同的反演权重.

1.3   非对称地震矩张量反演是否过度拟合

即使不涉及震源的理论模型，单从数据反演的角度，需要回答的一个重要问题是：非对称地震矩张量反演是否只是因为增加了3个模型参数而提高了对数据的拟合程度，因此增加的3个参数并不是必需的？或者说，从对称地震矩张量反演到非对称地震矩张量反演，是否存在过度拟合问题？

为回答上述问题，本文引入赤池信息准则 (Akaike Information Criterion, 简称AIC准则).按照AIC准则，衡量是否过度拟合的AIC值可表示为 (Akaike, 1974)

式中已略去常数项，N_c为波形数，N为每条波形的抽样点数，N_MT为反演地震矩张量的自由度 (对于对称地震矩张量N_MT=6，非对称地震矩张量N_MT=9)，E为反演残差，其定义为

式中，u_i为观测波形，s_i为拟合波形.利用AIC值可挑选同一反演条件下的最佳模型.所谓最佳模型是指该模型既能很好地拟合数据，又具有最少的模型参数. AIC值的大小不具有绝对意义，仅具有相对意义.在比较不同模型的AIC值时，较小的AIC值对应较优模型.

在非对称地震矩张量反演时，可利用同一组数据进行对称地震矩张量反演.比较两次反演结果的AIC值即可判定对称地震矩张量与非对称地震矩张量孰为更优的模型，从而判定非对称地震矩张量反演是否存在过度拟合的问题.

2.   地震矩张量的矢量表示法

在反演结果的误差分析中，需要定量地描述地震矩张量之间的差异，这涉及到地震矩张量的分解和分类方法.地震矩张量的分解是不唯一的，地震矩张量的分解与分类应视问题的物理内涵而定.在设定地震矩张量的分解分类标准时，往往希望这种分类方法既能直接服务于特定的研究内容，又能具有清晰直观的物理图像.

常见的地震矩张量分类表示方法通常具有以下一些不足：对地震矩张量间的差异无法完整统一地定量描述；用于描述地震矩张量的参数的权重难以定量确定；方法难以推广至非对称地震矩张量的分析中.在地理坐标系中，可通过定义走向、倾角和滑动角来表示地震矩张量中的直流分量 (“零频”分量) 所对应的断层面参数，但这3个参数间的权重难以定量地确定:倾角较小时，震源类型的相似性较大；倾角较大时，震源类型的相似性较小. Kagan (1991)提出利用两个地震矩张量的直流分量对齐所需的最小空间旋转角 (通常称为Kagan角) 来描述其间的差异；但是，考虑到两个完全相反的直流分量 (T轴与P轴互换)，绕B轴转动90°即可对齐，而另一些直流分量之间则需转动120°才能对齐，可见Kagan角并不是一个能够特别合理地描述地震矩张量之间差异的参数.利用地震矩张量的T-k分类法 (Hudson et al, 1989) 可以定量地分析各向同性分量 (isotropic component，简写为ISO)，双力偶分量 (double couple component, 简写为DC) 和补偿线性矢量偶极分量 (compensated linear vector dipole component, 简写为CLVD)，其中T表示震源中体积为常量 (即剪切) 分量的类型，k表示体积变化 (即膨胀) 分量的大小.虽然如此，但仍然难以确定T-k参数与本征轴参数之间的权重.尤为重要的是，上述方法均很难运用于非对称地震矩张量分析.

为克服上述困难，这里引入地震矩张量的矢量表示法 (Willemann, 1993)，利用地震矩张量内积的性质，将地震矩张量之间的距离定义为矢量之间的夹角.地震矩张量的矢量表示法还可以自然地推广至非对称地震矩张量的分析中.

在以e_i为基矢量的直角坐标系中，张量M=M_ije_ie_j和N=N_ije_ie_j的内积，即

对于迹为零 (ISO=0) 的地震矩张量，其自身的内积作为第三不变量，与坐标系选择无关 (Kagan, Knopoff, 1985).利用标量地震矩M₀可以定义归一化地震矩张量m_ij，使得对于任意地震矩张量M_ij，均满足

即归一化地震矩张量m_ij自身的内积为2.很自然地，可利用归一化地震矩张量的内积J₃描述m_ij与n_ij之间的差异，即

J₃在[-2, 2]上变化：当m_ij=n_ij时，J₃= 2；当m_ij=-n_ij时，J₃ =-2.

定义m_ij与n_ij之间的距离D为

这样定义的D具有角度的量纲，与地震矩张量在矢量空间中的旋转相联系. D在[0, π]上变化：当地震矩张量相同 (m_ij=n_ij) 时，距离最小, D=0；当地震矩张量相反 (m_ij=-n_ij) 时，距离最大, D=π.

选择这种定义距离的方法，使得D在数学上有以下优点：

1) 距离D与坐标系选择无关，是地震矩张量之间差异的客观量度；

2) D是非负的；

3) D满足三角不等式D_ij+D_jk≥D_ik；

4) 相同的地震矩张量之间的距离最小，为0；相反的地震矩张量之间的距离最大，为π；

5) 可自然地推广至非对称地震矩张量.

从距离D的定义还可以看出，与m_ij距离为D的地震矩张量n_ij是不唯一的，这一特点在地震矩张量的矢量表示中需要注意.

3.   格林函数对地震矩张量反演的影响

本文运用广义反射、透射系数矩阵法 (Kennett, Kerry, 1979; Kennett, 1980, 1983) 和相关程序 (李旭, 1993; 李旭等, 1994; 李旭, 陈运泰, 1996a, b) 计算格林函数.将原程序用Matlab重新编写并优化，以提高计算速度.利用该方法可较快地求出分层均匀介质中的全波理论地震图，因受限制因素较少，适于反演中格林函数的计算.

本节将通过分析格林函数与矩张量各分量之间的关系说明：在非对称矩张量反演时，若仅用垂直向地动位移数据，将无法区分M_xy与M_yx这两个分量；需要同时运用垂直向与水平向地动位移数据进行联合反演才能区分M_xy与M_yx.

利用地球展平变换，将地球的球形分层速度结构模型变换为水平分层速度结构模型.采用柱坐标系 (r, θ, z)，坐标原点取在震中，z轴取垂直向下为正，相应的基矢量为 ()，坐标系遵循右手法则.

引入一组面谐矢量 (R_k^m, S_k^m, T_k^m)

式中，Y_k^m(r, θ)=J_m(kr) exp (imθ)，，J_m(kr) 是第一类m阶贝塞尔函数，m是整数.位移u可表示为

式中，位移u各分量的表示式为

由上式可以看出，经过面谐矢量展开，在位移u的表示式中，u_z只与U有关，与V, W无关；u_r和u_θ则均与V, W有关，与U无关：

假设在水平分层均匀介质内部z=z_S处有一点源，它在自由表面z=z₀所产生的位移为 (Kennett, Kerry, 1979)

式中，R_EV为自由表面的接收系数矩阵，为自由表面的反射系数矩阵，R_D^RS为自由表面与震源之间下行波的广义反射系数矩阵，T_U^RS为自由表面与震源之间上行波的广义透射系数矩阵，R_D^SL为震源与底界面之间下行波的广义反射系数矩阵，R_U^FS为震源与自由表面之间上行波的广义反射系数矩阵 (包括自由表面的反射作用)，Σ_D和Σ_U为表示震源的波矢量间断矩阵.

由式 (23) 可看出，面谐矢量位移解W(z₀) 包含了自由表面的反射、层间的多次反射、透射以及PS--V波的相互转换等所有作用，是完全响应的位移解.为突出震源项与位移解的关系，可将其改写为

对于P波和SV波，

将其代入式 (23) 得

将上式改写为

对于SH波，

将其代入式 (23) 得

由式 (27) 和式 (29) 可以看出，面谐矢量位移解W中，U和V与震源项φ^D, ψ^D, φ^U, ψ^U有关，W与震源项χ^D和χ^U有关.利用式 (27) 和式 (29)，可将式 (22) 改写为

在直角坐标系中 (李旭, 1993；李旭等，1994; 李旭, 陈运泰，1996a):

当m=0时，

当m=±1时，

当m=±2时，

式中，

将含有地震矩张量的表示式代入震源波矢量间断的表示式，并略去与单力ε_x, ε_y和ε_z相关的项，得到：

当m=0时,

当m=±1时,

当m=±2时,

式 (35)—(37) 即为仅含矩张量项的震源波矢量间断的表示式.

现在来分析垂直向位移的表示式.由式 (30) 的第一式

可以看出：在φ^D, ψ^D, φ^U, ψ^U的表示式 (35)—(37) 中，式 (35) 和式 (36) 中不含有与M_xy和M_yx有关的项，式 (37) 中的M_xy和M_yx只以M_xy+M_yx的形式出现，所以垂直向位移u_z仅与M_xy+M_yx有关，不携带任何单独与M_xy或M_yx有关的信息；若单用垂直向位移u_z进行反演，则只能获得M_xy+M_yx的大小，而不能单独分辨出M_xy和M_yx的大小.在对称地震矩张量反演中，由于M_xy=M_yx，利用M_xy+M_yx可得M_xy=M_yx=(M_xy+M_yx)/2；但是，在非对称地震矩张量反演中，由于M_xy≠M_yx，所以仅利用垂直向的位移u_z无法单独反演出M_xy和M_yx.

若引入水平向位移u_r或u_θ与垂直向位移u_z进行联合反演，由式 (30) 的第二和第三式

可以看出：u_r和u_θ不仅与震源项φ^D, ψ^D, φ^U, ψ^U有关，还与震源项χ^D和χ^U有关.在χ^D和χ^U的表示式式 (35) 中，M_xy和M_yx以M_xy-M_yx的形式出现.这表明u_r和u_θ所携带的与M_xy和M_yx有关的信息不仅以M_xy+M_yx的形式出现，还以M_xy-M_yx的形式出现，这使得利用u_r和 (或)u_θ与u_z进行联合反演，有望单独反演出M_xy和M_yx.

利用u_z即可单独反演地震矩张量的其它分量，不存在上述分辨问题，分析方法与对M_xy分量和M_yx分量的分析类似.在下一节数值试验中将通过数值计算验证本节的结论，即仅利用u_z无法单独反演出M_xy分量和M_yx分量，需要利用u_r和 (或) u_θ与u_z进行联合反演, 才有可能单独反演出M_xy和M_yx.

4.   数值试验

4.1   用于数值试验的数据

为了检测上述非对称地震矩张量反演方法的可行性，本文采用合成数据进行了数值试验.数值试验中设地震震源深度为10 km，标量地震矩M₀=1×10¹⁹ N·m，相当于矩震级M_W6.6的地震.选择两组具有代表性的非对称地震矩张量用于数值试验：第一组中的非对称地震矩张量的对称部分，对应于左旋走滑为主的震源机制，有少量CLVD分量，用M₁表示；第二组中的非对称地震矩张量的对称部分，对应于逆冲为主的震源机制，含有较大的CLVD分量，用M₂表示.同一组矩张量，均具有相等的对称部分和幅值不同的反对称部分.

第一组M₁中的非对称地震矩张量的对称部分M_1sym为

利用与M_1sym有关的反对称地震矩张量M_{1anti_sym}

可以合成得到非对称地震矩张量M_1asym

式中，μ为标量，分别取0，0.05，0.1和0.21.由上式计算出M_1asym并归一化 (刘超, 2011)，即可得到M₁中的4个非对称地震矩张量，分别用M₁⁽¹⁾, M₁⁽²⁾, M₁⁽³⁾和M₁⁽⁴⁾表示：

式中，M₁⁽¹⁾为对称地震矩张量.若考虑具有厚度的滑动弱化区地震断层模式 (刘超，2011；刘超，陈运泰，2014), 那么与较大力偶相联系的节面即为断层面，其矩为T₁; 与较小力偶相联系的节面即为辅助面，其矩为T₂.为了表示力偶矩大小的差异，定义标量λ为

λ反映的是在两个节面上，沿位错矢量方向力偶矩大小的差异.当λ=0时，T₁=T₂, 矩张量退化为对称矩张量.在此情形下，无法分辨出断层面和辅助面.随着反对称部分比重的增加，T₁逐渐增大，T₂逐渐减小，λ也随之增大.由式 (42) 可计算出M₁⁽¹⁾, M₁⁽²⁾, M₁⁽³⁾和M₁⁽⁴⁾对应的λ值分别为0, 0.09, 0.19和0.43.

类似地，在生成第二组 (M₂) 非对称地震矩张量时，其对称部分M_2sym和反对称部分M_{2anti_sym}分别为

μ分别取0, 0.05, 0.1和0.21.由式 (40) 可计算出M_2asym并归一化，即得到M₂中的4个非对称地震矩张量，分别用M₂⁽¹⁾, M₂⁽²⁾, M₂⁽³⁾和M₂⁽⁴⁾表示：

式中，M₂⁽¹⁾为对称地震矩张量.由式 (42) 可计算出M₂⁽¹⁾, M₂⁽²⁾, M₂⁽³⁾和M₂⁽⁴⁾对应的λ值分别为0, 0.06, 0.11和0.24.

数值试验选用了37个震中距处于30°与90°之间的虚拟台站 (图 1).对每一个台站，均采用归一化的、半周期为15 s的正弦函数作震源时间函数，与格林函数卷积合成地震图.对地震图分别加上0, 20%, 40%和80%的高斯噪声，再用4阶巴特沃斯 (Butterworth) 滤波器滤波，滤波频率上限为0.03 Hz，下限为0.002 Hz.截取滤波后合成地震图P波或S波到时前10 s至到时后90 s的数据，组成非对称地震矩张量反演数值试验的合成地震图.

图  1  数值试验中的台站分布

八角星表示震中，三角形表示地震台

Figure  1.  Station distribution in numerical tests

The star denotes earthquake epicenter, and triangles denote seismic stations

下载: 全尺寸图片 幻灯片

4.2   结果

4.2.1   M₁组中地震矩张量反演结果

以矩张量M₁⁽¹⁾为震源，反演结果如图 2左上图和图 3所示.由图 2左上图可知，在不同噪声强度水平下，反演结果的误差D随噪声增大而增加，其中以单用P波反演 (以下简称P波反演) 的误差最大，P波和S波联合反演 (以下简称联合反演) 的误差较小.由图 3可知:在40%噪声下，P波反演结果的误差均超过0.03，最高可达0.07；在80%噪声下，单用S波反演 (以下简称S波反演) 结果的误差和联合反演结果的误差才超过0.03，分别达到0.04和0.05. M₁⁽¹⁾为对称地震矩张量，其λ=0.反演得到的λ最大为0.02，相对较小，反演方法对λ的分辨较高，没有产生较大的误差. P波反演和S波反演分别得出的震源时间函数S^P(t) 和S^S(t) 均与S(t) 基本上一致, 准确稳定.

图  2  M₁的数值试验结果

Figure  2.  Results of numerical test for M₁

下载: 全尺寸图片 幻灯片

图  3  M₁⁽¹⁾的数值试验结果

第1列为“真实”的地震矩张量和震源时间函数S(t), 第2列至第5列为不同噪声强度水平下的反演结果；第1行为P波解，第2行为S波解，第3行为P波和S波联合反演解，第4行为震源时间函数S(t)

Figure  3.  Results of numerical test for M₁⁽¹⁾

Column 1 shows the 'real' seismic moment tensor and source time function S(t). Columns 2 to 5 show inversion results with different noise level. Rows 1 to 4 show the inverted seismic moment tensor solutions solely from P wave, S wave, and both P and S waves, and the source time function S(t), respectively

下载: 全尺寸图片 幻灯片

以矩张量M₁⁽²⁾为震源，反演结果如图 2右上图和图 4所示.由图 2右上图可知:在不同噪声强度水平下，反演结果的误差D随噪声增大而增加, 其中以P波反演的误差最大，联合反演的误差较小.由图 4可知:无论噪声大小，P波反演的误差均超过0.03，最小为0.05，最大可达0.10；S波和联合反演结果的误差在40%噪声下，达到0.03. M₁⁽²⁾为非对称地震矩张量，其λ=0.09. P波反演结果均未能正确地分辨出断层面，但S波反演和联合反演结果均能正确地分辨出断层面，λ值的误差随噪声增大而增加，在0，20%和40%噪声下λ值误差较小，在80%噪声下λ值有一定的误差. P波反演和S波反演分别得出的震源时间函数S^P(t) 和S^S(t) 均与S(t) 基本上一致, 准确稳定.

图  4  M₁⁽²⁾的数值试验结果

第1列为“真实”的地震矩张量和震源时间函数S(t), 第2列至第5列为不同噪声强度水平下的反演结果；第1行为P波解，第2行为S波解，第3行为P波和S波联合反演解，第4行为震源时间函数S(t)

Figure  4.  Results of numerical test for M₁ ⁽²⁾

Column 1 shows the 'real' seismic moment tensor and source time function S(t). Columns 2 to 5 show inversion results with different noise level. Rows 1 to 4 show the inverted seismic moment tensor solutions solely from P wave, S wave, and both P and S waves, and the source time function S(t), respectively

下载: 全尺寸图片 幻灯片

以矩张量M₁⁽³⁾为震源，反演结果如图 2左下图和图 5所示.由图 2左下图可知，在不同噪声强度水平下，反演结果的误差D大致随噪声增大而增加，其中以P波反演的误差最大，联合反演的误差较小.由图 5可知:无论噪声大小，P波反演结果的误差均超过0.03，最大可达0.16；S波反演结果的误差在40%噪声下，可达0.04, 在80%噪声下，又降至0.03；联合反演结果的误差在80%噪声下达到0.03. M₁⁽³⁾为非对称地震矩张量，其λ=0.19. P波反演结果均未能正确地分辨出断层面；但S波和联合反演结果均能正确分辨出断层面，λ值误差随噪声增大而增加，在0, 20%和40%噪声下λ值误差较小，在80%噪声下λ值的误差较大. P波反演和S波反演分别得出的震源时间函数S^P(t) 和S^S(t) 均与S(t) 基本上一致, 准确稳定.

图  5  M₁⁽³⁾的数值试验结果

第1列为“真实”的地震矩张量和震源时间函数S(t), 第2列至第5列为不同噪声强度水平下的反演结果；第1行为P波解，第2行为S波解，第3行为P波和S波联合反演解，第4行为震源时间函数S(t)

Figure  5.  Results of numerical test for M₁⁽³⁾

Column 1 shows the 'real' seismic moment tensor and source time function S(t). Columns 2 to 5 show inversion results with different noise level. Rows 1 to 4 show the inverted seismic moment tensor solutions solely from P wave, S wave, and both P and S waves, and the source time function S(t), respectively

下载: 全尺寸图片 幻灯片

以矩张量M₁⁽⁴⁾为震源，反演结果如图 2右下图和图 6所示.由图 2右下图可知，在不同噪声强度水平下，反演结果的误差D基本上随噪声增大而增加，其中以P波反演的误差最大，联合反演的误差较小.由图 6可见:无论噪声大小，P波反演结果的误差均超过0.03，为0.19；S波和联合反演结果的误差在40%噪声下均达到0.05. M₁⁽⁴⁾为非对称矩张量，其λ=0.43. P波反演结果均未能正确分辨出断层面，但S波反演结果和联合反演结果均能正确分辨出断层面，λ值误差大致随噪声增大而增加. P波反演和S波反演分别得出的震源时间函数S^P(t) 和S^S(t) 均与S(t) 基本上一致, 准确稳定.

图  6  M₁⁽⁴⁾的数值试验结果

第1列为“真实”的地震矩张量和震源时间函数S(t).第2列至第5列为不同噪声强度水平下的反演结果；第1行为P波解，第2行为S波解，第3行为P波和S波联合反演解，第4行为震源时间函数S(t)

Figure  6.  Results of numerical test for M₁⁽⁴⁾

Column 1 shows the 'real' seismic moment tensor and source time function S(t). Columns 2 to 5 show inversion results with different noise level. Rows 1 to 4 show the inverted seismic moment tensor solutions solely from P wave, S wave, and both P and S waves, and the source time function S(t), respectively

下载: 全尺寸图片 幻灯片

4.2.2   M₂组中地震矩张量的反演结果

以矩张量M₂⁽¹⁾为震源，反演结果如图 7左上图和图 8所示.由图 7左上图可知，在不同噪声强度水平下，反演结果的误差D大致随噪声增大而增加，其中以P波反演的误差最大, 联合反演的误差较小.由图 8可知，在20%噪声下，P波反演结果的误差超过0.03，可达0.11；在40%噪声下，S波反演的结果误差超过0.03，可达0.06；在40%噪声下，联合反演结果的误差可达0.03. M₂⁽¹⁾为对称矩张量，其λ=0.反演得到的λ值随噪声的增大而增加，在40%和80%噪声下，λ值误差相对较大. P波反演和S波反演分别得出的震源时间函数S^P(t) 和S^S(t) 均与S(t) 基本上一致, 准确稳定.

图  7  M₂的数值试验结果

Figure  7.  Results of numerical test for M₂

下载: 全尺寸图片 幻灯片

图  8  M₂⁽¹⁾的数值试验结果

第1列为“真实”的地震矩张量和震源时间函数S(t), 第2列至第5列为不同噪声强度水平下的反演结果；第1行为P波解，第2行为S波解，第3行为P波和S波联合反演解，第4行为震源时间函数S(t)

Figure  8.  Results of numerical test for M₂⁽¹⁾

Column 1 shows the 'real' seismic moment tensor and source time function S(t). Columns 2 to 5 show inversion results with different noise level. Rows 1 to 4 show the inverted seismic moment tensor solutions solely from P wave, S wave, and both P and S waves, and the source time function S(t), respectively

下载: 全尺寸图片 幻灯片

以矩张量M₂⁽²⁾为震源，反演结果如图 7右上图和图 9所示.由图 7右上图可知，在不同噪声强度水平下，反演结果的误差D随噪声增大而增加，其中以P波反演的误差最大，联合反演的误差较小.由图 9可知:无论噪声大小，P波反演结果的误差均超过0.03，最小为0.04，最大达0.39; S波反演和联合反演结果的误差在40%噪声下，均超过0.03. M₂⁽²⁾为非对称矩张量，其λ=0.06. P波反演结果均未能正确分辨出断层面.在无噪声情况下，S波反演和联合反演结果均能正确分辨出断层面；在20%噪声下，只有联合反演的结果能正确分辨出断层面，λ值误差随噪声增大而增加. P波反演和S波反演分别得出的震源时间函数S^P(t) 和S^S(t) 均与S(t) 基本上一致, 准确稳定.

图  9  M₂⁽²⁾的数值试验结果

第1列为“真实”的地震矩张量和震源时间函数S(t), 第2列到第5列为不同噪声强度水平下的反演结果；第1行为P波解，第2行为S波解，第3行为P波和S波联合反演解，第4行为震源时间函数S(t)

Figure  9.  Results of numerical test for M₂⁽²⁾

Column 1 shows the 'real' seismic moment tensor and source time function S(t). Columns 2 to 5 show inversion results with different noise level. Rows 1 to 4 show the inverted seismic moment tensor solutions solely from P wave, S wave, and both P and S waves, and the source time function S(t), respectively

下载: 全尺寸图片 幻灯片

以矩张量M₂⁽³⁾为震源，反演结果见图 7左下图和图 10.由图 7左下图可知，在不同噪声强度水平下，反演结果的误差D随噪声增大而增加, 其中以P波反演的误差最大，联合反演的误差较小.由图 10可知:无论噪声大小，P波反演结果的误差均超过0.03，最小为0.07，最大达0.41; S波反演结果误差在20%噪声下达到0.03；联合反演结果误差在80%噪声下超过0.03，可达0.07. M₂⁽³⁾为非对称矩张量，其λ=0.11. P波反演结果均未能正确分辨出断层面，但S波反演和联合反演结果均能正确分辨出断层面. P波反演和S波反演分别得出的震源时间函数S^P(t) 和S^S(t) 均与S(t) 基本上一致, 准确稳定.

图  10  M₂⁽³⁾的数值试验结果

第1列为“真实”的地震矩张量和震源时间函数S(t), 第2列到第4列为不同噪声下的反演结果；第1行为P波解，第2行为S波解，第3行为P波和S波联合反演解，第4行为震源时间函数S(t)

Figure  10.  Results of numerical test for M₂⁽³⁾

Column 1 shows the 'real' seismic moment tensor and source time function S(t). Columns 2 to 5 show inversion results with different noise level. Rows 1 to 4 show the inverted seismic moment tensor solutions solely from P wave, S wave, and both P and S waves, and the source time function S(t), respectively

下载: 全尺寸图片 幻灯片

以矩张量M₂⁽⁴⁾为震源, 反演结果如图 7右下图和图 11所示.由图 7右下图可知，在不同噪声强度水平下，误差随噪声增大而增加, P波反演的误差最大，联合反演的误差较小.由图 11可知:无论噪声大小, P波反演结果的误差超过0.03，最小为0.14, 最大为0.38; S波反演结果误差在20%噪声下超过0.03，达到0.05;联合反演结果误差在40%噪声下超过0.03，达到0.05. M₂⁽⁴⁾为非对称矩张量，其λ=0.24. P波反演结果均未能正确分辨出断层面，但S波反演和联合反演结果在噪声小于等于40%时，均能正确地分辨出断层面. P波反演和S波反演分别得出的震源时间函数S^P(t) 和S^S(t) 均与S(t) 基本上一致, 准确稳定.

图  11  M₂⁽⁴⁾的数值试验结果

第1列为“真实”的地震矩张量和震源时间函数S(t), 第2列到第5列为不同噪声下的反演结果；第1行为P波解，第2行为S波解，第3行为P波和S波联合反演解，第4行为震源时间函数S(t)

Figure  11.  Results of numerical test for M₂⁽⁴⁾

Column 1 shows the 'real' seismic moment tensor and source time function S(t). Columns 2 to 5 show inversion results with different noise level. Rows 1 to 4 show the inverted seismic moment tensor solutions solely from P wave, S wave, and both P and S waves, and the source time function S(t), respectively

下载: 全尺寸图片 幻灯片

4.3   非对称地震矩张量反演方法的可行性

以上通过数值试验验证了非对称地震矩张量反演方法的可行性.数值试验涉及具有代表性的非对称地震矩张量，包括以走滑为主、含有少量CLVD分量，以及以倾滑为主、含有较大CLVD分量两种情形，并且这两种情形下的矩张量均具有相等的对称部分和幅值不同的反对称部分.数值试验涉及范围广泛，试验结果具有相当大的适用性.由P波反演和S波反演分别得出的震源时间函数S^P(t) 和S^S(t) 均与S(t) 基本上一致, 准确稳定; 同时, S^P(t) 与S^S(t) 之间又有主要由破裂传播效应产生的微小差别, 有力地说明了反演结果的可靠性与反演方法的可行性.

观察P波反演结果误差与λ值的关系即可知，随着λ值的增大，P波反演结果的误差也增大.这是由于P波无法分辨非对称矩张量中M_xy和M_yx这两个分量所造成 (刘超, 2011).由P波反演得到的矩张量，其M_xy分量与M_yx分量始终相等.随着λ增大，矩张量的非对称性越大，M_xy与M_yx的差异也越大.但因P波无法分辨M_xy与M_yx，若在反演中始终令M_xy与M_yx相等，会造成结果误差也随λ的增大而增大.在S波反演以及联合反演中, 便不存在类似的问题.

除了对称地震矩张量M₁⁽¹⁾和M₂⁽¹⁾以外，M₂⁽²⁾对应的λ值最小，为0.06.当λ=0.06时，联合反演能在噪声小于等于20%的条件下正确地分辨出断层面；当λ≥0.09时, 由M₁⁽²⁾的结果可知, 联合反演能在噪声小于等于40%的条件下正确分辨出断层面.

这些数值试验结果表明: P波反演的误差最大，S波反演和联合反演的误差较小; S波携带的震源信息与P波不同，P波与S波联合反演可提高结果的准确性; 与P波反演和S波反演的结果比较，联合反演结果有较高的断层面判定正确率，充分说明了在非对称地震矩张量反演中，引入S波进行P波与S波联合反演的必要性.

5.   讨论与结论

本文研究了非对称地震矩张量时间域反演的理论与方法.研究表明, 虽然非对称地震矩张量反演与对称矩张量反演类似，但存在诸如是否过度拟合、如何定量表示地震矩张量之间的差异以及地震矩张量各分量的分辨等问题.为此, 本文引入了赤池信息准则 (AIC准则) 以判断非对称地震矩张量反演相对于对称地震矩张量反演，是否存在过度拟合的问题；引入了地震矩张量的矢量表示法以定量地描述地震矩张量之间的差异；通过分析格林函数与地震矩张量各个分量之间的关系, 得出需要同时运用垂直向与水平向地动位移数据进行联合反演才能区分M_xy与M_yx.为了检验非对称地震矩张量反演方法的可行性, 本文还利用合成地震图进行了数值试验, 验证了引入S波进行P波与S波联合反演对于提高非对称地震矩张量反演的准确性以及提高对断层面判定能力的必要性.这些研究结果对于运用实际观测资料反演非对称地震矩张量的工作均是很有益的参考.