Dynamic Slip Rupture at a Layer-Half-space Interface: A Hybrid Finite Difference and Spectral Boundary Integral Equation Method for Numerical Simulation

The propagation of rupture at the interface between a layer and a half-space, arising from their relative slipping, is commonly observed in geophysics and in technological materials. Fault ruptures develop as a consequence of ongoing aseismic interfacial movement caused by geological forces; and accumulation of stress over time ultimately leads to a sudden stress release resulting in seismic phenomena. In mathematical terms, it is a highly nonlinear and multiscale phenomenon that necessitates the accurate solution of elastodynamic equations and interfacial fault friction across an extensive domain over an extended time period. The most numerically efficient algorithm for simulating this process is the Boundary Integral Equation Method (BIEM). It computes field quantities at the fracture plane and reduces domain dimensionality of the problem by one, but its current applicability is mostly limited to planar interfaces and to unbounded geometries. BIEM computes elastodynamic fields at the fracture plane from a space-time traction history by utilizing suitable convolution kernels (dependent on the geometry and the material/bi-material characteristics across the interface). Applying BIEM to a geometry involving a layer over a half space presents significant challenges, particularly due to the difficulty in analytically deriving convolution kernels for finite geometry and in-plane deformations. Here, we develop an integrated approach, combining the traction formulation of BIEM (Ranjith, 2015) with the Finite Difference Method (FDM), to overcome these challenges. A hybrid FDM-BIEM has previously been used by Hajarolasvadi and Elbanna (2017) using the velocity formulation of BIEM (Geubelle and Rice, 1995) for unbounded geometries. In the present work, a fourth-order staggered-grid FDM is employed to model dynamic rupturing in a layer. The layer is adjoined by an elastic half-space. The elastodynamic equations in the half-space are handled using a BIEM. The FDM region incorporates a slip-weakening friction law at the fault interface using a thick fault zone model consisting of two grid rows (Madariaga, Olsen and Archuleta, 1998). By applying traction boundary conditions, it computes updated velocities in the domain, and passes them to BIEM. BIEM, in turn, responds with updated tractions at the boundary, utilizing convolution kernels derived by Ranjith (2015). We utilized the hybrid FDM-BIEM numerical algorithm to study the effect of finite layer thickness on rupture propagation, considering both homogeneous and bi-material fault interfaces. The results suggest that the present method has the capacity to effectively deal with a wide array of problems in finite domains often countered in both geophysics and technological domains.