Spectral-element simulations of long-term fault slip: Effect of low-rigidity layers on earthquake-cycle dynamics

We develop a spectral element method for the simulation of long-term histories of spontaneous seismic and aseismic slip on faults subjected to tectonic loading. Our approach reproduces all stages of earthquake cycles: nucleation and propagation of earthquake rupture, postseismic slip and interseismic creep. We apply the developed methodology to study the effects of low-rigidity layers on the dynamics of the earthquake cycle in 2-D. We consider two cases: small (M ~ 1) earthquakes on a fault surrounded by a damaged fault zone and large (M ~ 7) earthquakes on a vertical strike-slip fault that cuts through shallow low-rigidity layers. Our results indicate how the source properties of repeating earthquakes are affected by the presence of a damaged fault zone with low rigidity. Compared to faults in homogeneous media, we find (1) reduction in the earthquake nucleation size, (2) amplification of slip rates during dynamic rupture propagation, (3) larger recurrence interval, and (4) smaller amount of aseismic slip. Based on linear stability analysis, we derive a theoretical estimate of the nucleation size as a function of the width and rigidity reduction of the fault zone layer, which is in good agreement with simulated nucleation sizes. We further examine the effects of vertically-stratified layers (e.g., sedimentary basins) on the nature of shallow coseismic slip deficit. Our results suggest that low-rigidity shallow layers alone do not lead to coseismic slip deficit. While the low-rigidity layers result in lower interseismic stress accumulation, they also cause dynamic amplification of slip rates, with the net effect on slip being nearly zero

Three-dimensional dynamic rupture simulation with a high-order discontinuous Galerkin method on unstructured tetrahedral meshes

Accurate and efficient numerical methods to simulate dynamic earthquake rupture and wave propagation in complex media and complex fault geometries are needed to address fundamental questions in earthquake dynamics, to integrate seismic and geodetic data into emerging approaches for dynamic source inversion, and to generate realistic physics-based earthquake scenarios for hazard assessment. Modeling of spontaneous earthquake rupture and seismic wave propagation by a high-order discontinuous Galerkin (DG) method combined with an arbitrarily high-order derivatives (ADER) time integration method was introduced in two dimensions by de la Puente et al. (2009). The ADER-DG method enables high accuracy in space and time and discretization by unstructured meshes. Here we extend this method to three-dimensional dynamic rupture problems. The high geometrical flexibility provided by the usage of tetrahedral elements and the lack of spurious mesh reflections in the ADER-DG method allows the refinement of the mesh close to the fault to model the rupture dynamics adequately while concentrating computational resources only where needed. Moreover, ADER-DG does not generate spurious high-frequency perturbations on the fault and hence does not require artificial Kelvin-Voigt damping. We verify our three-dimensional implementation by comparing results of the SCEC TPV3 test problem with two well-established numerical methods, finite differences, and spectral boundary integral. Furthermore, a convergence study is presented to demonstrate the systematic consistency of the method. To illustrate the capabilities of the high-order accurate ADER-DG scheme on unstructured meshes, we simulate an earthquake scenario, inspired by the 1992 Landers earthquake, that includes curved faults, fault branches, and surface topography

Spectral-element modeling of spontaneous earthquake rupture on rate and state faults: Effect of velocity-strengthening friction at shallow depths

We develop a spectral-element methodology (SEM) for simulating dynamic rupture on rate and state faults and use it to study how the rupture is affected by a shallow fault region of steady-state velocity-strengthening friction. Our comparison of the developed SEM and a spectral boundary-integral method (BIM) for an anti-plane (two-dimensional) test problem shows that the two methods produce virtually identical solutions for the finest resolution we use and that the convergence with grid reduction of the developed SEM methodology is comparable to that of BIM. We also use the test problem to compare numerical resolution required for different state evolution laws and for linear slip-weakening friction. Using our three-dimensional implementation of the methodology, we find that a shallow velocity-strengthening fault region can significantly alter dynamic rupture and ground motion. The velocity-strengthening region suppresses supershear propagation at the free surface occurring in the absence of such region, which could explain the lack of universally observed supershear rupture near the free surface. In addition, the velocity-strengthening region promotes faster fall-off of slip velocity behind the rupture front and decreases final slip throughout the entire fault, causing a smaller average stress drop. The slip decrease is largest in the shallow parts of the fault, resulting in the depth profile of slip qualitatively consistent with observations of shallow co-seismic slip deficit. The shallow velocity-strengthening region also reduces the amplification of strong ground motion due to a low-velocity bulk structure

Inelastic material response in multi-physics earthquake rupture simulations

Dynamic rupture models are able to shed light on earthquake source dynamics where direct observations are rare or non-existent. These multi-physics simulations incorporate earthquake rupture along a fault governed by frictional constitutive laws, which is coupled to seismic wave propagation described by the linear elastic wave equation. To accurately model the earthquake process, numerical models need to include realistic material properties such as the ability of rocks to deform plastically. This dissertation extends the Arbitrary High Order Derivative Discontinuous Galerkin (ADER-DG) framework of the dynamic rupture software SeisSol to account for non-linear off-fault plasticity. The impact of plasticity on rupture dynamics and the emitted seismic wave field is investigated in realistic simulations motivated by past earthquakes on geometrically complex faults. We first present the implementation of off-fault plasticity, which is verified in community benchmark problems and by three-dimensional numerical refinement studies. Motivated by the high efficiency of the implementation, we present a large-scale simulation of earthquake rupture along the segmented fault system of the 1992 Landers earthquake including plasticity. The results indicate that spatio-temporal rupture transfers are altered by plastic energy absorption, correlating with locations of geometrical fault complexity. In a next step, the model of the 1992 Landers earthquake is further extended to account for a new degree of realism among dynamic rupture models by incorporating high-resolution topography, 3D velocity structure, and viscoelastic attenuation in addition to off-fault plasticity. The simulation reproduces a broad range of observations including moment release rate, seismic waveform characteristics, mapped off-fault deformation patterns, and peak ground motions. We find that plasticity reduces the directivity effect and the spatial variability of peak ground velocities in comparison to the purely elastic simulation. In addition to this continental strike-slip earthquake, we investigate the effect of off-fault plasticity on source dynamics and seafloor deformation in a 3D subduction zone model of the 2004 Sumatra-Andaman earthquake. Simulated seafloor displacements are drastically altered by inelastic processes within the entire accretionary wedge, depending on fault- strike and the applied regional stress field, which potentially affects the tsunamigenesis. Finally, since these application scenarios show that rupture dynamics and the occurrence of off-fault plasticity are highly influenced by the assumed initial stresses and fault geometry, we propose a workflow to constrain dynamic rupture initial conditions with plasticity by long-term seismic cycling modelling. The exploited seismo-thermo-mechanical model provides a self-consistent slab geometry as well as initial stress and strength conditions that evolve according to the tectonic stress build-up and the temperature-dependent strength of the rocks. The geomechanically constrained subduction zone model suggests that the accretionary wedge is very close to plastic failure such that the occurrence of plastic strain hampers rupture to the trench, but locally increases the vertical seafloor uplift

Inelastic material response in multi-physics earthquake rupture simulations

Dynamic rupture models are able to shed light on earthquake source dynamics where direct observations are rare or non-existent. These multi-physics simulations incorporate earthquake rupture along a fault governed by frictional constitutive laws, which is coupled to seismic wave propagation described by the linear elastic wave equation. To accurately model the earthquake process, numerical models need to include realistic material properties such as the ability of rocks to deform plastically. This dissertation extends the Arbitrary High Order Derivative Discontinuous Galerkin (ADER-DG) framework of the dynamic rupture software SeisSol to account for non-linear off-fault plasticity. The impact of plasticity on rupture dynamics and the emitted seismic wave field is investigated in realistic simulations motivated by past earthquakes on geometrically complex faults. We first present the implementation of off-fault plasticity, which is verified in community benchmark problems and by three-dimensional numerical refinement studies. Motivated by the high efficiency of the implementation, we present a large-scale simulation of earthquake rupture along the segmented fault system of the 1992 Landers earthquake including plasticity. The results indicate that spatio-temporal rupture transfers are altered by plastic energy absorption, correlating with locations of geometrical fault complexity. In a next step, the model of the 1992 Landers earthquake is further extended to account for a new degree of realism among dynamic rupture models by incorporating high-resolution topography, 3D velocity structure, and viscoelastic attenuation in addition to off-fault plasticity. The simulation reproduces a broad range of observations including moment release rate, seismic waveform characteristics, mapped off-fault deformation patterns, and peak ground motions. We find that plasticity reduces the directivity effect and the spatial variability of peak ground velocities in comparison to the purely elastic simulation. In addition to this continental strike-slip earthquake, we investigate the effect of off-fault plasticity on source dynamics and seafloor deformation in a 3D subduction zone model of the 2004 Sumatra-Andaman earthquake. Simulated seafloor displacements are drastically altered by inelastic processes within the entire accretionary wedge, depending on fault- strike and the applied regional stress field, which potentially affects the tsunamigenesis. Finally, since these application scenarios show that rupture dynamics and the occurrence of off-fault plasticity are highly influenced by the assumed initial stresses and fault geometry, we propose a workflow to constrain dynamic rupture initial conditions with plasticity by long-term seismic cycling modelling. The exploited seismo-thermo-mechanical model provides a self-consistent slab geometry as well as initial stress and strength conditions that evolve according to the tectonic stress build-up and the temperature-dependent strength of the rocks. The geomechanically constrained subduction zone model suggests that the accretionary wedge is very close to plastic failure such that the occurrence of plastic strain hampers rupture to the trench, but locally increases the vertical seafloor uplift

Data-driven optimization of seismicity models using diverse data sets: generation, evaluation and ranking using inlabru

Recent developments in earthquake forecasting models have demonstrated the need for a robust method for identifying which model components are most beneficial to understanding spatial patterns of seismicity. Borrowing from ecology, we use Log‐Gaussian Cox process models to describe the spatially varying intensity of earthquake locations. These models are constructed using elements which may influence earthquake locations, including the underlying fault map and past seismicity models, and a random field to account for any excess spatial variation that cannot be explained by deterministic model components. Comparing the alternative models allows the assessment of the performance of models of varying complexity composed of different components and therefore identifies which elements are most useful for describing the distribution of earthquake locations. We demonstrate the effectiveness of this approach using synthetic data and by making use of the earthquake and fault information available for California, including an application to the 2019 Ridgecrest sequence. We show the flexibility of this modeling approach and how it might be applied in areas where we do not have the same abundance of detailed information. We find results consistent with existing literature on the performance of past seismicity models that slip rates are beneficial for describing the spatial locations of larger magnitude events and that strain rate maps can constrain the spatial limits of seismicity in California. We also demonstrate that maps of distance to the nearest fault can benefit spatial models of seismicity, even those that also include the primary fault geometry used to construct them.K. B. was funded during this work by an EPSRC PhD studentship (Grant 1519006) and during the writing of this paper by NERC‐NSF grant NE/R000794/1 and by the Real‐time Earthquake Risk Reduction for a Resilient Europe "RISE" project, which has received funding from the European Union's Horizon 2020 research and innovation program under grant Agreement 821115

Importance of earthquake rupture geometry on tsunami modelling: the Calabrian Arc subduction interface (Italy) case study

SUMMARY The behaviour of tsunami waves at any location depends on the local morphology of the coasts, the encountered bathymetric features, and the characteristics of the source. However, the importance of accurately modelling the geometric properties of the causative fault for simulations of seismically induced tsunamis is rarely addressed. In this work, we analyse the effects of using two different geometric models of the subduction interface of the Calabrian Arc (southern Italy, Ionian Sea) onto the simulated tsunamis: a detailed 3-D subduction interface obtained from the interpretation of a dense network of seismic reflection profiles, and a planar interface that roughly approximates the 3-D one. These models can be thought of as representing two end-members of the level of knowledge of fault geometry. We define three hypothetical earthquake ruptures of different magnitudes (Mw 7.5, 8.0, 8.5) on each geometry. The resulting tsunami impact is evaluated at the 50-m isobath in front of coastlines of the central and eastern Mediterranean. Our results show that the source geometry imprint is evident on the tsunami waveforms, as recorded at various distances and positions relative to the source. The absolute differences in maximum and minimum wave amplitudes locally exceed one metre, and the relative differences remain systematically above 20 per cent with peaks over 40 per cent. We also observe that tsunami energy directivity and focusing due to bathymetric waveguides take different paths depending on which fault is used. Although the differences increase with increasing earthquake magnitude, there is no simple rule to anticipate the different effects produced by these end-member models of the earthquake source. Our findings suggest that oversimplified source models may hinder our fundamental understanding of the tsunami impact and great care should be adopted when making simplistic assumptions regarding the appropriateness of the planar fault approximation in tsunami studies. We also remark that the geological and geophysical 3-D fault characterization remains a crucial and unavoidable step in tsunami hazard analyses

The role of fluid pressure in induced vs. triggered seismicity. Insights from rock deformation experiments on carbonates

Fluid overpressure is one of the primary mechanisms for tectonic fault slip, because fluids lubricate the fault and fluid pressure reduces the effective normal stress that holds the fault in place. However, current models of earthquake nucleation, based on rate- and state- friction laws, imply that stable sliding is favoured by the increase of pore fluid pressure. Despite this controversy, currently, there are only a few studies on the role of fluid pressure under controlled, laboratory conditions. Here, we use laboratory experiments, to show that the rate- and state- friction parameters do change with increasing fluid pressure. We tested carbonate gouges from sub hydrostatic to near lithostatic fluid pressure conditions, and show that the friction rate parameter (a−b) evolves from velocity strengthening to velocity neutral behaviour. Furthermore, the critical slip distance, Dc, decreases from about 90 to 10μm. Our data suggest that fluid overpressure plays an important role in controlling the mode of fault slip. Since fault rheology and fault stability parameters change with fluid pressure, we suggest that a comprehensive characterization of these parameters is fundamental for better assessing the role of fluid pressure in natural and human induced earthquakes

Linked 3-D modelling of megathrust earthquake-tsunami events: from subduction to tsunami run up

How does megathrust earthquake rupture govern tsunami behaviour? Recent modelling advances permit evaluation of the influence of 3-D earthquake dynamics on tsunami genesis, propagation, and coastal inundation. Here, we present and explore a virtual laboratory in which the tsunami source arises from 3-D coseismic seafloor displacements generated by a dynamic earthquake rupture model. This is achieved by linking open-source earthquake and tsunami computational models that follow discontinuous Galerkin schemes and are facilitated by highly optimized parallel algorithms and software. We present three scenarios demonstrating the flexibility and capabilities of linked modelling. In the first two scenarios, we use a dynamic earthquake source including time-dependent spontaneous failure along a 3-D planar fault surrounded by homogeneous rock and depth-dependent, near-lithostatic stresses. We investigate how slip to the trench influences tsunami behaviour by simulating one blind and one surface-breaching rupture. The blind rupture scenario exhibits distinct earthquake characteristics (lower slip, shorter rupture duration, lower stress drop, lower rupture speed), but the tsunami is similar to that from the surface-breaching rupture in run-up and length of impacted coastline. The higher tsunami-generating efficiency of the blind rupture may explain how there are differences in earthquake characteristics between the scenarios, but similarities in tsunami inundation patterns. However, the lower seafloor displacements in the blind rupture result in a smaller displaced volume of water leading to a narrower inundation corridor inland from the coast and a 15 per cent smaller inundation area overall. In the third scenario, the 3-D earthquake model is initialized using a seismo-thermo-mechanical geodynamic model simulating both subduction dynamics and seismic cycles. This ensures that the curved fault geometry, heterogeneous stresses and strength and material structure are consistent with each other and with millions of years of modelled deformation in the subduction channel. These conditions lead to a realistic rupture in terms of velocity and stress drop that is blind, but efficiently generates a tsunami. In all scenarios, comparison with the tsunamis sourced by the time-dependent seafloor displacements, using only the time-independent displacements alters tsunami temporal behaviour, resulting in later tsunami arrival at the coast, but faster coastal inundation. In the scenarios with the surface-breaching and subduction-initialized earthquakes, using the time-independent displacements also overpredicts run-up. In the future, the here presented scenarios may be useful for comparison of alternative dynamic earthquake-tsunami modelling approaches or linking choices, and can be readily developed into more complex applications to study how earthquake source dynamics influence tsunami genesis, propagation and inundation