Cataclasites are a characteristic rock type found in drill cores from active faults as well as in exposed fossil subduction faults. Here, cataclasites are commonly associated with evidence for pervasive pressure solution and abundant hydrofracturing. They host the principal slip of regular earthquakes and the family of so-called slow earthquakes (episodic slip and tremor, low to very low frequency earthquakes, etc.). Slip velocities associated with the formation of the different types of cataclasites and conditions controlling slip are poorly constrained both from direct observations in nature as well as from experimental research. In this study, we explore exposed sections of subduction faults and their dominant microstructures. We use recently proposed constitutive laws to estimate deformation rates, and we compare predicted rates with instrumental observations from subduction zones. By identifying the maximum strain rates using fault scaling relations to constrain the fault core thickness, we find that the instrumental shear strain rates identified for the family of “slow earthquakes” features range from 10−3s−1 to 10−5s−1. These values agree with estimated rates for stress corrosion creep or brittle creep possibly controlling cataclastic deformation rates near the failure threshold. Typically, pore-fluid pressures are suggested to be high in subduction zones triggering brittle deformation and fault slip. However, seismic slip events causing local dilatancy may reduce fluid pressures promoting pressure-solution creep (yielding rates of <10−8 to 10−12s−1) during the interseismic period in agreement with dominant fabrics in plate interface zones. Our observations suggest that cataclasis is controlled by stress corrosion creep and driven by fluid pressure fluctuations at near-lithostatic effective pressure and shear stresses close to failure. We posit that cataclastic flow is the dominant physical mechanism governing transient creep episodes such as slow slip events (SSEs), accelerating preparatory slip before seismic events, and early afterslip in the seismogenic zone.

Geodetic and seismological observations in the past two decades have led to the discovery of transient slip behavior at convergent plate boundaries, complementing regular subduction earthquakes (e.g., Dragert et al., 2001; Rogers and Dragert, 2003; Schwartz and Rokosky, 2007; Obara and Kato, 2016; Bürgmann 2018, and references therein). Distinct types of slip transients occur in characteristic depth intervals within and below the seismogenic zone with slip localized in the plate interface zone (e.g., Peng and Gomberg, 2010; Audet and Kim, 2016). Typical features include shallow, very low frequency earthquakes and tremor in the updip part of the plate interface, long-term silent slip events at and below the level of seismic locking, deep episodic tremor accompanied by very low frequency earthquakes, and short-term silent slip events beneath the seismogenic zone at the transition to inferred stable creep—all hereafter referred to as slip transients or family of slow earthquakes.

A plethora of studies exist showing that movement along plate bounding faults displays a multitude of transient slip events. Significant advance has been made in estimating the respective slip velocities, slip patch sizes, stress drops, and pore-pressure conditions through analysis of instrumental data and modeling (e.g., Liu and Rice, 2005; Rubin, 2008; Saffer and Wallace, 2015; Gao and Wang, 2017; Shapiro et al., 2018; Luo and Liu, 2019; Im et al., 2020). Moreover, a potential mechanical interaction between slow slip and low-frequency earthquakes has been repeatedly discussed (e.g., Frank et al., 2018), including their role in loading the region updip and triggering of large earthquakes. However, the physical processes governing this broad spectrum of deformation events are still largely unknown.

During the past decade, the drilling of active plate interfaces (e.g., Chester et al., 2013; Rowe et al., 2013; Fabbri et al., 2020; Brodsky et al., 2020) and field observations of rocks exhumed from subduction zones have provided crucial structural evidence that allows defining the major deformation processes. The field examples covered to date represent the entire depth range of a plate interface zone from the former near-trench domain to well below the downdip end of seismic locking and the adjoining transition zone that hosts many of the instrumentally observed transients (Rowe et al., 2005, 2011; Vannucchi et al., 2008; Bachmann et al., 2009a; Fagereng et al., 2010, 2014; Angiboust et al., 2015; Fagereng and den Hartog, 2017; Agard et al., 2018; Menant et al., 2018; French and Condit, 2019; Ioannidi et al., 2020; Behr and Bürgmann, 2020). Apart from pseudotachylites, representing devitrified frictional melts that are generally considered to result from earthquakes, key observations include the ubiquitous presence of scaly fabrics grading into cataclasite-lined localized faults, pressure-solution structures in phyllosilicate-rich rocks pervasively increasing with depth, and various vein types indicating hydrofracturing.

Based on these observations, a diversity of fabric types and inferred deformation mechanisms are considered to result from slip transients to date. These include (1) mineral veining exhibiting crack-seal type fibrous vein fillings in conjunction with structures indicating local dissolution-precipitation (Vannucchi et al., 2010; Fagereng et al., 2011; Fagereng and Harris, 2014; Fisher and Brantley, 2014; Dielforder et al., 2015; Ujiie et al., 2018; Cerchiari et al., 2020; Muñoz-Montecinos et al., 2020); (2) narrow principal slip zones showing foliated cataclasites exhibiting a combination of fracturing and subsequent and/or intervening solution precipitation creep (Vannucchi and Leoni, 2007; Brantut et al., 2011; Rowe et al., 2011; Angiboust et al., 2015; Fabbri et al., 2020); (3) evidence for coeval brittle-viscous flow in these rocks that may exhibit strain-rate sensitivity to fluctuations in stress leading to unstable slip (Fagereng and den Hartog, 2017); (4) mélanges indicating strain partitioning between brittle fracturing of blocks forming force bridges or enveloped by a weaker viscous matrix (Collettini et al., 2011; Kimura et al., 2012; Fagereng et al. 2014; Hayman and Lavier, 2014; Behr et al., 2018; Beall et al., 2019; Kotowski and Behr, 2019; Barnes et al., 2020; Phillips et al., 2020); (5) evidence for enhanced pressure solution and fault healing related to postseismic creep acceleration and afterslip (Gratier et al., 2013, 2014; Fisher et al., 2019a); (6) phyllosilicate-controlled change of frictional properties at low-slip velocities in mafic rocks (Ikari et al., 2020); (7) (reaction-controlled) fluid pressurization driving brittle failure at inferred slow slip velocities (Taetz et al., 2017; French and Condit, 2019; Tarling et al., 2019), and more. Based on a review of features found in the rock record of faults, Kirkpatrick et al. (2021) suggest further correlations to various aspects related to instrumentally observed properties of silent slip events drawing particular attention to material heterogeneity and geometric complexity of a fault system. Linking fabrics to slip velocities is a classic problem in structural geology of relating fabrics with slip velocities as first summarized by Cowan (1999) and recently reassessed by Rowe and Griffith (2015). An obvious challenge is posed by the massive difference of scales of the instrumentally observed dimensions of slow slip with the scales of fabric types suggested being responsible. Partly because of the limitations in successfully linking instrumental observation and fabric type, Kirkpatrick et al. (2021) question the fundamental possibility of finding a simple relationship.

The sole feature that matches the large rupture area of single slow slip events (several tens to hundreds of km2) with a displacement front propagating from several km/day to more than 100 km/h (e.g., Bürgmann, 2018) is the interconnected fault network constituting the plate interface zone itself as also emphasized by Kirkpatrick et al. (2021). Hence, it is this fault network, its materials, and fabrics that represent the one candidate storing the record of slip conditions. The dominant fabrics observed in these fault rocks suggest that deformation involves brittle cataclasis, pressure solution, and hydrofracturing partly manifested in mutually overprinting structures. Interestingly, these cataclasites have received less attention so far, although they are considered important in possibly controlling slip at sub-seismic velocities.

Here, we combine geodetic, seismological, and thermobarometric data on slip transients with microstructural observations from fault drilling of active margins and exposed shear zones constituting former plate interface faults with a continental upper plate. We first describe available structural evidence from the subduction faults providing information about fault core types and deformation mechanisms and the respective fault widths. To this end, we summarize key observations on cataclasites reported for plate interface zones along the entire depth range from the uppermost near-trench domain to a depth below the 500 °C isotherm. This includes estimates of effective pressure and flow stresses. We compare these data to instrumentally observed data from active subduction zones used to infer effective stresses, temperature conditions, displacements, and slip velocities. Using observed widths of exposed fault cores and established displacement-fault width scaling relationships, we then estimate shear strain rates during slip events. Finally, we compare this range of estimated shear strain rates to predictions from experimentally derived constitutive laws. We find an excellent match between shear strain rates predicted by stress corrosion creep (“brittle creep”) and shear strain rates calculated from instrumentally observed slow earthquake slip velocities at a stress level near the failure threshold.

2.1 Fault Zone Widths

Although the total width of the subduction fault network may increase with depth, the width of an individual subduction fault shows limited to no change with depth. Where a recent earthquake has broken the entire fault zone to the trench—for example, in the case of the 2011 Tohoku-Oki earthquake (Chester et al., 2013) —the updip end of the fault exhibits <5-m-thick scaly clay fabric with thin and sharp contacts (<1 cm) inferred to result from seismic slip (see also Vannucchi, 2019). Samples from the Japan Trench Fast Drilling Project indicate fault widths between 2 and 9 cm at a depth of ~830 m (Brodsky et al., 2020). Similarly, drilled active subduction faults at a depth of a few hundred meters (e.g., Vannucchi and Leoni, 2007; Fabbri et al., 2020) exhibit thin ultracataclasites (millimeters to centimeters thick) with no or little foliation. The ultracataclasite layers have sharp and straight boundaries forming principal slip surfaces embedded in more strongly foliated cataclasites richer in fragments. Cataclasites often display transitional boundaries toward a surrounding damage zone comprising fragmented blocks of all sizes separated by thin fractures and cataclasite bands—all embracing a width of only a few meters.

The cumulative total width of the faults in deeper parts of ancient, fossilized subduction faults amounts to a few tens of meters as summarized by Rowe et al. (2011, 2013). While their compilation only includes exhumed examples from somewhere within the plate interface zone with a focus on underplated rocks, similar observations were reported for an ancient subduction interface zone exposed in the European Alps (Bachmann et al., 2009a; Angiboust et al., 2015). Here, a completely preserved seismogenic zone is exposed embracing the base of the upper plate and the underlying mélanges and accreted lower-plate rocks (Fig. 1; see also global compilation by Agard et al., 2018). Localized shearing in this plate interface zone is mostly focused at the top of a mélange zone that is several hundred meters to up to a km wide containing shear zones forming an anastomosing network (Bachmann et al., 2009a, 2009b; Angiboust et al., 2015; Wakabayashi and Rowe, 2015; Ioannidi et al., 2020). In the Alps, a complex fault network separating large slivers of mélange was found to juxtapose upper-plate basement rocks and mélanges or block assemblages. Individual shear zone thicknesses were found between 5 and 35 m. A relatively constant fault zone width with depth in subduction zones reported here is in striking contrast to observations and conceptual models of deep-rooted continental thrust faults. In the continental crust, thrust faults are suggested to broaden significantly with depth over a similar or even smaller pressure and temperature range (Sibson, 1977, 1983; Lusk and Platt, 2020).

2.2 Fault Zone Fabrics

We focus on case studies displaying the fabrics from fault cores and surrounding hanging-wall and footwall rocks. Key structural observations of fault core fabrics and their surroundings are summarized in Table 1. The ages of these fabric types relative to each other and with respect to the stage of active subduction are inferred from mutual crosscutting relationships or age dating and assessment of the pressure-temperature conditions. Here, we only include fabric observations that are coeval to the subduction stage within the limits of uncertainty and reflect the conditions of active subduction (cf. Rowe et al., 2011, 2013). We have excluded fabrics due to reactivation during collision (Tonai et al., 2016; Fisher et al., 2019b). So far, these temporal relations could be established for the Alpine subduction zone and a few other cases (Bachmann et al., 2009a; Angiboust et al., 2015).

We use the Alpine case as a framework to relate observed structures and rock types to depth-temperature segments of the plate interface (Table 1). The analyzed fault zones reveal deformation temperatures ranging from ~100 °C in the upper parts to ~500 °C in the deepest exposed sections (Bachmann et al., 2009a; Angiboust et al., 2012, 2015; Ioannidi et al., 2020). We complement observations from the Alps with data from other exposed plate boundary systems. This allows us to estimate the range of structural variability and to extend our observations to shallower and deeper domains not preserved in the Alps.

Several types of cataclastic fabrics were found along the exposed subduction fault zone network in the Alps; these fabrics partly correlate with lithology and partly with depth (Table 1). A striking observation is that the cataclastic fabric types from exposed shear zones over the entire depth range are largely similar to those described from core material retrieved from shallow active fault zones (cf. Vannucchi and Leoni, 2007; Chester et al., 2013; Fabbri et al., 2020). However, fabrics are less well preserved at temperatures above ~300 °C and in non-basement rocks (see Table 1). At higher temperatures, pervasive pressure solution together with minor dislocation creep progressively overprint brittle fabrics (e.g., Wassman and Stöckhert, 2013; Angiboust et al., 2015).

Pelitic rocks abundant in shallow subduction faults display a characteristic scaly fabric with slight striations on the slip surfaces (cf. Vannucchi, 2019, and references therein; also called “localized deformation zones” by Bachmann et al., 2009a; Figs. 1A and 2A). Similar fabrics are found in more competent rocks (e.g., paragneiss basement and dolomites) at greater depth (Fig. 2C). Extremely fine-grained foliated ultracataclasites tend to occur with sharply bounded principal slip zones as described in many fault studies (Figs. 2B, 2G, 3A, 3B, and 3D; Chester and Logan, 1987; Chester and Chester, 1998; Kitamura et al., 2005; Kondo et al., 2005; Rowe et al., 2005, 2011; Kimura et al., 2013; Wakabayashi and Rowe, 2015). At greater depths, folding of these straight slip zones may occur alongside progressive overprint by pressure-solution foliation in a phyllosilicate-rich environment. Rowe et al. (2005) found this cataclasite type to eventually grade into pseudotachylites, showing injection veins (Figs. 2F, 3G, and 4C). In addition, this type is progressively reworked toward depth in a zone widening to form a very fine grained, foliated gouge, up to several meters thick, with a very low clast-matrix ratio (Figs. 2D and 3C). Individual straight to curved striated slip zones with a width of a few millimeters to centimeters are nested in this wider cataclasite zone, mostly at the boundaries of this main fault zone and along splay faults into the host rock. At temperatures>350 °C, this fault type becomes progressively overprinted in sedimentary protoliths by pressure solution and mylonitization (in quartz-rich rocks) with rare ultracataclasite bands preserved (Figs. 2G, 3B, and 3C).

In quartz-feldspar–dominated basement rocks at the base of the upper plate or in large entrained basement blocks, cataclastic shear zones are often less localized with a high clast-matrix ratio. These are generally foliated with stronger grain-size comminution occurring only locally. The cataclasites may attain a thickness of several tens of meters and mostly lack sharp and straight boundaries (Figs. 2E and 3H). For high quartz content, these basement rocks may be mylonitic (Figs. 2F, 3E, 4E, and 4F), containing devitrified pseudotachylite in some cases. This attests to significant structural and material heterogeneity along strike of the shear zones.

In general, subduction zones display a broad range of rock types with widely differing mineralogy. This is reflected in a broad range of fabric types. Overall, mafic lithologies are more prone to brittle fracture and cataclasis (or preservation thereof, cf. Fig. 4J) and constitute a large fraction of the blocks, while felsic lithologies (mainly granitoid and amphibolite-paragneiss series from the hanging wall and flysch from trench deposits) are subject to dissolution-precipitation creep (Figs. 3B, 3D, 3F, and 3H). Felsic rocks also dominate the main cataclastic shear zones. Mineral reactions or fluid-mineral reaction textures become more prominent with depth and temperature. This is evidenced by phyllosilicate growth in the cataclasites, breakdown reactions of feldspar, amphibole, biotite, etc., and the growth of pressure fringes around competent clasts.

At elevated temperatures (T = 400–550 °C), a competition between pressure solution and dislocation creep in shear zone rocks has been found. This is evidenced by the observation of pressure-solution fringes alternating with monomineralic segregates in blueschist-facies micaschists and metabasalts (e.g., calcite, quartz) with a marked lattice-preferred orientation (e.g., Stöckhert, 2002). Locally, foliated cataclasite layers are still visible, characterized by very strong grain-size reduction, a very dark groundmass, and finely comminuted host fragments floating in the fault zone. At eclogite-facies conditions (T>500 °C), a competition between several deformation mechanisms remains present for mafic lithologies where matrix clinopyroxene takes the bulk of the deformation via both diffusion and dislocation creep, while garnet rotates within clinopyroxene (e.g., Philippot and van Roermund, 1992; Godard and van Roermund, 1995; Keppler, 2018; Yamato et al., 2019). Locally, at nearly lithostatic pore-fluid pressures, some of the mafic layers may develop cataclasites or breccias within ductile shear zones (e.g., Pennacchioni, 1996; Angiboust et al., 2012, 2017; Hertgen et al., 2017; Locatelli et al., 2018). Above 500 °C, rocks with felsic lithologies deform mostly by dislocation creep of micas and quartz at relatively low stress levels. At nearly lithostatic pore-fluid pressures the formation of cataclasites is likely caused by increasing slip rates. In mafic-ultramafic lithologies, potential slip rates are likely controlled by serpentinization reactions and the corresponding evolution of fluid networks at the boundary to the mantle wedge (Mizukami et al., 2014). In contrast, at fluid-deficient conditions, seismic instabilities (as indicated by the finding of pseudotachylites) may develop both within felsic and mafic lithologies (e.g., Menant et al., 2018; Pennacchioni et al., 2020).

In most analyzed plate interface faults, pressure-solution fabrics are associated with veining and hydrofracturing in cataclasite and mélange zones (Wassmann and Stöckhert, 2013), starting close to the trench and increasing downdip (Figs. 2D, 3F, 4A, and 4B). In the Alps, veining observed in the continental basement rocks of the upper plate subducted to depths>30 km is restricted to the first tens of meters above the paleo-interface. Bachmann et al. (2009a) suggested that retrograde hydration reactions in upper-plate basement rocks reduce permeability focusing fluid flow into the underlying subduction shear zones and mélanges. Veining may cut the foliated fabric of the fault core cataclasites at both high angles as well as being nearly parallel to the foliation fabric (Bachmann et al., 2009a; Dielforder et al., 2015; Phillips et al., 2020). Moreover, re-worked vein fragments are ubiquitous among the clast components of the cataclasites (Figs. 2D and 4D). In some cases, fragmented vein material in the cataclasites crosses the fault damage zone with no to minor offsets, indicating synchronous cataclastic flow and hydrofracturing (Angiboust et al., 2015). In the Alpine case, the veins in the main shear zone and the immediately underlying mélange units exhibit little fibrous vein growth involving crack seal textures. From the frequently found preservation of crystal faces in the vein fill, Bachmann et al. (2009a) argue that blocky crystal growth mostly occurred into an open fluid-filled cavity (e.g., Bons et al., 2012; Fagereng et al., 2018; Epstein et al., 2020), requiring rapid opening and lithostatic pore pressures during mineral growth. In shear zones within sedimentary accretionary complexes and mélanges, fibrous crystal growth and crack seal textures are more common than blocky crystal growth forming shear and tensile veins (Fagereng et al., 2010; Fagereng and Harris, 2014; Dielforder et al., 2015; Raimbourg et al., 2018; Phillips et al., 2020). Finally, it has been shown that at all depths, vein geochemical composition can differ from host composition (e.g., Bachmann et al., 2009a, 2009b; Angiboust et al., 2014b, 2015; Jaeckel et al., 2018) thus indicating long-range transport of released fluids. This disequilibrium is also witnessed by the identification of warmer incoming fluids in ancient vein systems infiltrating shallower plate interface material (e.g., Vrolijk et al., 1988) and episodic, incremental precipitation of solutes on their way up along veined and hydrofracs pathways (e.g., Taetz et al., 2017; Muñoz-Montecinos et al., 2021a).

2.3 Stress Estimates from Microstructure Observations and Modeling

Stöckhert (2002) and Wassmann and Stöckhert (2013) report extensive evidence for the role of dissolution precipitation creep along the plate interface zone identifying characteristic fabrics. They report original pre-subduction textures preserved in blocks enclosed within a deformed matrix showing almost no microstructural evidence for crystal-plastic flow. Based on these observations, these authors posit very low shear and effective stresses ~1–10 MPa governing flow in the subduction mélange. Close to the trench, shear stresses from core analysis drilling the Tohoku-Oki earthquake were found to be as low as 1.5 MPa (Brodsky et al., 2020). At depths of 25–35 km, flow stresses of ~10 MPa (Platt et al., 2018) to 30 MPa (Angiboust et al., 2015) have been estimated for deeper parts of the subduction shear zone employing paleopiezometry on quartzite or chert blocks embedded in the mélange. These values likely represent an upper bound of former shear stresses as deformation is focused in the mélange matrix or at the margins of the blocks. In the absence of evidence for crystal-plastic deformation, these inferred very low stresses require very low frictional strength of the faults and cataclastic shear zones with an effective coefficient of friction of well below 0.1. These very low values along the entire plate interface fault zone are in agreement with friction coefficients from large continental plate boundaries such as the San Andreas fault (Carpenter et al., 2015). However, these low friction values are in contrast with much higher values estimated for other large faults in the continental crust (e.g., Sibson, 1983; Lusk and Platt, 2020). These overall low friction values reported for subduction faults indicate that pore-fluid pressures were likely close to lithostatic. This is in line with ubiquitous evidence from exposed sections for hydrofracture formation increasing downdip as well as with the growth of blocky crystals in open fluid-filled cavities (e.g., Fisher, 1996; Bachmann et al., 2009a). Moreover, fault fabrics show tensile veins have formed repeatedly during slip at acute as well as very high angles to faults and cataclasite fabric, indicating repeated changes in stress state (Vannucchi and Leoni, 2007; Bachmann et al., 2009a; Brodsky et al., 2020; Phillips et al., 2020).

3.1 Subduction Zone Slip Rates and Displacements of Deformation Transients

The increasing wealth of instrumental observations from geodesy and seismology suggests that deformation along subduction zones displays superposition of transient slip events within a wide range of characteristic duration times. Several key geometric and kinematic properties of this transient behavior have been observed or inferred such as slip velocities, the size of slip patches, the duration of the events, total displacement, and more (e.g., Ide et al., 2007; Rubinstein et al., 2010; Bedford et al., 2013; Nishimura, 2014; Bostock et al., 2015; Thomas et al., 2016; Jiang et al., 2017; Frank et al., 2018; see Bürgmann, 2018, and Behr and Bürgmann, 2020, for an extensive summary of findings). Here, we focus on kinematic and physical properties relevant to our analysis.

At shallow depth above the locked seismogenic zone, slow slip is often accompanied by seismic tremor. Slow slip displacements are on the order of several centimeters on patch sizes of several 10 km width at slip velocities typically ~1–2 × 10−8 m/s (Nishimura, 2014; Yamashita et al., 2015; Jiang et al., 2017). In the seismogenic zone and immediately below, long-term silent slip events may last days to months (Ide et al., 2007; Bostock et al., 2015; Frank et al., 2015). Toward even deeper parts of the subduction interface, deep episodic slip and tremor show similar slip velocities of ~10−8 m/s, and transient slip has been observed to last for several months. However, individual episodic tremor and slow slip (ETS) events may involve slip patches with a length up to 200 km and much larger than at shallow depth (Ide et al., 2007; Rubinstein et al., 2010; Wallace and Eberhart-Philipps, 2013). Bruhat and Segall (2016) and Gao and Wang (2017) note that a gap may separate the locked seismogenic zone from deeper zones displaying slip transients. Only in some cases, this gap seems to be filled by deep long-term silent-slip events. Finally, using the highly resolved afterslip evolution of the 2010 Mw 8.8 Maule earthquake in Chile, Bedford et al. (2013) found afterslip to be composed of slip bursts within and beneath the zone slipping coseismically. The Maule afterslip bursts exhibit the same characteristics as silent-slip events, both in terms of patch size, displacement, and duration of slip events. Overall, the slip velocities reported for this entire group of transients range in a narrow band between 10−8 to 10−6 m/s (e.g., Ikari et al., 2015, Montgomery-Brown and Syracuse, 2015; kinematic properties summarized in Table 2 and Fig. 5A). Although it has been noted that the duration of transient slip at any location on the fault—the rise time—is likely shorter than the duration of the entire slow slip event (e.g., Bürgmann, 2018), it is unlikely that slip velocities are underestimated (e.g., Behr and Bürgmann, 2020). Also, slip propagation and subsequent secondary slip in the updip, downdip, and along-strike directions may follow complex modes lasting for extended periods of time (see Behr and Bürgmann, 2020, for extensive description).

The range of slip velocities reported for slip transients is narrow (Fig. 5A) irrespective of vastly varying temperatures, pressures, and lithologies involved in active subduction zones. The slip velocities of low- and very low frequency earthquakes range between 10−5 to 10−3 m/s (Ito et al., 2007; Bostock et al., 2015; Ghosh et al., 2015; Thomas et al., 2016). This is several orders higher than slow slip but much lower than slip velocities of regular earthquakes. A potential transition between slow slip and long-term creep at the plate interface accommodating convergence is beyond current instrumental resolution (e.g., Page et al., 2009; Ide and Yabe, 2014; Bürgmann, 2018).

Low- and very low frequency earthquakes are now found to be accelerating bursts temporally and spatially linked to short-term slow slip events—often radiated from the same rock volume (e.g., Ito et al., 2007; Araki et al., 2017; Bürgmann, 2018; Frank et al., 2018; Nakano et al., 2018; Hall et al., 2019, and references therein; Hutchison and Ghosh, 2019). Similarly, aftershocks of the 2010 M8.8 Maule earthquake in Chile are closely linked to afterslip bursts in space and time (Bedford et al., 2013; Lange et al., 2014). This suggests that individual slip events may show a variation of slip velocities that, in some cases, may trigger earthquakes locally at asperities. Conversely, long-lasting slow slip events typically located around and below the downdip limit of seismic locking (Obara and Kato, 2016; Gao and Wang, 2017) appear to occur without (very) low-frequency earthquakes.

While early observations of slow slip transients were related to slow subduction of young oceanic crust, more recent observations indicate that slow slip occurrence and rates seem independent of slab age, rate of subduction, mode of mass transfer—i.e., subduction erosion versus accretion, etc. (e.g., Schwartz and Rokosky, 2007; Obara and Kato, 2016; Bürgmann, 2018). However, major topographic oceanic plate features, such as seamounts or ridges, appear to promote the occurrence and spatial propagation of aseismic creep and silent-slip events (Wang and Bilek, 2014; Yamashita et al., 2015).

3.2 Seismic Velocities Vp/Vs Ratios and Inferred Pore-Fluid Pressures

Dense local seismometer networks have provided high-resolution studies showing plate interface zones to correlate with elevated Vp/Vs ratios >2 (e.g., Kodaira et al., 2004; Audet et al., 2009; Haberland et al., 2009; Kato et al., 2010; Audet and Bürgmann, 2014; Audet and Kim, 2016; Audet and Schaeffer, 2018; Bürgmann, 2018, and references therein). Peacock et al. (2011) showed that high Vp/Vs ratios require near-lithostatic pore-fluid pressure, pf, in conjunction with significant porosity of up to 4% at the Cascadia margin. From Vp/Vs data, Moreno et al. (2014) predicted very high values for the pore-fluid pressure (0.85< λ ≤1) in the lower seismogenic zone at the Chile margin. Their study shows that the geodetic locking degree inversely scales with the pore-fluid pressure ratio. Also, Saffer and Tobin (2011), Kitajima and Saffer (2012), Saffer and Wallace (2015), and Arnulf et al. (2021) find similarly high pore-fluid pressure ratios for the shallow updip parts of plate interfaces. However, Arnulf et al. (2021) not only find significant small-scale heterogeneity in fluid pressure using reflection seismology, but also point out that the methods employed to obtain these values may have limits resolving the physical state within the thin domains hosting the main fault. On average, regions of expected high pore-fluid pressure correlate spatially with well-localized slow slip.

However, all data from the studies cited above were collected at varying stages of the seismic cycle and may not be representative of the entire seismic cycle. Recent findings at the Japan margin explored the seismological response above the subduction zone megathrust over time. Earthquakes are found to be related to the slow slip event cycle and are likely controlled by fluctuating pore-pressure cycles (Nakajima and Uchida, 2018). Such variation during the seismic cycle has also been suggested by Frank et al. (2017) from the analysis of the 2015, Mw 8.3 Illapel afterslip and aftershock evolution. Using local seismicity above slow slip patches from the Hikurangi margin, Warren-Smith et al. (2019) were able to demonstrate that stresses vary episodically over complete slow slip cycles. This is attributed to cyclic hydraulic loading of the system. Continuous dehydration of the subducting slab induces elevated pore pressures causing fracturing, dilation, and subsequent pore-pressure drop. Interseismic fault healing via precipitation drastically reduces permeability leading again to pore-pressure buildup (e.g., Miller and Nur, 2000). These findings underscore the transient nature of the hydraulic state of the plate interface in the seismic cycle (cf. also p-wave velocity fluctuations correlating with silent-slip events; Gosselin et al., 2020). Stress drops observed for some of the slow slip events including very low frequency earthquakes are very low (<10 kPa). This may support very low effective stresses and near-lithostatic, pore-fluid pressures in the slipping domains. During the loading and slip periods, pore-fluid pressures are expected to vary (Saffer and Tobin, 2011; Bostock et al., 2015; Ghosh et al., 2015; Saffer and Wallace, 2015; Thomas et al., 2016; Warren-Smith et al., 2019; Behr and Bürgmann, 2020, and references therein; Gosselin et al., 2020). This may be in line with the variations of the stress state from the pre- to the postseismic stage as observed by Brodsky et al. (2020) for the Tohoku-Oki earthquake. Finally, from a global compilation of force balance estimates at convergent margins, Dielforder et al. (2020) find effective friction to be below 0.08 everywhere requiring pore-fluid pressure ratios, λ = pfgz, to be beyond 0.9 on average.

3.3 Temperature Estimates

The seismogenic zone of the plate interface has been known to be approximately limited by the 100–150 °C isotherm at the updip end and by ~350 °C at the downdip end of locking and seismogenesis (Hyndman et al., 1997; Oleskevich et al., 1999; Peacock and Hyndman, 1999). Hence, the shallow-slip transients occur at temperatures below 150 °C, while those near and below the downdip end would occur at temperatures around and above 350 °C, and at pressures of ≥1–1.5 GPa. Employing thermal modeling of a range of convergent margins that are known to produce slip transients, Gao and Wang (2017) showed that deep long-term silent slip occurs at a temperature range of 300–500 °C, while episodic slip and non-volcanic tremor (NVT) is limited to the 350–550 °C range.

To convert instrumentally observed slip velocities to shear strain rates, we estimate the width TF of the fault-core during a slip event. This requires some assumptions: (1) reported displacements represent a single shear event; (2) observed slip was localized in the fault core; and (3) finite fault displacement scales with fault-core width. Fault width TF has been found empirically to scale with displacement D as TF ≈ 0.01D (see Scholz, 1987, 1997, 2019; Torabi and Berg, 2011). Because real faults have ruptured multiple times, their fault core is predicted to grow in width according to the above scaling relationship. Hence, calculating shear strain rates using a single event fault-core width likely underestimates fault width. This assumption results in an upper strain rate bound. Contributing to the uncertainty in estimating shear strain from slip data is the observation that fault cores often display multiple narrow slip bands, each of which may have resulted from a single or a few slip events. This is found in exposed shear zones and core material from fault zone drilling. Hence, individual slip events may occur along different slip bands forming in the fault core (e.g., Chester and Chester, 1998; Sibson, 2003; Childs et al., 2009; Rowe et al., 2011, 2013; Chester et al., 2013; Fabbri et al., 2020). Only a few studies indicate earthquake-related slip to possibly associate with wider fault zones (e.g., Rowe and Griffith, 2015).

Torabi and Berg (2011) suggest that displacement-thickness scaling of faults holds for many different lithologies, but a slight difference may exist between large displacement faults (D>10 m) and those showing small displacements at the m-scale and lower. The maximum variability in their relation of fault core thickness to displacement encompasses two orders of magnitude, with the vast majority of data falling within one order of magnitude (cf. similar results by Childs et al., 2009). Sibson (2003) finds coseismic fault thicknesses of individual events exhibiting principal slip surfaces, slickensides, pseudotachylite, etc. to rarely exceed a few cm, and to mostly be at the mm-level for single events. Field observations reported here for likely single slip events show similar fault thicknesses. However, in general, the width of a principal slip zone associated with a single seismic event is difficult to constrain (see Rowe et al., 2013; Rowe and Griffith, 2015). Hence, shear strain rates calculated here as based on the single-event assumption may be biased resulting in an upper bound of true strain rates. From the total thickness of individual long-lived fault cores reported (1–10 m), and assuming homogeneous distribution of shear within the entire fault core, shear strain rates may be smaller by up to 2–3 orders of magnitude in individual cases, thereby forming a lower bound.

Low-frequency earthquakes and slip transients are inferred to be closely related in space and time with (V)LFEs possibly occurring on the same fault utilized by long-lasting slip transients (e.g., Ito et al., 2007; Ide and Yabe, 2014; Ghosh et al., 2015). Interestingly, shear strain rates,graphic, from earthquakes (graphic ≈ 101−103/s) including low-frequency events are of the same order of magnitude within uncertainty. Earthquake strain rates are separated by several orders of magnitude from slow and aseismic slip transients (graphic ≈ 10−3−10−5/s) suggesting that different physical mechanisms govern slow slip events and dynamic earthquake ruptures, respectively (Fig. 5). We also observe that all aseismic slip transients within uncertainty display a rather narrow range of shear strain rates. This suggests that the effect of depth, temperature, and also of lithology on the localized aseismic slip is rather small.

5.1 Methods

We estimate maximum strain rates assuming that slip rates are controlled by specific deformation mechanisms inferred from microstructural observations of exposed cataclasite-bearing shear zones. This complements similar estimates published previously (Hayman and Lavier, 2014; Fagereng and den Hartog, 2017; French and Condit, 2019). We employ experimentally derived constitutive laws to explore the contribution of various deformation mechanisms to slow slip. We apply flow laws to the full spectrum of depth and temperature ranges drawing from observations of shear zone fabrics and composition. Where several flow laws are available, we chose relations yielding the highest strain rates at given conditions to conform to the instrumentally observed velocities and strain rates estimated from single event rates.

Motivated by microstructural observations, we assume that quartz dislocation creep controls the deformation of the felsic continental basement rocks just below the brittle-ductile transition. Since wet quartz is weak compared to the plastic flow of other silicates (Bürgmann and Dresen, 2008), a quartz flow law will predict an upper bound to the rate of plastic deformation. Observations from the exposed Alpine subduction zone indicate a partial contribution of dislocation creep in felsic as well as carbonate-rich rocks at temperatures above 300–350 °C. To estimate strain rates, we employ the approach adopted by Tokle et al. (2019) based on their combination of geologic and experimental observations. Predicted strain rates from extrapolated constitutive laws of Luan and Paterson (1992), Gleason and Tullis (1995), and Rutter and Brodie (2004) with a stress exponent of n»4 and assuming high water fugacity predict fairly similar strain rates (see Table 3 for flow laws used here).

The dominant matrix fabric observed at nearly all depths is a foliation formed by solution-precipitation creep (SPC; e.g., Bachmann et al., 2009a; Wassmann and Stöckhert, 2013). This fabric dominates coarser-grained sheared sediments, hydrated mafic rocks forming greenschists, and the fine-grained cataclasites of more localized faults and/or shear zones. Solution-precipitation creep has been investigated, but the controlling micromechanisms and constitutive laws are still a matter of considerable debate (see detailed overviews by Gratier et al., 2013, and Wassmann and Stöckhert, 2013). Although compositional variability and grain-size evolution play a strong role in a convergent plate interface setting, nearly all rocks forming the matrix, the shear zones, and a large fraction of the enclosed blocks are dominated by phyllosilicates. These are either of primary—i.e., mostly sedimentary—origin or formed as a hydration mineral reaction product in the more mafic or feldspar-dominated basement lithologies. Microstructural evidence for pressure solution is only found to a lesser degree in mélange blocks, such as carbonates, ultramafic blocks, and weakly altered metabasalts. Cataclastic shearing with the formation of fine-grained gouges favors phyllosilicate growth and a massive increase in grain surface area. This may increase the contribution from SPC or even permanently change the dominant deformation mechanism (cf. Angiboust et al., 2015 for evolution from quartz-controlled mylonites to fractured and foliated cataclasites dominated by brittle grain-size reduction and pressure solution).

To estimate strain rates of the complex subducted rock assemblages and the main shear zone, we chose constitutive laws that best represent the typical lithology of quartz-feldspar–rich sediments with a large phyllosilicate fraction. Also, we consider a range of grain sizes and porosities typically found in rock samples from a wide spectrum of subduction shear zones (1–100 µm). For porous sediments (~15% porosity), we use the flow law of Niemeijer et al. (2002). The flow law predicts compaction in quartzose rocks at temperatures between 150 and 600 °C where dissolution is the rate-limiting process. For pressure solution in granitic rocks with near-zero porosity, we use the flow law of Rutter (1983) suggesting that deformation is rate-limited by diffusion. Granitic rocks and quartz-feldspar-mica-sediments reflect the compositional range encountered in the subduction shear-zone matrix as well as within the base of the continental crust overlying the plate interface. Within the uncertainties involved in the pressure-solution flow law and mineralogical composition, the predicted strain rates yield an upper and lower bound on the contribution of SPC in accommodating creep on the plate boundary fault. A comparison shows that the constitutive law suggested by Gratier et al. (2009) is close to that of Rutter (1983), albeit showing a lower T-dependence. The rather porous materials studied in the compaction experiments by Niemeijer et al. (2002) exhibit strain rates higher by 1–3 orders of magnitude that strongly increase with temperature. These latter strain rates therefore most probably form an upper bound only valid during the short period of elevated fracture porosities following earthquake rupture.

Next to abundant fabrics showing evidence for SPC, cataclasite microstructures are most frequently observed. Explicit constitutive laws that quantitatively constrain strain rates in cataclasites are scarce (cf. Brantut et al., 2012; Chen and Spiers, 2016). Based on experimental work, Brantut et al. (2012) recently developed a model describing cataclastic flow (by fracturing, sliding, and rolling of grains) in rocks of granitic, sedimentary, or even basaltic composition—perfectly reflecting the dominant rock types described here. Sub-critical crack growth accommodated by stress corrosion—or “brittle creep”—likely is the rate-limiting process in all creep regimes close to the short-term failure strength of rocks. Crack growth at differential stresses below 70%–75% of the failure strength results in negligibly slow strain rates. However, accelerating crack growth results from accelerating crack interaction and coalescence at elevated stresses approaching failure. While no explicit temperature limit is given, Brantut et al. (2012) assume strength to change only moderately with temperature. Here, we chose their data for granite as the representative rock composition for the faults observed (see Table 3) and simply extrapolate their flow law to the explored maximum temperature of 500 °C. Brantut et al. (2012) also find that pressure-solution creep in porous rocks is enhanced with increasing depth and decreasing strain rates. We here use the constitutive law suggested by Brantut et al. to estimate potential brittle creep strain rates controlled by a combination of crack growth through stress corrosion, rolling, and sliding of grains in gouge at differential stresses immediately below the short-term brittle-frictional failure strength.

We first estimate the maximum differential stress operating at the respective depth section in the subduction zone, which depends on effective pressure. To this end, we assume friction is limited by Byerlee's law assuming an average friction coefficient and an average density for the forearc (for physical properties and model parameters used, see sources of constitutive laws given in Table 3). To estimate the effect of pore-fluid pressure on effective stress, we show both, the result for a pore pressure ratio of λ = 0.7 reflecting the lower end of reported plate interface pore pressures (Saffer and Wallace, 2015) and for a close to lithostatic pore pressure at the upper end (Brodsky et al., 2020). The latter is usually considered closer to the values estimated for the domains exhibiting slow earthquakes (see above summary).

For the estimated stresses, thermodynamic boundary conditions, and respective deformation mechanisms, we then calculate strain rates using the flow laws listed in Table 3. We finally convert axial strain rates obtained from the constitutive laws into shear strain rates using the formalism by Ramsay and Huber (1983; Equation 2.7):
formula

5.2 Results

To estimate the contribution of the different deformation mechanisms identified from exposed sections of plate interfaces, we present a modified deformation mechanism map. With this map, we compare the estimated strain rates to instrumental observations. We plot shear strain rates against a generic pressure and temperature gradient along a plate interface zone (Fig. 6) using an average of the temperature-depth models of Syracuse et al. (2010). Alongside, we include the depth domains that host the different instrumentally observed slip transients according to Obara and Kato (2016) and Gao and Wang (2017). In this graph, we also include an estimate of the typical strain rate required to fully accommodate the plate convergence rates. Here, we assume homogeneous distribution of convergence (1–10 cm/yr) in a network of faults with a minimum total core width of 1 m and a maximum of up to 100 m. This choice reflects findings in natural systems such as the Tohoku-Oki earthquake (Chester et al., 2013) or exposed plate boundary fault networks (Bachmann et al., 2009a; Rowe et al., 2013; Angiboust et al., 2015; Wakabayashi and Rowe, 2015). The focus here is on the fault cores since the enveloping matrix or damage zone forming most of the subduction shear zones likely experienced significantly lower strain rates (e.g., Beall et al., 2019; Ioannidi et al., 2020). Finally, Figure 6B also includes the calculated shear strain rates for the different types of slow earthquakes along the plate interface.

Uncertainties in extrapolated strain rates are due to errors of flow law parameters, such as activation energy, stress exponents, etc. In addition, flow laws are defined for steady-state experimental conditions (i.e., for secondary creep) and may not adequately quantify rates for slow slip events that by their nature are transients. For brittle creep, Brantut et al. (2012) additionally discuss values for individual physical properties reported from other experimental studies. These tend to reduce our strain-rate estimates by 1–2 orders of magnitude. Moreover, the constitutive law of Brantut et al. is based on experiments performed to maximum temperatures of 300 °C at a limited range of pore pressures. Hence, while recent observational evidence suggests that brittle creep is an important deformation mechanism at still higher temperatures and elevated pore-fluid pressures, extrapolation of this flow law to greater depths and near-lithostatic pore pressures invariably adds uncertainty. For example, eclogite-facies metagabbro cataclasites (or breccias) occur coeval with mylonites in the Alpine plate interface zone that equilibrated at 550 °C and 2.6 GPa. The rocks display multiple veins suggesting repeated hydrofracturing events (see above; Philippot and van Roermund, 1992; Angiboust et al., 2012; Locatelli et al., 2018). Second, recently observed slow slip events at even greater depths within the upper-mantle domain exhibit similar slip velocities of ~10−8 m/s as known to occur updip (Bedford et al., 2020; Fig. 5A).

Shear strain rate estimates of slip transients span a narrow range between 10−3–10−5 s−1 (Fig. 6). This is roughly five orders of magnitude slower than strain rates expected from seismic slip but more than five orders faster compared to constant shearing at plate convergence rates. The field data suggest that strain rates of slip transients remain similar irrespective of depth and temperatures. Assuming that these slip events are accommodated by a similar dominant deformation mechanism, this already suggests that any mechanism controlling viscous deformation will likely have rather small thermal activation energies, as have been suggested for brittle creep of granite (50 kJ/mol, Brantut et al., 2012). Comparing strain-rate estimates from field data and extrapolated flow laws strongly suggests that the full depth range from near trench to below the transition zone will be dominated by either pressure-solution creep (in fine-grained rocks) or by brittle creep (Fig. 6). This is in excellent agreement with microstructure observations of fault zone rocks from all exposed plate interface zones. In fact, at critical differential stresses close to failure—i.e., the stress levels assumed in Figure 6—resulting strain rates by SPC and brittle creep suffice to accommodate plate convergence. However, the extrapolated strain rates in Figure 6 are strongly sensitive to effective stress and reported (near-) lithostatic pore-fluid pressures—except for brittle creep for which available experimental data do not indicate a relevant sensitivity.

For peak differential stresses close to failure, the extrapolated flow law for brittle creep predicts shear strain rates that are in excellent agreement with instrumentally observed shear strain estimates of slow slip. In contrast, strain rates predicted using existing flow laws for SPC are more than six orders smaller. At best, extrapolated SPC rates agree roughly with long-term plate convergence rates in deeper parts of the plate interface zone. This indicates that slip transients are too fast to be accommodated by SPC except—possibly—for deep afterslip in highly fractured and porous rock at lower fluid pressures. However, deep afterslip is likely a short-lived process in the postseismic period.

In a second step, we sum shear strain rates for all key processes observed over the full temperature spectrum showing the potential evolution during the complete loading cycle from low differential stress (reflecting the immediate postseismic stage following a megathrust rupture; Brodsky et al., 2020) to loading at the effective strength limit. Stress is normalized to the stress at failure,
formula

where σ is the stress parameter used in the flow law and σf is the stress at the failure limit. For dislocation creep, differential stress at failure is σf = 2τf; for solution precipitation creep and brittle creep, the largest principal stress at failure is σf = σn (1 − λ) + τf (where τf is shear stress at failure and σn is normal stress). Shear stress at failure is derived from the Mohr-Coulomb law corrected for effective pressure (τf = C0 + µσn(1 − λ)) and using the above parameters for friction coefficient and rock density for the calculation of normal stress.

The resulting shear strain rates of key processes observed over the range of temperatures between <100 °C to ~500 °C and assumed pore-fluid pressures λ yield two dominant regimes (Fig. 7). Fine-grained fault rocks in a low-stress regime at average to slightly elevated pore pressures (Fig. 7A) are expected to deform by solution precipitation creep (with possibly some contribution of dislocation creep at high temperatures). At lower temperatures and differential stresses, creep rates are predicted to be too small to accommodate plate convergence, suggesting partial locking. Increasing pore pressure to near-lithostatic values will strongly enhance this shift toward locking in this low-stress state (Fig. 7B). At high normalized differential stresses exceeding 70% of the failure limit, brittle creep will strongly dominate deformation behavior at all pore pressures. Brittle creep and/or cataclastic flow is predicted to allow for creep strain rates sufficiently fast for slow slip transients. Hence, kinematic locking is shown to be controlled chiefly by the temperature and state of fluid pressure at differential stresses below some 60%–70% of the short-term strength. This relationship reverts with an advanced stage of the loading cycle, or with near-lithostatic pore pressures by reducing the effective strength of rocks. As a consequence, both the loading as well as the hydraulic cycles will promote brittle creep accommodating cataclastic flow toward later stages of the seismic cycle.

6.1 Brittle Creep and Slow Slip

Instrumental observations of active subduction zones show that fault slip in a temperature range of <100° to 500 °C occurs either by dynamic earthquakes, by a range of slow slip transients, or by steady-state creep with a varying degree of kinematic locking (e.g., Bürgmann, 2018). Interestingly, the dominant deformation mechanisms deduced from fabric types in exposed subduction faults and strain rates predicted from experimentally constrained constitutive laws also indicate three separate deformation regimes in excellent agreement with instrumental observations. Brittle creep controlling localized cataclastic flow at an effective speed range likely determines localized slip in a narrow fault network. Episodic deformation and stress changes during the seismic cycle correspond to intermittent changes between fast brittle creep and slow solution precipitation creep, as revealed by multiple crosscutting relations in exposed cross sections of the plate interface zone.

Accelerated slip is associated with strain localization into a network of localized cataclastic shear zones several meters to tens of meters wide and eventually into principal slip zones of just a few millimeters to centimeters thickness displaying ultracataclasites or pseudotachylite layers. Constitutive laws of dominant mechanisms (Figs. 6 and 7) show that deformation is primarily controlled by fault strength, pore pressures, and stress levels. In particular, toward an advanced stage of the seismic loading cycle at stress levels>70%–80% of peak strength, cataclastic flow is likely to control slip (Fig. 7). Estimating the interface strength and hydraulic properties at the terminal stage of the seismic cycle combining stress drop analysis of megathrust earthquakes and dynamic Coulomb wedge modeling, Wang et al. (2019, and references therein) find values as low as 5–10 MPa at an effective coefficient of friction of µb» 0.03 requiring fluid pressure ratios of>0.95. Recently, Brodsky et al. (2020) suggest even smaller coseismic shear stresses <1 MPa from temperature measurements at the updip edge of the Japan trench. In consequence, the same transition to brittle creep would also result from fluid-pressure pulses to quasi-lithostatic pore pressure. This has been suggested to occur from pulsed prograde dehydration reactions or boosts of updip fluid percolation (Brantut et al., 2011; Fagereng and Diener, 2011; Angiboust et al., 2014b; Muñoz-Montecinos et al., 2021a). Slow slip events preceding large subduction earthquakes (Bürgmann, 2018) may also result from rapid fluid pressure changes (see Fig. 7A). As a consequence, brittle creep pulses are not only an appropriate measure of the stress state but also an indicator of imminent rupture or cyclic instability of loading as well as hydraulic pulses (Fig. 7). In fact, the fault-valve behavior or percolation pulsing as suggested by Sibson (1992), Nakajima and Uchida (2018), Shapiro et al. (2018), and Gosselin et al. (2020) would concur with this conclusion.

Recently available high-resolution observations of seismic slip and changing physical properties suggest a rather complex picture of subduction zones with depth and along strike. This raises the question to what extent continuous creep of a plate interface—exhibiting no or a very low degree of locking—is expected to occur at all. Continuous creep might be controlled by homogeneous and constant near-lithostatic pore pressure throughout, continuously releasing stresses by brittle creep at plate convergence rate. Conversely, continuous creep at near-convergence rates through solution precipitation creep is a potential alternative mechanism in the case of an extremely fine-grained fault gouge of several tens of meters thickness at rather low pore pressures. However, this is less likely, given instrumental and fabric observations of the plate interface region.

The resulting kinematic complexity is also in agreement with the observed behavior of transients constituting the slow earthquake family. Stress states assumed for the fault network fluctuate around values close to failure as indicated by propagation of slip over several tens of kilometers, succession by secondary slip fronts, back-propagation, coalescence of diverse ETS nucleation sites to a large slip zone, periodic to irregular recurrence of ETS, of tremor streaking in the direction of slip (Ghosh et al., 2010), temporal and spatial association of silent slip and tremor as well as low and very low frequency earthquakes at the tips of slow slip zones, and more (e.g., Peacock et al., 2011; Bürgmann, 2018; Hall et al., 2019; Warren-Smith et al., 2019, and references therein). Likewise, the slip rates may also vary significantly within the fault network.

Brittle creep likely controls cataclastic fabric formation reflecting the kinematic and physical properties of slow slip. Extrapolation of existing constitutive relations predicts a narrow bound of strain rates and slip velocities as instrumentally observed for the family of slow slip features. However, it does not explain the acceleration toward slip velocities of jointly appearing regular or very low frequency earthquakes. Recent microphysical models that capture the frictional and healing behavior of rate and state friction (den Hartog and Spiers, 2014; Chen and Spiers, 2016) and the velocity-dependence of frictional strength (Aharonov and Scholz, 2018; Im et al., 2020) may yield a potential answer. The models suggest that the competition of fracture and pressure-solution–controlled dilatancy changes may take a key role. Alternatively, recent experimental studies provide clear evidence for decreasing nucleation patch size and critical slip displacement and a transition from slow slip to dynamic rupturing with increasing loading rate (e.g., Guérin-Marthe et al., 2019). This is in line with some of the fabric observations reported here; however, possible mechanisms controlling slip acceleration are still a matter of debate.

6.2 Brittle Creep and Cataclasites

Cataclastic flow by brittle creep is localized in shear zones along the entire extent of the plate interface zone up to temperatures exceeding 500 °C (Figs. 4J and 8). Progressive grain comminution and localization produce ultracataclasites and associated pseudotachylites that may represent either extensive comminution or frictional melts, respectively. Although microstructural evidence does not provide unique information on slip rates, different cataclasite types may be related to different strain rates or slip velocity regimes. Thin ultracataclasite bands partly showing pseudotachylites and lining curvy-planar principal slip zones may represent products of dynamic ruptures and high slip rates encountered during earthquakes and in agreement with experimental observations (e.g., Scuderi et al., 2017). This is supported by direct observations available from earthquakes triggered by mining activity, where the epicenter and resultant earthquake fault rock are directly accessible after the event. In the case of the TauTona Mine in South Africa, Heesakkers et al. (2011) find a network of sharply bounded faults of several mm to cm width cutting through a quartzite sequence with principal slip surfaces and fine-grained sintered ultracataclasite and/or rock flour lining the rupture zones. Coexistence and mixing of sharply bounded ultracataclasites and pseudotachylites also support their formation in dynamic ruptures as has been suggested for parts of the subduction fault network exposed on Kodiak Island (Rowe et al., 2005, 2011; see also Rowe and Griffith, 2015, for an exhaustive review of earthquake-related fault rock).

Linking foliated and clast-rich cataclasites to instrumentally observed shear strain rates and subseismic slip rates is non-unique to date. Experimental and field observations are limited and suggest that the same fault and fault rock may exhibit a range of slip velocities (Kaproth and Marone, 2013; Leeman et al., 2016; Scuderi et al., 2017; Ikari et al., 2013; Ikari and Kopf, 2017; Ikari, 2019). Here, multiple activation of slip bands over an extended period, multiple overprints during inter-event periods, and subsequent exhumation invariably blur the image. The above experimental work suggests that both, sharply bounded principal slip zones with nanoscale ultracataclasite and more diffuse shear bands, may exhibit a velocity spectrum reflecting silent-slip events (cf. Kirkpatrick et al., 2021, and references therein). A potential resolution of this dilemma is as follows: The narrow principal slip zones lined with ultracataclasite nested in wider foliated cataclasite zones are found both at the very shallow parts of a plate interface zone near the trench (Fabbri et al., 2020) as well as at temperatures ≥400 °C (Angiboust et al., 2015; Fig. 8). While the upper domain may experience seismic slip at earthquake speed from trench-breaking faults—such as in the case of the Tohoku earthquake, the drilling of which also exhibited the above cataclasite spectrum (e.g., Chester et al., 2013; Brodsky et al., 2020)—the domain below the seismogenic zone shows various types of silent slip only. There, the range of slip velocities and strain rates between very low frequency earthquakes and slow slip events may explain two aspects: (1) the cataclasite types at this depth domain are probably the preserved micro-structural record with the locally fastest slip rate potential. Hence, we hypothesize that they are the most likely candidate to reflect the equally fastest instrumentally recorded slip at these depths indicating that brittle deformation here indeed occurs at subseismic creep velocities. (2) The (V)LFEs with their narrow slip zones (Table 2 and Fig. 5A) and the slow slip events are spatially and temporally co-located in the same rock volume. Both possibly reflect the different cataclasite types observed here. A similar argument, linking subseismic slip velocities to foliated cataclasites with a higher clast-matrix ratio, may also hold for the seismogenic zone proper, where long-term silent slip and/or afterslip bursts may occur in addition to coseismic slip velocity. This is supported by observations of a wide range of seismic coupling coefficients at subduction zones (Scholz and Campos, 2012).

Often-observed diffuse and generally foliated cataclasite zones with high clast-matrix ratios possibly include the entire transition between earthquakes, slow earthquakes, and transient slow slip events (e.g., Rowe et al., 2011; Fabbri et al., 2020). Experimental evidence clearly shows that the formation of pressure solution seams and/or mineral indentation and grain fracturing may occur coevally in the same experiment, both in experiments designed to study pressure solution as well as in brittle shear experiments (Niemeijer et al., 2002, 2012; Brantut et al., 2012). The constitutive laws suggest this would potentially be possible in very fine grained, matrix-dominated rocks at lowered pore pressures and high porosity. However, porosities suggested from Vp/Vs analysis for subduction fault rocks range between 2.7%–4% at>30 km depth (Peacock et al., 2011; see also Rowe et al., 2011, finding similar values in exposed cataclasites). These numbers are nearly an order of magnitude lower than those reported for the high-porosity sandstone experiments by Niemeijer et al. (2002). Hence, the above conditions for coeval action at similar rates are more likely to occur in the postseismic period related to afterslip (and pulses thereof) and fault healing, when enhanced fracture-induced permeability lowers pore pressure. Similarly, because of its lower stress dependence, Gratier et al. (2009, 2014) suggest that solution precipitation creep will likely only dominate in the early parts of the loading cycle. Conversely, brittle creep will dominate the stages toward the end of the seismic cycle as increasing pore pressure to near lithostatic values and increasing stress to near the short-term failure strength will replace SPC (Fig. 7). However, a small contribution of solution precipitation creep and dislocation creep is also possible as indicated in Figure 7A. At higher temperatures, both mechanisms add to the summed higher rates even toward the late stage of the loading cycle with brittle creep dominating strain accumulation. Whether this implication from the analysis of the flow laws is securely identifiable in rock fabrics of fossilized systems is questionable. Coeval action as well as repeated sequential predominance of mechanisms with a long series of loading cycles is unlikely to be resolvable and has yet to be demonstrated experimentally. But the ubiquitous presence of pressure-solution features overprinting fracture fabrics, and the strong increase in this overprinting relationship toward depth—i.e., toward higher temperatures—is a predicted outcome from the relationships shown.

This argument also surmises that SPC does not have the potential to generate strain-rate bursts at the SSE velocity required, even at high temperatures and very small grain sizes. These bursts may be due to fluctuation of pore pressure at near-lithostatic values. Such cyclic changes of the hydraulic conditions have been suggested by Miller and Nur (2000) using numerical modeling and recently by Phillips et al. (2020) from observation of both tensile and shear veins occurring in the plate interface environment as well as from observation of variable vein orientations parallel and oblique and/or across the main cataclastic fabric and slip zones (see also Angiboust et al., 2015). Such switches in stress orientation have, for example, been shown from fabric analysis and modeling of the plate interface stress states of the Costa Rica margin (Vannucchi and Leoni, 2007) and the recent Tohoku-Oki earthquake by Brodsky et al. (2020) and Wang et al. (2019). Also, cyclic fluid conditions on the megathrust are seen from variations in vein geochemistry (Cerchiari et al., 2020) as well as from fault sealing through veins with fibrous growth and incremental crack-seal features (Fagereng et al., 2010; Vannucchi et al., 2010; Ujiie and Kimura, 2014; Dielforder et al., 2015; Fisher et al., 2019a; Muñoz-Montecinos et al., 2020). Altogether, these observations underscore the dominant role of the state of the hydraulic system on the switch of mechanisms and the acceleration bursts occurring in the cataclasite network. Adding the observations of fabric distribution with depth and of kinematic regimes, the above observations warrant an extension and adaptation of Sibson's (1983) classical fault scheme (see also Lusk and Platt, 2020) developed for the continental crust toward a subduction fault scheme showing the competition of brittle fracture and solution precipitation creep to great depth that is tuned by the very low geotherm of subduction systems and the fluctuating hydraulic state (Fig. 9).

The above arguments linking cataclasite types to slip velocities mostly hinge on indirect reasoning possibly except for the relationship between earthquake slip velocities and pseudotachylites as well as at least some of the thin principal slip zones lined by ultracataclasite (e.g., Rowe and Griffith, 2015). Reducing the remaining uncertainties and achieving a clear quantitative relationship between cataclasite types and subseismic slip velocities calls for additional efforts in closing the observational gaps. Two avenues of research stand out: experimental rheology will need to explore transient slip events, their velocity spectrum, and the related fabric styles and dominant microphysical mechanisms controlling these. Complementing this goal, observation in nature will require the monitoring of faults that exhibit silent slip events and their sampling, via drilling or direct sampling at the Earth's surface, in order to ascertain the related fabric types.

Exploring recent experimental data and instrumentally observed slip velocities across the entire range of slow earthquakes, we assess the maximum shear strain rates that can be assigned to observed slip events and deformation mechanisms and their fabrics. We find that the cataclasite types typically found at drilled active plate boundary faults, as well as at exposed plate interface faults exhumed from depths of more than 40 km and temperatures of some 500 °C, match the strain rates predicted using constitutive laws for brittle creep or stress corrosion creep (Brantut et al., 2012). The extrapolated strain rates from lab tests are in good agreement with shear strain rates deduced for instrumentally observed range of slow slip transients. Several key results stand out from the analysis of potential strain rates and their dependence on the key constraints:

  1. A single deformation mechanism—brittle creep or stress corrosion creep as the physical mechanism underlying cataclastic flow—can accommodate the full spectrum of instrumentally observed deformation rates of slow slip transients, such as SSEs, NVT, afterslip bursts, etc. As a consequence, the variety of instrumentally observed slow slip features likely produces similar rock fabrics.

  2. Brittle creep likely controls deformation at stresses >70% of the effective peak strength and at elevated pore-fluid pressures. The occurrence of brittle creep may hence serve as a qualitative stress meter with respect to the state of the loading cycle.

  3. A range of cataclasite types characterizes the main plate interface fault network, possibly reflecting different slip velocity regimes: seismic slip velocities—that also embrace very low and low-frequency earthquakes—are preferentially linked to sharply confined curvy-planar principal slip surfaces lined by ultracataclasite (and sometimes pseudotachylite); the subseismic shear strain rates of slow slip events are likely reflected by more diffuse and clast-richer cataclasites with well-developed pressure-solution foliation.

  4. Also, the ubiquitous presence of synkinematic, as well as variably oriented tensile veins near and within the main fault network, indicates cataclastic flow or brittle creep to occur at near-lithostatic fluid pressures and varying stress orientations. This suggests a cyclic evolution of the hydraulic state of the plate interface controlling kinematics. This cycle is also responsible for the apparent fabric-stress paradox from the coeval occurrence of solution precipitation creep and brittle creep alongside hydrofracturing and pseudotachylites to depths well below the “brittle-ductile” transition.

  5. Solution precipitation creep progressively overprints fabrics toward depth and likely dominates plate interface deformation over the earlier stages of the seismic loading cycle and in cases of lowered pore pressures.

  6. Geodetically observed locking would hence preferentially occur at a combination of early tectonic and/or hydraulic loading cycle with temperatures below some 400 °C (also depending on average grain size).

  7. Conversely, constant creep of the plate interface zone is facilitated also at lower temperatures when constant near-lithostatic fluid pressures allow combined brittle creep and pressure-solution creep to entirely accommodate plate convergence.

We gratefully acknowledge the very valuable reviews by an anonymous reviewer and Paola Vannucchi, both of which greatly helped to clarify the manuscript. S.A. acknowledges the A. von Humboldt Foundation for supporting the early stages of this research. This study contributes to the IdEx Université de Paris ANR-18-IDEX-0001.

Science Editor: Shanaka de Silva
Guest Associate Editor: Gray E. Bebout
Gold Open Access: This paper is published under the terms of the CC-BY license.