On the critical region of long-range depinning transitions

The depinning transition of elastic interfaces with an elastic interaction kernel decaying as $1/r^{d+\sigma}$ is characterized by critical exponents which continuously vary with $\sigma$. These exponents are expected to be unique and universal, except in the fully coupled ($-d<\sigma\le 0$) limit, where they depend on the"smooth"or"cuspy"nature of the microscopic pinning potential. By accurately comparing the depinning transition for cuspy and smooth potentials in a specially devised depinning model, we explain such peculiar limit in terms of the vanishing of the critical region for smooth potentials, as we decrease $\sigma$ from the short-range ($\sigma \geq 2$) to the fully coupled case. Our results have practical implications for the determination of critical depinning exponents and identification of depinning universality classes in concrete experimental depinning systems with non-local elasticity, such as contact lines of liquids and fractures.

The depinning transition of elastic interfaces with an elastic interaction kernel decaying as 1/r d+σ is characterized by critical exponents which continuously vary with σ. These exponents are expected to be unique and universal, except in the fully coupled (−d < σ ≤ 0) limit, where they depend on the "smooth" or "cuspy" nature of the microscopic pinning potential. By accurately comparing the depinning transition for cuspy and smooth potentials in a specially devised depinning model, we explain such peculiar limit in terms of the vanishing of the critical region for smooth potentials, as we decrease σ from the short-range (σ ≥ 2) to the fully coupled case. Our results have practical implications for the determination of critical depinning exponents and identification of depinning universality classes in concrete experimental depinning systems with non-local elasticity, such as contact lines of liquids and fractures.

I. INTRODUCTION
Many dissipative disordered systems display a collective depinning transition, from an almost static (or inactive) to a sliding (or active) regime at a threshold value of a driving force. Examples range from field driven domain walls in ferromagnetic [1][2][3] or ferroelectric materials [4,5], crack propagation under stress in heterogeneous materials [6,7], contact lines of liquids on a rough substrate [8,9], imbibition of fluids in porous and fractured media [10], reaction fronts in porous media [11], solid-solid friction [12], sheared amorphous solids or yield stress fluids [13], dislocation arrays in sheared crystals [14], current driven vortex lattices in superconductors [15][16][17], skyrmion lattices in ferromagnets [18], to even collective cell migration during wound healing or cancer invasion [19]. The collective nature of the depinning transition is often spectacularly manifested at low temperatures through some kind of "crackling noise" which, well beyond the lab scale, much resembles earthquakes, motivating also their study within the very same framework [20,21].
A very fruitful analogy of this problem with equilibrium phase transitions emerges when considering the driving force as the control parameter and the mean sliding velocity as the order parameter. This analogy has been useful in pointing out directions for seeking universal behaviour, and inspiring new methodologies [22,23]. The analogy has also been useful to point out the relevance of genuine non equilibrium effects [24], and to detect non standard features of the transition [25][26][27]. Among the very different models that can be proposed, the depinning transition of elastic manifolds in random media has become a paradigmatic basic problem as it presents the essential ingredients for a non trivial universal behaviour, together with an advantageous combination of analytical [28] and numerical [29] tractability. Moreover, it is directly relevant for predicting universality classes of various concrete systems where the elastic approximation can be justified, notably magnetic domain walls and contact lines of liquids menisci.
The depinning transition at zero temperature of an overdamped elastic interface in a random potential is continuous, non hysteretic, and occurs at a characteristic threshold force f c . Close enough and above the threshold the mean velocity v of the interface in the direction of the force is well described by the putative depinning law v ∼ (f − f c ) β , with β a non-trivial critical exponent. A divergent correlation length l ∼ (f − f c ) −ν and a divergent correlation time τ ∼ l z characterize the jerky motion as we approach f c from above. Concomitantly, the rough geometry of the interface becomes self-affine with the displacement field growing as u ∼ x ζ for length-scales x below l. Hence v ∼ l ζ−z and β = ν(z − ζ). In this regime, the spatio-temporal fluctuations also display universal behaviour and are controlled by avalanches with a broad distribution of sizes S (and durations T ∼ S z/(1+ζ) ), such that P (S) ∼ S −τ , with τ = 2 − (ζ + 1/ν)/(d + ζ) in the quasistatic limit. The critical exponents can be estimated analytically [30,31] and numerically [32][33][34][35][36][37] to determine the different universality classes. These are determined by d, the range [38][39][40][41][42] or nature [43] of the elastic interactions, the anisotropic [44] or isotropic correlations of the pinning force [45,46], and by the presence of additional (i.e. apart from the pinning force) non-linear terms [24,44,[47][48][49][50]. Boundary [51] or ac driven [52] depinning of elastic interfaces have been also studied. If the so called statistical tilt symmetry holds, only two exponents are needed to fully characterize the depinning universality class. In any case, it is very convenient to consider separately the purely geometric ζ or ν, which do not involve time scaling, from z or β which do.
In order to quantify the universal properties of the depinning transition for a concrete experimental system (or microscopic model with a yet unknown coarse grained dynamics) it is important to determine the critical exponents accurately enough so to be able, at least, to differentiate between different candidate universality classes. Unfortunately, testing the depinning law is in general a rather difficult task experimentally, and in many cases also numerically. On one hand [53], fitting accurately β certainly requires an accurate estimation of the nonuniversal threshold force f c . In that respect it is impor-tant to note that the depinning force f c displays important sample to sample fluctuations in finite systems [54], , where L 0 is the linear size of the system and ν F S ≥ 2/(d + ζ) [30]. The thermodynamic limit L 0 → ∞ is also delicate as the value of f c can be strongly affected by the anisotropic sample aspectratio scaling we keep in such limit [54,55]. On the other hand, even if we are able to get a sharply defined f c , we are faced to the fact that the depinning law is expected to work only in an unknown critical region of size ∼ ∆f crit , such that the asymptotic power law scaling for v fully develops only for (f − f c ) ∆f crit . Knowing roughly ∆f crit is thus fundamental for practical applications of the theory. Little is known however about ∆f crit for the depinning transition, except that it is non-universal. How does the critical region depend on the microscopic shape of the disorder, the range of the elastic interactions or the dimensionality d? Do scaling corrections produce intermediate power-laws with effective exponents? If so, are the effective exponents expected to be larger or smaller than the true ones? Do they violate the expected asymptotic scaling relations among exponents?
In this paper we try to answer some of the above mentioned practical questions by performing numerical simulations on different microscopic models. We study depinning models with isotropic uncorrelated disorder and harmonic long-range elasticity, with elasticity kernel decaying with distance as 1/r d+σ . We vary the range of the elastic interactions from the (σ ≤ 0) fully coupled case to the (σ ≥ 2) short-range cases and compare the critical behaviour of the velocity-force characteristics for two different forms of the microscopic disorder. They are termed the "smooth" case (in which the force originated in the disorder does not have any discontinuities), and the "cuspy" case (in which the force has an abrupt jump at the transition point between different potential basins). For the cuspy potential the extent of the critical region tends to be large and rather independent on the value of σ. For the smooth potential the critical region where universality holds (i.e. where we get the same exponents than in the the cuspy case), decreases by increasing the range of elastic interactions and strictly vanishes, ∆f crit → 0, in the fully coupled limit. In such limit scaling corrections are not anymore "corrections" but control the ultimate asymptotic scaling. This explains the peculiar non-universality of the fully coupled model in the strong pinning phase (i.e. with f c > 0), which displays two different exponents, β = 3/2 and β = 1, for the smooth and cuspy cases, respectively (as it is well known from the equivalent and exactly solvable Prandtl-Tomlinson model). For σ > 0, where a unique value of β is the "right" one [56] , our results show nevertheless the great importance of scaling corrections and the emergence of dangerous effective power-laws which particularly affect the obtention of the asymptotic dynamical exponents β or z as compared with the roughness exponent ζ which is found to be more robust. These corrections are particularly relevant for a successful experimen-tal (and also numerical) identification of depinning universality classes in elastic systems with long-range interactions (0 < σ < 2) and to explain quantitative discrepancies with theory. To arrive to these results we devise a convenient model for comparing the critical behaviour of cuspy and smooth microscopic disorders accurately.

II. GENERALITIES OF THE BASIC MODEL
We model a d-dimensional interface embedded on a d + 1 disordered material as a collection of blocks i = 1, . . . , N , located at sites of a d-dimensional regular lattice, and characterized by a continuous displacement u 1 , . . . , u N in the d + 1 transverse direction. We will assume an overdamped Equation of motion, where the terms on the right hand side represent the sum of the elastic couplings, the disorder, and the uniform and constant pulling force, respectively. The G term in Eq. (1) accounts for the harmonic elastic interactions, with G ij being the spring constant associated with the blocks i and j. In order to model longrange elastic interactions we use G ij = κ/|i − j| d+σ , with the normalizing constant κ used to obtain j G ij = 1 (note that the value of G ii does not influence Eq. (1), and is taken as zero). The G term just described implies the convex elastic energy ij G ij (u i − u j ) 2 /2.
When σ ≥ 2 the elastic kernel represents a short-range elastic interaction, while for −d < σ ≤ 0, it represents the fully coupled case that can be exactly solved using mean field techniques. Periodic boundary conditions can be taken into account by summing the elementary kernel over periodic images of the finite system and by using its Fourier representation to obtain the elastic forces at each step of time integration through a numerically efficient convolution.
The second term in Eq.(1) (the only non-linear term of the equation of motion) accounts for the pinning forces. For the moment we will just assume it is statistically where · · · stands for the average over disorder realizations and ∆(u) is a short-ranged function with ∆(0) measuring the strength of the disorder.
The motion described by Eq.(1), with its convex elastic energy is characterized by a unique critical force f c in the large-size limit [57]. This critical force is important in determining the fate of the system at very long times. If f < f c the system reaches a static solution such thatu i (t) = 0 ∀i. For f > f c it reaches a unique steady-state witḣ u i (t) ≥ 0 univocally defined up to a global time shift, with the mean velocity defined as v ≡ N −1 N j=1u i (t) or, by v ≡ lim t→∞ (N t) −1 N j=1 (u i (t) − u i (0)) thanks to self-averaging. The properties just described are conse- Schematic illustration of the peculiar, non-universal, fully-coupled limit in the strong pinning phase. The velocity exponent β is universal (i.e. independent of the shape of the microscopic potential), from its short-range value for σ ≥ 2, varying continuously to its mean-field value for σ ≤ d/2. In the region σ ≤ 0 the universality is suddenly broken: the smooth potential β (solid red line) does not coincide with its cuspy potential value (solid blue line). How does this affect the critical region for the long-range depinning transitions?
quences of the Middleton theorems [58] which assure that the dynamics described by Eq.(1) with its convex elastic energy converges for f ≥ f c to a unique "Middleton attractor". The properties of this attractor can be exploited to devise smart algorithms to target the critical force and critical configuration in finite samples without solving the true dynamics [59,60]. It also allows to cleanly visualize the convergence towards the steadystate by reparametrizing the time with the system center of mass u(t) ≡ N −1 L d j=1 u i (t) [61]. In this paper we will rely (apart from the unicity of the dynamical attractor) on the general propertyu i (t) ≥ 0, valid in the f ≥ f c steady-state. This will be particularly important in relation to the model discussed in Sec.IV B.
It is also worth remarking here that the so called statistical tilt symmetry (STS) holds for Eq.(1), so ν = 1/(σ − ζ) for d/2 ≤ σ ≤ 2, ν = 1/(2 − ζ) for σ ≥ 2, and ν = 1/2 for σ < d/2 [23]. Therefore only two exponents are needed to fully characterize the depinning universality class. [62] Using the STS, in this work we will consider separately ζ and β, to characterize the universality classes. As we will discuss later this arbitrary separation is nevertheless quite convenient, as ζ (and ν) has a purely geometric origin unlike β (and z) which are related with time-scaling and thus affected by local non-universal bottlenecks. This has important consequences from the practical point of view. For instance, it is easier to get much more accurate values for ζ (or ν) by exploiting different methods which do not involve a true temporal evolution, such as the variant Monte Carlo algorithm [59,60] or the metastable configurations obtained by relaxing a flat configuration below the depinning threshold [63].
For f f c the effect of the disorder can be treated as a perturbation. At first order disorder mimics an effective temperature proportional to v −1 . In this so called "fast flow" regime the interface can be fairly described by the forced Edwards-Wilkinson equation and v ≈ f [64]. For the low velocity critical regime f f c we are interested in, perturbation theory fails. Numerical simulations and the functional renormalization group (FRG) approach applied to Equation (1) teach us that, for σ > 0, the above description uniquely determines the critical depinning exponents of the model. Their values depend on d, and smoothly evolve with decreasing σ, from their short-range values for σ ≥ 2, to the mean field value for σ ≤ d/2 or equivalently d ≥ d c (σ), with d c = 2σ the upper critical dimension [65]. The exact "shape" of the microscopic pinning force F i is believed to be unimportant in many respects. Indeed, FRG tells us that for the model of Equation (1), the bare correlator of the pinning force ∆(u) flows, under coarse graining above the fundamental Larkin scale L c , towards a correlator with a "cuspy" singularity. The existence of such cusp nicely accounts for the existence of a critical force, and also for the existence of avalanches. The fixed point of the renormalization flow equations for the pinning correlator function gives us access to unique values of β ≡ β(d, σ) and ζ ≡ ζ(d, σ), which completely characterize the depinning universality class of our model. These FRG calculations are performed assuming in principle a small separation 1 from the upper critical dimension, = d c (σ) − d, with d c = 2σ. One may thus question their validity for the experimentally relevant case d = 1 for instance. Numerical simulations by Rosso et.al. [36] fairly confirm however the FRG picture for one-dimensional interfaces with short-range elastic interactions (σ ≥ 2), showing its validity for the extreme = 3 case.
The above picture, valid for σ > 0 and in principle any dimension d ≥ 1, sharply contrasts with the σ ≤ 0 fully coupled limit, where the depinning model becomes equivalent to the exactly solvable one-particle Prandtl-Tomlinson model, one of the most popular models in nanotribology [66]. For the strong pinning phase of this model, which has f c > 0 (see e.g. Ref. [67]), the critical behaviour of the velocity v ∼ (f − f c ) β becomes non-universal for different microscopic potentials. On one hand for a smooth random potential V i such that does not have jumps, one has β = 3/2. On the other hand, for a cuspy random potential with force discontinuities β = 1 is obtained, a value that coincides with the mean field value expected from FRG for d ≤ d c (σ). One may thus ask: what is exactly happening in the σ = 0+ limit of the smooth potential case? (See Fig.1)

III. ORGANIZATION OF THE PAPER
To clarify the last issue and discuss its practical consequences we will focus on the d = 1 model with long-range interactions (the d > 1 short range case will be discussed in the appendix). First, in Section IV A we will consider the standard model with continous displacements, and compare the critical behaviour of smooth and cuspy microscopic pinning potentials. This will allow us to illustrate the kind of effects that can be expected in the critical region for σ > 0 in both cases. In Section IV B we will propose an alternative model to compare the same two pinning cases but much more accurately, using discrete displacements and an effective microscopic potential described by traps and suitable transition rates. In Section VI we summarize our results and discuss their practical implications for the study of the depinning transition.

A. Numerical results: Continuous potentials
It is convenient to see the actual results of direct numerical simulations of Eq. 1 to have a first clear picture of the differences that appear between smooth and cuspy pinning potentials.
In concrete, the numerical potentials we used are defined as follows (see Fig. 2) For each site i a potential V i (u i ) is constructed. The generic potential V (u) is constructed piecewise, by dividing the u axis in segments through a set of values a n . In each interval a n -a n+1 (defining a ≡ (a n+1 + a n )/2, and ∆ ≡ a n+1 − a n ) the potential is defined as for the cuspy case, and for the smooth case. Note that even in the smooth case the potential is not analytic, but it has a continuous second derivative, which is enough for our purposes. The separation ∆ between a n and a n+1 is stochastically chosen from a flat distribution between ∆ min = 1 and ∆ max = 2. In Fig. 3(a) we see the value of v as a function of f for the case of cuspy potentials for a few values of σ going from nearest neighbor interaction (σ → ∞) to mean field interaction (σ → 0). The plot in logarithmic scale with respect to the critical force f c (Fig. 4) (fitting in each case the value of f c ) displays a robust critical region in which the β exponent can be defined. The value of β as a function of σ (reported also in Fig. 4) increases when σ moves from large values to σ = 0. Moreover, the actual values of β obtained for different σ accurately fit those known from the literature. [68] This represents the "standard" behavior that is compatible with the analysis using renormalization group techniques.
The results of simulations using smooth potentials are shown in Fig. 3(b). They apparently show larger values of β than the values for cuspy potentials at the same σ ( Fig. 3(a)). For instance, the curve for σ = 1 seems to have a slope close to one, instead of displaying β 0.62 as in the cuspy case.
As we have anticipated, the way out of this conundrum is to realize that the critical region at which the results for smooth potentials should coincide with those of cuspy potentials may be small, and we might not be observing it in Fig. 3(b). The critical region should become observable when plotting the results of Fig. 3(b) in logarithmic scale with respect to the critical force f c . However, the identification of the critical region and the true value of β relies on the accurate determination of the value of f c , and this has to be done at the same time that fitting the value of β, so it is very difficult to get reliable values of β if the critical region is expected to be very small. To overcome this difficulty and try to set this point, we have done simulations in a modified model, that is described in the next section.

B. Numerical results: Discrete pinning potential
The results for the flow curves contained in the previous section (Figs. 3 and 4) suggest that there are strong non-universal effects associated to the form of the pinning potential that is used. These non-universal effects can mask the true critical behavior (that must be independent of the form of the potential for σ > 0), and so they have a great practical importance. Yet to determine accurately the behavior close to the critical force f c we face the problem mentioned at the end of the previous section: The value of the critical force is not known in advance, and it has to be determined during the fitting process itself. This may be quite inconvenient when the extent of the critical region is very small, as a slight un-certainty in the critical force can completely alter the results in this critical region.
In this Section we present a modification of the model used in the previous Section IV A that does not have this drawback, and allows a more precise characterization of the effects we are studying. In addition, it also allows a very precise determination of other exponents of the transition, in particular the dynamical exponent z, something that we are also interested in.
We start by defining this modified model, that was already introduced in the context of thermal creep in [27]. Its main characteristic is that the position of the interface at each spatial location is not continuous (as in the previous section) but discrete. We first define the applied force f app at site i as The interface jumps between successive discrete positions (that are taken completely random, with an average separation of 0.1) when f app i exceeds the local critical force, f th i (for simplicity in the simulations we take the values of f th i as constant: f th = 2.5). The jump itself is considered to be instantaneous, but the transition between consecutive positions does not occur immediately after the critical force is exceeded. Instead, a transition rate is considered. In this version of the model, the cuspy and smooth cases of the model in the previous section can be simulated by considering different forms of the transition rate. In fact, as it will be further discussed below, the cuspy and smooth cases of the previous section differ mainly in the time it takes for a particle that has exceeded the stability limit of one potential well to reach the next equilibrium point at the next well. In the case of cuspy potentials this time is roughly constant, independent of the extent by which the threshold force has been exceeded. In the case of smooth potentials this time goes as ∼ (f − f th ) −1/2 , typical of saddle-node bifurcations. We can model this behavior by assigning a constant transition rate λ ∼ cte to mimic the effect of cuspy potentials (this will be referred to as the "constant rate" case), and a rate that depends on applied force as λ ∼ (f − f th ) 1/2 to simulate the case of smooth potentials (this will be referred to as the "variable rate" case). In the concrete implementation, we consider all unstable sites for which f app i > f th , and calculate an expected time τ i for each site to jump, taken from a Poisson distribution with the corresponding rate λ i . The lowest of all τ i is chosen and this is the site that is actually moved. Time is advanced by this minimum τ i , elastic forces are re-calculated and the process is continued. Average velocity is simply calculated as (u(t) − u(0))/t, in a run over a long time interval t.
The advantage of the discrete potential model is that its quasistatic properties are completely independent of the rate law that is used. In fact, let us consider for instance two avalanches in the system that start from the same configuration but that evolve according to two different rate laws. Since the effect of any forward jump of any portion of the interface is to increase the force over any other site in the system, the avalanches will be exactly the same whatever the rate law is, the only possible difference is in the order and time of activation of different sites, by virtue of the Middleton theorems. This means at once that all static critical exponents such as τ , ζ, ν... must be independent of the rate law. In addition, the critical force will also be independent of the rate law, which is very convenient from the point of view of accuracy of the simulations, as explained previously. The only exponents that can depend on the rate law are those that sense temporal properties of the dynamics. They are the flow exponent β and the dynamical exponent z. The flow curves obtained at constant rate for different values of σ are shown in Fig. 5. There is a clear difference in the flow curves between Figs. 5 and 3 for large values of f , which is a consequence of the details of the models (the discrete pinning model does not have a "fast flow" regime, but a velocity saturation at large forces. However, the values of the β exponent determined from the logarithmic plots (panel (b)) are in excellent agreement with those in Fig. 4, showing that the discrete model with constant rated in fact reproduces the behavior of the continuous cuspy pinning potential.  Fig. 5(b). The continuous straight lines have the same slope than in the previous Figure. case, namely λ ∼ (f −f th ) 1/2 , mimicking the smooth continuous potential case. Comparing Figs. 6(a) and 5(a) we observe again (as between Figs. 3(a) and 3(b)) a clear difference in the overall form of the curves with values of β that look larger in Fig. 6(a) than in Fig. 5(a). In order to quantify this difference in more detail it is necessary to look close to the critical force. It was already mentioned that the great advantage of the discrete pinning model is that the values of critical force for the curves in Fig.  6(a) are the same values as for the curves in Fig. 5(a), which were already fitted to construct Fig. 5(b). Using those values we construct the plot in Fig. 6(b). We note that as far as σ > 0 the curves eventually reach the same exponent than in Fig. 5 when approaching the critical force. However, the force range in which this limiting behavior is obtained shrinks as σ is reduced, and we are actually unable to observe it clearly once σ 0.5. Our conclusion is that the critical region with the same β as in the constant rate case remains finite for all σ > 0, but shrinks as σ is reduced, vanishing for σ ≤ 0, where the mean field value β = 3/2 is re-obtained.
An alternative way to look at the effect just described is to plot the ratio R between the velocities for variable and constant rates (note that this has sense only in the present case in which the critical force is the same for the two different rates). The results contained in Fig. 7 show a remarkable systematic trend. As a function of ∆f ≡ f − f c , the velocity ratio R behaves as R ∼ C + A∆f 1/2 . The value of A is almost independent of the long range interaction exponent σ, but C has a systematic dependence, reducing as σ decreases, and vanishing at σ = 0. For any σ > 0 the finite value of C implies the coincidence of the β values for constant and variable rate. However, as C decreases, the range to observe this coincidence decreases also, and for σ = 0 (i.e., C = 0) the value β = 3/2 is obtained for variable rate, instead of the β = 1 that is obtained for constant rate.
The estimation we have for R allows to quantify the extent of the critical region. We can say that the critical region extends roughly up to the point where C(σ) ∼ A∆f 1/2 . This provides ∆f crit ∼ C(σ) 2 . As C is observed to be roughly linear with σ, we finally obtain ∆f crit ∼ σ 2 .

V. IMPLICATIONS ON THE AVALANCHE DYNAMICS
The discrete pinning potential model is also useful to make an accurate determination of avalanche statistics. In order to analyze it, we did simulations using a quasistatic algorithm, consisting in progressively reducing the value of the external force f in Eq.(4) as dynamics proceeds. Namely, each time a site jump from u i to u i + δ, f is reduced to f − δ/N . In this way, f eventu- ally becomes lower than f c and the dynamic stops. At this point f is increased up to the point in which one site becomes unstable. This allows to analyze individual avalanches in the system at f c . In this way we generated large sequences of avalanches in the system, to which we can calculate duration T , size S (calculated as the difference between i u i after and before the avalanche), and spatial extent L (which is the number of sites that jumped during the avalanche). The statistical behavior of these three quantities allows to determine some important critical exponents of the transition. For instance, the (average) relation between S and L allows to determine the roughness exponent ζ through S ∼ L d+ζ . In addition the relation between T and L determines the dynamical exponent z: T ∼ L z .
In Fig. 8, S vs L is plotted for the "standard" long range case σ = 1, which is particularly relevant for propagating fractures, contact lines of liquids or "magnetically charged" domain walls [39]. By construction of the model, this plot is valid both for constant, and also for variable rates. Although the individual data points are quite scattered, averaging over small intervals of L allows to obtain a good estimation of ζ. The value obtained (ζ 0.39) perfectly coincides within the error bar with the value reported in the literature (see e.g. Ref. [ 40]).
We now focus on the duration T vs spatial extent L relation, determining the dynamical exponent z as T ∼ L z . This result depends on the form of the rates. For constant rate the result is shown in Fig. 9(a), and is consistent with the expected value of z, namely z 0.77 [38,41]. The results in Figs. 8 and 9(a) further support the claim that the constant rate discrete pinning potential is a realization of the cuspy continuous potential case.
The results for T vs L for the variable rate case are presented in Fig. 9(b). The first thing that is observed is that data for individual avalanches (small dots) are much more scattered compared to Fig. 9(a). This be- havior is clearly related to the variable rate: whereas for constant rate the avalanche duration is at most of the order of avalanche size, for variable rate even small avalanches can last for quite long, as a single site may take a very long time to be activated if it is only slightly above the local critical force. The next observation in Fig. 9(b) is that a relation T ∼ L z with z being the same exponent as in Fig. 9(a) is obtained, but only for sufficiently large avalanches. This is in fact consistent with our view that the critical region for the variable rate case is much smaller than for the constant rate case, and only large avalanches display the correct critical T vs L dependence. For small avalanches this dependence deviates towards larger values of T . Interestingly, cracks experiments [6] also show an"excess duration" for small avalanches, suggesting the non-universal effect of smooth microscopic disorder or "variable rate" local instabilities. Note also that for these small avalanches the typical duration depends also on the system size N , an effect that was already discussed in Ref. 69 in the context of the yielding transition.
Since the results for the flow curve support the idea that the size of the critical region shrinks to zero as σ → 0, the question naturally arises of what is the nature of the T vs L dependence and the value of z in this limit, for constant and variable rates. To address this point we generated avalanches in a mean field case (σ = 0), and plot the results in Fig. 10, both for constant and variable rate, for system of different sizes. The results for constant rate are consistent with the standard value of z in mean field, namely z = 1/2. However, results for variable rate sharply deviate of this behavior. We obtain a value of z = 1/4 instead, and also an overall   Fig. 1. The extent of the critical region for smooth potentials (measured by ∆f crit ) reduces continuously to zero as σ = 0 is approached. In the region σ ≤ 0 the non-universal critical exponent β = 3/2 is observed. This raises the possibility (in experiments or numerical simulations) to observe strong corrections to scaling in the smooth potential case, specially at low values of σ, that may induce to adjust effective values of the β exponent (roughly indicated by the hatched region), that are expected to be larger than the true values. dependence on the system size, in such a way that we can write T ∼ L 1/4 N 1/2 . This last relation can in fact be obtained analytically [70].

VI. CONCLUSIONS
In Fig.11 we summarize the picture that emerges from our results. The peculiar breakdown of universality in the σ → 0 limit is explained in terms of a vanishing critical region for smooth potentials (whenever disorder is strong enough to have a finite f c in such limit). We argue that the vanishing of the critical region has practical implications for the analysis of long-range depinning. Since the fully-coupled model has a β = 3/2 exponent for smooth potentials, larger than the universal β = 1 for cuspy potentials, effective values 3/2 > β eff > 1 are plausible to be observed for σ 0, but also for even larger σ we expect an excess, i.e. β eff (σ) > β(σ). The difference between the effective and the right exponents are found to be more important for the dynamical exponents, β and z, than for the geometric exponents, ζ and ν. From a simple microscopic model we have shown that this is related to the competition between the characteristic time τ 1 associated to single particle instabilities, with the time associated to collective instabilities, τ ∼ (f − f c ) −zν , which is roughly controlled by the number of active particles involved in the spreading of correlations at lengths l ∼ (f − f c ) −ν . This competition is made more clear when we analyze avalanche dynamics in the quasistatic limit.
To illustrate the kind of effects we can expect from the non-universal corrections to scaling arising from smooth microscopic pinning potentials, in Fig.12 we show estimations of β and f c from the raw data points of Fig.6 for σ = 1 corresponding to the variable jump rate case discussed in Section IV B. We use two fitting methods which are often used in the literature. The best pairs f c and β are obtained from least squares fits to log v = β log(f − f c ) + cte varying the number N p of points considered, starting from the three lowest values of f . For each case we slowly decrease f c from the lowest value of f . In Fig.12(a) we show that for each N p (labels for each N p are shared with Fig.12(b)) the fit parameter error ∆β displays a minimum for a given value of f c . If we choose the values corresponding to this minimum we obtain the fairly good fit of the data shown in Fig.12(c). If instead we choose the values corresponding to the minimum of the standard deviation χ of the fit we obtain the fit shown in Fig.12(c) (note that in this case, the optimum N p and f c are different than with the previous criteria). The obtained values must be compared with those obtained from the more robust constant rate simulations (Fig. 5) β = 0.62, and f c = 2.37670. It is worth noting that in either case the effective β is larger than its true asymptotic value, which is accurately obtained by using the constant rate discrete model.
Our results motivate a reexamination of the empirical experimental and numerical (smooth potential) long range depinning data analysis. In Ref. [7] the depinning exponent β ≈ 0.8 was directly measured for cracks propagating in an elastic inhomogeneous material. A less direct estimate, also for propagating cracks, can be obtained from the experimental results for the avalanche duration exponent γ = 1.67 reported in Ref. [42]. Using that γ ≈ β + ζ/(1 + ζ) and assuming ζ = 0.39 [40,41] we get β ≈ 0.72. Both values appear to be larger than the ones predicted for the universality class of one dimensional elastic interfaces with σ = 1 long-range elastic couplings and uncorrelated isotropic disorder, where β ≈ 0.63 [41,42] and β ≈ 0.68 [38] were found numerically using "cuspy" or cellular automata lattice models. One can thus argue that the excess in the effective value of β may arise, in part, from the strong corrections to scaling we expect for long-range depinning with smooth microscopic pinning potentials, due to the vanishing of the critical region approaching the fully-coupled limit. In addition to considering a long range interaction in a one dimensional system, there is a second standard way to move towards the mean field limit. This is to consider the short range depinning problem in progressively larger number of dimensions. For short range interactions, the critical dimension of the depinning problem is d c = 4, i.e, for d > d c we expect mean field critical exponents, in particular β = 1. The short range case in d c = 4 corresponds to the σ = 1/2 case of a 1D system. It is then natural to ask if increasing the dimension d we observe the same trends we observe by decreasing the exponent σ, particularly in the flow curves.
We have found that the answer to this question is affirmative. The relevant results are contained in Figs. 13,14, and 15 (which should be directly compared to Figs. 5, 6, and 7). There we present simulations of the discrete pinning potential model, in different number of spatial dimensions, namely d = 1, 2, 3, 4 and 6, with interactions only among nearest neighbor sites (corresponding to σ → ∞). In order to compare different dimensions more easily, and to have a well defined limit as d → ∞, here the elastic interaction G ij (in Eq. 1) is normalized differently, namely we take the value of G ij for neighbor sites as 2/d. The trend we observe as dimension is increased is equivalent to what we have obtained as σ is reduced, in one dimension: For constant rate a robust critical region is obtained and the value of β increases with the number of dimensions, reaching the mean field value β = 1 for d = 4. For variable rates the extent of the critical region is smaller, and is reduced as d increases.
Note however that the critical region does not vanish at the upper critical dimension d = 4: our data are consistent with the critical region vanishing only as d → ∞ where the mean field value (β = 3/2) is fully observed. As in the case of varying σ, here we can estimate the extent of the critical region as a function of d (see Fig.  15) and the result is ∆f crit ∼ d −1 . Ratio between the velocity for variable and constant transition rates. We observe the same trend as it was obtained in 1D long range interacting systems as a function of σ (Fig.  7). The value of C decreases with d as d −1/2