Earthquake nucleation is currently explained using rate and state stability analysis, which successfully models the behavior of laboratory simulated faults with constant thickness gouge layers. However, roughness is widely observed on natural faults and its influence on earthquake nucleation is little explored. Here we conduct frictional sliding experiments with different roughness on granite samples at upper crustal conditions (30–200 MPa). We observe a wide range of behaviors, from stable sliding to stick slip, depending on the combination of roughness parameters and normal stress. Stick slip is repeatedly observed in velocity-strengthening regimes, and increases in normal stress stabilize slip; these features are not fully predicted by current stability analysis. We derive a new instability criterion that matches our observations, based on fracture energy considerations and the size of weak patches created by fault roughness.

A central question regarding tectonic faults concerns the onset of earthquake-generating stick slip as opposed to aseismic stable sliding. This problem has been addressed in observational (Dodge et al., 1996), theoretical (Rice and Ruina, 1983), and experimental studies (Leeman et al., 2016; Scuderi et al., 2016) using the predictions of rate and state friction laws and stability analysis, where instability develops under velocity-weakening friction and low mechanical stiffness (Leeman et al., 2016; Marone, 1998; Scuderi et al., 2016). In general, most experiments have been conducted on homogeneous materials, either generating slip in a constant thickness gouge layer (e.g., Leeman et al., 2016; Scuderi et al., 2016) or on roughened cohesive rock surfaces (e.g., Passelegue et al., 2013). However, natural faults are highly heterogeneous features with variable composition, physical properties, and complex slip surface geometries (Bistacchi et al., 2011; Brodsky et al., 2016; Candela et al., 2012; Sagy et al., 2007).

In this study, we investigate the effects of heterogeneity due to the roughness of fault surfaces, and its influence on the onset of unstable sliding. Roughness is observed on faults at all scales (Bistacchi et al., 2011; Candela et al., 2012; Sagy et al., 2007), and plays a key role in faulting by controlling the size and distribution of asperities (Dieterich and Kilgore, 1994; Scholz, 1988) and the stress distribution on the fault surface (Persson, 2013; Selvadurai and Glaser, 2017). Therefore, roughness should have significant implications for both the static (Brodsky et al., 2016) and dynamic frictional strength of fault zones (Fang and Dunham, 2013), critical slip distances (Candela and Brodsky, 2016; Ohnaka and Shen, 1999; Okubo and Dieterich, 1984), and nucleation size (Ohnaka and Shen, 1999; Okubo and Dieterich, 1984). Currently only a narrow range of conditions has been investigated experimentally (Marone and Cox, 1994; Ohnaka and Shen, 1999; Okubo and Dieterich, 1984).

Here we report results of the first systematic experimental study investigating the occurrence of frictional instability under a range of roughness and normal stress conditions. We show that the combination of these parameters controls the onset of frictional instability of faults. We then provide a novel explanation based on the interaction between maximum weak patch scaling, roughness, and normal stress.

We performed 23 frictional sliding experiments with Westerly Granite (northeastern United States) subjected to stress conditions representative of Earth’s upper crust (30 MPa ≤ σn ≤ 200 MPa; σn is normal stress) and loaded in a direct shear configuration to force sliding. Samples were subjected to an initial period of 1–1.5 mm run-in before carrying out velocity steps between 0.1 and 10 µm s–1. To create different distributions of heterogeneity on the simulated faults, the sliding surfaces were axially precut and carefully polished to obtain variable degrees of roughness, characterized in terms of root mean square roughness (Zrms) using stylus profilometry measurements (see Sections DR1 and DR2 in the GSA Data Repository1).

All experiments show initial elastic loading followed by frictional rollover, where the contacting surfaces begin to slide (see Fig. DR1 and the summary of results in Sections DR3 and DR4). Once past this initial stage, the frictional strength remains relatively constant and a steady state is reached (typically requiring a displacement of 0.75–1.5 mm). The full spectrum of frictional sliding behaviors is observed, from stable sliding to seismic stick slip, across the range of experimental conditions. In several experiments, it was possible to determine the rate and state friction parameters a and b by modeling the frictional data to load-point velocity stepping during stable sliding episodes (Section DR3). Figure 1 shows examples of typical slip dynamics observed in different experiments (full plots are provided in Section DR3).

At lower normal stress (σn = 30 MPa), rougher faults (Zrms ≥ 8 µm) are observed to slide stably with velocity-neutral friction (Fig. 1A). Velocity-weakening friction and marginal stability are confined to the smoothest faults (Zrms ≤ 4.3 µm), manifested by fast stress drops during stepwise velocity increases (Fig. 1B).

When normal stress is increased to 100 MPa, smooth faults (Zrms ≤ 4.3 µm) are observed to become fully unstable with repetitive fast stick slip instabilities (Fig. 1C). Fast slip is confirmed by observations of frictional melting in scanning electron microscopy (SEM) imaging of slip surfaces (Fig. 2C; Section DR5). Intermediate roughness surfaces (Zrms = 8 µm) show marginal stability with velocity-weakening to neutral friction accompanied by slow stress drops upon increases in velocity (Fig. 1D). Rougher faults (Zrms ≥ 18.4 µm) are stable throughout the course of experimentation with velocity-strengthening friction, and abundant cataclasis observed in SEM imaging (Fig. 2D; Section DR5).

For σn > 100 MPa, sliding shows a wider spectrum of behaviors, with some unexpected results. At 150 MPa, smooth faults (Zrms ≤ 4.3 µm) remain unstable with repetitive fast stick slip cycles. All rougher faults (Zrms > 4.3 µm) are marginally stable, with evidence of fast stress drops nucleating spontaneously (without a velocity kick) or upon stepwise velocity increases, in spite of velocity-strengthening friction measurements (Fig. 1E; Sections DR3 and DR4), and evidence of frictional melt in SEM images (Fig. 1B; Section DR5). It is surprising that increasing the normal stress to 200 MPa results in the smoothest faults (Zrms ≤ 4.3 µm) becoming marginally stable (Fig. 1F). Similar behavior is also observed on intermediate roughness faults (4.3 < Zrms < 28.2 µm), which are stable with velocity-neutral to velocity-strengthening friction. Unexpectedly, given the consistently velocity-strengthening friction at lower normal stresses, the roughest fault (Zrms = 28.2 µm) is unstable with repetitive dynamic stick slip (Section DR3).

The variety of slip behaviors observed is summarized in Figure 2A, where points correspond to experimental conditions (Zrms σn), allowing approximate definition of differing frictional domains. Two characteristic trends emerge in the data (Fig. 2). First, there is a transition from stable to unstable and marginally stable slip as normal stress is increased, in accordance with the predictions of rate and state (Marone, 1998; Rice and Ruina, 1983), with the transition at increasingly higher normal stress as faults become rougher. Second, with further normal stress increase up to 200 MPa, instability is suppressed on all but the roughest fault, which becomes unstable (Fig. 2A). The occurrence of spontaneous rupture nucleation in a velocity-strengthening regime for several experiments and the second trend of the stabilizing effect of normal stress are not predicted using a standard stability analysis (Marone, 1998; Rice and Ruina, 1983).

We now discuss these results in light of rupture stability criteria and develop a theoretical model based on roughness-induced weak fault patches. Studies of natural fault surfaces show that faults have a characteristic self-affine roughness, described by a power density spectrum over nine orders of magnitude (from 10–4 to 105 m):
where α is the amplitude scaling (in m3), k is the inverse wavelength (in m–1), H is the Hurst exponent, and k0 is a normalizing factor (here 1 m–1) with α = 10–3 to 10–1 m3 and H = 0.6–0.8 (Bistacchi et al., 2011; Candela et al., 2012). At shorter length scales of <1–50 µm, this scaling diminishes and becomes isotropic as a result of plastic yielding at asperity contacts (Candela and Brodsky, 2016). From stylus profilometry measurements of our pre-experimental and post-experimental surfaces, we have a corner inverse wavelength, kmin, identified using Fourier analysis (see Section DR2), above which surfaces obey self-affine scaling. The power law parameters and kmin can be related to the root mean square of elevation such that forumla (Section DR6). Therefore, Zrms is selected as a single representative parameter of the surface statistics in Figure 2, and throughout the discussion.

Onset of rupture propagation can be interpreted either (1) in the context of rate and state dependent friction (Marone, 1998; Rice and Ruina, 1983), when a sliding patch reaches a critical size; or (2) as the consequence of stress concentration around a weak patch that may propagate unstably according to fracture energy considerations (Griffith, 1921), which have been adapted to the problem of shear cracks and earthquake faulting (Andrews, 1976; Ida, 1972).

According to criterion 1, stability is controlled by the ratio of the mechanical stiffness Kf to the frictional stiffness Kc, defined as forumla⁠, where a and b are rate and state friction dimensionless parameters and Dc is the critical slip distance. When the stiffness criterion forumla is satisfied, instability can develop; otherwise, sliding is conditionally stable (Marone, 1998; Rice and Ruina, 1983). In the case of tectonic faults embedded in an elastic medium, Kf represents the stiffness of the fault and can be expressed as forumla⁠, where G is the shear modulus, h is the linear fault dimension, and C is a dimensionless shape factor (Rubin and Ampuero, 2005). The stiffness criterion allows definition of the minimum dimension h*,
of a slip patch required for instability to develop. Rate and state friction laws and Equation 2 provide effective tools to model the full spectrum of fault slip behaviors observed across relatively homogeneous sliding interfaces such as gouge-dominated faults (Leeman et al., 2016; Scuderi et al., 2016).

However, stability criterion 2, based on fracture energy, surmises the presence of a pre-existing flaw or weak patch of finite size. Material flaws are inherent in Griffith’s (1921) original crack theory, and in the case of tectonic faults we may equate them to an elastic bridge between asperities (Figs. 3A and 3B). For earthquake nucleation, instability arises when the growth of the weak patch is energetically favorable, requiring that the strain energy release balances or exceeds the work done against residual frictional strength. According to a model of dynamic slip weakening, fault friction drops over a characteristic distance δc (Andrews, 1976; Ida, 1972); here we extend the significance of δc to low sliding velocities.

Criterion 2 allows the definition of a critical length, forumla at which a shear crack undergoes unstable failure (Andrews, 1976), where τ0 is the static shear stress on the fault, τp = µpσn is the peak stress, µp is the peak friction coefficient, and τr = µrσn is the shear stress after weakening where µr is the sliding or weak friction coefficient. To enhance similarity with h* (Equation 2), we derive a lower bound length estimate for Lc by assuming that the stress state on the fault is close to the peak stress during experiments (i.e., τ0 ≈ τp), yielding

Although the critical patch length h* and Lc share some scaling similarities, they may differ by orders of magnitude: a – b is usually small (typically −0.01) while we can expect µp – µr to be quite large. Here the frictional strength at asperities equates to µp, and within elastic bridges or zones of reduced asperity density it equates to µr. We adopt µp – µr = 0.2, as suggested by the observations of Selvadurai and Glaser (2017) on rough surfaces, which show that stress fluctuations can be as much as 40% of the peak stress.

Estimates of nucleation size for experiments showing stick-slip instability at moderate normal stress using Equation 2 for a fault of Zrms = 3.6 µm [a – b = −0.003, Dc = 5 µm, G = 50 GPa (estimated from loading curves), and σn = 100 MPa] yield h* values of ∼1 m, ∼2 orders of magnitude larger than the size of samples utilized. Following Rubin and Ampuero (2005), if we neglect the rate parameter a, we find that hb* reduces by a factor of 3 (a = 0.005), which is still an order of magnitude larger than the size of the sample. If we calculate the nucleation length using Equation 3, we obtain Lc = 1.25–3.75 cm (with G = 50 GPa, σn = 100 MPa, µp – µr = 0.2, and δc = 5–15 µm). Values of Lc obtained are in agreement with other studies that posit that the nucleation length should be smaller than the sample length (L0 = 4 cm) for laboratory stick slip occurrence (Ohnaka and Shen, 1999; Okubo and Dieterich, 1984; Passelegue et al., 2013). The values estimated here for δc = 0.05 kmin–1 (Ohnaka and Shen, 1999) are consistent with previous estimates using high-frequency strain gauges (Okubo and Dieterich, 1984) and those predicted by numerical modeling of elastic surface closure (see Section DR7). This is in contrast to values of Dc obtained during velocity steps that do not show a systematic dependence on roughness. The onset of instability observed at higher normal stress for increasing roughness is in accord with δc ∝ kmin–1.

While a stability criterion based on fracture energy (e.g., Equation 3) can explain the onset of stick slip during our experiments at low to moderate normal stress (30–150 MPa), the surprising observation that slip instability is suppressed at higher normal stress indicates the presence of some limiting process (Fig. 2A). As we discuss in the following, this behavior could be explained by considering the microphysical properties of contact asperity distribution in relation to fault zone roughness, and the associated stress heterogeneity (Scholz, 1988).

Previous studies imaged the distribution of frictional contacts with increasing normal stress and varying surface roughness in transparent materials (Dieterich and Kilgore, 1994; Selvadurai and Glaser, 2017). With increasing normal stress, contact asperities increase in number and grow, as shown in Figure 3 (Section DR6). Theory indicates that stress and asperity sizes will follow a power law distribution for a self-affine surface under load (Scholz, 1988). The asperity bridging length, λc, which is the maximum supportable elastic length or bridge between asperities, is also shown to decrease as λc ∝ σn–2 for a self-similar surface (Scholz, 1988). Here we consider a generalization of this result to any self-affine surface with 0 < H < 1 as
where forumla is a scaling factor of length dimension (Section DR6). From measurements of experimental fault surfaces, the Hurst exponent is typically 0.6–0.9 above kmin, yielding λc ∝ σn–2.5 – σn–10. In comparison, Equation 3 gives Lc ∝ σn–1, demonstrating that as normal stress increases, the bridging length will decrease at a faster rate than that of the nucleation length (Fig. 4). In extreme cases, bridges of length scale λc may represent voids, as shown in Figure 3, but more generally represent zones of reduced normal stress, or weak patches of low stiffness, filled with undercompacted gouge. These weak zones can act as stress concentrators and initiate rupture, provided that Lc < λc, as shown in Figure 4. However, with increasing normal stress, the bridges will gradually be closed and the maximum open patch will decrease until λc < Lc, and rupture nucleation is no longer possible in accordance with our experimental observations (Fig. 2). In general, instability leading to rupture nucleation vin our experiments is only observed when the nucleation length Lc satisfies the condition Lc < λc and Lc < L0. Conversely, the conditions Lc > L0 at lower normal stress inhibits stick slip; λc < Lc at higher normal stress leads to stable sliding (Fig. 4).

Our findings also have implications for the larger scale behavior of natural fault zones. In principle, the model shown in Figure 4 suggests that the transition from seismic to aseismic faulting may be controlled by the stabilizing influence of increasing normal stress upon asperities, as originally suggested by Brace and Kohlstedt (1980), in addition to currently accepted temperature-induced rheological changes (Scholz, 1988). In addition to this our results qualitatively support observations of subduction zone seismicity, where rough seafloor topography is observed to be related to creeping behavior, and smooth seafloor topography is related to seismicity and large earthquake nucleation (Wang and Bilek, 2014).

Our results highlight the key role of fault heterogeneity in earthquake nucleation. On larger scales such heterogeneity includes fault jogs, compositional contrasts, and fluid injection in addition to fault roughness, as suggested herein. These results complement rate and state friction stability analysis, providing a physical framework to include the complexity of roughness in earthquake nucleation models.

Harbord acknowledges Natural Environment Research Council (NERC) Case studentship under grant EP/M507854/1, and Nielsen acknowledges NERC capital grant CC019. We are grateful to Eric Dunham and two anonymous reviewers who significantly helped to clarify the manuscript. We thank Cyril Bourgenot, Ian Chaplain, and Leon Bowen for technical support.

1GSA Data Repository item 2017311, summary of experimental conditions (sections DR1 and DR4), surface topography measurements (section DR2), mechanical dataset (section DR3), microstructural images (section DR5), derivation of Equation 4 (section DR6), and details of elastic modelling (section DR7), is available online at http://www.geosociety.org/datarepository/2017/ or on request from [email protected].