Search for: All records

Earthquakes are rupture-like processes that propagate along tectonic faults and cause seismic waves. The propagation speed and final area of the rupture, which determine an earthquake’s potential impact, are directly related to the nature and quantity of the energy dissipation involved in the rupture process. Here, we present the challenges associated with defining and measuring the energy dissipation in laboratory and natural earthquakes across many scales. We discuss the importance and implications of distinguishing between energy dissipation that occurs close to and far behind the rupture tip, and we identify open scientific questions related to a consistent modeling framework for earthquake physics that extends beyond classical Linear Elastic Fracture Mechanics.

Megathrust geometric properties exhibit some of the strongest correlations with maximum earthquake magnitude in global surveys of large subduction zone earthquakes, but the mechanisms through which fault geometry influences subduction earthquake cycle dynamics remain unresolved. Here, we develop 39 models of sequences of earthquakes and aseismic slip (SEAS) on variably‐dipping planar and variably‐curved nonplanar megathrusts using the volumetric, high‐order accurate codetandemto account for fault curvature. We vary the dip, downdip curvature and width of the seismogenic zone to examine how slab geometry mechanically influences megathrust seismic cycles, including the size, variability, and interevent timing of earthquakes. Dip and curvature control characteristic slip styles primarily through their influence on seismogenic zone width: wider seismogenic zones allow shallowly‐dipping megathrusts to host larger earthquakes than steeply‐dipping ones. Under elevated pore pressure and less strongly velocity‐weakening friction, all modeled fault geometries host uniform periodic ruptures. In contrast, shallowly‐dipping and sharply‐curved megathrusts host multi‐period supercycles of slow‐to‐fast, small‐to‐large slip events under higher effective stresses and more strongly velocity‐weakening friction. We discuss how subduction zones' maximum earthquake magnitudes may be primarily controlled by the dip and dimensions of the seismogenic zone, while second‐order effects from structurally‐derived mechanical heterogeneity modulate the recurrence frequency and timing of these events. Our results suggest that enhanced co‐ and interseismic strength and stress variability along the megathrust, such as induced near areas of high or heterogeneous fault curvature, limits how frequently large ruptures occur and may explain curved faults' tendency to host more frequent, smaller earthquakes than flat faults.

Structural inversion of rifted basins is generally associated with surface uplift and denudation of the sedimentary infill, reflecting the active contractional deformation in the crust. However, worldwide examples of inverted rifts show contrasting basin-scale subsidence and widespread sedimentation patterns during basin inversion. By conducting a series of three-dimensional coupled geodynamic and surface processes models, we investigated the dynamic controls on these subsidence anomalies during the successive stages of rifting and basin inversion, and we propose a new evolutionary model for this process. Our models show that the inherited thermo-rheological properties of the lithosphere influence the initial strain localization and subsequent migration of crustal deformation during inversion. The sense of the vertical movements (i.e., uplift or subsidence), however, is not directly linked to the underlying crustal stress patterns; rather, it reflects the balance among contraction-induced tectonic uplift, postrift thermal subsidence of the inherited lithosphere, and sediment redistribution. Based on the interplay among the competing differential vertical movements with different amplitudes and wavelengths, inversion of rifted basins may lead to the growth of intraplate orogens, or the contraction-driven localized uplift may be hindered by the thermal sag effects of the inherited shallow lithosphere-asthenosphere boundary, resulting in basin-scale subsidence. In such basins, dating the first erosional surfaces and other unconformities may not provide accurate timing for the onset of inversion.

Slow slip events (SSEs) have been observed in spatial and temporal proximity to megathrust earthquakes in various subduction zones, including the 2014M_w7.3 Guerrero, Mexico earthquake which was preceded by aM_w7.6 SSE. However, the underlying physics connecting SSEs to earthquakes remains elusive. Here, we link 3D slow‐slip cycle models with dynamic rupture simulations across the geometrically complex flat‐slab Cocos plate boundary. Our physics‐based models reproduce key regional geodetic and teleseismic fault slip observations on timescales from decades to seconds. We find that accelerating SSE fronts transiently increase shear stress at the down‐dip end of the seismogenic zone, modulated by the complex geometry beneath the Guerrero segment. The shear stresses cast by the migrating fronts of the 2014M_w7.6 SSE are significantly larger than those during the three previous episodic SSEs that occurred along the same portion of the megathrust. We show that the SSE transient stresses are large enough to nucleate earthquake dynamic rupture and affect rupture dynamics. However, additional frictional asperities in the seismogenic part of the megathrust are required to explain the observed complexities in the coseismic energy release and static surface displacements of the Guerrero earthquake. We conclude that it is crucial to jointly analyze the long‐ and short‐term interactions and complexities of SSEs and megathrust earthquakes across several (a)seismic cycles accounting for megathrust geometry. Our study has important implications for identifying earthquake precursors and understanding the link between transient and sudden megathrust faulting processes.

Two‐phase flow, a system where Stokes flow and Darcy flow are coupled, is of great importance in the Earth's interior, such as in subduction zones, mid‐ocean ridges, and hotspots. However, it remains challenging to solve the two‐phase equations accurately in the zero‐porosity limit, for example, when melt is fully frozen below solidus temperature. Here we propose a new three‐field formulation of the two‐phase system, with solid velocity (v_s), total pressure (P_t), and fluid pressure (P_f) as unknowns, and present a robust finite‐element implementation, which can be used to solve problems in which domains of both zero porosity and non‐zero porosity are present. The reformulated equations include regularization to avoid singularities and exactly recover to the standard single‐phase incompressible Stokes problem at zero porosity. We verify the correctness of our implementation using the method of manufactured solutions and analytic solutions and demonstrate that we can obtain the expected convergence rates in both space and time. Example experiments, such as self‐compaction, falling block, and mid‐ocean ridge spreading show that this formulation can robustly resolve zero‐ and non‐zero‐porosity domains simultaneously, and can be used for a large range of applications in various geodynamic settings.

The 2023 Turkey earthquake sequence involved unexpected ruptures across numerous fault segments. We present 3D dynamic rupture simulations to illuminate the complex dynamics of the earthquake doublet. Our models are constrained by observations available within days of the sequence and deliver timely, mechanically consistent explanations of the unforeseen rupture paths, diverse rupture speeds, multiple slip episodes, heterogeneous fault offsets, locally strong shaking, and fault system interactions. Our simulations link both earthquakes, matching geodetic and seismic observations and reconciling regional seismotectonics, rupture dynamics, and ground motions of a fault system represented by 10 curved dipping segments and embedded in a heterogeneous stress field. The Mw 7.8 earthquake features delayed backward branching from a steeply branching splay fault, not requiring supershear speeds. The asymmetrical dynamics of the distinct, bilateral Mw 7.7 earthquake are explained by heterogeneous fault strength, prestress orientation, fracture energy, and static stress changes from the previous earthquake. Our models explain the northward deviation of its eastern rupture and the minimal slip observed on the Sürgü fault. 3D dynamic rupture scenarios can elucidate unexpected observations shortly after major earthquakes, providing timely insights for data-driven analysis and hazard assessment toward a comprehensive, physically consistent understanding of the mechanics of multifault systems.

Dynamic perturbations reveal unconventional nonlinear behavior in rocks, as evidenced by field and laboratory studies. During the passage of seismic waves, rocks exhibit a decrease in elastic moduli, slowly recovering after. Yet, comprehensive physical models describing these moduli alterations remain sparse and insufficiently validated against observations. Here, we demonstrate the applicability of two physical damage models—the internal variable model (IVM) and the continuum damage model (CDM)—to provide quantitative descriptions of nonlinear co‐seismic elastic wave propagation observations. While the IVM uses one internal variable to describe the evolution of elastic material moduli, the CDM damage variable is a mathematical representation of microscopic defects. We recast the IVM and CDM models as nonlinear hyperbolic partial differential equations and implement 1D and 2D numerical simulations using an arbitrary high‐order discontinuous Galerkin method. We verify the modeling results with co‐propagating acousto‐elastic experimental measurements. Subsequently, we infer the parameters for these nonlinear models from laboratory experiments using probabilistic Bayesian inversion and 2D simulations. By adopting the Adaptive Metropolis Markov chain Monte Carlo method, we quantify the uncertainties of inferred parameters for both physical models, investigating their interplay in 70,000 simulations. We find that the damage variables can trade off with the stress‐strain nonlinearity in discernible ways. We discuss physical interpretations of both damage models and that our CDM quantitatively captures an observed damage increase with perturbation frequency. Our results contribute to a more holistic understanding of co‐seismic damage and post‐seismic recovery after earthquakes bridging the worlds of theoretical analysis and laboratory findings.

Dynamic rupture simulations generate synthetic waveforms that account for nonlinear source and path complexity. Here, we analyze millions of spatially dense waveforms from 3D dynamic rupture simulations in a novel way to illuminate the spectral fingerprints of earthquake physics. We define a Brune-type equivalent near-field corner frequency (fc) to analyze the spatial variability of ground-motion spectra and unravel their link to source complexity. We first investigate a simple 3D strike-slip setup, including an asperity and a barrier, and illustrate basic relations between source properties and fc variations. Next, we analyze >13,000,000 synthetic near-field strong-motion waveforms generated in three high-resolution dynamic rupture simulations of real earthquakes, the 2019 Mw 7.1 Ridgecrest mainshock, the Mw 6.4 Searles Valley foreshock, and the 1992 Mw 7.3 Landers earthquake. All scenarios consider 3D fault geometries, topography, off-fault plasticity, viscoelastic attenuation, and 3D velocity structure and resolve frequencies up to 1–2 Hz. Our analysis reveals pronounced and localized patterns of elevated fc, specifically in the vertical components. We validate such fc variability with observed near-fault spectra. Using isochrone analysis, we identify the complex dynamic mechanisms that explain rays of elevated fc and cause unexpectedly impulsive, localized, vertical ground motions. Although the high vertical frequencies are also associated with path effects, rupture directivity, and coalescence of multiple rupture fronts, we show that they are dominantly caused by rake-rotated surface-breaking rupture fronts that decelerate due to fault heterogeneities or geometric complexity. Our findings highlight the potential of spatially dense ground-motion observations to further our understanding of earthquake physics directly from near-field data. Observed near-field fc variability may inform on directivity, surface rupture, and slip segmentation. Physics-based models can identify “what to look for,” for example, in the potentially vast amount of near-field large array or distributed acoustic sensing data.

In traditional modeling approaches, earthquakes are often depicted as displacement discontinuities across zero‐thickness surfaces embedded within a linear elastodynamic continuum. This simplification, however, overlooks the intricate nature of natural fault zones and may fail to capture key physical phenomena integral to fault processes. Here, we propose a diffuse interface description for dynamic earthquake rupture modeling to address these limitations and gain deeper insight into fault zones' multifaceted volumetric failure patterns, mechanics, and seismicity. Our model leverages a steady‐state phase‐field, implying time‐independent fault zone geometry, which is defined by the contours of a signed distance function relative to a virtual fault plane. Our approach extends the classical stress glut method, adept at approximating fault‐jump conditions through inelastic alterations to stress components. We remove the sharp discontinuities typically introduced by the stress glut approach via our spatially smooth, mesh‐independent fault representation while maintaining the method's inherent logical simplicity within the well‐established spectral element method framework. We verify our approach using 2D numerical experiments in an open‐source spectral element implementation, examining both a kinematically driven Kostrov‐like crack and spontaneous dynamic rupture in diffuse fault zones. The capabilities of our methodology are showcased through mesh‐independent planar and curved fault zone geometries. Moreover, we highlight that our phase‐field‐based diffuse rupture dynamics models contain fundamental variations within the fault zone. Dynamic stresses intertwined with a volumetrically applied friction law give rise to oblique plastic shear and fault reactivation, markedly impacting rupture front dynamics and seismic wave radiation. Our results encourage future applications of phase‐field‐based earthquake modeling.

To reach Earth’s surface, magma must ascend from the hot, ductile asthenosphere through cold and brittle rock in the lithosphere. It does so via fluid-filled fractures called dykes. While the continuum mechanics of ductile asthenosphere is well established, there has been little theoretical work on the cold and brittle regime where dyking and faulting occurs. Geodynamic models use plasticity to model fault-like behaviour; plasticity also shows promise for modelling dykes. Here we build on an existing model to develop a poro-viscoelastic–viscoplastic theory for two-phase flow across the lithosphere. Our theory addresses the deficiencies of previous work by incorporating (i) a hyperbolic yield surface, (ii) a plastic potential with control of dilatancy and (iii) a viscous regularization of plastic failure. We use analytical and numerical solutions to investigate the behaviour of this theory. Through idealized models and a comparison to linear elastic fracture mechanics, we demonstrate that this behaviour includes a continuum representation of dyking. Finally, we consider a model scenario reminiscent of continental rifting and demonstrate the consequences of dyke injection into the cold, upper lithosphere: a sharp reduction in the force required to rift.