††thanks<EMAIL_ADDRESS><EMAIL_ADDRESS>
# Velocity-driven frictional sliding: Coarsening and steady-state pulse trains

Thibault Roch1 Efim A. Brener2,3 Jean-François Molinari1 Eran Bouchbinder4
1Civil Engineering Institute, Materials Science and Engineering Institute,
Ecole Polytechnique Fédérale de Lausanne, Station 18, CH-1015 Lausanne,
Switzerland
2Peter Grünberg Institut, Forschungszentrum Jülich, D-52425 Jülich, Germany
3Institute for Energy and Climate Research, Forschungszentrum Jülich, D-52425
Jülich, Germany
4Chemical and Biological Physics Department, Weizmann Institute of Science,
Rehovot 7610001, Israel

###### Abstract

Frictional sliding is an intrinsically complex phenomenon, emerging from the
interplay between driving forces, elasto-frictional instabilities, interfacial
nonlinearity and dissipation, material inertia and bulk geometry. We show that
homogeneous rate-and-state dependent frictional systems, driven at a
prescribed boundary velocity — as opposed to a prescribed stress — in a range
where the frictional interface is rate-weakening, generically host self-
healing slip pulses, a sliding mode not yet fully understood. Such velocity-
driven frictional systems are then shown to exhibit coarsening dynamics
saturated at the system length in the sliding direction, independently of the
system’s height, leading to steadily propagating pulse trains. The latter may
be viewed as a propagating phase-separated state, where slip and stick
characterize the two phases. While pulse trains’ periodicity is coarsening-
limited by the system’s length, the single pulse width, characteristic slip
velocity and propagation speed exhibit rich properties, which are
comprehensively understood using theory and extensive numerics. Finally, we
show that for sufficiently small system heights, pulse trains are accompanied
by periodic elasto-frictional instabilities.

##  I Introduction

Frictional systems are composed of two bodies coupled at a contact interface,
formed by compressive forces that hold them together. Frictional motion is
typically driven by shear forces that are applied far from the frictional
interface. The frictional interface and its spatiotemporal dynamics are
generically characterized by strong nonlinearity and dissipation, where the
interfacial response depends on the local slip velocity $v$ and on the
structural state of the interface, carrying memory of its history [1, 2, 3, 4,
5, 6, 7, 8]. Different parts of the frictional interface are coupled through
long-range spatiotemporal interactions mediated by the bodies in contact. The
latter correspond to the bulk elastodynamics of the bodies, dependent on their
elastic response functions — which in turn depend also on their geometry — and
on material inertia. Consequently, frictional dynamics inherently emerge from
the coupled effects of interfacial and bulk physics [9, 10, 11].

The interplay between external driving forces, interfacial nonlinearity and
dissipation, material inertia and bulk geometry gives rise to very rich
spatiotemporal dynamics, which characterize a wide variety of natural and man-
made frictional systems, ranging from geological earthquake faults to a
multitude of engineering structures and devices [12, 13, 14, 15, 16, 17].
Understanding, predicting and controlling frictional dynamics remain major
scientific and technological challenges. An inseparable aspect of these
challenges is that frictional systems host various spatiotemporal
instabilities; among these, the most well-characterized instability is
associated with rate-weakening friction, i.e. with physical situations in
which the steady frictional resistance is a decreasing function of the slip
velocity $v$ [18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33,
34]. This interfacial destabilizing process is counteracted by stabilizing
bulk elastic interactions, giving rise to a critical elasto-frictional length
for the onset of instability (to be accurately defined below). Such
instabilities typically result in rapid slip propagation along frictional
interfaces, mediated by rupture modes [35, 36, 37, 38, 39, 40, 41, 42, 43, 34,
44, 45, 46].

Spatiotemporal rupture propagation modes can be generally classified into
expanding cracklike rupture fronts and compact self-healing slip pulses [47,
35, 36, 37, 38, 39, 48, 49, 50, 51]. In the former, $v$ at an interfacial
position behind the propagating mode remains finite as long as propagation
persists, while in the latter, $v$ vanishes over a finite time as propagation
persists. The conditions for the emergence of self-healing slip pulses and
their properties are not yet fully understood, though it is currently accepted
that this mode of frictional rupture propagation is prevalent [47, 52, 53, 54,
55, 51]. For example, it is not yet clear how and under what conditions the
very same frictional system can feature both cracklike rupture and slip pulses
[35, 36, 37, 38, 39, 40, 41, 42, 43, 34, 44, 45, 46]. In this paper, we show
that there exist generic, and in fact widely used, external driving forces
that may lead to the generation of slip pulses. In particular, we show that
velocity-driven frictional sliding — as opposed to stress-driven frictional
sliding — may give rise to propagating periodic slip pulse trains, whose
emerging properties are extensively studied below.

Figure 1: (a) A schematic representation of the class of frictional systems
under consideration. Two identical elastic bodies of width $W$ and height $H$
are in contact along an interface located at $y\\!=\\!0$ (brown line). The
bodies are under constant normal (compressive) stress $\sigma_{0}$ and are
driven anti-symmetrically at the upper and lower boundaries with a constant
shear velocity $v_{0}/2$ such that the overall applied slip rate is $v_{0}$.
Periodic boundary conditions are imposed in the $x$ direction. (b) The steady-
state normalized frictional strength $f_{\rm ss}$ (solid brown line) vs. the
steady-state slip velocity $v_{\rm ss}$ (normalized by $v_{\rm min}$, the
minimum of $f_{\rm ss}$) on a semi-logarithmic scale. The curve has a generic
N shape [56], with two rate-strengthening branches ($df_{\rm ss}/dv_{\rm
ss}\\!>\\!0$) separated by a rate-weakening branch ($df_{\rm ss}/dv_{\rm
ss}\\!<\\!0$), see text for additional discussion. The dashed horizontal
orange line represents an imposed driving stress
$f_{0}\\!=\\!\tau_{0}/\sigma_{0}$, which intersects the steady-state friction
curve at three points. Two are stable fixed-points on the rate-strengthening
branches of the friction law (black squares), while the third one is an
unstable fixed-point on the rate-weakening branch (black circle). The dashed-
dotted vertical green line corresponds to an imposed driving velocity $v_{0}$
— as in panel (a) — and intersects the friction law at a single point on a
rate-weakening branch. As is extensively discussed in the text, imposing a
slip velocity $v_{0}$ on the _unstable_ rate-weakening branch cannot result in
spatially-homogeneous sliding; rather, it leads to a dynamic and spatially-
inhomogeneous phase separation between low velocity (‘stick phase’, left-
pointing arrow) and high velocity (‘slip phase’, right-pointing arrow)
regions. The outcome is a propagating phase-separated state in the form of a
pulse train. (c) A generic pulse train observed under velocity-controlled
frictional sliding. The train travels at a velocity $c_{\rm p}$, and features
periodicity $W_{\rm p}$ and pulses of width $w_{\rm p}$. The pulses feature a
maximum slip velocity $v_{\rm m}$ and propagate into the stick phase of a
characteristic stick velocity $v_{\rm s}$, which is vanishingly small. The
dashed-dotted horizontal green line corresponds to the driving velocity
$v_{0}$.

To understand the qualitative differences between velocity-driven and stress-
driven frictional sliding, consider the frictional system illustrated in Fig.
1a, composed of two identical bodies of height $H$ and length $W$, and
characterized by elastic constants $\mu$ (shear modulus) and $\nu$ (Poisson’s
ratio), and mass density $\rho$. The bodies are held together by a homogeneous
normal (compressive) stress $\sigma_{0}$ and are subjected to some shear-
related boundary conditions. In the figure, the boundary conditions at the
upper and lower edges are denoted by $v_{0}/2$, implying that in this case the
applied velocity (of overall magnitude $v_{0}$) is kept fixed. However, one
could also consider a situation in which the applied shear stress $\tau_{0}$
is kept fixed.

To explain why these two types of driving forces may lead to qualitatively
different physical consequences, we plot in Fig. 1b the steady-state
frictional resistance, $f_{\rm ss}\\!=\\!\tau_{\rm ss}/\sigma_{0}$ (where
$\tau_{\rm ss}$ is the steady-state frictional strength/stress), as a function
of the logarithm of the steady-state slip velocity $v_{\rm ss}$. By steady-
state we mean that we focus on a point along the frictional interface that
experiences a slip velocity $v_{\rm ss}$ for a sufficiently long time, and
measure the resulting frictional strength/stress $\tau_{\rm ss}$; the steady-
state friction coefficient, which out of steady-state depends on the
instantaneous slip velocity $v$ and on a set of structural (internal)
variables (see below), is simply given by $f_{\rm ss}\\!=\\!\tau_{\rm
ss}/\sigma_{0}$, as stated above. The steady-state friction curve $f_{\rm
ss}(v_{\rm ss})$ features a generic N shape [56], where friction is rate-
strengthening at both low and high slip velocities (the latter occurs above
the minimum of the curve, denoted by $v_{\rm min}$), and rate-weakening at
intermediate velocities (typically spanning a few orders of magnitude, note
again the logarithmic $v_{\rm ss}$-axis).

To highlight the differences between velocity-driven and stress-driven
sliding, consider the frictional system illustrated in Fig. 1a under the
application of a fixed (total) slip velocity $v_{0}$, which resides in the
rate-weakening branch of the steady-state friction curve (marked by the
dashed-dotted vertical green line in Fig. 1b). Consider also the stress-driven
counterpart of this system, i.e. the case in which a fixed shear stress
$\tau_{0}$ is applied. The latter corresponds to
$f_{0}\\!=\\!\tau_{0}/\sigma_{0}$ (marked by the dashed horizontal orange line
in Fig. 1b) and is chosen such that $f_{0}\\!=\\!f_{\rm ss}(v_{0})$, i.e. the
dashed-dotted vertical and dashed horizontal lines intersect the steady-state
friction curve at the very same point, marked by a black circle in Fig. 1b.
The horizontal line, corresponding to stress-driven sliding, intersects the
steady-state friction curve also at two other points on the two rate-
strengthening branches, marked by black squares. Let us focus first on the
stress-driven case and ask whether spatially homogeneous and stable steady
sliding can emerge under these conditions. It is clear that steady sliding of
the whole system at $v_{0}$ is not possible because the rate-weakening branch
is unstable (for $W$ larger than the critical elasto-frictional length). On
the other hand, homogeneous and stable steady sliding at the velocities
corresponding to the black squares is possible, because the latter reside on
rate-strengthening branches, which are stable [9].

This situation is in sharp and qualitative contrast to the velocity-driven
case, where the sliding velocity $v_{0}$ is enforced on the outer boundaries.
In this case, spatially homogeneous sliding is _impossible_ because $v_{0}$
resides on the rate-weakening branch, implying instability, and no other
spatially homogeneous velocity solutions are possible either (unlike the
stress-driven case). Consequently, either steady sliding does not exist at all
or spatially inhomogeneous steady-state $v(x,t)$ emerges ($x$ here denotes the
spatial coordinate along the frictional interface and $t$ is the time) such
that $W^{-1}\\!\int_{0}^{W}\\!v(x,t)\,dx\\!=\\!v_{0}$ at any time. In the
latter situation, $v(x,t)$ must take the form of a steadily propagating pulse
train, illustrated in Fig. 1c. The pulse train is characterized by a spatial
periodicity $W_{\rm p}$ and a propagation velocity $c_{\rm p}$. Each pulse
within the train features a characteristic width $w_{\rm p}$ (to be accurately
defined below), a maximal velocity $v_{\rm m}\\!>\\!v_{0}$ and a minimal
velocity $v_{\rm s}\\!<\\!v_{0}$. Since $v_{\rm s}$ is typically much smaller
than $v_{0}$, it is termed the ‘stick velocity’, corresponding to a nearly
non-sliding state, which is termed the ‘stick phase’ (cf. Fig. 1b). The parts
of a pulse train that feature $v(x,t)\\!>\\!v_{0}$ can be regarded as ‘slip
phases’ (cf. Fig. 1b); consequently, a pulse train may be viewed as a
propagating mode composed of alternating stick and slip phases.

This physical picture of pulse trains as propagating modes composed of
alternating stick and slip phases, emerging under velocity-driven conditions,
may suggest an analogy to phase separation (and the associated Maxwell
construction) in equilibrium thermodynamics [57]. In the latter, phase
separation emerges when the pressure-volume isotherm of a system in
equilibrium features a non-monotonic behavior, corresponding to a non-convex
free energy [57]. If then the system is enforced to have a volume in the non-
monotonic region of the pressure-volume isotherm, where the pressure is an
increasing function of the volume (thermodynamic stability requires the
pressure to be a decreasing function of the volume), it cannot attain a
(single) homogeneous phase; rather, the system undergoes a phase transition
that leads to the co-existence of two phases of different densities [57].
While $f_{\rm ss}(v_{\rm ss})$ of Fig. 1b is by no means an equilibrium
pressure-volume isotherm, rather it corresponds to a strongly dissipative
interfacial response function of a driven open system, and while the imposed
slip velocity $v_{0}$ inside the rate-weakening (unstable) branch is by no
means the volume of an equilibrium system, there exists a clear and direct
analogy between the two physical situations. In the frictional case, the
result is not static phase separation, but rather a periodic and propagating
phase-separated state of alternating stick and slip phases, the train of
pulses illustrated in Fig. 1c.

It is important to stress that in practical terms velocity-driven frictional
sliding is the rule, rather than the exception. That is, velocity is almost
always what is actually being controlled; in order to maintain a fixed stress,
one needs to employ a feedback loop such that the velocity is precisely varied
so as to keep the stress fixed. This procedure depends on the dynamics of the
system, which may be fast, and is in general difficult to achieve. In most
cases, it is just the velocity which is prescribed.

Our goal in this paper is to understand the spatiotemporal dynamics of
velocity-driven frictional systems and in particular the emergence of steadily
propagating pulse trains. Moreover, we aim at understanding the selection of
the train properties, i.e. its spatial periodicity and propagation velocity,
as well as the properties of a single pulse within the train. These goals are
achieved using extensive numerical simulations — employing several
computational methodologies — and theoretical analysis performed within a
generic rate-and-state friction constitutive framework (to be detailed below).
In particular, we employ the Boundary Integral Method in the
$H\\!\to\\!\infty$ limit; while this method is formulated for stress-driven
sliding in terms of $\tau_{0}$ [58, 59, 60], we show that it can nevertheless
be employed to mimic velocity-driven frictional dynamics. For finite $H$, we
employ the Finite Element Method, as explained below.

We find that frictional systems under velocity-driven conditions, where the
imposed velocity $v_{0}$ resides on an unstable rate-weakening branch of the
steady-state friction curve, feature coarsening dynamics that lead to pulse
trains whose periodicity $W_{\rm p}$ is determined by the system length $W$
(where periodic boundary conditions in the sliding direction are employed, see
below), independently of $H$. We also show that in the small $H$ limit,
coarsening competes with elasto-frictional instabilities, giving rise to pulse
trains with $W_{\rm p}\\!=\\!W$ that experience repeated/periodic
instabilities. Furthermore, we show that the pulse train propagation velocity
$c_{\rm p}$, the average single pulse slip velocity and the single pulse width
$w_{\rm p}$ are related through an equation of motion inspired by fracture
mechanics [61]. The latter involves the pulse leading edge singularity and the
emergence of an effective fracture energy. Finally, we show that the single
pulse width $w_{\rm p}$ reveals non-trivial dependencies on $W_{\rm
p}\\!=\\!W$ and $v_{0}$, featuring properties that are qualitatively different
from those of ideal pulse solutions [62, 63, 64]. Taken together, we provide a
comprehensive physical picture of velocity-driven frictional sliding in
general, and of the properties of the emerging pulse trains in particular.

##  II Simulating velocity-driven frictional dynamics

The rich spatiotemporal dynamics featured by frictional systems, with the
multitude of physical factors at play (as discussed above), make purely
analytical treatments of this class of problems practically impossible.
Consequently, one needs to resort to numerical simulations, at least at the
initial stages of investigating a given set of questions. Independently of the
approach taken — either analytic, computational or hybrid — one should first
adopt bulk and interfacial constitutive relations, which is done in Sect.
II.1. To fully define the problem at hand, one should then specify the bulk
geometry and external driving forces (boundary conditions) — here following
Fig. 1a —, and finally one needs to select a solution method for the coupled
bulk-interface problem. Focusing first on computational methods, we explain in
Sect. II.2 how the Boundary Integral Method — conventionally formulated in
terms of stress boundary conditions — can be used to mimic velocity-driven
frictional dynamics in the $H\\!\to\\!\infty$ limit. Next, in Sect. II.3, we
discuss the usage of the Finite Element Method to address the finite $H$
regime.

###  II.1 Bulk and interfacial constitutive relations: Linear elastodynamics
and rate-and-state friction

The frictional system illustrated in Fig. 1a is formed by two symmetric
bodies, each satisfying its own continuum momentum balance equation
$\rho\,\ddot{\bm{u}}(x,y,t)\\!=\\!\nabla\\!\cdot{\bm{\sigma}}(x,y,t)$, where
$\rho$ is the mass density, $\bm{u}(x,y,t)$ is the displacement vector field,
$\bm{\sigma}(x,y,t)$ is the stress tensor field, $(x,y)$ is a two-dimensional
Cartesian coordinate system and $t$ is the time (a superposed dot represents a
partial time derivative). $\bm{\sigma}$ in each body is related to $\bm{u}$
through a bulk constitutive relation, which is taken here to be that of linear
elasticity, i.e. we adopt Hooke’s law [65] of the form
$(1+\nu)\mu\left[\nabla{\bm{u}}\\!+\\!(\nabla{\bm{u}})^{\mbox{\scriptsize
T}}\right]\\!=\\!{\bm{\sigma}}-\nu({\bm{I}}\operatorname{tr}{\bm{\sigma}}-{\bm{\sigma}})$.
Here ${\bm{I}}$ is the identity tensor, $\nu$ is Poisson’s ratio and $\mu$ is
the shear modulus of each body. Note that body forces are neglected in the
momentum balance equation and that the interface resides at $y\\!=\\!0$ (cf.
Fig. 1a). An interfacial constitutive law is an implicit boundary condition
for the two bulk problems defined above, formulated in terms of a functional
relation between the interfacial shear stress, the interfacial normal stress
and the slip velocity (and typically also additional interfacial state fields,
see below). In the problems considered in this paper, the interfacial normal
stress is constant (uncoupled to frictional sliding), i.e.
$\sigma_{yy}(x,y\\!=\\!0,t)\\!=\\!-\sigma_{0}$, where $\sigma_{0}$ is the
applied compressive stress (cf. Fig. 1a).

The slip velocity $v(x,t)$ is the time derivative of the slip displacement
$\delta(x,t)$, $v(x,t)\\!=\\!\dot{\delta}(x,t)$. Under in-plane shear (the so-
called mode-II) conditions, where
${\bm{u}}(x,y,t)\\!=\\!(u_{x}(x,y,t),u_{y}(x,y,t),0)$ (here the $z$ component
of the displacement vector field, $u_{z}(x,y,t)$, vanishes), one has
$\delta(x,t)\\!\equiv\\!{u}_{x}(x,y\\!\to\\!0^{+},t)-{u}_{x}(x,y\\!\to\\!0^{-},t)$,
where $+/-$ correspond to the upper/lower bodies, respectively. The relevant
interfacial shear stress in this case is $\sigma_{xy}(x,y\\!=\\!0,t)$. Under
anti-plane shear (the so-called mode-III) conditions, where
${\bm{u}}(x,y,t)\\!=\\!(0,0,u_{z}(x,y,t))$ and $z$ is the out-of-plane
direction (perpendicular to the $x\\!-\\!y$ plane), one has
$\delta(x,t)\\!\equiv\\!{u}_{z}(x,y\\!\to\\!0^{+},t)-{u}_{z}(x,y\\!\to\\!0^{-},t)$,
and relevant interfacial shear stress in this case is
$\sigma_{yz}(x,y\\!=\\!0,t)$. The interfacial shear stress, either
$\sigma_{xy}(x,y\\!=\\!0,t)$ (mode-II) or $\sigma_{yz}(x,y\\!=\\!0,t)$ (mode-
III), is continuous across the interface and equals the frictional
stress/strength $\tau(x,t)$. Below we present results for both mode-II and
mode-III, which are qualitatively and even semi-quantitatively similar. Mode-
III (corresponding to scalar elastodynamics) is, however, mathematically and
computationally simpler.

The frictional stress/strength $\tau(x,t)$ is related to $v(x,t)$,
$\sigma_{yy}(x,y\\!=\\!0,t)$ (which in our case equals $-\sigma_{0}$) and
additional interfacial state fields through the interfacial constitutive law.
The latter, at any position $x$ along the interface and at any time $t$, is
described by the following local relation

$\tau(\sigma_{0},v,\phi)=\sigma_{0}\,\operatorname{sgn}(v)\,f(|v|,\phi)\ ,$
(1)

where $\phi(x,t)$ is a non-equilibrium order parameter, sometimes termed an
internal-state field, which represents the structural state of the interface
and encodes its history [1, 2, 3, 4, 5, 6, 7, 8]. Extensive evidence indicates
that $\phi$ physically represents contact’s age/maturity [2, 3, 4, 5, 6, 7, 8]
and that its evolution follows

$\dot{\phi}=g\\!\left(\frac{\phi|v|}{D}\right)\ ,$ (2)

with $g(1)\\!=\\!0$. The characteristic slip displacement $D$ controls the
transition from a stick state $v\\!\approx\\!0$ to a steadily slipping/sliding
state $v_{\rm ss}\\!\neq\\!0$, with $\phi_{\rm ss}\\!=\\!D/v_{\rm ss}$ (the
latter corresponds to $\dot{\phi}\\!=\\!g(1)\\!=\\!0$). Under steady-state
sliding conditions and a controlled normal stress $\sigma_{0}$, the function
$f_{\rm ss}(v_{\rm ss})\\!=\\!f(|v_{\rm ss}|,\phi_{\rm ss}\\!=\\!D/v_{\rm
ss})\\!=\\!\tau_{\rm ss}(v_{\rm ss})/\sigma_{0}$ has been measured over a
broad range of slip rates $v$ for many materials [56]. $f_{\rm ss}(v_{\rm
ss})$ is generically N-shaped, as shown in Fig. 1b, where the precise
functional form of $f(\cdot)$ and $g(\cdot)$ are detailed in Appendix A.

###  II.2 Mimicking velocity-driven frictional dynamics in infinite systems
using the boundary integral method

With the bulk and interfacial constitutive relations at hand, the coupled
bulk-interface problem is fully defined once the bulk geometry is specified.
The system length is taken to be $W$, as shown in Fig. 1a; in order to avoid
lateral edge effects, we employ periodic boundary conditions along this
direction. The remaining geometric length scale in the problem is $H$, i.e.
the height of each body (cf. Fig. 1a). It makes a difference whether $H$ is
taken to be arbitrarily large or finite, as is explained below. Here we first
discuss the $H\\!\to\\!\infty$ case.

In the $H\\!\to\\!\infty$ limit, i.e. when the upper and lower boundaries are
located infinitely away from the frictional interface, and no wave reflections
from these boundaries take place, on can naturally invoke a Green’s function
approach. The latter allows to eliminate the two bulk problems altogether,
reducing the coupled bulk-interface problem to an interfacial integro-
differential equation of the form [58, 19, 59, 60, 66]

$\tau[v(x,t),\phi(x,t)]=\tau_{0}(t)-\frac{\mu}{2c_{s}}[v(x,t)-v_{0}]+s(x,t)\
.$ (3)

where $\phi(x,t)$ satisfies Eq. (2) and the right hand side is the interfacial
shear stress. The latter contains three physically distinct contributions, to
be discussed next.

The first contribution is the applied (spatially homogeneous) stress
$\tau_{0}(t)$, i.e. this approach is directly applicable to stress-controlled
conditions. The second contribution is the so-called radiation damping term
[67, 35, 38, 68], where $c_{s}$ is the shear wave speed and $v_{0}$ is set as
a reference slip velocity (to be identified with the applied slip velocity
$v_{0}$ in our velocity-controlled setting, see below). The radiation damping
term locally depends on $v(x,t)$ and physically represents plane-waves being
radiated away from the interface into the surrounding bulks, serving as
effective damping from the perspective of the interface. Finally, the third
contribution $s(x,t)$ is non-local in space and time, and physically
represents the spatiotemporal interaction of different points on the
interface, mediated by bulk elastodynamics. $s(x,t)$ generally does not admit
real-space representation, and is related to the gradient of $\delta(x,t)$ in
the spectral domain via a convolution integral, which is different for mode-II
and mode-III [58, 59, 60]. The spectral nature of the formulation fits the
choice of periodic boundary condition in the lateral/sliding direction (with
periodicity $W$), and is reflected in its common name, the spectral Boundary
Integral Method (BIM).

As explained above, and as evident in Eq. (3), the BIM is most suitable for
stress-controlled boundary conditions represented by $\tau_{0}(t)$. Yet, we
propose here an approach in which the BIM formulation can be nevertheless used
to mimic velocity-driven frictional dynamics that are of interest here. The
idea is the following: as explained above, under velocity-driven conditions
_and_ once the system reached steady-state, one has

$\frac{1}{W}\int_{0}^{W}\\!v(x,t)\,dx=v_{0}\ .$ (4)

Consequently, one can choose $\tau_{0}(t)$ — an a priori unknown function of
time $t$ — such that Eq. (4) is satisfied at any time $t$. That is, we propose
to treat $\tau_{0}(t)$ in Eq. (3) as unknown and instead to impose Eq. (4) for
any $t$. Thinking about a numerical implementation of the formulated problem,
where both the time $t$ and the spatial coordinate $x$ are discretized, it is
clear that the above suggestion leads to a well-defined problem; at each
discrete time $t_{i}$, we added a single unknown $\tau_{0}(t_{i})$ and a
single constraint (Eq. (4) at $t_{i}$). Moreover, as $v_{0}$ is the relevant
slip velocity in this modified BIM formulation, we used it as a reference
velocity in the radiation damping term in Eq. (3).

The modified BIM formulation, aimed at mimicking velocity-driven frictional
sliding, has both clear advantages and potential limitations. On the one hand,
it is a relatively computationally cheap and very robust approach, which is
expected to reveal velocity-driven steady-states in the $H\\!\to\\!\infty$
limit (if these exist). On the other hand, as Eq. (4) is strictly valid only
in steady-state — i.e. some deviations from it are expected in early-time, out
of steady-state dynamics — some dynamical aspects of the full velocity-driven
problem may not be accurately captured. We address this potential limitation
in two ways; first, we mainly focus on the long-time, steady-state behavior of
the system, where our modified BIM formulation is strictly valid. Second, we
verify through Finite Element Method (FEM) that the possible deviations of the
transient dynamics in the $H\\!\to\\!\infty$ BIM calculations from the exact
transient dynamics have no effect on the obtained steady-state solutions. More
importantly, the FEM formulation to be discussed next allows us to understand
the roles played by a finite $H$ on the physics of the problem at hand.

###  II.3 Simulating velocity-driven frictional dynamics in finite systems
using the finite element method

The length scale $H$ may play important roles in velocity-driven frictional
sliding. Consequently, it is essential to supplement the $H\\!\to\\!\infty$
BIM calculations with FEM ones, which allow to probe the finite $H$ regime.
Moreover, finite $H$ FEM calculations can be used to verify the validity of
the modified BIM calculations, as discussed above, for selected test cases in
which $H$ is chosen to be large enough (see below). The new physics introduced
by a finite $H$ is the wave interaction of the frictional interface with the
boundaries at $y\\!=\\!\pm H$, which is absent in the $H\\!\to\\!\infty$
limit, where wave are being radiated from the interface to $\pm\infty$ without
reflections.

FEM simulations of frictional systems, especially when rate-and-state friction
is taken into account, are significantly more challenging and computationally
demanding than their BIM counterparts [69]. They are also prone to numerical
and physical instabilities (associated with the finite $H$) that are absent in
their BIM counterparts. One such physical instability is encountered in our
large-$H$ FEM calculations. This issue, along with the full details of our FEM
simulations, are addressed in Appendix D. Agreement between our FEM and BIM
calculations for mode-II (simpler mode-III calculations are performed only
using BIM) is demonstrated in the next section, where we start discussing the
emergent physics of velocity-driven frictional systems, and pulse trains in
particular.

##  III Coarsening dynamics: The selection of pulse train periodicity

The first question we aim at addressing is whether steady-state pulse trains
indeed emerge under velocity-driven sliding conditions and if so, what
determines the train periodicity $W_{\rm p}$. As the latter is of length
dimension, one can a priori ask what quantities of length dimension existing
in the posed problem could possibly determine $W_{\rm p}$. The problem at
hand, as formulated above, features 3 length scales: the system height $H$,
the system length $W$ and in addition to these two geometric/extrinsic length
scales, we also have the intrinsic interfacial length $D$. As $D$ is a
mesoscopic length that characterizes the interfacial response, it does not
appear in itself in the macroscopic coupled bulk-interface problem. Rather, it
appears through the elasto-frictional length $L_{\rm c}$ that characterizes
the rate-weakening ($df_{\rm ss}/dv\\!<\\!0$) elasto-frictional instability
[1, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33, 34].

The elasto-frictional length $L_{\rm c}$ can be determined using a systematic
linear stability analysis [33, 10], which is very briefly reviewed in Appendix
B. Here we are mainly interested in the scaling structure of $L_{\rm c}$,
taking the form

$L_{\rm c}(H,D)=H\,{\cal
F}\\!\left(\\!\frac{\mu\,D}{-\sigma_{0}Hv_{0}\,df_{\rm ss}/dv_{0}}\\!\right)\
,$ (5)

where ${\cal F}({\cal X})\\!\sim\\!{\cal X}$ for ${\cal X}\\!\ll\\!1$ and
${\cal F}({\cal X})\\!\sim\\!\sqrt{\cal X}$ for ${\cal X}\\!\gg\\!1$ [10],
with ${\cal X}\\!\equiv\\!-\mu D/(\sigma_{0}Hv_{0}\,df_{\rm ss}/dv_{0})$.
${\cal X}$ clearly manifests the coupled bulk-interface nature of $L_{\rm c}$,
incorporating the shear modulus $\mu$ of the bodies in contact, their height
$H$ (a geometric/extrinsic length), the applied normal stress $\sigma_{0}$,
the mesoscopic interfacial length $D$ and the constitutive interfacial
property $v_{0}\,df_{\rm ss}/dv_{0}\\!=\\!\tfrac{df_{\rm
ss}(v_{0})}{d\\!\log{v}}\\!<\\!0$. Note that the limiting behaviors of ${\cal
F}({\cal X})$, stated above, imply that $L_{\rm c}$ is independent of $H$ in
the limit $H\\!\to\\!\infty$. The length $L_{\rm c}$ implies that
infinitesimal perturbations of wavelength larger than $L_{\rm c}$, on top of
homogeneous sliding with slip velocity $v_{0}$, are linearly unstable. That
is, for $L_{\rm c}\\!<\\!W$ the system is linearly unstable as unstable
perturbations can fit in (to ensure that this is the case, we use $L_{\rm
c}\\!\ll\\!W$ throughout this work). Moreover, the wavelength of the fastest
growing mode of instability features the same scaling properties as $L_{\rm
c}$ in Eq. (5).

To address the question of whether $W_{\rm p}$ exists, i.e whether steady-
state pulse trains emerge, and if so how $W_{\rm p}$ depends on $L_{\rm
c}(H,D)$, $H$ and $W$, we first explore the $H\\!\to\\!\infty$ limit using the
BIM formulation of Sect. II.2. That is, we first consider the case $L_{\rm
c}\\!\ll\\!W\\!\ll\\!H$. A representative result under mode-III conditions is
presented in Fig. 2a, where a space-time plot of $v(x,t)/v_{0}$ is shown. It
is observed that, as expected, at early times the velocity-driven frictional
system experiences instabilities, resulting in multiple interacting slip
pulses. As time progresses (increasing vertical direction in the space-time
plot), the number of pulses decreases until — in the long-time limit — a
single pulse survives. The latter corresponds to a steady-state pulse train of
periodicity $W$, i.e. $W_{\rm p}\\!=\\!W$ (recall that throughout this work we
employ periodic boundary conditions along the sliding direction). That is, the
frictional system undergoes coarsening dynamics that are saturated at the
system length $W$. To test the robustness of this observation, we varied $W$
and the interfacial parameters over a wide range (cf. Appendix A and C);
independently of these variations, we always observed $W_{\rm p}\\!=\\!W$ in
our $H\\!\to\\!\infty$ calculations, indicating coarsening dynamics that are
truncated at the system length.

Figure 2: (a) A space-time plot of the slip velocity field $v(x,t)$
(normalized by $v_{0}$) obtained using a BIM simulation ($H\\!\to\\!\infty$)
with $v_{0}\\!=\\!3\\!\times\\!10^{-3}\text{m/s}$, $W\\!=\\!28.14\text{m}$ and
the N-shaped steady-state friction law of Fig. 1b. Black regions correspond to
the stick phase, $v(x,t)\\!<\\!v_{0}$. At early times, the system host several
instabilities. The interface progressively coarsens with increasing time until
a single pulse remains in the periodic domain. The asymptotic pulse
propagation velocity $c_{\rm p}/c_{\rm s}$ is shown as the slope of the slip
phase $v(x,t)\\!>\\!v_{0}$. The normalized pulse width $w_{p}/W$ corresponds
to the extent of the slip phase along the $x$ axis at a given time $t$. A
movie of this simulation is available in [70]. (b) A snapshot of $v(x,t)$ of a
steady-state pulse train, comparing a BIM simulation (dashed orange line) and
a FEM simulation (solid brown line) with
$v_{0}\\!=\\!1\\!\times\\!10^{-3}\text{m/s}$, $W\\!=\\!6\text{m}$ and the
steady-state friction law of Fig. S1. The dashed-dotted horizontal green line
corresponds to the driving velocity $v_{0}$. Both pulses are propagating
towards the left at the same velocity $c_{\rm p}\\!\simeq\\!0.35c_{s}$ and
their shapes are almost identical.

Next, we aim at understanding whether a finite height $H$ can result in a
qualitative change in this physical picture, i.e. whether the pulse train
periodicity $W_{\rm p}$ can be affected/determined by $H$. To address this
question, one clearly needs to resort to FEM calculations. Due to numerical
stability considerations, we use a rate-and-state friction law that does not
feature the very low velocity rate-strengthening branch, see Appendix D.2 for
discussion and details. In particular, we focus on the regime where
$H\\!\simeq\\!W\\!\gg\\!L_{\rm c}$ in order to see whether $H$ affects $W_{\rm
p}$ in addition to $W$ (as they are comparable). In this regime, we
encountered in addition to the elasto-frictional instability associated with
$L_{\rm c}$ also a system-size instability. In this instability, the slip
velocity vanishes throughout the interface for some period of time, and then
interface starts sliding almost _homogeneously_. That is, this instability
appears to feature a space-independent stick-slip behavior (see Appendix D.3),
which is likely to be related to the finite-$H$ elasto-frictional
instabilities discussed in [71]. We have not studied this instability in depth
as we strongly suspect that it does not manifest itself in the limit of large
$W$ (which is not accessible in our FEM calculations) and large times, i.e. it
does not affect the long-time behavior of the system. Instead, we effectively
eliminated this space-independent instability by breaking translational
symmetry along the interface using a constraint similar to Eq. (4), as
explained in Appendix D.3.

Figure 3: (a) A space-time plot as in Fig. 2a for a FEM simulation with
$v_{0}\\!=\\!2\\!\times\\!10^{-3}\text{m/s}$, $W\\!=\\!4\text{m}$,
$H\\!=\\!0.2\text{m}$ and the dashed-dotted orange steady-state friction law
of Fig. S1. The observed behavior is similar to the BIM example
($H\\!\to\\!\infty$) of Fig. 2a, with the notable difference that the pulse
train that remains in the system in the long-time coarsening limit features
repeated oscillations (visible from the variation of the width and slip rate
of the pulse). A movie of this simulation is available in [70]. (b) Successive
snapshots of $v(x)/v_{0}$ are presented in the co-moving frame of the pulse
train (with $x_{\rm p}$ being the pulse’s leading edge). The pulse profiles
are ordered in time as follows (line colors): dashed yellow, dashed-dotted
orange, dotted brown, and solid dark brown. The dashed-dotted horizontal green
line corresponds to the driving velocity $v_{0}$. As time progresses, both the
pulse width $w_{\rm p}$ and the slip velocity increase, indicating an elasto-
frictional instability. A movie of the instability in the co-moving frame of
the pulse is available in [70]. (inset) The time evolution of the normalized
maximum velocity $v_{\rm m}$ (solid black curve) of the dynamics shown in
panel (b), presented on a semi-logarithmic scale. The green dashed line is the
best linear fit, corresponding to an exponential growth characterizing a
linear instability, with the growth rate being the slope of this line. $t_{\rm
c}$ corresponds to the time at which the growth becomes nonlinear. (c) The
time evolution of $v_{\rm m}/v_{0}$. The vertical lines indicate the time of
the snapshots in panel (b), with the same style and color codes. Two events of
this repeated behavior (growth and decay of the pulse) are shown.

In Fig. 2b, we present a snapshot of the long-time behavior of an FEM
simulation (as described above) of a frictional system with $H\\!=\\!2W$
(solid line), where yet again a steady-state pulse train with $W_{\rm
p}\\!=\\!W$ is observed. This result is further strengthened by additional FEM
calculations where the value of $H$ has been varied for a fixed $W$, still
having $H\\!\simeq\\!W\\!\gg\\!L_{\rm c}$, yielding $W_{\rm p}\\!=\\!W$
independently of $H$. The independence of $W_{\rm p}$ on $H$ also suggests
that for a fixed set of frictional parameters and $W$, the $H\\!=\\!2W$ FEM
calculation and the $H\\!\to\\!\infty$ BIM calculation should give rise to
very similar results. Such a quantitative agreement is possible if the
radiation damping term is not too large. This is the case if $c_{\rm p}$ is
not close to the relevant wave speed.

In Fig. 2b, we superimposed the $H\\!\to\\!\infty$ BIM results (dashed line,
for mode-II as in the FEM calculation) on top of the $H\\!=\\!2W$ FEM ones
(the two snapshots are shifted along $x$ for visual clarity). It is observed
that in both cases $W_{\rm p}\\!=\\!W$ and that the pulse shape is almost
identical (the same applies to the pulse train velocity $c_{\rm
p}\\!\simeq\\!0.35c_{s}$, see figure caption). These convergent results
demonstrate the robustness and validity of our FEM calculations. More
importantly, they show that velocity-driven frictional systems exhibit
coarsening dynamics that are saturated at the system length, independently of
the system’s height $H$ in the regime $L_{\rm c}\\!\ll\\!W\\!\leq\\!H$,
leading to steadily propagating pulse trains with $W_{\rm p}\\!=\\!W$. In
light of the independence of $W_{\rm p}$ on $H$, we would like next to
understand whether $L_{\rm c}$ can affect $W_{\rm p}$, which requires to
consider the small $H$ regime.

##  IV The small $H$ limit: The competition between coarsening and elasto-
frictional instabilities

The results of the previous section strongly suggest that we have $W_{\rm
p}\\!=\\!W$ independently of $H$ for $L_{\rm c}\\!\ll\\!W\\!\leq\\!H$, where
$L_{\rm c}$ does not seem to play a role in the long-time pulse train as well.
To address the possible role of $L_{\rm c}$ in the selection of the pulse
train periodicity $W_{\rm p}$, we explore here the small $H$ limit for which
$H\\!<\\!L_{\rm c}\\!\ll\\!W\\!$. Before discussing our FEM calculations in
this regime, let us invoke some theoretical considerations regarding their
possible outcomes. One possibility is that we find steady-state pulse trains
with $W_{\rm p}\\!\simeq\\!L_{\rm c}$. While we cannot a priori rule out this
possibility, we note that — if true — it implies that coarsening dynamics play
no role whatsoever in this regime. It is not easy to imagine this in light of
the fact that coarsening dynamics strongly dominate the physics for $L_{\rm
c}\\!\ll\\!W\\!\leq\\!H$. Another possibility is that coarsening dynamics
remain dominant, leading to $W_{\rm p}\\!=\\!W$ as before. This, however,
cannot lead to a strict steady-state.

To understand this statement, let us consider some of the properties of the
slip velocity field inside a pulse within a pulse train, for both large and
small $H$. In the former case, the slip velocity inside the pulse is
significantly amplified compared to $v_{0}$ over a large fraction of the pulse
width $w_{\rm p}$, such that most of the pulse is characterized by a slip
velocity out of the unstable rate-weakening branch of the friction curve. This
is demonstrated in the pulse shown in Fig. 2b (see Sect. V.1 for additional
discussion of the slip velocity amplification). This is qualitatively
different from the small $H$ case, $H\\!<\\!L_{\rm c}\\!\ll\\!W$; here, the
slip velocity is amplified over a scale $H$ that is significantly smaller than
$w_{\rm p}$. This situation is closely related to the finite height strip
problem of fracture mechanics [64], where the controlling length scale is the
strip height $H$. This implies that a significant fraction of the
pulse/interface is characterized by $v(x,t)\\!\simeq\\!v_{0}$, and since
$w_{\rm p}\\!\gg\\!L_{\rm c}$ (see Sect. V.2 for additional discussion of the
properties of $w_{\rm p}$), this extended region is characterized by a slip
velocity belonging to the rate-weakening branch of the friction curve and
hence is unstable. Consequently, if $W_{\rm p}\\!=\\!W$ also for
$H\\!<\\!L_{\rm c}\\!\ll\\!W$, then we expect to observe repeated
instabilities on top of the pulse train. That is, we expect the competition
between coarsening and elasto-frictional instabilities to prevent the system
from reaching a strict steady-state.

In Fig. 3a we present a space-time plot (in the same format of Fig. 2a)
obtained using FEM simulations in the $H\\!<\\!L_{\rm c}\\!\ll\\!W$ regime. It
is important to note that in this regime we do not observe the space-
independent stick-slip instability discussed in Sect. III and hence no
additional constraints that break translational symmetry along the interface
are introduced. The results in Fig. 3a indicate that coarsening dynamics
dominate in this regime as well, leading to a pulse train with $W_{\rm
p}\\!=\\!W$. This pulse train, however, is not in strict steady-state as it is
clearly observed to be repeatedly/periodically interrupted by bursts of large
slip velocities (compare to Fig. 2a, where these bursts are absent and the
system reaches a strict steady-state). These observations support the physical
picture discussed in the previous paragraph, where coarsening dynamics compete
and coexist with repeated elasto-frictional instabilities in the
$H\\!<\\!L_{\rm c}\\!\ll\\!W$ regime.

In order to better understand these complicated spatiotemporal dynamics, we
focus in Fig. 3b on a sequence of snapshots of $v(x)$ of the pulse within the
pulse train observed in the long-time dynamics of Fig. 3a, during one of the
bursts (the snapshots are shown in the pulse train co-moving frame, see figure
caption). The snapshots are ordered in time, where the smallest $v(x)$
corresponds to early time and the largest to late time. The time points in
which the snapshots were taken are marked by vertical lines in Fig. 3c (same
colors and line styles as in the curves of Fig. 3b), where the maximal slip
velocity $v_{\rm m}$ is plotted as a function of time. The early time field
reveals a long plateau, larger than $L_{\rm c}$ (recall that here $L_{\rm
c}\\!\ll\\!W$), featuring $v(x)\\!\simeq\\!v_{0}$, as predicted above.
Consequently, an elasto-frictional instability is indeed expected.

As is evident from both Fig. 3b and Fig. 3c, the slip velocity inside the
pulse grows significantly, indicating an elasto-frictional instability. A
clear signature of the linear elasto-frictional instability associated with
$L_{\rm c}$ is an exponential growth of the slip velocity at the early stages
of the instability development. This is indeed demonstrated in the inset of
Fig. 3b, supporting the physical picture discussed above. Moreover, the
observed exponential growth rate appears to be in the ballpark of the
theoretically estimated growth rate (not shown). Finally, the repeated nature
of the instability, already observed in Fig. 3a, is evident also in Fig. 3c.
There, the maximal slip velocity $v_{\rm m}$ is plotted as a function of time,
and it is observed that after the instability shown in Fig. 3b relaxes,
another instability with very similar properties occurs after some interval of
time.

The results of this section and of the previous one suggest that velocity-
driven frictional systems, at least the rather generic class of frictional
systems considered in this paper, are strongly dominated by coarsening
dynamics that lead to pulse trains characterized by periodicity determined by
system length in the sliding direction. For $H\\!<\\!L_{\rm c}\\!\ll\\!W$, the
train pulse experiences repeated elasto-frictional instabilities, but for
larger $H$’s the pulse train is a stable steady-state solution. It is
important to note that in this paper we do not study in depth and
quantitatively the coarsening dynamics themselves (this will be done
elsewhere), but rather focus on their long-time outcome. Next, we shift our
focus to other salient features of the observed pulse trains, i.e. their
propagation velocity $c_{\rm p}$ and the properties of single pulses within
the train.

##  V Single pulse properties in the large $H$ limit

Self-healing slip pulses, as explained in Sect. I, are frictional rupture
modes that are believed to be quite prevalent, yet they are not fully
understood [47, 38, 52, 34, 55, 49]. Our analysis above established that self-
healing slip pulses, as part of pulse trains, generically and robustly emerge
in velocity-driven frictional sliding dynamics. These results offer a rather
unusual opportunity to better understand the physical properties of single
pulses, which is the main goal of this section.

A steady-state pulse train is composed of single pulses that repeat themselves
with spatial periodicity $W_{\rm p}$, extensively discussed above. Once the
train periodicity $W_{\rm p}$ is known — it was shown above to be coarsening-
limited (i.e. equal to the system length, $W_{\rm p}\\!=\\!W$) —, one is
interested in the selection of the train velocity $c_{\rm p}$ (which is
obviously also the single pulse propagation velocity) and in the spatial
distributions of slip velocity $v(x)$ and stress $\tau(x)$ within the single
pulse. Obtaining closed-form solutions for the field distributions is a
difficult challenge not addressed here; instead, we focus below on the
behavior of $v(x)$ and $\tau(x)$ near the leading edge of the single pulse,
and on some characteristic properties of $v(x)$ and $\tau(x)$. Mort
importantly, as already introduced in Fig. 1, we are interested in the single
pulse width $w_{\rm p}$ and in its average slip velocity $v_{\rm p}$.

The single pulse width $w_{\rm p}$, as illustrated in Fig. 1c, is defined as
the size of the portion of the single pulse for which $v(x)\\!\geq\\!v_{0}$.
Accordingly, the average slip velocity of a single pulse is defined as $v_{\rm
p}\\!\equiv\\!w_{\rm p}^{-1}\\!\int_{0}^{w_{\rm p}}v(\tilde{x})d\tilde{x}$,
where $\tilde{x}\\!=\\!0$ in this context corresponds to the spatial point on
the pulse’s leading edge such that $v(\tilde{x}\\!=\\!0)\\!=\\!v_{0}$ (the
corresponding relation at the trailing edge reads $v(\tilde{x}\\!=\\!w_{\rm
p})\\!=\\!v_{0}$). Our goal in this section is to gain physical insight into
the single pulse quantities $c_{\rm p}$, $w_{\rm p}$ and $v_{\rm p}$ as a
function of the control parameters $W$ and $v_{0}$ for the adopted interfacial
constitutive law, in the $H\\!\to\\!\infty$ limit. Moreover, we would like to
understand the interrelations between $c_{\rm p}$, $w_{\rm p}$ and $v_{\rm
p}$, and clearly three such relations are needed in order to determine these
three quantities.

One such relation is readily obtained from the steady-state condition of Eq.
(4). The latter implies

$W^{-1}\\!\left[w_{\rm p}v_{\rm p}+(W-w_{\rm p})v_{\rm
s}\right]\\!=\\!v_{0}\quad\\!\Longrightarrow\\!\quad w_{\rm p}v_{\rm
p}\\!\simeq\\!Wv_{0}\ ,$ (6)

where the slip velocity in the ‘stick region’, $v_{\rm s}$ (cf. Fig. 1c), has
been assumed to be negligibly small (as is indeed the case). Note that Eq. (6)
does not explicitly include the propagation velocity $c_{\rm p}$. In Sect. V.1
we show that another relation between $c_{\rm p}$, $w_{\rm p}$ and $v_{\rm p}$
emerges from the pulse’s leading edge behavior, a relation that may be viewed
as a pulse equation of motion. A third relation is discussed in Sect. V.2,
where we focus on the pulse width $w_{\rm p}$.

###  V.1 Pulse equation of motion: Leading edge singularity, effective
fracture energy and propagation velocity

A slip pulse, like other rupture fronts, is expected to feature nearly
singular fields near its leading edge, where the negligibly small slip
velocity ahead of the pulse, $v_{\rm s}$, dramatically increases to the peak
velocity, $v_{\rm m}$ (cf. Fig. 1c), over a short length scale. It has been
recently shown [9, 11, 72, 73] that for the class of rate-and-state friction
constitutive laws illustrated in Fig. 1b, where the rate dependence of
friction is rather weak, the singular leading edge behavior is described by
the classical square root singularity of Linear Elastic Fracture Mechanics
(LEFM) [61] to a good approximation.

In the framework of LEFM, the slip velocity behind the leading pulse edge is
expected to take the following universal form [61]

$v(x)\simeq\frac{2\,c_{\rm p}\,K_{\rm III}}{\alpha_{\rm s}(c_{\rm
p})\,\mu\,\sqrt{2\pi(x-x_{\rm p})}}\ ,$ (7)

where $\alpha_{\rm s}(c_{\rm p})\\!=\\!\sqrt{1-(c_{\rm p}/c_{s})^{2}}$,
$x_{\rm p}$ is the leading edge position and $K_{\rm III}$ is the mode-III
stress intensity factor, quantifying the amplitude of the universal square
root singularity (Eq. (7) refers to mode-III conditions, which we focus on in
this section, but the same singularity applies to mode-II conditions [61]). In
Fig. 4 we plot $v(x)$ (properly normalized, see legend) over the entire pulse
train period $W_{\rm p}\\!=\\!W$ (solid brown line, left $y$-axis). To test
whether the observed amplified slip velocity near the leading edge indeed
follows the square root singularity of Eq. (7), we plot in the inset the same
normalized $v(x)$ behind the leading edge against $(x-x_{\rm p})/W$ on a
double-logarithmic scale. As is indicated by the dashed line of slope
$-\tfrac{1}{2}$, there exists a spatial range near the leading edge over which
$v(x)$ follows the square root singular behavior of Eq. (7).

Figure 4: The normalized (see legend) profiles of the velocity $v(x)$ (solid
brown line, left $y$-axis), the shear/frictional stress $\tau(x)$ (dashed
orange, left $y$-axis) and $v(x)\phi(x)/D$ (dashed-dotted green line, right
$y$-axis) of a steady-state pulse train propagating to the left at $c_{\rm
p}\\!=\\!0.85c_{\rm s}$, as obtained by a BIM simulation ($H\\!\to\\!\infty$)
with $v_{0}\\!=\\!3\\!\times\\!10^{-3}\text{m/s}$, $W\\!=\\!14.07\text{m}$ and
the N-shaped steady-state friction law of Fig. 1b. The leading edge is
indicated by $x_{\rm p}/W$ and the pulse width by $w_{\rm p}/W$ (and a double-
arrow). Both the stress and the velocity fields reveal a singular behavior
near the pulse leading edge, see text for discussion. (inset) The velocity
field $v(x)$ (solid brown line) behind the leading edge, on a double-
logarithmic scale, with the $x$ axis being $(x-x_{\rm p})/W$. The dashed black
line features a slope of $-1/2$, characteristic of the square root singularity
of LEFM, see text for additional discussion.

The singular slip velocity behavior of Eq. (7) behind the pulse’s leading edge
is accompanied by a shear/frictional stress field $\tau(x)$ that features the
same square root singularity _ahead_ of the leading edge and a predominantly
constant residual stress behind it [9, 11]. In Fig. 4 we present $\tau(x)$
(properly normalized, see legend) again over the entire pulse train period
$W_{\rm p}\\!=\\!W$ (dashed orange line, sharing the _same_ left $y$-axis with
the normalized $v(x)$ field). It is observed that indeed $\tau(x)$ features
significant amplification ahead of the leading edge (its quantification is not
presented here) and that it is approximately constant behind the leading edge
inside the pulse. Finally, we present in Fig. 4 $v(x)\phi(x)/D$ (dashed-dotted
green line, right $y$-axis), whose deviation from unity provides a measure of
the degree by which $\phi(x)$ is out of steady-state with respect to $v(x)$.
It is observed that inside the pulse $v(x)\phi(x)/D\\!\simeq\\!1$ (see
horizontal dashed line) almost everywhere. That is, inside the pulse away from
the leading edge transition region, $\phi(x)$ is “equilibrated” with $v(x)$.

The approximate validity of Eq. (7) has significant implications for the pulse
dynamics. Most notably, the square root singular fields are accompanied by a
finite (non-singular) flux of energy $G$ into the leading edge region. The
latter is dissipated near the leading edge, i.e. it is balanced by an
effective fracture energy $G_{\rm c}$, which quantifies the dissipation
involved in slip pulse propagation (on top of the background frictional
dissipation [74, 11]). That is, we expect the following relation [61]

$G_{\rm c}\simeq G=\frac{K^{2}_{\rm III}}{2\alpha_{\rm s}(c_{\rm p})\,\mu}\ ,$
(8)

to approximately hold. The effective fracture energy $G_{\rm c}$ can be
independently estimated using the interfacial constitutive law, following [24,
25] and as shown in Appendix E. With $G_{\rm c}$ at hand, one can eliminate
$K_{\rm III}$ between Eqs. (7)-(8), yielding an expression for the singular
part of $v(x)$ near the leading edge of the pulse.

The resulting $v(x)$ allows to derive a relation between $c_{\rm p}$, $w_{\rm
p}$ and $v_{\rm p}$. It is important to note that the latter two quantities
characterize the whole pulse, not just its singular part. To bridge over this
gap, we assume that the singular part of $v(x)$ describes reasonably well the
entire pulse, up to an overall shift that dominates the trailing edge
behavior. The shift can be determined by demanding $v(x_{\rm p}\\!+\\!w_{\rm
p})\\!=\\!v_{0}$, which follows from the definition of $w_{\rm p}$, leading to

$\frac{v(x)}{c_{s}}\\!\simeq\\!\sqrt{\frac{4G_{c}}{\pi\mu}}\frac{c_{p}}{c_{s}\sqrt{\alpha_{s}(c_{p})}}\left(\frac{1}{\sqrt{x-x_{\rm
p}}}\\!-\\!\frac{1}{\sqrt{w_{p}}}\right)+\frac{v_{0}}{c_{s}}\ .$ (9)

With this approximation at hand, we can then calculate $v_{\rm p}$ according
to its definition $v_{\rm p}\\!=\\!w_{\rm p}^{-1}\\!\int_{x_{\rm p}}^{x_{\rm
p}+w_{\rm p}}v(\tilde{x})d\tilde{x}$, leading to

$\frac{c_{\rm p}}{c_{s}\sqrt{\alpha_{s}(c_{\rm
p})}}\sim\sqrt{\frac{\pi\,\mu\,w_{\rm p}}{G_{\rm c}}}\,\frac{v_{\rm
p}-v_{0}}{2c_{s}}\ .$ (10)

The prediction in Eq. (10) is a relation between $c_{\rm p}$ (left hand side),
and $w_{\rm p}$ and $v_{\rm p}$ (right hand side), where the pre-factor is
expected to be ${\cal O}(1)$. It may be viewed as an equation of motion for
the pulse. Moreover, note that Eq. (6) (which has been independently verified)
can be used to eliminate $v_{\rm p}$ from the right hand side of Eq. (10),
making it a function of $w_{\rm p}$ alone. Equation (10) is tested in Fig. 5
over a broad range of $W$ and $v_{0}$ values. It is observed that the
theoretical prediction is strongly supported by the simulational data, where
the pre-factor is indeed ${\cal O}(1)$. Consequently, only one additional
relation between $c_{\rm p}$, $w_{\rm p}$ and $v_{\rm p}$ is needed in order
to fully determine the single pulse properties. In particular, predicting the
pulse width $w_{\rm p}$ in terms of $W$ and $v_{0}$ would be sufficient since
Eqs. (6) and (10) would then allow to calculate $v_{\rm p}$ and $c_{\rm p}$,
respectively. Therefore, we next consider $w_{\rm p}(W,v_{0})$.

Figure 5: Testing the theoretical prediction of Eq. (10) for a broad range of
parameters, including variations of the system length $W$, the driving
velocity $v_{0}$, the shear modulus $\mu$ and the characteristic slip
displacement $D$, see Appendix A and C for additional details and the exact
values used. The data correspond to BIM simulations with the N-shaped friction
law of Fig. 1b, revealing very good agreement with the prediction, with an
${\cal O}(1)$ pre-factor. $G_{\rm c}$ has been estimated following [24, 25],
see Appendix E for more details.

###  V.2 The pulse width

Our goal in this subsection is to discuss the pulse width $w_{\rm
p}(W,v_{0})$. It would be useful in this context, as a preparatory step, to
consider the “ideal” solution for an isolated pulse. The latter corresponds to
an isolated pulse in an infinite medium ($H\\!\to\\!\infty$) whose slip
velocity $v^{(i)}(x)$ (the superscript $(i)$ denotes hereafter “ideal pulse”)
vanishes for $x\\!<\\!x_{\rm p}$ and $x\\!>\\!2w^{(i)}_{\rm p}+x_{\rm p}$, and
the frictional strength is constant for $x_{\rm
p}\\!\leq\\!x\\!\leq\\!2w^{(i)}_{\rm p}+x_{\rm p}$. The solution takes the
form [62, 63, 64]

$v^{(i)}(x)=v_{0}\sqrt{\frac{2w^{(i)}_{\rm p}+x_{\rm p}-x}{x-x_{\rm p}}}$ (11)

for $x_{\rm p}\\!\leq\\!x\\!\leq\\!2w^{(i)}_{\rm p}+x_{\rm p}$. Note that the
solution in Eq. (11) features the square root singularity of Eq. (7) in the
$x\\!\to\\!x_{\rm p}^{+}$ limit and it is consistent with the single pulse
width definition adopted throughout this paper. That is, we have
$v^{(i)}(x\\!=\\!w^{(i)}_{\rm p}\\!+\\!x_{\rm p})\\!=\\!v_{0}$ such that
$w^{(i)}_{\rm p}$ retains its meaning as the pulse width.

The physics underlying the ideal pulse solution in Eq. (11) differs from the
physics underlying the train pulse solutions derived above in three major
respects: (i) Our solutions correspond to rate-and-state friction, while Eq.
(11) corresponds to a _constant_ frictional strength inside the pulse. (ii)
Eq. (11) corresponds to an isolated pulse, i.e. to a single pulse that does
not interact with other pulses in the train over the periodicity scale $W_{\rm
p}\\!=\\!W$. (iii) $w^{(i)}_{\rm p}$ in Eq. (11) is not selected, an issue
that is obviously related to points (i) and (ii).

An attempt to remedy points (ii) and (iii) above, which is inevitably
_superficial_ , would be to impose the steady-state pulse train condition of
Eq. (4) on Eq. (11), resulting in $w^{(i)}_{\rm p}/W\\!=\\!1/\pi$. That is, we
construct an approximate pulse train solution by concatenating/superimposing
the ideal pulse solutions of Eq. (11) with periodicity $W_{\rm p}\\!=\\!W$ and
a single pulse width $w^{(i)}_{\rm p}/W\\!=\\!1/\pi$. Below we show that this
artificially constructed pulse train qualitatively fails to describe our rate-
and-state friction pulse train solutions, where $w_{\rm p}/W$ is not a
constant, but rather depends on both $W$ and $v_{0}$. Moreover, while both the
ideal pulse solution of Eq. (11) and our pulses (cf. Figs. 2b and 4) do not
feature a singularity at the trailing edge, the pulse shape $v(x)$ near the
trailing edge differs quite significantly (not shown). Another difference
between the two solutions is that while the ideal pulse solution is truly
singular at the leading edge, the singularity in our solutions is self-
consistently regularized on small scales near the leading edge.

In light of the discussion of the ideal pulse above, which mainly served to
highlight the differences compared to our pulse train solutions, our next goal
is to better characterize the function $w_{\rm p}(W,v_{0})/W$. The first
question we need to address is how to properly nondimensionalize $W$ and
$v_{0}$, i.e. the dimensionless pulse width $w_{\rm p}/W$ should be expressed
in terms of dimensionless quantities. We suggest that the proper way to
nondimensionalize $W$ and $v_{0}$ is through their relation to the elasto-
frictional instability, which is a necessary condition for the emergence of
pulse trains.

Figure 6: The dimensionless pulse width $w_{\rm p}/W$ as a function of both
the normalized driving velocity $v_{0}/v_{\rm min}$ and system size $L_{\rm
c}/W$, see text for additional discussion. Each black circle corresponds to
one of the simulations shown in Fig. 5, spanning a broad range of system
parameters (including variations of the system size $W$, driving velocity
$v_{0}$, shear modulus $\mu$ and the characteristic slip displacement $D$).
Variations of $\mu$ alter the elasto-frictional length $L_{\rm c}$, while
variations of $D$ affect both $L_{\rm c}$ and $v_{\rm min}$. The surface is
interpolated from the data points, supporting the existence of a smooth and
monotonically increasing function $w_{\rm p}/W\\!=\\!{\cal G}(L_{\rm
c}/W,v_{0}/v_{\rm min})$, as predicted theoretically (see text for
discussion).

The elasto-frictional instability exists if the system is driven at a velocity
$v_{0}$ that resides on the rate-weakening branch of the steady-state friction
curve $f_{\rm ss}(v)$ (cf. Fig. 1b), i.e. when $v_{0}$ is larger than the slip
velocity at which $f_{\rm ss}(v)$ attains its maximum and smaller than $v_{\rm
min}$. In the limit $v_{0}\\!\to\\!v_{\rm min}$, the system will approach a
stable homogeneous sliding state; that is, as $v_{0}/v_{\rm min}$ is increased
toward unity, we expect the pulse width to approach $W$ (and obviously its
slip velocity to approach $v_{0}$). For $v_{0}\\!<\\!v_{\rm min}$, the elasto-
frictional instability manifests itself only for $L_{\rm c}\\!<\\!W$, i.e.
when an unstable perturbation can fit into the system. Consequently, in the
limit $W\\!\to\\!L_{\rm c}^{+}$ the system will approach a stable homogeneous
sliding state; that is, as $L_{\rm c}/W$ is increased toward unity, we expect
the pulse width to approach $W$.

Taken together, these physical considerations lead us to expect that there
exists a smooth function $w_{\rm p}/W\\!=\\!{\cal G}(L_{\rm c}/W,v_{0}/v_{\rm
min})$ that monotonically increases with its two arguments. This prediction is
tested in Fig. 6, where $w_{\rm p}/W$ is plotted against $L_{\rm c}/W$ and
$v_{0}/v_{\rm min}$. It is observed that a smooth $w_{\rm p}/W\\!=\\!{\cal
G}(L_{\rm c}/W,v_{0}/v_{\rm min})$ function seems to exist, and that $w_{\rm
p}/W$ indeed increases with both of its arguments. The theoretical derivation
of the precise functional form of $w_{\rm p}/W\\!=\\!{\cal G}(L_{\rm
c}/W,v_{0}/v_{\rm min})$ is left as a challenge for future research. With this
analysis, we conclude our discussion of the single pulse properties $c_{\rm
p}$, $w_{\rm p}$ and $v_{\rm p}$ in the large $H$ limit.

##  VI Summary and Discussion

In this work we extensively studied velocity-driven frictional systems, using
a combination of computational and theoretical approaches. Velocity-driven
frictional systems exist in a broad range of engineering and tribological
applications, as well as in geophysical contexts. We show that such frictional
systems, described within the experimentally motivated rate-and-state friction
constitutive framework, give rise to the emergence of pulse trains once driven
at a sliding velocity for which the frictional interface is rate-weakening.
The propagating pulse train is a non-equilibrium dissipative analog of
equilibrium phase separation in thermodynamic systems. We find that such
velocity-driven frictional systems undergo coarsening dynamics leading to
train periodicity set by the size of the system in the sliding direction,
independently of the height of the system. The properties of single pulses
within the train are quantitatively and comprehensively analyzed and
interpreted. In particular, the pulse propagation velocity, width and average
slip rate are shown to be related through an equation of motion that is
associated with the nearly singular fields in the vicinity of the leading
pulse edge.

Throughout this work we employed periodic boundary conditions in the sliding
direction. Similar periodic boundary conditions have been employed in the
velocity-driven frictional simulations of [51], where instability has been
induced by various physical processes leading to a reduction in the effective
normal stress, rather than by rate-weakening friction per se (note that the
possible emergence of pulse trains under rate-strengthening friction has been
discussed in [41]). While coarsening dynamics have not been discussed in [51],
the results therein (cf. Fig. 4 in [51]) do indicate similar coarsening
dynamics as found here. Such periodic boundary conditions naturally emerge in
annular/rotary shear geometry experiments [75, 76], which is also
characteristic of many engineering systems. Such an experimental setup offers
a natural test bed for our predictions. Moreover, we believe that the
coarsening dynamics revealed in this work are physically relevant for
sufficiently long systems (large $W$) that do not feature periodic boundary
conditions, as well as for geophysical fault dynamics, issues that should be
addressed in future work.

This work offers a rather robust framework for studying self-healing pulses,
which are of importance in many frictional systems. Our results also pose
several interesting questions. First, calculating the single pulse shape and
consequently its width $w_{\rm p}$ (cf. Fig. 6) remain open challenges.
Second, we have not studied in any detail the coarsening dynamics themselves,
but rather focused on their long-time outcome. Yet, earlier time coarsening
(cf. Fig. 2a) may be characterized by rather complex spatiotemporal dynamics,
possibly chaotic ones. Quantitatively analyzing these dynamics and their
acoustic emission signature may be important for understanding the non-steady
frictional response of various physical systems. Moreover, in an even broader
statistical physics context, such complex dynamics may reveal how anomalous
statistical properties (e.g. fat power-law statistical distributions [77, 78])
spontaneously emerge even in the absence of input disorder (quenched or
thermal). These interesting issues will be addressed in future investigations.

Acknowledgements E. B. acknowledges support from the Israel Science Foundation
(Grant no. 1085/20), the Ben May Center for Chemical Theory and Computation
and the Harold Perlman Family. J.-F. M. and T. R. acknowledge internal support
from EPFL.

Appendices

##  Appendix A The interfacial constitutive law

The rate-and-state friction constitutive framework [5, 1, 3, 4] has been
formulated in Eqs. (1)-(2). Here we specify the constitutive functions
$f(\cdot)$ and $g(\cdot)$, which appear in these two equations respectively,
taking the explicit forms [5, 1, 3, 4, 8, 79, 56, 80, 10, 33, 49, 9]

$\begin{split}f^{\rm
N}(|v|,\phi)&=\left[1+b\log\left(1+\frac{\phi}{\phi_{*}}\right)\right]\\\
&\times\left[\frac{f_{0}}{\sqrt{1+(v_{*}/|v|)^{2}}}+a\log\left(1+\frac{|v|}{v_{*}}\right)\right]\end{split}$
(1)

and

$g(|v|,\phi)=1-\frac{|v|\phi}{D}\sqrt{1+(v_{*}/v)^{2}}\ .$ (2)

These functions feature 6 parameters: $f_{0}$, $a$, $b$, $D$, $v_{*}$ and
$\phi_{*}$, whose values are specified in Table 1. The superscript ‘N’ in Eq.
(1) denotes the fact that the steady-state friction curve $f_{\rm ss}(v)$
corresponding to Eqs. (1)-(2) features an N shape, as shown in Fig. 1b and
again here in Fig. S1 (solid brown curve).

The N-shaped steady-state friction curve features rate-strengthening at very
low velocities, rate-weakening at intermediate velocities and again rate-
strengthening at high velocities, i.e. it follows a strengthening-weakening-
strengthening (SWS) sequence with increasing slip velocities. As highlighted
in the text, the most crucial branch for the emergence of pulse trains is the
rate-weakening one. Consequently, we also considered two variants of SWS-
related function $f^{\rm N}(|v|,\phi)$ of Eq. (1). One variant corresponds to
replacing Eq. (1) with

$f^{\rm
WS}(|v|,\phi)\\!=\\!f_{0}\\!\left[1\\!+\\!b\log\left(\\!1\\!+\\!\frac{\phi}{\phi_{*}}\\!\right)\\!\right]\\!+\\!a\log\left(\\!1\\!+\\!\frac{|v|}{v_{*}}\\!\right)\
.$ (3)

The resulting steady-state curve is rate-independent at very low velocities,
rate-weakening at intermediate velocities and rate-strengthening at high
velocities. Consequently, we use the superscript ‘WS’ for this function and
plot the corresponding steady-state curve in Fig. S1 (dotted-dashed orange
curve).

The second variant corresponds to omitting the ‘+1’ in the logarithm that
multiplies $b$ in the first square brackets on the right hand side of Eq. (1).
This modification eliminates the high velocities rate-strengthening branch of
the steady-state friction curve [56]. Consequently, we term the resulting
$f(\cdot)$ that replaces Eq. (1) $f^{\rm SW}(|v|,\phi)$ and plot the
corresponding steady-state curve in Fig. S1 (dashed green curve). All three
steady-state curves, corresponding to $f^{\rm N}(|v|,\phi)$, $f^{\rm
WS}(|v|,\phi)$ and $f^{\rm SW}(|v|,\phi)$, essentially share the same rate-
weakening behavior at intermediate slip velocities. Results for $f^{\rm
WS}(|v|,\phi)$ and $f^{\rm SW}(|v|,\phi)$ are presented below.

Parameter |  Value | Unit
---|---|---
$f_{0}$ |  0.28 | …
$a$ |  0.005 | …
$b$ |  0.075 | …
$D$ |  5 $\times$ 10-7 | m
$v_{*}$ |  1 $\times$ 10-7 | m/s
$\phi_{*}$ |  3.3 $\times$ 10-4 | s
Table 1: Typical values of rate-and-state parameters used in this work. In
addition, $D$ has been varied from $1.25\\!\times\\!10$-7 to
$1.5\\!\times\\!10$-6m in the results appearing in Figs. 5, 6, S3.  Figure S1:
The normalized steady-state frictional strength $f_{\rm ss}$ vs. the steady-
state normalized slip velocity $v_{\rm ss}$, presented on a semi-logarithmic
scale for three constitutive laws. The first (solid brown curve), which
features an N shape and is denoted by $f^{\rm N}_{\rm ss}(v_{\rm ss})$ (see
legend), corresponds to Eqs. (1)-(2) and has already been presented in Fig.
1b. The local minimum of $f^{\rm N}_{\rm ss}(v_{\rm ss})$, denoted by $v_{\rm
min}$, is used to normalize $v_{\rm ss}$. The second constitutive law (dotted-
dashed orange curve) corresponds to Eqs. (2)-(3) and is denoted by $f^{\rm
WS}_{\rm ss}(v_{\rm ss})$, see text for additional discussion. Finally, the
third constitutive law (dashed green curve) features no local minimum (see
text for discussion) and is denoted by $f^{\rm SW}_{\rm ss}(v_{\rm ss})$. All
three constitutive laws share the same rate-weakening behavior of the steady-
state frictional strength at intermediate slip velocities, which plays a major
role in this work.

##  Appendix B The elasto-frictional length $L_{\rm c}$

The calculation of the elasto-frictional length $L_{\rm c}$ associated with
the rate-weakening instability appears in previous works [2, 1, 66, 5, 8, 33,
10], so here we just very briefly highlight the structure of the calculation.
A linear stability analysis of an interface sliding at a constant velocity
$v_{0}$ is performed; that is, the starting point is a space and time
independent solution featuring a steady-state frictional stress
$\tau\\!=\\!\sigma_{0}f_{\rm ss}(v_{0})$. One then introduces small
perturbations to all fields, each assumed to be proportional to a Fourier mode
$e^{\Lambda t-ikx}$, where $\Lambda$ is the complex growth rate and $k$ is the
wavenumber. Obtaining expressions for the perturbation of the interfacial
shear stress, from the bulk elastodynamic equations corresponding to bodies of
height $H$, and for the perturbation of the frictional strength, one obtains
(by equating the two) the linear perturbation spectrum $\Lambda(k)$. The
elasto-frictional length is related to the critical wavenumber $k_{\rm c}$ for
the onset of instability according to $L_{\rm c}(H)\\!=\\!2\pi/k_{\rm c}$,
where $k_{\rm c}$ is determined from the zero crossing of the real part of
$\Lambda$, i.e. by $\Re[\Lambda(k_{\rm c},H)]\\!=\\!0$. The general structure
of $L_{\rm c}(H)$ appears in Eq. (5).

##  Appendix C The spectral boundary integral formulation

The infinite height ($H\\!\to\\!\infty$) calculations are performed using an
in-house open-source implementation (called cRacklet [81]) of the spectral
boundary integral formulation of the elastodynamic equations [58, 59, 60]. The
basic relation between the interfacial shear stress and the slip displacement
in this case is given in Eq. (3), where the Fourier representation of the
spatiotemporal integral term $s(x,t)$ for both mode-II and mode-III symmetries
can be found in [60]. In our numerical calculations, the interface is assumed
to be initially at steady-state with $v_{\rm ss}\\!=\\!v_{0}$ and
$\phi_{ss}\\!=\\!D/v_{0}$. The interface is then perturbed by adding spatial
Gaussian noise to the state variable $\phi(x,t)$, and the resulting slip
velocity is computed by combining Eq. (3) and the rate-and-state friction law
$\tau\\!=\\!\sigma_{0}\operatorname{sgn}(v)f(|v|,\phi)$. Note that
$\tau_{0}(t)$ in Eq. (3) is treated as unknown and that we impose the
constraint of Eq. (4). The slip displacement $u(x,t)$ is then integrated in
time using an explicit time-stepping scheme $u(x,t+\Delta
t)\\!=\\!u(x,t)+\tfrac{1}{2}v(x,t)\Delta{t}$, with the time step being
$\Delta{t}\\!=\\!\alpha_{{}_{\rm BIM}}\Delta{x}/c_{s}$, where $\Delta{x}$ is
the numerical grid spacing. The numerical parameter $\alpha_{{}_{\rm BIM}}$ is
chosen to ensure the stability and the convergence of the numerical scheme,
and is typically set equal to $0.1$.

The exact shape of the initial perturbation has no impact on the long-time
behavior of the system (i.e. it does not alter the properties of the emergent
steady pulse trains). The coarsening dynamics, which are not studied in detail
in this work, may depend on the initial perturbation. The bulk parameters in
our calculations have been set equal to $\mu\\!=\\!9\\!\times\\!10$9Pa,
$\nu\\!=\\!0.33$ and $\rho\\!=\\!1200$kg/m3. In addition, $\mu$ has been
varied from $2.25\\!\times\\!10$9 to $9\\!\times\\!10$10Pa in the results
appearing in Figs. 5, 6, S3.

##  Appendix D Finite element method formulation

###  D.1 The numerical method

The finite height $H$ calculations follow the system configuration sketched in
Fig. 1a. The top and bottom boundaries are loaded by a horizontal velocity
$v_{0}/2$ (in opposite directions, cf. Fig. 1a) and by a constant compressive
normal stress of magnitude $\sigma_{0}$. The bodies are initially moving
uniformly in opposite directions at a velocity $v_{0}/2$, and periodic
boundary conditions are enforced at the lateral edges, $x\\!=\\!0$ and
$x\\!=\\!W$. The interface is initially at steady-state with $v_{\rm
ss}\\!=\\!v_{0}$ and $\phi_{ss}\\!=\\!D/v_{0}$. Perturbations are introduced
by adding spatial Gaussian noise to the state variable, as in the BIM case.
The bulk parameters used for the FEM simulations and their BIM counterparts
(see, for example, Fig. 2b) are set equal to $\mu\\!=\\!3.1\\!\times\\!10$9Pa,
$\nu\\!=\\!0.33$ and $\rho\\!=\\!1200$kg/m3.

The FEM calculations are performed using the explicit dynamic finite element
framework, based on an in-house open-source finite element library called
Akantu [82]. The domain is discretized into a regular mesh composed of
bilinear quadrilateral elements (Q4). The sliding interface between the two
elastic bodies is modeled using a node to node contact algorithm, whose
details can be found in [69]. Time is integrated using the central difference
method and the time step is taken small enough to eliminate the numerical
instabilities associated with the explicit finite element modeling of rate-
and-state friction, as explained in [69]. In our simulations, we set the time
step to $\Delta{t}\\!=\\!\alpha_{{}_{\rm FEM}}\Delta{t}_{{}_{\rm CFL}}$, where
$\Delta{t}_{{}_{\rm CFL}}$ is set by the Courant-Friedrichs-Lewy condition and
$\alpha_{{}_{\rm FEM}}$ is typically taken to be ${\cal O}(0.01)$.

###  D.2 The WS friction law

As explained in the main text and in Appendix A, the most important aspect of
the interfacial constitutive relation for our results is the existence of a
rate-weakening branch of the steady-state friction curve and that the applied
velocity $v_{0}$ resides on it. To demonstrate this point, we performed
exactly the same BIM simulations except that in one case we used $f^{\rm
N}(|v|,\phi)$ of Eq. (1) and in the other $f^{\rm WS}(|v|,\phi)$ of Eq. (3),
which are mainly distinguished by the different steady-state behavior at very
low velocities, see Fig. S1. The long-time behavior of the system in the two
cases is presented in Fig. S2 (employing the same line type and color scheme
of Fig. S1). It is observed that the results are quantitatively similar, in
terms of the train periodicity ($W_{\rm p}\\!=\\!W$), single pulse shape and
propagation velocity (see figure legend), substantiating our claim.

Figure S2: Snapshots of the steady-state slip velocity field $v(x)/v_{0}$
obtained in two BIM simulations, one with the N-shaped friction law (Eq. (1),
solid brown curve, cf. Fig. S1) and the other with the WS friction law (Eq.
(3), dotted-dashed orange curve, cf. Fig. S1). The two simulations are
otherwise identical, employing $v_{0}\\!=\\!1\\!\times\\!10^{-3}\text{m/s}$
and $W\\!=\\!40\text{m}$. In both cases, the resulting pulse train periodicity
satisfies $W_{\rm p}\\!=\\!W$, the single pulse shape is quantitatively
similar (note that the two snapshots are shifted along $x$ for visual clarity)
and the propagation velocity is also quantitatively similar ($c_{\rm p}^{\rm
N}\\!\simeq\\!0.73c_{s}$ for the N-shaped law and $c_{\rm p}^{\rm
WS}\\!\simeq\\!0.70c_{s}$ for the WS law).

It turned out that this robustness of our results against changes in the low
velocity behavior of the friction law (while preserving a very similar rate-
weakening behavior) is of technical importance in relation to our FEM
calculations. In particular, we found that the FEM calculations suffered from
numerical instabilities associated with the low frictional strength featured
by the N-shaped interfacial constitutive relation at low velocities. These
numerical instabilities can be eliminated by using Eq. (3) instead of Eq. (1),
without affecting the emerging physics. Consequently, in our FEM calculations
we adopted the WS friction law of Eq. (3).

###  D.3 The space-independent stick-slip like behavior for
$H\\!\simeq\\!W\\!\gg\\!L_{\rm c}$

As explained in the main text, our FEM calculations with
$H\\!\simeq\\!W\\!\gg\\!L_{\rm c}$ featured time periods where the slip
velocity almost vanished homogeneously across the interface, followed by
nearly homogeneous large slip velocity periods. This space-independent stick-
slip like behavior appears to substantially deviate from the steady-average
velocity condition of Eq. (4), and hence is not expected to persist in the
large $W$ and long-time limits. Consequently, we slightly modified our
numerical scheme to eliminate this behavior by adding a constraint in the FEM
$H\\!\simeq\\!W\\!\gg\\!L_{\rm c}$ simulations. In particular, at each time
step we computed the average slip rate $\langle v\rangle$ of the interface and
compared it to the driving velocity $v_{0}$. We then introduced a small
spatially-homogeneous shift $\Delta{v}\\!=\\!v_{0}-\langle v\rangle$ to the
velocity field and $\Delta{a}\\!=\\!2(v_{0}-\langle v\rangle)/\Delta{t}$ to
its acceleration counterpart ($a\\!=\\!\dot{v}$ and $\Delta{a}$ is the small
shift in $a$) in computational nodes along the interface (with opposite signs
for nodes belonging to the top/bottom body, to preserve the shear symmetry).
This amounts to the application of a spatially-homogeneous external stress on
the top/bottom nodes along the interface. This stress is typically very small,
less than one percent of the initial frictional stress
$\tau_{0}\\!=\\!\sigma_{0}f_{\rm ss}(v_{0})$, yet it nevertheless ensures that
the average driving condition of Eq. (4) is nearly satisfied at all times.

##  Appendix E The effective fracture energy

The effective fracture energy $G_{\rm c}$ of interfaces obeying rate-and-state
friction law is self-selected by the interfacial dynamics. In order to test
the pulse equation of motion of Eq. (10), we computed $G_{\rm c}$ using two
different methods. The first method relies on an approximated mapping of slip-
weakening (rather than rate-weakening) friction laws — characterized by a
slip-weakening distance — to rate-and-state ones, as discussed in [40]. In
[40] it has been shown that an effective slip-weakening distance $\delta^{\rm
eff}_{\rm c}$ can be extracted in rate-and-state friction calculations, and
that it approximately follows the relation $\delta^{\rm eff}_{\rm
c}\\!\approx\\!D\ln(v_{\rm r}/v_{\rm bg})$, where $v_{\rm bg}$ is the
background slip velocity that a rupture mode propagates into and $v_{\rm r}$
is the residual slip velocity left behind it. This relation has been later
used in [24, 25] to obtain an approximate expression for the effective
fracture energy in the form $G_{\rm
c}\\!\approx\\!\tfrac{1}{2}Dbf_{0}\sigma_{0}\ln^{2}(v_{\rm r}/v_{\rm bg})$. We
used this estimate for $G_{\rm c}$, with $v_{\rm bg}\\!=\\!v_{\rm s}$ and
$v_{\rm r}\\!=\\!v_{0}$, in Fig. 5. The results strongly supported the
prediction in Eq. (10).

Figure S3: Testing the theoretical prediction of Eq. (10), using the same data
as in Fig. 5, but this time with $G_{\rm c}$ being estimated from the fit of
the singular fields near the rupture edge, which is denoted by $G^{\rm
fit}_{\rm c}$ (see text for discussion). The results reveal reasonably good
agreement with the prediction, with an ${\cal O}(1)$ pre-factor, similarly to
the results of Fig. 5.

The second method relies on the extraction of the nearly singular fields in
the vicinity of the leading edge of the propagating pulse, shown in Fig. 4 to
be reasonably well approximated by the classical square root singular fields
of LEFM. In particular, such a procedure allows to extract the stress
intensity factor $K_{\rm III}$ by fitting simulational data to Eq. (7).
Invoking then the leading edge energy balance of Eq. (8), one can obtain an
estimate of the effective fracture energy, which in this context is denoted by
$G^{\rm fit}_{\rm c}$. In Fig. S3, we present the very same data (and
combinations of physical quantities in the $x$ and $y$ axes) as in Fig. 5, but
this time using $G^{\rm fit}_{\rm c}$ for $G_{\rm c}$, estimated from fitting
the singular fields near the pulse leading edge, as just described. The
results yet again agree with the prediction in Eq. (10) reasonably well,
demonstrating the validity of the latter independently of the method used to
estimate the effective fracture energy.

##  Appendix F The SW friction law

Many studies available in the literature employ a rate-and-state friction law
that does not feature rate-strengthening at relatively high slip velocities
(and hence no local minimum as observed in Fig. 1b) [35, 38, 24, 25, 83]. Such
a behavior is exhibited by the strengthening-weakening (SW) friction law
discussed in Appendix A (which corresponds to using $f^{\rm SW}(|v|,\phi)$ in
Eq. (1)), whose steady-state curve is shown in Fig. S1 (dashed green curve).
For completeness, we performed mode-III BIM calculations ($H\\!\to\\!\infty$)
with the SW friction law. We find qualitatively similar results to those
obtained for the N-shaped law, i.e. coarsening dynamics towards a pulse train
with periodicity set by the system length, $W_{\rm p}\\!=\\!W$. A
representative movie of such calculations is available in [70], see Appendix G
for details. We do note that the SW friction law leads to sharper pulses and
typically features higher slip and propagation velocities compared to their
N-shaped counterparts (not shown).

##  Appendix G Movies

Movies are available at [70] as follows:

  * •

M1.mp4: Infinite domain with the N-shaped friction law, see Fig. 2a

  * •

M2.mp4: $H\\!<\\!L_{\rm c}<W$ with the WS friction law, corresponding to Fig.
3

  * •

M3.mp4: $H\\!<\\!L_{\rm c}\\!<\\!W$ with the WS friction law, corresponding to
Fig. 3 in the co-moving frame of the pulse.

  * •

M4.mp4: $H\\!\simeq\\!W\\!\gg\\!L_{\rm c}$ with the WS friction law.

  * •

M5.mp4: Infinite domain with the SW friction law, see Appendix F.

## References

  * Ruina [1983] A. L. Ruina, Slip instability and state variable friction laws, J. Geophys. Res. 88, 10359 (1983).
  * Rice and Ruina [1983] J. R. Rice and A. L. Ruina, Stability of Steady Frictional Slipping, J. Appl. Mech. 50, 343 (1983).
  * Marone [1998] C. Marone, Laboratoty-derived friction laws and their application to seismic faulting, Annu. Rev. Earth Planet. Sci. 26, 643 (1998).
  * Nakatani [2001] M. Nakatani, Conceptual and physical clarification of rate and state friction: Frictional sliding as a thermally activated rheology, J. Geophys. Res. Solid Earth 106, 13347 (2001).
  * Baumberger and Caroli [2006] T. Baumberger and C. Caroli, Solid friction from stick–slip down to pinning and aging, Adv. Phys. 55, 279 (2006).
  * Dieterich [2007] J. H. Dieterich, Applications of rate-and state-dependent friction to models of fault slip and earthquake occurrence, Treatise Geophys. 4, 107 (2007).
  * Nagata  _et al._ [2012] K. Nagata, M. Nakatani, and S. Yoshida, A revised rate- and state-dependent friction law obtained by constraining constitutive and evolution laws separately with laboratory data, J. Geophys. Res. Solid Earth 117, B02314 (2012).
  * Bhattacharya and Rubin [2014] P. Bhattacharya and A. M. Rubin, Frictional response to velocity steps and 1-D fault nucleation under a state evolution law with stressing-rate dependence, J. Geophys. Res. Solid Earth 119, 2272 (2014).
  * Barras  _et al._ [2019] F. Barras, M. Aldam, T. Roch, E. A. Brener, E. Bouchbinder, and J.-F. Molinari, Emergence of cracklike behavior of frictional rupture: The origin of stress drops, Phys. Rev. X 9, 041043 (2019).
  * Bar-Sinai  _et al._ [2019] Y. Bar-Sinai, M. Aldam, R. Spatschek, E. A. Brener, and E. Bouchbinder, Spatiotemporal Dynamics of Frictional Systems: The Interplay of Interfacial Friction and Bulk Elasticity, Lubricants 7, 91 (2019).
  * Barras  _et al._ [2020] F. Barras, M. Aldam, T. Roch, E. A. Brener, E. Bouchbinder, and J.-F. Molinari, The emergence of crack-like behavior of frictional rupture: Edge singularity and energy balance, Earth Planet. Sci. Lett. 531, 115978 (2020).
  * Ben-Zion [2008] Y. Ben-Zion, Collective behavior of earthquakes and faults: Continuum-discrete transitions, progressive evolutionary changes, and different dynamic regimes, Rev. Geophys. 46, RG4006 (2008).
  * Vanossi  _et al._ [2013] A. Vanossi, N. Manini, M. Urbakh, S. Zapperi, and E. Tosatti, Colloquium: Modeling friction: From nanoscale to mesoscale, Rev. Mod. Phys. 85, 529 (2013).
  * Armstrong-Hélouvry  _et al._ [1994] B. Armstrong-Hélouvry, P. Dupont, and C. C. De Wit, A survey of models, analysis tools and compensation methods for the control of machines with friction, Automatica 30, 1083 (1994).
  * Wojewoda  _et al._ [2008] J. Wojewoda, A. Stefański, M. Wiercigroch, and T. Kapitaniak, Hysteretic effects of dry friction: modelling and experimental studies, Philos. Trans. R. Soc. A Math. Phys. Eng. Sci. 366, 747 (2008).
  * Massi  _et al._ [2007] F. Massi, L. Baillet, O. Giannini, and A. Sestieri, Brake squeal: Linear and nonlinear numerical approaches, Mech. Syst. Signal Process. 21, 2374 (2007).
  * Tonazzi  _et al._ [2013] D. Tonazzi, F. Massi, A. Culla, L. Baillet, A. Fregolent, and Y. Berthier, Instability scenarios between elastic media under frictional contact, Mech. Syst. Signal Process. 40, 754 (2013).
  * Yamashita and Ohnaka [1991] T. Yamashita and M. Ohnaka, Nucleation process of unstable rupture in the brittle regime: A theoretical approach based on experimentally inferred relations, J. Geophys. Res. Solid Earth 96, 8351 (1991).
  * Ben-Zion and Rice [1997] Y. Ben-Zion and J. R. Rice, Dynamic simulations of slip on a smooth fault in an elastic solid, J. Geophys. Res. Solid Earth 102, 17771 (1997).
  * Ben-Zion [2001] Y. Ben-Zion, Dynamic ruptures in recent models of earthquake faults, J. Mech. Phys. Solids 49, 2209 (2001).
  * Ampuero  _et al._ [2002] J.-P. Ampuero, J.-P. Vilotte, and F. J. Sánchez-Sesma, Nucleation of rupture under slip dependent friction law: Simple models of fault zone, J. Geophys. Res. Solid Earth 107, 2324 (2002).
  * Lapusta and Rice [2003] N. Lapusta and J. R. Rice, Nucleation and early seismic propagation of small and large events in a crustal earthquake model, J. Geophys. Res. Solid Earth 108, 2205 (2003).
  * Uenishi and Rice [2003] K. Uenishi and J. R. Rice, Universal nucleation length for slip-weakening rupture instability under nonuniform fault loading, J. Geophys. Res. Solid Earth 108, 2042 (2003).
  * Rubin and Ampuero [2005] A. M. Rubin and J.-P. Ampuero, Earthquake nucleation on (aging) rate and state faults, J. Geophys. Res. Solid Earth 110, B11312 (2005).
  * Ampuero and Rubin [2008] J.-P. Ampuero and A. M. Rubin, Earthquake nucleation on rate and state faults – Aging and slip laws, J. Geophys. Res. Solid Earth 113, B01302 (2008).
  * Kaneko and Lapusta [2008] Y. Kaneko and N. Lapusta, Variability of earthquake nucleation in continuum models of rate-and-state faults and implications for aftershock rates, J. Geophys. Res. Solid Earth 113, B12312 (2008).
  * McLaskey and Kilgore [2013] G. C. McLaskey and B. D. Kilgore, Foreshocks during the nucleation of stick-slip instability, J. Geophys. Res. Solid Earth 118, 2982 (2013).
  * Latour  _et al._ [2013] S. Latour, A. Schubnel, S. B. Nielsen, R. Madariaga, and S. Vinciguerra, Characterization of nucleation during laboratory earthquakes, Geophys. Res. Lett. 40, 5064 (2013).
  * Viesca [2016a] R. C. Viesca, Stable and unstable development of an interfacial sliding instability, Phys. Rev. E 93, 060202 (2016a).
  * Viesca [2016b] R. C. Viesca, Self-similar slip instability on interfaces with rate- and state-dependent friction, Proc. R. Soc. A Math. Phys. Eng. Sci. 472, 20160254 (2016b).
  * Kaneko  _et al._ [2016] Y. Kaneko, S. B. Nielsen, and B. M. Carpenter, The onset of laboratory earthquakes explained by nucleating rupture on a rate-and-state fault, J. Geophys. Res. Solid Earth 121, 6071 (2016).
  * Kaneko  _et al._ [2017] Y. Kaneko, B. M. Carpenter, and S. B. Nielsen, Nucleation process of magnitude 2 repeating earthquakes on the San Andreas Fault predicted by rate-and-state fault models with SAFOD drill core data, Geophys. Res. Lett. 44, 162 (2017).
  * Aldam  _et al._ [2017] M. Aldam, M. Weikamp, R. Spatschek, E. A. Brener, and E. Bouchbinder, Critical Nucleation Length for Accelerating Frictional Slip, Geophys. Res. Lett. 44, 11,390 (2017).
  * Gabriel  _et al._ [2012] A.-A. Gabriel, J.-P. Ampuero, L. A. Dalguer, and P. M. Mai, The transition of dynamic rupture styles in elastic media under velocity-weakening friction, J. Geophys. Res. Solid Earth 117, B09311 (2012).
  * Perrin  _et al._ [1995] G. Perrin, J. R. Rice, and G. Zheng, Self-healing slip pulse on a frictional surface, J. Mech. Phys. Solids 43, 1461 (1995).
  * Beeler and Tullis [1996] N. M. Beeler and T. E. Tullis, Self-healing slip pulses in dynamic rupture models due to velocity-dependent strength, Bull. Seismol. Soc. Am. 86, 1130 (1996).
  * Cochard and Madariaga [1996] A. Cochard and R. Madariaga, Complexity of seismicity due to highly rate-dependent friction, J. Geophys. Res. Solid Earth 101, 25321 (1996).
  * Zheng and Rice [1998] G. Zheng and J. R. Rice, Conditions under which velocity-weakening friction allows a self-healing versus a cracklike mode of rupture, Bull. Seismol. Soc. Am. 88, 1466 (1998).
  * Nielsen  _et al._ [2000] S. B. Nielsen, J. M. Carlson, and K. B. Olsen, Influence of friction and fault geometry on earthquake rupture, J. Geophys. Res. Solid Earth 105, 6069 (2000).
  * Bizzarri and Cocco [2003] A. Bizzarri and M. Cocco, Slip-weakening behavior during the propagation of dynamic ruptures obeying rate- and state-dependent friction laws, J. Geophys. Res. 108, 2373 (2003).
  * Brener  _et al._ [2005] E. A. Brener, S. V. Malinin, and V. I. Marchenko, Fracture and friction: Stick-slip motion, Eur. Phys. J. E 17, 101 (2005).
  * Rubin and Ampuero [2009] A. M. Rubin and J.-P. Ampuero, Self-similar slip pulses during rate-and-state earthquake nucleation, J. Geophys. Res. 114, B11305 (2009).
  * Ben-David  _et al._ [2010] O. Ben-David, G. Cohen, and J. Fineberg, The dynamics of the onset of frictional slip, Science 330, 211 (2010).
  * Bar-Sinai  _et al._ [2013] Y. Bar-Sinai, R. Spatschek, E. A. Brener, and E. Bouchbinder, Instabilities at frictional interfaces: Creep patches, nucleation, and rupture fronts, Phys. Rev. E 88, 060403 (2013).
  * Svetlizky and Fineberg [2014] I. Svetlizky and J. Fineberg, Classical shear cracks drive the onset of dry frictional motion, Nature 509, 205 (2014).
  * Putelat  _et al._ [2017] T. Putelat, J. H. Dawes, and A. R. Champneys, A phase-plane analysis of localized frictional waves, Proc. R. Soc. A Math. Phys. Eng. Sci. 473, 20160606 (2017).
  * Heaton [1990] T. H. Heaton, Evidence for and implications of self-healing pulses of slip in earthquake rupture, Phys. Earth Planet. Inter. 64, 1 (1990).
  * Nielsen and Madariaga [2003] S. Nielsen and R. Madariaga, On the self-healing fracture mode, Bull. Seismol. Soc. Am. 93, 2375 (2003).
  * Brener  _et al._ [2018] E. A. Brener, M. Aldam, F. Barras, J.-F. Molinari, and E. Bouchbinder, Unstable Slip Pulses and Earthquake Nucleation as a Nonequilibrium First-Order Phase Transition, Phys. Rev. Lett. 121, 234302 (2018).
  * Brantut  _et al._ [2019] N. Brantut, D. I. Garagash, and H. Noda, Stability of pulse-like earthquake ruptures, J. Geophys. Res. Solid Earth 124, 8998 (2019).
  * Heimisson  _et al._ [2019] E. R. Heimisson, E. M. Dunham, and M. Almquist, Poroelastic effects destabilize mildly rate-strengthening friction to generate stable slow slip pulses, J. Mech. Phys. Solids 130, 262 (2019).
  * Somerville  _et al._ [1999] P. Somerville, K. Irikura, R. Graves, S. Sawada, D. Wald, N. Abrahamson, Y. Iwasaki, T. Kagawa, N. Smith, and A. Kowada, Characterizing Crustal Earthquake Slip Models for the Prediction of Strong Ground Motion, Seismol. Res. Lett. 70, 59 (1999).
  * Lu  _et al._ [2007] X. Lu, N. Lapusta, and A. J. Rosakis, Pulse-like and crack-like ruptures in experiments mimicking crustal earthquakes, Proc. Natl. Acad. Sci. 104, 18931 (2007).
  * Shlomai and Fineberg [2016] H. Shlomai and J. Fineberg, The structure of slip-pulses and supershear ruptures driving slip in bimaterial friction, Nat. Commun. 7, 11787 (2016).
  * Melgar and Hayes [2017] D. Melgar and G. P. Hayes, Systematic Observations of the Slip Pulse Properties of Large Earthquake Ruptures, Geophys. Res. Lett. 44, 9691 (2017).
  * Bar-Sinai  _et al._ [2014] Y. Bar-Sinai, R. Spatschek, E. A. Brener, and E. Bouchbinder, On the velocity-strengthening behavior of dry friction, J. Geophys. Res. Solid Earth 119, 1738 (2014).
  * Baus and Tejero [2007] M. Baus and C. F. Tejero,  _Equilibrium statistical physics: phases of matter and phase transitions_ (Springer Science & Business Media, 2007).
  * Geubelle and Rice [1995] P. H. Geubelle and J. R. Rice, A spectral method for three-dimensional elastodynamic fracture problems, J. Mech. Phys. Solids 43, 1791 (1995).
  * Morrissey and Geubelle [1997] J. W. Morrissey and P. H. Geubelle, A numerical scheme for mode III dynamic fracture problems, Int. J. Numer. Methods Eng. 40, 1181 (1997).
  * Breitenfeld and Geubelle [1998] M. S. Breitenfeld and P. H. Geubelle, Numerical analysis of dynamic debonding under 2D in-plane and 3D loading, Int. J. Fract. 93, 13 (1998).
  * Freund [1998] L. B. Freund,  _Dynamic Fracture Mechanics_ (Cambridge university press, Cambridge, 1998).
  * Kostrov [1964] B. V. Kostrov, Self-similar problems of propagation of shear cracks, J. Appl. Math. Mech 28, 1077 (1964).
  * Freund [1979] L. Freund, The mechanics of dynamic shear crack propagation, J. Geophys. Res. Solid Earth 84, 2199 (1979).
  * Broberg [1999] K. B. Broberg,  _Cracks and fracture_ (Elsevier, 1999).
  * Landau and Lifshitz [1986] L. Landau and E. Lifshitz,  _Theory of Elasticity, Third Edition: Volume 7 (Course of Theoretical Physics)_ (Butterworth-Heinemann, 1986).
  * Lapusta  _et al._ [2000] N. Lapusta, J. R. Rice, Y. Ben-Zion, and G. Zheng, Elastodynamic analysis for slow tectonic loading with spontaneous rupture episodes on faults with rate- and state-dependent friction, J. Geophys. Res. Solid Earth 105, 23765 (2000).
  * Ben-Zion and Rice [1995] Y. Ben-Zion and J. R. Rice, Slip patterns and earthquake populations along different classes of faults in elastic solids, J. Geophys. Res. Solid Earth 100, 12959 (1995).
  * Crupi and Bizzarri [2013] P. Crupi and A. Bizzarri, The role of radiation damping in the modeling of repeated earthquake events, Ann. Geophys. 56, R0111 (2013).
  * Rezakhani  _et al._ [2020] R. Rezakhani, F. Barras, M. Brun, and J.-F. Molinari, Finite element modeling of dynamic frictional rupture with rate and state friction, J. Mech. Phys. Solids 141, 103967 (2020).
  * [70] See [url will be inserted by publisher] for movies of dynamic simulations conducted using BIM and FEM, as detailed in the manuscript. The movies are available here: https://www.youtube.com/watch?v=UVD4VeJHnRM&list=PLT7c4IN61XQM8dk-hq9ZWtwGCvhT54nt7
  * Brener  _et al._ [2016] E. A. Brener, M. Weikamp, R. Spatschek, Y. Bar-Sinai, and E. Bouchbinder, Dynamic instabilities of frictional sliding at a bimaterial interface, J. Mech. Phys. Solids 89, 149 (2016).
  * Brener and Bouchbinder [2020] E. A. Brener and E. Bouchbinder, Unconventional singularities and energy balance in frictional rupture (2020), arXiv:2008.04697 [cond-mat.soft] .
  * Brener and Bouchbinder [2021] E. A. Brener and E. Bouchbinder, Theory of unconventional singularities of frictional shear cracks (2021), arXiv:2103.10826 [cond-mat.soft] .
  * Kanamori and Heaton [2000] H. Kanamori and T. H. Heaton, Microscopic and macroscopic physics of earthquakes, Geophys. Monogr. Am. Geophys. Union 120, 147 (2000).
  * Galeano  _et al._ [2000] J. Galeano, P. Espanol, and M. Rubio, Experimental and theoretical results of stress relaxations in a model of earthquake dynamics, EPL (Europhysics Letters) 49, 410 (2000).
  * Ma  _et al._ [2014] S. Ma, T. Shimamoto, L. Yao, T. Togo, and H. Kitajima, A rotary-shear low to high-velocity friction apparatus in beijing to study rock friction at plate to seismic slip rates, Earthq. Sci. 27, 469 (2014).
  * Rundle  _et al._ [2003] J. B. Rundle, D. L. Turcotte, R. Shcherbakov, W. Klein, and C. Sammis, Statistical physics approach to understanding the multiscale dynamics of earthquake fault systems, Rev. of Geophys. 41 (2003).
  * Kawamura  _et al._ [2012] H. Kawamura, T. Hatano, N. Kato, S. Biswas, and B. K. Chakrabarti, Statistical physics of fracture, friction, and earthquakes, Rev. Mod. Phys. 84, 839 (2012).
  * Bar-Sinai  _et al._ [2012] Y. Bar Sinai, E. A. Brener, and E. Bouchbinder, Slow rupture of frictional interfaces, Geophys. Res. Lett. 39, L03308 (2012).
  * Bar-Sinai  _et al._ [2015] Y. Bar-Sinai, R. Spatschek, E. a. Brener, and E. Bouchbinder, Velocity-strengthening friction significantly affects interfacial dynamics, strength and dissipation, Sci. Rep. 5, 7841 (2015).
  * Roch  _et al._ [2021] T. Roch, F. Barras, P. H. Geubelle, and J.-F. Molinari, cRacklet simulations for ”Velocity-driven frictional-sliding, Coarsening and steady-state pulse trains” (2021).
  * Richart and Molinari [2015] N. Richart and J. F. Molinari, Implementation of a parallel finite-element library: Test case on a non-local continuum damage model, Finite Elem. Anal. Des. 100, 41 (2015).
  * Nielsen [2017] S. B. Nielsen, From slow to fast faulting: recent challenges in earthquake fault mechanics, Philos. Trans. R. Soc. A Math. Phys. Eng. Sci. 375, 20160016 (2017).