Finding Point J figure by figure in O’Hern et al. (2003)

Start with the protocol: where the states come from

Section II A and Eq. (1) define the model: two particles repel only when they overlap, with energy $\varepsilon(1-r_{ij}/\sigma_{ij})^{\alpha}/\alpha$ . In 2D they use a 50–50 mixture of disks with size ratio 1.4, which prevents crystallisation, and particle numbers from 4 to 4096.

States start at infinite temperature: centres placed completely at random in a periodic box, then conjugate gradients walk downhill in energy. Minimisation stops when the energy per particle falls below $10^{-16}$ or changes by less than $10^{-15}$ between iterations. The first condition is the test for “unjammed”, equivalent to a pressure below $10^{-10}$ .

Keep in mind that “jammed” in this paper is a statement about this protocol. Another way of making states, such as slow cooling, gives other thresholds. Every distribution and extrapolation later in the paper is a statistic over these random starts.

Source: C. S. O’Hern

How one state’s threshold is measured (Section II B)

For each zero-temperature state: if the pressure is zero, every particle is grown by the same factor; if not, every particle is shrunk; after each step the energy is minimised again, until the pressure just becomes non-zero. Steps range from $10^{-6}$ to $10^{-4}$ , smaller near the threshold, and are refined until the final configuration no longer depends on the step.

The key finding of this section is that the bulk and shear moduli vanish at the same packing fraction. Resisting compression and resisting shear begin together, so each state has a well-defined threshold. The authors also checked the dynamical matrix: at least 1e-6 above the threshold, the only zero-frequency modes belong to rattlers and to uniform translations.

Watch out at Fig. 3. In the arXiv v1 we read, the caption gives reference slopes of 2.0 and 2.5, while Table I gives pressure exponents of 1.01 and 1.50, the same values Chapter 4 derived from $p\propto\Delta\phi^{\alpha-1}$ . The values 2.0 and 2.5 are the energy exponent $\alpha$ . We checked pressure against Table I; the caption of the published version is one we have not seen.

Source: C. S. O’Hern

From the jammed fraction to the threshold distribution (Figs. 5 to 8)

Fig. 5 shows $f_j(\phi)$ , the fraction of random starts already jammed at a given packing fraction: a step rising from 0 to 1. Fig. 6 is its derivative, the threshold distribution $P_j(\phi_c)$ . The distributions for three contact exponents lie on top of each other within error: the threshold is set by geometry alone, not by how stiff the contacts are.

Fig. 7 plots the full width at half maximum against particle number. Read it in two parts: for very few particles the width first grows with size, and only past its maximum does it fall as Eq. (3). The values $w_0=0.16\pm0.04$ and $\Omega=0.55\pm0.03$ apply to the second part only.

Fig. 8 deals with the position of the peak, $\phi_0$ . It approaches a limit $\phi^{*}$ as a power of the system length $L=N^{1/d}$ , Eq. (4), with fitted $\nu=0.71\pm0.08$ and $\delta_0=0.12\pm0.03$ . So $\phi^{*}$ is an extrapolation, not something any single simulation measured.

w=w_0\,N^{-\Omega},\qquad \phi_0-\phi^{*}=\delta_0\,L^{-1/\nu},\qquad L=N^{1/d}

Eqs. (3) and (4) of the paper. The width going to zero and the peak going to

\phi^{*}

together say that in an infinite system nearly every random start jams at the same packing fraction. Fig. 8 plots the absolute value of the difference.

Source: C. S. O’Hern

Contacts and Table I: which numbers rest on an assumption

Section II E fits Fig. 9 with Eq. (6), $Z-Z_c=Z_0(\phi-\phi_c)^{\zeta}$ , and the text gives $\zeta=0.50\pm0.03$ . The $Z_c=3.98$ in Table I is a different matter: the authors first remove rattlers (about 5%) and then assume $\zeta=0.5$ to obtain this “precise” value. So $Z_c$ agreeing with 4 does not independently prove that the exponent is one half.

Goodrich, Liu and Nagel (2012) refine the count: N disks in a periodic box need at least $2N-1$ contacts for a positive bulk modulus, and the square-root law for the contact number fails when $p^{1/2}N$ approaches one. Read Table I of the 2003 paper with that finite-size limit in mind.

Our recomputation: on the load-bearing network at the threshold, the contacts beyond the minimum count are 0 for 32 and 64 disks and at most 1 for 128 and 256 disks. Rattlers average 4.0% to 5.0%, in line with the paper’s roughly 5%. The contact-number exponent at 128 disks came out as 0.577 ± 0.009, against 0.49 in the paper.

Source: C. S. O’Hern · C. P. Goodrich

Worked example: checking our threshold statistics against Eqs. (3) and (4)

Question: in periodic boxes we took 24 random starts each for 32, 64, 128 and 256 disks. The full widths of the threshold distributions are 0.0259, 0.0148, 0.0129 and 0.00654, and the mean thresholds 0.827, 0.835, 0.838 and 0.839. Are they consistent with Eqs. (3) and (4), and roughly where would an infinite system jam?

Step one, the widths. The central values of Eq. (3), $0.16\times N^{-0.55}$ , give 0.0238, 0.0162, 0.0111 and 0.00758. Taking the ends of the uncertainties on $w_0$ and $\Omega$ , the range allowed at 256 disks is 0.0048 to 0.0112, and our 0.00654 lies inside; the other three sizes also lie inside their ranges. Note that we estimate the full width as the standard deviation times $2\sqrt{2\ln2}$ , which assumes a Gaussian shape.

Step two, the extrapolation. In 2D $L=\sqrt{N}$ , so 256 disks give $L=16$ , and the central values of Eq. (4) give $0.12\times16^{-1/0.71}=0.0024$ . The mean threshold rises with size, so the peak approaches the limit from below and $\phi^{*}\approx0.839+0.0024=0.841$ . Within the uncertainties of $\delta_0$ and $\nu$ , the correction lies between 0.0011 and 0.0045.

Conclusion: the widths are consistent with Eq. (3), and the extrapolated value of about 0.841 sits near the large-system value of 0.84 in Fig. 5(a). Two approximations remain: Eq. (4) concerns the peak of the distribution while we used the mean of 24 samples, and 256 disks is itself only a moderate size.

Source: C. S. O’Hern

A checklist for papers like this one

First, find the protocol: how states are made, what decides “jammed” and how tight the convergence is. Second, separate measurement from fitting: the location of Point J and the limiting contact number come from extrapolation or assumption. Third, find the range of validity: Eq. (3) holds for large systems only, and by the captions of Figs. 3 and 9 the 2D scaling data come from 1024 disks. Fourth, check captions against the text: when they disagree, as at Fig. 3, go with the table and the derivation, and note what you did not verify.

Fifth, when recomputing, state what counts as agreement. We required the widths to fall within the range allowed by the uncertainties of Eq. (3), and the contact exponent to fall within the $\pm0.04$ of Table I. The first passed; the second, 0.577 against 0.49, did not. Both results are on record.

Source: C. S. O’Hern · C. P. Goodrich

Try it in the Lab

Find the Table I quantities in a barely jammed packing

Predict first: just past the threshold, are there many or few contacts beyond isostatic? Roughly what share of the disks are rattlers?
Choose the Barely jammed preset, play to the deepest packing fraction and read “Contacts beyond isostatic”; then read the numbers of rattlers and load-bearing disks under “More numbers and how they are computed”.
Repeat with the Deep squeeze preset. Which quantity changes clearly with compression? Match it to Z − Zc in Eq. (6), and note that the Lab uses a walled box, so its onset differs from the paper’s periodic box.

Exercises with explained answers

With the central values of Eq. (3), how many disks are needed for the full width of the threshold distribution to fall below 0.001?
Show the answer
$N=(0.16/0.001)^{1/0.55}=160^{1.82}\approx1.0\times10^{4}$ . That is an extrapolation: the paper’s simulations go up to 4096 particles.
Table I gives Zc = 3.98, within its uncertainty of the isostatic value 4. Why does that not independently prove ζ = 1/2?
Show the answer
The authors obtained that Zc after assuming ζ = 0.5 and removing rattlers. A different exponent would extrapolate to a different Zc. An independent test would fit Zc and ζ together, or count contacts directly as Goodrich and co-workers did.
Our 128-disk data give a contact-number exponent of 0.577; the paper gives 0.49. Name two possible causes and say which one we have ruled out.
Show the answer
One is finite size: the square-root law fails when $p^{1/2}N$ approaches one. We used only data with $p^{1/2}N\ge3$ and the exponent is still too large, so this cause is not fully ruled out but cannot explain the gap alone. The other is the reference count: the paper removes rattlers and assumes ζ = 0.5 to find Zc, while we use the minimum count of Goodrich and co-workers, and a different reference point changes the fit. Separating the two needs larger systems.

Sources and what we recomputed

C. S. O’Hern, L. E. Silbert, A. J. Liu and S. R. Nagel, Jamming at zero temperature and zero applied stress: the epitome of disorder, arXiv:cond-mat/0304421v1; Physical Review E 68, 011306 (2003)
Version 1, 17 April 2003. Eq. (1): the finite-range repulsion between soft spheres. Section II: random infinite-temperature starts, conjugate-gradient minimisation in a periodic cell, and the energy test for an unjammed state. Figs. 5 to 7 and Eq. (3): the fraction of jammed states, the distribution of thresholds and its width, which falls as the number of particles to the power −0.55. Table I: pressure and contact-number exponents and prefactors. Fig. 20: a 2D force network. We recomputed the threshold width, the pressure exponent and the contact-number exponent with our own code; the contact-number exponent did not reproduce at 128 disks.
C. P. Goodrich, A. J. Liu and S. R. Nagel, Finite-size scaling at the jamming transition, arXiv:1204.4915v3; Physical Review Letters 109, 095704 (2012)
Version 3, 6 June 2012. Eq. (2): the contact count needed to remove every zero mode except global translations. Eq. (3): one more contact is needed for a positive bulk modulus, so a periodic packing of N disks in two dimensions needs at least 2N − 1 contacts. Eqs. (4) and (5): the corresponding minimum contact number. The finite-size crossover where the usual power laws break down at low pressure. We checked the 2N − 1 count on our own threshold packings.