format: text

text:
# Depinning with dynamic stress overshoots: A hybrid of critical and pseudohysteretic behavior

## I Introduction

Extended elastic manifolds pulled through a quenched random medium by an applied force $`F`$ exhibit, in the absence of thermal fluctuations, a sharp transition from a pinned phase to a moving phase as $`F`$ is increased through a critical value $`F_c`$ . Examples include interfaces between two fluids in porous media  or between oppositely magnetized ferromagnetic domains, vortex lines and lattices in type II superconductors , , charge density waves , and planar crack fronts in solids .

Although the depinning transitions of interest are driven non-equilibrium transitions, it is instructive to draw an analogy with equilibrium phase transitions with the average velocity $`\overline{v}`$ playing the role of an order parameter, $`F`$ a tuning parameter, and the quenched variations of the random potential loosely analogous — at least as giving rise to an ensemble — to thermal fluctuations. The character of depinning transitions can, one might expect, be either discontinuous transitions with hysteresis — loosely like first-order transitions —, or critical — analogous to second-order — transitions depending on the system and, perhaps, on its history; such history dependence is an effect that cannot occur in true equilibrium. Theoretical analysis has shown that a broad class of realistic models undergo a critical depinning transition with a unique, history independent critical force in the limit of a large system and

$$\overline{v}(FF_c)^\beta $$
(1)
for $`F`$ just above $`F_c`$. Several different universality classes have been studied, including both short and long-range interactions, and random forces with or without periodicity — the former arising for manifolds with a periodic structure in the direction in which they move ,, , , .

But most of the theoretical analysis has focused on dissipative dynamics for which both inertia and any wave or other non-local stress propagation effects are ignored. The purpose of this paper is to study some of the consequences of these and other effects which we shall generically refer to, for reasons to be explained shortly, as stress overshoot effects.

In order to understand the potential of these effects to significantly change the nature of depinning, it is necessary to consider the nature of the irregular local motion that underlies the critical depinning phenomena. The elasticity of the manifold mediates between two competing types of forces: the applied driving force and the random local pinning forces. For small $`F`$, the pinning dominates and the system relaxes to one of many static locally-stable configurations. But as the force is slowly increased, there will be local instabilities when the driving force exceeds the random pinning in some small region. A segment of the manifold will then move forward rapidly and there will be some transient motion, limited in spatial extent, until a new static configuration is reached. For large $`F`$, in contrast, the applied force will dominate and the system will approach a nonequilibrium statistically steady state with a non-zero mean velocity. Nevertheless, especially if the force is not too large, the motion on short length and time scales will be very irregular with instantaneous local velocities that far exceed the macroscopic average velocity $`\overline{v}`$.

In the absence of inertia or wave propagation, each segment of the manifold will move in response to the total force applied to it: from the applied drive, from the random pinning, and from the other segments via elasticity. As long as the applied force is non-decreasing in time, the motion (at least after initial transients have decayed away) will be only in the “forward” direction in which the system is driven. This, combined with the convexity properties of elasticity, means that the static configuration the system will settle into after it is disturbed by an increase — global or local — in the applied force does not depend on the details of its dynamics. The motion can thus be considered as quasistatic, in spite of the rapid local motion that occurs.

Now consider what can happen in the presence of either inertia or elastic waves that carry stress from one region to another in response to motion of one segment. The local dynamics would then appear to be crucial: If the local motion is rapid enough that the relaxation to a new static configuration is underdamped, a moving segment can overshoot one or more potential static configurations before settling, if at all, in another. Even if the inertia is small enough that such local prolonged jumps do not occur, as long as the motion is sometimes underdamped, a segment can temporarily overshoot a static configuration before relaxing back into it; this will produce a temporary overshoot of the stress — above its eventual static value — that this motion induces on neighboring segments. Any arbitrary small overshoot in the stress has the potential to dislodge another segment if there is one nearby that was sufficiently close to being destabilized in the absence of the overshoot; again, the effects of this will be to cause the system to skip through a potentially static configuration without stopping. Elastic waves, just like their electrodynamic cousins, carry with them pulses in stress that are larger than the eventual static stress that will obtain long after the waves have passed by. These stress overshoots, like those from the inertia of local motion, have the potential to cause overjumping.

Very generally, overjumping of any kind means that which configuration a pinned manifold stops in depends on details of its local dynamics and on its history, in a way that cannot occur in the absence of inertial effects. One particularly interesting consequence of this is the coexistence in two identical samples at the same value of the driving force of a static locally stable configuration, and a moving configuration that will “overtake” the static configuration. What the macroscopic consequences of this are is the primary subject of this paper. As stress overshoots can occur more readily than local overjumps caused by inertia, we will generally refer to both types of effects as stress overshoots, although, with slight inconsistency that we trust will not be confusing, we will characterize their strength by a parameter that we denote $`M`$.

What is the nature of the depinning transition as $`M`$ is increased from zero? There are various scenarios one can readily envisage. For large enough $`M`$, in the model we introduce in the next section, an infinitesimal increase in $`F`$ from any pinned state will result in a non-zero average velocity at long times. This is because sufficiently large stress overshoots always induce other segments to move when triggered by an initial segment that moves in response to an increase in $`F`$. The increased triggering will spawn further motion, despite the fact that the stress overshoot is only temporary. This process will run away; and the manifold will acquire a non-zero average velocity. While it might thus seem likely that the transition will become “first-order” for large enough $`M`$, does it do so for arbitrary small $`M`$? If there are regimes in which the depinning is indeed discontinuous in some way, is there macroscopic hysteresis? In what sense?

More generally, what happens to the depinning transition beyond the dissipative limit? If it remains critical — at least in some respects — for a range of $`M`$, what is its nature? Can the quasistatic behavior persist macroscopically for small $`M`$ in spite of the presence of additional microscopic hysteresis? If so, what is the size dependence of the hysteresis and related phenomena? The system size dependence is particularly relevant for geological faults for which the statistics of the earthquakes are affected both by the nature of the drive and the distribution of the “sizes” of faults.

Several recent papers have undertaken some preliminary studies of the effects of stress overshoots in models of depinning of elastic manifolds. Reference  studied several one-dimensional models with long-range elasticity and stress overshoots motivated by planar crack fronts driven by applied loads. Reference  introduced a particularly simple model with short-range static elasticity and analyzed its infinite-range limit in which all segments of the manifold are equally coupled to each other. In this limit, the spatial properties of the manifold are averaged away; only time dependence remains and mean field theory becomes exact. Such mean field models were the starting point for theoretical understanding of the finite-dimensional physics in the quasistatic limit . Whether or not they provide a useful starting point beyond the dissipative limit, is one of the questions that we must address.

In this paper, we investigate numerically and phenomenologically the finite-range version of the stress overshoot model introduced in .

### A Outline

The remainder of this paper is organized as follows: In the next section we introduce the basic lattice model on which we focus. In Section III the general scaling picture is introduced and known results for the dissipative case are summarized. In Section IV the critical behavior in the moving phase is studied, initially for the dissipative case, and then in the presence of stress overshoots. We summarize a variety of evidence that the critical behavior which occurs as the driving force is decreased until the system stops is in the same universality class as the dissipative case.

Section V turns to the key aspect of overshoots: hysteresis. We analyze the hysteresis loops that occur when the system is stopped from the moving phase and restarted by gradually increasing the force. Various puzzling aspects of the data are discussed and some understanding of the hysteresis loops in terms of a low density of segments that can be readily triggered by an increase in $`F`$ is reached. In the following section, studies of the dynamics and statistics of avalanches that occur as the force is gradually increased are presented. These again lead to puzzling dependence on overshoot and pinning strengths although some aspects of the avalanches appear very similar to those in the dissipative limit. The data suggest various subtle crossovers may be occurring. In Section VII the dynamics of the nucleation of restarting after a system has been stopped from the moving phase are analyzed. It is found that over a substantial range of sizes, bubble-like nucleation can occur.

In Section VIII the puzzling aspects of the various sets of data are tentatively resolved in terms of a crossover as a function of length scale and system size that manifests itself in different ways as the various parameters are varied. Finally, in Section IX, the conclusions are summarized and applications to the dynamics of earthquakes are discussed briefly,

In the main body of this paper, we restrict consideration to weak enough overshoots that they do not totally change the local dynamics. But for sufficiently large $`M`$, the overshoots cause dramatic changes in the macroscopic behavior. Although these are interesting, they are probably peculiar to certain aspects of the model; results on these will be presented elsewhere .

## II Model

Near the depinning transition, the dynamics is very jerky with segments of the manifold spending most of their time stationary or almost so, but occasionally getting unpinned by the forces from other segments and moving forward only to get pinned again by a combination of the newly explored random forces and the elasticity. The inherent discreteness of these local jumps suggests that we model the manifold as a large number of segments that can jump discontinuously from one pinning position to another; this is also convenient for numerical studies. We define $`h(𝐱,t)`$ to be the single-valued scalar displacement of the manifold from some undeformed reference configuration with both the position, $`𝐱`$ and the time, $`t`$, taken to be discrete. Note that by constraining the displacement field to be single-valued, we exclude “overhangs” as well as defects such as dislocations that could otherwise occur in periodic systems. The forces on a segment of the manifold consist of three terms: the applied force $`F`$, a static random pinning force $`\eta (𝐱,h(𝐱,𝐭))`$, and the stress caused by the elasticity $`\sigma (𝐱,t)`$.

The stress depends linearly on the displacements of other parts of the manifold via

$$\sigma (𝐱,t)=\underset{𝐲}{}\underset{\tau >0}{}J_{\mathrm{𝐱𝐲}}(\tau )h(𝐲,t\tau )\widehat{J}h(𝐱,t),$$
(2)
where

$$\widehat{J}=\underset{𝐲,\tau }{}J_{\mathrm{𝐱𝐲}}(\tau )$$
(3)
and the sum is over nearest neighbors of $`𝐱`$. To model stress overshoots, we assume the simplest possible form: that the overshoot only applies to neighbors and only lasts for one time step so that

$$J_{\mathrm{𝐱𝐲}}(\tau =0)=\frac{1}{Z}(1+M)$$
(4)
and

$$J_{\mathrm{𝐱𝐲}}(1)=\frac{M}{Z},$$
(5)
with $`Z`$ the number of nearest neighbors. With this stress transfer, the jump of any nearest neighbor of $`𝐱`$ induces an extra temporary stress on the $`𝐱`$th segment. When $`M=0`$, the stress transferred to the $`𝐱`$th segment is simply proportional to the static curvature at $`𝐱`$. This stress will not decrease with time as long as the $`𝐱`$th segment does not move and the other segments only move forward or remain at rest; this limit is thus the dissipative dynamics already studied extensively. However, for positive $`M`$ the stress on the $`𝐱`$th segment caused by a jump forward of one of its neighboring segments will first increase by a larger amount than in the absence of $`M`$ and then decrease at the next time step to reach its quasistatic value.

Modeling of the local pinning forces also involves substantial arbitrariness. In Ref.  we chose randomly spaced pinning positions for each segment with uniform pinning strengths, but this choice did not affect substantially the mean field behavior . For our present purposes, it is more convenient to choose the pinning positions for each segment to be uniformly spaced but with their yield strengths, the maximum force they can sustain, randomly distributed. In particular, we take the distribution of these yield strengths to be uniformly distributed from $`[0,f_p^{max}]`$.

Because these random forces pin, or hold back, the manifold the corresponding forces $`\eta (𝐱,h(𝐱))`$, take on negative values in the range $`(f_p(𝐱,h(𝐱)),0)`$. With this form of the pinning forces, the equation of motion is simply given by

$$h(𝐱,t+1)=h(𝐱,t)+\mathrm{\Theta }[\sigma (𝐱,t)+Ff_p(𝐱,h(𝐱,t))],$$
(6)
where $`\mathrm{\Theta }`$ is the unit step function. The theta function is imposed so that a segment can move only forward and does so when the net force $`f(𝐱)`$ on it (the argument of the theta function), is positive; otherwise it remains stationary. With this dynamics, when a segment jumps its displacement always increases by one. As long as $`\overline{v}1`$, the (artificial) upper limit on the velocity, this automaton dynamics mimics the continuous time motion reasonably well.

Note that in the absence of elasticity, with $`M=0`$, $`F_c=f_p^{max}`$ with each segment becoming stuck on an anomalously strong pinning site. In the presence of the elasticity, not all of the segments can be simultaneously pinned on strong pinning sites, and the critical force will decrease. For weak pinning forces, however, the definition of the applied force in this model is somewhat pathological: when released from a pinning segment, a segment can jump forward far enough that the total force on it becomes negative and it can then pull forward other segments resulting in overall motion even if $`F`$ is negative. To make it more realistic, one could replace $`F`$ by $`F+1`$; so that there is always enough force to make the forces at pinning segments non-negative. With non-zero $`M`$, another adjustment should really be made as more realistic forms of stress overshoot involve a concomitant negative force on the segment that has moved. Since the zero of the applied force is entirely a convention, we will not make these adjustments; this will mean, however, that in some regimes the critical “force” will be negative.

Several additional aspects of the model need to be specified: the initial conditions and the order of the updating. To avoid lock-step or other “faceting” like behavior , we choose the pinning positions of each segment to be offset from one another by random amounts in the interval $`[0,1]`$. Because of the integer character of the jumps, the fractional part of the displacement of a given segment does not evolve with time. It might be thought that the displacement-independent randomness induced by this constraint could dominate over the randomness of interest, especially as far as determining the variations of critical forces, etc. in finite size systems. \[Indeed, just such an effect does occur for systems with periodic randomness such as charge density waves . But in the present case it can be shown that the additional randomness is analogous to a spatially random force that is the derivative of a random function. This, combined with the statistical tilt symmetry, mean that its effects are subdominant for large systems (although they could give rise to additional corrections to scaling.) The fact that the variations in the critical forces in finite size systems decrease substantially faster than the inverse square root of their area supports this assertion.

The updating of the displacements are done in parallel after computing all of the stresses. While there are alternate sequential methods of updating, this parallel method requires the least amount of computation and does not appear to introduce any troublesome artifacts in the regime of smallish $`M`$ of primary interest here.

Finally, to limit boundary effects, we impose periodic boundary conditions on $`h(𝐱,t)`$. This is especially important as we will use finite size scaling to analyze much of the data; this is substantially more straightforward with periodic boundary conditions. Our simulations are restricted to two dimensions which we chose because of the availability of the widest range of system sizes without running into the complications associated with very large stresses that arise in one dimensional depinning with short range interactions. We study systems of size $`L\times L`$ up to $`256\times 256`$ with most of our “large” system data on $`128\times 128`$ samples.

## III Scaling and Dissipative Dynamics

Before presenting new results for the systems of interest with stress overshoots, we briefly summarize the scaling behavior that obtains near the critical force in the absence of stress overshoots; i.e. for $`M=0`$. As the force is adiabatically increased from zero, local instabilities lead to a succession of avalanches, most of which will be small, but which can occasionally become large as the unique critical force $`F_{c0}`$ is approached. Above $`F_{c0}`$ the mean velocity in the statistical steady state rises continuously with an exponent $`\beta `$. The motion is jerky out to length scales of order the velocity correlation length $`\xi _v`$ which diverges at the critical force as

$$\xi _v\frac{1}{(FF_{c0})^\nu }.$$
(7)
The characteristic time for relaxation on scales of order $`\xi _v`$ is

$$\tau \xi _v^z$$
(8)
and in this time the manifold typically moves forward by an amount

$$\mathrm{\Delta }h\xi _v^\zeta .$$
(9)
These three exponents, $`\nu `$, $`z`$ and $`\zeta `$ characterize the scaling behavior near the transition. The velocity exponent is related to these via the observation that the mean velocity is of order the characteristic displacement per characteristic time so that

$$\overline{v}\xi _v^{\zeta z}(FF_{c0})^\beta $$
(10)
with

$$\beta =(z\zeta )\nu .$$
(11)
In the pinned phase, the critical behavior as the force is adiabatically increased can, in the absence of stress overshoots, be related to that in the moving phase reviewed above. In particular, the scaling of the dynamics and shape of the avalanches, the probability that they will be large, the divergence as $`F_{c0}`$ is approached of the cutoff size in their distribution, and the “roughness” of the manifold at the critical point are all given in terms of the same three exponents. In Section VI we will discuss the avalanches in detail, but for now we focus on the macroscopic behavior such as the velocity in the moving phase and the mean displacements in the pinned phase.

The mean displacement in response to a spatially varying applied force yields, via a statistical symmetry of the system, a scaling law that relates two of the exponents. This relation can be derived from the average static polarizability

$$\chi (q,\omega =0)\frac{\overline{h(q,\omega =0)}}{ϵ(q)}$$
(12)
to a perturbing force $`ϵ(q)\mathrm{cos}(𝐪𝐱)`$. A change of variables to $`h(𝐱,t)=h^{}(𝐱,t)^2ϵ(x)`$, yields an equation of motion for $`h^{}`$ that is statistically identical to the original one for $`h`$, independent of the perturbing force. Therefore

$$\frac{\overline{h(q,\omega =0)}}{ϵ(q)}=\frac{1}{1J(q)}\frac{1}{q^2}.$$
(13)
Since the polarizability should scale as $`\mathrm{\Delta }h/\mathrm{\Delta }F\xi ^{\zeta +1/\nu }X(q\xi )`$, with $`X`$ a scaling function, this yields the scaling law

$$\zeta +\frac{1}{\nu }=2.$$
(14)
### A Dissipative exponents

The critical exponents for $`M=0`$ take simple mean field values of $`\nu =1/2`$, $`z=2`$, characteristic of diffusive dynamics, and $`\zeta =0`$ above the critical dimension of $`d_c=4`$ for short range elasticity. The velocity scales with $`\beta =1`$ .

Below four dimensions, renormalization group expansions have been performed that justify the scaling laws and claims of universality as well as yielding results for the exponents as expansions in powers of

$$ϵ4d:$$
(15)
$$\zeta \frac{ϵ}{3}$$
(16)
and

$$z2\frac{2}{9}ϵ.$$
(17)
These yield

$$\beta 1\frac{1}{9}ϵ.$$
(18)
Recently Chauve et al  have computed these exponents to second-order in $`ϵ`$, obtaining

$$\zeta \frac{ϵ}{3}(1+0.14331ϵ)$$
(19)
and

$$z2\frac{2}{9}ϵ0.04321ϵ^2,$$
(20)
although there are some doubts about the validity of these second-order results .

## IV Critical Behavior in Moving Phase

We will shortly turn to presentation of our numerical results for the critical behavior in the moving phase. But first, it is instructive to summarize the behavior found in the mean field limit for small $`M`$ and to consider several possible scenarios that might obtain in short-range systems. We can then determine which scenario is most consistent with the data.

### A Mean field limit and scaling scenarios

In the mean field limit all of the sites are coupled to all of the others so that the number of “nearest neighbors” $`Z`$ is equal to the number of segments, $`N`$ (more precisely, $`Z=N1`$). In this limit, the critical force is found to be unchanged for $`M`$ less than a critical value, $`M_c=1`$. The velocity versus force curve is modified for any non-zero $`M`$, however, but the exponent $`\beta `$ remains at its quasistatic value of $`\beta =1`$ for $`M<M_c`$. The other universal properties of the transition are also unchanged for small $`M`$, including the lack of hysteresis in steady state and the asymptotics of the distribution of large avalanches as the critical force is approached from below.

The simplest scenario for the short range systems of interest would be like that of the mean field limit: unchanged critical behavior and no macroscopic hysteresis for small $`M`$. But previous work has shown that this cannot be the case: As shown in references , , any stress overshoot will cause the critical force to be shifted downwards, in the sense that dynamic behavior can persist for some ($`M`$-dependent) range of forces below the quasistatic critical force $`F_{c0}`$. This implies that some form of macroscopic hysteresis can exist since locally stable — at least linearly stable — static configurations exist up to $`F_{c0}`$. But whether such configurations are non-linearly stable to, for example, an arbitrarily small increase in $`F`$, is a question of substantial importance to which we will return later. For now, we focus on the moving states and how they stop as the force is lowered.

The simplest scenario that cannot immediately be ruled out is a modified version of the mean field scenario: a velocity versus force curve with a well defined critical force, $`F_c^{}`$, that is non-hysteretic as the force is decreased; history independent steady states above $`F_c^{}`$; and critical exponents, $`\beta ,\nu ^{},z^{}`$ unchanged from their quasistatic values. In renormalization group language, this would correspond to $`M`$ being an irrelevant perturbation at least as far as behavior in the moving phase.

If this mean-field-like scenario indeed applies for small $`M`$ to the finite-range model, we expect the following scaling behavior for both $`M`$ and the proximity to the $`M`$-dependent critical force,

$$fFF_c^{}(M).$$
(21)
small:

$$\overline{v}(M,F)f^\beta B(\frac{M}{f^\varphi })$$
(22)
with the crossover exponent $`\varphi <0`$ indicating the irrelevance of $`M`$, and $`B`$ a scaling function. In analogous situations in equilibrium statistical mechanics, if a parameter such as $`M`$ does not change the nature of the transition, the effect of $`M`$ on $`F_c^{}`$ can be taken into account perturbatively. Because of the singular nature of the critical fixed point that describes the quasistatic depinning — resulting, in part, from the absence of thermal fluctuations but, more essentially, from the jerky nature of the motion — the “analytic parts” might not be smooth functions of $`M`$, but under the assumption that they are, we expect that

$$F_c^{}(M)=F_{c0}aMbM^2+\mathrm{}.$$
(23)
Nevertheless, there can be singular corrections to scaling determined by the form of $`B`$ and the crossover exponent $`\varphi `$.

This scenario, in which $`M`$ is irrelevant when it is small, we call the dissipative scenario. We note, however, that this scenario is compatible with a change in behavior at a critical value of $`M`$ as occurs in mean field theory. This would give rise, for $`M`$ close to its critical value, to a crossover to some kind of multicritical behavior emerging at larger velocities.

More interesting behavior would occur if the qualitative mean field results on the effects of small $`M`$ do not simply carry over to the finite-range case. This would be the case if $`M`$ is a relevant perturbation and would correspond to a crossover scaling function, such as $`B`$ in Eq. (22) with $`\varphi `$ positive. A relevant perturbation $`M`$ would yield a singular correction to the critical force of the form

$$[F_{c0}F_c^{}(M)]_{\mathrm{sing}}M^{1/\varphi },$$
(24)
which would dominate over the leading (or subdominant) analytic shift if $`\varphi >1`$ (or $`\varphi >\frac{1}{2}`$). Earlier numerical results by Ramanathan and Fisher suggested that this might be the case .

There are two simple scenarios for the velocity versus force if $`M`$ is relevant. One possibility is that the depinning transition is driven discontinuous immediately for any non-zero $`M`$. The average velocity would then have a discontinuity of

$$\overline{v}_{min}M^{\frac{\beta }{\varphi }}$$
(25)
We call this the first-order scenario.

If $`M`$ is relevant but the transition is still continuous, one would expect it to be in a different critical universality class. In this case, the scaling function $`B(y\mathrm{})y^\rho `$ so that $`\overline{v}(f)`$ remains continuous but with a new exponent

$$\beta =\beta _0\varphi \rho $$
(26)
with $`\beta _0`$ the quasistatic value. Asymptotically close to the critical force the new critical behavior would obtain but for $`f>>M^{\frac{1}{\varphi }}`$, the average velocity curve would crossover to the dissipative behavior. We refer to this as the new universality class scenario.

While more exotic scenarios may be possible, we will limit our consideration to the three scenarios enumerated above.

Note that we have explicitly not considered scenarios in which the velocity in the moving phase is hysteretic. Although we cannot rule this out entirely, the fact that the random environment through which the manifold moves acts, to some extent, like thermal noise, suggests that if there were more than one possible moving phase for the same applied force, there would be some stochastic process by which the system could jump from one to the other. The result of this would be that, as in equilibrium transitions, true “coexistence” would not be possible over a range of parameters. This argument does not apply to coexistence between static and moving phases, as the former are not subject to time dependent “noise”.

### B Numerical results in moving phase

Our primary numerical results in the moving phase were carried out on square two-dimensional samples of linear dimension $`L=128`$. The maximum pinning force $`f_p^{max}`$, was chosen to be either $`0.5`$ or $`1.0`$ so that the typical pinning force is comparable to the change in elastic forces caused one neighbor of a segment jumping.

We focus on the results for $`M`$ small enough that sublattice effects are not too important — for $`f_p^{max}=0.5`$, we study $`M0.8`$. Fig. 1(a) shows $`\overline{v}(F)`$ for $`0M0.6`$ in increments of 0.2. To generate these curves we start, for each $`M`$ at $`F>F_c^{}`$, with almost flat initial conditions — segment displacements random in $`(0,1)`$ — and then decrease the applied force very slowly until $`F_c^{}(M)`$ is reached.

For our finite systems with $`L=128`$, $`F_c^{}(M)`$ is defined as being the force below which the system halts after $`12,000`$ time steps, a value chosen so as to be long enough for transients to decay, but not so long that rare configurations of the randomness with anomalously strong pinning forces can dominate.

For $`M=0`$, there is a finite-size crossover regime in which the system may stop due to an anomalously strongly pinned region and the infinite system behavior will no longer be observed. This tends to occur when the average velocity while it is still moving is of order $`L^{\frac{\beta }{\nu }}`$. For $`L=128`$, this value is about $`0.013`$ and with this average velocity, the manifold will typically travel a distance several times the characteristic displacement $`L^\zeta `$ within $`12,000`$ time steps. As we shall see, our estimates are that at least one measure of the characteristic time scale, which grows as $`L^z`$, is only weakly dependent on $`M`$. Thus, we use this same criterion for non-zero $`M`$ while being aware that it may bias our scaling results slightly such that if time scales do change substantially with $`M`$ we may be observing more of either finite-size effects or non-equilibrium effects, the former if the time scale decreases with $`M`$ and the latter if it increases with $`M`$. For other system sizes, the equilibration time is decreased correspondingly, roughly with $`L^{z_0}`$, i.e., according to the dynamic scaling found in the dissipative limit.

We first study the dependence on $`M`$ of the critical force, $`F_c^{}`$, below which the steady state motion ceases, in particular to test whether $`F_c^{}(M)`$ is a singular or smooth function of $`M`$. See Fig. 2 The results, along with a quadratic fit,

$$F_c^{}(0)F_c^{}(M)=0.27(\pm 0.01)M+0.012(\pm 0.01)M^2$$
(27)
are shown in Fig. 2. \[Note that if we had included a constant in the fit, the constant would have vanished within one standard deviation, as it should, and the linear coefficient would have been only slightly modified to $`0.26\pm 0.01`$.\] A natural expectation — although overly naive — is a linear decrease of $`F_c^{}`$ by an amount $`M/Z`$. For small $`M`$, this appears to work rather well. The reason for this linear shift and the corrections to it will be discussed later.

The analytic fit should be compared to a fit — with the same number of parameters — to an arbitrary power-law $`M`$ dependence of the shift in $`F_c^{}`$ such as would obtain if $`M`$ were a relevant perturbation. The best fit to the $`F_c^{}`$ data yields an exponent $`\frac{1}{\varphi }=1.18\pm 0.02`$; note, however, that by eye the quadratic fit looks slightly better than the power-law fit. Although with a weakly relevant $`M`$ with a crossover exponent less than unity, as the power law fit suggests, one would presumably have a linear analytic term as well and thus the inferred $`\varphi 0.85`$ should not be taken too seriously, In any case, one must ask whether it is consistent with our other data. It does not appear to be: if we use this value of $`\varphi `$ to try and find a scaling function $`B(\frac{M}{f^\varphi })`$ for the velocity data of Fig. 1(a), the curves do not collapse. This suggests that either $`M`$ is irrelevant, or that it is sufficiently weakly relevant that the crossover exponent $`\varphi `$ is small enough that it would not dominate the shift in $`F_c^{}`$. Other data, as summarized below, suggests that, in fact, $`M`$ is irrelevant, at least for the steady state moving phase.

The mean velocity data can be used to obtain the critical exponent $`\beta `$ and see whether it depends on $`M`$. Figure 3 shows a log-log plot of the $`\overline{v}(F)`$ curves with the best $`F_c^{}(M)`$ value, that which makes the curve the most linear, determined by hand for each $`M`$. The error bars indicated are the rms variations in the average velocity over $`10`$ samples. The values of the applied forces used are separated by an interval of $`510^4`$, one-tenth of those used in Fig. 1(a).

The velocity critical exponent $`\beta `$ inferred from these data is, for $`M=0`$,

$$\beta _0=0.66\pm 0.03.$$
(28)
Surprisingly, this value appears to be consistent within one standard deviation with the data for all the $`M`$ values shown. Even for $`M=0.6`$, we find

$$\beta =0.67\pm 0.1,$$
(29)
although the straight line fit is only over one and a half decades in the reduced force, $`f`$, substantially less than the three decades of the fit for $`M=0`$. From Fig. 1(a), it is apparent that this reduced range of scaling is primarily due to a larger amplitude for the singular velocity for larger $`M`$. See Fig. 1(b).

The data suggest that the evidence is at least consistent with the dissipative critical behavior obtaining asymptotically for small $`M`$; i.e., with $`M`$ being an irrelevant perturbation.

We must be careful, however, especially as Figures 1(a) and (b), and 3 indicate that the minimum average velocity is increasing with increasing $`M`$. Might this observation suggest that there is a discontinuity opening up as a power of $`M`$, perhaps suggesting that the transition is driven discontinuous immediately? We can check the data against the finite-size crossover behavior expected in the $`M=0`$ quasistatic limit; this would yield

$$\overline{v}_{\mathrm{min}}L^{\frac{\beta }{\nu }}.$$
(30)
Figure 4 tests for this scaling and finds it to be consistent with the data: even for $`M=0.4`$, the slope of the log-log plot is $`0.91\pm 0.01`$, which agrees within one standard deviation with the $`M=0.0`$ slope. While we do observe an increase in the minimum average velocity with increasing $`M`$ at fixed size, it is not the size-independent power-law increase with $`M`$ that would have been expected if the transition became discontinuous. Instead, there is an $`M`$-dependent coefficient

associated with the average velocity, just as occurs for small $`M`$ in mean field theory. Similar minimum velocity data was also obtained for stronger randomness, with $`f_p^{max}=1`$. We thus see that the finite-size data are consistent with the dissipative scenario, with small stress overshoots being irrelevant for the velocity versus force curves.

But it is still possible that a new critical universality class is emerging for small $`M`$ if $`\rho <<1`$, so that the emergence of the new universality class would be difficult to detect by simply measuring the velocity exponent as this would be little changed from its dissipative-limit value. We therefore look more closely at the finite-size crossover regime to investigate whether other aspects of the behavior really look similar to the $`M=0`$ quasistatic depinning.

Anticipating that it might be the stopping behavior that would distinguish between quasistatic and overshoot dynamics, we have explored some of the dynamics of the stopping process. Given that there is a distribution of $`F_c^{}`$’s, we let the manifold equilibrate at an $`F`$ three standard deviations above the average of $`F_c^{}(M)`$. We then lowered $`F`$ to the average of $`F_c^{}(M)`$ and waited for the manifold to come to a stop. Figure 5 plots the instantaneous spatially averaged velocity $`v(t)`$ averaged over many samples as they come to a stop at time $`t_{stop}`$. Here we explicitly see that the various $`M`$ samples come to a stop in a similar gradual manner; only the amplitudes vary. The differences between the curves are presumably due to the $`M`$ dependence of the amplitude of the steady state $`\overline{v}(F)`$ curves. The final stages of the stopping process can be analyzed — at least for the dissipative ($`M=0`$) case — by a simple scaling argument. Once only a small fraction of the system is still moving, the average velocity will be inversely proportional to the area $`L^d`$. Since the velocity scales as $`\xi ^{\zeta z}`$ and lengths as $`t^{\frac{1}{z}}`$, we expect

$$v(t)\frac{(t_{stop}t)^{(d+\zeta z)/z}}{L^d};$$
(31)
the data on a log-log plot, Fig. 5, are reasonably consistent with this for all values of $`M`$ tested.

We can also probe the finite-size crossover regime in terms of the roughness. We define the maximum width of the manifold $`W_{\mathrm{max}}`$ as the absolute value of the maximum deviation of the displacement from its spatially averaged value, $`\overline{h}`$: $`W_{\mathrm{max}}=<\mathrm{max}|h(𝐱)\overline{h}|>`$. This should scale as $`L^\zeta `$ with $`\zeta 2/3`$ for $`M=0`$. If $`M`$ is irrelevant, then the same should hold true for all small $`M`$. Figure 6 demonstrates that the maximum width does indeed obey this scaling with system size, but with apparent values of $`\zeta `$ that are somewhat larger than $`\frac{2}{3}`$ even for $`M=0.0`$ for which $`\zeta =0.75\pm 0.05`$ is inferred from the data with $`f_p^{max}=0.5`$. Although the exponent appears to be $`M`$-independent, there is an $`M`$-dependent coefficient with the overall width increasing — albeit only slightly — with increasing $`M`$. For the stronger randomness data, $`f_p^{max}=1`$, this tendency is less strong but that is because the value of $`M`$ at which the behavior changes character is larger for stronger pinning. Overall it appears that, at least up until $`M=0.6`$, the finite-size crossover regime looks similar to the dissipative $`M=0`$ case.

We can also determine $`\zeta `$ from the power spectrum of the displacement just after the motion has stopped from the moving phase; we use the same equilibration time and applied force increments as in Figs. 7(a) and 7 (b). In two dimensions,

$$<|\widehat{h}(k)|^2>=d^2xe^{i𝐤𝐱}<h(x)h(0)>k^{2\zeta 2},$$
(32)
where the brackets denote averaging over samples. From Fig. 8(a) we see that fitting to this form over a range of one and a half decades for samples with $`L=256`$ yields $`\zeta =0.72\pm 0.02`$ and $`0.74\pm 0.01`$ for $`M=0.0`$ and $`M=0.5`$ respectively, again in mutual agreement and similar to the values from the maximum width discussed above. But again, apparently slightly larger than the $`\zeta =\frac{2}{3}`$ value from the first-order epsilon expansion. Again, the cumulative evidence suggests that at least up until $`M>0.7`$, it is not likely that a new universality class is emerging in the “equilibrium” moving or stopped phase as $`\overline{v}0`$.

It is useful to compare our values of the exponent $`\zeta `$ with those previously obtained. Leschhorn et. al. found, in numerical simulations, $`\zeta =0.75\pm 0.02`$ for $`M=0.0`$.  The most solid theoretical result is that $`\frac{ϵ}{3}`$ is a lower bound for $`\zeta `$. This comes from application of finite size scaling to the connections between the variations of the critical force and the correlation length exponent — see discussion below. While Narayan and Fisher  had argued that the value of $`\frac{ϵ}{3}`$ is exact to all orders in $`ϵ`$, Chauve et. al.  have computed the exponents to second-order in $`ϵ`$ and found that the roughness exponent is increased to $`\zeta =\frac{ϵ}{3}(1+0.14331ϵ)+𝒪(ϵ^3)`$, which naively extrapolated to two dimensions yields $`\zeta 0.86`$ substantially higher than the value inferred numerically. Whether or not this discrepancy is due to the neglect of terms higher order in $`ϵ`$, to some problem with the expansion, or to corrections to scaling remains in doubt. It is worth noting in this context that for one-dimensional systems with the long range interactions appropriate for crack fronts, Ramanathan and Fisher  found that to obtain reliable and universal values of $`\zeta `$ analysis of corrections to scaling were needed. With these included, they found a value of $`\zeta `$ very close to that predicted from the first-order $`ϵ`$-expansion without any higher order corrections.

For completeness, we have also studied the spatial power spectrum of the displacements just above the depinning transition. Figure 7(a) and (b) shows the square root of the averaged power spectrum for $`\overline{v}<<1`$. The roughness exponents obtained for $`M=0.0,0.2,0.5`$ are again roughly in agreement with each other — $`\zeta =0.65\pm 0.01`$, $`0.65\pm 0.02`$, and $`0.69\pm 0.01`$ respectively — but slightly smaller than those found for the stopped manifold. In the moving phase, however, there should be a crossover at long wavelengths, $`k\xi _v1`$, to the Edwards-Wilkinson universality class with only logarithmic roughness . This arises because on scales larger than $`\xi _v`$, the motion makes the randomness appear like white noise in both space and time and the displacement correlation function becomes

$$<\widehat{h}(k,\omega )\widehat{h}(k,\omega )>\frac{1}{|i\omega +Dk^2|^2},$$
(33)
with $`D(F)`$ an effective diffusion constant. This crossover is observed for the smallest $`k`$ in Fig. 7 (a) where the slope of the power spectrum decreases in contrast to the data taken after the motion has stopped. We note that this difference in slope between the moving and stopped spectra is not as prominent when the randomness is stronger as shown in Fig. 7 (b).

For the moving configurations we observe a peak in the power spectrum at $`𝐤=(\pi ,\pi )`$. This peak is caused by the tendency of one segment’s motion to trigger jumps of its neighbors at the next time step. The structure and amplitude of the peak looks similar for all of the $`M`$’s shown, although the wavevector dependence indicates that there are somewhat more segments participating in the sublattice behavior at $`M=0.5`$ than at $`M=0`$; but only about a third more.

The dynamic exponent in the moving phase, $`z`$, can be determined in various ways from data near to the critical force. We first study the non-equilibrium roughening of the manifold starting with almost flat initial conditions. We define

$$w^2(t)=\overline{(h(𝐱,t)\overline{h(t)})^2}t^{\frac{2\zeta }{z}},$$
(34)
with the overbar denoting spatial averaging over the sample. The scaling behavior is expected for times short compared to the critical correlation time $`\tau `$ which diverges at $`F_c^{}`$. Once the roughness exponent has been calculated independently from the same set of configurations in Fig. 7(a), $`z`$ can be extracted from the log-log plot of $`w^2(t)`$ vs. $`t`$ as is done in Fig. 8(a)(b) For $`M=0.0`$ and $`M=0.2`$ we see that the dynamic exponents $`z`$ are very similar: $`z=1.51\pm 0.03`$ for both. However, for $`M=0.5`$, the data are somewhat further above the transition as they correspond to an equilibrium average velocity of $`\overline{v}=0.08`$. While these data could be fit with the same $`z`$ over a limited range of times, there is clearly some new physics emerging: an upward curvature on the log-log plot and a substantial — more than a factor of two — overshoot in the velocity before it settles down to its steady state value. This effect is related to the change in the dynamical onset of the motion to which we will turn in Section VII.

Aside from the transient effects associated with approach to steady state, all our measurements in the moving phase suggest that the critical behavior is most consistent with the dissipative universality class obtaining for all sufficiently small $`M`$ as the transition is approached from above; this despite the $`M`$-dependent shift in $`F_c^{}`$ and the concomitant hysteresis that is possible because of the existence of linearly stable static configurations up to $`F_{c0}`$ which is larger than $`F_c^{}`$. We will analyze this paradox later. For now, it appears that the mean-field-like scenario has won out over the first-order scenario and the new universality class scenario.

### C Amplitude ratios

If we accept that $`M`$ is irrelevant for the steady state moving critical behavior, the observed changes as $`M`$ is increased in the moving phase for $`F`$ just above $`F_c^{}`$ appear to be primarily attributable to the increased amplitude $`A_v(M)`$ of the velocity near $`F_c^{}`$, $`\overline{v}A_v(M)f^\beta `$. As seen in Figs. 1(b) and 3: For larger $`M`$, the velocity rises more rapidly with increasing $`F`$ until it becomes close to $`\frac{1}{2}`$ at which point sublattice effects set in. A similar increase in the amplitude of the velocity is observed in mean field theory, with the exponent $`\beta `$ unchanged for small $`M`$ but the amplitude of the velocity growing with $`M`$. In both cases, this implies that the width in $`F`$ of the depinning transition narrows with increasing $`M`$ — naively just by the narrowing of the range of $`FF_c^{}`$ over which $`\overline{v}`$ is small.

If the universality class of the critical behavior in the moving phase is independent of $`M`$ over a range of $`M`$, then there should be various universal relations between non-universal coefficients: universal amplitude ratios. A priori, we would expect three non-universal scale factors associated with the scaling relationships between length and, respectively, deviation from criticality, $`FF_c^{}`$; displacement, and time. We can define these, for example, by the scaling of the correlation length,

$$\xi _v\left(\frac{A_F}{FF_c^{}}\right)^\nu ;$$
(35)
mean square displacements at separations smaller than $`\xi `$,

$$[h(𝐱,t)h(0,t)]^2A_h^2|𝐱|^{2\zeta };$$
(36)
and velocity

$`\overline{v}`$ $``$ $`{\displaystyle \frac{A_h}{A_t}}\xi _v^{\beta /\nu }`$ (37)
$``$ $`{\displaystyle \frac{A_h}{A_t}}\left({\displaystyle \frac{FF_c^{}}{A_F}}\right)^\beta `$ (38)
$`=`$ $`A_v(FF_c^{})^\beta `$ (39)
where we measure all lengths in units of the lattice constant. But in the absence of stress overshoots the statistical “tilt” symmetry of the system that relates the exponents $`\zeta `$ and $`\nu `$ via the triviality of the averaged response to a static spatially varying additional applied force, also relates the associated coefficients:

$$\frac{A_F}{A_h}=C_KK$$
(40)
where $`C_K`$ is a universal dimensionless coefficient and $`K`$ is the long-wavelength elastic constant, which, in our model, is simply the inverse of the coordination number $`Z`$.

For $`M=0`$, we thus expect that amplitudes of scaling laws that only involve the three exponents $`z,\nu `$ and $`\zeta `$ should be expressible in terms of the two amplitudes $`A_h`$ and $`A_t`$ only. For example, the rms variations in the critical force $`F_c^{}`$, in finite size systems should be expressible in terms of $`A_F`$ or, via Eq. (36), $`A_h`$:

$$\sqrt{\mathrm{var}(F_c^{}(L))}C_\mathrm{\Delta }C_FL^{\frac{1}{\nu }}$$
(41)
with $`C_\mathrm{\Delta }`$ universal.

It is not clear, a priori, whether the statistical tilt symmetry argument can be applied in the presence of stress overshoots and the concomitant local hysteresis. This is because, in essence, it relies on the history independence of linear response of at least some quantities. In the moving phase, on which we are currently focusing, it seems reasonable that the argument should apply and the amplitude ratios hence be related as in the dissipative case. But we should remain alert to the possibility that apparent failure of expected scaling laws may be due to this assumption.

Qualitative examination of the numerical data for the roughness (e.g. Fig. 7(a) and the variations in the critical force $`F_c^{}`$, (see Fig. 9)) suggest that neither of the amplitudes $`C_h`$ and $`C_F`$ are strongly dependent on $`M`$ in the range studied. One would then guess that the dependence of the amplitude of the velocity, $`A_v(M)`$, on $`M`$ is primarily caused by a decrease in the characteristic time scale as $`M`$ increases. This will have consequences for the behavior of other quantities as will now be discussed.

A useful quantity to study is the temporal fluctuations in the instantaneous spatially averaged velocity $`v(t)`$ in finite-size samples: this has information about both length and time scales. For a fixed average velocity $`\overline{v}`$, the magnitude of the fluctuations should increase with increasing $`M`$ because the system is effectively closer to the transition. Since regions of size of order the velocity correlation length will fluctuate roughly independently, the variance of the instantaneous velocity should, in a system much larger than $`\xi `$, be

$$\mathrm{var}(v)(v(t)\overline{v})^2)\overline{v}^2\frac{\xi _v^d}{L^d}.$$
(42)
Since we have not yet defined the correlation length precisely, this could well serve as its definition, thereby fixing the definition of the amplitude $`A_F`$. The universality of the amplitude ratios can then be checked by comparing the $`M`$ dependence of the amplitude of $`\mathrm{var}[v(t)]`$ with those of the variations in $`F_c^{}`$ and the roughness. These results are presented in Table I.

The correlation time $`\tau _v`$ can also be obtained from the truncated velocity fluctuations,

$$\mathrm{\Gamma }(t)v(t)v(0)\overline{v}^2.$$
(43)
Integrating over time yields,

$$\mathrm{\Gamma }(t)𝑑tC_{vv}\tau _v\overline{v}^2\frac{\xi ^d}{L^d}$$
(44)
with $`C_{vv}`$ a universal coefficient. The resulting $`\tau _v`$ can, together with $`\xi _v`$ be used to check other scaling relations. See Table I. Note however that there are difficulties associated with the subtraction needed to obtain the truncated correlations in the most interesting regime in which the fluctuations are large and an accurate extraction of $`\overline{v}`$ problematic.

The velocity-velocity correlations can also be used to probe the nature of the dynamics; in particular by studying the power spectrum which is the Fourier transform, $`\widehat{\mathrm{\Gamma }}(\omega )`$ of $`\mathrm{\Gamma }(t)`$. To understand how this is expected to behave in the dissipative limit, it is useful to consider the local velocity-velocity correlation function, $`v(𝐱,t)v(0,0)`$. At long distances, $`|𝐱|\xi _v`$ or time separations, $`t\tau _v`$, this will approach $`\overline{v}^2`$. But within a correlation space-time volume, the local velocities will be characteristic of avalanche events and hence be fractal. The correlations will be proportional to $`\overline{v}`$ times a conditional expectation of $`v(𝐱,t)`$ given that there is motion — i.e. a jump — at $`(0,0)`$. These conditional correlations within the space-time correlation volume will reflect the fractal structure, being of order $`1/|𝐱|^{z\zeta }`$ or $`1/t^{1\frac{\zeta }{z}}`$ whichever is smaller. Integrating the associated scaling forms over $`𝐱`$ and Fourier transforming in time, one finds that for $`\omega \frac{1}{\tau }`$,

$$\widehat{\mathrm{\Gamma }}(\omega )C_{PS}\frac{\overline{v}^2}{\omega ^{\frac{d+\zeta }{z}}}\frac{\xi ^d}{L^d}\frac{1}{\xi ^{dz+\zeta }A_t^{\frac{dz+\zeta }{z}}},$$
(45)
with $`C_{PS}`$ a universal coefficient.

The log-log plot of the square root of the velocity power spectrum (Fig. 10(b)) appears to exhibit power-law behavior for large $`\omega `$ for both the values of $`M`$ shown. For $`M=0`$, the observed exponent is close to the expected value of $`\frac{d+\zeta }{2z}0.86`$ over two and a half decades in frequency. For $`M=0.6`$, the best fit slope is somewhat larger, but some curvature is evident and consistency with the dissipative result is not ruled out. In both of these sets of data, there is a crossover at low frequencies to a flat spectrum. This is more pronounced in the $`M=0`$ data which are effectively further from the critical regime — $`\xi /L`$ smaller — than the $`M=0.6`$ data because the two sets of data were taken at the approximately the same $`\overline{v}`$ which is closer to the corresponding $`\overline{v}_{min}(L,M)`$ for $`M=0.6`$.

Not only does the value of $`F_c^{}(M)`$ provide us useful information, but its variations do as well. In the finite-size limited scaling regime in which $`L\xi `$, $`\xi `$ is replaced with $`L`$ in scaling laws and we expect

$$[\mathrm{var}(F_c^{}(M))]^{1/2}L^{1/\nu }.$$
(46)
Figure 9 demonstrates this scaling with the system length, the obtained correlation length exponents being $`0.78\pm 0.02`$, for $`M=0`$ and for $`M=0.6`$, $`\nu =0.81\pm 0.02`$. Calculations of Narayan  yield the leading irrelevant eigenvalue at the quasistatic fixed point as approximately $`\frac{ϵ}{3}`$. This suggests a fit of the data of Fig. 9 with the form $`\mathrm{var}(F_c^{}(M))=CL^{\frac{1}{\nu }}/(1+B_\nu L^{\frac{2}{3}})`$ yielding $`\nu =0.77\pm 0.03`$ with $`B_\nu =0.78\pm 0.01.`$

We have found that a variety of properties of the steady state moving phase as well as the variations of the critical force at which the system stops on decreasing the drive are all consistent with critical behavior that is independent of the magnitude of the stress overshoots. Given the local hysteresis that is intrinsic with stress overshoots, this universality is more than a little surprising. In the next section, we consider to what extent this applies more generally, in particular as far as macroscopic hysteresis.

## V Hysteresis

We now turn to an analysis of the hysteretic phenomena that are implied by the coexistence of moving and stationary solutions at the same force in the presence of stress overshoots. A crucial question which we must address is whether hysteresis persists in macroscopic systems that are not prepared in special ways. In particular, are there hysteresis loops with a width that is non-zero in the limit of large systems? If not, as we shall see is the case, how does the hysteresis depend on system size? Can one understand this in terms of the purely dissipative dynamics that appear to control the properties of the steady state moving phase? Or is new physics needed?

In Fig. 1(a), hysteresis loops are shown for typical samples of size $`128^2`$ with $`f_p^{max}=0.5`$ and $`M`$ from $`0`$ to $`0.6`$. An upwards arrow indicates the force $`F_{c_{\mathrm{hyst}}}^{}`$, at which the system starts moving again after it has been stopped at $`F_c^{}`$ by a gradual decrease in the force. For $`M=0`$ no hysteresis is apparent while for positive $`M`$, the difference $`F_{c_{\mathrm{hyst}}}^{}F_c^{}`$ appears to be close to $`M/4=M/Z`$, the magnitude of the stress overshoot. On the basis of these data, it would appear that the situation is rather simple: on decreasing the force the steady state moving phase and the stopping process are not qualitatively dependent on $`M`$, but once the system has stopped, an increase of the force by $`M/Z`$ is required to start it up again. Once restarted, the velocity rapidly increases to that of the apparently unique moving “state”. The reason for this macroscopic hysteresis would appear to be simple: If the force on each of the segments caused by the last motion of its neighbors before the system stopped was not enough to cause it to move, then at later times the force will be less than that needed to make a segment move by at least $`M/Z`$ as the stress overshoots from its neighbors jumping will no longer be in effect. If this applies to all of the segments, it should be necessary to increase the force back up again by at least $`M/Z`$ before anything can start to move. This would imply truly macroscopic hysteresis that is independent of size for large systems.

A more careful examination of both the data and the argument above shows that it is fallacious: even with $`f_p^{max}=0.5`$, a small fraction of samples have substantially narrower hysteresis loops. There must thus be some segments that can be restarted by an increase in the force from the stopped state by less than $`M/Z`$. In the next subsection, we discuss the origin of this effect, but first, we present and analyze the numerical data.

### A Distributions of $`F_{c_{\mathrm{hyst}}}^{}`$

For reasons that will become clear later, we can obtain more useful data on the hysteresis by increasing the strength of the randomness. Most of our detailed hysteresis data is for $`f_p^{max}=1.0`$ and $`M=0.8`$, the latter being sufficiently large that the effects of overshoots are strong, but not so large (for this larger value of $`f_p^{max}`$) that sublattice effects start to play a role. For these parameter values, the mean force at which the system stops on decreasing $`F`$ under the procedure discussed in Section IV is

$$F_c^{}0.0665$$
(47)
with the rms variations about this of

$$\sqrt{\mathrm{var}(F_c^{}(f_p^{max}=1.0,M=0.8))}0.0006.$$
(48)
In Fig. 11(a), the distributions of $`F_{c_{\mathrm{hyst}}}^{}`$ are shown for various system sizes. It can be seen that they are much broader — by almost two orders of magnitude for $`L=128`$ — than the distributions of $`F_c^{}`$. The shapes of the distributions at first appears rather strange: For the smaller system sizes, a substantial fraction of the weight is in a narrow peak that has similar width to that of the distribution of $`F_c^{}`$, but is shifted up from this by an amount $`F_{c_{\mathrm{hyst}}}^{}F_c^{}M/Z`$. In the largest samples, this peak has completely disappeared and we see that the width of the hysteresis loop has narrowed considerably. The narrowing with size of the median width of the hysteresis loops is shown in a log-log plot in Figure 12 . As was evident in the shape of the distributions, a crossover length of around $`L=20`$ is seen in these data. For small sizes the hysteresis loops have width that is typically close to $`M/Z`$. But for the large sizes, the typical width appears to decrease as a power of $`L`$:

$$F_{c_{\mathrm{hyst}}}^{}F_c^{}\frac{1}{L^{\frac{1}{\mu }}}$$
(49)
with

$$\mu 1.15\pm 0.03$$
(50)
for $`M=0.8`$ and

$$\mu 1.35\pm 0.02.$$
(51)
for $`M=1.0`$. Note that these exponents are obtained from a narrow range of length scales and crossover behavior is likely to be playing a role; we will return to the issue of crossovers later.

Nevertheless, in spite of uncertainties in the asymptotic size dependence of the hysteresis loops, the overall trend is clear: in the limit of large systems, the width of the hysteresis loops vanish! This in spite of the fact that there are many linearly stable static configurations that coexist with the moving state up to the force $`F_{c0}`$ which is substantially greater than $`F_{c_{\mathrm{hyst}}}^{}`$.

Before discussing this result, it is instructive to consider what happens in the dissipative limit, $`M=0`$. Although it might appear that there would be no hysteresis in this case, this is not strictly correct for finite size systems that have been stopped from a moving state. At the force at which the system stops, it gets stuck in a somewhat anomalously strong pinning region (how anomalous depends on the rate of decrease of the force). Getting it unstuck from such a region, in the sense that all parts of the system move for at least some distance, can require an increase in the force that is comparable to the width of the distribution of $`F_c^{}`$. Thus we expect “hysteresis loops” in the dissipative limit to have a width of order $`L^{1/\nu }`$.

If the systems with stress overshoots behaved like the dissipative case in all universal aspects, one would expect the asymptotic large system-size dependence of the width of the hysteresis loops to have the same exponent as the dissipative case, i.e. that $`\mu =\nu `$. Up to questions about estimation of uncertainties in the presence of complicated crossovers, it appears that this is not the case: we seem to find

$$\mu >\nu $$
(52)
corresponding to system size dependence of $`F_{c_{\mathrm{hyst}}}^{}`$ being slower than that of $`F_c^{}`$. If this were indeed the case, we would expect that it would most likely hold asymptotically for any $`M>0`$. Unfortunately, the range of data is not so large as to conclusively rule out equality rather than inequality, although if $`\mu =\nu `$ one would probably need either a large dimensionless amplitude ratio between the coefficients of the size dependence of the two critical forces, or strongly nonmonotonic behavior; we will later explore such scenarios. But for now we focus on the scaling behavior that seems to be emerging for system sizes larger than of order 20 or so for ($`f_p^{max}=1.0,M=0.8`$).

It appears that the larger system sizes do exhibit scaling behavior of the distributions, a more stringent test than exponents. Indeed, the size dependence of the distributions provides a useful way to understand the causes of the size dependence of the hysteresis.

Let us assume that restarting on increasing the force after stopping occurs via some kind of nucleation process whose occurrence is dominated by scales that are much smaller than the system size. Then we expect a density of nucleation segments — or “seeds” — with a distribution of values of the local critical forces $`F_s`$ needed to restart. If this distribution extends down to $`F_c^{}`$, larger systems are more likely than smaller ones to have a seed with a small $`F_s`$. As the lowest $`F_s`$ in a given stopped configuration will be the one that determines $`F_{c_{\mathrm{hyst}}}^{}`$, this will yield a distribution of critical forces that becomes squeezed down to $`F_c^{}`$ as $`L\mathrm{}`$. A simple check on the assumption of locality of the seeds is accomplished by estimating the distribution of $`F_{c_{\mathrm{hyst}}}^{}`$ for a sample of size $`L`$, by considering it as being made of $`b^d`$ independent samples of size $`L/b`$ whose $`F_c^{}`$s are drawn independently from the observed distribution for these smaller size samples. In Figs. 13(a) and 13(b), this is carried out for $`L=128`$ and $`b=2,4`$, and $`8`$. As can be seen, the distributions obtained agree very well with those measured for the $`L=128`$ samples directly. This agreement is particularly striking given that the data for $`F_{c_{\mathrm{hyst}}}^{}`$ for $`L=16`$ in Figs. 11(a) and 11(b) have a very different form than those for the large samples. (This difference in form for the distributions is the source of the crossover observed in the size dependence of the width of the hysteresis loop.)

The agreement of the actual distribution of $`F_{c_{\mathrm{hyst}}}^{}`$ for $`L=128`$ with that obtained from the distribution with $`L=16`$ suggests that the nucleation process in a sample of size 128 is typically dominated by regions whose diameter is less than 16. Once such a small nucleation region gets going, it will typically expand to make the whole system restart, independent of the existence or lack thereof of seeds in other regions, or of other stochastic properties of the rest of the system. Before analyzing the consequences of this, we must caution that a crucial question is whether such a relation between the distributions of $`F_{c_{\mathrm{hyst}}}^{}`$ for systems of size $`L`$ and $`L/b`$ holds, in the limit of large $`L`$, for only a limited range of $`b`$, for any $`bL`$, or for $`b`$ up to some (subdominant) power of $`L`$.

For now, we will extract the shape of the distribution from the observation that there is at least a substantial range of $`b`$ over which the relationship between the distributions of $`F_{c_{\mathrm{hyst}}}^{}`$ for size $`L`$ and size $`L/b`$ does hold.

The basic picture of the restarting being controlled by the least pinned of many independent seeds enables one to relate the size dependence of the median $`F_{c_{\mathrm{hyst}}}^{}`$ to the form of the distribution. Picking the minimum of the $`F_{c_{\mathrm{hyst}}}^{}`$’s from $`S=b^d`$ subsystems each of which has a distribution of $`F_{c_{\mathrm{hyst}}}^{}`$ that vanishes as $`(F_{c_{\mathrm{hyst}}}^{}F_c^{})^{d\mu 1}`$ for small $`F_{c_{\mathrm{hyst}}}^{}F_c^{}`$, yields a power law decrease with $`L`$ of the width of the distribution, and, indeed, the actual form of the distribution: We expect a Weibull distribution with one non-universal scale parameter, $`f_h`$, .

$`\mathrm{Prob}[dF_{c_{\mathrm{hyst}}}^{}]dF_{c_{\mathrm{hyst}}}^{}{\displaystyle \frac{d\mu }{f_h}}\left({\displaystyle \frac{F_{c_{\mathrm{hyst}}}^{}F_c^{}}{f_h}}\right)^{d\mu 1}`$ (53)
$`\mathrm{exp}\left(\left({\displaystyle \frac{F_{c_{\mathrm{hyst}}}^{}F_c^{}}{f_h}}\right)^{d\mu }\right).`$ (54)
As can be seen in Fig. 14, this yields a rather good fit to the data for $`L=128`$ with the value of $`\mu `$ extracted from the size dependence of the median. If, instead, we do a best fit to the shape of the distribution for the largest size, we find, $`\mu 1.08`$ for $`M=0.8`$ and for $`M=1.0,f_p^{max}=1.0`$, $`\mu 1.25`$. Note that these values are slightly smaller than those obtained from the size dependence and thus somewhat closer to $`\nu `$; this may well be a sign of slow crossover to asymptotic behavior that is like the dissipative limit.

### B Origins of seeds for restarting

We next develop an understanding of the origins of the unusual hysteretic behavior that was found in the numerical studies: in particular, the origins of the seeds for nucleation of motion when the force is increased back up again after the motion has stopped. To do this, we need to understand how the manifold stops moving as $`F`$ is decreased to $`F_c^{}(M)`$ as it is this that sets up the configurations in which the seeds exist. We first analyze the basic role of the stress overshoots in the steady state moving phase.

A very crude approximation to the effects of the stress overshoots is to ignore their local and transient natures. We thus consider an artificial model in which if any segment has moved on the previous time step, the force on all the segments is increased by $`M/Z`$ above what it would be with purely dissipative dynamics. As long as something is always moving, this is identical to merely increasing the applied force by $`M/Z`$. But once the system has stopped — because of a decrease in the applied force or because of running into a strongly pinned region — no segment can move again until the force is increased by $`M/Z`$; this is because any segment that could move with less of an increase in $`F`$, should, a fortiori have moved already because of the stress overshoot that was present before the motion stopped. In terms of the distribution of the total force, $`\{f(𝐱)\}`$, on a segment, which must be positive for it to jump, the stopped configuration will have a depletion layer: no segments with $`f(𝐱)`$ in the interval $`(M/Z,0)`$ (and of course none with $`f(𝐱)>0`$). The behavior of this crude model is thus very simple: the $`v(F)`$ curve is shifted down (in $`F`$) by $`M/Z`$, the steadily moving states are identical to those at $`M=0`$ with $`FF+M/Z`$. When the force is increased after stopping, no motion will occur until the depletion layer disappears, therefore $`F_{c_{\mathrm{hyst}}}^{}=F_c^{}+M/Z`$ for each sample.

Before we return to the model or primary interest, it is worth noting that a model which is much less pathological than the crude model discussed above nevertheless has much of the same behavior. This non-additive stress overshoot model has nearest neighbor stress overshoots that last for one time step and a “self overshoot” (analogous to inertia) that likewise lasts for one time step. But the stress overshoots are non-additive so that any site has a stress overshoot that is either zero or $`\frac{M}{Z}`$. For example, if two nearest neighbors of a segment jump, that segment feels a stress overshoot of $`\frac{M}{Z}`$, in contrast to the $`\frac{2M}{Z}`$ overshoot it would feel in our primary model. The partial equivalence between this non-additive stress overshoot model and the crude model can be simply understood: If the total force on a segment does not change from the previous time step it cannot jump at the next time step. Therefore changes in the total force on a segment are what determines whether a segment jumps or not at the next time step. At fixed $`F`$, changes in the total force on a segment arise from nearest neighbors jumping at the previous time step and from the segment itself jumping at the previous time step. In the non-additive model, each of these will involve an extra $`\frac{M}{Z}`$ stress. Thus as long as motion has existed somewhere in the systems for more than one time step, given a configuration of the crude model and which segments have jumped on the previous time step, there is an exactly equivalent configuration of the non-additive model which will have the same dynamics at all future times as long as $`F`$ remains fixed. But the dynamics are not fully equivalent: when the first segment moves in the crude model, it can trigger others far away; this cannot happen in the non-additive model, Nevertheless, the steady state velocity as a function of force will be the same in these two models, with the critical force in the infinite system limit shifted down by exactly $`\frac{M}{Z}`$ from the dissipative case. The hysteresis loops will also be similar, but not identical: in both cases there will be a depletion layer of width $`\frac{M}{Z}`$ after the system has stopped and the force will have to be increased by this for motion to start again, But in finite size samples the behavior will be slightly different as it is much more likely in the non-additive model that motion could start in one region but die out: the actual critical force for restarting would then be slightly higher. The dynamics of the transient motion on restarting would also differ due to the locality of the non-additive model.

The decrease of $`F_c^{}`$ with $`M`$ in the crude and non-additive models is the underlying cause for the linear decrease of $`F_c^{}(M)`$ in the primary model: As long as segments only move in response to their neighbors moving, what is crucial in determining $`F_c^{}`$ is how the system gets through potential sticking points. Some of these are likely to involve only one neighbor of a segment moving at the previous time step; if they do, then the critical force at which they can proceed will, in the absence of other changes of the dynamics due to $`M`$, be just $`M/Z`$ lower than it would be with $`M=0`$. In the limit of small $`M`$, we expect that the sequence of jumps at $`F`$ will be very close to that at $`F+M/Z`$ in the absence of overshoots.

The non-linear part of the dependence of $`F_c^{}`$ on $`M`$ is of a different origin. As $`M`$ grows, the equivalence between the sequences of jumps at different values of $`M`$ no longer obtains because of, for example, the effects of two neighbors jumping at the same time which increases the stress on a segment by $`2M/Z`$. This will tend to make stopping less likely as a region can be restarted by motion in other regions that is caused by such multiple-neighbor jumps. As one would thus expect, $`F_c^{}`$ decreases faster than linearly as $`M`$ is increased.

The focus on particular sites and whether they can be retriggered by a given increase in $`F`$ is also useful for understanding the hysteresis in the model of primary interest. A crucial question about the local dynamics is: How close can a segment be to moving without one of its neighbors having moved on the previous time step? We must consider the most recent time in the past, say time one, at which a neighbor of the segment $`𝐱`$ of interest moved. For simplicity, let us assume that none of the neighbors moved at time zero. At time zero, the total force on $`𝐱`$ is then

$$f(𝐱,t=0)=Ff_p(0)+\sigma (0)<0$$
(55)
with $`\sigma (0)=\sigma (𝐱,0)`$ and $`f_p(0)=f_p(𝐱(0),h(𝐱(0)))`$ the initial elastic force and initial pinning strength, respectively, at $`𝐱`$. If $`n`$ out of the $`Z`$ neighbors jump at time one,

$$f(𝐱,1)=Ff_p(0)+\sigma (0)+n(1+M)/Z.$$
(56)
If this is negative, then $`𝐱`$ will not move and, at later times, the force $`f(𝐱,t>1)<nM/Z`$ as the stress overshoot will no longer apply; thus segment $`𝐱`$ will not be in the depletion layer. If, however, $`𝐱`$ does jump at time two in response to its neighbors jumping, i.e., if $`f(𝐱,1)>0`$, then the total force on it at later times will be

$$f(𝐱,2)=Ff_p(1)+\sigma (0)1+n/Z$$
(57)
with a new random pinning force $`f_p(1)`$. As long as $`f(𝐱,2)<0`$, then this segment will not move further unless one of its neighbors does. Thus this represents a possible local configuration when the system has just stopped.

The condition $`f(𝐱,0)<0`$ implies, from Eq. (55) and Eq. (56) that $`f(𝐱,1)<n(1+M)/Z`$. The maximum of $`f(𝐱,2)`$ is then obtained when $`f_p(1)`$ is minimal (i.e. zero) and $`f_p(0)`$ maximal (i.e. $`f_p^{max}`$); this yields $`f(𝐱,2)<f_p^{max}1+n/Z`$. If the force is now increased back up by $`\mathrm{\Delta }F`$ so that $`f(𝐱,2)+\mathrm{\Delta }F=0`$, the segment $`𝐱`$ will jump and could trigger restarting of the overall motion. We must thus ask how close to zero $`f(𝐱,2)`$ can be.

If the last motion in any region were always via single neighbors triggering each other, then there would be a depletion layer on stopping and macroscopic hysteresis for $`f_p^{max}<1\frac{1}{Z}=\frac{3}{4}`$ in our case with coordination number $`Z=4`$. \[This would obtain if the stress overshoots were not additive but instead such that any number of neighbors jumping yielded the same value of the stress overshoot as if only one had.\] But in reality, it is possible that any number of neighboring segment on a segment could jump at one time and then not again; some results on the simultaneous hopping of a number of neighbors are presented in Table II. As such sets of simultaneous jumps can occur for any $`nZ`$ — including for $`n=Z`$ — even as the system is stopping, there will be no depletion layer even for arbitrarily weak randomness and no hysteresis in the infinite system limit. Nevertheless, for weak pinning, the depletion layer will only be filled by simultaneous jumps of multiple neighbors followed by a jump of the central segment that does not trigger further jumps of any of its neighbors, Although we expect that a local condition such as this will always occur for some finite fraction of the segments, it appears likely that the potential seeds for restarting with small increases in $`F`$ will be very rare for weak pinning; this, indeed, turns out to be the case. It is the proximate cause of the long crossover lengths apparent in the size-dependent distributions of $`F_{c_{\mathrm{hyst}}}^{}`$.

The above analysis gives a qualitative explanation for the lack of macroscopic hysteresis. But to explain the observed dependence of the widths of the hysteresis loops on system size, we must understand the density of states of segments with small negative total force on them in the stopped configurations, and what happens after motion is triggered by one of these seeds as the force is increased. Quite generally, the continuous nature of the distributions of yield strengths and the discrete nature of the stress transfer means that the density of states for local properties should either be zero or be strictly positive. Specifically, from the above we expect that the density of states of the local forces will be positive at zero in the stopped state. To check this, we have computed the probability density per site, $`r(F)`$, that some motion is triggered with a small increase in $`F`$ to $`F+dF`$ in stopped systems; this is normalized so as to include only those samples that have not yet restarted macroscopically (although they might have already had some transient local motion). This rate of triggering of jumps appears to go to a constant as $`F`$ decreases to $`F_c^{}`$ and exhibits relatively weak dependence on $`F`$ as the force is increased; we will present the data later in the paper.

From the data for the distribution of $`F_{c_{\mathrm{hyst}}}^{}`$, which vanishes approximately linearly at $`F_c^{}`$, a constant density of states for triggering motion is perhaps surprising. The reason for this must lie not in the seeds themselves, but in how they grow. In particular, very close to $`F_c^{}`$ local triggering must be less likely to induce restarting than it does at higher forces. In order to understand this, it is necessary to investigate the avalanche dynamics, i.e., the transient motion in response to triggering of one segment. Before considering this in the context of the hysteresis loops and restarting, we analyze avalanches that occur in the approach to depinning from below.

## VI Avalanche Dynamics

In the previous section we have seen that to understand the hysteretic phenomena observed on cycling the force up and down, we need to understand how macroscopic motion starts once it has been triggered by a local instability that leads to one segment jumping. Before studying the case of interest for hysteresis loops, which involves initial conditions that are set by the stopping process, we analyze the behavior as the force is slowly increased from far below the depinning transition starting from more generic initial conditions. We will call this initial depinning. In particular, we are interested in the behavior as the depinning transition is approached from below.

Even though there is no steady state motion in this regime, there can be local, transient motion in response to small increases in $`F`$. Such avalanches will not persist indefinitely for small $`F`$ because the pinning forces in other regions will eventually dominate as long as $`F<F_{c_{\mathrm{aval}}}^{}(M)`$, the — possibly history dependent — critical force on increasing $`F`$.

In this section we investigate the dynamics that result when $`F`$ is increased adiabatically: initially by just enough that one segment moves. This can then trigger other segments to hop forward, while $`F`$ is held fixed until the avalanche stops. The same procedure is then repeated until $`F=F_{c_{\mathrm{aval}}}^{}`$. For an infinite system, this is defined as the force above which the motion persists indefinitely in the absence of any further increase. In finite systems, there are some ambiguities in how it is defined; we choose to define it as the lowest force at which all of the segments move during a single avalanche. The primary quantities of interest below the depinning are the sequence of avalanches and their statistical properties: numbers, sizes, durations, etc. More macroscopic quantities, such as macroscopic responses to a small but non-infinitesimal increase in $`F`$ can be determined by integrating over the properties of the avalanches.

There are various measures of the size of an avalanche. Three of these will be of particular interest. The moment, $`m`$, of an avalanche is defined as the total motion that occurs:

$$m=\underset{𝐱}{}\left[h(𝐱)_{\mathrm{after}}h(𝐱)_{\mathrm{before}}\right];$$
(58)
this is the quantity of primary interest for earthquakes. Alternatively, one can consider the area, $`a`$, (in the two dimensional case of interest): the total number of segments that move at least once during the avalanche. Lastly, is the linear size, $`\mathrm{}`$, of an avalanche, one measure of which is its diameter defined, for example, either via some weighted sum of distances of moving segments from its center, or as the diameter of the smallest circle that will enclose the avalanche.

### A Scaling

For the purely dissipative case, $`M=0`$, the scaling of the avalanches is related to that of the various quantities — $`h`$, $`x`$, $`FF_{c_{\mathrm{aval}}}^{}`$ and $`t`$ — discussed in the context of the moving phase. In particular, if an avalanche has a diameter $`\mathrm{}`$, its area will scale as $`\mathrm{}^{d_f}`$, its duration as $`\mathrm{}^z`$, the typical maximum displacement — change in $`h`$ — as $`\mathrm{}^\zeta `$, and its moment as

$$m\mathrm{}^{d_f+\zeta }.$$
(59)
As long as the dimension is less than the upper critical dimension, $`d_c=4`$ for short range interactions, the avalanches will not be fractal and hence

$$d_f=d$$
(60)
so that the area is a good surrogate, which we will use, for the length scale of an avalanche:

$$\mathrm{}a^{\frac{1}{d}}.$$
(61)
Well below the depinning transition, most avalanches are small. But as the transition is approached from below, larger ones become possible, although the distribution of their diameters is cutoff by a correlation length $`\xi `$ that, in the dissipative limit, diverges at $`F_{c0}`$ as

$$\xi \frac{1}{(F_{c_{\mathrm{aval}}}^{}F)^\nu }$$
(62)
with the same exponent $`\nu `$ as determines the scaling of the physically different characteristic length in the moving phase, the velocity correlation length. At the critical point in the dissipative limit, the distribution of avalanche sizes is a power law:

$$\mathrm{Prob}[\mathrm{diameter}>\mathrm{}]\frac{1}{\mathrm{}^\kappa }$$
(63)
and a similar relation applies for other measures of size; for example, for the area, the exponent is simply changed to $`\kappa /d`$. Near to the critical force, the distribution of the areas of large avalanches has the scaling form

$$\mathrm{Prob}[da]=p(a;F)da\frac{1}{a^{\frac{\kappa }{d}}}𝒫(a/\xi ^d)\frac{da}{a}$$
(64)
where, of the avalanches that occur within a small force interval around $`F`$, $`p(a,F)da`$, is the fraction that of these that whose area is between $`a`$ and $`a+da`$. The scaling function $`𝒫(y\mathrm{})`$ decays rapidly while for $`y0`$, it goes to a constant.

The statistical “tilt” symmetry of the system that was used earlier to yield the scaling law, $`\frac{1}{\nu }=2\zeta `$, can also be used, via relating the polarizability to the avalanche production rate and the distribution of their sizes, to show that for $`M=0`$,

$$\kappa =d\frac{1}{\nu }=d2+\zeta $$
(65)
as derived in Ref. . The one crucial assumption is that the rate, $`r(F)`$, of avalanche production, defined as $`1/\mathrm{\Delta }F`$ times the number of avalanches per unit area of the system as the force is increased by a small amount from $`F`$ to $`F+\mathrm{\Delta }F`$, tends to a finite non-zero constant at the critical force.

We now turn to an analysis of the data for avalanche statistics and properties of the avalanches: for the dissipationless limit, to check the theoretical predictions outlined above; and for non-zero $`M`$, to investigate the effects of stress overshoots.

An easy quantity to measure is the cumulative distribution of all the avalanches as the applied force is increased to $`F_{c_{\mathrm{aval}}}^{}(M)`$. This is given by

$$p_{cum}(a)da_{\mathrm{}}^{F_{c_{\mathrm{aval}}}^{}}𝑑Fp(a;F)r(F)\frac{da}{a}\frac{1}{a^{K_{cum}}}da.$$
(66)
Assuming that $`r(F)`$ approaches a constant as $`FF_{c_{\mathrm{aval}}}^{}`$, the scaling laws for the dissipationless case yield a power law cumulative avalanche area distribution with an exponent of unity,

$$K_{cum}(M=0)=1$$
(67)
independent of the values of the other exponents. This thus provides a good test of the general scaling theory that does not depend on particular predictions for exponents. Although we do not expect the universal aspects of the avalanche statistics to depend on details of the initial conditions, to avoid effects that might arise from smoothening out rough initial conditions, we take the initial configuration to be approximately flat: specifically, the initial $`\{h(𝐱)\}`$ uniformly distributed in the interval $`[0,1]`$.

### B Dissipative limit

We first analyze the data for $`M=0.0`$. In Fig. 15(a), the cumulative avalanche area statistics are shown; a fit to the data on a log-log plot yields an exponent somewhat less than the theoretical expectation of $`K_{cum}=1`$, If we restrict consideration to those avalanches that occur in the region close to the critical force in which most of the activity occurs, specifically, within the applied force region of $`[F_{c_{\mathrm{aval}}}^{}0.1,F_{c_{\mathrm{aval}}}^{}]`$, the apparent exponent is roughly the same, $`K_{cum}=0.85\pm 0.03`$ as shown in Fig. 15(a). Before trying to understand the apparent discrepancy of this with the scaling prediction, we consider the statistics of avalanches that occur in the critical regime, specifically, only those that occur for $`F>F_c^{}(M)`$. The distribution of these also decays as a power of the area, as shown in Fig. 15(b). But the power is much smaller: $`0.37\pm 0.02`$. Because the correlation length that would cutoff the avalanche distribution is of order the system size in this regime, the distribution should essentially be that of critical avalanches with an exponent $`\kappa /d`$. We thus obtain an estimate

$$\kappa (M=0)0.74\pm 0.04$$
(68)
to be compared with the theoretical expectation of $`\kappa =d2+\zeta =d\frac{1}{\nu }`$ equal to $`\zeta `$ in two dimensions. We see that the agreement with our data for $`\zeta `$ is quite good, suggesting that the basic scaling scenario is correct.

There are various possible sources for the substantial discrepancy of the apparent cumulative exponent from unity. As can be seen in Fig. 15(a),(b), the rate, $`r(F)`$ of avalanche production for $`M=0`$ increases sharply as the critical force is approached. Although it does appear to go to a finite constant, as it should on general grounds since $`r(F)`$ is a locally determined property, the precursor increase will bias the cumulative statistics as there are more avalanches produced near the critical point than further away, and these are the ones that have possibilities of being large. This will tend to put more weight in the large avalanche part of the cumulative distribution thereby decreasing the apparent $`K_{cum}`$. \[The extreme limit of this weighting of those near the critical force would just yield the critical exponent $`\kappa /d`$ instead as in Fig. 15(b).\] In particular, if there is a cusp singularity in the rate of avalanche production as $`FF_{c_{\mathrm{aval}}}^{}`$ with an exponent $`\alpha `$, then this would induce a multiplicative correction to the cumulative avalanche size distribution of the form $`1C/a^{\frac{\alpha }{d\nu }}`$ which could complicate interpretation of data.

In general we expect there to be corrections to scaling arising from weakly irrelevant operators at the RG fixed point that governs the depinning critical behavior. In terms of areas, the leading irrelevant eigenvalue $`\theta \frac{4d}{3}`$  would give corrections of the form $`a^{K_{cum}}(1+Ba^{\theta /d})`$ with $`\theta /d\frac{1}{3}`$ in two dimensions. Such corrections would dominate over those from a cusp in the avalanche production rate unless $`\alpha <0.5`$. Indeed the simplest expectation is that the cusp in the avalanche production rate is controlled by exactly this correction exponent,

$$r(F)r(F_{c0})C_r(F_{c0}F)^{\theta \nu }$$
(69)
corresponding to $`\alpha 0.5`$. A fit to this form with $`r(F_{c0})`$ and $`F_{c0}`$ fixed, and $`C_r`$ and $`\theta \nu `$ as the free parameters, yields a lower value of $`\alpha `$. If we take into account a correction to scaling by fitting the log-log plot of the cumulative avalanche area distribution with a form that includes the leading correction (shown in Fig. 15(a)), the cumulative avalanche size exponent is changed from $`0.89\pm 0.02`$ to $`0.83\pm 0.02`$, i.e. in the wrong direction. But if the data are fit with the expected $`K_{cum}=1`$ and an $`a^{1/3}`$ correction, — a fit with the same number of parameters as an undetermined power law with no correction — the inferred correction to scaling amplitude is $`B1.2\pm 0.1`$, not unreasonably large. If this were the actual form of the distribution, some downward curvature at the largest sizes would be expected as seen in the dashed fitting line in Fig. 15(a). A competing tendency, however, is the flattening of the distribution at areas of order the total system area; this is evident in Fig. 15(a) in which data for $`L=64`$ and $`L=128`$ are shown. These effects combine to make the real uncertainties in the exponents here, and probably for other quantities, substantially larger than the apparent uncertainties.

Note, however, that the useful range of length scales available for avalanche data and other quantities for which one is trying to extract infinite system results from finite system data are, because of the crossover when the length scales approach the system size, less than is available for quantities, such as the variance of the critical forces, that are intrinsically properties of finite size systems. Thus we might hope to have more confidence in exponent estimates extracted from such intrinsic finite size properties.

#### 1 Durations, moments, and roughness

From the avalanche data, specifically for the durations of avalanches, one can determine the dynamic exponent $`z`$: with non-fractal avalanches, the duration $`\tau `$ will scale as $`a^{\frac{z}{d}}`$. In Fig. 17(a) we plot the mean duration of avalanches whose area is within a factor of $`\sqrt{2}`$ of $`a`$ as function of $`a`$. The slope of the log-log plot yields

$$z(M=0)1.40\pm 0.02.$$
(70)
As for the avalanche statistics, we can attempt to take into account corrections to scaling by fitting to the form $`\tau (a)\frac{Ca^{\frac{z}{2}}}{1+B_za^{1/3}}`$. This three parameter fit yields $`z=1.52\pm 0.02`$ and $`B_z=0.58\pm 0.05`$; this value of $`z`$ is closer to the value of $`\frac{14}{9}`$ from a naive extrapolation of the $`4ϵ`$ expansion, although $`z`$ is likely to have higher order corrections in $`ϵ`$ even if $`\zeta `$ does not. Ref.  shows that $`z=2\frac{2}{9}ϵ0.04321ϵ^21.38`$. The $`z`$ extracted from the avalanche dynamics agrees quite well with this value without any corrections to scaling.

The roughness exponent $`\zeta `$ can be extracted from the moment as a function of the area: $`ma^{\frac{d+\zeta }{d}}`$. Factor-of-two logarithmically binned histograms are shown in Fig. 21; for the dissipative case these appear to yield a rather small value of

$$\zeta _{aval}(M=0)0.44\pm 0.02$$
(71)
however there is definitely upward curvature observable at the larger sizes and a fit with $`\zeta =\frac{2}{3}`$ and an $`a^{\frac{1}{3}}`$ correction to scaling (same number of parameters) is somewhat better and yields a correction to scaling amplitude of order unity, more precisely $`0.91`$. Note, however, that, as for the distribution of avalanche sizes, such fits suffer from finite size effects.

A roughness exponent can be more directly extracted from the statistics of the roughness of the manifold exactly at the critical force, $`F_{c_{\mathrm{aval}}}^{}`$. In Figs. 19(a)(b), fits to the spatial power spectrum, $`|\widehat{h}(k)|^2`$, yield $`\zeta (M=0)0.65\pm 0.02`$ and $`0.69\pm 0.02`$ for $`f_p^{max}=0.5`$ and $`f_p^{max}=1.0`$ respectively. These values are close to those observed for systems that have stopped after a decrease of $`F`$ from the moving phase.

The apparent discrepancy between the scaling of the typical displacement of avalanches and the roughness at the critical point is somewhat troubling, although, as we have seen, it can readily be accounted for by corrections to scaling. This sort of discrepancy has been seen previously in simulations of manifold depinning . In that case, the roughness exponent at the critical force was similar to ours, but the displacement of the avalanches scaled with a larger exponent. The fact that the discrepancy can be in either direction supports the belief that it is due to corrections to scaling rather than a difference between the two exponents, as had been conjectured in .

As was done for the critical force defined from the moving phase, we can also extract the correlation length exponent from the finite-size scaling of the variance of the critical force. In Fig. 20, the data are shown and exponent estimates of $`\frac{1}{\nu (M=0)}1.22\pm 0.04`$ and $`1.28\pm 0.02`$ extracted from the data for, respectively, $`f_p^{max}=0.5`$ and $`f_p^{max}=1.0`$. These yield

$$\nu (M=0)0.82\pm 0.04$$
(72)
for $`f_p^{max}=0.5`$ and $`\nu (M=0)=0.78\pm 0.02`$ for the stronger pinning. The scaling law $`\zeta +\frac{1}{\nu }=2`$ is consistent with the inferred exponents at the level of a few times the apparent error bars.

In spite of the discrepancies noted above for the typical displacement of avalanches and their cumulative statistics, overall the data for avalanches in the dissipative limit appear to be consistent with those obtained from the moving phase and with theoretical predictions. In particular, if one includes corrections to scaling of the anticipated form, all the data is consistent.

### C Avalanches with stress overshoots

We now turn to avalanche properties for $`M>0`$. Because the distribution of avalanche sizes will turn out to raise the most questions, we first study the duration versus area and moment versus area.

#### 1 Durations and moments

In Figs. 17(a) and (b), the mean and median durations of avalanches whose area is within a factor of $`\sqrt{2}`$ of $`a`$ are plotted for various $`M`$. Except possibly at the largest sizes, the data look remarkably similar to the data for $`M=0`$ and an exponent of $`z1.4`$ would be inferred from each. In Fig. 21 the mean moments are plotted; again there seems to be remarkably little dependence on $`M`$. For both the duration and the moment data, not only do the exponents look similar to those at $`M=0`$, but the amplitudes do also. Given the relatively strong dependence of the velocity amplitude on $`M`$ (Sec. IV), this is particularly surprising.

An amplitude that would be expected to be proportional to the time-length scaling factor, $`A_t`$, can be extracted from the mean duration of avalanches as a function of their area. Although these would be expected to differ by an unknown numerical factor from those extracted from measurements in the moving phase, the ratios might have been expected to be universal. The amplitude $`A_t^{\mathrm{aval}}1.33`$ for $`M<0.4`$ and $`A_t^{\mathrm{aval}}1.44`$ for $`M=0.6,f_p^{max}=0.5`$. It is apparent that these amplitudes exhibit substantially less dependence on parameters than do the $`A_t`$ inferred from the moving phase. Although at first this seems troubling, it should, perhaps, not be: Even in the absence of overshoots, the properties of avalanches can depend on how the system is prepared. Nevertheless, the differences here are probably related to some of the more puzzling differences that occur in distributions of avalanche areas and distributions of critical forces; these we discuss further below.

#### 2 Distributions of avalanche area

In contrast to the properties of avalanches of a given area, the distributions of avalanche sizes depend strongly on $`M`$. In Figs. 15 (a) and 16 (a) the cumulative distributions of all avalanches are shown for $`f_p^{max}=0.5`$ and $`f_p^{max}=1.0`$ respectively: the data are straight on a log-log plot over a substantial range of areas but the slope appears to vary continuously with $`M`$, for $`f_p^{max}=0.5`$ from $`0.89\pm 0.02`$ for $`M=0`$ to $`1.53\pm 0.02`$ for $`M=0.6`$ and for $`f_p^{max}=1.0`$, from $`0.90\pm 0.01`$ for $`M=0`$ to $`1.14\pm 0.02`$ for $`M=1.0`$. At large areas, of order 10% or so of the system area, the distributions of avalanches falls off increasingly more rapidly as $`M`$ increases. Whether this is an intrinsic effect or a finite size effect is an important issue to understand.

Data taken for stronger pinning, $`f_p^{max}=1.0`$, are shown in Fig. 16(a). Surprisingly, these cumulative data depend far less on $`M`$ over the range shown. Although the apparent exponent varies slightly, all are probably consistent with $`K_{\mathrm{cum}}=1`$ with corrections to scaling. More strikingly, there is no dropoff seen for large sizes even with $`M=1.0`$.

The dependence of the avalanche distributions on stress overshoot are more marked if we restrict consideration to avalanches that occur for forces above $`F_c^{}`$, Fig. 15 and 16. As mentioned earlier, for $`M=0`$ these yield an exponent, expected to be equal to $`\kappa /d=\zeta /2`$, of $`0.38\pm .03`$. But as $`M`$ increases the distributions appear to follow roughly the dissipative behavior for small avalanches only to fall substantially below it for large ones. This difference is most pronounced for large $`M`$: the largest avalanches are about an order of magnitude rarer for $`M=0.4`$ than for $`M=0`$ for the weaker pinning and similarly for $`M=1.0`$ versus $`M=0`$ for the stronger pinning. In contrast to the data for the cumulative numbers of avalanches, these data are suggestive of typical crossover behavior: from dissipative behavior for small sizes to something else for large sizes, with a crossover length that is long for small $`M`$ and shrinks as $`M`$ grows. Note that this form of $`M`$ dependence of a crossover length is what would be expected if $`M`$ were a relevant perturbation about the dissipative depinning fixed point. But we must exercise caution before drawing such a conclusion.

From the data of Fig. 15 for $`f_p^{max}=0.5,M=0.6`$ it appears that there may be some interesting size dependence to the crossover: these data for $`L=64`$ fall somewhat below those for $`L=128`$ for large avalanches. One can test whether the crossover is a finite size effect (in contrast to one arising from putative relevance of $`M`$) by plotting the data versus the scaled area, $`a/L^2`$: this is done in Fig. 22(a) for a wide range of system sizes with $`M=0.6`$. It can be seen that the rollover at large areas does not scale simply with system size. But neither is it consistent with a system-size independent crossover: there is much less than the expected factor of 64 difference in the crossover value of $`a/L^2`$ between $`L=16`$ and $`L=64`$. A plausible intermediate conjecture would be that the crossover scales as a power $`1\mathrm{{\rm Y}}`$ of system area: this could arise from dangerous irrelevancy of $`M`$. The plot in Fig. 22(b) versus $`a/L`$ appears roughly consistent with such a conjecture with $`\mathrm{{\rm Y}}0.5`$. But this requires the introduction of an extra scale and exponent into the interpretation of the data. Another interpretation, one that does not require such additional hypotheses. is that systems with $`f_p^{max}=0.5,M=0,6`$ are in the midst of some type of crossover for the range of system sizes investigated. As discussed later, we believe that this is probably the case.

#### 3 Roughness from avalanches

The absence of large avalanches in the presence of substantial stress overshoots for $`f_p^{max}=0.5`$ suggests that, right up to $`F_{c_{\mathrm{aval}}}^{}`$, the manifold should be less rough on long length scales than in the dissipative limit. This is indeed found to be the case: In Figures 19(a) and (b), the spatial power spectrum of the roughness is shown as a function of wavevector just below $`F_{c_{\mathrm{aval}}}^{}`$. In the dissipative limit, a power law that is consistent with a roughness exponent of $`\zeta =0.65\pm 0.02`$ for $`f_p^{max}=0,5`$ and $`\zeta =0.71\pm 0.05`$ for $`f_p^{max}=1.0`$ obtains down to the smallest wavevectors. In contrast, for $`f_p^{max}=0.5`$, the spectra with stress overshoots are found to become flat at long length scales; this is particularly pronounced at the larger $`M`$. Although we have not studied the system size dependence of these spectra in detail, for $`M=0.6`$ the rough magnitude of the wavevector at which the crossover from power-law rough to flat occurs — as seen in the blowups of the small $`|𝐤|`$ regime in Fig. 19(b) — is similar to the wavevector inferred from the inverse square-root of the system area at which the drop-off in the avalanche area distributions occurs. But to again confuse the interpretation, different behavior is once again seen for the stronger pinning, $`f_p^{max}=1,M=0.8`$, samples: these show apparent roughness exponents that are only slightly smaller than in the absence of overshoots and the roughness has less of a tendency to saturate at small wavevectors than the weaker pinning data.

Information can also be garnered from the system size dependence of the average maximum width, $`W_{\mathrm{max}}=<\mathrm{max}(|h(𝐱)\overline{h(𝐱)})|>`$, of the manifold at $`F_{c_{\mathrm{aval}}}^{}`$. In the dissipative case, $`W_{\mathrm{max}}L^\zeta `$ as expected for simple scaling. The stronger pinning samples with $`f_p^{max}=1.0,M=0.8`$ show behavior for $`W_{\mathrm{max}}`$ that is relatively similar to the dissipative case, consistent with the spatial power spectra $`|\widehat{h}(k)^2|`$. For the samples with $`f_p^{max}=0.5,M=0.6`$, however, much weaker size dependence is observed: Fig. 23. We consider several possible reasons that this peculiar behavior for these parameter values might occur.

If, as appears to be the case for samples with $`f_p^{max}=0.5,M=0.6`$, the avalanche production were sharply cutoff at a diameter of order $`L^{1\mathrm{{\rm Y}}}`$ from its distribution in the dissipative case, the largest avalanches that occur would change $`h`$ locally by of order $`L^{(1\mathrm{{\rm Y}})\zeta }`$, yielding a $`W_{\mathrm{max}}`$ of this order. If, instead, the crossover were to a distribution decaying as a larger power-law of the area for large avalanches, the apparent roughness exponent extracted from $`W_{\mathrm{max}}`$ would be somewhat larger. But if the largest avalanches were qualitatively different from dissipative ones —- perhaps like failed nucleation bubbles with $`WL`$ — this could lead to an increased apparent roughness exponent for $`W_{\mathrm{max}}`$. Unfortunately, at this point, it is hard to conclude much from this set of data except that something peculiar seems to be going on for $`f_p^{max}=0.5,M=0.6`$. But peculiarity in this range of parameters appears also in other quantities.

#### 4 History dependence of roughness

The data at $`F_{c_{\mathrm{aval}}}^{}`$ for both the spatial power spectrum of the roughness and for $`W_{\mathrm{max}}`$ should be contrasted with those in the steady state moving phase , Fig. 7(a), and with those in the “just started” moving state that has not reached steady state: Fig. 19(a). The latter two both show power-law roughness with exponents that are consistent with being universal independent of $`M`$ and independent of whether the manifold has been slowly stopped from the moving phase or has barely started moving on initial depinning — albeit with a velocity dependent cutoff at long scales in the latter case. In contrast, the data just below $`F_{c_{aval}}^{}`$ exhibit roughness that depends substantially on $`M`$.

The trends with changing $`M`$ are also different for the stopped samples (at $`F_c^{}`$) and the samples at initial depinning ($`F_{c_{\mathrm{aval}}}^{}`$): When the motion has been stopped from the moving phase, the roughness is slightly larger for larger $`M`$ at least for the weaker pinning. In contrast, the long length scale roughness in systems that started flat and had $`F`$ increased to $`F_{c_{\mathrm{aval}}}^{}`$ tends to decrease with increasing $`M`$ because of the rarity of large system-roughening avalanches with this history.

#### 5 Variations of $`F_{c_{\mathrm{aval}}}^{}`$

We have seen that there appear to be interesting differences between the physics of avalanches in the dissipative limit and those with stress overshoots. In the dissipative case, the critical force, $`F_{c_{\mathrm{aval}}}^{}`$, defined from the first system-sized avalanche, and that, $`F_c^{}`$, defined by when steady-state motion ceases, are essentially the same; the latter being slightly larger and with a somewhat smaller variance because of the way it is defined (the system stops in a somewhat anomalously strongly pinned region). See Table III. With substantial stress overshoots, in contrast, while the lowest force at which a finite size system can have a system-sized avalanche is roughly $`F_c^{}(M)`$, the distribution of $`F_{c_{\mathrm{aval}}}^{}`$ extends far further up, into the region in which, once the motion starts, it is extremely unlikely to stop again. From Fig. 20. it is seen that for $`f_p^{max}=0.5`$ and $`L=128`$, the rms variations of $`F_{c_{\mathrm{aval}}}^{}`$ are about a factor of four larger for $`M=0.4`$ than those for $`M=0`$, and for $`f_p^{max}=1.0`$, they are a factor of ten larger for $`M=0.8`$ than for $`M=0`$, even though the corresponding rms variations of $`F_c^{}`$ hardly differ. This factor of ten difference in $`F=F_c^{}`$ corresponds crudely to a factor of $`10^\nu 6`$ difference in length scale suggesting that some scale of order a sixth of the system size might appear for the initial depinning history with $`f_p^{max}=1.0,M=0.8`$. Recall that such a scale was also apparent in the distribution of the force, $`F_{c_{\mathrm{hyst}}}^{}`$, at which the system restarts after being stopped from the moving phase.

The data for the rms variations of $`F_{c_{\mathrm{aval}}}^{}`$ are shown in Fig. 20 and histograms of the distributions in Figs. 24(a),(b),(c). Note the strong deviations from simple power law behavior of the intermediate $`M`$ data for the variances. Surprisingly, the large $`M`$ samples with $`f_p^{max}=1.0`$ have size dependence of the variance of $`F_{c_{\mathrm{aval}}}^{}`$ that is, except for the overall amplitude, quite like that of the dissipative limit. This again recalls the puzzling behavior that was seen for the distribution of avalanche areas: the strong pinning large $`M`$ systems appear more similar to the dissipative limit than do intermediate $`M`$ ones with weak pinning.

#### 6 Crossover behavior

The data we have presented on the statistics of avalanches and the roughness they induce show complicated dependence on the system size, on the magnitude of the overshoots, and on the strength of the random pinning. But collectively they do suggest a plausible scenario involving a subtle crossover as a function of length scale: For small $`M`$ one would expect dissipative-like behavior for small samples as the effects of small overshoots will need a substantial range of length scales to build up. At the other extreme, the data for large $`M`$ and strong pinning suggest that very large samples show behavior characterized by the same exponents as the dissipative limit, but with very different amplitudes. In between, at an $`M`$ dependent scale, there must then be a crossover from one to the other amplitude. On, for example, a log-log plot of the rms variations of $`F_{c_{\mathrm{aval}}}^{}`$ versus system size this would show up as a regime of slope $`\frac{1}{\nu }`$, followed by a crossover regime of lower slope to an asymptotic regime of slope again $`\frac{1}{\nu }`$. The intermediate $`M`$ samples could then appear to have substantially different scaling by virtue of being in the crossover regime over the range of length scales studied, In the next subsection we see how such crossover behavior might arise from the production and development of avalanches and the effects each has on subsequent ones.

#### 7 Avalanche production and development

We have seen that with stress overshoots, static configurations exist in finite size samples at forces that are substantially above (on the scale of variations of $`F_c^{}(L)`$) the force, $`F_c^{}`$, at which the whole system can move. This suggests that in this regime there may be a change in character of the avalanches as they become large. In contrast to below $`F_c^{}`$ where all avalanches will stop, if an avalanche above $`F_c^{}`$ becomes larger than some characteristic size, naively, perhaps the velocity correlation length, $`\xi _v(F)`$, in the moving phase, it will runaway and the whole system will move. This will strongly affect the distribution of avalanche sizes in this regime, yet it might not have much effect on the durations or moments of those that do not runaway: the evidence discussed above is that the durations and moments of these are indeed rather similar to the large ones that occur with $`M=0`$.

The rates, $`r(F)`$, of avalanche production as $`F`$ is increased adiabatically are shown in Figs. 24(a),(b) and (c) for a variety of values of $`M`$ and $`f_p^{max}`$. These are normalized by the number of samples that are still pinned, i.e. still below their $`F_{c_{\mathrm{aval}}}^{}`$ at the given $`F`$. In contrast to the dissipative limit for which $`r(F)`$ increases sharply with $`F`$ as $`FF_{c_{\mathrm{aval}}}^{}F_c^{}`$, the data for non-zero $`M`$ show a peak in the avalanche production rate somewhat below $`F_c^{}`$ and then a marked decrease as $`F`$ increases into the region in which some samples have self sustaining avalanches and become depinned. For the strongly pinned systems this suppression is strong; the behavior in the region in which $`r(F)`$ is small is shown blownup in the figure. Note that in this regime the avalanche production rate is two orders of magnitude lower for $`M=0.8`$ and $`M=1.0`$ than for $`M=0`$.

Although we have not studied it in detail, this decrease in avalanche production rate is presumably related to the depletion of close-to-unstable sites in regions in which moderately sized avalanches have already occurred; this is analogous to the depletion that occurs after the whole system is stopped from the moving phase. Such depletion will suppress both the density of “seeds” of avalanches and the probability of them becoming large. Most new avalanches would be expected to occur in regions that have not yet had substantial sized avalanches. The exception to this are avalanches that get up enough steam that they can run through the depleted regions. Once an avalanche does this, it will continue growing until it sweeps through the whole system and becomes self sustaining. As the chances of this happening will depend on the properties of various regions of the system, how large it has to become to runaway is likely to be subtle.

There are three possible behaviors for the runaway process of large avalanches in the regime above $`F_c^{}`$: (i) that there is a finite sample-size-independent length scale, $`\xi _R(M)`$ above which avalanches typically run away, (ii) that the scale at which runaway occurs grows with system size, but as a decreasing fraction of system size, or (iii) that the runaway scale is proportional to system size in large samples. In the dissipative limit, the last of these obtains with a proportionality constant close to unity. In contrast, in the artificial model discussed earlier with overshoot stresses felt for the duration of an avalanche, once all segments have moved at least once, an increase of the force by $`M/Z`$ is needed to cause another avalanche and when one does occur it is likely to run away after its diameter becomes of order the velocity correlation length $`\xi _v`$ of the moving phase at that force. Why this is the controlling length scale is explained in the next section.

As discussed in the previous subsection, the peculiar data for avalanche distributions etc, for $`f_p^{max}=0.5,M=0.6`$ appear to indicate that in the regime $`F>F_c^{}`$ the cutoff for the distribution of avalanche areas grows as a power of the area. This would suggest the second scenario as such a crossover in the distribution is presumably associated with a tendency of larger avalanches to runaway. In contrast, the avalanche distributions observed in this regime with strong pinning suggest that there may be a crossover from dissipative-critical-like behavior for smaller sizes to dissipative-cumulative-like behavior for larger sizes, with the crossover occurring for $`\sqrt{a}`$ of order $`0.10.3L`$; if this fixed fraction of $`L`$ persists for larger sizes, it would support the third scenario.

In the next section we develop these ideas further and test them in the context of the restarting of the manifold after it has come to a stop from the moving phase.

## VII Dynamics of Nucleation

In the previous section we studied some aspects of how the system becomes depinned when it has not previously undergone macroscopic motion. We now return to the related question of how the manifold restarts after it has been stopped by a slow decrease of the force from the moving phase.

As discussed earlier, the stopping process leaves behind an at-least-somewhat depleted region of small local forces. This means that as the force is increased back up, the rate of avalanche production, $`r(F)`$, (which in an infinite system is just the probability density of the local forces at zero) will be small. In Fig. 25(d), $`r(F)`$ is plotted as a function of $`F`$ for manifolds that have previously been moving. Even with relatively strong pinning, $`f_p^{max}=1.0`$, it is seen that for $`M=0.8`$ and $`M=1.0`$ the avalanche production rate is almost three orders of magnitude lower than its typical magnitude in a never-moved dissipative system. Yet $`r(F)`$ is relatively flat over the regime in which almost all systems will restart: from $`F_c^{}`$ to $`F_c^{}+M/Z`$. Note, however, that in the upper end of this range, at least for the data we have taken for $`L=128`$, there are no samples still pinned. For smaller samples, the range over which the distribution of $`F_{c_{\mathrm{hyst}}}^{}`$ extends is broader, i.e. from $`F_c^{}`$ to $`F_c^{}+M/Z`$. Note that for all these strong pinning samples, the range over which the restarting occurs is two orders of magnitude larger than the width of the distribution of $`F_c^{}`$.

### A Distribution of restarting avalanches

Although some fraction of the samples of size $`L=128`$ will restart the first time a site is triggered, by no means all will and in general there will be a distribution of avalanche sizes before restarting. In Fig. 16(b) the distribution of the areas of these avalanches are shown for the same parameters as above. The data are somewhat sparse due to the small number of avalanches that typically occur in a given sample, yet some are observed out to about a tenth of the system area. Over the observed range, these avalanche distributions are, perhaps surprisingly, quite close to power laws, with exponents close to $`0.9`$ This should be compared with the much smaller avalanche area distribution exponent $`\kappa /20.39`$ of the critical region avalanches in the dissipative limit. But the measured exponent it is close to that of the cumulative distribution of dissipative avalanches.

### B Bubble nucleation

We conjecture that in the stopped states the length scale above which an avalanche is likely to runaway and restart the overall motion is of order the velocity correlation length, $`\xi _v(F)`$, of the moving phase at the same force. This can be rationalized from the behavior of the moving phase.

At a force that is only a small amount $`f`$ above $`F_c^{}`$,

$$fFF_c^{},$$
(73)
the moving phase is strongly fluctuating on scales smaller than its correlation length, $`\xi _v\frac{1}{f^\nu }`$. On scales smaller than this, there is substantial starting and stopping and the motion is fractal. Given the $`M`$ independence of many of the qualitative properties of the moving phase, the motion on scales smaller than $`\xi _v`$ is probably qualitatively indistinguishable from a finite dissipative avalanche. Thus a reasonable guess is that an avalanche from the stopped state at $`f`$ will be fractal on small scales and is not likely to runaway unless it becomes of order $`\xi _v(f)`$. But if it does become bigger than this, it will “think” that it is part of the moving phase and its local velocity will be unlikely to fluctuate to zero. As its local motion continues, it will impart higher stresses to neighboring regions and its size will grow. With its interior moving at approximately the steady state velocity $`\overline{v}(F)`$, we would expect the diameter of the resulting “bubble” to grow linearly in time.

In order to test this hypothesis, we investigate in our simulations how the manifold begins to move when the applied force is increased after it has been stopped. For $`f_p^{max}=0.5`$ and $`M=0.2`$, it can be seen from Fig. 26 that restarting will typically occur at $`F_{c_{\mathrm{hyst}}}^{}=F_c^{}+M/Z`$ by which point the steady state velocity of the moving phase is already substantial and its correlation length short. In Figs. 28 (a-l) the temporal evolution of the restarting process is shown; we see that it is *not* fractal, but looks more like an expanding nucleation bubble as occurs after supercooling through an equilibrium first-order phase transition. After a short initial transient. there is a front that propagates outward from an initial seed in a roughly deterministic manner. We observe that the radius of the nucleation bubble $`R(t)`$ grows approximately linearly with time so that $`R(t)=ct`$, $`c`$ thus being an expansion rate. The interior of the bubble is, as expected, moving at roughly the steady state velocity, $`\overline{v}(F)`$ of the moving phase. \[Note that $`c`$ is a speed in the spatial rather than that in the displacement $`(h)`$ direction which is $`v`$.\]

The expansion rate can be estimated from the time dependence of the spatially averaged velocity, $`v(t)`$ by assuming a circular bubble of radius $`r(t)=ct`$ so that

$$v(t)\frac{\pi (ct)^2}{L^2}\overline{v}.$$
(74)
In Fig. 27 we thus plot

$$c(t)\left(\frac{v(t)L^2}{\overline{v}\pi t^2}\right)^{\frac{1}{2}}$$
(75)
as a function of the time $`t`$ after the restarting was triggered, for several different $`M`$’s and $`f_p^{max}`$’s. In general, at short times after nucleation, fractal growth that is like that in a dissipative avalanche occurs within an erratically growing region of typical diameter $`\mathrm{}(t)t^{\frac{1}{z}}`$ with resulting cumulative moment of the avalanche: $`m(t)\mathrm{}^{d+\zeta }`$. This results in an apparent $`c(t)\left(\frac{1}{t^2\overline{v}}\frac{dm}{dt}\right)^{\frac{1}{2}}t^{\frac{d+\zeta }{2z}\frac{3}{2}}\frac{1}{t^{0.64}}`$. See Fig. 27 (b). When the size of the avalanche has reached $`\xi _v(F)`$, $`c(t)`$ stops decreasing and the flattish parts of the curves at intermediate times are indicative of approximately uniform expansion; the inferred $`c(t)`$ in this regime is roughly the expansion rate, $`c`$, of the bubble. If the expansion rate is small, as will occur if the steady state velocity is low, then we expect $`c\frac{\xi _v}{\tau _v}\xi _v^{1z}`$ with $`\tau _v`$ the characteristic relaxation time in the moving phase; this is also of order the time at which the crossover from fractal to bubble growth occurs.

At long times, the whole system is moving and $`c(t)\frac{L}{\sqrt{\pi }t}`$ as shown. This saturation occurs at $`tL/c`$. Because the region near the origin of the bubble has already been moving for some time when the regions far away start moving, at the time when the furthest regions start moving, the roughness of the manifold will be greater than it is in steady state The elasticity will diffusively smooth out this roughness at longer times leaving behind only the logarithmic roughness characteristic of the moving phase. This excess roughness near the crossover time from bubble nucleation to steady state motion is the cause of the peak in the width, $`w^2(t)`$, of the manifold shown in Fig. 8. Naively, we expect the peak $`w^2(t)`$ to scale as system diameter, $`L`$, but subtleties associated with elastic slowing down of the velocity inside the bubble might need to be taken into account to understand the process more fully.

For $`f_p^{max}=0.5`$, the forces at which the restarting occurs are substantially above $`F_c^{}`$ and bubble growth appears to occur with an expansion rate for $`M=0.6`$ of $`c0.7`$, close to the maximum possible rate of unity, while for $`M=0.2`$, it is somewhat slower: $`c0.5`$. In the former case, the motion inside of the bubble is essentially alternating sublattices, behavior characteristic of the plateau in the steady state velocity at $`\overline{v}\frac{1}{2}`$ that occurs with substantial overshoots; see Fig. 1(a). For the stronger randomness, $`f_p^{max}=1.0`$, the hysteresis loops are narrower and the nucleation occurs to a lower velocity state. For $`M=0.2`$, the expansion rate is quite small, $`c0.13`$ but the flat constant $`c`$ regime in the figure is observable over a factor of five in time. For these parameters, at early stages the evolution appears fractal, consistent with expectations. For $`M=0.8,f_p^{max}=1.0`$, there is a much narrower flat region and the evolution is not obviously bubble-like. This can be seen in the snapshots of the local motion as the bubble expands shown in Figs. 28. At long times, all of the samples show, as they must, approach of the spatially averaged velocity to its steady state value: this asymptotic behavior is indicated in Fig. 27(a) by the solid line. We see that the time at which $`v(t)`$ becomes close to $`\overline{v}`$ is roughly that that would be expected for a bubble to encompass the whole system and run into its periodic “images”.

More investigation of the dynamics of nucleation, in particular of the transition from fractal on small scales to bubble-like on larger, and the dependence of the crossover scale on system size, on $`FF_c^{}`$, and on the other parameters, is clearly needed.

### C Nucleation avalanche statistics

Armed with the picture of avalanches in the hysteresis loops as fractal on scales out to of order $`\xi _v(F)`$ and bubble-like on larger scales, we can understand the key aspects of the avalanche statistics in this regime. Avalanches that are triggered very close to $`F_c^{}`$ — at forces within the distribution of $`F_c^{}`$ for that system size — will be fractal on all scales and very unlikely to run away. We expect their statistics to be similar to those that occur in the regime near $`F_c^{}`$ when the force is increased from an initially flat configuration: the probability of their area being larger than $`a`$ will decay only as $`a^{\frac{\kappa }{d}}`$ out to the largest sizes. But as the force is increased, the small rate of avalanche production means that there is still a substantial probability that no system sized avalanche will have occurred until a force well outside the distribution of $`F_c^{}`$. In this regime, avalanches that do occur will still have a probability of reaching area $`a`$ that decays as $`a^{\frac{\kappa }{d}}`$ out to scales of order $`\xi _v(F)`$, but this is now considerably smaller than the system size. Those that do reach areas of order $`\xi _v^d`$ will tend to runaway and so the distribution of areas, $`p(a,F)da`$, will be cutoff at this $`F`$-dependent crossover scale. The probability of runaway of an avalanche triggered at $`F`$ is roughly the probability that its area become larger than $`a_v(F)\xi _v(F)^d`$:

$$P_{\mathrm{run}}(F)1_0^{a_v(F)}p(a,F)𝑑a\frac{1}{\xi _v(F)^\kappa }$$
(76)
where $`p(a,F)da`$ is the fraction of avalanches that occur near $`F`$ which have areas within $`da`$ of $`a`$. In order to obtain the cumulative distribution of avalanches in the hysteresis loop we must multiply the distribution $`p(a,F)`$ by the probability that runaway has not yet occurred,

$$P_{\mathrm{pinned}}(F,L)=\mathrm{exp}\left(_{F_c^{}}^Fr(F^{})L^dP_{\mathrm{run}}(F^{})𝑑F^{}\right)$$
(77)
and then integrate this over $`F`$. If the rate of avalanche production per unit area were of order unity, then $`P_{\mathrm{pinned}}`$ would be small for $`FF_c^{}>L^{\frac{1}{\nu }}`$ and the cumulative avalanche distribution would be dominated by avalanches in the regime very near $`F_c^{}`$ for which the distribution is characteristic of the critical point; the cumulative distribution of avalanche areas in the hysteresis loop would then decay with the exponent $`\kappa /d`$ consistent with what we observe for the dissipative case. But with $`r(F)`$ very small for strong pinning and $`M=0.8`$ or $`M=1.0`$ — almost three orders of magnitude smaller than for the dissipative case, (see Fig. 25(d)) $`P_{\mathrm{pinned}}`$ will only become small when $`\xi _v^d<rL^d`$ corresponding to $`FF_c^{}>r^{\frac{1}{d\nu }}L^{\frac{1}{\nu }}`$ a factor of fifty or so further from $`F_c^{}`$ than in the dissipative case, roughly consistent with what is observed.

## VIII Crossovers and Asymptopia

### A Hysteresis loops

The form of the avalanche distribution in hysteresis loops with a very low density of nucleation segments should be apparent from the discussion at the end of the previous section: Most avalanches will take place when $`F`$ is close to the $`F_{c_{\mathrm{hyst}}}^{}`$ at which the system restarts, thus it is the distribution $`p(a,F_{c_{\mathrm{hyst}}}^{})da`$ in this regime that will dominate the cumulative statistics. This results in a distribution of avalanche areas that for small avalanches is critical with exponent $`\frac{\kappa }{d}`$, but for large avalanches decays with the larger cumulative exponent $`K_{\mathrm{cum}}=1`$. With small $`r`$, as occurs in the strong pinning systems with $`M=0.8`$ and $`M=1.0`$, the crossover in a system of linear size $`L`$ occurs at a length scale of order $`Lr^{\frac{1}{d}}`$ or corresponding crossover area of order $`rL^2`$. For the size ranges investigated here — up to $`L=128`$ — this means that in practice almost the full range of avalanche areas will be in the large avalanche tail region with exponent $`K_{\mathrm{cum}}=1`$. But extrapolating to very large systems from these data would be highly misleading: in the limit of large systems, the distribution should scale like the critical dissipative limit over most of the range of avalanche sizes with the rapid decay observed in our data obtaining only for the largest few orders of magnitude of avalanche areas.

The dynamics of a typical runaway event will reflect this crossover. Initially, it will start as a fractal avalanche similar to those in the dissipative limit. But when — and if — an avalanche reaches the crossover length of order $`Lr^{\frac{1}{d}}`$, it will crossover to bubble like expansion with its interior behaving like the moving phase and its radius expanding linearly in time at a rate, $`c`$, that scales as $`(Lr^{\frac{1}{d}})^{1z}`$. As the distribution of $`F_{c_{\mathrm{hyst}}}^{}F_c^{}`$ is broad, however, there will be substantial run-to-run variations in both the crossover size and the rate of growth of the bubble.

The behavior of the size dependence of the distribution of widths of the hysteresis loops, as shown in Figs. 11(a) (b), can also be better understood in this framework. In small systems, the force will typically have to be increased by $`M/Z`$ from $`F_c^{}`$ in order to trigger restarting. In contrast, in very large samples, we expect that the restarting will typically occur at forces $`F=F_c^{}+gL^{\frac{1}{\nu }}`$ which are above $`F_c^{}`$ by of order the width of the distribution of $`F_c^{}`$. Yet because of the strong suppression of the density, $`r(F)`$, of seeds for nucleation of avalanches in the stopped system, in practice the numerical factors, $`g`$ can be very large. The distribution of $`g`$ in the small $`r`$ limit is broad, and has a shape close to a Weibull distribution with characteristic scale $`1/r^{\frac{1}{d\nu }}`$ and shape parameter $`\frac{1}{d\nu }`$ but we do not expect it to be exactly of this form. The reason is due to the breakdown of the argument given in Section V: In sufficiently small subregions of size $`\left(\frac{L}{b}\right)^d`$ of a system of size $`L^d`$, the possible nucleation processes are not independent from one subregion to the next. This is because the diameter of a critical nucleation avalanche when it crosses over to approximately deterministic bubble-like growth is of order $`\xi _v(F)`$ which is of order a fixed — albeit $`r`$ dependent — fraction of $`L`$ at the typical forces at which the first runaway avalanche occurs in a large system. Thus the division into approximately independent subregions will not work for $`b>1/r^{\frac{1}{d}}`$.

Because the hysteresis loops vanish in the limit of large systems, the total number of avalanches per unit area that will occur in a hysteresis loop decreases as a power of system size. This will also occur in the dissipative limit, however in that case the “runaway” avalanches that occur for $`FF_c^{}`$ will not have all that much larger moment than ones that only involve, say, half the system area. But due to the strong suppression of avalanche production after a large event, in the presence of substantial stress overshoots the total number of avalanches will be further reduced in spite of the hysteresis loop becoming wider. The above analysis implies that the total number of avalanches per unit area that occur in the hysteresis loop will be of order:

$$n_{aval}r^{1\frac{1}{d\nu }}\frac{1}{L^{\frac{1}{\nu }}}.$$
(78)
Furthermore, the runaway avalanche will typically run for a very long way — probably exponentially far in $`1/r`$ — before stopping at a displacement, $`h`$, for which the total pinning is anomalously strong.

At this point, it is not clear that the numerical data are consistent with this picture. In particular, it appears that the size dependence of the widths of the hysteresis loops is somewhat weaker than would be expected asymptotically. Nevertheless, we believe that the deviation of the scaling of the width of the hysteresis loops from the $`L^{\frac{1}{\nu }}`$ scaling expected asymptotically is most likely associated with the crossover from the small system to large system behavior. In sufficiently large systems, we conjecture that

$$\mu =\nu .$$
(79)
### B Initial depinning

The complicated avalanche behavior observed when depinning is approached for the first time starting from almost flat initial conditions is probably due to a hybrid of effects. For the strong pinning samples with substantial-sized overshoots, in systems of size $`128^2`$ enough avalanches typically occur before runaway that a strong suppression of the avalanche production rate occurs over the whole system. This suppression enables the force to be increased well past the range in which runaway would occur in the dissipative limit, although not as far as it can for samples that have been stopped from the moving phase. The statistics of the avalanches in this regime will have a similar form to that of avalanches in the hysteresis loops discussed above, with a crossover from critical-like to cumulative-like at a length scale that depends on how strongly $`r(F)`$ has been suppressed.

For the weaker pinning samples, for example with $`f_p^{max}=0.5`$ and $`M=0.6`$, the rate of avalanche production is not suppressed nearly as strongly and it is most likely that the crossover regime is right in the middle of the range of sample sizes studied.

Further systematic tests to investigate this overall scenario would be very useful. In particular, one would like to understand on what aspects of the model, parameters, and history, crucial quantities such as the suppression of avalanche production rate after a large avalanche, depend.

## IX Discussion and Conclusions

In this paper we have introduced and studied numerically and by scaling analyses a simple model of depinning transitions in which there are dynamic stress overshoots caused by rapid motion of segments of the manifold.

### A Macroscopic behavior

We have produced a substantial amount of evidence that the macroscopic behavior of this model with various kinds of generic initial conditions is not hysteretic for a wide range of stress overshoots in spite of the coexistence over a range of the driving force of static configurations and steady state moving configurations. In particular, the critical force, $`F_c(M)`$ is uniquely defined and history independent in the limit of large systems; the distributions of the critical force, $`F_c^{}`$, at which steady state motion ceases, that at which it restarts when the force is increased from such a stopped state, $`F_{c_{\mathrm{hyst}}}^{}`$, and that at which macroscopic motion commences when the force is increased from a flat initial condition for the first time, $`F_{c_{\mathrm{aval}}}^{}`$, are all narrow and converge to $`F_c(M)`$ in the limit of large systems. Furthermore, the steady state velocity, $`\overline{v}(F)`$, is a unique function of the driving force, $`F`$ that vanishes as a universal power of $`FF_c`$. Nevertheless, in principle there are static configurations that could be reached by careful enough control of the dynamics of all the system for all forces up to $`F_{c0}=F_c(M=0)>F_c(M)`$. But if the system were to be started in such a configuration and the force increased by an amount that is arbitrarily small in the limit of a large system, an avalanche would be triggered that would runaway and cause the whole system to start moving.

The universality class of the depinning transition for all $`M`$ less than a multicritical value $`M_c`$ appears to be that of the purely dissipative limit without stress overshoots that has been analyzed by renormalization group methods previously. This universality class is thus much broader than had been conjectured.

The crucial feature of the dynamics on which this macroscopic uniqueness and universality relies is the existence of a non-zero density of nucleation sites for new avalanches after an avalanche has run through a region. These must exist for arbitrarily small increases of the applied force even if their density, $`r`$, is much lower than it would have been in the absence of stress overshoots. If these did not exist — more precisely if there were a depletion layer of local forces after an avalanche— then there would be macroscopic hysteresis with the force needing to be increased by the width of the depletion layer to restart the system once it has stopped after having moved previously.

### B Finite-size effects

In the remainder of this last section we summarize the behavior that occurs in our model for a range of $`M`$ in which there is a strongly suppressed density of nucleation seeds. We also consider briefly under what circumstances similar behavior will obtain; and discuss some of the consequences for a system in which understanding the finite size effects is crucial: earthquakes on a geological fault.

In the dissipative limit of our model, all finite size effects in the vicinity of the depinning transition occur within a range of forces around $`F_c`$ of width of order $`1/L^{\frac{1}{\nu }}`$. In contrast, in the presence of stress overshoots finite-size effects can be important over a much wider range of forces. These effects will be particularly pronounced if the suppression of nucleation seeds for further avalanches is very strong, i.e. if after a system spanning avalanche the density of nucleation segments, $`r`$, is very small.

With $`r`$ small, unless samples are sufficiently large, $`L^d>Z/(Mr)`$ there will usually be no nucleation until the force has been increased by the magnitude, $`M/Z`$, of the overshoots. Unless $`M`$ is very small, the nucleation will then be bubble-like with rapid motion in a linearly expanding bubble occurring as soon as one segment is triggered. For larger systems, the width of hysteresis loops will start to decrease and eventually, if the scenario presented above is correct, decrease as $`1/L^{\frac{1}{\nu }}`$ but with a large coefficient that is a randomly varying multiple of $`1/r^{\frac{1}{d\nu }}`$. Sites that are triggered as $`F`$ is increased after the whole system has been stopped may result in finite avalanches with a distribution that crosses over from one power law to another as a function of their size. Eventually one of these avalanches will runaway: initially it will grow in a fractal manner, but when its diameter reaches of order $`L/r^{\frac{1}{d}}`$, it will began to expand in a more deterministic manner at an approximately constant rate until it covers the whole system and the macroscopic motion restarts. In the limit of small $`r`$ this regime of bubble growth will obtain over a wide range of length scales.

### C Model earthquake dynamics and statistics

The picture presented here for effects of stress overshoots on the dynamics of finite size systems that occur as the driving force is varied have particularly interesting consequences for models of earthquakes on disordered faults. The appropriate driving to model a geological fault is not a constant force, by rather driving by a weak spring. This is roughly equivalent for our case with short range elastic interactions to replacing the applied force by a “pulling spring” with

$$F(t)=\frac{G}{L^2}\left(v_st\overline{h}(t)\right)$$
(80)
where $`G`$ is an effective elastic constant of order unity, the $`1/L^2`$ factor arises from the $`^2`$ elasticity, and

$$\overline{h}(t)=\frac{1}{L^d}\underset{𝐱}{}h(𝐱,t)$$
(81)
is the spatially averaged displacement.

As has been discussed previously,  in the absence of stress overshoots a system driven with an infinitesimal “shearing” velocity, $`v_s`$, “self-organizes” into a statistically steady state with a power law distribution of avalanche diameters that falls off as $`\frac{1}{\mathrm{}^\kappa }\frac{d\mathrm{}}{\mathrm{}}`$ out to length scales of order $`L`$, (the cutoff being affected by the magnitude of $`G`$). There is no particular qualitative distinction between the character of “earthquakes” that have substantial $`\mathrm{}L`$ and those with $`\mathrm{}L`$.

With stress overshoots that give rise to a substantial reduction of the density of seeds for nucleating new events, the behavior is quite different. As the pulling spring advances there will be a string of avalanches until one runs away and the whole system moves. This motion will only stop when the force $`F(t)`$ has decreased as a result of the motion to $`F_c^{}`$, a force that will have variations of order $`1/L^{\frac{1}{\nu }}`$. As the spring slowly restretches, the pulling force will increase again and a series of quakes of various sizes will occur at a rate given by $`rL^d/v_s`$. These events will not relieve much of the accumulating stress and so the force will continue to increase until a force $`F_{c_{\mathrm{hyst}}}^{}`$ at which a runaway event occurs. The typical magnitude of $`F_{c_{\mathrm{hyst}}}^{}F_c^{}`$ and the width of its distribution are both of order $`\frac{1}{r^{\frac{1}{d\nu }}}\frac{1}{L^{\frac{1}{\nu }}}`$, The evolution of the runaway event will at first be fractal, but then a crossover to bubble-like growth will occur and it will grow with constant expansion rate until the whole system is moving.

The initial slip velocity of the large earthquake will be of order $`\overline{v}(F_{c_{\mathrm{hyst}}}^{})`$ which will vary substantially, but have typical magnitude

$$v_R\frac{1}{(r^{\frac{1}{d}}L)^{z\zeta }}$$
(82)
As it runs, the driving force will gradually decrease until it gets down to another $`F_c^{}`$ at which point it will stop. The total displacement in such an event will be of order $`\mathrm{\Delta }hL^2\left(F_{c_{\mathrm{hyst}}}^{}F_c^{}\right)`$ and its moment hence of order

$$m_RL^d\mathrm{\Delta }h\frac{L^{d+\zeta }}{r^{\frac{1}{d\nu }}}$$
(83)
bigger by the $`r`$ dependent factor than the typical largest events in the dissipative limit. \[In obtaining this result, the scaling law $`\zeta +\frac{1}{\nu }=2`$ was used.\]

The distribution of earthquake moments in this model will be a composite of characteristic earthquake and Gutenberg-Richter-like statistics , . Small events will have a power law distribution with exponent $`\kappa /(d+\zeta )`$ that will obtain out to moments of order

$$m_X\left(Lr^{\frac{1}{d}}\right)^{d+\zeta },$$
(84)
much smaller than $`m_R`$ for small $`r`$, The distribution of larger events will fall off more rapidly with exponent $`d/(d+\zeta )`$ out to events with moments of order $`m_0L^{d+\zeta }`$. Events larger than this will not occur unless they are runaway events. As these run for a considerably larger $`\mathrm{\Delta }h`$ before stopping, there will be a strong suppression of the distribution of moments from $`m_0`$ out to of order $`m_R`$, but the runaway events will have a non-trivial distribution.

This picture implies interesting anticorrelations between the occurrences of the largest of the intermediate size and the runaway events: the former are more likely to occur relatively soon after a runaway event rather than just before one because at later stages in the “earthquake cycle”, once an event gets to a small fraction of the system size, it will almost certainly runaway. Some information on the statistics of the largest events can also be inferred from the approximate Weibull form of the distribution of $`F_{c_{\mathrm{hyst}}}^{}`$ which will primarily determine the moments of the largest events.

In order to translate these results to models more appropriate for geological faults, the elasticity needs to be long-ranged, with static stress transfer falling off as

$$J(𝐱,𝐲)\frac{\mathrm{𝟏}}{|𝐱𝐲|^\mathrm{𝟑}}$$
(85)
in the two dimensional case of interest for all but the biggest earthquakes. This modifies scaling laws to $`\zeta +\frac{1}{\nu }=1`$ and $`\kappa =d1+\zeta `$. In the absence of long range correlations in the random properties of the fault, in the dissipative limit the predicted value of $`\zeta `$ is zero, up to logarithms. But with long-range correlated randomness it could be substantially larger.

The consequences of the scenario suggested by our results for geological fault dynamics certainly merits further investigation along these and other lines. But to conclude this paper, we turn to another important issue: how much can be carried over from the stress overshoot model to other types of local dynamic effects? We should remark that this scenario is valid within the context of $`\mu =\nu `$ for very large systems.

### D Other dynamic overshoot effects

At the beginning of this paper we argued that some of the effects of local inertia of a driven, pinned elastic manifold would be qualitatively similar to the effects of stress overshoots. A key question to raise at this point is whether the finite density of seeds for future avalanche events will be left behind in inertial systems. We conjecture that this will be the case if the inertia is not too large, as we found here for stress overshoots, but this bears more careful thought.

In the context of earthquakes, a more important issue is dynamic frictional weakening. An extreme form of frictional weakening was considered in : if once a segment has moved it is always easier for it to move again at any time in the same earthquake, this will certainly give rise to macroscopic hysteresis. As long as there is healing back up to a higher strength between events, a system spanning event will leave behind a configuration that cannot start moving again until the driving force has increased again by enough to overcome the difference between the static and dynamic friction. In reality, the history dependence of frictional forces is complicated and some healing will occur already during an event on those segments that slow down or stop while other parts are still moving. But whether this, combined with dynamic stress overshoots, will leave behind a finite density of nucleation sites for easily triggerable events is a difficult but important question.

Although this paper has perhaps raised as many questions as it has answered, the progress it represents in understanding of the dynamics of driven elastic manifolds in the presence of both randomness and non-dissipative dynamics should help frame and address some of the key questions that remain.

## X Acknowledgments

The authors would like to thank Jim Sethna and Ron Maimon for very useful discussion. This work was supported in part by the National Science Foundation via grants DMR-9976621, 9809363. J.M.S. was also supported by NSF via grants DMR-0109164, 9805818.