i
PHYSlCA ELSEVIER
Physica A 230 (I 996) 329-335
A sandpile model with dip Jan H e m m i n g s s o n 1 P.M.M.H., E.S.PC.L, 10 rue Vauquelin, 75231 Paris Cedex 05, France
Received 8 December 1995; revised 20 December 1995
Abstract
We present a simple model of forces inside a granular medium. The model is applied to simulate the forces within a sandpile, and has further the feature of generating the sandpile itself. At the bottom of the sandpile the distribution of vertical forces shows a local minimum under the apex of the pile, in coherence with experimental observations. PACS: 03.20
I. Introduction
Among the many-facetted problems of granular media, we will focus here on one special aspect of the forces inside a medium. The observations by Smid and Novosad [ 1] of the dip in the distribution of vertical forces under the apex of a pile of granular material have yet not been completely understood. Several different models have been developed recently; those that arrive at a uniform force field [2-4] and a distribution that becomes constant within a certain distance of the center of the pile [5]. Edwards and Mounfield proposed a model [6] along the lines of an earlier work [7]. The arching inside the sandpile is modeled by plates, resulting in a double tent-shaped force distribution that takes the value zero under the apex of the pile, which is not the case in the experiments. An analytic work by Bouchaud et al. [8] does produce a dip after a reinterpretation of A in their model [9]. Another, also very promising model by Liu et al. and Coppersmith et al. [10] link together experiments, numerical simulations and exact results. In the following, we will show how a simple model of particles might reproduce the effect of a finite depression under the apex of the pile. 1 Permanent address: IFM, Link6ping University, S-58183 Link0ping, Sweden. 0378-4371/96/$15.00 Copyright (~ 1996 Elsevier Science B.V. All rights reserved Pll S 0 3 7 8 - 4 3 7 1 ( 9 6 ) 0 0 0 2 6 - X
330
J. Hemmingsson/Physica A 230 (1996) 329-335
2. A model of forces in granular media We have developed a simple two-dimensional model for simulating the forces in granular media. We will first explain the model in general, and then we will describe how it may be used to simulate a sandpile on a free surface. Before we go any further, we should state clearly that the proposed model is an extreme simplification of a real system. The approach is phenomenological and not based on some profound theory. In any case, since we are interested in statistical properties of relatively large systems, we hope that the details of the model do not change the main characteristics in a dramatic way. On a triangular lattice, we let each lattice point represent a grain of some granular medium. The idea is to 'relax' the forces within one layer before moving on to the next one. Since forces not are allowed to go upward in the system, every layer only has to be relaxed once. This makes this model simple. Starting by the uppermost layer, each grain experiences a downward force due to gravity. With some probability p, the grain looses contact with one of its two downward neighbors due to inhomogeneities in the medium. This will be referred to as a bridge site. When a site is marked as a bridge site, the grain on that site will with equal probability loose the contact to the grain to the down left or to the down right (but only with one of them). As a result, the force acting on the particle will be decomposed into the horizontal direction and the direction of the particle with which it is still in contact. The horizontal forces are transmitted to the nearest neighbors in the same layer. When a horizontal force arrives from a neighbor, it is added to the sum of forces acting on that grain. If this resulting force is pointing between the two sites underneath, there is no horizontal component left, like for site A in Fig. 1. Otherwise, it is decomposed in the same manner as just described, and the horizontal part is transferred to the next site. When the wall is reached, the arriving force is simply absorbed. This scheme is followed once in each horizontal lattice direction. Thereafter, all force resultants are pointing in a direction between bonds to the two supporting sites, or in the extreme case, in the direction of one of the bonds. Now the forces are decomposed in these two directions and transferred to the next layer where they are added. For the next layer, a force corresponding to gravity is added, and the procedure is repeated. The result does depend on the order in which the sites are updated. When the disorder parameter p is small, however, the bridge sites can be considered as independent, and the updating order is not very important. For p -- 1, i.e. when all sites are bridge sites, the model would be identical to the Scheidegger model if the horizontal forces were discarded. This is a thoroughly studied model for river flows and aggregation [11,12]. If a force is very large, it might deform the grains and reestablish a bond that was previously broken. We introduce this mechanism in a simple way. If the resulting force at a certain site exceeds a value Fc, p is zero for that site. This Fc corresponds to some ratio between the hardness of the particles and the gravity.
331
J. Hemmingsson/Physica A 230 (1996) 329-335
{
D
/,,
E
,,
F
/,,
G
/
Fig. 1. In this example, site B has no contact with its neighbor E Therefore, the horizontal force component is transferred to its right neighbor C. The resultant of site C will thereafter point in the direction of its bond to site G. Site F does not receive force from neither B nor C, and an arch is thus created.
3. Model of a pile When simulating a sandpile using the described model, one might proceed as follows: Start with one grain, which will be the apex of the pile. In the next layer, one grain has to be added in order to support the previous layer. So the second layer will contain two grains, the third contains three, and so on. This works well as long as p = 0, that is when the system is completely ordered. As soon as p is larger than zero, however, forces traveling horizontally in the system will appear and reach either end of the pile. (For simplicity, p -- 0 for the border sites.) Here there are no wails to absorb the forces, so one has to proceed in some other way. Since the system that is simulated is in equilibrium, there must be another force that balances such a force. Consider the rightmost site y of some layer x. If there are horizontal forces acting on this grain, it means that there must be some grain(s) to the right of it, since otherwise it would start moving. One grain is thus added to the position right of y. We also assume that there are grains in the next layer x + 1 to support this grain. The horizontal force is added to the gravity force of this new site, and if this resulting force points in the direction between the bonds to the grains in the next layer, the grain is in equilibrium. If not, the force is decomposed (as described above), and another grain to the right of it has to be added. This is repeated until the horizontal force has vanished. In fact, all sites that are added will be supported by only the outmost site of the two sites in the next layer that could support it, except for the last one. Therefore, the value of p will not be important for these added sites. The system is thus produced during the simulation by adding sites in order to remain in a steady state. In the ensemble of all possible systems in equilibrium with certain parameter choices, we find one. As the algorithm proceeds, the actual shape of the pile is found (see, e.g., Fig. 2). The procedure of adding a grain is just one way of choosing a representative of all the possible equilibrium configurations. The sandpile found here will be represented by an angle which we will refer to as the angle of stability. This angle will vary with both p and -Pc, as one would expect.
J. Hemmingsson/Physica A 230 (1996) 329-335
332 2000
,
i
A
i
1800
1600
1400
1200
1000
800
600
400
200 0
7000
I
I
I
I
I
8000
9000
10000
11000
12000
13000
Fig. 2. The physicalshapeaveragedover50 sandpiles,p = 0.01, Fc = 850 (in units wherethe massof one grain times the gravitationalconstantg is 1).
4. Results We simulated the described model for different values of p, and made histograms of the vertical forces as a function of position. The granular material resides in a box of some width (in this case 20000 sites), and when the granular material touches the walls of the box, the forces that reach the walls are simply absorbed. This is made in order to control the size of the growing pile. We found that for a range of values of p and Fc, a dip could be found as shown in Fig. 3. It does not depend on the walls; the same result is observed before the pile reaches the walls. For very low values for p (about 0.0001) or very low values of F~, the angle of stability approaches 60 degrees, and the force distribution under the apex is almost uniform. If, on the other hand, p is very large, there will be a peak under the apex, since the angle of stability is very low. If _Pc is infinite, a small dip can still be observed, but the effect is very weak. The angle of stability decreases with growing p, see Fig. 4. This is also what one would expect, since more disorder will cause the heap to broaden. The influence of F~ on the angle of stability depends on the height of the system. For the typical values of p and _Pc given here, systems of 3000 to 5000 layers were considered, and typical values of angle of stability were 30 to 40 degrees. In our model, the angle remains constant only for small sandpiles even when _Pc goes to infinity. For high piles (typically a few thousand layers), very large horizontal forces will emerge. What is the origin of the observed dip? As suggested by others [7], arching seems to
J. Hemmingsson/Physica A 230 (1996) 329-335 2500
,
,
,
333 ,
2000
1500
~
1000
500
0
7000
--
8000
9000
10000 position
11000
12000
13000
Fig. 3. Vertical force distribution under a pile of 2000 layers. Units are chosen such that the mass of one grain times the gravitational constant 9 is 1. The bridge parameter is p = 0.01, compression parameter Fc = 850. This is the average from 50 runs, the angle of stability being 33.8 degrees. The position of the "horns" are decided by the underlying lattice, corresponding to the arch angle/3 in [6].
be crucial. In our model, we found this effect also without introducing the compression mechanism, but the dip is much more pronounced when compression is included. In the presence of compression, the forces can be thought of as (random) walkers. The walkers representing forces that exceed Fc are not random but will always go straight. Since the medium has some order, the large forces will tend to recede from the center of the pile. This reminds a little of the piles suggested in [5] by Liffman, Chan and Hughes, where two slopes are imposed on the pile; one for the upper part and one for the lower part of the pile. In this way they actually get a small dip in the distribution.
5. Conclusions We have found that a simple model is able to reproduce a local minimum (or a dip) of the vertical force field found under the apex of a sandpile. The underlying properties taken into consideration are acting and reacting forces, inhomogeneities and enhancement of contact due to compression. In the experiments by Smid and Novosad [ 1], two kinds of granular materials were used; sand and granulated fertilizer NPK- 1. For both experiments, the main characteristics of the force distribution were the same. Therefore material constants do not necessarily have to influence the overall behavior considerably. We do hope that the approach taken here will
J. Hemmingsson/Physica A 230 (1996) 329-335
334 50
48"
46"
44
¢,
42
40 .'t8 36 34 32 0.0001
0.001
0.01
P
Fig. 4. The angle of stability as a function of p. The height was 2000 layers. Fc was 850 (in units where the mass of one grain times the gravitational constant g is 1). For even higher values of p, the bridge sites can not be considered as independent, and the model is not valid.
be useful for a better understanding of granular packing, and that it might help to refine the theory. Further investigations of a similar model are presented elsewhere [ 13].
Acknowledgements I would like to express my thanks to Hans He,Lmann and Stephane Roux for useful, helpful and inspiring discussions, and to Harald Puhl for support during debugging sessions. This work was supported by a grant from the French government and by the Swedish Natural Research Council.
References [1] [2] [3] [4] [5] [6] [7l [8] [9]
J. Smid and J. Novosad, Int. Chem. Eng. Symp. 63 (1981) D3/V/1. J.M. Hunfley, Phys. Rev. E 48 (1993) 4099. D.C. Hong, Phys. Rev. E 47 (1993) 760. M.K. Haft Menon, in: Non-linearity and Breakdown in Soft Condensed Matter, K.K. Bardhan, B.K. Chakrabarti and A. Hansen, eds. (Springer, Berlin, 1994). K. Liffman, C.Y. Chart and B.D. Hughes, Powder Technol. 72 (1992) 255. S.F. Edwards and C.C. Mounfield, Physica A 226 (1996) 1, 12, 25. S.E Edwards and R.B.S. Oakeshott, Physica D 38 (1989) 88. J.-P. Bouchand, M.E. Cates and P. Claudin, J. Phys. 1 (France) 5 (1995) 639. J.-P. Bouchaud, private communication.
J. Hemmingsson/Physica A 230 (1996) 329-335
335
[10] C.-h. Liu, S.R. Nagel, D.A. Schecter, S.N. Coppersmith, S. Majumdar, O. Narayan and T.A. Witten, Science 269 (1995) 513; S.N. Coppersmith, C.-h. Liu, S. Majumdar, O. Narayan and T.A. Witten, A model for force fluctuations in bead packs, preprint (1995). [11] A.E. Scheidegger, Bull. lASH 12 (1967) 15. [12] H. Takayasu, M. Takayasu, A. Provata and G. Huber, J. Stat. Phys. 65 (1991) 725. [13] J. Hemmingsson, H.J. Herrmann and S. Roux, preprint.