Ergodic-nonergodic transition in tapped granular systems: The role of persistent contacts

 

Paula A. Gago,^1-3 Diego Maza,⁴ Luis A. Pugnaloni¹'²*

* E-mail: luis.pugnaloni@frlp.utn.edu.ar

¹ Dpto. Ingeniería Mecánica, Facultad Regional La Plata, Universidad Tecnológica Nacional, Av. 60 Esq. 124, 1900 La Plata, Argentina.
² Consejo Nacional de Investigaciones Científicas y Táecnicas, Argentina.
³ Current address: Department of Earth Science and Engineering, Imperial College London, South Kensington Campus, London SW7 2AZ, UK.
⁴ Departamento de Física y Matematica Aplicada, Facultad de Ciencias, Universidad de Navarra, Navarra, Spain.

Received: 2 November 2015,  Accepted: 22 December 2015;
DOI: http://dx.doi.org/10.4279/PIP.080001

---

Abstract

Static granular packs have been studied in the last three decades in the frame of a modified equilibrium statistical mechanics that assumes ergodicity as a basic postulate. The canonical example on which this framework is tested consists in the series of static configurations visited by a granular column subjected to taps. By analyzing the response of a realistic model of grains, we demonstrate that volume and stress variables visit different regions of the phase space at low tap intensities in different realizations of the experiment. We show that the tap intensity beyond which sampling by tapping becomes ergodic coincides with the forcing necessary to break all particle-particle contacts during each tap. These results imply that the well-known "reversible" branch of tapped granular columns is only valid at relatively high tap intensities.

Keywords Ergodicity breaking; Tapped Granular Packs

---

 

I. Introduction

Granular matter is ubiquitous in nature. However, due to the complexity of the real particle-particle interactions, the standard approaches of continuum mechanics and thermodynamics are still limited in providing meaningful descriptions of the states in which these systems can be. Edwards and Oakeshot introduced a tentative approach inspired by the ideas of equilibrium statistical mechanics to formally describe the global properties of a static granular pack. Since the introduction of this theory —where the entropy of the systems is governed by the spatial disorder of the grains [1], a number of studies have used it to frame the interpretation of the results of specific experiments. The most relevant case is the so-called "Chicago experiment", where a column of grains was repeatedly tapped following an annealing-type protocol [2,3]. The main outcome of this experiment is that a stationary state can be reached, where the mean volume fraction,Φ, is a well defined function of the tap amplitude, r. Others have also obtained seemingly reproducible states without the need of annealing [4]. However, it has been shown recently, by simulation of frictionless grains, that these stationary states are not necessarily ergodic [5]. At low Г, different members of an ensemble of steady-states prepared with a well defined protocol may sample a different region of the phase space, as the fluctuations of Φ indicate.

In this paper, we demonstrate that not only the volume but also the force moment tensor, ∑, are sampled in a non-ergodic fashion and that ergodicity seems to be recovered if all particle-particle contacts are lost during each tap. This sets a clear limit to the range of driving forces able to generate a sequence of configurations for which the Edwards framework can be applied.

II. Numerical protocol

We simulated using the LAMMPS package [6] a quasi-two-dimensional cell containing N = 1000 spherical particles of diameter d. The cell is 1.1 d thick and 27.8 d wide (the granular column is about 35 layers deep) to have a one to one representation of a previously introduced experimental device [7, 8]. We use a model for soft frictional spheres described in Refs. [9, 10]. The normal component, F_n, of the contact interaction is given by an elastic repulsive force proportional to the overlap of the interacting spheres and a dissipative term proportional to the normal component of the relative velocity. The tangential term, F_t, implements an elastic shear force and a damping force. The shear force takes into account the cumulated tangential displacement between the particles while they remain in contact. Whenever F_t >F_n is the friction coefficient), this lower dynamic friction force is used. In this work we use the same interaction parameters as in Ref. [11, 12]. The wall-particle interaction is defined with the same parameters as the particle-particle force. Tapping is simulated by imposing an external vertical motion to the cell. This pulse is a single sinusoidal cycle A sin(wt). We fix x 33 Hz and use the tap amplitude A as control parameter. The tap intensity is characterized by The mechanical equilibrium after each tap is deemed achieved if the kinetic energy of each particle has fallen (in average) below 10^-6mgd. Where m is the mass of one particle and g the acceleration of gravity.

We study 20 independent realizations of a decreasing ramp of the tap amplitude. We initially fill the cell by placing the spheres at random positions before letting them deposit under the action of gravity. In each realization, we decrease r in small steps, from r = 20.0 down to Г = 0.8, and apply 200 taps for each Г. Note that for Г < 1.0, the column of grains does not detach from the base during a tap. The 200 taps at each value of r are enough to reach a steady-state. We do not observe any drift of the mean values of Φ or ∑ after the initial 100 taps, which we will discard later in our analysis. Finally, we also study a cyclic annealing protocol: starting from the final configuration at Г = 0.8 for each of the former realizations; the tap amplitude is cyclically increased and decreased every 200 taps in the range 0.8 <Г< 5.5 in order to compare the steady-states reached with an alternative method.

III. Data analysis

To measure the packing fraction we use the 2D Voronoi tessellation (implemented in [13]) of the x-z plane projection of the particle positions, disregarding the third coordinate on the thin direction of the cell. Then, we associate to each particle a "local volume" fraction by dividing the particle area by the corresponding Voronoi area. In order to avoid boundary effects, we disregard particles closer than 2d to the lateral walls. Following the recommendations in Ref. [14], we analyze horizontal slices of the granular column 15 d thick measured at approximately the same depth with respect to the free surface in order to retrieve unbiased results for the force moment tensor due to the uneven free surface. Averaging over the N particles contained in the slice of interest, we obtain the volume fraction of each static configuration. To obtain the steady-state Φ corresponding to a given r, we averaged this quantity over the last 100 configurations obtained for each tap intensity. We also obtain the force moment tensor of each particle in the slice of interest. Here, the sum runs over all the contacts c of the particle,is the vector from the center of the grain to the contact c and is the corresponding contact force. We apply the same averaging protocol used for ∑ to obtain the force moment tensor for the configuration and ∑ for the steady-state of a given realization and Г.

IV. Results

In Fig. 1(a) the ensemble average of the steady-state Φ (i.e., averaged over the 20 independent realizations) is displayed as a function of Г. The error bars correspond to the standard deviation over the 20 realizations. As observed by a number of authors, the curve seems to be very well defined with independent realizations falling within a very narrow range of Φ values for any given r. In the past, this led to the conclusion that this was a truly reversible process, where lowering or raising r would lead to the same steady-state Φ . In the inset of Fig. 1(a) two of the independent realizations are shown for the low-tap-intensity region. From this picture, it is clear that steady-states corresponding to a given r can differ from one realization to another. Notice that in the inset the error bars for the two isolated realizations correspond to the standard error of the mean (SEM), which gives an estimate of the uncertainty of the mean value reported rather than the size of the Φ fluctuations. For these two realizations, although the mean Φ seems to agree within the estimated error for intensities Г > 1.5, it is clear that they are different for low Г. This is consistent with the findings of Pail-lusson and Frenkel [5] for frictionless spheres under event-driven simulations. However, in our simulations we are able to extract the stress state of the system as well as the history of the contacts. These reveal valuable information, as we discuss below.

In Fig. 1(b), we show the trace Tr(∑) of ∑ averaged over all 20 realizations as a function of Г. As before, the error bars indicate the standard deviation over the 20 realizations. As we can see, the variability of the mean stress is significantly large at low r. This is not due to the large fluctuations during a given series of taps but to the variations observed from one realization of the protocol to the other. Indeed, the inset in Fig. 1(b) shows that, for low r, the mean values of Tr(∑) for two of the 20 realizations have a relatively small SEM (i.e., the fluctuations in each realization are small). However, realizations differ from each other. The difference between results corresponding to different realizations become much more evident here than in the case of the Φ-Г plot. We suggest that the stress tensor may be more sensitive and then more suitable to sense if ergodicity is fulfilled in experimental data. Overall, as r is decreased, different realizations explore non-overlapping ranges of volume/stress. Therefore, temporal averages (on a single time series) do not match with ensemble averages (over the realizations).

Since the number of taps we have explored for each r may be small to assume that the steady-state has been properly sampled, we carried out 2000 additional taps at Г = 0.8 for each realization. Since some of the signals are not normally distributed, we confirmed the stationarity of these states by using a non-parametric test at a level of significance of 5% [15, 16]. Two of the normally distributed realizations are displayed in Fig. 2 as a function of the tap number. Comparing the corresponding notched boxes and "violin" diagrams of both signals, it is clear that the states do not match. Hence, the reversible branch found in Ref. [3] is not so at low r in our case since truly stationary states with distinguishable Φ and Tr(∑) may be obtained on different realizations of the annealing protocol. Of course, this is harder to detect in Φ since the dispersion between realizations is much smaller compared with the range of Φ values obtained at different Г.

The former results also confirm the hypothesis that non-ergodicity is present in a typical tapping protocol beyond the special case reported in Ref. [5], where the steady-states obtained did not follow any annealing-type protocol. Hence, we observe this non-ergodic behavior even after annealing the system from high tapping strengths. To stress this point and assess if the speed of the annealing may prevent the system from reaching a unique steady-state on each realization, we apply a slower cyclic annealing protocol (similar to the one introduced in [2]) to each of the 20 final states at low r in order to reproduce the "reversible branch". In Fig. 3 we display a sequence of two successive up and down ramps applied to one of the 20 initial realizations using T-steps about one half of those used in Fig. 1. Although in the scale used for < the steady-state packing fraction seems reversible, a close inspection shows that the states have distinguishable < at low Г [see Fig. 3(a) and the corresponding inset]. This is much more evident when the stress is analyzed [see Fig. 3(b)].

In order to set a criterion to decide if the steady-states are not ergodic for a given Г, we show in Fig. 4 the p-values for the Kruskal-Wallis [17] one-way analysis of variance performed on the 20 realizations at each Г. This simple non-parametric test allows for the rejection of the null hypothesis that all 20 data series are drawn from a unique distribution (which does not need to be normal), hence that they correspond to a unique steady-state. If the p-value is significant (in our case p > 0.01), then we cannot rule out the possibility that the 20 series come from the same steady-state. As we can see from Fig. 4(a), the test run on the data for Φ indicates that for Г < 5.0, the null hypothesis must be rejected and therefore there exist at least two out of the 20 steady-states that are not the same. However, for higher T, the test is significant and then the 20 realizations may correspond to the same steady-state. Interestingly, when the test is run on ∑ [see Fig. 4(b)], the steady-state seems to be unique for all 20 realizations if Г > 3. 75. Although differences between realizations are simpler to detect on visual inspection of Tr(∑), it is actually that sets a higher threshold for the Г values needed to ensure an ergodic steady-state Г (i.e.,> 5.0).

The previous results indicate that the steady-states sampled at low tap intensities do not only depend on Г but on the particular history of each realization. Notice that this goes beyond the history dependent out-of-equilibrium trajectories already reported in tapped systems [18] since here we are focusing on the steady-states. One may hypothesize that the constraints imposed by the contacts is one of the reasons for the non-ergodic behavior at low tap. If a contact persists from one tap to the next, the contact force after coming back to rest will depend on the history of all contacts of that particular grain. In order to test this idea, we analyze the evolution of all contacts during each tap and identify those that persist (i.e., contacts that did not break at any time during the pulse of energy). Figure 5 shows the average ratio C/C₀ of persistent contacts, C, to the total number of contacts, C₀, as a function of Г. For this calculation, each contact was tracked during the final 5 taps for each Г on 6 of the independent realizations and only grains that fall within the layer of interest, as discussed above, where included in the analysis. The percentage of persistent contacts is very small but non-zero up to T « 5.0. As it is expected, when T is increased sufficiently all the contacts are broken and new ones are made during each tap (resulting in C/C₀ = 0). This transition coincides with the value of Г where the realizations seem to sample the same steady-state (see the p-values included in Fig. 5). Therefore, when small taps are applied, the aging of some of the contacts seem to lead the system to sample different regions of the phase space during independent realizations. However, if all contacts are made anew at each tap, the sampling becomes compatible with the idea of ergodicity introduced in Fig. 4. In order to generalize this result, we also simulated frictionless grains. Interestingly, the same conclusion drawn for frictional grains is true for frictionless ones: different realizations seem to sample the same steady-state only if all contacts are made anew upon each tap [see Fig. 5 (b)].

V. Conclusions

Our analysis of the steady-states of tapped granular systems indicate that these states are history-dependent for tap intensities below a certain threshold. This is in contradiction with the general assumption that macroscopic time averages —such as the volume fraction— can be recovered when the amplitude of the perturbation applied to the system is tuned back and forth. The differences between independent realizations become particularly noticeable in the stress distribution. These findings show that the postulates of the equilibrium statistical thermodynamics may not be always fulfilled to describe the steady state of static granular systems (see also Ref. [19] for a discussion on the Boltzmann distribution failure for an analytically solvable model). Focusing on tap intensities that warrant that all contacts are made anew after each tap may allow exploring the available phase space in agreement with the ergodic hypothesis. However, gentle perturbations deserve an approach that includes memory effects to suitably describe the states. In that sense, non-equilibrium thermodynamic approaches may be a suitable alternative [20]. Further research on such alternative formalisms, the effect of other types of forcing mechanisms (e.g., shear), and possible extensions to other complex systems (e.g., active matter) become necessary.

 

Acknowledgements - This work has been partially supported by projects PICT-2012-2155 ANPCyT (Argentina), FIS2011-26675 MINECO (Spain) and PIUNA (Universidad de Navarra).

 

1 S F Edwards, R B Oakeshot, Theory of powders, Physica 157D, 1091 (1989).         [ Links ]

2 E R Nowak, J B Knight, E Ben-Naim, H M Jaeger, S R Nagel, Density fluctuations in vibrated granular materials, Phys. Rev. E 57, 1971 (1998).         [ Links ]

3 E R Nowak, J B Knight, M L Povinelli, H M Jeager, S R Nagel, Reversibility and irre-versibility in the packing of vibrated granular material, Powder Technol. 94, 79 (1997).         [ Links ]

4 Ph Ribiere, P Richard, P Philippe, D Bideau, R Delannay, On the existence of stationary states during granular compaction, Eur. Phys. J. E 22, 249 (2007).         [ Links ]

5 F Paillusson, D Frenkel, Probing ergodicity in granular matter, Phys. Rev. Lett. 109, 208001 (2012).         [ Links ]

6 S Plimpton, P Crozier, A Thompson, LAMMPS-large-scale atomic/molecular massively parallel simulator, Sandia National Laboratories (2007).

7 L A Pugnaloni, I Snchez, P A Gago, J Damas, I Zuriguel, D Maza, Towards a relevant set of state variables to describe static granular packings, Phys. Rev. E 82, 050301 (2010).         [ Links ]

8 S Ardanza-Trevijano, I Zuriguel, R Arevalo, D Maza, Topological analysis of tapped granular media using persistent homology, Phys. Rev. E 89, 052212 (2014).         [ Links ]

9 N V Brilliantov, F Spahn, J-M Hertzsch, T Poschel, Model for collisions in granular gases, Phys. Rev. E. 53, 5382 (1996).         [ Links ]

10 L E Silbert, D Ertas, G S Grest, T C Halsey, D Levine, S Plimpton, Granular flow down an inclined plane: Bagnold scaling and rheology, Phys. Rev. E 64, 051302 (2001).         [ Links ]

11 L A Pugnaloni, J Damas, I Zuriguel, D Maza, Master curves for the stress tensor invariants in stationary states ofstatic granular beds. Implications for the thermodynamic phase space, Papers in Physics 3, 030004 (2011).         [ Links ]

12 L A Pugnaloni, M Mizrahi, C M Carlevaro, F Vericat, Nonmonotonic reversible branch in four model granular beds subjected to vertical vibration, Phys. Rev. E. 78, 051305 (2008).         [ Links ]

13 C Rycroft, VORO++: A three-dimensional Voronoi cell library in C++, Chaos 19, 041111 (2009).         [ Links ]

14 P A Gago, L A Pugnaloni, D Maza, Relevance of system size to the steady-state properties of tapped granular systems, Phys. Rev. E 91, 032207 (2015).         [ Links ]

15 W Constantine, D Percival, Fractal: Fractal Time Series Modeling and Analysis. R package version 2.0-0, http://CRAN.R-project.org/package=fractal (2014).         [ Links ]

16 R Core Team, R: A language and environment for statistical computing, R Foundation for Statistical Computing, Vienna, Austria (2013).         [ Links ]

17 W H Kruskal, W A Wallis, Use of ranks in one-criterion variance analysis, J. Am. Stat. Assoc. 47, 583 (1952).         [ Links ]

18 Ch Josserand, A V Tkachenko, D M Mueth, H M Jaeger, Memory effects in granular material, Phys. Rev. Lett. 85, 3632 (2000).         [ Links ]

19 R M Irastorza, C M Carlevaro, L A Pugnaloni, Exact predictions from the Edwards ensemble versus realistic simulations of tapped narrow two-dimensional granular columns, J. Stat. Mech. P12012 (2013).         [ Links ]

20 G Lebon, D Jou, JC Casas-Vazquez, Understanding Non-equilibrium Thermodinamics, Springer-Verlag, Berlin Heidelberg (2008).         [ Links ]