Finite size effects in the microscopic critical properties of jammed configurations: a comprehensive study of the effects of different types of disorder
Patrick Charbonneau, Eric Corwin, Cameron Dennis, Rafael Díaz Hernández Rojas, Harukuni Ikeda, Giorgio Parisi, Federico Ricci-Tersenghi
FFinite size effects in the microscopic critical properties of jammed configurations: acomprehensive study of the effects of different types of disorder
Patrick Charbonneau,
1, 2
Eric Corwin, Cameron Dennis, Rafael D´ıaz Hern´andezRojas, Harukuni Ikeda, Giorgio Parisi,
4, 6 and Federico Ricci-Tersenghi
4, 6 Department of Chemistry, Duke University, Durham, North Carolina 27708, USA Department of Physics, Duke University, Durham, North Carolina 27708, USA Department of Physics and Material Science Institute,University of Oregon, Eugene, Oregon 97403, USA Dipartimento di Fisica, Sapienza Universit`a di Roma, Piazzale Aldo Moro 5, 00185 Rome, Italy Graduate School of Arts and Sciences, The University of Tokyo 153-8902, Japan INFN, Sezione di Roma1, and CNR-Nanotec, unit`a di Roma, Piazzale Aldo Moro 5, 00185, Rome, Italy
Jamming criticality defines a universality class that includes systems as diverse as glasses, col-loids, foams, amorphous solids, constraint satisfaction problems, neural networks, etc. A peculiarlyinteresting feature of this class is that small interparticle forces ( f ) and gaps ( h ) are distributedaccording to non-trivial power laws. A recently developed mean-field (MF) theory predicts the char-acteristic exponents of these distributions in the limit of very high spatial dimension, d → ∞ and,remarkably, their values seemingly agree with numerical estimates in physically relevant dimensions, d = 2 and 3. These exponents are further connected through a pair of inequalities derived fromstability conditions, and both theoretical predictions and previous numerical investigations suggestthat these inequalities are saturated. Systems at the jamming point are thus only marginally stable.Despite the key physical role played by these exponents, their systematic evaluation has remainedelusive. Here, we carefully test their value by analyzing the finite-size scaling of the distributionsof f and h for various particle-based models for jamming. Both dimension and the direction ofapproach to the jamming point are also considered. We show that, in all models, finite-size effectsare much more pronounced in the distribution of h than in that of f . We thus conclude that gapsare correlated over considerably longer scales than forces. Additionally, remarkable agreement withMF predictions is obtained in all but one model, near-crystalline packings. Our results thus help tobetter delineate the domain of the jamming universality class. We furthermore uncover a secondarylinear regime in the distribution tails of both f and h . This surprisingly robust feature is thoughtto follow from the (near) isostaticity of our configurations. I. INTRODUCTION
Jammed systems may lack dynamics, but their study isfar from motionless. A surge of physical interest over thelast couple of decades has indeed led marked advances [1–8]. This interest stems partly from jamming occurring insystems as varied as grains, foams, and emulsions, andpartly from jamming exhibiting universal features. Themix of ubiquity and unversality has motivated the searchfor a common framework to explain the pervasiveness ofjammed systems and their properties, starting with theseminal works of Liu, Nagel and co-workers [9, 10]. Ithas since become clear that although different systemsreach jamming by tuning different physical variables, sev-eral properties near and at the onset of jamming areshared by all of them. In other words, the same un-derlying physics should be responsible for the jammingphenomenology. Even though a fully comprehensive the-ory remains to be formulated, a major step forward hasbeen the discovery that this jamming point is a criticalpoint that gives rise to a phase transition, albeit an out-of-equilibrium one [11].It is therefore without loss of generality that attemptsto better understand jamming [2] have focused on sys-tems of frictionless spherical particles[5]. An outstandingexample of the theoretical analysis that can be achieved by such geometric simplification is the recently devel-oped mean-field (MF) theory [1, 6, 12–15] that describes–exactly, in the infinite-dimensional limit– the behaviorof glass-forming liquids from the point they fall out ofequilibrium up to jamming. Even though one mightexpect this theory only to be valid in high spatial di-mensions, near jamming it describes many of the criticalproperties observed in dimensions as low as d = 2 and d = 3[1, 6, 16, 17]. (A different criticality is observed inquasi-one-dimensional systems [18, 19].) Jamming criti-cality is peculiar because not only thermodynamic vari-ables, such as the pressure or bulk and shear moduli,but collective quantities, such as the mean squared dis-placement and the average contact number, scale criti-cally with the distance from the jamming point. Morespecifically, denoting the configuration density (or pack-ing fraction) φ and its value at jamming φ J , several quan-tities either jump discontinuously or scale as power laws, | φ − φ J | µ , as the jamming point is approached [3, 10, 20].Naturally, the exponent µ is specific for each quantity,but surprisingly in the T → a r X i v : . [ c ond - m a t . d i s - nn ] N ov portantly, once a jammed state is reached for a givenpotential, the resulting configuration is an equally validjammed state for any other potential [16] .However, this broad universality does not prevent µ from depending on whether the jamming point isapproached either from below ( i.e. from the under-compressed (UC) phase, φ → φ − J ) or from above (over-compressed (OC) phase, φ → φ + J ). A salient exampleis pressure, P , which scales as P ∼ | φ − φ J | ± [10, 21],i.e., µ ± = ± φ → φ ± J . In the UC case, pressure thusdiverges as density approaches φ J , as found in granu-lar materials or glass-formers made out of infinitely hardparticles [4]. Conversely, in the OC case, pressure van-ishes linearly as the packing fraction is brought down to φ J , as found in soft-harmonic particles [10]. Another im-portant example is the average contact number, z . Sim-ulations of harmonic soft spheres, for instance, show that z exhibits a discontinuity exactly as φ → φ − J , and thengrows as z ( φ ) − z ( φ J ) ∼ ( φ − φ J ) / for φ > φ J [10].This discontinuity can be related to the condition thatthe number of contacts in a configuration should exactlymatch its number of degrees of freedom, i.e. , the onsetof isostaticity [3, 7, 22]. Recent studies have further ver-ified the expected finite-size scaling of P , z , and the bulkand shear moduli for a wide variety of potentials in d = 2and 3 [23, 24]. A Widom-like scaling function has furtherbeen derived for these variables as well as for the con-figurational energy and shear stress [25]. Furthermore,various studies have identified correlation lengths asso-ciated to the characteristic length scales of vibrationalresponse to perturbations [26, 27], the fluctuations in thenumber of contacts [28, 29], and the fluctuations of par-ticle mobility [20], all of which diverge at the jammingpoint. These observations for thermodynamic variablesand bulk properties provide some of the strongest evi-dence in support of the critical nature of the jammingtransition.Remarkably, some of the microscopic structural prop-erties of jammed configurations, such as the distributionsof contact forces and interparticle gaps, are also expectedto exhibit non-trivial critical scalings. In particular, in ajammed configuration of N spherical particles with centerpositions r i and diameters σ i , one can define a dimension-less gap between any pair of particles, h ij = | r i − r j | σ ij − σ ij = ( σ i + σ j ) /
2. Because jammed packings are dis-ordered, gap values are randomly distributed, but theo-retical predictions [15] state that the distribution of smallgaps should scale as g ( h ) ∼ h − γ , with γ = 0 . . . . (1)Similarly, the distribution of small contact forces is pre-dicted to scale algebraically, p ( f ) ∼ f θ , but initial re-ports found a strong dependence of θ on dimensionalityand jamming protocol, in apparent contradiction withthe theoretical expectation [30]. This paradox was re-solved by recognizing that two different types of forcescontribute in this regime [16, 30]. Opening the contact between a pair of particles can indeed give rise to two dis-tinct responses: (i) a localized rearrangement of neigh-boring particles; or (ii) a displacement field that extendsover the whole configuration, without decaying with dis-tance. The former is associated with a buckling motion,and hence remains localized; the latter is associated witha correlation length of the same order as the system size,and hence is a clear example of the criticality of jammedpackings. Considering these two types of forces sepa-rately yields two power laws with different exponents, p ( f (cid:96) ) ∼ f θ (cid:96) (cid:96) , with θ (cid:96) (cid:39) . , (2a) p ( f e ) ∼ f θ e e , with θ e = 0 . . . . ; (2b)for localized and extended excitations, respectively. Theability of MF theory [1, 6, 15] to predict the non-trivialvalues of γ and θ e is considered a major analytical suc-cess. MF theory, however, does not directly predict θ (cid:96) ,because bucklers are an intrinsically low-dimensional fea-ture [16], and are therefore absent from the d → ∞ de-scription. The critical exponents of gaps and contactforces are also of utmost importance because they areassociated with the mechanical stability of jammed pack-ings. By considering the displacement field that followsopening one of two types of contacts as well as the ensu-ing closure of gaps to form stabilizing contacts, a pair ofinequalities between γ , θ (cid:96) , and θ e can be derived [17, 31], γ ≥ − θ (cid:96) , (3a) γ ≥
12 + θ e . (3b)MF theory values as well as numerical simulations showthat both inequalities are in fact saturated, implying thatjammed packings are marginally stable [31, 32]. Thisresult is consistent with the MF description, which al-ways locates the jamming point within a critical Gard-ner phase that emerges deep in the glass phase andis characterized by the emergence of marginally stablestates [1, 6, 13, 15, 33, 34].The picture that coalesces from putting together theexact MF description with the critical scalings for ther-modynamic and other variables, and from considering therobustness of numerical experiments for several dimen-sions and for different protocols [4, 6, 11, 16], suggeststhat the jamming transition of spherical particles prop-erly defines a universality class. We now know that thisclass should encompass a broad range of problems andmodels beyond spherical particles, including the percep-tron [35, 36], neural networks [37–39], statistical infer-ence [40], and the SAT-UNSAT transition in continuousconstraint satisfaction problems [41, 42]. Recent workshave shown that universality persists even when the in-teractions are non-analytic, for instance, due to discon-tinuous forces [36, 43].Yet, a careful analysis of the values of θ (cid:96) , θ e , and γ inferred from numerical simulations has not system-atically been carried out. Conducting such an analy-sis is especially important considering that packings ofslightly polydisperse crystals are reported to exhibit a mi-crostructure characterized by exponents that differ con-siderably from those of Eqs. (1) and (2) [44, 45]. Ad-ditionally, recent works have shown that many of thesalient features of spherical packings depend sensitivelyon particle shape. For instance, introducing even an in-finitesimal amount of asphericity changes the universalityclass [46, 47], in which the isostatic condition no longerholds. An assessment of the extent of the jamming uni-versality class and an accurate test of its many theoreticalpredictions are therefore in order [34].In this work we systematically analyze the finite-sizescaling of the distributions of interparticle gaps and con-tact forces. These distributions are one of the funda-mental consequences of the presumed non-trivial criticalbehavior of jammed packings, hence their testing is a keystep toward rigorously validating a whole set of criticalproperties. Although a similar analysis has been carriedout for the perceptron [48] and for the gaps distributionof a two-dimensional binary mixture [46], no systematicresult exists for jammed packings of spherical particlesnor for amorphous packings with other sources of disor-der. Here, in addition to analyzing the most commoncases of jammed configurations, i.e. d polydisperse and3 d monodisperse packings, we consider two additionalsets of jammed packings: (i) polydisperse spheres in acrystalline FCC structure; and (ii) Mari-Kurchan (MK)hard spheres with random shifts distributed uniformlyover space [49]. By examining the impact of differentsources of disorder, we attempt to define precisely whichare the most robust features of jamming criticality, andthus better demarcate its physical universality. The restof this paper is organized as follows. In Sec. II we de-scribe the models used and the algorithms employed toproduce jammed configurations and extract the relevantstructural information, i.e. the interparticle gaps, h , andcontact forces associated with extended, f e , and local-ized, f (cid:96) , displacement fields. We also explain how finite-size effects in the distributions of these structural vari-ables are considered. In Sec. III we present a detailedanalysis of the finite size effects in jammed configurationsof monodisperse spherical particles in 3 d , where we revealthe striking contrast of such effects on the distributions f e and h . Then, in Sec. IV we present a similar analysis forthe other types of systems considered, finding importantdifferences with the results for d = 3 spherical systems.We nevertheless argue that most of these differences canbe explained from the other scaling corrections describedin Sec. II D. Because theory and previous numerical stud-ies suggest that f e and h are critically correlated acrossthe whole system, we first consider these two quantities.The distribution of localized forces, f (cid:96) , associated withbuckling effects is expected to be independent of systemsize, hence its analysis is postponed to Sec. V. A discus-sion and brief conclusion are given in Sec. VI. II. NUMERICAL METHODS, MODELSSYSTEMS, AND FINITE-SIZE SCALING
In this section, we describe the numerical techniquesused to produce jammed sphere packings, coming fromeither the OC or the UC phase. Studying independentlythese two regimes is useful because –as for other criticalpoints– there is no reason a priori to assume that thescalings from above and below φ J are the same. Becauseeach of these two phases is identified with different mate-rials, namely granular matter (from the UC regime) andglasses, foams, and colloids (from the OC phase), thisverification is an important test of materials universal-ity. We also describe the other models considered, whichare chosen to better appraise the extent of the jamminguniversality class. The methodology employed to ana-lyze the system size dependence on the distributions ofthe microstructural variables, Eqs. (1) and (2), is alsodetailed. A. Jammed states from the OC phase
We first consider three-dimensional configurations of N spheres of equal diameter σ in a cubic box under periodicboundary conditions. In a certain sense, this choice is theminimal model with which to produce jammed packings.Lower-dimensionality systems inevitably crystallize un-less polydisperse mixtures are used, but ordering can beavoided for monodisperse spheres in d ≥
3. Sphere po-sitions then serve as the only source of disorder. Giventhe set of vectors of positions { r i } Ni =1 , the jamming pointstarting from the OC phase is obtained for the harmoniccontact potential, U (cid:0) { r i } Ni =1 (cid:1) = (cid:15) (cid:88) i,j ( σ − | r i − r j | ) Θ( σ − | r i − r j | ) , (4)where (cid:15) is a constant that defines the energy scale, and Θis the Heaviside step function. Hence, a pair of particlesonly interacts if there is an overlap between them. Start-ing in the OC phase with φ > φ J and a Poisson randomdistribution of spheres, a series of energy minimizationsteps and packing fraction reduction steps are performeduntil the system has just a single state of self stress, whichis where jamming criticality occurs [17, 31, 50, 51]. Sucha state is characterized for having one contact above iso-staticity, i.e. when the total number of constraints ina system, N c , matches its number of degrees of free-dom, N dof . A single state of self stress is required forcritical jamming in order to achieve a finite bulk mod-ulus [23, 52]. Put differently, the system density is anadditional variable that needs to be fixed, and thus re-quires one additional contact above isostaticity [53]. Ata given density the FIRE algorithm, a damped dynamicsmethod, is used to achieve force balance in the configura-tion [54]. The energy of the configuration is then calcu-lated and the known scaling relation, U ∝ ( φ − φ J ) [16],is used to determine by how much the sphere radii shouldbe uniformly decreased to reduce the system energy bya fixed fraction. After several iterations of this proce-dure, the packing has precisely N c = N s d − d + 1 con-tacts where N s is the number of stable particles and thus N dof = d ( N s −
1) corresponds to the number of degreesof freedom in a system under periodic boundary condi-tions. A small fraction of particles, termed rattlers, re-main unconstrained at jamming and do not contributeto the overall rigidity of the packing[16, 24, 53]. In a d -dimensional system, these rattlers can be identified asparticles with fewer than d + 1 contacts. Although thefraction of rattlers changes from one configuration to an-other, it always lies within a small range of ∼ − N , is thus reported. After removing rattlers, the dy-namical matrix [7] is used to ensure that the packing isjammed. This algorithm is implemented in the pyCud-aPacking software using general purpose graphical pro-cessing units and quad-precision calculations [55–57]. B. Jammed states from the UC phase
For configurations initially in the UC regime, a hard-sphere potential is used and a combination of moleculardynamics (MD) and linear optimization algorithms areemployed to approach φ J from below. More precisely,we start from a low-density configuration of particleswith random positions and use event-driven MD witha Lubachevsky–Stillinger (MD-LS) growth protocol [58]to increase the (reduced) pressure up to P = 500. Thisfirst step is performed with a fast compression rate inorder to avoid any partial crystallization and is then fol-lowed by a second, much slower, growth protocol until P (cid:38) . At this point the configuration is used as inputfor the sequential Linear Programming (LP) algorithmused in Refs. 59 and 60 to produce jammed packings. Ateach step, the LP algorithm finds the optimal rearrange-ment of particles positions that maximizes their radius,considering a linearized version of the non-overlappingconstraint between any pair of particles. Upon conver-gence, this algorithm produces a jammed configuration,because neither particle displacements nor size increasesare possible. This approach also allows to easily buildthe full network of contacts at jamming, because genuinecontact forces can be identified, up to a proportionalityfactor, from the active dual variables associated to thenon-overlapping constraints.As with the OC phase, rattlers are removed and onlysystems at a single state of self stress are considered. C. Other models of jammed packings
We also investigate the jamming point of three othermodels.
Polydisperse disks:
Previous studies strongly suggestthat the upper critical dimension of the exact MF theoryis d = 2 [16, 17, 24]. However, as mentioned above, parti-cles of different sizes must then be utilized to inhibit crys-tallization. An additional source of disorder is thus intro-duced by extracting particle radii from a log-normal dis-tribution with parameters µ = 0 and σ = (cid:112) ln (0 . + 1)to achieve a polydispersity–defined as the ratio of stan-dard deviation to mean–of 20%. (Note that the radiidistribution parameters should not be confused with theparticle diameter used in monodisperse systems.) Thesesoft harmonic spheres are initially in the OC regime, andthus the FIRE-based algorithm is used to bring configu-rations to their jamming point. Crystalline polydisperse spheres:
Removing random-ness from particle positions while keeping size polydis-persity as the main source of disorder is achieved bygenerating jammed packings on the sites of a regularface-centered cubic (FCC) lattice. Radii are drawn froma log-normal distribution with a polydispersity of 3%.These nearly crystalline packings are brought to criticaljamming using the FIRE-based protocols for soft spheresinitially in the OC phase. Although this type of sys-tem displays many of the features associated with tra-ditional glasses [44], its distributions of forces and gapsoften markedly differ from those predicted by MF the-ory [44, 45]. By using a system with a different crystallinesymmetry we aim to quantify such discrepancy.
Monodisperse Mari-Kurchan (MK) spheres:
The MKmodel is a MF reference given that, by construction, theproperties of MK configurations are roughly independentof dimension. Specifically, we consider d = 3 systems ofmonodisperse spheres that interact according to a ran-domly shifted distance, d ( r i , r j ) = | r i − r j + A ij | , where A ij is a quenched random vector drawn uniformly fromthe total system volume. Introducing random shifts, A ij , suppresses almost completely correlations due toshort loops on the interaction graph. Even if d ( r i , r j ) = d ( r j , r k ) = σ it is very unlikely that d ( r i , r k ) (cid:39) σ . Inother words, while for particles interacting via the usualEuclidean distance neighbours of a given particle arelikely also neighbours, in the MK model, almost cer-tainly, they are not. Because this property is also the casefor systems using the Euclidean distance in the d → ∞ limit, it is expected that the microscopic structural prop-erties of MK jammed configurations should follow the MFtheory predictions closely. Besides, it has already beenverified that the MK model exhibits several features ofmore usual glass formers [61], that a Gardner transitionalso occurs deep in the glass phase [62], and that con-tact number fluctuations are critically correlated at jam-ming [29]. Consequently, any deviation from MF predic-tions observed for this system can safely be attributed tofinite-size corrections, which makes the MK model a par-ticularly useful reference to explain the contrasting scal-ing effects in the distributions of gaps and contact forces(Sec. VI). For this model, we consider hard sphere con-figurations initially in the UC phase, and use the MD-LSand LP algorithms to reach their corresponding jammingpoint. D. Expected finite-size scalings
For each system several independent jammed configu-rations are obtained, and interparticle gaps, h , and con-tact forces, f are computed. For the latter, contributionsassociated with localized buckling displacements, f (cid:96) , areseparated from those that produce extended excitations, f e , using the fact that with very high probability buck-lers are associated to particles with d + 1 contacts andlocalized excitations [16]. Hence, the set { f e } is taken asthe set of forces applied on particles with more than d +1contacts.To ensure that we sampled all the systems of a giventype with the same accuracy, M N independent configu-rations are produced for a fixed value of N , such thatdata of N × M N (cid:38) particles is obtained. (Specificvalues for each system are given below.) Forces and gapscan then be studied across many orders of magnitude,and finite-size corrections can be systematically identi-fied. Because testing for power-law distributions usinglogarithmic binning of the probability density function(pdf) leads to poor comparisons (due to the loss of res-olution when grouping data in a single bin to producea smooth trend [63]), the cumulative distribution func-tion (cdf) is considered instead. Note that if a randomvariable x is distributed according to a pdf of the form p ( x ) ∼ x α for α > −
1, then its cdf follows c ( x ) ∼ x α .When fitting a distribution to empirical data it shouldbe considered that even if x ideally follows such a dis-tribution all the way down to x →
0, finite sampling in-evitably leads to deviations. Here, the situation is furthercomplicated by our consideration of marginals of corre-lated variables. Gaps and forces distributions of finite N configurations are indeed prone to exhibit deviationsfrom their expected form due to both finite samplingand system-wide correlations. Fortunately, introducinga scaling function, as is usually done in the study of crit-ical phenomena [64, 65], can account for both effects, andhence the dependence of the cdf on system size can becarefully teased out.To derive the size scaling of the distributions of x , wefirst note that in a sample of size N (cid:29)
1, we can estimatethe order of the smallest value observed in the data, x min : (cid:90) x min p ( x ) d x ∼ x α min ∼ N . (5)It then follows that x min ∼ N − / (1+ α ) . Note that strictlyspeaking in this last equation N should be replaced by N c when analyzing, for instance, the distribution of con-tact forces. However, given that N c ∼ dN and that weare mostly concerned with the scaling exponent, we cansafely neglect the associated proportionality constants.The behavior of the gap distribution is expected to be similar, in that the amount of particles almost in contactshould be self-averaging. Next, we follow the traditionalpath for analyzing size scaling and write the pdf as p ( x ) ∼ N β ˜ p (cid:16) xN α (cid:17) (6)where ˜ p is the scaling function of the pdf such that˜ p ( x ) ∼ x α for x (cid:38)
1. The exponent β can be easily deter-mined by requiring that p ( x ) exhibit no N dependence.We thus get that β = − α α , whence the expressions usedfor the scalings studied in Ref. 48 are recovered. For thecumulative distributions, repeating the above analysis for c ( x ) ∼ N − δ ˜ c (cid:16) xN α (cid:17) gives ˜ c ( x ) ∼ x α , and it imme-diately follows that δ = 1, whence the relevant scalingrelation is c ( x ) ∼ N − ˜ c (cid:16) xN α (cid:17) . (7)Using the correct α should remove any dependence on N .Data for different system sizes should then be rescaledsuch that they follow a common curve, ˜ c . Finding agood collapse of the curves for different N thus indicatesthat deviations from the expected power laws fall out-side the thermodynamic limit, but are not caused by thevariables following a different power-law scaling. Addi-tionally, showing that the system size influences the cdfof a given variable strongly evinces that such a variableis correlated across the whole system. Hence, an upperbound to the correlation length can then be estimated.(Note that for microscopic variables of jammed config-urations, the situation is conceptually different from thatof standard critical phenomena, because the systems arealready at the critical point. We here do not investigatehow the distributions of contact forces and gaps convergeto their expected distributions as we move away from φ J ,but instead analyze how the system size affects the rangeover which power-law scalings are followed. Equation (7)can nevertheless be used to estimate the scaling func-tions of the cdf of gaps and forces obtained by integratingEqs. (1) and (2), respectively.)At the upper critical dimension d = 2, we expect alogarithmic correction to the size scaling law [24, 66–68], p ( x ) ∼ x α ( − log x ) ξ , for x (cid:28) . (8)We can then estimate x min as (cid:90) x min p ( x ) dx ∼ x α +1min ( − log x min ) ξ ∼ N , (9)leading to x min ∼ N − α ( − log x min ) − ξ α ∼ N − α (log N ) − ξ α . (10)Repeating the same argument as above, we get c ( x ) ∼ N − (log N ) − ξ ˜ c (cid:16) xN α (log N ) ξ α (cid:17) , (11) − − − − − − f e − − − c d f ( f e ) (a) N − − − − − − f e (b) N − − − − N / (1+ θ e ) f e − − − N c d f ( f e ) (c) y ( x ) ∼ x θ e y ( x ) ∝ x Figure 1. Cumulative distributions of extended contact forces associated with extensive excitations of monodisperse configura-tions of frictionless spheres for different system sizes N , as their jamming point is reached (a) from below (UC) and (b) fromabove (OC). To better distinguish between the two different regimes, results belonging to the UC (OC) phase are identifiedby circular markers (solid lines). (c) Rescaling (a) and (b) according to Eq. (7) clearly collapses the data. The red dashedline corresponds to the power-law scaling of Eq. (3b), and shows an excellent agreement between the MF predictions and ournumerical results. The coincidence of results from the UC phase and OC phase for various N confirms that θ e is the samewhen jamming is reached from either direction. In the left tail of the distributions of panel (c) we also include a comparisonwith the linear scaling expected for very small values, following Eq. (13). where the pre-factor is chosen such that c ( x ) does notdepend on N for x (cid:29) x min . For the cases considered inthis work, no theoretical prediction exists for the valueof ξ , and hence it here serves as a fitting parameter.We consider yet another correction to Eq. (6) thatcan also be derived from MF theory. Given thatjammed configurations have one extra contact than N dof (see Sec. II A), the power laws of the microstructuralcritical variables should be cut off at very small val-ues [41, 46, 69]. MF theory predicts that interparticlegaps are distributed as h − γ only for values larger than acut-off h (cid:63) ∼ δz − γ , where δz is the excess of contacts ina system. In our case, δz ∼ /N , so instead of Eq. (1)the pdf describing the distribution of h reads, g ( h ) ∼ N γ − γ g (cid:16) hN − γ (cid:17) , hN − γ (cid:28) h − γ , hN − γ (cid:38) g ( x ) ∼ x (cid:28) p ( f ) ∼ N − θe θe p (cid:16) f N θe (cid:17) , f N θe (cid:28) f θ e , f N θe (cid:38) , (13)where p ( x ) ∼ III. FINITE-SIZE EFFECTS IN d = 3 SYSTEMS
We first consider systems of monodisperse particles in d = 3 by generating, for each N , M N independent pack-ings, such that N × M N (cid:39) . × (5 . × ) particlesare considered when the jamming point is approachedfrom the UC (OC) phase. Figure 1 shows the distri-butions of f e obtained coming from below (UC, panel(a)) and from above (OC, panel (b)). Comparing theresults with the theoretical prediction for the power-law − − − − − h − − − − c d f ( h ) (a) N − − − − − h (b) N − − − − − − N / (1 − γ ) h − − − N c d f ( h ) (c) y ( x ) ∝ x − γ y ( x ) ∝ x Figure 2. Cumulative distributions of interparticle gaps for the same configurations as in Fig. 1, as their jamming point isreached (a) from below (UC) and (b) from above (OC). (c) Rescaling (a) and (b) according to Eq. (7) shows that finite-size corrections can be accounted for in all cases. For comparison, the power-law scaling derived from MF theory, Eq. (1),is also shown (red dashed line). Once again, the fact that datasets from both phases, i.e.
UC (markers) and OC (lines),neatly superimpose confirms that the exponents at the jamming point are the same, independently of how φ J is approached.Additionally, the secondary scaling regime g ( h ) ∼ scaling reveals an outstanding agreement over at leastthree decades. More importantly, no visible signature offinite-size corrections can be detected over the range of N considered. To verify more stringently the absence offinite-size effects, we attempted to collapse the differentcurves by rescaling the extended forces and their cdf fol-lowing Eq. (7), obtaining the curves reported in panel(c). This last figure evinces that the same critical dis-tribution of forces is found independently of whether thejamming point is generated from the UC or OC regimes.Yet, it is clear that our packings exhibit an excess of verysmall forces (an effect more noticeable when jamming isreached from below; see Fig. 1a), echoing earlier obser-vations [16, 17, 44, 48]. Note that the scaling of Eq. (7)does not remove these deviations from the predicted dis-tribution. Note also that these deviations roughly occurfor the same scaled force, N θe f e (cid:46)
1. It is thereforelikely that forces are subject to size effects caused by theonset of a second power law, p ( f ) ∼ h are strongly depen-dent on system size. In contrast to p ( f e ), the scalingcorrection given in Eq. (7) using the MF value of γ pre-cisely corrects for such effects over almost seven ordersof magnitude (Fig. 2c). The growing deficit of very smallgaps as the system size decreases is another manifestationof the cut-off of the main power law of g ( h ). It leads toa secondary linear regime, as given in Eq. (12), that is inagreement with the numerical results (Fig. 2c). This in-dicates that distances between nearby spheres are signif-icantly modified in finite-size configurations and, conse-quently, so is the distribution of gaps. This phenomenonis physically interesting. Heuristically, the finite N in-fluence on g ( h ) can be understood by relying on themarginal stability of jammed packings. In the thermody-namic limit, a system has always enough space to relaxany perturbation caused by a contact opening, and henceis always able to re-accommodate particle positions–evenif this requires bringing many of them infinitesimally − − − − − f e − − − − − − − c d f ( f e ) (a) N − − − − − − f e − − − − − − − c d f ( f e ) ∝ f θ e e ∝ f e (b) N − − − − − − f e − − − − − − c d f ( f e ) (c) N − − − N / θ e f e − − − − N c d f ( f e ) (d) Figure 3. Cumulative distributions of f e for jammed configurations of (a) d = 2 polydisperse disks packings, (b) polydispersespheres with a FCC crystalline structure, and (c) packings using the d = 3 MK model. Panel (d) depicts the same data fromthe MK model, rescaled according to Eq. (7); see text for details. Data in the upper (resp. lower) panels were produced asjamming was approached from above (resp. below). The expected power law, Eq. (2b) is shown (red dashed lines), as is thesecondary linear regime, see Eq. (13) (cyan dotted lines). close to each other–in order to guarantee stability. In afinite system, by contrast, no such unconstrained relax-ation can take place. Rearranging an extensive fractionof particles necessarily influences the pair of spheres in-volved in the contact just opened. There is therefore acertain scale, below which the occurrence of small gapsis disfavoured. If the system were further relaxed, thenat least one extra contact would form.At this point, we wish to stress that our results demon-strate the existence of two different types of finite-sizecorrections to the distributions of extended forces andgaps. The first is a consequence of large scale correla-tions and can thus be readily taken into account by thescaling of the cdf given in Eq. (7). Although this correc-tion is practically absent in the forces distribution, for g ( h ) it is the main source of deviation from the theoreti-cal prediction. The second is a consequence of the criticalscalings of Eqs. (1) and (2b) being cut off at very smallvalues. This effect, which is very likely related to the ex-cess contact with respect to N dof (see Sec. II D), affectsboth microstructural variables. We get back to this pointin Sec. IV, after having considered its signature in othermodels.Before concluding this section, it is worth emphasiz-ing that our numerical results are in excellent agreementwith the MF, d → ∞ predictions for the power-law scal- ing of the distributions of both the extended forces andthe interparticle gaps. These results confirm that thejamming criticality of these microstructural variables isrobust with respect to changes in the systems dimen-sionality, all the way down to d = 3, in agreement withearlier albeit less accurate studies [11, 17, 30]. Becauseresults from both OC and UC phases superimpose ontoeach other, we further conclude that the critical behavioris controlled by the same exponents on both sides of thejamming point. IV. FINITE-SIZE EFFECTS IN OTHERDISORDERED SYSTEMS
We next consider the finite-size scaling of the forceand gap distributions at jamming for the three othermodels mentioned above: (i) polydisperse disks, (ii) crys-talline polydisperse spheres, and (iii) monodisperse MKspheres. From Sec. III, we understand that the directionof approach to the jamming point does not influence onthe criticality of microstructural variables, so only onesuch direction is considered for each mode. The first twoapproach the jamming point from the OC phase with N × M N (cid:39) × particles, and the third from the UCphase with N × M N (cid:39) . − − − − − h − − − − c d f ( h ) (a) N − − − N / (1 − γ ) h − − − N c d f ( h ) (b) − − − − N − γ (log N ) ξ − γ h − − − − N ( l og N ) ξ c d f ( h ) (c) y ( x ) ∝ x − γ y ( x ) ∝ x Figure 4. (a) Cumulative distributions of h of jammed con-figurations of d = 2 polydisperse disks and different size N .(b) Scaling of the different curves following Eq. (7) using theMF value of γ . (c) Same scaling as in (b) but including a log-arithmic correction as in Eq. (11). Choosing ξ = − . Despite the marked differences between the three mod-els, their distributions of f e all follow the MF predictionsvery closely (Fig. 3). (a) The d = 2 packings show a verygood agreement with the cdf derived from Eq. (2b) overmost of the accessible range. (b) Results for the FCCsymmetry also follow the expected scaling, but becauseits onset takes place at smaller forces, the range of consis-tency with the MF power-law scaling is correspondinglyreduced. (c) Jammed configurations produced using theMK model exhibit a noticeable–albeit small–dependenceon N , but this dependence can be removed by rescalingthe cdfs according to Eq. (7) using the MF value of θ e (see Fig. 3d). Interestingly, all three systems display anexcess of very small contact forces for f e (cid:46) − , sim-ilarly to what was found for d = 3 configurations (seeSec. III). Our results suggest that this effect is due toa crossover to a second regime, in which forces are dis-tributed uniformly, as given by Eq. (13). A comparisonwith the corresponding linear behavior in each panel ofFigure 3 presents a reasonably good agreement, in sup-port of this hypothesis. A more careful analysis wouldnonetheless be needed to single out the true form of theleft tails of p ( f e ).We next consider the finite-size effects on the distri-bution of gaps of these three systems. From the spacingbetween different curves in d = 2 packings, it is clearthat such effects are pronounced (Fig. 4a). Rescalingthese distributions following Eq. (7) with MF value for γ yields a collapse (Fig. 4b) that is not as good as fortheir d = 3 counterparts. Section II D anticipated thisdiscrepancy on the basis that d = 2 is the upper criti-cal dimension for jamming [24], and hence a logarithmiccorrection should be included, as in Eq. (11). As shownin Fig. 4c, with such correction the data can be robustlycollapsed using the MF value of γ .By contrast, gaps distributions in the FCC jammedconfigurations are best described by a completely differ-ent exponent. Figure 5 clearly shows that (a) finite-sizecorrections are important, but that (b) a poor collapse isobtained when curves are rescaled using Eq. (7) with theMF value of γ . Using (c) a different γ F CC (cid:39) .
33, how-ever, satisfactorily captures the N dependence. This con-firms previous reports that γ is changed in presence of anunderlying crystalline structure [44, 45]. Reference [45]even found that γ depends on the system polydispersity,through the variance of the particle sizes. It is importantto stress that finding a smaller γ is not merely a mat-ter of scrupulous curve fitting. It also positively violatesthe marginal stability relations, Eqs. (3). We commentfurther on this point in Sec. VI.Figure 6 presents the gap distributions for the MKmodel. Here again, finite-size corrections to g ( h ) are sig-nificant, but now taking the MF value of γ in Eq. (7)yields a very good collapse, as expected from the MF na-ture of the model. It is important to note that althoughindividual distributions of h suggest that a smaller expo-nent would better fit the curves in Fig. 6a, doing so wors-ens significantly the quality of the scaling collapse. This0 − − − − − h − − − − c d f ( h ) (a) N − − − N / (1 − γ ) h − − − N c d f ( h ) (b) y ( x ) ∝ x − γ − − − N / (1 − . h − − − N c d f ( h ) (c) y ( x ) ∝ xy ( x ) ∝ x − . Figure 5. (a) Cumulative distributions of h for jammed con-figurations of polydisperse spheres with an FCC structure anddifferent N . Scaling the different curves according to Eq. (7)using (b) the MF value of γ and and (c) γ FCC = 0 .
33. For aclearer comparison, the trend for the expected power-law ex-ponent (red dashed line) and for γ FCC (pink dashed-dottedcurve) are shown. For FCC configurations, unlike for d = 2systems, the collapse obtained with the MF value of γ is poorover the whole interval considered of the scaled variables (seeFig. 4b). Note that when γ FCC is used, a linear scaling atvery small arguments is recovered (cyan dotted line). − − − − − h − − − c d f ( h ) (a) N − − N / (1 − γ ) h − − N c d f ( h ) (b) y ( x ) ∝ x − γ y ( x ) ∝ x Figure 6. (a) Cumulative distributions of h for jammed con-figurations of d = 3 MK systems of different size N . (b) Samedata but collapsed using the scaling in Eq. (7). Such scalingindicates that γ MK = γ , in agreement with MF theory. Theexpected power-law scaling for the cdf, Eq. (1) (red dashedline), shows that for this type of systems finite-size correctionsare particularly important (see main text for discussion). Theexpected linear scaling for the left tail of the distributions isalso shown (cyan dotted line). situation is typical of many critical scalings in finite- N systems [64, 65]. The most reliable way to determine truecritical exponents remains the finite-size scaling analy-sis. The size dependence of the (apparent) exponent, i.e. γ MK = γ MK ( N ), and of cdf ( f e ) are nonetheless substan-tial, especially compared to that of the other models atsimilar sizes N . This discrepancy likely results from theMK system being fully connected. In contrast with theirsparse counterparts, fully connected models indeed re-quire much larger system sizes for thermodynamic power-law scalings to be visible [70–72]. This feature can be1 − − − f ‘ − − − − − − c d f ( f ‘ ) (a) ∝ f θ ‘ ‘ N − − − f ‘ (b) N − − − f ‘ (c) N − − − f ‘ ∝ f ‘ (d) N Figure 7. Cumulative distributions of f (cid:96) for jammed packings of (a) d = 3 monodisperse spheres, (b) d = 2 polydisperse disks,(c) polydisperse spheres with FCC structure, (d) d = 3 MK model. Solid lines (circular markers) denote data obtained fromconfigurations from the OC (UC) phase. For reference, the expected power law, cdf ( f ) ∼ f θ (cid:96) , with θ (cid:96) = 0 .
17, is shown (reddashed lines), and in panel (d) the power-law fit found by inspection for the MK model, cdf ( f ) ∼ f , i.e., θ (cid:96),MK = 0, is alsoshown (pink dotted line). See text for more details. physically understood by recalling that the introductionof random shifts results in neighbors of a given particle(very likely) not being neighbors themselves. A particlecan thus have many more contacts than normally allowedin Euclidean space. For instance, it is not uncommon( ∼ d = 3 kissing number) or more. In general,particles are thus surrounded by many more particles–both actual and near contacts–than usual hard spheres.Additionally, jamming densities in this model are muchhigher than can be achieved with hard spheres. Usingour MD-LS+LP algorithm, as well as planting [62] tospeed up the growing protocol, results in jamming pack-ing fractions φ J,MK (cid:38) . φ J, d (cid:39) . φ ∼ σ /d , our MK configurations are made out ofparticles nearly twice as big as those of standard hardspheres. The combination of these two effects is thatparticles in MK packings are surrounded by a cluster ofmany relatively large neighbors. The effective size of thesystem being drastically reduced, finite-size correctionsare correspondingly more pronounced.Looking at the whole set of gap distributions, an in-teresting feature is the robust emergence of a regime ofuniform distribution at very small gaps, in a way entirelyanalogous to the distributions of extended contact forces.We argued in Sec. II D that this truncation of the leadingpower-law scaling in the distributions likely follows fromthe combined effect of the additional state of self stressand the system sizes being finite. All the models consis-tently exhibit this behavior and show very good agree-ment with the associated linear scaling (see Figs. 2c, 4c,5c and 6b). The invariance of this secondary power-lawscaling with dimensionality, inherent order or other sys-tem properties is reassuring, albeit somewhat surprising, given that the leading power-law scaling is more stronglyaffected by these same effects. The universality of thissecondary scaling has been previously predicted [41] forall models that can be mapped to jamming of spheri-cal particles, and it has been shown to occur even fornon-spherical particles [46], provided that their jammedstates remain sufficiently close to isostaticity. Such ro-bustness can be understood in part by considering thatisostaticity is a global property of the system related toa matching between constraints and degrees of freedom,and not to the specific distributions of its microstruc-tural variables. Because we have restricted our analysisto packings with exactly N c = N dof + 1, the ubiquity ofthe linear left tails in our distributions supports the hy-pothesis that the form of g ( x ) (see Eq. (12)) and p ( x )(see Eq. (13)) is determined by the single state of selfstress alone, and not by the inherent structure. V. CUMULATIVE DISTRIBUTIONS OF f (cid:96) The last microstructural variable we consider is theset of localized forces. Figure 7 presents the probabil-ity distributions for all our results. As expected, thisquantity exhibits no clear finite N signature for any ofthe models, even though some dispersion around the ex-pected behavior is observed in the left tails of d = 3monodisperse and MK configurations, (panels (a) and(d), respectively). This behavior is expected because theset { f (cid:96) } corresponds to contact forces acting on bucklers,for which opening a weak contact mostly results in local-ized displacement field [16, 17]. Because opening any ofthe contacts associated with a buckler only has a non-negligible effect over a few particle layers away from its2origin, it is reasonable to assume that their propertiesshould be insensitive to N , or to any border or periodiceffects.An intriguing finding is that only the cdf of d = 3monodisperse and d = 2 polydisperse particles follow theknown value of θ (cid:96) (cid:39) .
17 (see Fig. 7a-d). By contrast,FCC structures give rise to no obvious power-law scaling.The FCC arrangement induces strong spatial correlationsthat seem to suppress the appearance of localized forces,as seen from the smaller slope of the cdf. Observing adistribution with an exponent different from θ (cid:96) , or actu-ally failing to scale as a power law, is in striking contrastwith many other models, and even other crystalline struc-tures [44]. It nonetheless echoes very recent reports of adependence of θ (cid:96) on geometry for other near-crystals [45].These considerations highlight the need for further as-sessment of which aspects of jamming criticality are in-deed universal, which are more generically conserved [73],and which disappear in the presence of long-range spatialconstraints.Although a power-law scaling is also obtained for MKconfigurations, the best fit to the data is achieved with aunit slope, i.e., θ (cid:96),MK = 0 (see Fig. 7d). Localized forcesare thus distributed uniformly in this model. A care-ful analysis suggests that this unexpected distributionis in tune with the spatial properties of MK packings.First, note that even though bucklers follow a differentpdf, selecting particles with d + 1 contacts is still a validselection criterion. (If their contribution had not beenisolated, the remaining forces would not follow the MFpower-law scaling given in Eq. (2b), as it does in Fig. 3c,whereas if both kinds of forces are considered together,their joint pdf scales as ≈ .
1, which differs from theanalogous quantity for standard hard spheres [16].) Sec-ond, analyzing the distribution of dot products betweencontact vectors as in Ref. 60 reveals that particles with d + 1 contacts in MK packings have a very similar distri-bution as those in standard hard sphere packings. Buck-lers thus mainly give rise to a localized response thanksto them having three nearly coplanar contacts and onenearly orthogonal force. In order to understand why lo-calized forces are uniformly distributed, we follow Ref. 17,which showed that the two types of contact forces are re-lated to two types of floppy modes: extended forces arerelated to floppy modes that can couple strongly to ex-ternal perturbations, and hence their response is bulkdominated; and buckling forces are associated to floppymodes of a rapidly decaying displacement field. (Thevalue of θ (cid:96) ≈ .
17 was estimated from the statistics of displacements in the latter.) There is therefore a strongconnection between the distribution of forces in bucklersand the particle displacements their floppy modes pro-duce. Now, let us assume that in an MK packing weopen a buckling contact, ( ij ), between particles i and j ,in order to describe the associated displacement field. Inparticular, let us focus on the remaining contacts of anyof these particles, say i . Because of the random shifts,the other particles touching i are (very likely) not con- strained by each other nor by the other particles near i .Instead, the displacement of each neighbor of i is limitedby its own contacts, which are not neighbors themselves,and are typically far apart. By the same token, the effecton the rest of particles in contact with j is determined bysecondary contacts that–with high probability–are dis-tant from each other and from ( ij ). As a result, openinga buckling contact produces a small series of uncorre-lated displacements. No particular length scale is hencefavored over any other. Because of the close relation be-tween localized forces and displacements just mentioned,it is natural for f (cid:96) to be uniformly distributed.Before closing this section, we note that the distribu-tions of f (cid:96) for the FCC and MK packings violate thestability condition related to local excitations given byEq. (3a). We comment further on this point in Sec. VI.For now, we simply note that broader classes of disorderneed to be considered when studying the criticality as-sociated with localized contact forces, even though theirfinite-size effects are unimportant. VI. DISCUSSION
For clarity, we synthesize our results in Table I. Thefirst three rows, which consider the power-law scaling ofthe pdfs in Eqs. (1) and (2), assess the jamming critical-ity associated with microstructural variables for differenttypes of systems. Recall that not only were differentmodels considered, but so was the direction of approachto the jamming point. The systematic corroboration ofthe non-trivial distributions of forces and gaps for fullydisordered systems at jamming fully support the descrip-tion derived from the exact MF theory. Systems with anunderlying FCC symmetry, however, exhibit marked dis-crepancies. Our result thus validate earlier reports thatcrystalline structures fall outside the jamming universal-ity [44, 45] , even though some of its critical features areconserved [73].Our main finding is the contrasted system-size depen-dence of the distribution of gaps and contact forces, assummarized in the last two rows of Table I. Size effectsin p ( f e ) are practically nonexistent for all models, di-mensionality, and interaction type, while g ( h ) exhibitsclear and systematic signatures of finite- N deviationsfrom the expected power-law scaling. Logarithmic cor-rections to g ( h ) are further observed in two-dimensionalsystems. We emphasize that testing for such size scal-ings not only rigorously assesses the critical scaling andits exponents [64, 65], but also provides key insight intothe length scale of their correlations. Hence, we concludethat the MF exponents for all gap and f e distributionsin the MK systems are correct. Yet–leaving aside forthe moment the MK results–a second and more infor-mative conclusion is that the distribution p ( f e ) reachesits thermodynamic limit behavior at smaller values of N than g ( h ). Two different correlation lengths, ξ f e and ξ h ,therefore characterize the relevant length scales of corre-3 Property d = 3 Monodisperse d = 2 Polydisperse FCC MKUC and OC OC OC UC p ( f e ) with θ e = 0 . (cid:51) (cid:51) (cid:51) (but small range) (cid:51) g ( h ) with γ = 0 . (cid:51) (cid:51) (cid:55) : γ FCC (cid:39) . (cid:51) p ( f (cid:96) ) with θ (cid:96) = 0 . (cid:51) (cid:51) (cid:55) : no power law (cid:55) : θ (cid:96),MK = 0Eq. (7) scaling for forces (cid:55) (cid:55) (cid:55) (cid:51) (but small effect)Eq. (7) scaling for gaps (cid:51) (cid:51) (Eq. (11)) (cid:51) (using γ FCC ) (cid:51) Table I. Summary of our main results for the various properties and models considered. In the heading we also indicate ifthe respective jamming point was reached from the under- (UC) or over-compressed (OC) phase. In the first three rows acheck-mark ( (cid:51) ) denotes that the corresponding theoretical prediction was verified, and a cross ( (cid:55) ) that it was not. In the lasttwo rows symbols denote whether the size scaling was verified or not. Results that contradict MF predictions, or results fromprevious studies, are highlighted in red. lations of contact forces and gaps, respectively. For theformer, it must be that N /d (cid:29) ξ f e , and thus finite sizecorrections are negligible, while for the latter it must bethat N /d (cid:46) ξ h . In other words, gaps are correlated oversignificantly larger distances than forces, i.e. ξ h (cid:29) ξ f e .Such disparity in the correlation lengths is a surprisingconsequence of our results, considering that both quan-tities are usually treated on an equal footing from theperspective of the SAT-UNSAT transition in the percep-tron [35, 36], constraint satisfaction problems [41, 69],and neural networks [39] as well as from the point ofview of marginal stability in amorphous solids [31, 32].MK results also fit into this description if we considerthat their very high densities and connectivity reduce theeffective system size, as discussed in Sec. IV. Observingthe scaling of Eq. (7) for the cdf of f e is thus a manifes-tation of the smaller effective volume (for a similar N ),which confirms that finite-size corrections for p ( f e ) arepresent at jamming, but disappear for relatively smallsystem sizes. The significantly more pronounced N de-pendence of the distributions of h (Figs. 3c and 6b) thussupports our finding that ξ h (cid:29) ξ f e .Interestingly, our results further suggest that themarginal stability bounds for the exponents, as expressedin Eqs. (3), should be modified when different types ofdisorder are present. For instance, our findings alongwith other works [44, 45] evince that these inequalitiesare prone to be violated when crystalline lattices are usedto generate the jammed packings. The inherent geome-try of jammed configurations therefore plays a significantrole in formulating general stability criteria. Because thebounds in Eq. (3) were derived [17, 31, 32] assuming,implicitly, that particles positions are uncorrelated, itshould not be overly surprising that γ F CC violates bothrelations.The linear growth of cdf( f (cid:96) ) in the MK model is alsoat odds with the stability condition of Eq. (3a). Thisfinding is more surprising because there is no long-rangeorder in this type of system. At the end of Sec. V weused the peculiar geometry of these packings to suggesta physical explanation for the uniform distribution of f (cid:96) , but this reasoning does not explain why the stability con-dition between γ and θ (cid:96) is apparently violated. Giventhe drastic difference in the inherent structures of theFCC and MK packings, they highlight the need for morestudies to better understand the role played by disorderin determining how the response to external perturba-tions is related to spatial correlations between particlesin jammed systems.The most persistent observations was that all cumula-tive distributions of both gaps and extended forces be-have in a seemingly linear fashion at very small argu-ments, in agreement with the MF predictions, p and g in Eqs. (12) and (13), respectively. Such a cut-off ofthe main power-law scaling is due to the extra contactof isostatic configurations (see Sec. II D) and has beenpreviously reported for the gaps distributions of diskspackings [46], but we are not aware of analogous find-ings in any other model or for the f e distributions. Asdiscussed at the end of Sec. IV, our results suggest thatscalings caused by the additional contact with respect toisostaticity are more robust against changes in the typeof disorder and have a similar characteristic scale in bothtypes of microstructural variables. However, because ofundersampling of the left tails of these distributions, amore stringent analysis would need to be carried out toverify that p ( x ) ∼ g ( x ) ∼ x (cid:28)
1. A previ-ous work on the perceptron [48] also reported a similartransition to a uniform distribution of contact forces thatdepended on the type of algorithm used to reach the jam-ming point, but given that we have used two different al-gorithms to produce our packings, it is unlikely that bothcould produce the same systematic effect. This questionis particularly interesting because it would directly af-fect the robustness of jamming universality, albeit onlyfor the very smallest forces and gaps. Yet, given thatthe left tails of g ( h ) and p ( f e ) determine the smallestgaps and contact forces, accurately describing their truedistribution is key to assessing the stability of jammedpackings away from the thermodynamic limit. We nev-ertheless leave this and other related issues as topics forfuture consideration.4 ACKNOWLEDGMENTS
We want to thank Franceso Zamponi for insightful com-ments and suggestions to our work. RDHR thanks Geor-gios Tsekenis for very useful discussions during the initialstage of this work and Beatriz Seoane for helpful sug-gestions regarding the molecular dynamics simulations of the MK model. This work was supported by the Si-mons Foundation grant ( [1] Giorgio Parisi, Pierfrancesco Urbani, and FrancescoZamponi,
Theory of Simple Glasses: Exact Solutions inInfinite Dimensions (Cambridge University Press, Cam-bridge, 2020).[2] Ludovic Berthier and Giulio Biroli, “Theoretical perspec-tive on the glass transition and amorphous materials,”Reviews of Modern Physics , 587–645 (2011).[3] Andrea J. Liu and Sidney R. Nagel, “The Jamming Tran-sition and the Marginally Jammed Solid,” Annual Re-view of Condensed Matter Physics , 347–369 (2010).[4] Giorgio Parisi and Francesco Zamponi, “Mean-field the-ory of hard sphere glasses and jamming,” Reviews ofModern Physics , 789–845 (2010).[5] S. Torquato and F. H. Stillinger, “Jammed hard-particlepackings: From Kepler to Bernal and beyond,” Reviewsof Modern Physics , 2633–2672 (2010).[6] Patrick Charbonneau, Jorge Kurchan, Giorgio Parisi,Pierfrancesco Urbani, and Francesco Zamponi, “Glassand Jamming Transitions: From Exact Results to Finite-Dimensional Descriptions,” Annual Review of CondensedMatter Physics , 265–288 (2017).[7] Martin van Hecke, “Jamming of soft particles: geometry,mechanics, scaling and isostaticity,” Journal of Physics:Condensed Matter , 033101 (2010).[8] Adrian Baule, Flaviano Morone, Hans J. Herrmann, andHern´an A. Makse, “Edwards statistical mechanics forjammed granular matter,” Reviews of Modern Physics , 015006 (2018).[9] Andrea J. Liu and Sidney R. Nagel, “Nonlinear dynam-ics: Jamming is not just cool any more,” Nature ,21–22 (1998).[10] Corey S. O’Hern, Leonardo E. Silbert, Andrea J. Liu,and Sidney R. Nagel, “Jamming at zero temperature andzero applied stress: The epitome of disorder,” PhysicalReview E , 011306 (2003).[11] Patrick Charbonneau, Eric I. Corwin, Giorgio Parisi,and Francesco Zamponi, “Universal Microstructure andMechanical Stability of Jammed Packings,” Physical Re-view Letters , 205501 (2012).[12] Jorge Kurchan, Giorgio Parisi, and Francesco Zamponi,“Exact theory of dense amorphous hard spheres in highdimension I. The free energy,” Journal of Statistical Me-chanics: Theory and Experiment , P10012 (2012).[13] Jorge Kurchan, Giorgio Parisi, Pierfrancesco Urbani,and Francesco Zamponi, “Exact Theory of Dense Amor-phous Hard Spheres in High Dimension. II. The HighDensity Regime and the Gardner Transition,” The Jour-nal of Physical Chemistry B , 12979–12994 (2013).[14] Patrick Charbonneau, Jorge Kurchan, Giorgio Parisi,Pierfrancesco Urbani, and Francesco Zamponi, “Exacttheory of dense amorphous hard spheres in high dimen-sion. III. The full replica symmetry breaking solution,” Journal of Statistical Mechanics: Theory and Experi-ment , P10009 (2014).[15] Patrick Charbonneau, Jorge Kurchan, Giorgio Parisi,Pierfrancesco Urbani, and Francesco Zamponi, “Frac-tal free energy landscapes in structural glasses,” NatureCommunications , 3725 (2014).[16] Patrick Charbonneau, Eric I. Corwin, Giorgio Parisi,and Francesco Zamponi, “Jamming Criticality Revealedby Removing Localized Buckling Excitations,” PhysicalReview Letters , 125504 (2015).[17] Edan Lerner, Gustavo D¨uring, and Matthieu Wyart,“Low-energy non-linear excitations in sphere packings,”Soft Matter , 8252–8263 (2013).[18] Harukuni Ikeda, “Jamming below upper critical dimen-sion,” Physical Review Letters , 038001 (2020).[19] Yuxiao Zhang, M. J. Godfrey, and M. A. Moore,“Marginally jammed states of hard disks in aone-dimensional channel,” arXiv:2005.02901 [cond-mat](2020).[20] Atsushi Ikeda, Ludovic Berthier, and Giulio Biroli, “Dy-namic criticality at the jamming transition,” The Journalof Chemical Physics , 12A507 (2013).[21] E. DeGiuli, E. Lerner, and M. Wyart, “Theory of thejamming transition at finite temperature,” The Journalof Chemical Physics , 164503 (2015).[22] Cristian F. Moukarzel, “Isostatic Phase Transition andInstability in Stiff Granular Materials,” Physical ReviewLetters , 1634–1637 (1998).[23] Carl P. Goodrich, Andrea J. Liu, and Sidney R. Nagel,“Finite-Size Scaling at the Jamming Transition,” Physi-cal Review Letters , 095704 (2012).[24] Carl P. Goodrich, Simon Dagois-Bohy, Brian P. Tighe,Martin van Hecke, Andrea J. Liu, and Sidney R. Nagel,“Jamming in finite systems: Stability, anisotropy, fluc-tuations, and scaling,” Physical Review E , 022138(2014).[25] Carl P. Goodrich, Andrea J. Liu, and James P. Sethna,“Scaling ansatz for the jamming transition,” Proceed-ings of the National Academy of Sciences , 9745–9750(2016).[26] Carl P. Goodrich, Wouter G. Ellenbroek, and AndreaJ. Liu, “Stability of jammed packings I: the rigiditylength scale,” Soft Matter , 10993–10999 (2013).[27] Samuel S. Schoenholz, Carl P. Goodrich, Oleg Kogan,Andrea J. Liu, and Sidney R. Nagel, “Stability ofjammed packings II: the transverse length scale,” SoftMatter , 11000–11006 (2013).[28] Daniel Hexner, Andrea J. Liu, and Sidney R. Nagel,“Two Diverging Length Scales in the Structure ofJammed Packings,” Physical Review Letters , 115501(2018). [29] Daniel Hexner, Pierfrancesco Urbani, and FrancescoZamponi, “Can a Large Packing be Assembled fromSmaller Ones?” Physical Review Letters , 068003(2019).[30] Eric DeGiuli, Edan Lerner, Carolina Brito, andMatthieu Wyart, “Force distribution affects vibrationalproperties in hard-sphere glasses,” Proceedings of the Na-tional Academy of Sciences , 17054–17059 (2014).[31] Matthieu Wyart, “Marginal Stability Constrains Forceand Pair Distributions at Random Close Packing,” Phys-ical Review Letters , 125502 (2012).[32] Markus M¨uller and Matthieu Wyart, “Marginal Stabil-ity in Structural, Spin, and Electron Glasses,” AnnualReview of Condensed Matter Physics , 177–200 (2015).[33] Ludovic Berthier, Patrick Charbonneau, Yuliang Jin,Giorgio Parisi, Beatriz Seoane, and Francesco Zamponi,“Growing timescales and lengthscales characterizing vi-brations of amorphous solids,” Proceedings of the Na-tional Academy of Sciences , 8397–8401 (2016).[34] Ludovic Berthier, Giulio Biroli, Patrick Charbonneau,Eric I. Corwin, Silvio Franz, and Francesco Zamponi,“Gardner physics in amorphous solids and beyond,” TheJournal of Chemical Physics , 010901 (2019).[35] Silvio Franz and Giorgio Parisi, “The simplest model ofjamming,” Journal of Physics A: Mathematical and The-oretical , 145001 (2016).[36] Silvio Franz, Antonio Sclocchi, and Pierfrancesco Ur-bani, “Critical jammed phase of the linear perceptron,”Physical Review Letters , 115702 (2019).[37] Stefano Spigler, Mario Geiger, St´ephane d’Ascoli, LeventSagun, Giulio Biroli, and Matthieu Wyart, “A jammingtransition from under- to over-parametrization affectsloss landscape and generalization,” arXiv:1810.09665[cond-mat, stat] (2018).[38] Mario Geiger, Stefano Spigler, St´ephane d’Ascoli, LeventSagun, Marco Baity-Jesi, Giulio Biroli, and MatthieuWyart, “Jamming transition as a paradigm to under-stand the loss landscape of deep neural networks,” Phys-ical Review E , 012115 (2019).[39] Silvio Franz, Sungmin Hwang, and Pierfrancesco Urbani,“Jamming in Multilayer Supervised Learning Models,”Physical Review Letters , 160602 (2019).[40] Fabrizio Antenucci, Silvio Franz, Pierfrancesco Urbani,and Lenka Zdeborov´a, “Glassy Nature of the Hard Phasein Inference Problems,” Physical Review X , 011020(2019).[41] Silvio Franz, Giorgio Parisi, Maxim Sevelev, Pier-francesco Urbani, and Francesco Zamponi, “Univer-sality of the SAT-UNSAT (jamming) threshold in non-convex continuous constraint satisfaction problems,” Sci-Post Physics , 019 (2017).[42] Florent Krzakala and Jorge Kurchan, “Landscape analy-sis of constraint satisfaction problems,” Physical ReviewE , 021122 (2007).[43] Silvio Franz, Antonio Sclocchi, and Pierfrancesco Ur-bani, “Critical energy landscape of linear soft spheres,”SciPost Physics , 012 (2020).[44] Patrick Charbonneau, Eric I. Corwin, Lin Fu, GeorgiosTsekenis, and Michael van der Naald, “Glassy, Gardner-like phenomenology in minimally polydisperse crystallinesystems,” Physical Review E , 020901 (2019).[45] Georgios Tsekenis, “Jamming Criticality of Near-Crystals,” arXiv:2006.07373 [cond-mat] (2020). [46] Harukuni Ikeda, Carolina Brito, and Matthieu Wyart,“Infinitesimal asphericity changes the universality of thejamming transition,” Journal of Statistical Mechanics:Theory and Experiment , 033302 (2020).[47] Carolina Brito, Harukuni Ikeda, Pierfrancesco Urbani,Matthieu Wyart, and Francesco Zamponi, “Universalityof jamming of nonspherical particles,” Proceedings of theNational Academy of Sciences , 11736–11741 (2018).[48] Yoav Kallus, “Scaling collapse at the jamming transi-tion,” Physical Review E , 012902 (2016).[49] Romain Mari and Jorge Kurchan, “Dynamical transitionof glasses: From exact to approximate,” The Journal ofChemical Physics , 124504 (2011).[50] Corey S. O’Hern, Stephen A. Langer, Andrea J. Liu, andSidney R. Nagel, “Random Packings of Frictionless Par-ticles,” Physical Review Letters , 075507 (2002).[51] Adam B. Hopkins, Frank H. Stillinger, and SalvatoreTorquato, “Disordered strictly jammed binary spherepackings attain an anomalously large range of densities,”Physical Review E , 022205 (2013).[52] Varda F. Hagh, Eric I. Corwin, Kenneth Stephenson,and M. F. Thorpe, “A broader view on jamming: Fromspring networks to circle packings,” Soft Matter ,3076–3084 (2019).[53] Aleksandar Donev, Salvatore Torquato, and Frank H.Stillinger, “Pair correlation function characteristics ofnearly jammed disordered and ordered hard-sphere pack-ings,” Physical Review E , 011105 (2005).[54] Erik Bitzek, Pekka Koskinen, Franz G¨ahler, MichaelMoseler, and Peter Gumbsch, “Structural RelaxationMade Simple,” Physical Review Letters , 170201(2006).[55] Peter K. Morse and Eric I. Corwin, “Geometric Signa-tures of Jamming in the Mechanical Vacuum,” PhysicalReview Letters , 115701 (2014).[56] Patrick Charbonneau, Eric I. Corwin, Giorgio Parisi,Alexis Poncet, and Francesco Zamponi, “Universal Non-Debye Scaling in the Density of States of AmorphousSolids,” Physical Review Letters , 045503 (2016).[57] Peter K. Morse and Eric I. Corwin, “Echoes of the GlassTransition in Athermal Soft Spheres,” Physical ReviewLetters , 118003 (2017).[58] Monica Skoge, Aleksandar Donev, Frank H. Stillinger,and Salvatore Torquato, “Packing hyperspheres in high-dimensional Euclidean spaces,” Physical Review E ,041127 (2006).[59] Claudia Artiaco, Paolo Baldan, and Giorgio Parisi, “Ex-ploratory study of the glassy landscape near jamming,”Physical Review E , 052605 (2020).[60] Rafael D´ıaz Hern´andez Rojas, Giorgio Parisi, andFederico Ricci-Tersenghi, “Inferring the particle-wisedynamics of amorphous solids from the local struc-ture at the jamming point,” Soft Matter (2020),10.1039/C9SM02283J.[61] Patrick Charbonneau, Yuliang Jin, Giorgio Parisi, andFrancesco Zamponi, “Hopping and the Stokes–Einsteinrelation breakdown in simple glass formers,” Proceedingsof the National Academy of Sciences , 15025–15030(2014).[62] Patrick Charbonneau, Yuliang Jin, Giorgio Parisi, Cor-rado Rainone, Beatriz Seoane, and Francesco Zam-poni, “Numerical detection of the Gardner transition in amean-field glass former,” Physical Review E , 012316(2015). [63] M. E. J. Newman, “Power laws, Pareto distributions andZipf’s law,” Contemporary Physics , 323–351 (2005).[64] Daniel J Amit and Victor Martin-Mayor, Field the-ory, the renormalization group, and critical phenomena:graphs to computers (World Scientific Publishing Com-pany, 2005).[65] M Newman and G Barkema,
Monte Carlo methods instatistical physics (Oxford University Press: New York,USA, 1999).[66] J. Wang and A. P. Young, “Monte Carlo study of thesix-dimensional Ising spin glass,” Journal of Physics A:Mathematical and General , 1063 (1993).[67] Juan J. Ruiz-Lorenzo, “Logarithmic corrections for spinglasses, percolation and Lee-Yang singularities in six di-mensions,” Journal of Physics A: Mathematical and Gen-eral , 8773–8787 (1998).[68] R Kenna, “Finite size scaling for O(N) Φ4-theory at theupper critical dimension,” Nuclear Physics B , 292–304 (2004). [69] Harukuni Ikeda, Pierfrancesco Urbani, and FrancescoZamponi, “Mean field theory of jamming of nonspheri-cal particles,” Journal of Physics A: Mathematical andTheoretical , 344001 (2019).[70] C Lucibello, F Morone, G Parisi, F Ricci-Tersenghi, andTommaso Rizzo, “Anomalous finite size corrections inrandom field models,” Journal of Statistical Mechanics:Theory and Experiment , P10025 (2014).[71] C. Lucibello, F. Morone, G. Parisi, F. Ricci-Tersenghi,and Tommaso Rizzo, “Finite-size corrections to disor-dered Ising models on random regular graphs,” PhysicalReview E , 012146 (2014).[72] U. Ferrari, C. Lucibello, F. Morone, G. Parisi, F. Ricci-Tersenghi, and T. Rizzo, “Finite-size corrections to dis-ordered systems on Erd¨os-R´enyi random graphs,” Phys-ical Review B , 184201 (2013).[73] Harukuni Ikeda, “Jamming and replica symmetry break-ing of weakly disordered crystals,” Physical Review Re-search2