Relaxation in yield stress systems through elastically interacting activated events

We study consequences of long-range elasticity in thermally assisted dynamics of yield stress materials. Within a two-dimensinal mesoscopic model we calculate the mean-square displacement and the dynamical structure factor for tracer particle trajectories. The ballistic regime at short time scales is associated with a compressed exponential decay in the dynamical structure factor, followed by a subdiffusive crossover prior to the onset of diffusion. We relate this crossover to spatiotemporal correlations and thus go beyond established mean field predictions.

Relaxation of the microscopic structure in glasses, and more generally in soft yield stress materials, is a topic of long standing interest and great complexity. Broad ranges of time, energy and length scales are involved, together with non-equilibrium aspects such as aging and a strong dependence on the preparation scheme. As a result, no unique scenario has emerged to describe the relaxation of density fluctuations in systems that, quenched from a liquid into a glassy (solid) state, still display internal dynamics strong enough to produce structural relaxation on a measurable time scale. The complexity of the relaxation is usually quantified by the manner in which it deviates from exponential. In many cases, stretching, corresponding to a broad distribution of relaxation times, is observed. However, the opposite situation of compressed relaxation (i.e., faster than exponential) has emerged in the last years as a new paradigm. In this work, we confirm through the numerical study of a simplified model that this behavior can result from interacting elastic events akin to the shear transformations observed in yield stress solids undergoing external deformation, but now, when such events are thermally activated. The importance of the space-time correlations between events for the overall relaxation picture is highlighted.
A milestone in the experimental analysis of the relaxations processes at hand has been achieved by a series of dynamic light scattering experiments on fractal colloid gels [1][2][3][4][5]. More recently, x-ray photon correlation spectroscopy has been used to study slow dynamics, not only in supercooled liquids [6], colloidal suspensions [7] and gels [7][8][9], but also in hard amorphous materials as metallic glasses [10,11]. A common denominator of these experiments is the decay of the dynamical structure factor as a compressed exponential in time t and scattering vector q: f (q, t) ∼ exp[−(t/τ f ) γ ], with τ f ∼ q −n , n 1 and γ > 1. The observed dynamics was a priori unexpected, not only because of the faster than exponential decay of the correlations but also for the ballistic dynamics contrasting the usual diffusive behavior (τ f ∼ q −2 ) in molecular dynamics simulations of glassy systems (see for example [12]). Regarding the last point, the only exception we know are simulations on a gelformer model [13,14] showing compressed exponential relaxation with a q-dependent γ(q) > 1.
Originally, a heuristic explanation for the observed phenomena was based on the syneresis of a gel: the gel shrinks locally and the inhomogeneity acts as a dipole force with a long range elastic effect [1]. A simple meanfield model approach [15][16][17] further encouraged the view that anomalous relaxation has its origin in elasticity effects. On the other hand, the mean-field approach was reported to fail to describe a q-dependence of γ observed in experiments [5,6,9]. Independently, a phenomenological continuous time Lévy flight model was introduced [5]. It accounts for the discrete nature of the rearrangements and it is able to describe not only the compressed behavior but also the stretched one observed in soft [6,9,18] and hard [10] amorphous materials. However, the model assumes a Poissonian distribution for the waiting times between events to fit the data, which has not been confirmed in measurements and could be questionable. In this letter we propose a novel minimalistic model description for thermally activated relaxation dynamics in yield stress materials. The crucial assumption we make is an elastic response of the surrounding material to the occurrence of local yield events. In a first part, we present the mesoscopic dynamics that we study numerically using a GPU-based parallel implementation [19]. Then, we validate our results against mean-field predictions for elasticity effects in the relaxation. Later, we go beyond meanfield and reveal effects due to correlations in the dynamics that give rise to new interesting phenomena. We conclude with a discussion of our results and their impact on the understanding of recent experimental findings.
The model -Our model for thermally activated yield stress dynamics is inspired by former works on steadily sheared disordered material [20,21]. A bi-dimensional amorphous system is described by a coarse-grained scalar stress field σ(r, t) obeying the over-damped dynamics where µ is the elastic modulus and G(r, r ) the solution for the far field elastic response to a deformed inclusion. The propagator for an infinite system reads in polar coordinates [20,22] G ∞ (r, θ) = 2cos(4θ)/πr 2 (see Fig.1). Note that we make here the simplifying assumption that rearrangements occur in a preferred direction, which should not alter the qualitative features discussed in this work. We discretize Eq. (1) in time and space, on a square lattice (typical linear sizes 2 8 , 2 9 ) with periodic boundary conditions. We set the discrete propagator in Fourier space, where we calculate the convolution [20,21]. The yielding of a site leads to a rearrangement with a local deformation rate given by ∂ t pl (r, t) = n(r, t)ε(r, t)/(2τ ), where τ = 1 is a mechanical relaxation time defining our time scale, ε(r, t) = ±ε 0 integrated over the average duration of an event gives the typical strain caused by a rearrangement, depending only on a sign according to the yielding direction, as in [23]. Here n(r, t) is a binary local "state variable" indicating whether a site is active (n = 1) or not (n = 0).
The local stochastic rules for yielding depend both on the local stress and the global temperature T , through activation (on) and deactivation (off) probabilities: where Γ 0 is an attempt frequency and τ res the typical duration of a restructuring event. All characteristic times are chosen to be equal to τ = 1/Γ 0 = τ res = 1, unless otherwise specified. We implement a fixed absolute value of the local yield stress, σ Y = ±1 (units of stress) [24]. The lack of disorder in the local yield stress leads to steady state dynamics without aging effects. This simplifies enormously the analysis of the dynamics and we ensure that our results are by construction independent of non-equilibrium features and thus, for instance, different from the reported origin of the compressed exponential relaxation on Refs. [13,14].
To establish an analogy with experiments we need to introduce particles in the model. We calculate at each time step the vectorial displacement field u(r, t) associated to the discretized stress field [20] and implement tracer particles, that follow this displacement field with no further interactions, mimicking the underlying particle dynamics [25] (see Fig. 1). For each temperature T , with an appropriate protocol, we reach a unique steady state characterized by stress fluctuations around zero with an approximately Gaussian distribution [26]. The typical strain change due to a rearrangement will be material dependent. To enhance spatial correlations in the system and make the resulting effects more visible we choose here a large value of ε 0 = 1. We checked that qualitatively all phenomena persist for a range of model parameters, although possibly less pronounced.
The observed dynamical features can essentially be separated in two categories: On one hand, part of the phenomenology is purely due to the characteristic spatial decay of the elastic propagator and can easily be captured in mean-field descriptions. On the other hand, the presence of spatio-temporal correlations in the system leads to the prediction of a new dynamical regime. Yet, before going into the discussion of the correlation effects we will show the compatibility of our model with the mean-field predictions of Bouchaud and Pitard [15][16][17].
Elasticity effects -In order to compare our thermal model with the mean-field results, we implement independently an analogous system with the sole difference of having randomly activated sites; a Poissonian rule for activation instead of the stress-dependent rule (2). In this second random model, we control the activity, that is, the number of events per unit time. We measure the mean activity a = 1 N i n i (time average) in the thermal model and plug it in the random model as a parameter to compare equivalent systems [27]. While the thermal model activation rule depends on the evolution of the local stresses and therefore of their distances and relative location to other active sites, the random model does not. This produces a striking difference which is at the background of all the observed phenomenology.
In Fig. 2a we compare for both models the evolution of the diffusion coefficient D(t) = ∆r 2 /(4t), where the mean square displacement of the tracers ∆r 2 on a time window t is averaged both over number of tracers and sliding time t 0 , for distances ∆r i = |r i (t 0 + t) − r i (t 0 )| traveled by each tracer i in the time-lapse (t 0 , t 0 + t]. We observe that the initial ballistic regime (regime I), whose duration is correlated with the persistence τ res of the events, does not depend on the model choice. Also, for both dynamics we find a long time diffusive behavior, only with different values of the diffusion coefficient (regime III). Although there exists a third intermediate sub-diffusive regime in the spatial model (regime II), we will focus first on the two regimes that can be predicted by mean-field considerations.
In analogy with experimental measurements, we calculate the dynamical structure factor, where M is the total number of tracers, and the brackets average over the sliding time window t 0 and the different discretized wave vectors that share the same modulus. From mean-field considerations [15] we expect for the ballistic regime a decay of S(q, t) as a compressed exponential with a dimensionality dependent shape parameter γ 2d = 2 (γ 3d = 3/2) [28]. If we search for the best fit of the data for τ f ∝ q −1 we find indeed γ ≈ 2. But the best collapse of the data is achieved for τ f ∝ q −1.125 yielding a fit with γ = 1.8 (see Fig. 2b). In the diffusive regime, we can collapse the data by plotting S(q, t) as a function of q 2 t and we obtain a pure exponential decay with γ = 1 as expected from the mean-field considerations (see Fig. 2d). Even when we show these results for a particular temperature, they hold for all the range of analyzed temperatures (and further, also in the equivalent random model); only the prefactors α b , α d are T -dependent. Assuming a long range elastic response to the local re-laxation processes, we expect that the displacement field decays as 1/r d−1 , where r is the distance to the event and d the dimensionality of the system. From a mean-field analysis it is expected that the distribution of large particle displacements should decay as P (u) ∝ u −(2d−1)/(d−1) , this is a result of the short range behavior of the elastic propagator, yielding for our 2d study P (u) ∝ u −3 . The probability for small displacements on the other hand should grow as u d−1 , due to the far field effect or one rearrangement when there are no other interfering events in between. The crossover between these two regimes should depend on the density of events, that is the activity a. We confirm these scalings within our simulations for low temperatures (see Fig. 3a). For high temperatures we expect the assumptions of the mean-field description to break down, due to the high density of events that leads to a screening of the large displacements.
Correlations effects -One of the main differences of the thermal model compared to random dynamics is the appearance of sub-diffusion. While the random model changes from a pure ballistic to a pure diffusive behavior for all activation probabilities, a comparable (same activity) thermal model, where spatial correlations are allowed to arise, develops an intermediate sub-diffusive regime for low enough temperatures (regime II in Fig. 2a). To determine the origin of this effect, we first check if the tracer displacements are essentially changed when considering systems with and without spatial correlations. We find that the distribution for the absolute displacements is not altered (see Fig.3a). The change of the dynamics is rather due to negative correlations in the twotime auto-correlation function of the vectorial displacements (Fig.3c). Note that the resolution of the correlation measurement is not sufficient to determine the extension of the sub-diffusive regime, which instead is seen in the intermittent dynamics of the local rearrangement (Fig. 3b). Defining τ ev as the elapsed time between two consecutive activations of the same site, the resulting distribution Ψ(τ ev ) shows a trivial exponential shape in the random case (expected by construction) but a power-law form with an exponential cutoff in the thermal model, well scalable with the mean activity in a master curve Ψ(aτ ev ) ∼ (aτ ev ) −2/3 . We devise a simplified continuous time random walk model [29], where we assume that the tracer particles only move when there is an event close by. In this picture the mean square displacement is given by the number N of events in a time interval t times the typical displacement during jumps. The Laplace transform of N (t) is related to the waiting time distribution asÑ (s) =ψ(s)/(s(1 −ψ(s)) and is shown in the inset of Fig. 3c. We observe two regimes: one proportional to s −2 , this is, the long time diffusion, and a second one proportional to s −5/3 , sub-diffusive. This leads to a prediction of ∆r 2 (t) ∝ t 1/3 that overestimates the exponent for the anomalous diffusion, we expect, due to the rough assumption of dynamical arrest between large jumps. Still, the qualitative picture and the duration of the sub-diffusive regime are captured. We deduce that  the new dynamical regime results from the negative correlations in the displacements combined with the intermittent dynamics for the activity, a feature that we like to call statistical caging. It is the local activity intermittency (power-law distributed) which allows the emergence of a correlated dynamical regime, like the sub-diffusive one observed here, and gives rise to crossovers in time scale that impact any measurement involving them. For instance, coming back to the analysis of S(q, t), we notice that during the sub-diffusive regime, the relaxation time τ r (such that S(q, τ r ) = e −1 ) scales as τ r ∼ q −n with n > 2.
When rescaling time and wavelengths as q n t with the appropriate n, curves corresponding to a time window where D(t) decreases, collapse onto a new master curve S(q, t) = exp(α s (q n t) γ ), now with γ < 1. This anomalous diffusion with a stretched behavior of S(q, t) is not accessible in a mean-field approximation, and can be a relevant observation to interpret experimental results.
Discussion -Despite the strong simplifications we made to derive our model description, we expect the qualitative features to be relevant in experimental systems. We tested our model reproducing mean-field predictions for the distribution of the absolute values of the tracer displacements and related values of the dimensiondependent shape parameters in the decay of the dynamical structure factor S(q, t). We perfectly fit S(q, t) in the ballistic regime with a compressed exponential of shape parameter γ ≈ 1.8 (given τ f ∝ q −1.125 ) a value close to expected mean-field value in 2d γ MF = 2. We have observed (not shown) that we further approach γ ≈ 2 when we address smaller time scales compared to the event duration by increasing τ res , and that this is accompanied by a clearer ballistic (τ f ∝ q −1 ) scaling of the curves. Note that the anomalous structural relaxation does very well coexist with a stretched exponential decay of two time auto-correlations in the local stresses ( Fig. 3 (d)).
We insist that the commonly referenced γ MF = 3/2 is valid in 3d only, and for times smaller than the typical rearrangement duration. We observe that the comparison between experiments and the mean-field prediction is frequently inaccurate in the literature, failing in basic aspects as dimensionality mismatching and/or overlooking the limits of validity of the theory. Interestingly, the observed value γ ≈ 1.8 coincides with experimental results on effectively 2d systems in the high density and small q regime [9]. In that work, a q-dependence of the shape parameter is also reported; our model does not resolve it, since it considers only the far field effect of the rearrangements. Annexed experimental estimations of rearrangement typical size and duration are fundamental to interpret the q-dependence of the measured shape factor and compare with theoretical predictions.
Beyond the mean field results, we predict a new phenomenon evidenced by the behavior of tracer particles, that we call statistical caging. This is different from the traditional atomic caging effect, that happens at much smaller length and time scales and is beyond the scope of mesoscopic approaches. Statistical caging due to correlation between events, leads to a sub-diffusive behavior of the tracer particles, reinforced as temperature decreases, that should result in a second shoulder in the evolution of their mean-square displacement. In this regime, we observe a slightly stretched exponential behavior of S(q, t) with an anomalous relaxation time. We attribute this effect to anti-correlations in the displacement field and power-law intermittency in the relaxation dynamics, shown to be a feature of spatial interactions.