Statistical properties of cell stresses in 2D granular solids

We investigate numerically the statistics of the stresses, defined for cells the irreducible loops surrounded by particles in contact in two-dimensional granular solids and the relation between these cell stresses and cell structures. We prepared a number of mechanically stable granular assemblies at various packing fractions and mean coordination numbers, and measured the stress statistics on each cell. We find the following. (1) The mean cell stress agrees with the bulk external pressure; (2) the probability density function of the normalised cell stress ratio collapses into a master curve; (3) the mean stress ratio increases with increasing intergranular friction; (4) the mean stress ratio per cell order decreases with increasing cell order because of their reduced stability; (5) the maximum cell stress direction correlates well with the cell’s long axis. These results are strong evidence that the cell structure and the stress field self-organise together and reveal the co-evolving correlations between the two.


Introduction
Granular materials are ubiquitous in our daily life, but their mechanics are still not well understood. In particular, granular solids support stresses rigidly as regular solids, but their stiffness, shear strength and other mechanical properties depend strongly on both the grain properties and the packing structure. Since Seelig and Wulff [1] and Oda and Konishi [2] showed evidence of force chains in two-dimensional (2D) granular assemblies, much effort went into the development of stress field models, based on this inhomogeneous structure [3][4][5][6]. However, this approach encounters serious difficulties in relating the stress along a chain of particles with the relevant strain increment. A different set of studies has focused in recent years on cell structures, where a 'cell' is defined as an irreducible loop surrounded by particles in contact. Ball and Blumenfeld [7] introduced the concept of loop forces and derived a stress-structure constitutive equation for 2D isostatic systems. Tordesillas et al. [8] showed that the minimal 3-cell structure in 2D granular system stabilizes the load-bearing force chains. Nicot and Darve [9] proposed a constitutive model, focusing on 6-cell structures. Matsushima and Blumenfeld [10,11] studied the structural properties of cells in details and showed that the cell shape is elongated the most when the cell order is around six, because elongated cells of higher cell order are less stable mechanically. Recently, Wanjura et al. [12] developed a model for cell structure evolution, based on the statistics of cell making and breaking during quasi-static dynamics. A significant advantage of relating granular mechanics to cell structural characteristics is that cell stress deformation is expected to be directly linked to the coarse-grained and * Corresponding author: tmatsu@kz.tsukuba.ac.jp upscaled bulk stress and strain. Several fundamental questions remain in this approach: (1) What are the relations between the local cell stress and the bulk stress? (2) How does the local cell structure affect the stress and vice versa? Providing even partial answers to these is very important. Here, we investigate in detail the statistics of cell stresses, defined below, in 2D circular grains assemblies in mechanical equilibrium under isotropic bulk stress field. We first compute the statistical distribution of each component of the cell stress, and then discuss the correlations between the cell stress ratio and its geometrical characteristics.

Cell structure and cell stress
In 2D disc systems, a cell can be defined as the smallest (irreducible) closed loop of disc centres around a void in the granular solid (Fig. 1). With this definition, the total volume of the system is a sum of all the cell volumes, . This definition differs from that in Matsushima and Blumenfeld [10,11], where a cell is defined by the shortest (irreducible) closed loop of contact points around a void. Nevertheless, this difference is of no major consequence in disc systems. The stress of grain , , is defined via contact forces on it: in which ′ is the vector extending from the grain centre to the contact point with the neighboring grain ′ , and ′ is the contact force that ′ , applies to (Fig.  2). The normalising area, , is defined as that enclosed within the centres of surrounding grains, as shown in Fig.  1. Hence, the system's area is the sum of over all grains. The cell stress is then defined by the area-weighted average of the grain stress of the grains around the cell: which is different from . † A little manipulation leads to the following equation showing that the cell stress can be determined from the boundary forces, i.e. the contact forces around the cell.  For our numerical experiments we used the Discrete Element Method (DEM, [13,14]). We analysed the large set of systems used in Matsushima and Blumenfeld [10,11] to study in detail the geometric properties of † Using the concept of loop forces [6], we may define a cell stress corresponding to the area in the case of ideal isostatic state, but it is not straightforward in general case. packing structure. Each system consisted of about 21,000 grains, with an identical slightly-dispersed grain size distribution. Starting from various initial unjammed states of different packing fraction ( =0.72, 0.76, 0.82, denoted as LIS, IIS and DIS, respectively) and intergranular friction coefficient ( =0.01, 0.1, 0.2, 0.5 and 10), an isotropic compressive stress 0 was applied by changing the periodic length in the two directions. The state with 0 was defined as that in which the average overlap of grains was = 0 ⁄ = 10 −5 , where is the normal contact stiffness. The compression continued until the dynamic fluctuations of both grain positions (per mean grain diameter) and intergranular forces (per mean average contact force) were below 10 −9 and 10 −6 , respectively. The obtained mechanically equilibrated systems covered a wide range of combinations of the packing fraction, , and the mean coordination number, ̅, as shown in Fig. 3. Note that ̅ in the figure is computed by discarding rattlers, i,e. grains with one or no force-carrying contact, and its upper and lower bounds =4 and =3 correspond generically to isostatic states of smooth ( =0) and highly frictional ( → ∞ ) grains, respectively.

Statistics of cell stresses
Figs. 4-6 show the probability density of the mean cell stress = ( 1 + 2 )/2, the deviatoric cell stress = ( 1 − 2 )/2 , and their ratio, ℎ ≡ / , respectively, for various systems shown in Fig. 3. In these expressions, 1 and 2 are the maximum and minimum principal stresses, respectively. is widely distributed, ranging from almost 0 to twice the bulk stress 0 , and while the standard deviation increases slightly with the intergranular friction , the mean value ����������� is almost identical to the bulk pressure 0 for all the systems, within an accuracy of 3.0-5.5%. The distribution of appears to have a slightly longer tail for larger , and accordingly ����������� is also larger for larger . The distribution of the ratio ℎ also has a slight but clear dependence on , irrespective of the initial condition. However, and significantly, on normalising ℎ by its mean value ℎ ̅ , all the distributions collapse onto a unique curve, as shown in Fig. 7. We established that this distribution could not be fitted well by a gamma distribution. To pursue the origin of the PDF of ℎ, we considered the decomposition where ( ) is the occurrence probability of cells of order , and (ℎ| ) is the conditional PDF given that it corresponds to a cell of order . Since ( ) has already been reported in Matsushima and Blumenfeld [11], we focused here on (ℎ| ) . The mean value of the conditional PDF as a function of , ℎ ̅ ( ) = ∑ ℎ (ℎ| ) ∈ , is shown in Fig. 8. Interestingly, ℎ ̅ ( ) decreases roughly linearly with at a rate that is almost identical in all the systems. Moreover, for a given value of , ℎ ̅ ( ) increases with the intergranular friction . No effect of the initial condition was observed, which is consistent with the findings of Matsushima and Blumenfeld [10], that the difference of the packing fraction shown in Fig. 3 for the same friction disappears if rattlers are removed. We also observed that normalising the conditional probability of ℎ by the mean, ℎ ̅ ( ) , (ℎ/ℎ ̅ ( )| ) , collapses again all the distributions onto a single curve, regardless of the initial condition, intergranular friction and cell order . This result is striking and encourages pursuing for understanding the relations between micromechanics of granular solids and cell structures.

Co-axiality between cell stress and cell geometry
Next, we explore the correlation between cell stresses and cell geometries. To study the statistical properties of cell shape, Matsushima and Blumenfeld [11] approximated cells to lowest order as best-fitted ellipses and showed that there is a relation between cell elongation and cell order. They argued that this relation can be understood on mechanical stability grounds. Here we find that the cell long axis direction is correlated well with the direction of the maximum principal stress, as we show in Fig. 10. This co-axiality can also be explained by mechanical stability: an elliptical cell structure is more stable against loading from its long axis than its short axis. Taken together with the observation in Fig. 8, this is evidence that: (1) stability decreases with cell order e and the higher the cell order the less elongated the cell shapes; (2) high-order cells support lower cell stresses than low-order ones; (3) the more elongated the cell the larger the cell stress ratio it can support; (4) the distribution of the cell stress ratio of a given cell order is determined by the distribution of cell elongation and granular friction. Fig. 11 displays an example of this interplay between cell geometry and cell stress. In the figure, cells are drawn in black and the direction of the maximum principal stress and the magnitude of cell stress ratio ℎ are indicated by the orientation and the length of red segments, respectively. It demonstrates clearly points (1)-(3) above. Our results suggest that cell-based mechanics of granular solids is a promising direction to explore.

Conclusions
We defined cell stresses and studied numerically, using DEM, their statistics and their relation to cell structures in 2D granular solids. We have found: (i) that the conditional PDFs of the cell stress ratios collapse onto a master curve and (ii) a co-axiality between the large principal cell stress and the orientation of the cell's major axis. These results not only suggest that it would be useful to analyse cell shape and geometry statistics together with cell stresses but they also show that it is essential to understand the structural organization on the cell-scale for to modelling the bulk behaviour of granular solids. A more detailed discussion will be presented elsewhere.
This work was supported by JSPS KAKENHI Grant Number 18KK0115 and 20K20434.