Shear response of granular packings compressed above jamming onset

We investigate the mechanical response of jammed packings of repulsive, frictionless spherical particles undergoing isotropic compression. Prior simulations of the soft-particle model, where the repulsive interactions scale as a power law in the interparticle overlap with exponent $\alpha$, have found that the ensemble-averaged shear modulus $\langle G\left(P\right)\rangle$ increases with pressure $P$ as $\sim P^{(\alpha -3/2)/(\alpha -1)}$ at large pressures. $\langle G\rangle$ has two key contributions: (1) continuous variations as a function of pressure along geometrical families, for which the interparticle contact network does not change, and (2) discontinuous jumps during compression that arise from changes in the contact network. Using numerical simulations, we show that the form of the shear modulus $G^f$ for jammed packings within near-isostatic geometrical families is largely determined by the affine response $G^f\sim G_a^f$, where $G_a^f/G_{a0}={(P/P_0)}^{(\alpha -2)/(\alpha -1)}-P/P_0, P_0\sim N^{-2(\alpha -1)}$ is the characteristic pressure at which $G_a^f=0, G_{a0}$ is a constant that sets the scale of the shear modulus, and $N$ is the number of particles. For near-isostatic geometrical families that persist to large pressures, deviations from this form are caused by significant nonaffine particle motion. We further show that the ensemble-averaged shear modulus $\langle G\left(P\right)\rangle$ is not simply a sum of two power laws, but $\langle G\left(P\right)\rangle \sim {(P/P_c)}^a$, where $a\approx (\alpha -2)/(\alpha -1)$ in the $P\to 0$ limit and $\langle G\left(P\right)\rangle \sim {(P/P_c)}^b$, where $b\gtrsim (\alpha -3/2)/(\alpha -1)$, above a characteristic pressure that scales as $P_c\sim N^{-2(\alpha -1)}$.

Figure 1

Snapshots of $N=32$ jammed disk packings with repulsive linear spring interactions [$\alpha =2$ in Eq. (1)] undergoing isotropic compression before and after a point change where a new contact is added to the contact network (indicated by thin blue lines). In (a), we show an isostatic packing at $P={10}^{-7}$ with $N_c=N_c^0=59$ contacts. (b) The packing in (a) has been compressed to $P=4.60\times {10}^{-6}$ without any changes in the contact network. The vectors indicate the displacements of the particle centers relative to the packing in (a) after multiplying by 100. (c) The packing in (c) has been compressed to $P=4.65\times {10}^{-6}$, which results in the addition of one contact, $N_c=N_c^0+1$, indicated by the thick black line. The shaded particles indicate rattlers.

Figure 2

(a) Shear modulus $G^{\left(1\right)}$ vs pressure $P$ (black points) for isostatic packings with $N=32$ disks within 50 geometrical families that maintain their interparticle contact networks for purely repulsive linear ($\alpha =2$) spring interactions. The blue lines give the affine contribution $G_a^{\left(1\right)}$ [Eq. (2)] for each geometrical family. Blue open circles (red crosses) at the end of $G^{\left(1\right)}$ indicate that the packing experienced a point (jump) change at that particular pressure. In this panel, we do not show packings with $G^{\left(1\right)}<0$. (b) $(G^{\left(1\right)}/G_0){(P/P_0)}^{-1/2}$ for isostatic geometrical families plotted vs $P/P_0$ for disk packings with repulsive linear ($\alpha =2$) spring interactions, for several system sizes: $N=32$ (black upper triangles), 64 (blue circles), 128 (red diamonds), 256 (green downward triangles), and 512 (magenta asterisks). The dashed line is Eq. (4) for $\alpha =2$ and the solid lines are fits to Eq. (5). In the inset, we show $\mathrm{\Delta }{G}^{\left(1\right)}=G^{\left(1\right)}(P/P_0)-G^{\left(1\right)}\left(0\right)$ (black upward triangles), $\mathrm{\Delta }{G}_n^{\left(1\right)}=G_n^{\left(1\right)}(P/P_0)-G_n^{\left(1\right)}\left(0\right)$ (red dots), and $\mathrm{\Delta }{G}_a^{\left(1\right)}=G_a^{\left(1\right)}(P/P_0)-G_a^{\left(1\right)}\left(0\right)$ (blue $\times$) with best fits to Eq. (2) (black, red, and blue solid lines, respectively) for an example $N=32$ packing with $\alpha =2$. (c) $G^{\left(1\right)}$ vs $P$ for the packings in (a) plotted on linear scales (with $G^{\left(1\right)}<0$). The blue lines indicate best fits to Eq. (5) for each geometrical family.

Figure 3

(a) Shear modulus $G^{\left(1\right)}$ vs pressure $P$ (black points) for isostatic packings with $N=32$ disks within 50 geometrical families that maintain their interparticle contact networks for Hertzian ($\alpha =5/2$) spring interactions. The blue lines give the affine contribution $G_a^{\left(1\right)}$ [Eq. (2)] for each geometrical family. Blue open circles (red crosses) at the end of each $G^{\left(1\right)}$ indicate that the packing experienced a point (jump) change at that particular pressure. We do not show packings with $G^{\left(1\right)}<0$. (b) $(G^{\left(1\right)}/G_0){(P/P_0)}^{-2/3}$ for isostatic geometrical families plotted vs $P/P_0$ for disk packings with Hertzian ($\alpha =5/2$) spring interactions, for several system sizes: $N=32$ (black upper triangles), 64 (blue circles), 128 (red diamonds), 256 (green downward triangles), and 512 (magenta asterisks). The dashed line is Eq. (4) for $\alpha =5/2$ and the solid lines are fits to Eq. (5) for $\alpha =5/2$. In the inset, we show $G^{\left(1\right)}$ (black upward triangles), $G_n^{\left(1\right)}$ (red dots), and $G_a^{\left(1\right)}$ (blue $\times$) with best fits to Eq. (2) (black, red, and blue solid lines, respectively) for an example $N=32$ packing with $\alpha =5/2$. (c) $G^{\left(1\right)}$ vs $P$ for the packings in (a) plotted on a linear scale (and including packings with $G^{\left(1\right)}<0$). The blue lines indicate best fits to Eq. (5) for each geometrical family.

Figure 4

Probability distributions of the characteristic pressure $P_0\left(\alpha \right)$ for (a) $\alpha =2$ and (b) $5/2$ for $N=32$ (black upward triangles), 128 (blue circles), 512 (red diamonds). [$P_0\left(\alpha \right)$ was obtained using a best fit of the data for $G^{\left(1\right)}$ to Eq. (4).] The insets to (a) and (b) display $\langle P_0\rangle$ (averaged over geometrical families) vs system size $N$ for $\alpha =2$ and $5/2$, respectively. The dashed lines have slopes equal to $-2$ and $-3$ in the insets to (a) and (b).

Figure 5

Probability distributions of the characteristic shear modulus $G_{a0}\left(\alpha \right)$ for (a) $\alpha =2$ and (b) $5/2$ for $N=32$ (black upward triangles), 128 (blue circles), 512 (red diamonds). [$G_{a0}\left(\alpha \right)$ was obtained using a best fit of the data for $G^{\left(1\right)}$ to Eq. (4).] The insets to (a) and (b) display $\langle G_{a0}\rangle$ (averaged over geometrical families) vs system size $N$ for $\alpha =2$ and $5/2$, respectively. The dashed lines have slopes equal to $-1$ and $-2$ in the insets to (a) and (b).

Figure 6

Probability distributions of the characteristic pressure $P_1\left(\alpha \right)$ for (a) $\alpha =2$ and (b) $5/2$ for $N=32$ (black upward triangles), 128 (blue circles), 512 (red diamonds). [See Eq. (5) for the definition of $P_1\left(\alpha \right)$.] The insets to (a) and (b) display $\langle P_1\rangle$ (averaged over geometrical families) vs system size $N$ for $\alpha =2$ and $5/2$, respectively. The dashed lines have slopes equal to $-2$ and $-3$ in the insets to (a) and (b).

Figure 7

Shear modulus $G^{\left(i\right)}$ of a series of $N=32$ disk packings with repulsive (a) linear and (b) Hertzian spring interactions as the system undergoes a point change during isotropic compression (at $P\approx 1.29\times {10}^{-5}$ for $\alpha =2$ and $P=1.37\times {10}^{-6}$ for $\alpha =5/2$ indicated by vertical dashed lines) from the isostatic (black upward triangles) geometrical family with $N_c^0$ contacts to the second geometrical family (blue circles) with $N_c^0+1$ contacts.

Figure 8

Shear modulus $G^{\left(i\right)}$ of a series of $N=32$ disk packings with repulsive (a) linear and (b) Hertzian spring interactions as the system undergoes a jump change during isotropic compression (at $P\approx 1.23\times {10}^{-4}$ for $\alpha =2$ and $P=4.57\times {10}^{-6}$ for $\alpha =5/2$ indicated by vertical dashed lines) from the isostatic (black upward triangles) geometrical family with $N_c^0$ contacts to a second geometrical family (blue circles) with $N_c^0+1$ contacts.

Figure 9

(a) $(G^{\left(2\right)}/G_0){(P/P_0)}^{-1/2}$ vs $P/P_0$ for packings with $\alpha =2$ in the second geometrical family following a point or jump change in the contact network. $(G^{\left(2\right)}/G_0){(P/P_0)}^{-2/3}$ is plotted vs $P/P_0$ for packings with $\alpha =5/2$ in the second geometrical family following (b) a point change or (c) a jump change in the contact network. In (a)–(c), several system sizes are shown: $N=32$ (black upper triangles), 64 (blue circles), 128 (red diamonds), 256 (green downward triangles), and 512 (magenta asterisks). The dashed lines in (a)–(c) give Eq. (4) for $\alpha =2$ in (a) and $5/2$ in (b) and (c). In the insets to (a)–(c), we show $\mathrm{\Delta }{G}^{\left(1\right)}=G^{\left(1\right)}(P/P_0)-G^{\left(1\right)}\left(0\right)$ or $G^{\left(1\right)}(P/P_0)$ (black upward triangles), $\mathrm{\Delta }{G}_n^{\left(1\right)}=G_n^{\left(1\right)}(P/P_0)-G_n^{\left(1\right)}\left(0\right)$ or $G_n^{\left(1\right)}(P/P_0)$ (red dots), and $\mathrm{\Delta }{G}_a^{\left(1\right)}=G_a^{\left(1\right)}(P/P_0)-G_a^{\left(1\right)}\left(0\right)$ or $G_a^{\left(1\right)}(P/P_0)$ (blue $\times$) with best fits to Eq. (2) (black, red, and blue solid lines, respectively) for an example $N=32$ packing with $\alpha =2$ [inset to (a)] and $5/2$ [insets to (b) and (c)].

Figure 10

Shear modulus $G^{\lambda }$ for initial condition $\lambda$ at $P=0$ undergoing isotropic compression as a function of pressure $P$ for an $N=32$ packing with repulsive linear (black upward triangles) and Hertzian spring interactions (blue circles). $G^{\lambda }=G^{f,\lambda }+G^{r,\lambda }$ can be decomposed into the contributions from the continuous geometrical families $G^{f,\lambda }$ and discontinuities $G^{r,\lambda }$ caused by point and jump changes in the contact network.

Figure 11

Ensemble-averaged shear modulus $\langle G\rangle$ (black upward triangles) as a function of pressure $P$ for $N=128$ packings with (a) $\alpha =2$ and (b) $5/2$ decomposed into contributions from geometrical families $|\langle G^f\rangle |$ (blue circles) and changes in the contact network $\langle G^r\rangle$ (red diamonds). In (a), the dashed line has slope equal to $1/2$, and in (b), the dashed and dotted lines have slopes equal to $1/3$ and $2/3$, respectively.

Figure 12

(a) Ensemble-averaged shear modulus $(\langle G\rangle /G_c){(P/P_c)}^{-(a+b)/2}$ vs $P/P_c$ for jammed disk packings with $\alpha =2$ and system sizes $N=32$ (black upward triangles), 64 (blue circles), 128 (red diamonds), 256 (green downward triangles), and 512 (magenta asterisks). The dashed lines have slopes equal to $-1/4$ and $1/4$. The dotted line gives Eq. (6). (b) $\langle G\rangle /G_c$ plotted vs ${(P/P_c)}^{(a+b)/2}$ for the same data in (a). The solid lines in (a) and (b) are fits to Eq. (7) with $p=2$–5 for system sizes $N=32$–512. The upper left inset shows $P_c$ (black asterisks) and $G_c$ (blue circles) vs $N$. The dotted and dashed lines have slopes equal to $-1$ and $-2$, respectively. The lower right inset gives the exponents, $a$ (black upper triangles) and $b$ (blue diamonds), used in fits to Eq. (7) vs $N$. The horizontal dotted and dashed lines indicate $a=0$ and $b=0.5$, respectively.

Figure 13

(a) Ensemble-averaged shear modulus $(\langle G\rangle /G_c){(P/P_c)}^{-(a+b)/2}$ vs $P/P_c$ for jammed disk packings with $\alpha =5/2$ and system sizes $N=32$ (black upward triangles), 64 (blue circles), 128 (red diamonds), 256 (green downward triangles), and 512 (magenta asterisks). The dashed lines have slopes equal to $-1/6$ and $1/6$. The dotted line gives Eq. (6). (b) $\langle G\rangle /G_c$ plotted vs ${(P/P_c)}^{(a+b)/2}$ for the same data in (a). The solid lines in (a) and (b) are fits to Eq. (7) with $p=2$–15 for system sizes $N=32$–512. The upper left inset shows $P_c$ (black asterisks) and $G_c$ (blue circles) vs $N$. The dotted and dashed lines have slopes equal to $-2$ and $-3$, respectively. The lower right inset gives the exponents, $a$ (black upper triangles) and $b$ (blue diamonds), used in fits to Eq. (7) vs $N$. The horizontal dotted and dashed lines indicate $a=1/3$ and $b=2/3$, respectively.

Figure 14

(a) Shear modulus $G^{\left(1\right)}$ within isostatic geometrical families vs pressure $P$ for individual $N=64$ sphere packings in 3D with $\alpha =2$ repulsive interactions. The solid blue lines are fits to Eq. (2). (b) $(G^{\left(1\right)}/G_0){(P/P_0)}^{-1/2}$ plotted vs $P/P_0$ for the data in (a) and isostatic geometrical families with $G^{\left(1\right)}<0$. The dashed line gives Eq. (4) for $\alpha =2$, and the solid black lines are fits to Eq. (5).

Figure 15

Shear modulus $G^{\left(1\right)}$ within isostatic geometrical families vs pressure $P$ for individual $N=64$ sphere packings with $\alpha =5/2$ repulsive interactions. The solid blue lines are fits to Eq. (2). The dashed line in (a) has a slope equal to $1/3$. (b) $(G^{\left(1\right)}/G_0){(P/P_0)}^{-2/3}$ plotted vs $P/P_0$ for the data in (a) and isostatic geometrical families with $G^{\left(1\right)}<0$. The dashed line gives Eq. (4) for $\alpha =5/2$, and the solid black lines are fits to Eq. (5).

Figure 16

Shear modulus $G^{\left(1\right)}$ within isostatic geometrical families vs pressure $P$ for individual $N=32$ jammed disk packings with repulsive interactions with power-law exponent $\alpha =3$ in Eq. (1). The dashed line has a slope equal to $1/2$. The solid blue lines are fits to Eq. (2) for each of the 50 geometrical families. (b) $(G^{\left(1\right)}/G_0){(P/P_0)}^{-3/4}$ plotted vs $P/P_0$ for the data in (a) and isostatic geometrical families with $G^{\left(1\right)}<0$. The dashed line gives Eq. (4) for $\alpha =3$, and the solid black lines are fits to Eq. (5).

Figure 17

The probability distribution $p\left(P_{-}\right)$ of the pressure $P_{-}$ at which the shear modulus for the isostatic geometrical family first becomes negative $G^{\left(1\right)}<0$ for disk packings with (a) $\alpha =2$ and (b) $5/2$ and system sizes $N=32$ (black upward triangles), 64 (blue circles), 128 (red diamonds), 256 (green downward triangles), and 512 (magenta asterisks).

Figure 18

Fraction of packings $F\left(P\right)$ at each pressure that possess a negative shear modulus for disk packings with (a) $\alpha =2$ and (b) $5/2$ and system sizes $N=32$ (black upward triangles), 64 (blue circles), 128 (red diamonds), 256 (green downward triangles), and 512 (magenta asterisks).