第三类模式反平面地震动力学模拟系统

资 源 简 介

该项目包含一系列MATLAB脚本，专门用于求解嵌入在二维均匀弹性无限介质中的一维断层上的反平面（第三类模式/Mode III）地震动力学问题。系统采用了滑移弱化摩擦定律（slip-weakening friction）来描述断层面的摩擦行为。在核心算法上，该脚本将物理问题公式化为边界积分方程（Boundary Integral Equation），其弹性动力学核函数是在谱域（空间波数域）中进行解析推导的，从而实现了高精度的数值求解。该仿真方法基于Morrissey和Geubelle（1997）的理论基础，并结合了Nadia Lapusta、Alain Cochard等学者后续改进的谱边界积分法（SBI）技术。项目主要功能包括模拟地震破裂的成核、动态传播及停止过程，精确计算断层上的应力演化、滑移速率以及弹性波在无限介质中的辐射效应，适用于地震震源物理、断层力学以及摩擦滑动稳定性等领域的科学研究与教学演示。

详 情 说 明

第三类模式(Mode III)反平面地震动力学模拟系统

项目简介

该项目是一个基于MATLAB开发的高精度地震动力学仿真系统，专注于解决二维无限弹性介质中一维断层上的反平面（Mode III）剪切破裂问题。该模拟器结合了谱边界积分法（Spectral Boundary Integral Method, SBI）的空间高精度特性与显式时间积分方案，能够模拟地震破裂从成核、动态传播到最终停止的全过程。

核心物理机制采用了经典的线性滑移弱化摩擦定律（Linear Slip-Weakening Friction Law），通过在频域（谱域）计算弹性应力传递，并结合辐射阻尼项（Radiation Damping）来近似处理动力学效应，实现了高效且物理意义明确的断层力学行为数值模拟。

功能特性

• 反平面（Mode III）地震模拟：模拟剪切波（SH波）主导的断层破裂过程。
• 谱域高精度算法：利用快速傅里叶变换（FFT）和逆变换（IFFT）在波数域进行空间导数和应力卷积计算，相比传统有限差分或有限元方法具有更高的空间精度和计算效率。
• 滑移弱化摩擦模型：实现了依赖于滑移量的本构关系，精确描述从静态摩擦系数过渡到动态摩擦系数的过程，捕捉断层的脆性破裂特征。
• 人工成核机制：通过叠加高斯型应力扰动，在设定的背景应力场中人为诱发地震成核。
• 吸收边界条件：实现了“海绵层”（Sponge Layer）吸收边界，有效抑制了由于FFT周期性边界条件导致的数值边界反射波，模拟无限介质环境。
• 全过程物理量记录：记录并可视化断层滑移量（Slip）、滑移速率（Slip Rate）以及剪切应力（Shear Stress）的时空演化。

算法原理与实现细节

该系统完全基于MATLAB脚本实现，其核心算法流程如下：

1. 物理模型与参数化

系统定义了弹性介质的基本属性（剪切模量、密度、波速）以及断层的几何尺度。空间上将断层离散化为均匀网格，利用FFT所需的点数（2的幂次）进行配置。时间积分遵循Courant-Friedrichs-Lewy (CFL) 稳定性条件，确保显式时间步进的收敛性。

2. 谱域应力计算 (Spectral Stress Transfer)

利用边界积分方程，断层上的弹性相互作用通过卷积形式表示。在实现上，这部分通过谱域乘法完成：

• 对当前的断层滑移量进行FFT变换。
• 将变换结果乘以静态弹性核函数 $K(k) = -0.5 mu |k|$。该核函数描述了Mode III裂纹在波数域的静态应力响应。
• 通过IFFT将结果变换回空间域，得到由断层非均匀滑移引起的弹性静力传递应力。

3. 摩擦求解器与辐射阻尼 (Friction Logic)

系统采用准动态（Quasi-dynamic）公式，将动力学效应显式包含在辐射阻尼项中：

• 总负载应力：由初始背景应力与弹性传递应力叠加而成。
• 强度计算：根据当前累积滑移量，线性插值计算当前的摩擦强度。当滑移量超过临界距离 $D_c$ 后，强度维持在残余摩擦水平。
• 速率求解：基于动力学平衡方程 $tau_{load} - eta V = tau_{strength}$ 求解滑移速率 $V$。其中 $eta = mu / (2c_s)$ 为辐射阻尼系数，代表弹性波向周围介质辐射能量产生的能量损失。若计算出的速率为负，则强制设为0（断层锁定状态）。

4. 时间积分与状态更新

采用显式时间积分策略（如欧拉法），根据当前计算出的滑移速率更新下一时刻的滑移量。

5. 边界处理 (Sponge Layer)

为了处理FFT隐含的周期性边界条件（即模拟域左端的波会从右端进入），在空间网格的两端设置了吸收层。通过构建一个空间掩膜（Mask），对边缘区域的滑移速率和增量滑移进行二次方衰减，从而模拟无限大介质，防止虚假反射波干扰破裂区。

系统要求

• MATLAB: 版本 R2018a 或更高（推荐），需包含基础工具箱。
• 内存: 取决于离散点数 $N$ 和记录的时间步数，通常 8GB RAM 足以应对 $N=512$ 的演示级模拟。

使用方法

1. 确保MATLAB环境已正确安装。
2. 打开包含脚本的文件夹。
3. 直接运行主程序脚本。
4. 程序将初始化物理场，并在计算过程中弹出图形窗口（视代码具体绘图实现而定），并在控制台输出模拟进度信息（如点数、时间步长、总步数）。
5. 模拟结束后，工作区将保留 Store_Slip, Store_Rate, Store_Stress 等矩阵，可用于后续的数据分析或动态视频制作。

关键参数说明

在代码中可以通过修改以下变量来调整模拟设置：

• N: 空间离散点数（建议保持为2的幂，如256, 512, 1024）。
• L: 模拟域的物理总长度。
• TotalTime: 模拟的总物理时间。
• sigma_n: 断层上的有效正应力，直接控制摩擦强度的大小。
• fs / fd: 静态与动态摩擦系数，决定了应力降的大小。
• Dc: 临界滑移距离，控制破裂尖端的内聚区（Cohesive Zone）尺寸。
• nucleation_width: 人工成核区域的宽度。