Table Of Content

Temporal and spatial heterogeneity in aging colloids: a mesoscopic model Nikolaj Becker and Paolo Sibani FKF, Syddansk Universitet, DK5230 Odense M, Denmark Stefan Boettcher and Skanda Vivek Department of Physics, Emory University, Atlanta, GA, USA A model of dense hard sphere colloids building on simple notions of particle mobility and spatial coherence is presented and shown to reproduce results of experiments and simulations for key quantities such as the intermediate scattering function, the particle mean-square displacement and the χ mobility correlation function. All results are explained by two emerging and interrelated 4 dynamicalproperties: i) arateofintermittentevents, quakes,whichdecreasesastheinverseofthe 4 system age t, leading to µq(tw,t) ∝ log(t/tw) as the average number of quakes occurring between 1 the ‘waiting time’ tw and the current time t; ii) a length scale characterizing correlated domains, 0 whichincreaseslinearlyinlogt. Thisleadstosimpleandaccuratescalingformsexpressedinterms 2 of the single variable, t/t , preferable to the established use of t and of the lag time τ = t−t w w w as variables in two-point correlation functions. Finally, we propose to use χ (t ,t) experimentally n 4 w to extract the growing length scale of an aging colloid and suggest that a suitable scaling of the a J probability density function of particle displacement can experimentally reveal the rate of quakes. 5 PACSnumbers: 82.70.Dd,05.40.-a,64.70.pv 2 ] t Aging is a spontaneous off-equilibrium relaxation pro- thespatial formoftheclustersonalatticeinanydimen- f o cess,whichentailsaslowchangeofthermodynamicaver- sion. Thisenablesustocalculatetheinternalenergy,the s ages. Inamorphousmaterialswithquencheddisorder[1– intermediate scattering function, the mean square dis- . t 6], measurable quantities such as the thermo-remanent placement,themobilitycorrelationfunctionχ andother a 4 m magnetization [7] and the thermal energy [5, 8–10] de- measures of spatial complexity. We can hence compare crease, on average, at a decelerating rate during the ag- with simulational [20] and experimental [21, 22] data, - d ing process. In dense colloidal suspensions, light scatter- and suggest a number of experiments to probe colloidal n ing[11,12]andparticletrackingtechniques[13–16]have dynamics. o uncovered intermittent dynamics and a gradual slowing Our particles reside on a lattice with periodic bound- c down of the rate at which particles move. Intermittency ary conditions, each lattice site occupied by exactly one [ suggest a hierarchical dynamics, instead of coarsening, particle. Particles are either mobile singletons (cluster- 1 as the origin of this process. However, changes in spa- size h = 1) or form immobile contiguous clusters of size v tiallyaveragedquantitiessuchasenergyandparticleden- h > 1. When picked for an update, mobile particles ex- 1 2 sity are difficult to measure and the question of which change position with a randomly selected neighbor and 5 physical properties are actually evolving in an aging col- join that neighbor’s cluster. If the particle is not mobile, 6 loid [17, 18] lacks a definite answer. either its entire cluster “shatters” into h newly mobile . 1 particles with probability In a recent paper [19], two of the authors proposed 0 that kinetic constraints bind colloidal particles together 4 P (h)=e−h, (1) 1 in ‘clusters’. As long as a cluster persists, its center : of mass position remains fixed, on average, but once it v ornoactionistaken. Whenstartingwithaninitialstate breaksdowntheparticleswhichbelongtoitcanmove in- i consisting of singletons, i.e., without any structure, the X dependentlyinspace,andareabletojoin otherclusters. model develops spatially heterogeneous clusters with a r The dynamics is controlled by the probability per unit a length scale growing logarithmically in time. The state of time, P(h), that a cluster of size h collapses through ofthesystemonasquarelatticewithL=256after1015 a quake. Specifically, if the cluster-collapse probability is sweeps is depicted in Fig. 1(a) while Fig. 1(b) shows the exponential, as in Eq. (1) below, quakes follow a Poisson logarithmic growth of the average cluster size. process whose average is proportional to the logarithm As previously shown [19], the rate of events deceler- of time. Log-Poisson processes describe the aging phe- ates as 1/t, see inset of Fig. 1(b), which makes random- nomenology of a wide class of glassy systems [5, 7–10] sequential updates inefficient. In our simulations, we and, specifically in our case, imply that particle motion thereforeusetheWaiting Time Method [23,24], wherea is(nearly)diffusiveonalogarithmictimescale, asfound random“lifetime”isassignedtoeachclusterbasedonthe in our analysis [19] of tracking experiments [14]. geometric distribution associated with P(h); the cluster InthisLetter, weintroduceamodelofclusterdynam- withtheshortestremaininglifetimeisshatteredandlife- ics based on these principles that explicitly accounts for times for other pre-existing or newly formed clusters are 2 40 0.7 30 1.5 i t0.6 hhi20 )1100−04 i heIn0.5 tλ(10−8 nt 10 10−12 eI 1 0.4 1001041081012 h 0.2 0.25 0.3 0.35 t 0 1/plog(t) 100 103 106 109 1012 t FIG. 1: (a) Snapshot of a 256×256 system after t = 1015 0.5 sweeps. Random colors are assigned to different clusters for 0 2 4 6 8 10 12 visibility. Cluster sizes apparently are peaked around some 10 10 10 10 10 10 10 average value (cid:104)h(cid:105) ∼ logt, as shown in (b). The decelerating t rate λ(t)∼1/t of cluster break-up events that emerges from the probability in Eq. (1) is shown in the inset. FIG. 2: Average interface energy per particle, (cid:104)eInt(cid:105), for L= 64. For large t, the interface energy follows the form (cid:104)e (cid:105)∼ Int (cid:112) 1/ log(t) derived in Eq. (2), as confirmed by the inset. adjusted or newly assigned, following the Poisson statistics. Withthisevent-drivenalgorithm,wehavebeenable sesses two-time correlations used to resolve dynamical tofollowourmodelevolutionover15decadesintime,far characteristics of non-equilibrium systems. Formally, it exceeding current experimental time windows. is defined as the spatial Fourier transform, Important aspects of aging dynamics are described by ˆ observablequantitieswithtwotimearguments. Here,we f ((cid:126)q,t ,t)= d(cid:126)rG ((cid:126)r,t ,t)exp[−i(cid:126)q·(cid:126)r], (3) denote by t the current time and by t the waiting time s w s w w before measurements are taken for a system initialized of the self-part of the van Hove distribution function, at t = 0. To conform to common usage, the lag time τ ≡t−tw isusedasabscissainthemainplotofrelevant 1 (cid:88) figures. However, we also provide a collapse of the data, Gs((cid:126)r,tw,t)= N δ[(cid:126)r−∆(cid:126)rj], (4) which is best accomplished when global time t is scaled j by t as independent time variable. w with ∆(cid:126)r (t ,t)=(cid:126)r (t)−(cid:126)r (t ) as displacement of par- Using the details of the Lennard-Jones potential in j w j j w ticlesj inthetimeintervalbetweent andt. Ingeneral, theirmoleculardynamicssimulation,El-Masrietal. [20] w SISF can be interpreted as a measure of the “reciprocal were able to determine the evolution of the internal en- ofmovement”,meaningtheaveragetendencyofparticles ergy of a colloidal system in terms of its pressure. We to stay confined in cages whose size scales with inverse simply monitor the interface between clusters as a proxy magnitude of the wave vector (cid:126)q. Using symmetry and of the internal energy, assuming that a shrinking inter- theintegervaluesofthepositions,thediscreteversionof face indicates a decline in free volume which allows par- SISF reduces to ticles within clusters to relieve their mutual repulsion. The average number of clusters (cid:104)n(cid:105) can be written in (cid:42) N (cid:43) 1 (cid:88) terms of average cluster size (cid:104)h(cid:105) as (cid:104)n(cid:105) = L2/(cid:104)h(cid:105). Since fs(q,tw,t)= N cos((cid:126)q·∆(cid:126)rj) . (5) the average cluster size increases with (cid:104)h(cid:105) ∼ log(t), see j=1 Fig. 1(b), and since for compact clusters in two dimen- (cid:112) Duetospatialisotropy,theSISFisonlyafunctionofthe sions the interface-length scales as S(h) ∝ (cid:104)h(cid:105), the √ magnitude q, with q =2π/L≤q ≤π 2=q . average energy per particle (cid:104)e (cid:105) is estimated as min max Int Figure 3 shows the results of simulating an L = 64 (cid:104)n(cid:105) 1 1 system using 2000 instances and waiting times varying (cid:104)e (cid:105)=S(h) ∼ ∼ √ . (2) Int L2 (cid:112)(cid:104)h(cid:105) logt from 210 to 218 in powers of two. Panel (a) depicts the behavior of SISF as a function of lag time. For large Fig. 2 shows that the approximation holds after a more t ,toagoodapproximationthedatacanberepresented w rapid initial decay. The slow decay matches that of the by a power law of the scaling variable tˆ = log(t/t ), w Lennard-Jones simulations in Fig. 1 of Ref. [20], and it f ∼ Ctˆ−A, where C is a constant and A a positive, s is reminiscent of granular compactification [25], where non-universal exponent, see panel (b) of Fig. 3. Some noisy tapping slowly anneals away excess free volume. curvature remains, nevertheless, and data do not com- The same process drives our cluster growth, although pletely collapse. In contrast, panel (c) achieves excellent densitychangesarenotexplicitlyexpressedinthemodel. scaling and collapse using the form (discussed later) Readilyavailablethroughlightscatteringexperiments, the self-intermediate scattering function (SISF) f as- f ∝exp(cid:8)−Atˆ(cid:0)1−(cid:15)tˆ(cid:1)(cid:9), (6) s s 3 1 (a) 0.18 (b) tw: 210 0.8 0.6 t : 211 0.4 w 0.6 t : 212 w ) 0.2 tw: 213 x ma0.4 100 102 104 106 108 t : 214 f(qs 0.18 t/tw (c) tww: 215 0.6 tw: 216 0.2 0.4 t : 217 w 0.2 tw: 218 10 0 103 10τ6 109 1012 0 Atˆ(11−ǫtˆ) 2 FIG. 4: Persistence on a 64×64 lattice. For times t from √ t =212 to 4t activations are recorded. Within a dominant w w FIG.3: DecayofSISFatq = 2πandsystemsizeL=64 max inert background (white), domains are encircled with darker for t =2k, k =10,11,...,18, using three different forms of w to lighter contours to mark one to ten activations. This mo- independent variable: (a) the lag time τ, (b) the total run bility pattern demonstrates dynamic heterogeneity and the time t scaled by waiting time t , and (c) the correction to w associatedthepreferentialreturnofactivitytothesamesites. tˆ=log(t/t )inEq.(6)(seetextfordetails). Thefittedvalues w A and (cid:15) remain approximately constant across the range of differentt ;withA≈0.1andacorrectionofmerely(cid:15)≈1%. w then over the ensemble. Using (cid:104)·(cid:105) for the ensemble average and |·| for the Euclidean norm, the MSD is written as with (cid:15)≈0.012. Thepower-lawexponent A≈0.1weakly depends on tw and changes systematically by about 40% (cid:42) 1 (cid:88)N (cid:43) overtwodecadesoftw,likelyreflectingtheeffectofhigher ∆r2(tw,t)= N |(cid:126)rj(t)−(cid:126)rj(tw)|2 . (8) order corrections in Eq. (6). Note that A is comparable j=1 to the same exponent found in expensive Lennard-Jones simulations, see Fig. 2 of Ref. [20]. Figure5showstheMSDforasystemofsizeL=64with waiting times t = 2k for k = 10,11...18. In analogy An alternative characterization of immobility, persis- w with Fig. 3, the MSD is plotted vs. three different vari- tence, measures the average fraction of particles that ables: Panel (a) uses the lag time, panel (b) the scaling never move [19, 26]. Conceptually simple and easily ac- variablet/t ,andpanel(c)usesthesametypeofcorrec- cessible in simulations, persistence must be deduced in- w tion as in Eq. (6), that is, ∆r2 ∝ A(cid:48)tˆ(cid:0)1−(cid:15)tˆ(cid:1) with the directly in experiments from the SISF at its peak wave- same(cid:15)=0.012andthe“log-diffusion”constantA(cid:48) ≈0.2. vector. In terms of Eq. (4), persistence is defined as Note that a system aged up to time t has a “plateau” w P(tw,t)=Gs(|(cid:126)r|≤d,tw,t), (7) of inactivity for lag times up to τ ∼tw. These plateaus, often associated with the “caging” of particles [13], are i.e., the fraction of particles that have coordinates at timest andtwith|(cid:126)r(t)−(cid:126)r(t )|≤d,wheredisathresh- old repwresenting the largest dwistance a particle can move (a) 3(b) tw: 210 3 without detection. For small d, the SISF in Eq. (5) re- 2 tw: 211 ducestopersistencefor|(cid:126)q|→∞. Toavoidover-counting t : 212 1 w particles that return to their original position, in simu- t : 213 slainticoentsww.eOounrlyrecsouultnstfpoarrptiecrlseissttehnactehwaivtehbdee=n1acatrievavtierd- ∆2ri2 1000 102 1t0/4tw106 108 tww: 214 tually indistinguishable from those for SISF at qmax in h 3(c) tw: 215 Fig. 3. Figure 4 illustrates the recurrance of quake ac- 1 t : 216 w 2 tivity to sites on a L = 64 lattice for a time interval t : 217 w [t ,4t ] with t = 212. While most particles persist in 1 w w w t : 218 their position, mobility concentrates in areas scattered 0 0 w about the system (“dynamic heterogeneities”[13]) with 10 0 103 10τ6 109 1012 0 A1′tˆ(1−2ǫtˆ) 3 future activity favoring previously mobile sites. The positional variance or mean square displacement (MSD) between times t and t is computed by averag- FIG.5: MSDforarangeofdifferentwaitingtimesanddiffer- w ent choices of independent variable. The system parameters ingthesquaredisplacement,firstoverofallparticlesand are the same as in Fig. 3. 4 Deviations from log(t/t )-scaling are visible at long w 1.5 tw: 20 times in both the MSD and the SISF (or, equivalently, t : 22 the persistence). Interestingly, experimental data show w t : 24 a similar behavior, see Fig. 1 in Ref. [19] and Fig. 4 in w t : 26 Ref. [20]. As shown in the insets, these deviations can 1 w be eliminated by a new scaling variable with a small t : 28 χ4 tw: 210 ((cid:15) ≈ 1%) correction in log(t/tw), whose origin is as fol- w lows: Consider first the fraction of persistent particles 0.5 tw: 212 pn after n quakes (marked white in Fig. 4), and neglect tw: 214 that the size of a quake slowly increases with cluster size t : 216 and that particle hits are not uniformly spread through- w tw: 218 out the system. That persistence curve pn then decays 01 00 102 104 106 108 1010 1012 exponentiallywithn. Averagingpn overthePoissondis- t tribution of n gets FIG. 6: Plot of χ4 as a function of time t for various tw at P(tw,t)=exp(−Aµq(tw,t))=exp(−Atˆ), (10) L=64. Note that the peak-height increases with logt . w see Eq. 6, where A is a small constant and µ (t ,t) ∝ easily removed with t/t as independent variable, see q w w log(t/t ) as the average number of quakes occurring be- Fig. 5(b), leading to the approximate scaling behavior w ∆r2 ∼ log(t/t ) and a reasonable data collapse. The tween tw and t. In reality, spatial heterogeneity means w that quakes increasingly hit the same parts of the sys- residual curvature is removed altogether in panel (c) us- tem,seeFig.4,andthattheireffectonpersistencehence ing the same correction as in Fig. (3), (cid:15) = 0.012, but a constant A(cid:48) ≈0.2 that is unrelated to A. gradually decreases. Heuristically, this effect is accom- modated by correcting the exponent with the O((cid:15))-term To gauge correlated spatial fluctuation and dynamical in µ , as in Eq. 6. More precisely, since all moments of heterogeneity, we consider the 4-point susceptibility [22] q the quaking process can be expressed in terms of µ , the q χ (t ,t)=(cid:104)M(t ,t)2(cid:105)−(cid:104)M(t ,t)(cid:105)2, (9) correction is the first term of a Taylor expansion of the 4 w w w actualexponent. Furthermore,thedependenceofcluster with mobility measure M(tw,t) = (cid:80)jcj(tw,t), size distribution on tw leads to a similar dependence of c (t ,t)=exp{−||(cid:126)r (t)−(cid:126)r (t )|| }, where ||·|| denotes A. The downward curvature of the MSD plotted vs. tˆ j w j j w 1 1 the Manhattan norm. Figure 6 shows χ vs. t for a is analogously explained, thus, the weak curvature seen 4 range of waiting times t . While the data does not al- in our data appears to be a direct consequence of spa- w low a global collapse, the logarithmic increase in peak- tial heterogeneity. The mobility correlation function in height vs. t is clearly visible. This is expected since, by Eq. (9) reveals the presence of a growing length scale in w model construction, heterogeneous dynamics is mainly colloidal systems, here, the average linear cluster size. due to the collapse of clusters, locked in at a typical size Finally, wesuggestthreemeasurementstofurtherelu- h ∼ logt , that re-mobilizes a corresponding number of cidatecolloidaldynamics. Thefirstusestheχ suscepti- w 4 particles at a time t > t , which is reflected in the bilityfunction, aswepresentlydo, toinvestigateagrow- peak w height of the peak. Fig. 6 is reminiscent of results in ing correlation length as a function of the age. The sec- Ref. [21] for experiments and simulations, where χ is ond collects the PDF of fluctuations in particle positions 4 measured for a sequence of equilibrium states prepared over short time intervals of length ∆t = at , where a is w ever closer to jamming, either by increasing density or a small constant, uniformly covering the longer interval decreasing temperature. (t ,2t ). If the rate of intermittent quakes decreases as w w Insummary,ourmodelcoarse-grainsawaythe“in-cage 1/t , the log-Poisson statistic in Eq. (10) predicts PDFs w rattling” of particles while incorporating time intermit- whichareindependent oft ,asshowninRef.[6]forspin- w tencyandspatialheterogeneity. Itsbehavior,whichqual- glasssimulations. Finally,wesuggestthattheslightcur- itativelyaccountsforrelevantexperimentalfindings, can vatureseeninexperimentsforMSDvs.log(t/t )[19]and w be described analytically using the log-Poisson statistics in the tail of the SISF [20] reflects spatial heterogeneity, of cluster collapses [19]. asdiscussedabove. Moredetailedexperimentsmightas- Two-point averages have been plotted versus the lag certain if the correction producing our data collapse has time, τ =t−t , to adhere to the established usage and, a similar effect on experimental data. w ininsets,versusthescaledvariablet/t . Thefirstchoice NB thanks the Physics Department at Emory Univer- w lacks a theoretical basis in the absence of time transla- sityforitshospitality. TheauthorsareindebtedtotheV. tional invariance. The second indicates that the distinc- Kann Rasmussen Foundation for financial support. SB tion between an early dynamical regime τ < t and an is further supported by the NSF through grant DMR- w asymptotic aging regime τ >t is moot. 1207431, and thanks SDU for its hospitality. w 5 [16] R.Candelier,O.Dauchot,andG.Biroli,Phys.Rev.Lett. 102, 088001 (2009). [17] H. G. E. Hentschel, V. Ilyin, N. Makedonska, [1] L. Struik, Physical aging in amorphous polymers and I. Procaccia, and N. Schupper, Phys. Rev. E 75, other materials (Elsevier Science Ltd, New York, 1978). 050404 (2007), URL http://link.aps.org/doi/10. [2] P. Nordblad, P. Svedlindh, L. Lundgren, and L. Sand- 1103/PhysRevE.75.050404. lund, Phys. Rev. B 33, 645 (1986). [18] G. C. Cianci, R. E. Courtland, and E. R. Weeks, Solid [3] H. Rieger, J. Phys. A 26, L615 (1993). State Comm. 139, 599 (2006). [4] W. Kob, F. Sciortino, and P. Tartaglia, Europhys. Lett. [19] S. Boettcher and P. Sibani, J. Phys.: Condens. Matter 49, 590 (2000). 23, 065103 (2011). [5] A.CrisantiandF.Ritort,Europhys.Lett.66,253(2004). [20] D.E.Masri,L.Berthier,andL.Cipelletti,Phys.Rev.E [6] P. Sibani and H. J. Jensen, Europhys. Lett. 69, 563 82, 031503 (2010). (2005). [21] L.Berthier,G.Biroli,J.-P.Bouchaud,L.Cipelletti,D.E. [7] G.G. Kenning, G.F. Rodriguez and R. Orbach, Phys. Masri,D.L’Hote,F.Ladieu,andM.Pierno,Science310, Rev. Lett. 97, 057201 (2006). 1797 (2005). [8] P. Sibani, Eur. Phys. J. B 58, 483 (2007). [22] L. Berthier, Physics 4, 42 (2011). [9] P. Sibani and S. Christiansen, Phys. Rev. E 77 (2008). [23] A.B.Bortz,M.H.Kalos,andJ.L.Lebowitz,Journalof [10] S. Christiansen and P. Sibani, New Journal of Physics Computational Physics 17, 10 (1975). 10, 033013 (2008). [24] J. Dall and P. Sibani, Comp. Phys. Comm. 141, 260 [11] L. Cipelletti, S. Manley, R. C. Ball, and D. A. Weitz, (2001). Phys. Rev. Lett. 84, 2275 (2000). [25] E. Ben-Naim, J. Knight, E. Nowak, H. Jaeger, [12] D.E.Masri,M.Pierno,L.Berthier,andL.Cipelletti,J. and S. Nagel, Physica D: Nonlinear Phenomena Phys.: Condens. Matter 17, S3543 (2005). 123, 380 (1998), ISSN 0167-2789, annual Interna- [13] E. R. Weeks, J. Crocker, A. C. Levitt, A. Schofield, and tional Conference of the Center for Nonlinear Stud- D. Weitz, Science 287, 627 (2000). ies, URL http://www.sciencedirect.com/science/ [14] R. E. Courtland and E. R. Weeks, J. Phys.: Condens. article/pii/S0167278998001365. Matter 15, S359 (2003). [26] R. Pastore, M. P. Ciamarra, A. de Candia, and [15] J.M.Lynch,G.C.Cianci,andE.R.Weeks,Phys.Rev. A. Coniglio, Phys. Rev. Lett. 107, 065703 (2011). E 78, 031410 (2008).