以克罗地亚MS6.4 地震为例检验震前应力集中现象

孔 华， 万永革*

（1. 防灾科技学院， 河北 三河 065201； 2. 河北省地震动力学重点实验室， 河北 三河 065201）

0 引言

自上个世纪七八十年代开始， 国内外的一些地震学家便试图在震前小震与主震之间找到某种联系， 进而推动人们对发震机理的深入认识以及在地震预测等研究领域的实质性进展。震源机制解是地球物理学研究过程中重要的工具， 一方面它可以用来描述发震机制类型和发震构造[1-5]， 另一方面通过研究地震序列间震源机制的空间变化关系， 我们还能刻画出震源区应力场形态， 这些对我们加深发震机理的理解以及进行地震预测问题的研究探索均具有重要意义。 在研究海城地震序列时， 陈颙[6]最早提出了震源机制的一致性的概念， 希望可以将其当作判别大震发生的新参数。 但是因为当时大震前的震源机制资料有限以及作者忽略了断层的倾角和滑动角等影响因素而直接采用比较断层标准差， 因此该方法的正确性需要做进一步探讨。 在研究夏威夷Kaoiki 地区的震源机制变化时， 刁桂苓等[7]将大震前构造应力场的最大、 中等和最小主压应力轴与震源机制的P、 T 和B轴之间的夹角之和作为震源机制的一致性参数，然后采用千岛岛弧地区地震的矩心矩张量解（CMT） 分析了本地区震源机制的一致性特征，结果表明， MW≥7.5 地震之前都有一致性参数降低的现象[8]。 但是因为震源机制P、 T、 B 三轴的二义性导致他们得出的震源机制一致性参数具有不可忽视的振荡现象， 所以不可避免地给甄别震源机制参数一致性的工作带来了困难。 在研究1992 年LandersMW7.3 地震和1999 年Hector Mine MW7.1 地震时， 万永革[9]利用Kagan[10-11]提出的双力偶地震震源机制空间旋转的方法， 采用Hauksson[12]研究得到的美国南加州详尽的小震震源机制资料， 计算出目标地区10 年内地震震源机制相对主震震源机制平均空间旋转角的变化， 结果表明， 在两个大震发生前的一段时间里， 震前小震都出现了趋于主震震源机制的现象。 2020 年12 月29 日克罗地亚地区发生了MS6.4 地震， 该地震具有前震和余震， 并且均可以得到这些地震的震源机制， 这为研究本次地震前的应力集中及地震之后的应力变化提供了绝好的机会。 本研究以该地震序列为例进行研究， 为避免不同机构求解的震源机制解对分析结果的影响， 本研究利用万永革[13]提出的震源机制中心解思想， 采用Levenberg-Marquardt 方法求解出各个地震的震源机制中心解， 然后再利用通过计算两组震源机制解最小空间旋转角的方法， 比较得到目标地区震前、 主震和余震震源机制解的最小空间旋转角的差别， 以此来得到大震前后应力变化过程， 为大震孕育过程提供新的资料。

1 震源机制中心解的求解

北京时间2020 年12 月29 日19 时19 分（当地时间12 月29 日12 时19 分） 在克罗地亚发生MS6.4 地震。 震中位于（45.45°N， 16.21°E），震源深度10km（数据来自中国地震台网）。 该地震的震源机制结果表明， 破裂发生在一个接近垂直的断层上。 不论是从震中位置还是震源深度来看， 此次地震都属于发生在欧亚板块内的块内地震， 而不是沿着主要的构造板块边界的构造地震。 震中地处地中海沿岸， 是欧亚大陆与非洲的边界交汇地区， 涉及许多微板块和区域尺度结构的运动， 地质构造背景较为复杂。特别是该地震具有前震和余震， 网站（如EMSC 网站， https://www.emsc-csem.org/Earthqua ke/tensors.php?view=1） 给出了本次地震的主震震源机制、 震前两组地震的震源机制和两组余震的震源机制。 文中所选用的地震序列目录见表1， 其中包含MS6.0 以上的地震事件1 次、MS5.0～MS6.0 地震事件2 次、 MS5.0 以下的地震事件2 次。

表1 克罗地亚MS6.4 地震序列目录

1.1 震源机制中心解模型的建立与求解

震源机制在地壳应力场的分析中起到重要作用， 但是由于不同机构采用的计算方法以及参考资料的不同， 导致最终得到的同一个地震的震源机制之间具有相当的离散度。 为了在多个震源机制中找到一个最优的震源机制进行后续研究分析， 本文利用了万永革[13]提出的震源机制中心解模型进行数据处理。 该模型先是将某个震源机制定为初始解， 在该解的邻域内进行一阶Taylor 展开， 然后采用Levenberg-Marquardt方法迭代求解并同时渐进调整初始解， 最终与各备选震源机制的三维空间旋转角平方和最小的解作为目标地震的中心解用来进行后续的研究分析。 目前此模型已被多位学者应用到地震学研究当中。

1.2 主震震源机制中心解的求解

主震数据如表2 所示。 分别将各个机构（表2 第2 列） 测定的克罗地亚MS6.4 地震的震源机制作为初始解， 计算得到该地震的震源机制中心解（表2 第4 列） 和对应初始解的标准差（表2 第5 列）。 通过对比可以发现： 不论将哪个机构得到的震源机制作为初始解， 最终得到的震源机制中心解之间的差别是很小的， 并且对应得到的标准差只在小数点后3 位有一定的涨落。 由此表明， 以该方法得到的震源机制中心解是稳健可靠的。 即便如此， 我们依然以标准差最小的震源机制解作为最终解。 对于本次地震， 我们发现以INGV 测得的震源机制解为初始解得到的标准差最小， 将此中心解（第一个节面走向131.84°， 倾角78.01°， 滑动角178.53°； 第二个节面走向222.14， 倾角88.56°，滑动角11.99°） 当做最终结果（表2 第7 行）。根据走向、 倾角和滑动角与P、 T、 B 轴间的关系[17]， 我们便可以得到P 轴走向为356.37°， 倾伏角为7.42°， 不确定范围为346.96°～365.96°和0.68°～14.18°； T 轴走向为87.61°， 倾伏角为9.48°， 不确定范围为78.20°～97.20°和2.84°～16.21°； B 轴走向为228.88°， 倾伏角为77.92°，不确定范围为187.37°～264.63°和68.66°～87.30°。

表2 2020 年12 月29 日11∶20 克罗地亚MS6.4 地震震源机制解

不同机构和作者测定的震源机制和震源机制中心解的最小空间旋转角如表2 第6 列所示。通过比对可以发现， 该震源机制与中心解的最小空间旋转角的最小值为3.82°， 最大值为15.46°， 最小三维空间旋转角的标准差最小为9.7067°。 需要强调的一点是， 虽然以INGV 提供的震源机制解为初始解得到的中心解的标准差最小， 但是这并不意味着此解与实际震源机制中心解的空间距离最近， 此次地震距离实际中心解最近的应为USGS 机构。 造成此现象的根本原因是求解中心解是一个非线性问题， 采用雅克比矩阵求解各个参数的偏量有可能略微超过或者未到最优解， 表2 第4 列便说明了该现象。 这也是本文用不同机构测得的震源机制解作为初始解逐一反演的原因。 将最终确定的中心震源机制解及其不确定性和空间三维辐射花样绘为图1。

图1 2020 年12 月29 日克罗地亚MS6.4 地震的中心震源机制解（a）和空间三维辐射花样（b）Fig.1 Central focal mechanism（a）and 3-D radiation pattern（b）of the 11∶20 Croatian earthquake MS6.4 on December 29， 2020

1.3 两组前震震源机制中心解的结果展示

1.3.1 12 月28 日05∶28 地震震源机制中心解的结果展示

将发震时间为05∶28 的地震的震源机制解整理为表3， 对于本次地震， 我们发现以GFZ测得的震源机制解为初始解得到的标准差最小，将此中心解（第一个节面走向63.97°， 倾角61.31°， 滑动角17.28°； 第二个节面走向325.47°，倾角74.90°， 滑动角150.19°） 当做最终结果（表3 第7 行）。 根据走向、 倾角和滑动角与P、T、 B 轴间的关系， 便可以得到P 轴走向为17.08°, 倾伏角为8.93°， 不确定范围为-0.92°～35.08°和-3.89°～22.17°； T 轴走向为281.54°, 倾伏角为31.58°， 不确定范围为263.54°～299.54°和18.92°～43.79°； B 轴走向为121.02°， 倾伏角为56.90°， 不确定范围为94.77°～155.47°和44.02°～72.87°。

不同机构和作者测定的震源机制和震源机制中心解的最小空间旋转角如表3 第6 列所示。通过比对可以发现， 该震源机制与中心解的最小空间旋转角的最小值为15.41°， 最大值为20.07°， 最小三维空间旋转角的标准差最小为18.5008°。 将最终确定的中心震源机制解及其不确定性和空间三维辐射花样绘为图2。

图2 2020 年12 月28 日05∶28 克罗地亚MS5.2 地震的中心震源机制解（a）和空间三维辐射花样（b）Fig.2 Central focal mechanism（a）and 3-D radiation pattern（b）of the 05∶28 Croatian earthquake MS5.2 on December 28， 2020

表3 2020 年12 月28 日05∶28 克罗地亚MS5.2 地震震源机制解

1.3.2 12 月28 日06∶49 地震震源机制中心解的结果展示

将发震时间为06∶49 的地震的震源机制解整理为表4， 对于本次地震,我们发现以INGV测得的震源机制解为初始解得到的标准差最小，将此中心解（第一个节面走向144.02°， 倾角86.49°， 滑动角-167.00°； 第二个节面走向53.21°，倾角77.02°， 滑动角-3.60°） 当做最终结果（表4 第7 行）。 根据走向、 倾角和滑动角与P、 T、B 轴间的关系， 我们便可以得到P 轴走向为9.31°, 倾伏角为11.66°， 不确定范围为7.23°～12.23°和9.61°～13.74°； T 轴走向为277.94°, 倾伏角为6.64°， 不确定范围为275.85°～280.85°和4.54°～8.14°； B 轴走向为158.85°, 倾伏角为76.54°， 不确定范围为146.57°～166.89°和74.54°～79.37°。

表4 2020 年12 月28 日06∶49 克罗地亚MS5.0 地震震源机制解

不同机构和作者测定的震源机制和震源机制中心解的最小空间旋转角如表4 第6 列所示。通过比对可以发现， 该震源机制与中心解的最小空间旋转角的最小值与最大值均为3.04°， 最小三维空间旋转角的标准差最小为3.0414°。 将最终确定的中心震源机制解及其不确定性和空间三维辐射花样绘如图3 所示。

图3 2020 年12 月28 日06∶49 克罗地亚MS5.0 地震的中心震源机制解（a）和空间三维辐射花样（b）Fig.3 Central focal mechanism（a）and 3-D radiation pattern（b）of the 06∶49 Croatian earthquake MS5.0 on December 28， 2020

1.4 两组余震震源机制中心解的结果展示

1.4.1 12 月30 日05∶15 地震震源机制中心解结果展示

将发震时间为05∶15 的地震的震源机制解整理为表5， 对于本次地震,我们发现以GFZ 测得的震源机制解为初始解得到的标准差最小，将此中心解（第一个节面走向54.14°， 倾角52.75°， 滑动角9.11°； 第二个节面走向318.60°，倾角82.76°， 滑动角142.40°） 当做最终结果（表5 第7 行）。 根据走向、 倾角和滑动角与P、T、 B 轴间的关系， 我们便可以得到P 轴走向为12.05°， 倾伏角为19.81°， 不确定范围为-23.27°～47.73°和-2.96°～44.93°； T 轴走向为269.48°，倾伏角为31.14°， 不确定范围为234.15°～305.15°和6.74°～54.73°； B 轴走向为129.31°， 倾伏角为51.81°， 不确定范围为63.89°～183.27°和22.86～86.42°。

不同机构和作者测定的震源机制和震源机制中心解的最小空间旋转角如表5 第6 列所示。通过比对可以发现， 该震源机制与中心解的最小空间旋转角的最小值为36.10°， 最大值为36.22°， 最小三维空间旋转角的标准差最小为36.1627°。 将最终确定的中心震源机制解及其不确定性和空间三维辐射花样绘如图4 所示。

表5 2020 年12 月30 日05∶15 克罗地亚MS4.8 地震震源机制解

图4 2020 年12 月30 日05∶15 克罗地亚MS4.8 地震的中心震源机制解（a） 和空间三维辐射花样（b）Fig.4 Central focal mechanism（a） and 3-D radiation pattern（b） of the 05∶15 Croatian earthquake MS4.8 on December 30， 2020

1.4.2 12 月30 日05∶26 地震震源机制中心解结果展示

将发震时间为05∶26 的地震的震源机制解整理为表6， 对于本次地震， 我们发现以GFZ测得的震源机制解为初始解得到的标准差最小，将此中心解（第一个节面走向344.32°， 倾角83.01°，滑动角176.93°； 第二个节面走向74.70°， 倾角86.95°， 滑动角7.00°） 当做最终结果（表6 第7行）。 根据走向、 倾角和滑动角与P、 T、 B 轴间的关系， 我们便可以得到P 轴走向为209.34°， 倾伏角为2.78°， 不确定范围为191.49°～227.49°和-10.03°～15.58°； T 轴走向为299.68°， 倾伏角为7.11°， 不确定范围为281.84°～317.84°和-5.46°～19.93°； B 轴走向为98.13°，倾伏角为82.36°， 不确定范围为8.65°～238.81°和77.81°～99.75°。

不同机构和作者测定的震源机制和震源机制中心解的最小空间旋转角如表6 第6 列所示。通过比对可以发现， 该地震震源机制解与中心解的最小空间旋转角的最小值为36.10°， 最大值为36.22°， 最小三维空间旋转角的标准差最小为36.1627°。 将最终确定的中心震源机制解及其不确定性和空间三维辐射花样绘为图5。

表6 2020 年12 月30 日05∶26 克罗地亚MS4.7 地震震源机制解

图5 2020 年12 月30 日05∶26 克罗地亚MS4.7 地震的中心震源机制解（a）和空间三维辐射花样（b）Fig.5 Central focal mechanism（a）and 3-D radiation pattern（b）of the 05∶26 Croatian earthquake MS4.7 on December 30， 2020

2 地震序列间震源机制中心解最小空间旋转角的求解

我们把因震源机制P、 T、 B 轴二义性而产生的四个待定解中的最小值称为最小空间旋转角[13]。 通过比较两个震源机制中心解间的最小空间旋转角便可以定量得出震源机制中心解间的差别——其值越小越表明两个地震间的差别越小， 反之亦然。

2.1 震前-主震/震后-主震的震源机制最小空间旋转角的求解

按各地震的发震时间（表7 第1 列） 的先后顺序， 将第二部分得到的两组前震和两组余震的震源机制中心解（表7 第2 列） 以及主震震源机制中心解（表7 第3 列） 整理成表7 并计算得到各个地震的中心解与主震中心解的最小空间旋转角（表7 第4 列）。

表7 震前-主震/余震-主震的震源机制最小空间旋转角

通过研究我们可以发现： 震前一段时间内，该地发生的前震震源机制中心解有靠近主震震源机制中心解的趋势， 但是余震与主震的中心解的最小空间旋转角较前者又有增大的现象。

2.2 震前-震前/震后-震后的震源机制最小空间旋转角的求解

以主震发生时间为界， 将震前的两组中心解（表8 第2 行第2 列、 第3 列） 与两组余震中心解（表8 第3 行第2 列、 第3 列） 整理成表8 并计算得到最小空间旋转角（表8 第4 列）。

表8 震前-震前/震后-震后的震源机制最小空间旋转角

通过对比数据可以发现， 震前两组中心解的最小空间旋转角明显小于两组余震的中心解最小空间旋转角。

3 结论与讨论

本文先是求解了2020 年克罗地亚MS6.4 主震及震前、 余震序列的震源机制中心解， 然后在此基础上又计算了震前/震后-主震和震前-震前/震后-震后两组中心解的最小空间旋转角。当前学者普遍认为， 一个地震发生前震源区会先后处于地应力积累、 集中和加强的阶段， 应变能在应力集中区的突然释放是地震发生的直接原因[18]。 我们通过比较最小空间旋转角发现， 在主震发生前的一段时间内， 该地区发生的地震之间的震源机制中心解有趋近主震震源机制中心解的趋势； 主震过后， 余震的震源机制中心解和主震震源机制中心解则离散开来。这一结果说明了震前因应力集中导致该地区出现了靠近主震机制的优势应力主轴， 震后由于地下积累应变能的突然释放， 震源区的应力水平明显下降， 进而使得应力分布趋向暂时的平衡而进入到下一个应力集中和应变能积累期，本文的研究结果也证明了该假说的正确性。