Accounting for Thermal and Poroelastic Processes in Rupture Dynamics (invited talk)
Eric M. Dunham (Harvard)
Studies of fault zones suggest that slip localizes to a narrow shear zone, less than a few hundred microns in width that is embedded within a low permeability ultracataclastic fault core. This is surrounded by a damage zone that can be several orders of magnitude more permeable than the fault core. Given this structure, and the fact that fault zones are fluid-saturated at seismogenic depths, what are the relevant physical processes that influence the evolution of fault strength during earthquakes?
We first consider poroelastic effects in the absence of thermal processes. In-plane slip compresses one side of the fault and extends the other, generating a gradient in fluid pressure that drives flow across the fault. Any asymmetry across the fault, either in elastic properties measured over the scale of the rupture or in poroelastic properties measured over the hydraulic diffusion length (as would occur if slip localizes at the boundary between the fault core and damage zone), leads to effective normal stress changes on the fault [Rudnicki and Rice, 2006]. The sign of the effective normal stress change reverses if the rupture propagates in the opposite direction, in a manner similar to the well-known bimaterial effect of normal stress changes during slip between dissimilar elastic solids. The sign of the effective normal stress change cannot always be predicted solely from the contrast in elastic properties across the fault. In numerical models with opposing elastic and poroelastic effects, we observe, as the rupture accelerates, a reversal in the sign of effective normal stress change from that predicted by the poroelastic mismatch to that predicted by the elastic mismatch, provided that the wave-speed contrast exceeds about 5-10% (the precise value depends on the poroelastic contrast and Skempton's coefficient). For faults separating more elastically similar materials, there exists a minimum poroelastic contrast above which the poroelastic effect always determines the sign of the effective normal stress change, no matter the rupture speed.
A complementary study explores coupled thermal and poroelastic effects while neglecting the undrained pore pressure changes from stresses within the fault zone (which formed the basis of the previously considered bimaterial effect). Shear heating of fluids within the shear zone locally increases pore pressure. If diffusion of heat and fluid cannot occur sufficiently rapidly, then the effective normal stress is reduced. We combine this model for thermal pressurization with a rate-and-state friction law that includes strong velocity weakening at coseismic slip rates arising from flash heating of asperity contacts. Numerical modelingreveals both crack-like and pulse-like rupture modes. The rupture mode is determined by the background stress level on the fault; below a critical level, no crack-like solutions can exist. The critical stress level is a function of both the velocity-weakening characteristics of the friction law and a "hydrothermal diffusivity factor" that quantifies the efficiency of thermal pressurization; low diffusivities favor crack-like ruptures. The width of the shear zone, even if less than a few hundred microns, also influences the rupture mode; wide shear zones favor pulse-like ruptures.
SCEC/SURE: A Finite Differences Model for Dynamic Ruptures along Rough Faults (poster with undergraduate SCEC intern)
David Belanger (Harvard) and Eric M. Dunham (Harvard)
Many earthquake source models assume that slip occurs across a planar interface between identical linear elastic materials, for which the assumption that the walls of the fault do not open implies that the normal stress remains unaltered by the rupture. However, when we examine ruptures along rough, nonplanar faults, there will certainly be changes in normal stress. Both the no-opening assumption and the assumption of an ideally elastic material response seem less reasonable. To better understand the validity of these assumption in rough fault models, it is useful to develop a model that assumes the no-opening condition and elastic response, and then to investigate the magnitude of stress changes as a function of fault roughness. If normal stress on the fault ever becomes tensile, it would contradict the no-opening assumption. Additionally, we would like to quantify the degree of roughness necessary to invoke an inelastic material response.
We are developing a dynamic rupture model in which we can parameterize and adjust the roughness of the fault as desired. First we specify a smooth curve describing the shape of the fault. Next, we generate a curvilinear coordinate system that conforms to the fault and the other boundaries of the physical domain. We then construct a mapping that transforms between curvilinear coordinates in the irregular physical domain and Cartesian coordinates in a rectangular logical domain. A regularly spaced distribution of grid points in the logical domain maps to an irregularly spaced mesh that conforms to the boundaries of the physical domain; this mapping is obtained by numerically solving an elliptic partial differential equation given a control function that specifies grid spacing at every point in the domain, a feature that allows us to "zoom in" on certain regions.
We transform the elastodynamic equation into the coordinate system of the logical domain, and numerically solve it there using second-order finite differences. The particular numerical method used to integrate the system relies on generating an orthogonal grid, one in which all coordinate lines meet at right angles. The use of an orthogonal grid allows us to avoid troublesome cross terms that would couple adjacent points on the fault.
The code exhibits second-order convergence for specified traction boundary conditions on the fault, and we are implementing rate- and state-dependent friction laws that will permit us to model dynamic ruptures along rough faults.
Attenuation of Radiated Ground Motion and Stresses from Three-Dimensional Supershear Ruptures (poster)
Harsha S. Bhat (Harvard) and Eric M. Dunham (Harvard)
Radiating shear and Rayleigh waves from supershear ruptures form Mach fronts that transmit largeamplitude ground motion and stresses to locations far from the fault. We simulate bilateral ruptures on a finite-width vertical strike-slip fault (of width W and half-length L with L >> W) breaking the surface of an elastic half-space, and focus on the wavefield out to distances comparable to L. At distances much smaller than W, two-dimensional plane-strain slip-pulse models (i.e., models in which the lateral extent of the slip zone is unbounded) [Dunham and Archuleta, 2005; Bhat et al., 2007] accurately predict the subsurface wavefield. Amplitudes in the shear Mach wedges of those models are undiminished with distance from the fault. When viewed from distances far greater than W, the fault is accurately modeled as a line source that produces a shear Mach cone and, on the free surface, a Rayleigh Mach wedge. Geometrical spreading of the shear Mach cone occurs radially and amplitudes there decrease with the inverse squareroot of distance [Ben-Menahem and Singh, 1987]. The transition between these two asymptotic limits occurs at distances comparable to W. Similar considerations suggest that Rayleigh Mach waves suffer no attenuation in the ideally elastic medium studied here. The rate at which fault strength weakens at the rupture front exerts a strong influence on the off-fault fields only in the immediate vicinity of the fault (for both sub-Rayleigh and supershear ruptures) and at the Mach fronts of supershear ruptures. More rapid weakening generates larger amplitudes at the Mach fronts.
Self-healing vs. crack-like rupture propagation in presence of thermal weakening processes: The effect of small, but finite width of shear zone (poster)
Hiroyuki Noda (Kyoto), Eric M. Dunham (Harvard), and James R. Rice (Harvard)
We have conducted rupture propagation simulations incorporating the combined effects of thermal pressurization of pore fluid by distributed heating within a finite width shear zone, and flash heating of microscopic contacts. These are probably the primary weakening mechanisms at high coseismic slip rates. For flash heating, we use a rate- and state-dependent friction law in the slip law formulation, accounting for extreme velocity weakening above a slip rate of ~0.1 m/s that depends on the background temperature, and a very short state evolution distance of ~10 microns, which is comparable to the asperity length.
Under the assumption that shear heating is confined to a mathematical plane, Noda et al. [2006 SCEC Annual Meeting] proposed a criterion for the type of solution (crack-like and self-healing pulse-like) by using an analytical solution by Rice  and modifying a theorem proved by Zheng and Rice , which state that there is a critical background shear stress, taupulse, below which an expanding crack-like solution does not exist for mode III rupture propagation. The model of slip on a plane may be valid if width of shear zone is negligibly thin compared to the other length scales perpendicular to the fault (for example, the diffusion lengths of temperature and pore pressure). In our calculations, even without thermal pressurization, the slip rate at a node at the rupture tip reaches its maximum value ~1 microsec after its rupturing, probably due to the short evolution distance. Over this short time, fluid and heat diffuses over ~1 micron, which is much shorter than the typically observed width of a shear zone of hundreds of microns [Chester et al., 2003 for Punchbowl fault, Beeler et al., 1996 for low velocity friction experiments of granite gouge, Mizoguchi et al., 1994 for high velocity friction experiments of Nojima fault gouge]. Therefore, the width of the shear zone is not negligible even below 100 microns. Our calculations show that the width of the shear zone strongly influences the type of rupture propagation, with a wide shear zone favoring a selfhealing pulse-like solution. Noda et al. [2006 SCEC Annual Meeting] reported that, using a model of slip on a mathematical plane, crack-like solutions were obtained even if the background shear stress, taub, is below taupulse defined by the steady state shear stress as a function of slip rate evaluated with the initial values of temperature and pore pressure. In the case of a heat source distributed uniformly or in a Gaussian shape over a finite width zone, crack-like solutions, including a case with taub just below taupulse, become self-healing pulse-like as the width of the shear zone is increased up to 200 microns.
Finite Difference Modeling of Rupture Propagation with Strong Velocity-Weakening Friction (poster)
Otilio Rojas (SDSU), Eric Dunham (Harvard), Steven Day (SDSU), Luis A. Dalguer (SDSU), and Jose Castillo (SDSU)
We present second- and fourth-order finite difference (FD) implementations for simulating dynamic rupture on faults with rate- and state-dependent frictional resistance, extending a FD method previously verified only for slip-dependent friction (Rojas et al, 2007). The methods are tested by modeling ruptures when the frictional stress follows a slip evolution law, with thermal weakening due to flash heating of asperities (Rice 2006) represented through a strong velocity-weakening behavior at high slip rates. We also consider friction laws with moderate velocity weakening, such as the classic aging and and (low-velocity) slip laws. In some cases, we find stiff ODE integration techniques (based on backward differentiation formulae and Rosenbrock methods) to be necessary to resolve different time scales of fault-variables (slip rate, traction, and state) in their transition from initial interseismic conditions to dynamic values. A convergence analysis is carried out using highly-resolved reference solutions from a Boundary Integral Equation Method and error metrics based upon RMS differences in rupture time, final slip, and peak slip rate, as well as L-infinity misfits of time histories of fault-variables.
Last updated: September 14, 2007