Solution of the dynamics of liquids in the large-dimensional limit

Solution of the dynamics of liquids in the large-dimensional limit

Thibaud Maimbourg LPT, École Normale Supérieure, UMR 8549 CNRS, 24 rue Lhomond, 75005 Paris, France    Jorge Kurchan LPS, École Normale Supérieure, UMR 8550 CNRS, 24 rue Lhomond, 75005 Paris, France    Francesco Zamponi LPT, École Normale Supérieure, UMR 8549 CNRS, 24 rue Lhomond, 75005 Paris, France

We obtain analytic expressions for the time correlation functions of a liquid of spherical particles, exact in the limit of high dimensions . The derivation is long but straightforward: a dynamic virial expansion for which only the first two terms survive, followed by a change to generalized spherical coordinates in the dynamic variables leading to saddle-point evaluation of integrals for large . The problem is thus mapped onto a one-dimensional diffusion in a perturbed harmonic potential with colored noise. At high density, an ergodicity-breaking glass transition is found. In this regime, our results agree with thermodynamics, consistently with the general Random First Order Transition scenario. The glass transition density is higher than the best known lower bound for hard sphere packings in large . Because our calculation is, if not rigorous, elementary, an improvement in the bound for sphere packings in large dimensions is at hand.

Introduction –

The physics of liquids and glasses belongs to the group of fields that are victims of the lack of a small parameter. Many approximations have been proposed over the years, but they suffer from the uncertainty about what is the limit in which they are supposed to become exact. This has been true both for equilibrium and for dynamic properties. From the point of view of dynamics, an extreme case is that of Mode-Coupling Theory (MCT) Götze (2009, 1998); Reichman and Charbonneau (2005): it may be introduced by an (uncontrolled) resummation of an infinite subset of diagrams. The Mode-Coupling approximation yields Mode-Coupling dynamics: the phenomenology depends on the approximation itself Kob and Schilling (1991), somewhat like a harmonic approximation is expected to predict oscillations.

An often used remedy to the absence of a small parameter Witten (1980) is to promote the system to dimensions, solve the large limit, and (eventually) expand around. This strategy has been used with success for liquids Wyler et al. (1987); Elskens and Frisch (1988); Frisch and Percus (1999), strongly coupled electrons Georges et al. (1996), atomic physics Svidzinsky et al. (2014), gauge field theory Drouffe et al. (1979), and most recently the thermodynamics of amorphous systems Parisi and Zamponi (2010); Charbonneau et al. (2014a). In this paper we extend this procedure to the dynamics of liquids made of spherical particles. We restrict ourselves to equilibrium, although the extension to the glassy off-equilibrium “aging” regime is at hand. It has been a long-standing question whether MCT becomes exact in infinite dimension Kirkpatrick and Wolynes (1987); Kob and Schilling (1991); Schmid and Schilling (2010); Ikeda and Miyazaki (2010); Jacquin and van Wijland (2011); Charbonneau et al. (2011), and the present computation gives an answer.

Statement of the main result –

We shall consider a system of identical particles, interacting via a spherical potential of typical interaction length in dimensions, and obtain a solution for the equilibrium time correlations of the resulting liquid, that becomes exact in the limit . We need to confine the particles in a finite volume . It is very convenient to do this in such a way that the “box” does not break rotational and translational invariances, which are crucial in our developments. A practical way to do this is to consider particles living on points () on the -dimensional surface of a hypersphere with , (. The thermodynamic limit with constant density , in which the flat space is recovered, will be taken before . Rotations and translations in dimensions – with dimensions and , respectively – are encoded in the rotations in dimensions, with dimension . We consider a Langevin dynamics


where is a white noise with , and is the temperature. Here, and in what follows, denotes average over noise and/or initial conditions. The are Lagrange multipliers, imposing the spherical constraints . For , they do not fluctuate and their value is proportional to the equilibrium reduced pressure  Hansen and McDonald (1986). We shall in what follows treat the overdamped case , but the inertial term may be reinstalled at any stage (and in that case the thermal bath may be disconnected setting , to recover a purely Newtonian dynamics).

We define the adimensional scaled correlation and response functions (see Fig.1 in the Appendix):


where is an external field conjugated to , and we have also introduced the single particle mean square displacement  Hansen and McDonald (1986). As we shall see, these quantities have a finite and non-fluctuating limit for . In the following and in the Appendix, we shall derive exact equations for the correlations (2). For simplicity, we will restrict to equilibrium where , and ; here is the Heaviside step function.

Our main result is that correlation functions are obtained in terms of a memory kernel through the following equations:


the second relation being valid for .

To give the expression of the memory kernel , let us define a scaled density where is the volume of a -dimensional sphere of radius , a scaled , a scaled potential , and the inter-particle force . is then obtained in terms of a one-dimensional effective dynamics with a colored noise ,


from the self-consistent equation:


The average in (5) is over the process (4) with . The total effective potential has the form (see Fig.4 in the Appendix):


The quadratic part of the potential plays the role of the confining “box”. In fact, it is negligible for finite times and large , and the probability distribution is exponential in near (expressing the growth of entropy as a function of distance along the -dimensional sphere), the region relevant for those times. Reassuringly, box details are irrelevant at short times.

Finally, Eq. (5) has a microscopic counterpart (cf. Appendix); in fact, the function is, in the large limit, the autocorrelation of the inter-particle forces ,


a result that provides a physical interpretation of .

Relation with MCT –

If were a simple function of , , then Eq. (3) would be in the schematic MCT form. As is well known, schematic MCT is obtained as the exact dynamics of a system of spherical spins with -spin random interactions Cugliandolo (2003); Bouchaud et al. (World Scientific, 1997); Castellani and Cavagna (2005), for which


However, as soon as one considers non-spherical variables, e.g. soft-spins with a potential , one obtains an equation like (4) with this  Bouchaud et al. (World Scientific, 1997); Sompolinsky and Zippelius (1982); Cugliandolo and Kurchan (1994):


Here again, Eq. (8) holds and the system is closed by . Within the liquid phase, this more general form of dynamic equation has essentially the same phenomenology as schematic MCT. Our system of equations belongs to this more general class, and they thus show exactly the same MCT phenomenology for what concerns universal quantities that are independent of details of the memory kernel (e.g. the dynamical scaling forms and the relations between critical exponents) Andreanov et al. (2009); Janssen et al. (2014); Reichman and Charbonneau (2005).

However, important quantitative differences are observed with respect to applying the MCT approximation to the intermediate scattering function, which leads to the “standard” formulation of MCT for liquids Götze (1998, 2009). Standard MCT has the same qualitative structure as schematic MCT, but also provides quantitative results for the self and collective scattering functions in all dimensions, in particular in  Kob and Andersen (1995a, b); its limit was discussed in Refs. Ikeda and Miyazaki (2010); Schmid and Schilling (2010). Our result in is formulated in terms of , and most of the other natural observables are functionals of . For example, for , we have for the self intermediate scattering function (cf. Appendix):


in contrast to the non-Gaussianity in one finds within MCT close to the plateau Ikeda and Miyazaki (2010); Schmid and Schilling (2010). One could then write our equations in terms of . The result, however, is different from standard MCT and in particular our kernel is not an analytic function of .

Because our equations fall in the same universality class as schematic MCT, but provide different quantitative results with respect to standard MCT in , it remains a matter of taste if one wishes to call them with the same name or, more generally, “dynamic Random First Order Transition (RFOT)” Kirkpatrick and Wolynes (1987).

Sketch of the derivation –

Let us outline the main steps in the derivation; more details are given in the Appendix. In order to construct the high-dimensional limit, a virial expansion is a reliable method. Following Mari and Kurchan (2011), we exploit the well-known analogy between trajectories and polymers. The dynamics are generated by a sum of trajectories in -dimensional space, in our case, the surface of a -dimensional sphere. To each trajectory is associated an Onsager-Machlup probability weight which is the exponential of an action. The sum over all trajectories of this quantity is analogous to a partition function and is thus used to generate averages over the Langevin process (1). The action is expressed in terms of trajectories of the and auxiliary “response” variables (Martin-Siggia-Rose-De Dominicis-Janssen generating path integral) Cugliandolo (2003); Castellani and Cavagna (2005). The result reads, in the It convention:


Following standard liquid theory Hansen and McDonald (1986); Mari and Kurchan (2011), we introduce the density function for trajectories


where is the functional Dirac (a product of deltas over all times) and we construct a virial (Mayer) expansion as a power series in . One can show that all terms involving a product of more than two density fields are subleading for  Frisch and Percus (1999): they are exponentially suppressed because of the requirement that the three trajectories overlap (see Fig.1 in the Appendix), which is exponentially unlikely in large . Truncating the virial expansion accordingly, we get (cf. Appendix):


where and . The physical is determined by and the normalization . The first term in Eq. (13) is an ideal gas contribution and the second accounts for interactions.

Following the thermodynamic treatment Kurchan et al. (2012), we may now argue that due to rotational invariance on the hypersphere, where and . We can thus make a change of variables in the functional integration over to . The change of variables gives for density averages (cf. Appendix):


where is the Jacobian of the transformation (cf. Appendix) and we defined . The appearance of the dimension in the exponent leads to a narrowing of fluctuations of correlations, and saddle-point evaluation becomes exact (cf. Appendix). In this way we can compute the ideal gas term in Eq. (13). For the interaction term, that involves two functions, we need the variables corresponding to , and also , . The Jacobian may be calculated with the same methods (cf. Appendix), and the crucial result is that at the saddle-point level with determined by the same saddle-point as in Eq. (14), while the remaining integration over is effectively one-dimensional. Changing variables with leads to a finite integration over that eventually gives Eq. (4). Hence, the typical distance between two trajectories turns out to be . This scaling physically tells us that a particle vibrate inside a cage with amplitude and interacts with neighbors. Finally, equilibrium and causality at the saddle-point level imply and the Fluctuation-Dissipation relation (FDT) . The saddle-point evaluation for the two-time variables is the mathematical justification of the above-mentioned fact that do not fluctuate. These saddle-point equations give Eq. (5) and (3) (cf. Appendix).

In this paper we are treating an equilibrium situation. Within the liquid phase, this may be achieved by starting from any configuration in the distant past. A more practical way, however, is to assume equilibrium at a convenient time (e.g. ). How does one deal with a non-Markovian equation of motion like (4)? The answer is simple: either one makes the memory kernel extend to the remote past, or, alternatively, one may assume equilibrium at , in other words summing all the past histories passing through at . It turns out Hänggi (1997) that this is implemented simply by cutting the memory at a lower limit , as in Eq. (4) (cf. Appendix). This completes the derivation of our basic dynamical equations.

Dynamic transition and timescale separation –

We now apply the standard MCT methodology to locate the density or temperature at which a dynamic transition occurs with freezing in a cage, corresponding to the development of a plateau in  Götze (2009, 1998); Reichman and Charbonneau (2005); Castellani and Cavagna (2005).

Consider the case when falls from to a plateau value , and then, at much larger times, to zero. Concomitantly, grows to a plateau value , and then continues to grow at a slower (diffusive) pace. Denote the fast part . In the limit in which the plateau times are much larger than the microscopic times, and much smaller than the final relaxation times, the noise breaks into a fast variable and a slow random variable , as does the friction term. Their sum acts as an adiabatically slow field at those times. We may thus split the equilibration in two steps Cugliandolo and Kurchan (2000): and . When is in the plateau region, we may write the expectations:


where , , and


Obviously, . We therefore obtain the self-consistent equation for :


From Eq. (3) we obtain


The dynamical transition point, at which the plateau becomes infinite, happens when Eq. (17) first has a non-zero solution for . This point happens as a bifurcation, and may be quickly obtained by solving Eq. (17) together with (cf. Appendix).

For hard spheres, the result is . This result is fully consistent with the one based on thermodynamics Parisi and Zamponi (2010), and in fact Eq. (17) is exactly identical to the one that can be derived using the replica method (cf. Appendix), consistently with the general RFOT scenario Kirkpatrick and Wolynes (1987); Kirkpatrick and Thirumalai (1988, 1989); Kirkpatrick et al. (1989); Lubchenko and Wolynes (2007); Wolynes and Lubchenko (2012); Kirkpatrick and Thirumalai (2014). This is a particular instance of a general correspondence between thermodynamic and dynamic results that is verified by the infinite- solution, and can be extended to critical MCT exponents Caltagirone et al. (2012) and to correlation functions Franz et al. (2011); Parisi and Rizzo (2013); Franz et al. (2013). In fact, expanding around  Götze (2009, 1998); Reichman and Charbonneau (2005), one can show that upon approaching the plateau, while upon leaving the plateau. The exponents satisfy the famous relation Götze (2009, 1998); Reichman and Charbonneau (2005); Caltagirone et al. (2012); Parisi and Rizzo (2013):


For hard spheres, we obtain which implies and  Kurchan et al. (2013). Finally, from Eq. (10) one can show that the factorization property of MCT Götze (1998, 2009) still holds in , namely that close to the plateau, factorizes in a function of and a function of (cf. Appendix).

Diffusion, viscosity, Stokes-Einstein relation –

At long times, in the liquid phase , the motion is diffusive and , where is the diffusion coefficient. Plugging this form in Eq. (3), and recalling that decays to zero over a finite time, we obtain an exact result for :


At low density and we recover the bare diffusion coefficient . Upon increasing density, increases and decreases. For , where a finite plateau of emerges, diverges and the diffusion coefficient vanishes as with the exponent , which is consistent with the numerical results of Charbonneau et al. (2014b).

Within linear response theory, the shear viscosity of the liquid is given by Hansen and McDonald (1986); Yoshino (2012)


where are two arbitrary components of the stress tensor . In the Appendix we show that , and thus


This relation shows that for , , as it is found in MCT.

Combining Eqs. (20) and (22) we obtain a relation similar to the Stokes-Einstein relation (SER)


where the second expression holds close to . This relation is interesting because it predicts that the SER is not exactly satisfied in dense liquids: the quantity , which is constant in SER, has instead a small variation proportional to . This is in agreement with results of (Charbonneau et al., 2013, Fig.7b) that show a linear variation of with in the dense regime for large enough dimension.

Conclusions –

In this work we obtained an exact solution of the dynamics of liquids in the limit of infinite spatial dimension. The picture that emerges has a relevant “caging” lengthscale that is smaller than the particle radius by (it is about for real colloids Weeks and Weitz (2002)). The physics of diffusion stems from interactions at that small scale, while all particle motion beyond that scale consists of uncorrelated steps of displacement, and memory of what happens at distances is lost. Diffusion coefficients and viscosity are thus decided at distances much smaller than the particle radius. This strongly suggests that the wavevectors that matter for this transition are not the ones associated with the first neighbor distance , but rather the much larger ones corresponding to the cage size. Note that at the transition cages are correlated over large distances Kirkpatrick and Thirumalai (1988); Franz and Parisi (2000); Berthier et al. (2005); Biroli et al. (2006). Also, an expansion in the number of collisions Elskens and Frisch (1988) seems difficult to reconcile with our results, because we expect multiple collisions within a cage.

Is the high-dimensional dynamics related to MCT? The answer is that the result is not the one obtained from the usual procedure for building up a MCT equation, which for example gives a different scaling of with dimension Schmid and Schilling (2010); Ikeda and Miyazaki (2010). Instead, the system we obtain is formally quite close to the slightly more general case of soft-spin mean-field dynamics, Eq. (9), because we have mapped the system into a one-dimensional dynamics in the presence of a colored noise and friction, that have to be determined self-consistently through Eqs. (4) and (5).

In the dense glassy regime we obtain predictions for the scaling of the cage radius, of the dynamical transition density , and of the parameter , that differ from the ones of usual MCT Ikeda and Miyazaki (2010); Schmid and Schilling (2010); Charbonneau et al. (2011). Our results are fully consistent with those obtained from the thermodynamic approach Parisi and Zamponi (2010); Kurchan et al. (2013); Charbonneau et al. (2014a), which proves the exactness of the RFOT scenario Kirkpatrick and Thirumalai (1988, 1989); Kirkpatrick et al. (1989); Lubchenko and Wolynes (2007); Wolynes and Lubchenko (2012); Kirkpatrick and Thirumalai (2014) for statics and dynamics in , as conjectured in Kirkpatrick and Wolynes (1987). Our results are also in agreement with numerical simulations of hard spheres in large spatial dimension Charbonneau et al. (2014b); Charbonneau et al. (2013).

Interestingly, we find that an ergodic liquid phase of hard spheres exists for densities . This implies that hard sphere packings exist (at least) up to , and they can be constructed easily through a sufficiently slow compression of the liquid Lubachevsky and Stillinger (1990); Montanari and Semerjian (2006). Note that the value of is larger than the best known lower bound for the existence of sphere packings,  Vance (2011), and that it took 20 years to improve the previous best lower bound  Ball (1992) by a small factor . Our calculation is simple enough that there is hope to transform it into a rigorous proof, along the lines of Ben Arous et al. (2001). This would result in an improved constructive lower bound for sphere packings.

Future extensions of this work include the investigation of the effect of dissipation Berthier et al. (2000); Ikeda and Berthier (2013); Fuchs and Cates (2002); Szamel and Flenner (2011), the study of out-of-equilibrium aging dynamics Cugliandolo and Kurchan (1993), and the study of non-perturbative processes in through an instantonic expansion. The thermodynamic partition function of quantum systems is formally very similar to Eq. (27) and could also be studied along these lines.

Acknowledgements –

We warmly thank G. Biroli, J. P. Bouchaud, P. Charbonneau, M. Fuchs, H. Jacquin, Y. Jin, T. R. Kirkpatrick, C. Rainone, D. Reichman, R. Schilling, G. Tarjus and P. Urbani for very useful discussions. We specially thank G. Szamel for very useful comments and for pointing out an inconsistency in an earlier version of this work. T. M. acknowledges funding from a fondation CFM grant.



I Formulation of the dynamics

i.1 Definition of the model

We consider an assembly of spheres of mass interacting through a repulsive finite-ranged radial pair potential . The interaction Hamiltonian is then with the positions of the particles. is the diameter of the spheres, we note and the rescaled force . In the following, we define the usual inverse temperature and take . In order to exploit efficiently simplifications given by the translational and rotational symmetries of the system, we constrain each particle to live on the surface of a dimensional hypersphere of radius . A particle is thus represented by a point with . The volume of this space is . and are respectively the surface and the solid angle. We also denote by


the volume of a -dimensional ball of radius (i.e. the volume bounded by is ). Adding this extra dimension, rotation and translation invariances in usual dimensional space are transposed to rotational invariance only on the sphere . We recover in the limit the original definition in a dimensional periodic cubic volume. However, note that in subsections I.2 and I.3, we consider the original model in dimensional Euclidean space, so as to introduce the spherical model only when needed.

i.2 The dynamical action

A possible choice of the dynamics of the particles in the -dimensional space is the following Langevin process:


being a centered white Gaussian noise with . In the following, we will consider the overdamped111Either Newtonian or Brownian dynamics can be treated by changing the kernel associated to in the following sections. case where . Given a set of initial conditions, the Martin-Siggia-Rose-De Dominicis-Janssen (MSRDDJ) Castellani and Cavagna (2005); Cugliandolo (2003) generating path integral reads, with It convention:


where the action is:


The sum in is done over all possible trajectories, with the measure given by when discretizing the trajectory in time steps. Time integrals are taken over an interval noted in section IV where the are fixed (initial conditions, which would correspond to in the latter discretization) and the are summed over (which would correspond to ). In section IV we will need averages of the type with , and it is not needed to consider times larger than in the action owing to causality, as they will give no contribution to the average by probability conservation (i.e. for the same reason as ).
The action is rotation invariant (for both position and response fields using the same global rotation) and translation invariant along positions222The invariances for response fields follow if one studies the system modulo ‘rigid’ rotations and translations of the entire system.. Because of the so-called kinetic term, it is not translation invariant along the response fields (which are already centered at the origin owing to the white noise), though the interaction term is.

i.3 Derivation of the generating functional using a virial expansion

i.3.1 Dynamic virial

We define in order to use the Mayer expansion333This is why is cast in a symmetric form so that links in the product are not directed. as in liquid theory Hansen and McDonald (1986). This ‘grand canonical’ form is handy as a generating functional, but note that we assume there is no exchange of particles with a reservoir ; and are related by a Legendre transform and virial expansions are more fruitful with the former. We have:


with a generalized fugacity and a Mayer function . We Legendre transform with respect to since one has


where the mean is generated by the functional . Next, the usual Mayer expansion can be carried out in this dynamical case, and inverting the Legendre transform, can be written as an ideal gas contribution and the sum of all connected 1-irreducible Mayer diagrams4441-irreducible means that they do not disconnect upon removal of a node. with nodes and  bonds Hansen and McDonald (1986):


In infinite dimension, the Mayer expansion reduces to its first term. However, this is strictly true if we assume that we are in a regime where the trajectories have the time to wander away only a finite fraction of the box. Because we expect (and confirm) that all interesting dynamics (namely, the relaxation with the formation of a plateau in correlations and the onset of a relaxation towards equilibrium - relaxation -) occur on such scales, where the fluctuation around the initial position is of amplitude (a scaling consistent with the statics Parisi and Zamponi (2010); Kurchan et al. (2012, 2013); Charbonneau et al. (2014a)). We will show that one even gets the diffusive behaviour which is already decided at the scale (see figure 1(a)).

(a) (b)
Figure 1: (a) The different dynamical regimes, described by equation (107): most of the dynamics is determined at the cage size scaling . We assume that the diffusive regime found already at this scale extends trivially at longer times. Concerning the slope at very short times, the mean-square displacement is if inertia is not neglected; for purely overdamped motions it is . Finally, the particles feel the ‘box’ for distances of order , hence a saturation at this scale. (b) The dynamic triangle diagram with three trajectories.

To justify this truncation, let us consider for example the second term of this expansion (triangle diagram):


For a finite-support potential , (Heaviside’s step function555In the main text it is noted but this notation will be used for Grassmann variables from subsection I.5 on.), hence


Essentially, is just a numerical value depending on the details of (with ) and the physical content of lies in its finite support, i.e.


For finite times, a trajectory stays in a bounded region of space, represented as a ball of diameter the typical size of the trajectory in figure 1(b). The Mayer functions require that each couple of trajectories get closer than at some time. To each set of three trajectories one can associate a corresponding static diagram with three overlapping balls (figure 1(b)). One can see this static diagram as the equivalence class of all dynamic diagrams with trajectories contained inside these balls. Actually there are lots of trajectories contained by these bounding balls that do not contribute because they do not get close enough. Then the sum over trajectories of the value of the integrand is at most of the same order of the static diagram value due to the normalization which accounts for the huge number of equivalent dynamic diagrams compared to the static one. The same program applies to all terms in the expansion. We conclude, from Wyler et al. (1987), that the first term dominates the series666The critical packing fraction later found in (147) is of order where it is known that the virial series does not converge anymore; nevertheless, Frisch & Percus Frisch and Percus (1999) have shown that in high the truncation used here is valid largely above this bound. in the limit. A more careful but maybe less intuitive argument is given in I.3.2. Note that we did note specialize the discussion to hard spheres but to any finite-support interaction potential. The result is then777In the thermodynamic limit (as if formally the chemical potential was zero in the static partition functions analogy). That is why for simplicity we used in equation (13) of the main text.


with from the Legendre transform and the normalization . We neglected purely additive constants irrelevant for the dynamics.

i.3.2 The Mayer expansion in infinite dimension

We discuss here the truncation of the Mayer expansion when in a more pedestrian way. As in I.3.1, we focus on the triangle diagram in (31) with the Mayer function given by (33).
Consider two typical trajectories and which we expect from the static case to dominate the dynamics at large , and let us focus on the positions and (any other couple would do). These typical trajectories scale as follows: the other positions of the trajectory (respectively ) fluctuate around (respectively ) over a neighborhood of size , see figures 1(a),(b). Besides, the typical distance between the two trajectories is the diameter (plus