# Dynamical response functions in models of vibrated granular media

## Abstract

In recently introduced schematic lattice gas models for vibrated dry granular media, we study the dynamical response of the system to small perturbations of shaking amplitudes and its relations with the characteristic fluctuations. Strong off equilibrium features appear and a generalized version of the fluctuation dissipation theorem is introduced. The relations with thermal glassy systems and the role of Edwards’ compactivity are discussed.

Few years ago, Edwards formulated the hypothesis that it is possible to extend to powders the methods of standard statistical mechanics , an important challenge since powders are by definition “non thermal” systems . Actually, in granular media the role of a “temperature”, linked to the concept of Edwards’ compactivity , is played by the amplitude of external vibrations , and it is possible to properly individuate the corresponding “equilibrium” states . However, these systems are typically in off equilibrium configurations, as the presence of aging phenomena shows , and the above extension of statistical mechanics thus runs across its off equilibrium version.

A very important issue for both theoretical and practical reasons, is the understanding of relations between the response of a system to external perturbations and its characteristic fluctuations. Here, we study such a problem in the context of schematic lattice gas models  recently introduced to describe gently shaken granular systems. We show how in granular media it is possible to formulate a fluctuation dissipation theorem (FDT) where the shaking amplitudes play the role of usual temperatures, but in typical off equilibrium situations it coincides neither with its usual version at equilibrium nor with the extensions valid in off equilibrium thermal systems in the so called “small entropy production” limit . We also discuss the important relations with Edwards’ theory.

The lattice gas model we consider here, the Tetris, schematically describes the effects of steric hindrance and geometric frustration on grains in granular materials. Interestingly it shows many phenomena typical of granular media under shaking, as logarithmic compaction or segregation . It consists of a system of elongated grains, with two possible orientations, which occupy the sites of a square lattice tilted by $`45^{}`$. To avoid overlaps, two particles can be nearest neighbor only if they have the right reciprocal orientation. The motion of grains in absence of “vibrations” is subject only to gravity so grains only move downwards (without overlapping). The effect of vibration is introduced by allowing the particles to diffuse upwards with a probability $`p_{up}`$ and downwards with $`p_{down}=1p_{up}`$. The effects of geometric frustration on the microscopic dynamics are introduced with a kinetic constraint: particles can flip their orientation only if at least three of their neighbors are empty.

The parameter governing the dynamics is the adimensional quantity $`x_0=p_{up}/p_{down}`$, or $`\mathrm{\Gamma }1/\mathrm{ln}(x_0^{1/2})`$. $`\mathrm{\Gamma }`$ plays the same role as the amplitude of the vibrations in real granular matter , and, for not too high amplitudes, a good agreement is found between the model and experiments by posing $`\mathrm{\Gamma }^ba/g`$ with $`b`$ about 1 (here $`g`$ is the gravity and $`a`$ the peak acceleration of the shakes in the experiments) . In analogy with experiments, in our Monte Carlo simulations the shaking of the system is as follows: we prepare the system, confined in a box, in a given initial configuration at $`t=0`$ (see below), then we start to “shake” it continuously and indefinitely with a given “amplitude” $`x_0`$. We expect very similar results by considering, instead of a single long tap, a series of short taps, as experimentally more convenient (see Ref. ).

The Tetris can be described with the following Ising Hamiltonian in the limit $`J\mathrm{}`$ (see ): $`\beta H=J_{ij}f_{ij}(S_i,S_j)n_in_j\frac{g}{T}_iy_i`$ . Here $`g`$ is the gravity constant, $`y_i`$ is the $`i`$-th particle height, $`n_i=0,1`$ are occupancy variables, $`S_i=\pm 1`$ are Ising spin variables for the two possible orientations of grains. $`f_{ij}`$ is a function describing the hard core repulsion ($`J=\mathrm{}`$): $`f_{ij}(S_i,S_j)=0`$ if the orientation of neighbors, $`(S_i,S_j)`$, is allowed and $`f_{ij}=1`$ if it is not allowed ( $`f_{ij}(S_i,S_j)=1/2\left[S_iS_jϵ_{ij}(S_i+S_j)+1\right]`$, where $`ϵ_{ij}=+1`$ for bonds along one direction of the lattice and $`ϵ_{ij}=1`$ for bonds along the other). The temperature, $`T`$, of the above Hamiltonian system is related to the ratio $`x_0=p_{up}/p_{down}`$ via $`e^{2g/T}=x_0`$ (i.e., $`\mathrm{\Gamma }=T/g`$).

This Hamiltonian mapping shows that in our model the field coupled to grains height is the adimensional gravity: $`g/T=\mathrm{ln}(x_0^{1/2})`$. Thus a perturbation to the system, coupled to a easy observable, may be introduced by varying $`x_0`$. We record, with Monte Carlo simulations, the dynamical correlation functions and the response of the system to such a small perturbation in the “shaking amplitude”.

Inside a box of fixed size $`30\times 60`$, with periodic boundary conditions along the x-axis and rigid walls at bottom and top, we prepare the system in a uniform density initial configuration ($`\rho =0.5`$), corresponding to an highly fluidized state, and then we “shake” it at a given amplitude $`x_0`$. During this process we record the average height $`h_0(t)=_iy_i(t)`$ of the grains and the two times correlation function $`C_0(t,t^{})h_0(t)h_0(t^{})`$, or, the “mean square displacement” $`B(t,t^{})C_0(t,t)2C_0(t,t^{})+C_0(t^{},t^{})`$ which is a quantity relevant to test the FDT theorem. Time is measured in such a way that unity corresponds to an average update of all the degrees of freedom in the system, and the statistical averages run over 2048 noise and initial configuration realizations. In order to measure the response function, we also record the average height of grains, $`h_1(t,t_w)`$, in an identical copy of the system (a “replica”), which, after a fixed time $`t_w`$ (typically below we fix $`t_w=370,3700`$), is shaken at $`x_1=x_0+\mathrm{\Delta }x_0`$, i.e., is perturbed by a small increase of the shaking amplitude $`\mathrm{\Delta }x_0`$ (in what follows $`\mathrm{\Delta }x_0=0.002`$ ).

The quantity $`\mathrm{\Delta }h(t,t_w)=h_1(t,t_w)h_0(t)`$, i.e., the difference in heights between the perturbed and unperturbed systems, is by definition the integrated response. FDT or its generalizations concern the relation between $`\mathrm{\Delta }h(t,t_w)`$ and the displacement, $`B(t,t_w)`$, which is linked to the correlations in the unperturbed system. The simplest version of a possible generalization to off equilibrium granular matter of the FDT may be argued from thermal systems . The integrated response, $`\mathrm{\Delta }h`$, should be approximately proportional to $`B`$:

$$\mathrm{\Delta }h(t,t_w)\frac{X}{2}\mathrm{\Delta }(\mathrm{\Gamma }^1)B(t,t_w)$$
(1)   
Here, $`\mathrm{\Delta }(\mathrm{\Gamma }^1)\mathrm{\Gamma }_0^1\mathrm{\Gamma }_1^1`$, is the variation between the inverse shaking amplitudes of the reference and the perturbed systems and the prefactor, $`X`$, is a quantity to be determined, in principle function of $`t_w`$ and $`t`$ themselves. In thermal system, in the limit $`t,t_w\mathrm{}`$, $`X`$ is a piecewise constant, depending just on $`B(t,t_w)`$ and not on both $`t`$ and $`t_w`$ .

In the specific case of our model we have $`\mathrm{\Delta }(\mathrm{\Gamma }^1)\mathrm{\Delta }(\frac{g}{T})=\mathrm{ln}\left[(\frac{x_0+\mathrm{\Delta }x_0}{x_0})^{\frac{1}{2}}\right]`$, and if $`X=1`$ we recover the usual well known equilibrium version of FDT. In the study of glassy thermal systems the quantity $`\mathrm{\Gamma }/gXT/X`$ has the meaning of an “effective temperature” of the sample, which only for $`X=1`$ coincide with the equilibrium bath temperature .

In the present model, eq. (1) seems to be approximately valid, also if in typical off-equilibrium situations $`X`$ slowly depends on $`t`$ and $`t_w`$. To outline this, as first we explore the high shaking amplitude regime. In Fig.1 (top left), as a function of $`tt_w`$ are shown, for $`t_w=370`$ (analogous results are found for $`t_w=3700`$), the height $`h_0(t)`$ of a system “shaken” at $`x_0=0.8`$, and $`h_1(t,t_w)`$ of the replica at $`x_1=x_0+\mathrm{\Delta }x_0=0.802`$. $`h_0`$ and $`h_1`$ decrease in time signaling that the system is slowly compactifying while approaching the equilibrium density profile. In Fig.1 (top right), is also shown the “displacement”, $`B(t,t_w)`$, of the same unperturbed system. We show in Fig.2 (top left) the difference $`\mathrm{\Delta }h(t,t_w)h_1(t,t_w)h_0(t)`$ as a function of $`tt_w`$. In agreement with the expected simple behavior of a gas subject to gravity (which reacts to an increase of temperature, or decrease of gravity, by increasing its average height), $`\mathrm{\Delta }h`$ is a positive quantity: the “colder” system (with $`x_0=0.8`$) actually must have a lower equilibrium height respect the “hotter” one (with $`x_1=0.8+\mathrm{\Delta }x_0`$). In order to check the generalization of the FDT proposed in eq. (1), in Fig.2 (top right), we plot $`\mathrm{\Delta }h`$ as a function of the displacement $`B`$. Actually, $`\mathrm{\Delta }h`$ is approximately piecewise linear in $`B`$: after an early transient with $`X>1`$, in the long time regime $`X`$ is independent on $`t`$ and $`t_w`$ and equal to 1, showing that, in the high “temperature” and low density region, the usual equilibrium version of FDT is obeyed.

In the low $`x_0`$ region, we already know the presence of strong off equilibrium phenomena . Thus it is reasonable to expect that the above simple picture drastically changes. In Fig.2 (middle left), we plot the integrated response, $`\mathrm{\Delta }h(t,t_w)`$, as a function of $`tt_w`$ for a system “shaken” at $`x_0=0.5`$ and a replica shaken at $`x_0+\mathrm{\Delta }x_0=.502`$ after a $`t_w=370`$ (analogous results are found for $`t_w=3700`$). It is apparent that after an early transient, up to the time scales we explored we find a negative response, $`\mathrm{\Delta }h`$, in complete contrast with the above equilibrium scenario. Actually, both replicas start from a strongly off equilibrium state, but while the “colder” is “frozen” and remains longly trapped in metastable states, the “hotter” is able to more rapidly escape to approach its asymptote.

Notice that $`\mathrm{\Delta }h`$ is negative even after five decades in time and no change in its trend is observed. The approach to equilibrium is actually extremely slow : if the region with negative response would extend up to infinity (which is hard to say with computer simulations), this should results in a dynamical breaking of ergodicity introduced by gravity in our system of particles interacting just with excluded volume effects. Interestingly these results are partially supported by some recent experiments , which showed that the asymptotic density of a granular system compactified at a low shaking amplitude from a random initial configuration, is lower than the analog quantity in a system shaken at a slightly higher amplitude. This result is in contrast with the “equilibrium” measures (the “reversible branch”) in the experiments of Novak et al. . The above phenomenon has strong repercussions on the FDT. As plotted in Fig.2 (middle right), $`\mathrm{\Delta }h`$ may be approximately plotted as a piecewise linear function in $`B`$, as from eq. (1), but after an initial transient with $`X>1`$, the system enters a region with negative $`X`$ ($`X5`$), corresponding to a negative “effective temperature” $`T/X`$. The specific value of $`X`$, in this region, also slowly depends on $`t_w`$ (see below), and we are far from the small entropy production limit known in thermal systems where $`X=X(B)`$.

A similarly anomalous picture is found in a system shaken at $`x_0=0.001`$ and a replica at $`x_0+\mathrm{\Delta }x_0=.003`$. To show that the above general results do not depend on the details of the initial state, we discus the latter Monte Carlo shaking experiment in a differently prepared system, closer to experiments. Now the starting particle configuration is prepared by randomly inserting particles into the box from its top and then letting them fall down, with the above dynamics, until the box is half filled. Thus, the system starts from a more compact static state. In Fig.1, the heights of the original system, $`h_0`$, of the replica $`h_1`$ (bottom left), along with the displacement $`B(t,t_w)`$ (bottom right) are plotted as a function of $`tt_w`$ for $`t_w=333`$. This figure shows again a negative response $`\mathrm{\Delta }h`$ which is plotted as a function of $`tt_w`$ in Fig.2 (bottom left), over six order of magnitudes, for three different values of $`t_w`$ ($`t_w=33,333,3330`$). Three different regions are observed for $`\mathrm{\Delta }h`$, but the short times transient inflection zone seems to disappear for $`t_w`$ long enough. In agreement with the generalized FDT of eq. (1), $`\mathrm{\Delta }h`$ is at long times ($`10^4tt_w10^6`$) again approximately linear in $`B`$, as shown for $`t_w=333`$ in Fig.2 (bottom right) . The negative proportionality coefficient, $`X`$, slowly depends on $`t_w`$ (from $`X7`$ for $`t_w=33`$ to $`X15`$ for $`t_w=3330`$). For systems prepared in the previous fluidized initial state, discussed above, analogous results are found. It is important to stress that the above behaviors with negative responses are not found if the system starts from an equilibrium initial configuration. 

Interestingly, in the high density or low $`x_0`$ region, we observe an analogous scenario also in a similar lattice model, the IFLG , where, in order to describe the motion of grains in a disordered environment, quenched disorder is introduced in the Hamiltonian described above, but no kinetic constraints are present in the dynamics.

As stated, the presence of negative $`X`$ in the generalized FDT, is a feature of the very far from equilibrium dynamics. Actually, this is not expected for instance in off equilibrium thermal systems as glasses or spin glasses, at least in the small entropy production limit . Moreover, in spin glasses the parameter $`X`$, which is a seemingly dynamical quantity, may be asymptotically related to equilibrium static properties from replica theory . In granular media, as long as in fragile glasses , the static correspondent of $`X`$, if any, is still missing. It is useful to underline that while in glasses an experimental measure of the generalized FDT may be non trivial, in granular media, as shown above, this should be reasonably simple.

The interesting finding in our models for granular matter of a negative $`X`$, corresponding to negative responses and “effective temperatures”, may be clarified by stressing the relations with Edwards’ theory .

Actually, the present approach to granular media based on a standard Hamiltonian formalism, is very close in spirit to the more general statistical mechanics approach proposed by Edwards . The hard core repulsion term in the above Hamiltonian $`H`$ has a role similar to Edwards’ volume function $`W`$ (where a constraint of mechanical stability is explicitly present). The fundamental control parameter in Edwards’ theory which corresponds to our adimensional temperature, $`T/g`$, is the compactivity $`\lambda X_E`$. Edwards and Grinev  guess that the compactivity is related to the experimental “temperature” of a shaken granular medium according to the following “fluctuation-dissipation” relation: $`\lambda X_E(a/g)^2`$ (as above, $`a`$ is the shake and $`g`$ is the gravity acceleration). This statement is the analog of our suggestion, $`(T/g)^ba/g`$, which is based on the comparison of Monte Carlo and experimental data.

In Ref. is proposed that the state of a shaken granular medium in the equilibrium regime, corresponding for instance to configurations on the so called “reversible branch” in the experiments of Ref. , is characterized by a positive compactivity, and, by extrapolation, the off equilibrium dynamics, corresponding to the experimental “irreversible branch” of Ref. , has a negative compactivity, $`X_E<0`$. This observation is in agreement with our discover of different regions with positive as long as negative effective temperatures, $`T/X`$, in the study of the generalized FDT in granular media. A negative off equilibrium compactivity may correspond to the necessity of the system to cross states with “higher entropy” in order to lower its volume under shaking. Actually, the measure of integrated responses v.s. correlation functions in granular materials may open the way to a direct experimental access to Edwards’ compactivity for a clear settlement of a statistical mechanics for such systems.

In summary, in order to understand the off equilibrium statistical mechanics of powders, we have studied, in schematic lattice gas models for vibrated dry granular media, the dynamical response functions to small shaking amplitude perturbations and their relations to characteristic dynamical fluctuations. Strong off equilibrium features appeared, as long time regions with negative response functions, different from those observed in the small entropy production limit of thermal glassy systems , along with the necessity to introduce a generalized version of the fluctuation dissipation theorem. The novel properties we have found, as the presence of negative effective temperatures, are consistent with Edwards’ theory of powders and demand important experimental check.

We thank INFM-CINECA for computer time on Cray-T3D/E. Work partially supported by the TMR Network Contract ERBFMRXCT980183 and MURST-PRIN 97.