How to get the pressure from contact forces

Stress and pressure in two dimensions

In 2D a force acts across a line, so stress is force per length, $\varepsilon/d^2$ . $\sigma_{xx}$ is the push per unit length on a wall whose normal points along $x$ , and likewise $\sigma_{yy}$ . This course takes compressive stress as positive, because disks can only push.

The pressure is the average of the two normal stresses, $p=(\sigma_{xx}+\sigma_{yy})/2$ . If $\sigma_{xx}\ne\sigma_{yy}$ or $\sigma_{xy}\ne0$ , the packing is loaded differently in different directions; shear in Chapter 5 relies on that anisotropy.

Source: C. S. O’Hern

Stress computed inside the packing

For each contact take the branch vector $\mathbf b$ : from centre $j$ to centre $i$ for a disk–disk contact, and from the contact point on the wall to the centre for a wall contact. $\mathbf f$ is the force of that contact on disk $i$ . Contacts only push, so $\mathbf f$ points along $\mathbf b$ , $b_xf_x>0$ , and compression comes out positive.

Eq. (12) of Silbert et al. adds a kinetic part from velocity fluctuations. A static packing has none, so only the contact part remains. Eq. (2) of O’Hern et al. defines the pressure tensor from the same sum.

\sigma_{ab}=\frac{1}{A}\sum_{\text{contacts}} b_a\,f_b,\qquad p=\frac{\sigma_{xx}+\sigma_{yy}}{2}=\frac{1}{2A}\sum_{\text{contacts}}\mathbf b\cdot\mathbf f

A

is the box area and

a,b

stand for

x

y

. Each contact is counted once.

Source: L. E. Silbert · C. S. O’Hern

Why it equals the force on a wall

Let the box have width $L$ and height $H$ , with the left wall at $x=0$ and the right wall at $x=L$ . Every disk is balanced. Multiply each disk’s $x$ balance equation by its coordinate $x_i$ and add over all disks: the result is still zero.

Now regroup the sum by contacts. A disk–disk contact appears on both of its disks with opposite forces, giving $(x_i-x_j)f_x=b_xf_x$ . A left-wall contact gives $x_if_x$ , which is exactly $b_xf_x$ . A right-wall contact gives $x_if_x=(x_i-L)f_x+Lf_x=b_xf_x+Lf_x$ .

The $f_x$ of all right-wall contacts add up to the wall’s total push on the disks, $-F_R$ , where $F_R$ is the disks’ total push on the right wall. So the sum of $b_xf_x$ over all contacts equals $LF_R$ . Dividing by the area $LH$ gives $F_R/H$ : the force per unit length on the right wall.

\sum_i x_i F^{\,\mathrm{net}}_{i,x}=0\;\Longrightarrow\;\sum_{\text{contacts}}b_xf_x=L\,F_R\;\Longrightarrow\;\sigma_{xx}=\frac{F_R}{H}

Only the balance of each disk was used, not the contact law. While disks are still out of balance the two numbers disagree, so their difference tests whether a relaxation has finished.

Pressure is also an energy derivative

Scale the box and every centre by the same factor $1+e$ . The area becomes $A(1+e)^2$ , an increase $dA=2Ae$ . Every centre distance is also multiplied by $1+e$ , so the energy changes by $dU=\sum (dV/dr)\,r\,e=-\sum f\,r\,e$ .

Dividing, $-dU/dA=\sum f r/(2A)$ , which is the $p$ above. That lets us check the pressure from the energy alone. On a packing of 128 disks, a central difference of $-dU/dA$ agrees with the contact stress to a relative difference of at most $1.5\times10^{-7}$ ; the contact stress agrees with the forces on all four walls to about $10^{-12}$ .

p=-\frac{dU}{dA}\bigg|_{\text{affine}}=\frac{1}{2A}\sum_{\text{contacts}} f\,r

“Affine” means every coordinate scales by the same factor, with no extra rearrangement of the disks. If the disks are allowed to adjust, the packing relaxes further and the energy released by that relaxation must be subtracted; Chapter 6 uses this to find the energy lost over a compression and unloading cycle.

Source: C. S. O’Hern

Worked example: match all three

Take the two disks of Chapter 1: $R_1=0.5$ , $R_2=0.7$ , walls $L=2.39$ apart, contact force $f=4.59\times10^{-3}\,\varepsilon/d$ . Make the box $H=1.5$ tall so the large disk does not touch the top or bottom.

Force on the wall: $\sigma_{xx}=f/H=3.06\times10^{-3}\,\varepsilon/d^2$ and $\sigma_{yy}=0$ , so $p=1.53\times10^{-3}\,\varepsilon/d^2$ . Squeezed along one direction only, the pressure is half of $\sigma_{xx}$ .

Contact stress: the $x$ components of the three branch vectors are $R_1-\delta_{1w}=0.4989$ , $\sigma-\delta_{12}=1.1934$ and $R_2-\delta_{2w}=0.6978$ , which add to exactly $L=2.39$ . All three contacts carry the same force, so $\sum b_xf_x=fL$ , and dividing by $LH$ gives $f/H$ again.

Energy derivative: as $L$ shrinks, all three overlaps grow together. With the compliance of Chapter 1, $U=\varepsilon(L_0-L)^2/(2\times2.18)$ with $L_0=2.4$ . Along $x$ , $-dU/dL=f$ , and dividing by $H$ gives $\sigma_{xx}$ .

Try it in the Lab

Read one pressure two ways

Predict first: in a step that has not fully relaxed, would the wall pressure and the contact stress agree?
Play the default preset and, at several packing fractions, open “More numbers” and read the ratio of contact stress to wall pressure.
Switch the contact law to Hertzian, repeat, and explain why the ratio does not depend on the contact law.

Exercises with explained answers

A square box of side $L=10$ measures a total push of $0.05\,\varepsilon/d$ on each of its four walls. Find $\sigma_{xx}$ , $\sigma_{yy}$ and $p$ .

Show the answer

Force per length on each wall is $0.05/10=5\times10^{-3}$ , so $\sigma_{xx}=\sigma_{yy}=p=5\times10^{-3}\,\varepsilon/d^2$ .

If you use the convention that compressive stress is negative, how does the contact stress formula change?

Show the answer

Flip the overall sign: $\sigma_{ab}=-\frac1A\sum b_af_b$ . The physics is the same; compression simply reads negative. Check which convention a paper uses before comparing numbers.

A packing has energy $2.000\times10^{-3}\,\varepsilon$ at $A=100$ and $2.052\times10^{-3}\,\varepsilon$ after an affine compression to $A=99.9$ . Estimate the pressure.

Show the answer

$p\approx-\Delta U/\Delta A=-(0.052\times10^{-3})/(-0.1)=5.2\times10^{-4}\,\varepsilon/d^2$ . A one-sided difference gives the average over that small change of area; a central difference is more accurate.

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.

L. E. Silbert, D. Ertaş, G. S. Grest, T. C. Halsey, D. Levine and S. J. Plimpton, Granular flow down an inclined plane: Bagnold scaling and rheology, arXiv:cond-mat/0105071v1; Physical Review E 64, 051302 (2001)

Version 1, 3 May 2001. Section II: the frictional contact model, a normal spring with damping plus a tangential spring cut off by the Coulomb condition, with the ratio of tangential to normal stiffness of 2/7. Eq. (12): the stress tensor as a sum of branch vector times contact force divided by volume, plus a kinetic part. We use the contact part of Eq. (12) for static packings.