Velocity-strengthening friction significantly affects interfacial dynamics, strength and dissipation.

Frictional interfaces abound in natural and man-made systems, yet their dynamics are not well-understood. Recent extensive experimental data have revealed that velocity-strengthening friction, where the steady-state frictional resistance increases with sliding velocity over some range, is a generic feature of such interfaces. This physical behavior has very recently been linked to slow stick-slip motion. Here we elucidate the importance of velocity-strengthening friction by theoretically studying three variants of a realistic friction model, all featuring identical logarithmic velocity-weakening friction at small sliding velocities, but differ in their higher velocity behaviors. By quantifying energy partition (e.g. radiation and dissipation), the selection of interfacial rupture fronts and rupture arrest, we show that the presence or absence of strengthening significantly affects the global interfacial resistance and the energy release during frictional instabilities. Furthermore, we show that different forms of strengthening may result in events of similar magnitude, yet with dramatically different dissipation and radiation rates. This happens because the events are mediated by rupture fronts with vastly different propagation velocities, where stronger velocity-strengthening friction promotes slower rupture. These theoretical results may have significant implications on our understanding of frictional dynamics.

Frictional interfaces are abundant in biological (e.g. adherent cells and cell locomotion), engineering (e.g. microelectro-mechanical devices) and geophysical (e.g. earthquake faults) systems around us, and are of fundamental and practical importance. Consequently, understanding the dynamics of dry frictional interfaces has been the focus of intense scientific activity in the last few decades123456. It has been established that under steady-state sliding conditions, the frictional resistance features a non-trivial velocity dependence, and that this dependence has dramatic consequences on the dynamic response of frictional interfaces789101112. Specifically, it has been shown that for a broad range of materials friction is velocity-weakening – that is, the steady frictional resistance is a decreasing function of the sliding velocity – at least in the regime of low velocities, up to a few hundreds of microns per second. This feature favors various instabilities and stick-slip motion5131415.

A very recent compilation of a large set of experimental data for a broad range of materials, however, has revealed that for higher slip velocities, friction generically becomes velocity-strengthening over some range of slip velocities16. The existence of velocity-strengthening behavior might have significant effects on various aspects of frictional dynamics. In particular, recent laboratory experiments on fault-zone materials have documented slow slip interfacial events – an intensely debated issue – and have linked it to a crossover in the frictional response, from velocity-weakening to velocity-strengthening friction, with increasing slip velocity12. While the possible implications of the existence of velocity-strengthening friction have been rather sporadically discussed in the literature111217181920212223242526, to the best of our knowledge a comprehensive and systematic theoretical exploration of these important issues is currently missing.

As a step in closing this gap, we study here the effect of velocity-strengthening friction on spatiotemporal interfacial dynamics, energy dissipation and radiation, and the global interfacial strength, with a special focus on the nucleation, propagation and arrest of rupture fronts. We explore three variants of a realistic rate-and-state friction law, one which is purely velocity-weakening, one which crosses over at higher velocities to logarithmic velocity-strengthening friction, and one which crosses over to linear velocity-strengthening friction.

We show that the presence or absence of velocity-strengthening friction at relatively high slip velocities can significantly affect the global interfacial resistance (strength) and the energy released during frictional instabilities (“event magnitude”), even under quasi-static loading conditions. Different forms of velocity-strengthening friction, in our case logarithmic and linear, give rise to events of similar magnitude, yet with dramatically different dissipation and radiation rates. The difference stems from the broad range of the underlying rupture propagation velocities, where stronger velocity-strengthening friction promotes slower rupture, possibly orders of magnitude slower than elastic wave-speeds. This result is related to the recent experimental observations of Ref. 12. All in all, our results show that velocity-strengthening friction should be properly quantified and incorporated into friction theory as it appears to affect many basic properties of spatially extended frictional interfaces.

We start by discussing theoretical issues relevant to what follows, influenced by the works in Refs. 25,26,27 and then further developed in Refs. 16,20,21,23. The rate-and-state friction model we study has been introduced recently in Refs. 21,23, and is reviewed here briefly. Consider a multi-contact interface and write the ratio A of the real contact area (the area of all contact asperities) to the nominal one, in terms of a state parameter ϕ(x, t) (of time dimensions) as where σ is the normal (compressive) stress at the interface, σH is the material hardness, b is a dimensionless material parameter of order 10−2, and ϕ* is a short time cutoff16282930. ϕ is usually interpreted as the interface's effective age, and its evolution is given by where v is the local interfacial slip velocity and D is a lengthscale related to the contact asperities geometry. , with an extremely small v0 = 1 nm/s, is a regularization function that plays no important role, and is actually omitted in all of the analytic results that follow. As we focus here on unidirectional motion, we do not distinguish between v and |v|.

The frictional stress τ(x, t) is written as a sum of an elastic contribution, τ, and a rheological contribution, τ (corresponding to irreversible deformation of the bulk of the contact asperities), The rheological contribution takes the form τ = A(ϕ) w(v), where at least at low velocities, the rheological part w(v) corresponds to a stress-biased thermally-activated process111630 Here, k is Boltzmann's constant, T is the absolute temperature, Ω is an activation volume and v* is a very small velocity scale. A higher velocity variant of Eq. (4) will be discussed below.

The elastic stress follows the evolution equation where G0 is the interfacial shear modulus and h is the effective width of the interfacial region. G0 may differ from the shear modulus of the bulk (to be denoted below by G), e.g. due to the presence of free surfaces between contact asperities, near interface damage zone and gouge formation. The latter define a lengthscale h — the effective width of the interfacial region — over which G0 applies. The quantity that actually appears in Eq. (5) is the ratio G0/h.

Equations (1) – (5) describe the first variant of the friction model we study below. We begin by describing its behavior under steady sliding at a velocity v. The steady solution of Eq. (2) is ϕ(v) = D/v, from which it follows that the contact area is a logarithmically decreasing function of v13132. The fixed point of Eq. (5) reads , and hence the overall frictional resistance is given by where a higher order logarithmic term was omitted and the following definitions were used In the low velocity regime, i.e. , f is a logarithmic function of v, with ∂f/∂ log v ≈ α − β. Therefore, if α < β, friction is logarithmic velocity-weakening.

Physically, friction is velocity-weakening because the real contact area is a decreasing function of the sliding velocity, and its velocity dependence is stronger than the rheological dependence of τ. However, as discussed at length in Refs. 16,21,23, when the contact area saturates, and friction becomes logarithmically velocity-strengthening. We term this model the logarithmic velocity-strengthening (LS) friction model. The resulting steady-state friction curve is shown in Fig. 1. In case the contact area continues to decrease indefinitely with increasing v, friction remains velocity-weakening for arbitrarily high velocities. This is formally achieved by removing the “1” in the argument of the logarithm in Eq. (1), that is, replacing Eq. (1) by Although this is somewhat unphysical, this choice was widely used in the literature35113334, and we term it the pure velocity-weakening (PW) friction model. The resulting steady-state friction curve is also shown in Fig. 1.

Figure 1

The steady sliding friction coefficient τ/σ vs. the slip velocity for the three model variants (different colors, also marked by labels).

Note that all of the curves coincide at low velocities and that the driving velocity v is marked.

A third variant of the model is obtained by modifying the rheological function w(v), cf. Eq. (4). As discussed extensively in Ref. 16, and to some extent in Ref. 6, the simple picture of a single barrier, linearly biased, thermally-activated process is expected to break down when asperity-level stresses become sufficiently large. When this happens, a different dissipation mechanism is expected to dominate friction. While at the moment there is no general quantitative theory for the velocity dependence of friction in this regime, it is not expected to be logarithmic, but rather to exhibit a significantly stronger dependence on the slip velocity. In this work, we consider a simple model in which the logarithmic dependence crosses over continuously (but not smoothly) to a linear viscous rheology, directly motivated by experiments and theoretical considerations16. Explicitly, we replace Eq. (4) by where m is a dimensionless parameter. We term this model the stronger-than-logarithmic (STL) velocity-strengthening model. The resulting steady-state friction curve is shown in Fig. 1.

We stress that all three variants coincide in the low velocity regime, where they feature logarithmic velocity-weakening friction. At higher slip velocities, the LS variant, which is described by Eqs. (1) – (5), features a crossover to logarithmic velocity-strengthening friction. The PW variant does not feature any strengthening at all (i.e. it remains velocity-weakening), and is obtained from the LS model by using Eq. (8) instead of (1). The STL variant features linear velocity-strengthening friction, and is obtained from the LS model by using Eq. (9) instead of Eq. (4).

In order to investigate the implications of the different constitutive laws on frictional dynamics, we need to consider a spatially-extended interface under inhomogeneous sliding conditions. To this end, we consider a long elastic block of height H (in the y-direction) and length (in the x-direction), in frictional contact (at y = 0) with a rigid substrate (i.e. no deformation of the substrate is considered), see Fig. 2. The trailing edge of the elastic block (at x = 0) is moved at a constant velocity v in the positive x-direction, while the leading edge (at x = L) is stress-free. The block is driven quasi-statically with v = 10 µm/s, which is representative of typical laboratory experiments3536 and generically belongs to the steady-state velocity-weakening friction branch (cf. Fig. 1). The upper edge of the elastic block (at y = H) experiences a constant normal stress σ, σ(x, y = H, t) = σ, but no shear stress, i.e. σ(x, y = H, t) = 0.

Figure 2

A sketch of the spatially-extended frictional system.

An elastic block, which is in frictional contact with a rigid substrate, is loaded by a space- and time-independent normal stress σ(x, y = H, t) = σ (H is the block's height) and driven by a velocity v at its trailing edge (x = 0). The leading edge is at x = L. The shear stress at the interface, σ(x, y = 0, t), is equal to the frictional stress τ(x, t).

We focus on plane-strain deformation conditions and furthermore assume that H is smaller than the smallest lengthscale ℓ characterizing the spatial variation of various fields in the x-direction. Under the stated conditions, the momentum balance equation reduces to (see Ref. 23 for a detailed derivation)

Here ρ is the mass density, u are the components of the displacement vector and σ of Cauchy's stress tensor. Note that the plane-strain Hooke's law was used. In addition, (where G is the shear modulus of the bulk and ν is Poisson's ratio) and the shear stress at y = 0 simply equals the frictional stress, σ(x, y = 0, t) = τ(x, t). Note also that v(x, t) = ∂(x, t). Corrections to Eqs. (11) – (12) appear only to order (H/ℓ)2, a situation reminiscent of the shallow water approximation in fluid mechanics.

Finally, note that the lateral force required to maintain the velocity boundary condition at the trailing edge, u(x = 0, t) = v, reads and the traction-free boundary condition at the leading edge implies ∂(x = L, t) = 0. The geometry of the sliding body () and the sideways loading (tangential driving forces are localized at x = 0) may be more relevant to some systems (e.g. the laboratory experiments of Refs. 35,42) and less so to others (e.g. earthquake faults under tectonic loading). Yet, we believe that our analysis elucidates basic aspects of frictional dynamics and may be useful for understanding rupture dynamics in regions of macroscopic stress concentration even in geophysical contexts.

Equation (11), with the stated boundary conditions and with τ(x, t) corresponding to one of the three friction laws described above, has been solved numerically (see Methods). The model parameters for poly(methyl-methacrylate) (abbreviated PMMA), a polymeric glass that is widely used in laboratory experiments35363738, were extracted from a large set of experimental data (see Methods). The initial conditions are u(x, t) = 0, v(x, t) = 0, τ(x, t) = τ(x, t) = 0, and ϕ(x, t) = 1 s, the latter is typical of laboratory scale experiments. The results presented here are largely insensitive to the choice of the initial value of ϕ.

Results

Global frictional resistance

We begin by studying the macroscopic response of the system. Figure 3 shows the total frictional force exerted by the loading machine as a function of time, f(t). It is seen that the friction force increases gradually until it experiences an abrupt drop, followed by repeated cycles of gradual increases and abrupt drops, typical of frictional systems353940. The drops in the friction force, which appear as vertical lines in this figure, occur when sliding becomes unstable, and involve nucleation and propagation of rupture fronts, as will be discussed below.

Figure 3

The loading force f vs. time for the three models (color code as in Fig. 1).

It is seen that all of the curves coincide for short times, and then begin to diverge. The LS and STL models maintain the same “envelope”, while the PW model features more pronounced stress drops, larger inter-event times and a lower overall resistance. (inset) The same data as in the main panel, but this time is plotted vs. time. The red lines are linear fits to the values of at the rupture arrest times t (i.e. right after the force drops), cf. the prediction in Eq. (17).

Before the first drop, the friction force corresponding to the three variants is identical, as can be expected because the dynamics in this regime are slow and governed by the loading velocity v. In this range of velocities, the three variants coincide and consequently the first drop occurs almost exactly at the same point in time for all of the variants, suggesting that the instability mechanism is insensitive to the high velocity behavior (as predicted in Ref. 23). However, since the instabilities are accompanied by much larger velocities, the high velocity behavior of the friction law becomes important.

Figure 3 demonstrates that while the LS and STL models give rise to almost identical force profiles, the PW model results in significantly larger force drops, and a lower overall interfacial resistance. This suggests, and will be further substantiated in what follows, that while the total energy dissipated during these drops is similar in the LS and STL models, the energy dissipated in the PW model is significantly larger. Other features of the global friction curves shown in Fig. 3, such as the lower envelope of f(t) (corresponding to the values of f(t) after each drop), will be discussed and explained theoretically below.

Spatiotemporal interfacial dynamics

In order to understand the origin of these differences, one must examine the complex spatiotemporal dynamics that give rise to the “force drop events”, which are described at length in Ref. 23. As stated above, the instabilities result in the nucleation, propagation and arrest of rupture fronts, a scenario reported by many experimental, numerical and analytical works274142434445. Most of the remainder of this paper will be focused on the first event, which is marked in Fig. 3 by t. The rationale for focussing on the first event (rather than some later event) is that it ensures that the state of the interface is the same for all three model variants at the onset of instability (with no history effects), cleanly isolating the effects of the existence and form of the velocity-strengthening branch. Having said that, we note that it is clear from Fig. 3 that the differences between the three variants persist to any event. Furthermore, multiple-event properties will be explicitly discussed in relation to Eqs. (15) – (17) and the inset of Fig. 3.

Figure 4 shows the propagation and arrest of rupture fronts during the first event. First, we note the vast difference in the timescales involved: while rupture fronts in the LS and PW models arrest after a few 10 µs, in the STL model they last for a few ms. It is observed, however, that while the penetration depth of the front into the interface in the LS and STL models is comparable, for the PW model it is an order of magnitude larger. Furthermore, the rupture propagation velocity in the LS model is an order of magnitude smaller than in the PW model (the latter is of the order of the elastic wave-speed), and the propagation velocity in the STL model is yet two orders of magnitude smaller.

Figure 4

Propagation of rupture fronts in the first event for (a) All three models and (b) the LS and STL models.

x is the spatial location of the front tip, cf. Fig. 5. t is the time where the front starts to propagate, cf. Fig. 3. Note also the vast difference in timescales between the panels. The wave speed in this system is .

Both the LS and STL models give rise to rupture fronts that are much slower than the elastic wave-speed. These remarkably low rupture propagation velocities, as low as three orders of magnitude slower than the elastic wave-speed in the STL model, might be related to the important, and rather intensely debated, issue of slow rupture2246474849. Our calculations suggest that the emergence of slow rupture might be directly related to the existence and form of velocity-strengthening friction. This is in accord with recent laboratory experiments on fault-zone materials, which documented slow slip events together with a clear crossover from velocity-weakening to velocity-strengthening friction with increasing slip velocity12.

A lot can be learned from the state of the interface after the rupture front has passed. In Fig. 5 we plot the spatial distribution of the (normalized) friction stress just before the first rupture event and immediately after it for the three variant models. In both of these states, the higher slip rates associated with the rupture fronts are not present (before the event they have not yet been generated and after the event they have died off), and the mechanical state is quasi-static. In line with the previous results, prior to the inception of the first event the stress profiles in the three models essentially coincide. When the fronts propagate and eventually arrest, they leave behind them a residual stress profile, which is much smaller in the PW model compared to the LS and STL models. This residual stress is approximately homogeneous in space and is lower than the stress prior to the event. The elastic energy release during this stress relaxation process is the driving force to frictional dissipation.

Figure 5

The frictional stress 1s prior to the first event (dashed lines) and 1s after (solid lines).

The color code is as in Fig. 1. It is seen that the stress left at the tail of the rupture fronts, τ, is roughly homogeneous in space, and that it is much lower in the PW model than in the LS and STL models. The location of the fronts after the event is marked by x. The deeper penetration of the PW model, also shown in Fig. 4, is clearly visible.

The approximate spatial homogeneity of τ left behind any rupture front when it arrests, allows us to estimate the loading force at the discrete arrest times t (that is, there is t corresponding to each rupture event). For that aim, we neglect the contribution to the integral in the region x > x(t), where x(t) is the location of the peak of τ slightly after a rupture front arrested (cf. Fig. 5), and then assume that τ(x, t) can be replaced by a constant residual stress τ, obtaining To calculate x(t), we note that at the arrest times t Eq. (11) takes the form , valid in the range 0 < x < x(t) and neglecting inertia. With the approximate boundary conditions , this equation can be readily solved as This can be substituted in Eq. (15) to give where u(x = 0, t) = v was used (which is, of course, valid at any time, not only at the discrete arrest times t = t).

The prediction in Eq. (17), i.e. f(t)2 ~ t, is tested in the inset of Fig. 3 for all three models over many events (i.e. this is a multiple-event property, not only a property of the first event, which was the focus of the discussion up to now). The analytic prediction is observed to be in favorable agreement with the numerical data for all three models, where the prefactor (slope) in the relation f(t)2 ~ t is the same for the LS and STL models, but is significantly smaller for the PW model. These results show that τ is the same for every rupture event and lend direct support to the assumption that spatial variations of the residual stress left behind any rupture front can be neglected, consistent with the explicit stress profiles shown in Fig. 5 (for the first event in the three different models).

The latter observation allows us to extract τ, the only unknown quantity in Eq. (17) (all other quantities are known parameters, which are the same for all three models), yielding for the LS and STL models and for the PW model. The fact that the models that feature a nonmonotonic velocity dependence, i.e. the LS and STL models, give rise to an essentially identical residual stress τ is intimately related to the value of the steady state stress at the minimum of the friction curve (cf. Fig. 1), which is the same for both. Equation (17) then shows that the fact that the PW model produces a lower overall frictional resistance (and deeper force drops) compared to the LS and STL models is intimately related to the fact that the residual stress left behind the rupture fronts in the PW model is significantly lower than that of the LS and STL models. Furthermore, Eq. (16) suggests an explanation for why the penetration depth, i.e. x(t), is significantly larger in the PW model than in the other two models.

The “static friction coefficient” µ is ordinarily defined as the tangential force, normalized by the normal force, needed to initiate global motion of the block. This force also corresponds to the peak of the loading curve. From this perspective, all of the spatiotemporal dynamics discussed up to now are precursory274243, as they precede global motion which sets in only when a rupture front reaches the leading edge of the block (i.e. when x = L). Hence, we can estimate µ, which quantifies the global frictional resistance, as where Eq. (15) was used. This shows that the “static” frictional resistance of the interface, measured at slow loading velocities (here v = 10 µm/s), is influenced by dynamic processes at much higher slip rates and furthermore that the existence of velocity-strengthening friction behavior strongly affects µ through τ405051.

The results discussed above highlight two important points. First, an effectively constant residual stress τ is left behind rupture fronts in all of the models studied here. This property emerges spontaneously, unlike conventional slip-weakening models in which it is assumed a priori (see, for example Refs. 41,44,52, and the discussion in Ref. 53). The value of τ depends on the existence of velocity-strengthening friction, which in turn has significant implications on the strength of the interface, as evident from Fig. 3 and Eq. (18). Note also that the constancy of the residual stress τ implies that the mechanical fields associated with frictional shear cracks in 2D are well described by the classical theory of fracture54. Second, once τ is known, the arrest of rupture fronts is determined by global equilibrium conditions, rather than by dynamic considerations (cf. Eq. (15)).

Energy partition: Dissipation and radiation

As energy dissipation is at the heart of frictional phenomena, it will be interesting and instructive to consider the energy budget in the system. As a starting point, we briefly remind the reader that the linear momentum conservation law of Eq. (10) can be transformed into a continuity equation for the energy density (using Hooke's law and integration by parts). The result reads The first term is the rate of variation of the energy density (both kinetic and elastic), and the second term is the divergence of the energy flux vector. Their sum vanishes when energy is conserved.

Following the same procedure, one can derive the energy continuity equation for our model by combining Eqs. (5) and (11), obtaining where we defined Here ε is the kinetic energy density, ε is the (bulk) linear elastic strain energy density, ε is the interfacial elastic energy density and J is the energy flux. The interfacial energy density, ε, is dissipated during sliding due to the rupture of asperities, resulting in a dissipation rate p, in addition to the standard dissipation rate p = τ (p, defined in Eq. (21), in fact contains also a contribution of the form . As this contribution is negligibly small in our calculations, we omitted it).

Equation (20) has the same structure as Eq. (19), except for the non-vanishing dissipation rate p, which exists because frictional dynamics are dissipative, and the existence of an interfacial elastic contribution ε (both in the stored energy and in the dissipation power p).

The quantities defined in Eq. (21) are densities that exhibit complex spatiotemporal behaviors during frictional instabilities (which result in rupture events). In order to gain some insight into these complex energy-exchange processes, it will be useful to consider the corresponding space-integrated quantities and .

The interplay between these various quantities during frictional instabilities (“events”), shown for all three models in Fig. 6, is an essential feature of interfacial dynamics. Our goal is to quantify generic energy-exchange processes during frictional instabilities55 and in particular to understand the differences between the three models in this respect. As the dynamics during frictional instabilities are much faster than typical loading rates, we expect them to be exclusively driven by the already stored elastic energy. That is, we expect the rate of change of the sum of bulk and interfacial elastic energies, ∂(E + E), to be negative during an event. Figure 6 clearly demonstrates this, and that ∂ is negligible compared to ∂ (hence we neglect the former compared to the latter in what follows).

Figure 6

The rate of change of energies and dissipation rates (ε (x, t) and p(x, t) are defined in Eqs. (21)) during the first event for the three models.

The time integral of ∂ over the event duration is the total energy released, which is a natural measure of the magnitude of the event (other measures exist as well). The elastic energy released is either being dissipated directly or is being first transformed into kinetic energy (“radiation”). Eventually, the kinetic energy is also dissipated. This generic picture is demonstrated in Fig. 6 for all three models. In particular, it is observed that the dissipation contributions P and P are comparable, where the former is typically larger than the latter. Kinetic energy generation (“radiation”), ∂ > 0, is observed in the first part of the event. In the second part of the event ∂ < 0, when the kinetic energy decays and is being dissipated.

While this generic qualitative picture is similar in all three models, there are large quantitative differences that we wish to discuss now. The main characteristics of the first rupture event in the LS, STL and PW models are summarized in Table 1. As we already know from Fig. 4, the events are mediated by rupture fronts of vastly different velocities in the three models (~ 103 m/s in the PW model, ~ 102 m/s in the LS model and ~1 m/s in the STL model). The event duration is about 40% larger in the PW model as compared to the LS model, both in the few 10 µs range, while it is two orders of magnitude larger in the STL model (~ms). Despite the vast differences in the rupture propagation velocity and event duration, the total dissipated energy (which equals the amount of elastic energy released during the event) in the LS and STL models is essentially identical. This is in line with Fig. 3, which shows that the two models feature nearly identical stress drops and frictional resistance, and with Fig. 5 and the inset of Fig. 3, which show that the residual stress τ in the two models is essentially identical. This result clearly demonstrates that depending on the form of the velocity-strengthening friction branch (e.g. logarithmic vs. linear) one can observe events of the same magnitude (i.e. integrated dissipation/energy release) accompanied by very different dissipation rates (see Table 1). This feature might be related to geophysical observations indicating that slow rupture does not necessarily imply smaller integrated slip and energy release [ref. 46, especially Figure 5]. It is worth noting, though, that in both the LS and STL models the event under consideration is significantly slower than ordinary, wave-speed fast, rupture. This might be related to the fact that we focussed on the first event (for the reasons explained above), eliminating any history dependence and in particular residual stresses associated with previous events. Indeed, a recent study23 demonstrated that a properly aligned pre-stress can increase the velocity of rupture in the LS model to be of the order of the elastic wave-speed, while rupture in the STL model remains much slower (cf. Fig. 4 in Ref. 23).

Table 1

Summary of the main characteristics of the first rupture event in the LS, STL and PW models

	PW	LS	STL
Velocity strengthening	Absent	Logarithmic	Linear
Rupture propagation speed	1540 m/s	166 m/s	~ 3 m/s
Event's duration Δt	~ 50 µs	~ 35 µs	~ 2000 µs
Total dissipated energy	1.4 J/m	0.26 J/m	0.26 J/m
Maximal dissipation rate	52 kW/m	16 kW/m	0.19 kW/m
The maximum of Ek	0.27 J/m	4.5 mJ/m	0.54 µJ/m
Penetration length xtip(ta)	11.3 cm	5.52 cm	5.46 cm

The total dissipation in the PW model is about 5.4 times larger than the total dissipation in the LS and STL models, consistent with the much larger stress drops and the significantly reduced interfacial resistance observed in Fig. 3. Moreover, the amount of kinetic energy generated during the event is much larger in the PW model as compared to the other two models, and is about 19% of the total energy released (though eventually it is also dissipated). In systems of larger heights H, this radiated kinetic energy will decay on longer timescales, allowing it to interact with remote boundaries. The kinetic energy generated in the STL model is negligibly small, while in the LS it makes about 1.7% of the released energy (a similar value was reported in Ref. 55, although direct comparison is precarious). All in all, these results provide strong evidence that the existence and form of velocity-strengthening friction has significant implications on frictional dynamics and strength.

Discussion

In this work we studied the spatiotemporal dynamics in three variants of a realistic rate-and-state friction model under quasi-static side-loading conditions, and showed that the existence and form of velocity-strengthening friction may significantly affect various aspects of the frictional response of interfaces. These include the propagation velocity of coherent fronts that mediate interfacial rupture events, the emergence of slow rupture, the elastic energy released during events (i.e. their magnitude), the dissipation and radiation rates, and the global frictional resistance (strength). The clear connection between the existence of velocity-strengthening friction and slow rupture appears to be related to the recent experimental results of Ref. 12. It is also shown that events of similar magnitude (and hence stress drops) can be accompanied by substantially different dissipation and kinetic energy radiation rates.

It is important to note that while our analysis addressed the role of a crossover to velocity-strengthening friction in homogeneous interfaces (homogeneous in terms of the constitutive law, not the stress distribution), one should bear in mind that velocity-strengthening friction may have other important implications. For example, earthquake faults are typically spatially heterogeneous, featuring variation in the constitutive properties as a function of depth. Indeed, models in which velocity-weakening and velocity-strengthening friction segments coexist in spatially different parts of the fault have shown that the latter can play a role in rupture arrest and after-slip, cf. Ref. 56.

Our theoretical results, together with extensive experimental evidence16, highlight the need to quantitatively characterize the velocity-strengthening frictional response of interfaces, both experimentally and theoretically, and to systematically incorporate it into friction theory. Since frictional instabilities spontaneously lead to accelerated slip that probes relatively high-velocity properties of frictional interfaces, the latter – which include velocity-strengthening friction – affect the frictional response even under quasi-static loading conditions. This understanding may offer new ways to interpret existing observations in a broad range of frictional systems and to develop predictive theories of the dynamics of spatially extended frictional interfaces.

Methods

Governing equations

The evolution equations for each of the three models is a set of coupled partial differential equations for the fields ϕ(x, t), τ(x, t) and u(x, t) which are governed by Eqs. (2), (5) and (11), respectively (with τ in Eq. (11) replaced by Eq. (3)). The LS and PW models both use Eq. (4) for w(v), and the only difference is that the LS model uses Eq. (1) for A(ϕ), while the PW model uses Eq. (8). The STL model uses Eq. (1) for A(ϕ), and (9) for w(v).

Numerical integration

The governing equations were numerically integrated using the Method of Lines, by discretizing the spatial derivative in Eq. (11) and then employing a standard adaptive differential solver (NDSolve, Mathematica 9) for integrating the resulting set of ordinary differential equations. The spatial mesh was chosen to be small enough, such that numerical convergence was achieved.

Parameters

All the material parameters for PMMA, except for m and v, were extracted from available data in the literature as described in Ref. 21. They are summarized in Table 2. The parameters m and v, pertinent only for the STL model, were not directly measured for PMMA (for such measurements in other materials, see Fig. 1 in Ref. 16). For the sake of concreteness, we chose the values m = 25 and v = 7.5 mm/s. The results are qualitatively insensitive to this choice.

Table 2

Material parameters for PMMA

	9.3 GPa	σ	1 MPa
G0/h	300 MPa/µm	b	0.075
ρ	1,200 Kg/m3	D	0.5 µm
v*	0.1 µm/s	D/ϕ*	1.5 mm/s
	2.7 MPa	σH	540 MPa

Author Contributions

E.B. and Y.B.S. conceived research. Y.B.S. performed the calculations, wrote the numerical code, extracted the experimental parameters from available literature and generated all figures. Y.B.S., R.S., E.A.B. and E.B. discussed the results, contributed to their analysis and were involved in writing the manuscript.