(1)北京大学学报(自然科学版) 第62卷 第3期 2026年5月
Acta Scientiarum Naturalium Universitatis Pekinensis, Vol. 62, No. 3 (May 2026)
doi: 10.13209/j.0479-8023.2025.111
国家重点研发计划(2024YFF0506800)和国家自然科学基金(T2122012)资助
收稿日期: 2025–03–25;
修回日期:2025–07–24
摘要 海洋地震发生后可通过在海底布置压强传感器来记录地震引发的压强变化, 目前针对地震引发的压强变化还没有成熟的应用。本文采用 SPECFEM3D 程序构建海底地震模型, 模拟地震破裂过程以及由此引发的海底压强变化, 并分析海底压强信号的分布特征。研究结果表明, 在距离震源区较远时, 地震引起的低频压强信号变化 p 与海底垂直加速度 a 满足理论关系 p=rha, 距离震源较近时, 则不满足理论关系, 两者的边界与震源深度相关。在固定其他震源参数(如震级和滑动角)的情况下, 震源深度决定海底永久形变的大小, 因此低频压强变化特征对预测地震引发的海啸有指示意义。另外, 在震源区一定范围外可以记录到瑞利波, 并且瑞利波引发的海底压强变化相对于体波存在放大现象; 出现瑞利波放大现象的台站位置也与震源深度相关, 可以为海啸预警提供参考。
关键词 地震海啸生成过程; 海底永久形变; 海底压强; 瑞利波
作为全球极具有破坏性的自然灾害, 海啸时刻对沿海城市构成巨大的威胁, 而海底地震是引发海啸的主要原因, 占全球海啸的 90%以上[1]。为满足海洋地球物理学的研究需求, 海底地震观测系统已经在全球广泛部署, 例如位于日本东部的 DONET (Dense Ocean floor Network System for Earthquake and Tsunamis)[2–3]和 SNET (Seafloor Observation Net-work for Earthquakes and Tsunamis)[4–5]等。同时, 海底压强计(Bottom Pressure Recorder, BPR)在监测海啸生成和传播中发挥越来越重要的作用。美国国家海洋和大气管理局(National Oceanic and Atmosphe-ric Administration, NOAA)开发的用于全球海啸预警的 DART (Deep-ocean Assessment and Reporting of Tsunami)系统中的传感器正是位于海底的压强计。日本的 DONET 系统中, 每个台站除布置地震仪外, 也布置海底压强计。传统的海啸预警策略使用地震信号来反演约束震源参数, 再预测引发的海啸 波[6–7]。然而, 使用传统的地震学观测数据反演得到的震源参数有一定的误差: 远场地震波数据对断层上的永久位错敏感性不足[8], 而观测地表位移的GNSS 等台站往往位于陆地上, 导致对海洋下面断层位错的约束欠佳, 因此传统的震源模型预测海底的永久形变时存在明显的误差。例如, An 等[9]研究2014 年智利 Iquique 地震海啸, 发现只使用地震波反演得到的震源模型在浅部和深部同时存在位错, 预测出双峰的海啸波, 与实际观测不符。海洋中的压强观测仪器往往布放在潜在震源区附近的海底, 能够弥补传统地震学观测对海洋覆盖欠缺的问题。因此, 研究地震引发海底压强变化的特征, 建立海底压强与震源参数之间的关联, 对地震震源及海啸预警研究有积极的意义。
海底压强计不仅能记录潮汐和海啸波等水波引起的海面高度变化导致的海底压强变化, 还能记录地震波引起的海底运动导致的海底压强变化。这意味着海底地震发生后, 压强计能提供整个地震期间的动态变化和静态压力信号。最早的海底压强观测仪器 1979 年由 Fillloux 设计, 称为压差式压强计(di-fferential pressure gauge, DPG)。1979 年3 月 14 日墨西哥西南侧海岸发生 Mw 7.6 级地震, 引起微小的海啸波, Filloux[10]从震中距 1000km 外的海底压强计观测数据中观察到地震波之后传播的海啸波, 首次证实海底压强计能够记录并识别远距离传播的海啸波, 为海啸监测提供了新的观测手段; 同时指出, 如果瑞利波周期较长, 则瑞利波引起的海底压强 p与海底运动的垂直加速度 a 之间满足 p=rha (𝜌 为海水密度, h 为水深)。然而, Filloux 未对该关系式进行严格的推导。An 等[11]在忽略海水可压缩性并假设海底是刚性边界条件的情况下推导海底压强的显式理论解。在地震波频率远低于海水的共振频率时, 该理论解可以简化为 p=rha。An 等[11]还采用 3 个地震事件的观测数据对 p=rha 进行验证, 发现在低频区间 0.02～0.2Hz, 这一关系式能够很好地满足, 但在高频区间无法与观测数据吻合, 他们认为这是因为在推导时没有考虑海水的可压缩性造成的。胡晨彤[12]进一步考虑海水的可压缩性, 推导海底压强的理论解, 并与实际观测数据相比较, 发现在低频时两者能够很好地吻合, 但在高频区间上仍然存在较大的相位差和振幅差, 认为是没有考虑地球结构的弹性以及地震波不完全是垂直入射海底造成的。Deng 等[13]考虑海水可压缩性, 并假设海底运动是由入射地震波引起, 同时考虑地震波在水中的一系列反透射声波对海底的共同作用, 得到海底压强与海底垂直加速度之间的理论关系。他们采用实际观测数据验证该理论关系, 发现在共振频率以下, 该理论关系能够与观测数据很好地吻合。Saito[14]在考虑海水可压缩性和地球弹性后, 出现低频海底压强所满足的 p=rha 关系式仅在海底光滑时成立, 在海底存在剧烈变形时不成立。还有一些研究指出, 在海水的共振频率以上, 海底压强与海底垂直加速度之间满足p=rvc0 (c0 为海水中声波波速, v 为海底运动的速度)[15–17]。然而, 这一关系式仅考虑地震波在海底的初次入射和透射, 无法在后续时间段的海底运动中成立, 在实际观测数据中也无法复现地震波首波之后的海底运动[14]。
其他研究中也推导过海底压强与海底运动加速度之间的理论关系[5,18–19], 该理论关系一方面可用于修正压差式压强计记录的压强信号[11,13], 另一方面可以通过压强信号得到海底运动的加速度或位移, 进而反演震源并推测地震所致海底永久形变的幅度, 从而实现更快速准确的海啸预警。例如, 针对地震引发的海底压强变化, Sun 等[20]从海底压强计和海底地震仪数据中提取2011 年日本 Tohoku-Oki 地震浅层断层面上的滑动以及海底永久形变信息。Lay 等[21]利用 NOAA 深海压强计和智利 DART系统, 并结合断层滑移模型, 分析 2014 年 4 月 1 日智利 Mw 8.1 地震所致海底永久形变的分布特征。Nosov 等[22]对 2003 年日本Tokachi-Oki 海啸区域建立三维有限差分数值模型, 模拟震源区海啸生成的动态过程, 通过与震源区两个海底压强计观测数据的对比, 预测海底变形的幅度、变形速度和持续时间。Song 等[23]利用地震波与海啸波频率分布的差别以及低频海底压强满足 p=rha 关系式, 从震源区地震波和海啸波耦合的海底压强信号中提取纯净的海啸波信号, 通过海啸波波幅判断地震生成的海啸危害程度, 即使对近岸地震所致海啸, 可以提供有效的快速预警手段。
目前对震源约束的方法有很多, 仅对震源深度的判断就有传统的依据地震波走时的震源深度定位方法[24–27]、震相分析法(depth phase analysis)[28]和sPL 深度震相法[29], 以及近几年提出的虚拟震源地震探测方法[30]、CAP (Cut-and-Paste)方法[31–32]和TDMT 方法(time domain moment tensor inversion te-chnique)[33–34]等。然而, 地震学方法对震源参数的约束有一定的误差[35–36], 制约震后进行精确的海啸预警。2004 年印尼海啸后, 为满足快速海啸预警的需求, 海底压强计在全球广泛布设, 如 DART 海啸预警系统、DONET 和 SNET 等。DART系统通过观测海底压强, 实现对海啸波的实时监测, DONET 和SNET 等系统在潜在震源区附近布设海底压强计, 弥补了传统观测缺少海洋中观测数据的不足。
为了从海底压强计记录的压强变化中得到更多关于震源和海底永久形变的信息, 本文采用 SPEC-FEM3D 程序建立三维地震海啸生成模型, 模拟不同震源工况下地震海啸生成过程, 统计分析海底各台站的压强、位移和加速度分布特征。同时, 探讨低频时海底压强与地震所致海底永久形变的关系及压强信号出现瑞利波放大现象的原因, 以期对海底地震引发的海啸进行及时的准确预警提供参考。
为了模拟真实的地震海啸生成过程, 我们在直角坐标系中建立三维地震海啸生成模型(图 1)。假设海水为微可压缩流体, 海水初始时刻的密度为r0, 初始时刻的压强p0=-r0gz(g 为重力加速度), 微可压海水具有如下形式的状态方程:
(1)
其中, K 是海水的体积模量, p'为地震发生后引起的海水动压, r'表示海水的密度变化, p'和r'相比 p0 和r0 是一阶小量。这样, 海水的压强 p 和密度r分解为静压 p0 和动压 p'以及初始密度r0 和密度变化r'两部分, 即 p=p0+p'=-r0gz + p', r=r0+r'。将其带入流体的运动方程(式(2))并忽略高阶小量, 可将流体中重力项抵消, 得到动压的显式表达(式(3)):
图1 地震海啸生成模型
Fig. 1 Numerical model of earthquake and tsunami generation
(2)
(3)
其中, f是流体的速度势(即 v=Ñf)。计算得到速度势f后, 通过式(3)计算得到流体中的压强 p'。将式(3)带入流体的连续性方程, 可得
(4)联立式(1), (3)和(4), 可得流体中求解的波动方程:
(5)
其中, 
表示流体中声波(P 波)的波速。
模型假设下层固体为均匀的各向同性弹性体。在海啸生成过程中, 地震波是弹性波, 地震波的回复力是地球或海水的弹性力, 而海啸波是重力波, 其回复力是海水的重力。两者不仅回复力不同, 而且波速、周期和波长等物理特征也有显著差异。通常认为, 地震过程和海啸过程是解耦的, 因此可以将地震波和海啸波分开考虑[37–40]。海啸波的求解是水动力学中的常见问题, 通常采用浅水波方程或Boussinesq 方程求解, 而地震波通过波动方程求解。同理, 弹性固体中重力场不会影响地震波波动场, 地震波波动场也不会改变重力场。因此, 我们不考虑固体中的重力作用, 在下层地球结构中求解如下两个地震波方程:
(6)
(7)
其中, Φ 和 Ψ 分别是地震波 P 波和 S 波的势函数, λ 和μ 是两个弹性常数, rs 是弹性固体的密度。
在自由表面上不进行海啸波模拟, 边界条件设为表面压强等于大气压(设参考值为 0), 则由式(3)可得
(8)在模型的流固界面上, 需满足两侧法向速度连续和应力连续。在模型前、后、左、右以及下边界面上, 分别施加消波边界条件, 以便尽可能地吸收边界反射波。
采用 SPECFEM3D 程序[41–42]求解图 1 所示地震海啸生成模型的定解问题。在 SPECFEM3D 程序中设置模型 x 方向长度为 1380km, y 方向长度为 1080km, 上层海水水深为 3km, 上层海水密度r0=1.0× 103kg/m3, 下层固体密度rs=3.2×103kg/m3, 海水中 P波波速为 1500m/s, 固体中 P 波波速为 6500m/s, S波波速为 3000m/s。模拟采用的网格大小为 3km, 时间步长为 0.04s, 每个算例模拟总时长为 320s。采用的震源为单一点源, 震源位置在海底的投影位于海底面中心, 在图 1 中的坐标为 x= 690km, y=540km。将台站分别取在海底和海面上, 在海底, x 方向从 120 到 1260km, y 方向从 270 到810km, 每隔 15km 设置一个台站。在海底中心线y= 540km 上, x 方向从 120 到 1260km, 每隔 3km 设置一个台站。同时, 对震源区内的海底台站进行加密, 从 x 方向 645到 735km, y 方向 495 到 585km, 每隔 3km 设置一个台站。为了俘获自由表面的永久形变, 在自由表面上, x 方向从 240 到 1050km, y 方向从 360 到 720km, 每隔 3km 设置一个台站。
为了确保得到正确的数值模拟结果, 首先将采用 SPECFEM3D 程序模拟出的海底永久形变和 Oka-da 理论解[43]进行对比。以单一点源模拟的震级为7.3 级逆冲断层地震为例, 震源深度为 10.5km, 走向为 180º, 倾角为 15º, 滑动角为 90º。图 2 展示地震引发的中心线上垂直海底永久形变, 可见数值模拟结果与 Okada 理论解吻合得很好, 由此验证了数值模型的准确性。
图2 震级为 Mw 7.3 的点源地震引起的中心线上的垂直海底永久形变与 Okada 理论解的对比
Fig. 2 Comparison of the vertical seafloor permanent deformation along the centerline induced by a point source of magnitude Mw 7.3 with the Okada theoretical solution
地震发生后, 在远离震源区之外, 地震波引起的海底压强信号在低频时满足 p=rha[11,13,38]。通过数值模拟, 发现当震源深度较小时, 在离震源较近的区域, 低频海底压强不满足 p=rha。低频指远低于海水水层共振的频率(即 c/4h, 其中 c 为水中声波波速)[44]。对于本研究使用的数值模型, 水深为 3km, 水层的共振频率为 0.125Hz。我们以一个 Mw 7.3 级, 倾角为 15º, 滑动角为 90º, 深度为 4.5km 的点源地震为例, 给出满足和不满足 p=rha 的台站观测示例。图 3(a)展示位于震源正上方海底中心线上震源左侧 9km 处台站的海底压强记录, 滤波区间为0.01～0.02Hz, 使用海底垂直加速度 a 预测的海底压强。可以发现, 由于台站离震源区较远, 两者吻合得很好。图 3(b)展示位于震源正上方海底中心线上震源左侧 3km 处台站的海底压强记录, 滤波区间同样为 0.01～0.02Hz, 因为离震源区较近, 海底压强变化不满足 p=rha。
我们模拟并统计不同震源深度、不同震级和不同滑动角的地震工况下, 海底低频压强是否满足 p= rha 的分布特征。首先, 针对点源地震, 固定震级为 Mw 7.3、倾角为 15º 以及滑动角为 90º, 改变点源的深度分别为 4.5, 10.5, 19.5, 30.5 和 60.5km, 统计海底低频(0.01～0.02Hz)压强满足 p=rha 的范围, 结果如图 4 所示, 其中黑色圆点为海底低频压强不满足 p=rha 的位置, 其他位置则满足 p=rha。可以看到, 随着震源深度增加, 海底低频压强不满足 p= rha 的台站逐渐减少并趋于消失。点源深度为 60.5km 时, 海底已不存在低频压强不满足 p=rha 的台站。因此, 海底低频压强不满足 p=rha 的台站数量或范围可以大致反映震源深度。由于震源深度与地震引起的海底永久形变幅度直接相关, 所以海底低频压强不满足 p=rha 的台站数量或范围也可以反映海底(或海面)永久形变的分布。从 5 个地震工况对应的海面永久形变分布(图 4), 可以看到, 随着震源深度增加, 具有较大波高的海面永久形变分布趋于减小, 与海底低频压强不满足 p=rha 的台站数量或范围变化趋势一致。
图3 地震引发海底压强变化与p=rha的波形对比
Fig. 3 Comparison of waveforms of seafloor pressure variations induced by the earthquake with p=rha
为了进一步探讨海底低频压强是否满足关系式p=rha 及其与其他震源参数的关系, 我们统计在同一震源深度, 不同地震工况下海底低频压强不满足p=rha 的台站分布范围。固定震源深度为 10.5km, 倾角为 15º, 分别改变震级和滑动角, 统计结果如图5 所示。统计结果表明, 在其他震源参数固定的情况下, 仅改变地震震级, 海底低频压强不满足 p= rha 的台站分布范围不会发生变化, 如图 5 中同一行 3 个图的对比。仅改变滑动角时, 也不会明显改变海底低频压强不满足 p=rha 的台站分布范围, 如图 5 中同一列 3 个图的对比。只有当滑动角为 0º 时, 在垂直于断层滑动方向上存在较长距离的低频压强不满足 p=rha 的台站分布, 这排台站从震源所在的位置向左右两侧各延伸约 200km。在不同滑动角情况下, 海底低频压强不满足 p=rha 的台站分布位置的统计结果同样表明, 只有在滑动角为 0º(表示震源为走滑断层)时才出现如图 5(g), (h)和(i)所示在垂直于断层面滑动方向上存在的较长距离的海底低频压强不满足 p=rha 的台站分布。
综上所述, 根据海底低频压强是否满足 p=rha, 一方面可以大致判断震源区在海底的位置, 另一方面, 海底低频压强不满足 p=rha 的台站数目或范围可以用于指示震源深度。由于震源深度与地震所致海底永久形变的分布直接相关, 在震源深度较小时, 海底低频压强不满足 p=rha 的台站分布范围大致与海底永久形变的区域对应。地震所致海底永久形变范围越大, 潜在海啸的可能性越大, 因此海底低频压强关系式 p=rha 能够在海啸预警中发挥指导作用。在实际应用中, 需要在同一个海底观测点安装压强和加速度观测仪器。例如, 日本的 SNET 台站都包含两个海底压强计以及 4 个海底地震仪[4–5], DONET 系统的台站包括海底压强计、温度计、水听器以及地震仪等[2–3]。每个台站的海底压强计和海底地震仪可以分别测得该位置处的海底压强变化和海底运动的垂向加速度, 用于检验是否满足 p= rha, 从而使其发挥对海啸预警的指导作用。
利用 SPECFEM3D 程序模拟地震海啸生成模型, 当海底台站距离震源区一定范围后, 相较于垂直方向的加速度, 压强信号中的瑞利波出现放大的信号。以震级为 Mw 7.3, 倾角为 15º, 滑动角为90º, 震源深度为 10.5km 点源地震为例, 位于点源左侧, 与点源水平距离为 30, 150 和 270km的 3 个海底台站(S181, S141 和 S101)的垂直加速度和海底压强的模拟结果如图 6 所示, 可以看出, 瑞利波到时晚于 S波并且出现频散, 相比于加速度信号, 压强信号中被放大的部分为 S 波之后逐渐形成的瑞利波。距离震源较远的 S141 台站的压强信号在约 80s 之后出现放大的瑞利波波形(图 6(e)); 距离震源更远的S101 台站的压强信号在约 100s 之后出现被放大的瑞利波波形, 这部分瑞利波波形随着到震源距离的增加而持续延长(图 6(f)); 在距离点源较近的 S181台站, 相比于垂直加速度信号, 压强信号全程未出现被放大的瑞利波信号(图 6(d))。
黑色圆点代表海底低频压强不满足p=rha的台站位置
图4 震源深度分别为 4.5, 10.5, 19.5, 30.5 和 60.5km 时, 海底低频压强不满足p=rha的台站分布
Fig. 4 Station locations where the seafloor low-frequency pressure does not satisfy p=rha for source depths of 4.5, 10.5, 19.5, 30.5, and 60.5 km
黑色圆点代表海底低频压强不满足p=rha的台站位置
图5 震源深度为 10.5km, 倾角为 15º 时, 不同震级以及不同滑动角的海底低频压强不满足p=rha的台站分布
Fig. 5 Station locations where the seafloor low-frequency pressure does not satisfy p=rhafor a source depth of 10.5 km and a dip angle of 15º, with different magnitudes and rakes
我们统计了震级为 Mw 7.3, 倾角为 15º, 滑动角为 90º 的点源地震在不同震源深度下, 海底压强中出现放大的瑞利波的台站分布区域, 发现随着震源深度增加, 海底压强开始出现放大的瑞利波的台站范围逐渐从震源所在位置向外延展。当震源深度增至约 20km 时, 台站压强信号中已看不到被放大的瑞利波。图 7 给出点源深度分别为 4.5, 10.5 和19.5km 时海底压强中出现放大的瑞利波的台站分布。可以发现, 压强信号中存在放大的瑞利波的范围与震源深度相关, 震源深度越小, 离震源越近的位置越可能出现放大的瑞利波。
图6 海底台站压强信号出现相比于垂直加速度信号被放大的瑞利波的波形对比
Fig. 6 Comparison of waveforms showing the amplification of Rayleigh waves in seafloor pressure signals compared with the vertical acceleration signals
黑色点代表压强信号出现相比于垂直加速度信号被放大的瑞利波的台站位置
图7 震源深度分别为 4.5, 10.5 和 19.5km 时, 海底台站压强信号出现相比于垂直加速度信号被放大的瑞利波台站分布
Fig. 7 Station locations where Rayleigh waves in seafloor pressure signals are amplified compared with the vertical acceleration signals for source depths of 4.5, 10.5 and 19.5 km
下面分析当震源深度固定时, 压强信号出现被放大的瑞利波与其他震源参数的关系。当固定震源深度为 10.5km, 倾角为 15º, 滑动角为 90º 时, 分别改变点源的震级为 Mw 7.3, Mw 7.6 和 Mw 8.0, 发现压强信号中存在被放大的瑞利波在海底的分布位置不变。当固定震源深度为 10.5km, 倾角为 15º, 震级为 Mw 7.3 时, 分别改变点源的滑动角为 90º, 60º, 30º, 15º和 0º, 发现海底压强信号中存在放大的瑞利波的位置基本上一致, 只有在滑动角为 0º 时, 在海底垂直截面中心线上(即垂直于断层面滑动方向)台站的压强信号中没有出现放大的瑞利波。当固定震级为 Mw 7.3, 倾角为 15º, 滑动角为 0º 时, 逐渐增加震源的深度, 仍然可以发现海底压强中开始出现放大的瑞利波的台站范围逐渐从震源所在位置向外延展, 直至震源深度增加至约 20km 时, 压强信号上才看不到被放大的瑞利波。
综上所述, 当震源深度较小时, 离开震源一定距离后, 海底压强信号会出现相比于垂直加速度信号被放大的瑞利波, 随着震源深度增加, 被放大的瑞利波的出现范围逐渐从震源区向外延展, 直至消失, 因此压强信号中出现的放大的瑞利波可以用于指示震源深度。震源深度越小, 近源区海底压强信号中越早出现放大的瑞利波, 意味着海底可能存在较大的永久形变, 地震激发海啸的可能性较大。
真实地震的震源破裂情形往往为有限断层面上的一部分区域存在较明显的滑移量, 只有当震源到观测点的距离远大于有限断层的尺度时, 震源才被视为点源。因此, 有必要研究震源为有限断层面时海底低频压强不满足 p=rha 区域的特征。
为此, 我们设计两个 Mw 7.4 和两个 Mw 8.6 地震的有限断层面(图 8), 其中有限断层面的走向均为 180º, 沿图 1 中模型 y轴的负方向, 倾角均为 20º, 滑动角均为 90º。断层面总长度和宽度分别为 140和 200km, 整个断层面被划分为 16940 个子断层面, 每个子断层面的大小为 1.8km×1.8km。每个子断层面由位于中心的点源近似, 同时保证整个断层面上每个子断层面矩震级的加和满足给定的矩震级。为了预测海啸, 往往需要对断层面上复杂的位错分布进行简化, 通过不同的标度关系, 建立断层面上的破裂面积、破裂区域的长宽比与地震震级之间的定量关系[45–48]。本文采用 An 等[45]提出的标度关系式来确定断层面上的破裂面积:
(9)其中, S 表示断层面上能够激发海啸的具有显著滑移量的破裂面积(m2), 而不是真实的具有滑移量的破裂区域; M0为地震矩(Nm)。图 8 中 Mw 7.4 和 Mw 8.6 地震断层面上的破裂面积分别为 30km×30km 和115km×115km。我们假设断层面上的破裂从破裂面积的中心开始向外以 1.5km/s 的速度破裂。对于每个震级, 分别将破裂面积设置在较浅的延伸到距离海底 2km 处和较深的 60km 处, 使得每个地震能在海底引起不同大小和分布范围的永久形变区域。
位于不同震源深度的 Mw 7.4 和 Mw 8.6 有限断层地震引起的海面初始波高分布如图 9 所示。可以看到, 对于 Mw 7.4 地震, 当断层面破裂区域中心深度为 7.5km 时, 海底出现低频压强不满足 p=rha 的区域, 当断层面破裂区域中心深度为 60km 时, 海底没有出现不满足 p=rha 的台站; 对于 Mw 8.6 地震, 当断层面破裂区域中心深度为 20km 时, 海底出现较大范围的低频压强不满足 p=rha 的区域, 当断层面破裂区域中心为 60km 时, 海底没有出现不满足 p=rha 的台站。由此证明了前面的结论, 即低频海底压强不满足 p=rha 的区域既可以反映震源区在海底的位置, 也可以指示震源深度。
低频海底压强不满足 p=rha 的区域范围与海底(海面)永久形变的范围存在一定的对应关系。然而, 对比图 9(c)和(d)两个 Mw 8.6 地震可以看到, 当断层面破裂区域中心深度为 20km 时, 海底永久形变较大, 海底出现低频压强不满足 p=rha 的区域, 当断层面破裂区域中心深度为 60km 时, 海底永久形变依然很大, 海底却没有出现不满足 p=rha 的台站。再对比图 9(b)和(d)破裂中心深度为 60km 的 Mw 7.4和 Mw 8.6 地震, 可以看到, 无论海底永久形变是否较大, 海底都不存在低频压强不满足 p=rha 的台站。针对点源震源的研究结果(图 5)也表明, 大震级引发大的海底永久形变, 但低频压强不满足 p= rha 的区域与小震级相同。因此我们推测, 海底永久形变并不是造成海底出现低频压强不满足 p=rha区域的原因。一种可能原因是, 推导海底低频压强满足 p=rha 的关系式时, 做了地震波以平面波入射海底的假设[13], 而实际震源辐射出的地震波在近场区域是球面波, 随着震源深度的增加, 才逐渐趋近于以平面波入射海底。如果近场球面波为主因, 又很难解释走滑地震中在远场也出现不满足理论公式的现象(图 5(g)、(h)和(i))。对于走滑地震, 沿着经过震中且垂直于走向的水平线、存在较为明显的位移梯度, 并且显著的位移梯度可以延伸到较远的距离, 这也可能是造成不满足理论公式的原因。同时在近场区, 位移极大的地方梯度很小, 也不满足理论公式, 其内在机理需进一步研究。
图8 Mw 7.4 和 Mw 8.6 地震对应的破裂区域位于不同深度时的断层面位错分布
Fig. 8 Slip distribution for earthquakes of Mw 7.4 and Mw 8.6 with rupture areas located at different depths
黑色点代表海底低频压强不满足p=rha的台站位置
图9 Mw 7.4 和 Mw 8.6 地震海底低频压强不满足p=rha的台站分布
Fig. 9 Station locations where the seafloor low-frequency pressure does not satisfy p=rha for earthquakes of Mw 7.4 and Mw 8.6
在远离震源区时, 压强信号中出现瑞利波被放大的现象。也就是说, 在加速度信号中, 直达体波(S 波)的振幅高于瑞利波, 但在压强信号中, 瑞利波的振幅则与 S 波相当。图 6(c)和(f)显示, 在台站S101 处, 加速度中 S 波的最高振幅约为 0.3m/s2, 瑞利波的振幅约为 0.03m/s2, 体波振幅为瑞利波振幅的 10 倍; 在压强信号中的振幅分别为 3×105和 2× 105Pa, 两者相当。为解释这种现象, 我们对加速度信号进行短时傅立叶变换, 分析其主导频率, 结果如图 10 所示。图 10(b)显示, 加速度信号中 S 波的主导频率大约为 0.25Hz, 瑞利波的主导频率大约为0.12Hz。根据 Deng 等[13]推导的海底压强 p 与垂直加速度 a 的关系, 海底压强与 rha 之间存在依赖于频率 f的放大系数 F:
其中, 
为入射波频率, c0 为水内声波波速, q 为地震波的射线参数。针对该台站, S 波大约为水平入射海底, 可以假设其射线参数 q为 0, 估算其放大系数为 0.15, 即 S 波引起的海底压强变化振幅约为 0.15×rha=1.35×105Pa; 当入射波为瑞利波时, q 为其相速度的倒数, 可以得到瑞利波的放大系数约为 2.2, 即瑞利波引起的海底压强变化振幅约为 2.2×rha=1.98×105Pa, 与数值模拟结果一致。因此, 在压强信号中观察到的瑞利波放大现象是因为瑞利波引起压强变化的放大系数远大于体波引起压强变化的放大系数。
当震源深度较大时, 在震源区地震波体波(P 波和 S 波)成为主要的地震波信号, 瑞利波在距离震源较远处才会逐渐形成, 使得距离震源较近处海底压强信号上的瑞利波不会像浅震源地震时较早地出现放大现象。可见, 相比垂直加速度信号被放大的瑞利波, 压强信号也可以用于指示震源深度。当我们在震源区附近观测到这一现象时, 说明震源深度较小, 海底可能存在较大的永久形变, 地震有可能引发海啸。
本研究采用 SPECFEM3D 程序构建三维地震海啸生成模型, 模拟不同震源工况下地震海啸的生成过程及其所引起的海底压强变化。通过统计整个海底不同区域的压强变化规律, 得到如下结论。
1)地震海啸生成过程中, 海底存在低频压强不满足 p=rha 关系式的区域, 可能与震源辐射的球面波相关。这部分区域的位置可以反映震源区在海底的位置, 区域的大小与震源深度密切相关。
2)地震海啸生成过程中, 海底压强信号上存在相比于垂直加速度信号被放大的瑞利波, 其原因是瑞利波引起压强变化的放大系数远大于体波引起压强变化的放大系数。这部分区域的范围与震源深度相关。
图10 图6 中 S101 台站的加速度信号及其短时傅里叶变换结果
Fig. 10 Acceleration signal of station S101 from Fig. 6 and its short-time Fourier transform results
3)海底存在低频海底压强不满足 p=rha 的区域, 与地震引起的海底永久形变并无直接关联。但是, 在其他震源参数确定的情况下, 这部分区域的范围可以大致反映海底永久形变的分布范围, 对于准确的海啸预警具有参考价值。
参考文献
[1] 王培涛, 于福江, 赵联大, 等. 2011 年 3 月 11 日日本地震海啸越洋传播及对中国影响的数值分析. 地球物理学报, 2012, 55(9): 3088–3096
[2] Nakano M, Nakamura T, Kamiya S I, et al. Seismic activity beneath the Nankai trough revealed by DONET ocean-bottom observations. Marine Geophysical Re-search, 2014, 35: 271–284
[3] Matsumoto H, Araki E. Drift characteristics of DONET pressure sensors determined from in-situ and expe-rimental measurements. Frontiers in Earth Science, 2021, 8: 600966
[4] Mochizuki M, Uehira K, Kanazawa T, et al. S-net pro-ject: performance of a large-scale seafloor observa- tion network for preventing and reducing seismic and tsunami disasters // 2018 OCEANS-MTS/IEEE Kobe Techno-Oceans (OTO): Additional Information. Kobe, 2018: 1–4
[5] Saito T. Dynamic tsunami generation due to sea-bottom deformation: analytical representation based on linear potential theory. Earth, Planets and Space, 2013, 65: 1411–1423
[6] 温瑞智, 周正华, 谢礼立. 基于强震台网的我国沿海海啸走时预警. 地震工程与工程振动, 2006, 26 (2): 20–24
[7] 温瑞智, 公茂盛, 谢礼立. 海啸预警系统及我国海啸减灾任务. 自然灾害学报, 2006, 15(3): 1–7
[8] 岳汉, 张勇, 盖增喜, 等. 大地震震源破裂模型: 从快速响应到联合反演的技术进展及展望. 中国科学: 地球科学, 2020, 50(4): 515–537
[9] An C, Sepúlveda I, Liu P L F. Tsunami source and its validation of the 2014 Iquique, Chile, earthquake. Geo-physical Research Letters, 2014, 41(11): 3988–3994
[10] Filloux J. Tsunami recorded on the open ocean floor. Geophysical Research Letters, 1982, 9(1): 25–28
[11] An C, Cai C, Zheng Y, et al. Theoretical solution and applications of ocean bottom pressure induced by seis-mic seafloor motion. Geophysical Research Letters, 2017, 44(20): 10272–10281
[12] 胡晨彤. 考虑海水可压性的海底震动与海水运动关系及其在海啸生成中的应用[D]. 上海: 上海交通大学, 2020
[13] Deng H, An C, Cai C, et al. Theoretical solution and applications of ocean bottom pressure induced by seis-mic waves at high frequencies. Geophysical Research Letters, 2022, 49(9): e2021GL096952
[14] Saito T. Tsunami generation: validity and limitations of conventional theories. Geophysical Journal Interna-tional, 2017, 210(3): 1888–1900
[15] Bolshakova A, Inoue S, Kolesov S, et al. Hydroacou-stic effects in the 2003 Tokachi-oki tsunami source. Russian Journal of Earth Sciences, 2011, 12(2): 1–14
[16] Matsumoto H, Inoue S, Ohmachi T. Dynamic response of bottom water pressure due to the 2011 Tohoku earth-quake. Journal of Disaster Research, 2012, 7(suppl 1): 468–475
[17] Saito T, Tsushima H. Synthesizing ocean bottom pres-sure records including seismic wave and tsunami con-tributions: toward realistic tests of monitoring systems. Journal of Geophysical Research: Solid Earth, 2016, 121(11): 8175–8195
[18] Zha Y, Webb S C. Crustal shear velocity structure in the Southern Lau Basin constrained by seafloor com-pliance. Journal of Geophysical Research: Solid Earth, 2016, 121(5): 3220–3237
[19] Zhou Y, Ni S, Chu R, et al. Accuracy of the water column approximation in numerically simulating pro-pagation of teleseismic PP waves and Rayleigh wa- ves. Geophysical Journal International, 2016, 206(2): 1315–1326
[20] Sun T, Wang K, Fujiwara T, et al. Large fault slip peaking at trench in the 2011 Tohoku-oki earthquake. Nature Communications, 2017, 8(1): 14044
[21] Lay T, Yue H, Brodsky E E, et al. The 1 April 2014 Iquique, Chile, Mw 8.1 earthquake rupture sequence. Geophysical Research Letters, 2014, 41(11): 3818–3825
[22] Nosov M, Kolesov S. Elastic oscillations of water column in the 2003 Tokachi-oki tsunami source: in-situ measurements and 3-D numerical modelling. Natural Hazards and Earth System Sciences, 2007, 7(2): 243–249
[23] Song L, An C. Extraction of Tsunami Signals from Coupled Seismic and Tsunami Waves. Journal of Ma-rine Science and Engineering, 2025, 13(3): no. 419
[24] 张风雪, 李昱, 陈泆平. 浅析震源位置准确度及其影响因素. 地球与行星物理论评, 2025, 56(2): 182–192
[25] Chen W P, Nábělek J L, Fitch T J, et al. An intermediate depth earthquake beneath Tibet: source characteristics of the event of September 14, 1976. Journal of Geophy-sical Research: Solid Earth, 1981, 86(B4): 2863–2876
[26] 高原, 周蕙兰, 郑斯华, 等. 测定震源深度的意义的初步讨论. 中国地震, 1997, 13(4): 321–329
[27] 张国民, 李丽, 马宏生, 等. 中国大陆地震震源深度及其构造含义. 科学通报, 2002, 47(9): 663–668
[28] Spence W, Sipkin S, Choy G. Determining the depth of an earthquake. Earthquakes and Volcanoes, 1989, 21 (1): 58–63
[29] 项月文, 陈浩, 肖孟仁. sPL 震相在江西地区中小地震震源深度测定中的应用. 地震科学进展, 2019, (11): 14–19
[30] 张明辉, 徐涛, 田小波, 等. 虚拟震源地震探测方法及其应用. 地球与行星物理论评, 2025, 56(2): 215–224
[31] Templeton D, Rodgers A, Helmberger D, et al. Com-parison of the cut-and-paste and full moment tensor methods for estimating earthquake source parameters // AGU Fall Meeting Abstracts: Additional Informa-tion. San Francisco, 2008: S41C–1864
[32] Neto G S L, Julià J. Determination of intraplate focal mechanisms with the Brazilian seismic network: a simplified cut-and-paste approach. Journal of South American Earth Sciences, 2023, 121: 104149
[33] Ching J Y, Glaser S. A time domain moment tensor inversion technique and its verification. Trends in Rock Mechanics, 2000: 140–151
[34] Dreger D S, Helmberger D V. Determination of source parameters at regional distances with three‐component sparse network data. Journal of Geophysical Research: Solid Earth, 1993, 98(B5): 8107–8125
[35] Dettmer J, Hawkins R, Cummins P R, et al. Tsunami source uncertainty estimation: the 2011 Japan tsunami. Journal of Geophysical Research: Solid Earth, 2016, 121(6): 4483–4505
[36] Wei Y, Newman A V, Hayes G P, et al. Tsunami forecast by joint inversion of real-time tsunami waveforms and seismic or GPS data: application to the Tohoku 2011 tsunami. Pure and Applied Geophysics, 2014, 171: 3281–3305
[37] Kanamori H. Mechanism of tsunami earthquakes. Phy-sics of the earth and planetary interiors, 1972, 6(5): 346–359
[38] Saito T. Tsunami generation and propagation. Tokyo: Springer, 2019
[39] Lotto G C, Nava G, Dunham E M. Should tsunami simulations include a nonzero initial horizontal velo-city?. Earth, Planets and Space, 2017, 69: 1–14
[40] Hu C, Wu Y, An C, et al. A numerical study of tsunami generation by horizontal displacement of sloping sea-floor. Journal of Earthquake and Tsunami, 2020, 14(4): 2050018
[41] Komatitsch D, Tromp J. Spectral-element simulations of global seismic wave propagation — I. Validation. Geophysical Journal International, 2002, 149(2): 390–412
[42] Komatitsch D, Tromp J. Spectral-Element Simulations of Global Seismic Wave Propagation: 3-D Models, Oceans, Rotation, and Self-Gravitation // AGU Fall Meeting Abstracts: Additional Information. San Fran-cisco, 2001: U42B-09
[43] Okada Y. Surface deformation due to shear and tensile faults in a half-space. Bulletin of the Seismological Society of America, 1985, 75(4): 1135–1154
[44] Nosov M. Tsunami generation in compressible ocean. Physics and Chemistry of the Earth, Part B: Hydrology, Oceans and Atmosphere, 1999, 24(5): 437–441
[45] An C, Liu H, Ren Z, et al. Prediction of tsunami waves by uniform slip models. Journal of Geophysical Re-search: Oceans, 2018, 123(11): 8366–8382
[46] Blaser L, Krüger F, Ohrnberger M, et al. Scaling rela-tions of earthquake source parameter estimates with special focus on subduction environment. Bulletin of the Seismological Society of America, 2010, 100(6): 2914–2926
[47] Murotani S, Miyake H, Koketsu K. Scaling of charac-terized slip models for plate-boundary earthquakes. Earth, Planets and Space, 2008, 60: 987–991
[48] Wells D L, Coppersmith K J. New empirical relation-ships among magnitude, rupture length, rupture width, rupture area, and surface displacement. Bulletin of the Seismological Society of America, 1994, 84(4): 974–1002
Characteristic Analysis of Seafloor Pressure Variations Generated by Earthquakes Based on Numerical Simulation
Abstract When a submarine earthquake occurs, pressure sensors deployed on the seafloor can record the pressure variations induced by the earthquake. Currently there is no quantitative application for utilizing earthquake-induced pressure variations. In this study, SPECFEM3D program was employed to construct a numerical model of a submarine earthquake, and simulate the rupture process and the associated seafloor pressure variations. The spatial distribution characteristics of the seafloor pressure signals were analyzed. The results indicate that at stations far from the source region, the low-frequency pressure p induced by the earthquake and vertical acceleration a of the seafloor satisfies the theoretical relationship p=rha. In areas near the source, this theoretical relationship no longer holds, and the transition boundary between these two regimes is related to the earthquake source depth. If other source parameters are the same, such as magnitude and rake, the source depth determines the permanent seafloor deformation, which provides valuable insights for assessing tsunami threat potential. Additionally, it is found that Rayleigh waves can be recorded outside a certain range of the source region, and the pressure variations caused by Rayleigh waves exhibit an amplification effect compared with seismic body waves. The locations of stations where Rayleigh wave amplification occurs are also related to earthquake source depth, providing useful insights for tsunami early warning.
Key words earthquake-induced tsunami generation process; permanent seafloor deformation; seafloor pressure; Rayleigh waves