Distribution of dissipated energy in a multi-size granular system under vertical vibration

Powder Technology 260 (2014) 1–6

Contents lists available at ScienceDirect

Powder Technology journal homepage: www.elsevier.com/locate/powtec

Distribution of dissipated energy in a multi-size granular system under vertical vibration Zi-Ang Xie a, Ping Wu a,⁎, Wuhao Yang a, Jingbing Zhao a, Shiping Zhang a, Li Li a, Sen Chen a, Chao Jia b,c, Chuanping Liu b,c, Li Wang b,c a b c

School of Mathematics and Physics, University of Science and Technology Beijing, 100083, China School of Mechanical Engineering, University of Science and Technology Beijing, 100083, China Beijing Engineering Research Centre of Energy Saving and Environmental Protection, Beijing 100083, China

a r t i c l e

i n f o

Article history: Received 24 October 2013 Received in revised form 15 February 2014 Accepted 5 March 2014 Available online 21 March 2014 Keywords: Granular particle Multi-size system Vertical segregation Dissipated energy Distribution model

a b s t r a c t We studied the segregation behavior of multi-size granular particles enclosed in a vertically vibrated container. In the experiment, we applied three kinds of 13X molecular sieve particles. The diameters of these particles are approximately 1:3:6 in the multi-size granular system. Under vibration, the system finally forms a stable hybrid convection mode, and the apparent tendency of the Brazilian nut (BN) effect occurs. The dissipated energy distribution model was introduced, in which the multi-size system was compared with the species competition system. We demonstrated that each kind of particle's top free surface area occupation percentage approximately equals its dissipated energy power percentage in the entire system. The results calculated according to the model well matched the experimental results. © 2014 Elsevier B.V. All rights reserved.

1. Introduction The behavior of strong dissipative systems, such as the granular particle system that is enclosed in a vertically vibrated container, has been studied extensively [1–4]. One of the hottest topics in this field is the vibrated binary particle system [5–7]. In most situations, larger particles tend to rise to the top layer, a phenomenon which is known as Brazilian nut (BN) separation [8]. Conversely, reversed Brazilian nut (RBN) separation occurs in some situations [9–11]. The detailed behavior of the binary particle system appears to be very complex [11]. Parameters such as frequency f, dimensionless acceleration Γ = Aω2/g, material density ratio, and diameter ratio are commonly considered the main influential factors [12]. In the majority of previous literature, the experimental results show clearly identifiable separation patterns of either complete BN (CBN) or complete RBN (CRBN) separations [13, 14]. However, daily experience indicates that in granular piles with an initial random packing state, not all vibration parameter settings will finally form a definite separation result like CBN or CRBN [15–17]. A universal convection is often chaotic specifically in multi-size particle systems, whereas BN or RBN can only reveal themselves as separation tendencies. Therefore, numerous researchers treat BN and RBN as two separation phases and study the transitional states between them [18, ⁎ Corresponding author. E-mail address: [email protected] (P. Wu).

http://dx.doi.org/10.1016/j.powtec.2014.03.028 0032-5910/© 2014 Elsevier B.V. All rights reserved.

19]. A parameter λ between −1 and 1 was proposed to mark the extent to which the transitional states between CRBN and CBN mix in binary particle systems [18,19]. Generally, λ increases when f increases at constant Γ or when Γ decreases at constant f [13,19]. However, the detailed behavior of vibrating particle systems is very complex, and a satisfactory theory or model that describes this behavior completely has not yet been developed, especially for multi-size systems [1]. A study of the multi-size system's behavior mechanism has long been detained as a result of a lack of proper research perspective [20,21]. Recently, an energy conservation relationship has been widely acknowledged to exist between the external input of the vibrator and the dissipation losses within the particle system. The dissipation losses are caused by collision, shearing, and friction [6]. The total sum of the dissipation energies in the particle system can be calculated when f and Γ are determined and when the system stabilizes under vibration [22]. Few studies have been conducted regarding the distribution of the dissipated energy on each particle type involved in the process. In this paper, we propose a dissipated energy distribution model, and the corresponding calculation results are in accordance with the experimental results. 2. Experimental settings We paid particular attention to the distribution of the dissipated energy in each particle type. The granular system that contains three kinds of particles was investigated. We used 13X molecular sieve particles

2

Z.-A. Xie et al. / Powder Technology 260 (2014) 1–6

with the following diameters in the experiment: d1 = 1.12 ± 0.05 mm for the small particles (granular 1); d2 = 3.15 ± 0.05 mm for the medium-sized particles (granular 2); and d3 = 6.60 ± 0.05 mm for the large particles (granular 3). The ratio of the diameters of the three kinds of particles was approximately 1:3:6. The small, medium-sized, and large particles were colored blue, green, and yellow, respectively. Subscripts numbered from one to three were applied to each particle kind as the diameters increased. Initially, the three kinds of particles had the same total mass: m1 = m2 = m3 = 100 g. We eventually increased m1 to 200 g for comparison. The multi-size particle mixture has a bulk density of ρ0 = 0.69 ± 0.05 g/cm3, whereas the 13X molecular sieve particle's material density was ρ = 1.16 ± 0.05 g/cm3. The particle mixture was initially poured randomly into a columnshaped plexiglass container with a bottom radius of R0 = 5.0 cm and a height of 18.0 cm. The upper rim of the container was covered with carmine mask paper for optical detection. Shadowless lamp arrays cast uniform light on the granular piles from above to facilitate the recording of the CCD images. The container was then tightly fixed onto a standard Longdate LD-100XTL vibrator. The vibration was in the z = Asin(ωt) sinusoidal form, and a terminal controlling software determined the output value of f (10 Hz to 100 Hz) and Γ (1 to 8). The basic settings of the experiment are presented in Fig. 1. The height of the mixture's top free surface, i.e. top surface of the particle system, was Hp ≈ 5.5 cm. In our experiment we observed that, in all f and Γ parameter ranges we tested, Hp did not vary greatly, and the top free surface was almost flat. When the whole multi-size granular system was accelerated upwards by vibrations, the utmost speed the system acquired is Aω, neglecting the influence of the container wall's frictions. According to the energy conservation law, the system's mass center could rise (Aω)2/(2 g) = A·(Γ/2) = (g/(8π2))·(Γ2/f2). For a common case in our experiment, in which f = 50 Hz, and Γ = 3, (g/(8π2))·(Γ2/f2) ≈ 0.04 cm. Therefore, Hp ≈ 5.5 cm would not change much, and ρ0 remained almost constant.

3. Experimental results A final stable separation result was achieved within 30 s, and we observed no definite CBN or CRBN separation in this multi-size system. The particles mixed under vibration in all f and Γ parameter ranges and their motions constituted the double-cell cells throughout the entire system. A similar situation has long been observed in single-type particle systems, which exhibit two clear symmetric convection cells [23,24], or in binary systems, in which double-cell convection helps to form different separation phases [19]. However, the behavior of the multi-size particle system convection has not been studied thoroughly in previous literature. The situation in our multi-size particle system experiments is termed as the hybrid convection mode, in which all three kinds of particles were involved in the convection and were mixed. The features of BN separation were clear, although they were incomplete. The large particles still had a tendency to segregate in the upper layer, and only a small fraction of these particles sank to the bottom layer. Correspondingly, the small particles occupied much of the lower space, and only a few particles emerged in the top layer through convection. The medium-sized particles were mainly found between the large and small particle layers. As a result of convection, the spatial distributions of the three kinds of particles gradually stabilized. The spatial distributions of the three kinds of particles were actually based on a dynamic equilibrium, which is similar to that of binary particle systems [19]. In previous studies, similar behaviors in binary particle systems were studied by an extension of the usual statistical-mechanics concepts [25]. It was clarified that the steady probability of the granular system does not depend on the details of the configuration, but only on the number of particles of each of the two species [25]. When f increased at constant Γ, the convection motions of small and medium-sized particles were less agitated, and the appearance of these particles in the top layer contracted because of the percolation effect [1], i.e., the BN effect was reduced and vice versa. During vibration, the sum of the dissipated energy in the system should be equal to the external energy input when vibration is exerted for a long time and the system becomes stable, according to the energy conservation relationship [26,27]. In Ref. [22], a model based on both viscoelastic and two-phase theories was developed to describe the dissipation energy of the binary mixtures. In this model, only particle– particle collisions were considered because they were the main constitutions of energy dissipation during vibration. Particle–container collisions and frictions were therefore neglected [22,28]. We reorganized and modified the mathematical expressions of the dissipated energy power P, as proposed in Ref. [22], to conform to the multi-size system. First, the multi-size granular system was viewed as a whole, with a total mass of m, in our hybrid convection mode instead of considering only the interaction between two kinds of particles. Second, we focused on the relationships of P with f and with Γ and determined the distribution of P to the three kinds of particles. The expression of P is as follows:   pffiffiffiffiffi 2 3 2 mg2 Γ2 2η F l F l  W−W þ F l  W η p ffiffiffiffi ffi   P¼    2 m 16π2 f 2 F l −W 2 þ F l  Wη2

Fl ¼

qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi Yπ2 R40 ρ0 =m

2

W ¼ 4π mf

Fig. 1. (a) Schematic of the experimental facilities: 1) CCD camera, 2) shadowless lamps, 3) carmine mask paper, 4) Plexiglas container, 5) mixed particles inside the container, 6) standard vibrator. (b) Hybrid convection mode. (c) Image obtained by the CCD camera placed above the container (m1 = m2 = m3 = 100 g, f = 30 Hz, Γ = 2.5). (d) Same image as (c) in less detail for pixel calculation.

2

ð1Þ

ð2Þ

ð3Þ

where P represents the dissipated energy power; m is the mass of the entire particle system; and η is the loss factor of the mixture, which refers to the energy loss ratio of the multi-size granular system per vibration cycle, due to the internal damping and dissipation. Y is the Young's modulus of the multi-size granular system.

Z.-A. Xie et al. / Powder Technology 260 (2014) 1–6

3

Fig. 2. Relationships between ki(z), i = 1, 2, 3, and z. The horizontal axis is the section height, and the vertical axis marks the area coverage ratios. The BN effect is apparent.

Considering the orders of magnitudes of ordinary materials in Ref. [22], in the multi-size system we gave an estimation that η ≈ 0.1, and Y ≈ 1.0 × 105 Pa. The increase of f at constant Γ can be calculated by using Eq. (1), and P reaches its peak when fm

pffiffiffiffiffiffiffiffiffiffiffiffi R2 ¼ Y  ρ0  0 2m

Pm ¼

m2 gΓ 2 pffiffiffiffiffiffiffiffiffiffiffiffi : 2πηR20  Y  ρ0

ð4Þ

ð5Þ

Accordingly, we can obtain ωm  H p ¼

pffiffiffiffiffiffiffiffiffiffiffi Y=ρ0 ¼ v

ð6Þ

where ωm = 2πfm, and v represents the traveling speed of the solid wave if the granular pile is considered solid bulk, with a bulk density of ρ0. According to Eq. (4), fm = 34.6 Hz when m1 = m2 = m3 = 100 g. fm = 26.0 Hz when m1 = 200 g and m2 = m3 = 100 g. In Eq. (1), P decreases when f increased above fm at constant Γ. P increases when Γ increases at constant f above fm. Therefore, fm can be considered as an inherent frequency in the multi-size granular system.

Accordingly, P generally decreases when the system's BN tendency increases. 4. Dissipated energy distribution model After vibration was exerted on the system, the three kinds of particles began to interact (mainly colliding) with one another. The hybrid convection mode was finally formed. When the system and vibration parameters are determined according to Eq. (1), P can be directly determined as a constant. Previous literature determined that the dissipated energy within the granular system is mainly generated via collisions between particles, whereas frictions play only a very small part in this process [28]. In this study, all collisions in the granular system were considered direct binary collisions. Studies on the dissipative force between viscoelastic spheres demonstrated that when two particles collide, their corresponding dissipation energies are proportional to the local deformation rate of each particle [26,29]. Basically, the collision energy is absorbed by the harmonic lattice vibrations in the particles [27], and the dissipated energy is approximately proportional to the diameter ratio of the two particles [26,30]. In our experiments, P was in fact distributed to the three kinds of particles. This process is similar to the species competition process in nature. We viewed the three kinds of particles as different species ‘competing’ for the dissipated energy ‘resource’ P. Each particle has an initial

4

Z.-A. Xie et al. / Powder Technology 260 (2014) 1–6

occupation of Pi(0) = 0, i = 1, 2, 3. In a limited time, i.e., from the start of vibration to the final formation of the hybrid convection mode, a balance between these particles was gradually formed. The Lotka–Volterra model, which is concerned with the distribution of resources to the ecological population, was applied. Thus, the mathematical solution of this model remains stable [31,32]. After the vibration began, the three particle ‘species’ obtained the dissipated energy at various increase rates ri. The values of Pi(t) (i = 1, 2, 3) stabilized gradually as the sum of Pi(t) gradually approached P. Based on the Lotka–Volterra model, the dynamical equation was derived as follows: ! X dP i ðt Þ ¼ r i P i ðt Þ 1− P i ðt Þ=P : dt i

coverage ratios ki(z) are presented in Fig. 2. k1(z) + k2(z) + k3(z) = 1 is evident. The dissipated energy power in the upper layers is slightly greater than that in the lower layers because more large particles are present in the upper layers [26,30]. In the layer section from z to z + dz, the dissipating energy power is assumed to contain Ceξzdz. C is a constant slightly smaller than P/Hp, and ξ is a positive small quantity. We then obtain Z P¼

Hp 0

Z ξz Ce dz≈

Hp 0

  1 2 C ð1 þ ξzÞdz ¼ C  H p þ ξH p : 2

ð10Þ

ð7Þ

For ri, i = 1, 2, 3, which should be related to di, mi, and P. Prior to reaching the balance, in which P1(t) + P2(t) + P3(t) = P, the values of Pi(t) increase more quickly with a greater di in direct particle i–j collisions. A kind of particle becomes more sensitive to a change in P when the particle is smaller and has a greater total mass [19]. Therefore, we first set a factor of C1di2 as the particle–particle direct collision cross section. C1 is a constant related to the material's elasticity and frictional properties. Second, we set another factor of P Ai, in which Ai represents the functions of all three kinds of particles in terms of geometric size and total mass. For Ai, i = 1, 2, 3, which represent the sensitivity of a kind of particle to a change in P. Ai should be qualitatively proportional to mi, i.e., a greater total mass indicates that a particle type becomes more sensitive to a change in P. This sensitivity should also be related to the particle i–particle j direct collision probability τij, which is the reciprocal of the contact duration time tij [23]. In particle i–particle j direct collisions, the contact duration time tij is proportional to (meff2/Reff)1/5, where meff = mimj / (mi + mj) and Reff = RiRj / (Ri + Rj), i, j = 1, 2, 3, with Ri as each particle kind's radius. We then obtain τij = 1/tij [23]. Given that we are studying the sensitivity of particle i at present, we multiplied τij by particle j's total mass percentage mj/m and obtained the sum of j from 1 to 3. We then multiplied this value by C2mi, in which C2 is also a constant. The three kinds of particles have the same material density; therefore, Ai is apparently greater when mi is greater or when di is smaller, and vice versa. The expressions of ri and Ai are summarized as follows: 2

r i ¼ C 1  di  P

Ai ¼ C 2  mi 

Ai

ð8Þ

3 X mj j¼1

0

di d j @ m di þ d j

,0 3 3 12 11=5 d d @ i j A A : d3i þ d3j

ð9Þ

The numerical values of C1 and C2 in Eqs. (8) and (9) can be determined by using the least squares method on the experimental results (presented in subsequent sections). Originally, Pi(0) = 0 and i = 1, 2, 3. After vibration was initiated, Pi(t) began to increase with t at different rates, which were determined by ri. Finally, a balance, i.e., the hybrid convection mode, was established, where P1(t) + P2(t) + P3(t) = P. In the following section, we continued to express the dissipated energy powers ultimately acquired by the three kinds of particles, i.e., Pi(t)|t → ∞, as Pi, i = 1, 2, 3. To minimize all possible deviations, we carefully vacuumed the particle layers to avoid influencing the authentic packing structure of the final separation results. We then photographed sections of the separation result at heights of 5.5 (i.e., the free surface), 5, 4, 3, 2, and 1 cm by using the CCD camera. Each pixel within the scope of the carmine ring on the picture was then analyzed by using MATLAB 2011b programs to obtain the area coverage ratio of each kind of particle ki(z), i = 1, 2, 3 in these sections. The calculated results of the area

Fig. 3. Comparison between the experimental results (dots) and the dissipated energy distribution model results (lines), in ki(Hp) when the f, Γ, and the total mass of the small particle change. The diameters of the three kinds of particles are almost at 1:3:6 ratio. (a1) to (a3): m1 = m2 = m3 = 100 g. (b1) to (b3): m1 = 200 g, m2 = m3 = 100 g.

Z.-A. Xie et al. / Powder Technology 260 (2014) 1–6

5

Fig. 4. Comparison among the dissipated energy distribution model results calculated for the three kinds of granular particles. (a1) to (a3): m1 = m2 = m3 = 100 g. (b1) to (b3): m1 = 150 g, m2 = m3 = 100 g. (c1) to (c3): m1 = 200 g, m2 = m3 = 100 g. The horizontal and vertical axes represent f and Γ, whereas the z-axis represents the value of Pi/P, as marked by color. P1/P evidently increases when more small particles are present.

In Eq. (10), ξ is a positive small quantity, and z is typically less than H p = 5.5 cm, therefore, e ξz ≈ 1 + ξz. Based on Eq. (10), we obtain C¼

  P P 1  ≈  1− ξH p : 1 Hp 2 Hp  1 þ ξHp 2

ð11Þ

Z

ð12Þ

Hp

Z

ξz

0

Z

ki ðzÞ  Ce dz

Hp

ξz

ðai z þ bi Þ  e dz  Hp 1  ξHp ξz e −1 Ce dz ξ 0    a  ξH ai H p þ bi eξHp −bi þ i 1−e p ξ ¼ eξHp −1        ai H p þ bi 1 þ ξH p − ai H p þ bi ¼ ai H p þ bi ¼ ki Hp : ≈ ξHp

Pi ¼ P

Therefore, C is slightly smaller than P/H p , at approximately (1/2)ξP. As shown in Fig. 2, ki(z) has an approximate linear relationship with z in all of the cases studied. Therefore, we assume that: ki ðzÞ ¼ ai z þ bi ; i ¼ 1; 2; 3:

In Eq. (12), ai and bi are constants. When the dissipated energy power ratio, or Pi/P, of a kind of particle in the entire system is studied, the following is obtained:

¼

0

ð13Þ

6

Z.-A. Xie et al. / Powder Technology 260 (2014) 1–6

Intriguingly, Pi/P = ki(Hp) in the multi-size granular system with a hybrid convection mode. With this important relationship taken into consideration, experimental data were applied, and the values of C1 and C2 in the dissipated energy distribution model were determined in Eqs. (8) and (9) by using the least squares method, with C1 = 1.05 × 10− 2 and C2 = 1.20 × 105. Fig. 3 shows comparisons between the results calculated according to the dissipated energy distribution model in Eqs. (8) and (9) and the experimental results. The calculation results in Fig. 3 are acquired via numerical calculations of Eqs. (7)–(9), with the Runge–Kutta method. Generally, fewer small particles appeared on the top free surface when f was increased from 30 Hz to 70 Hz at constants Γ and mi and i = 1, 2, 3. Therefore, the system's BN tendency was enhanced, and P decreased. More small particles appeared on the top free surface when Γ was increased from 1 to 7 at constants f and mi and i = 1, 2, 3. The system's BN tendency was reduced, and P increased. Even more small particles appeared on the top free surface when m1 increased from 100 g to 200 g, as compared in Fig. 3(a1) to (a3) and Fig. 3(b1) to (b3). As m increased, P also increased, as calculated by using Eq. (1). However, the calculated results deviated slightly from the experimental data when the value of Γ is small, especially when Γ b 3, because the hybrid convection mode cannot stay robust when the vibration amplitude is weak. As a result, the BN separation effect becomes dominant. The system is in fact close to a multi-size CBN. The particles do not gain enough energy to form convection, and the movements of the small and middle-size particles are suppressed. Thus, the energy dissipation model mainly adapts to systems with an apparent overall convection to form the hybrid convection mode. The Pi/P–f–Γ relationships are presented in Fig. 4, with i = 1, 2, 3. As shown in Fig. 4(b1) to (b3), we also considered a situation in which m1 = 150 g, m2 = m3 = 100 g and fm = 29.6 Hz. As shown in Fig. 4, the small particles apparently obtain more dissipated energy from the system when total mass increases. Thus, more of these particles can appear on the top free surface. The introduction of the dissipated energy distribution model is based on the stable separation results in this study. However, additional challenges arise when the multi-size particle system exhibits more complicated behavior, such as oscillation and chaos [21]. 5. Conclusions A model for the distribution of the dissipated energy in a multi-size granular system is established in this study. The dissipated energy percentage of a certain kind of particle in the system is approximately equal to its coverage ratio on the top free surface. The calculated results agree with the experiment results. Our model linked the properties of granular dissipated energy distribution to theories with regard to the continuous process of species distribution to obtain resources. The method is not only applicable to the three kinds of particles that form a mingled convection, but it can also be used in future research to describe more complex segregation behavior of multi-size particle systems. Acknowledgments This article is sponsored by the National Natural Science Foundation of China (No. 51076010) and National Basic Research Program of China (973 Program 2012CB720406).

References [1] H.M. Jaeger, S.R. Nagel, Physics of the Granular State, Science 255 (1992) 1523–1531. [2] N. Burtally, P.J. King, M.R. Swift, Spontaneous air-driven separation in vertically vibrated fine granular mixtures, Science 295 (2002) 1877–1879. [3] C. Liu, L. Wang, P. Wu, M. Jia, Effects of gas flow on granular size separation, Phys. Rev. Lett. 104 (2010) 188001. [4] X. Cheng, Packing structure of a two-dimensional granular system through the jamming transition, Soft Matter 6 (2010) 2931–2934. [5] P. Eshuis, D. van der Meer, M. Alam, H.J. van Gerner, K. van der Weele, D. Lohse, Onset of convection in strongly shaken granular matter, Phys. Rev. Lett. 104 (2010) 038001. [6] D. Bi, J. Zhang, B. Chakraborty, R.P. Behringer, Jamming by shear, Nature 480 (2011) 355–358. [7] E.J.R. Parteli, J. José, S. Andrade, H.J. Herrmann, Transverse instability of dunes, Phys. Rev. Lett. 107 (2011) 188001. [8] A. Rosato, K.J. Strandburg, F. Prinz, R.H. Swendsen, Why the Brazil nuts are on top: size segregation of particulate matter by shaking, Phys. Rev. Lett. 58 (1987) 1038. [9] T. Shinbrot, F.J. Muzzio, Reverse buoyancy in shaken granular beds, Phys. Rev. Lett. 81 (1998) 4365–4368. [10] D.C. Hong, P.V. Quinn, Reverse Brazil nut problem: competition between percolation and condensation, Phys. Rev. Lett. 86 (2001) 3423–3426. [11] A. Kudroll, Size separation in vibrated granular matter, Rep. Prog. Phys. 67 (2004) 209–247. [12] M. Pica Ciamarra, A. Coniglio, M. Nicodemi, Thermodynamics and statistical mechanics of dense granular media, Phys. Rev. Lett. 97 (2006) 158001. [13] A. Breu, H.M. Ensner, C. Kruelle, I. Rehberg, Reversing the Brazil-nut effect: competition between percolation and condensation, Phys. Rev. Lett. 90 (2003) 014302. [14] D. Sanders, M. Swift, R. Bowley, P. King, Are Brazil nuts attractive? Phys. Rev. Lett. 93 (2004) 208002. [15] W.A.M. Morgado, E.R. Mucciolo, Numerical simulation of vibrated granular gases under realistic boundary conditions, Phys. A Stat. Mech. Appl. 311 (2002) 150–168. [16] L. Staron, J.P. Vilotte, F. Radjai, Preavalanche instabilities in a granular pile, Phys. Rev. Lett. 89 (2002) 204302. [17] A. Kabla, G. Debrégeas, Contact dynamics in a gently vibrated granular pile, Phys. Rev. Lett. 92 (2004) 035501. [18] Z.A. Xie, P. Wu, S. Wang, Y. Huang, S. Zhang, S. Chen, C. Jia, C. Liu, L. Wang, Behaviour of a binary particle system under the effects of simultaneous vertical vibration and rotation, Soft Matter 9 (2013) 5074–5086. [19] Z.A. Xie, P. Wu, S. Zhang, S. Chen, C. Jia, C. Liu, L. Wang, Separation patterns between Brazilian nut and reversed Brazilian nut of a binary granular system, Phys. Rev. E 85 (2012) 061302. [20] M.J. Metzger, B. Remy, B.J. Glasser, All the Brazil nuts are not on top: vibration induced granular size segregation of binary, ternary and multi-sized mixtures, Powder Technol. 205 (2011) 42–51. [21] P. Porion, N. Sommier, A.-M. Faugère, P. Evesque, Dynamics of size segregation and mixing of granular materials in a 3D-blender by NMR imaging investigation, Powder Technol. 141 (2004) 55–68. [22] T. Yanagida, A.J. Matchett, J.M. Coulthard, Energy dissipation of binary powder mixtures subject to vibration, Adv. Powder Technol. 12 (2001) 227–254. [23] A. Alexeev, A. Goldshtein, M. Shapiro, The liquid and solid states of highly dissipative vibrated granular columns: one-dimensional computer simulations, Powder Technol. 123 (2002) 83–104. [24] S.S. Hsiau, C.H. Chen, Granular convection cells in a vertical shaker, Powder Technol. 111 (2000) 210–217. [25] A. Prados, J.J. Brey, Dynamics and steady state of a vibrated granular binary mixture model, Europhys. Lett. 64 (2003) 29–35. [26] N.V. Brilliantov, F. Spahn, J.M. Hertzsch, T. Poschel, Model for collisions in granular gases, Phys. Rev. E 53 (1996) 5382. [27] W.A.M. Morgado, I. Oppenheim, Energy dissipation for quasielastic granular particle collisions, Phys. Rev. E 55 (1997) 1940. [28] Y. Srebro, D. Levine, Exactly solvable model for driven dissipative systems, Phys. Rev. Lett. 93 (2004) 240601. [29] D. Antypov, J.A. Elliott, B.C. Hancock, Effect of particle size on energy dissipation in viscoelastic granular collisions, Phys. Rev. E 84 (2011) 021303. [30] Y. Tsuji, T. Kawaguchi, T. Tanaka, Discrete particle simulation of two-dimensional fluidized bed original research article, Powder Technol. 77 (1993) 79–87. [31] C. Zhu, G. Yin, On competitive Lotka–Volterra model in random environments, J. Math. Anal. Appl. 357 (2009) 154–170. [32] X.A. Zhang, L. Chen, The global dynamic behavior of the competition model of three species, J. Math. Anal. Appl. 245 (2000) 124–141.