Granular Matter (2011) 13:149–158 DOI 10.1007/s10035-010-0204-9 Discrete element modelling of grain flow in a planar silo: influence of simulation parameters C. González-Montellano · F. Ayuga · J. Y. Ooi Received: 15 March 2010 / Published online: 23 October 2010 © Springer-Verlag 2010 Abstract There is extensive engineering literature concerning the prediction of pressure and flow in a silo. The great majority of them are based on continuum theories. The friction between the stored material and the silo wall as well as the inclination of the hopper at its base are considered to be the most influential parameters for the flow pattern within the silo. In this paper, the filling and discharge of a planar silo with a hopper at its base has been modelled using DEM. The aim is to investigate the influence of DEM model parameters on the predicted flow pattern in the silo. The parametric investigation particularly focused on the hopper angle of inclination and the contact friction between particles and walls. The shape of the particles was also considered by comparing spherical and non-spherical particles, thus providing an insight into how particle interlocking might influence solids flow behaviour in silos. The DEM computations were analysed to evaluate the velocity profiles at different levels as well as the wall pressure distribution at different stages during filling and discharge. A detailed comparison reveals several key observations including the importance of particle interlocking to predict a flow pattern that is similar to the ones observed in real silos. Keywords Discrete element method · Granular flow · Silo design · Mass flow · Funnel fow C. González-Montellano (B) · F. Ayuga BIPREE research group. ETSI Agrónomos, Universidad Politécnica de Madrid, Ciudad Universitaria S/N, 28040 Madrid, Spain e-mail: carlos.gonzalez.montellano@upm.es J. Y. Ooi Institute for Infrastructure and Environment, University of Edinburgh, The King’s Buildings, Edinburgh EH9 3JL, UK 1 Introduction and objectives Silos are containment structures widely used in industry for storing granular solids. The type of flow pattern developed in a silo as it empties under gravity is one of the main considerations in the design of a silo. In general, the type of flow can be classified as mass flow, funnel flow or mixed flow [1]. Mass flow occurs when all the stored particles are simultaneously in motion during discharge. In funnel flow, only particles within an internal channel above the hopper outlet are in motion, whilst the rest of the particles surrounding the channel remain at rest. Mixed flow is an intermediate situation where the flow channel reaches the vertical wall of the silo at a point below the solid surface. In silo design, the flow pattern occurring during discharge is taken as a function of the variables characterising the hopper geometry (inclination, shape, outlet size) and those related to the friction within the particles and between the particles and the silo walls [2]. Based on some of these variables and using continuum models, Jenike [4,5] developed design charts to estimate the expected flow pattern. Figure 4b shows an example of these charts, taken from Eurocode EN 1991 4 Standard [1] for the case of wedge hoppers. The pressure distribution exerted by the stored material on the silo wall is another important calculation in silo design. The pressure distributions are usually assessed in practice using analytical procedures based on Jansen’s or Reimbert’s equations [6,7], although the use of numerical models (such us the finite element method [8–11]) to evaluate these values is also possible. However, existing numerical and analytical models do not often predict accurately the flow pattern within the silo or the pressure exerted by the stored material on the silo walls [12]. The discrete element method (DEM) is a numerical technique that can be used to study the behaviour of granular 123 150 materials at particle interaction level and its use is becoming increasingly popular [13,14]. Its capability has been demonstrated in many cases provided that parameters involved in the simulation are established appropriately [15]. However simulation of a real handling situation is not always feasible, mainly due to the high compute time consumed in a DEM simulation. The aim of this work is to evaluate the effectiveness of a DEM model of a planar silo to predict the flow pattern and the pressures exerted by the stored material. The model parameters that most influence the flow pattern are studied: friction between particles and walls and hopper inclination. The influence of the particle shape was also considered by comparing spherical and non-spherical particles. 2 Methodology The discrete element method, based on the work developed by Cundall and Strack [16], is a numerical technique used to model the mechanical behaviour of granular assemblies. It is based on the use of an explicit calculation scheme in which the interaction between particles is monitored contact by contact (by means of a force-displacement law) and where the motion of the particles is modelled particle by particle (by applying Newton’s second law). The simulations in this study were conducted using EDEM version 2.1 software (DEM Solutions, 2009). DEM methodology has been covered extensively in the literature, so only the specific implementation and model parameters are described below. 2.1 Model and particle geometry In this study, the DEM model represents a slice of a silo with a hopper at its base. This slice is considered to be a simplification of a 3D planar silo and will reduce the computational time significantly. The DEM model is summarized in Fig. 1. The vertical section of the silo model has a width D of 0.4 m and was filled to a height of H = 0.5–0.55 m. The vertical section has a hopper at its base with an angle of inclination β and an outlet of 0.06 m. The angle of inclination can take two values: β = 45◦ and β = 60◦ . The slice is set to be 5 mm thick so that only a single layer of particles is allowed within the model: the model thus resembles a 2D planar representation. The model thus ignores any out of plane movements which is deemed an adequate approximation for silos of planar cross-section. The model is bounded by physical walls in all its faces (Fig. 1) except in the upper boundary. Lateral walls in the vertical section and the hopper section are each divided into 10 equal segments. For each segment, the mean contact force is calculated to evaluate the pressure distribution acting on the silo walls. The wall representing the hopper outlet (O) 123 C. González-Montellano et al. consists of a single section which is removed to initiate the discharge process in an unsteady-state condition. Frictionless frontal walls (W ′ ) were added to constrain out of plane displacements of the particles. Two types of particles were considered: single sphere (SP) or paired-sphere (PP). SP particles consist of a single sphere of 5 mm diameter whilst PP particles are formed by clustering two spheres of 5 mm diameter to give an aspect ratio of 1.5. 2.2 Contact model and baseline parameters used in simulations The force-displacement law used is based on the HertzMindlin no-slip contact model with viscous damping and a frictional slider in the tangential direction of the contact [17]. The model parameters for this study are shown in Table 1. These values have been taken from the literature ([18–20]) and will be considered as the baseline parameters for the parametric study which is defined later. 2.3 Particle generation procedure The model silo is initially filled by generating all particles at one instant. All particles are generated in such a way that the centre of each sphere (SP) or the line joining the centre of the combined spheres (PP) is contained in the mid-plane of the slice considered. The SP particles are generated randomly within this plane whilst the PP particles are generated in a rectangular lattice with a random orientation of the particle longitudinal axis. The particles are then allowed to settle down under gravity until the system reaches a static situation. 2.4 Parametric analysis A parametric study has been undertaken to investigate the influence of key factors on flow pattern and wall pressure in silos. These factors are: hopper angle of inclination (β = 45◦ /β = 60◦ ), static friction coefficient between particles and walls (W ) (µ = 0.13/µ = 0.9) and type of particle (SP/PP). All possible combinations of these variables are shown in Table 2 and lead to eight different model configurations to be calculated. 2.5 Flow pattern analysis Flow pattern analysis is sometimes a difficult task since different methods do not always reach the same and definitive conclusion. For that reason, different methods have been considered in this paper and are described below. Discrete element modelling of grain flow in a planar silo: influence of simulation parameters 151 Fig. 1 General characteristics of the DEM silo model Table 1 Baseline model parameters used in simulations Density (kg/m3 ) Young’s modulus (Pa) 2,500 7.1 · 1010 0.22 2,700 6.9 · 1010 0.33 Restitution coefficient (e) Static friction coefficient (µ) Rolling friction coefficient (µr ) Particle–particle 0.95 0.1 0.045 Particle–walls (W or O) 0.83 0.13 0.045 Particle–walls (W ′ ) 0.83 0 0 Material Poisson’s ratio Material properties Particle (SP and PP) Walls (W, W′ and O) Interacting materials Interaction properties Table 2 Adopted values for the parametric analysis Model β Particles Particle-wall friction (µpw ) M1 45◦ SP 0.13 PP 0.13 SP 0.13 PP 0.13 M2 0.9 M3 M4 M5 0.9 60◦ M6 M7 M8 0.9 0.9 123 152 C. González-Montellano et al. 2.5.1 Velocity profile and Mass Flow Index (MFI) 3 Results and discussion For each model defined in Table 2, three control volumes (CV1, CV2, CV3) at different heights (HCV1 , HCV2 and HCV3 ) of the model have been considered (see Fig. 1). Each control volume consists of several “bins” where the vertical velocities of the particles in each bin are spatially averaged at each output stage during the simulation. Thus, the velocity profiles at these locations at a particular simulation time can be obtained by plotting the vertical velocity in each bin against the horizontal coordinate (x) associated to each bin’s central point. Additionally, each velocity profile can be summarized through the Mass Flow Index (MFI), with is defined as [3]: 3.1 Analysis of velocity profile and Mass Flow Index (MFI) MFI = vwall vcentreline (1) For every control volume, vwall and vcentreline are, respectively, the particle velocity at the extreme bins and at the central bin. According to Johanson and Jenike [21], values of MFI > 0.3 are indicative of mass flow whereas values of MFI < 0.3 are indicative of funnel flow. In the case of the control volume CV3, the value of the velocity Vwall was calculated from the component parallel to the hopper wall. 2.5.2 Direct observation of discharge It is possible to visually assess the predicted flow pattern within the silo. To facilitate this task, all particles in the silo were divided into different horizontal layers at the beginning of the discharge process. These layers were coloured alternatively with two contrasting colours (Fig. 5) so that particle movement can be easily identified. 2.5.3 Residence time curves The residence time is the time a particle takes to exit the silo from the beginning of discharge. In this paper, the residence time curve is defined as the curve representing the residence time of all particles located at a particular vertical coordinate at the beginning of discharge against their initial horizontal coordinates. The shape of the residence time curves is indicative of the homogeneity of the particle velocities at this level, which is related to the flow pattern. In addition, it is possible to assess the presence of an internal flow channel by comparing residence time curves associated with different levels [22]. 123 The velocity profile predicted for each model in Table 2 has been determined. Figure 2 shows an example of these results for control volumes CV1 and CV3 in models M6 and M8 for 10 discrete instants during the initial period of discharge. Additionally, the mean velocity profile for this period is shown with a black line. Figure 2 shows how models where SP particles were used (M6) have a more fluctuating velocity profile, whereas those models using PP particles (M8) show a more stable velocity profile. This observation shows that even with spatial averaging, SP particles give rise to much more fluctuation of the velocity field. Figure 3 compares the mean velocity profiles during the discharge period adopted in Fig. 2 for all models in this study and for each of the three control volumes (CV1, CV2, CV3) considered. Additionally MFI value has been determined for every mean velocity profile and is shown in Fig. 3. It is indicative of the flow pattern as well as of the level of velocity variation at the height level considered. Using MFI as a first indication of flow pattern, Fig. 4a shows the MFI values, calculated from the mean velocity profiles shown in Fig. 3, against the associated control volumes, which represent the heights in the silo where the MFI values are evaluated. As a general trend, the value of MFI decreases (flow pattern approaches funnel flow) as the friction coefficient increases or as the hopper becomes shallower. Furthermore, the MFI decreases with depth for all cases which indicates a more uniform mass flow towards the upper region of the silo. From Fig. 4a, since the index MFI is always higher than 0.3 in CV1, no internal flow channel is predicted to develop at this height level in any of the models. Looking further at CV2 and CV3: if MFI > 0.3 in both control volumes, the flow pattern in the silo can be considered as mass flow. On the contrary, the flow can be considered as mixed flow when MFI < 0.3 in CV3 or in both CV2 and CV3. In the particular case of models M7/M5, flow pattern should be characterized strictly as mixed flow since MFI < 0.3 for the volume CV3. However, both models predicted MFI well above 0.3 in the vertical section which demonstrates a predominantly mass flow situation and have thus been classified as such. Figure 4b shows the flow pattern expected for each model according to the prediction of the Eurocode [EN 1991 4 [1]]. In order to determine the expected flow pattern, the microscopic values given in Table 2 for the friction between particles and walls (µpw ) have been converted into macroscopic wall friction (µw ). An accurate conversion between these values would require the simulation of a shear test for different µpw values [14] to produce a calibration curve for each Discrete element modelling of grain flow in a planar silo: influence of simulation parameters 153 Time within discharge 0.24 (a) 0.05 s 0.16 s 0.51 s 0.62 s 0.97 s 1.08 s 0.28 s 0.39 s 0.74 s 0.85 s 1.20 s MEAN 0.48 M6-CV 0.2 0.4 0.16 0.32 0.12 0.24 0.08 0.16 0.04 0.08 0 -0.2 0.24 -0.15 (c) -0.1 -0.05 0 0.05 0.1 0.15 0.2 0 -0.2 0.48 M8-CV1 0.2 0.4 0.16 0.32 0.12 0.24 0.08 0.16 0.04 0.08 0 -0.2 -0.15 -0.1 -0.05 0 0.05 0.1 0.15 0.2 (b) -0.15 (d) 0 -0.2 M6-CV3 -0.1 -0.05 0 0.05 0.1 0.15 0.2 -0.05 0 0.05 0.1 0.15 0.2 M8-CV3 -0.15 -0.1 Fig. 2 Velocity profiles at Levels CV1 and CV3 in models M6 (β = 60◦ , µ = 0.9, SP) and M8 (β = 60◦ , µ = 0.9, PP) during initial discharge. (Horizontal axis: Horizontal coordinate (x)/Vertical axis: Vertical velocity (m/s)) M1 (0.9) M5 (0.91) M2 (0.66) M6 (0.47) M3 (0.98) M7 (0.96) M4 (0.5) M8 (0.55) M1 (0.57) M5 (0.63) 0.16 0.16 0.14 0.14 0.12 0.12 0.10 0.10 0.08 0.08 0.06 M2 (0.56) M6 (0.35) M3 (0.93) M7 (0.91) M4 (0.23) M8 (0.29) 0.06 - CV1 - 0.04 0.04 0.02 - CV2 - 0.02 0.00 -0.2 -0.15 -0.1 -0.05 M1 (0.4) M5 (0.2) 0 M2 (0.33) M6 (0.09) 0.05 0.1 M3 (0.54) M7 (0.24) 0.15 0.2 0.00 -0.2 -0.15 -0.1 -0.05 0 0.05 0.1 0.15 0.2 M4 (0.03) M8 (0.01) 0.40 0.35 0.30 0.25 0.20 0.15 0.10 - CV3 - 0.05 0.00 -0.2 -0.15 -0.1 -0.05 0 0.05 0.1 0.15 0.2 Fig. 3 Mean velocity profiles for all models (Horizontal coordinate (x) vs. vertical velocity (m/s), with MFI values in brackets) 123 154 C. González-Montellano et al. 1.00 0.90 M1 (Mass Flow) 0.80 M2 (Mass Flow) 0.70 M3 (Mass Flow) MFI 0.60 M4 (Mixed Flow) 0.50 M5 (Mass Flow) 0.40 0.30 M6 (Mixed Flow) M7 (Mass Flow) MFI = 0,3 0.20 M8 (Mixed Flow) 0.10 0.00 CV1 CV2 CV3 (a) MFI values predicted in the silo (b) Expected flow pattern (EN 1991 4 [1]) as a function of hopper wall friction µ coefficient and hopper half angle β Fig. 4 Comparison of DEM predicted flow pattern with Eurocode (EN 1991 4 [1]) particle type. However, these simulations were not carried out and an approximate procedure was adopted instead: the equivalence between these two variables was approximated using the relationship given for a similar case in [2], considering the “with rotation” case when SP were used and the “no rotation” when PP were considered. Using this relationship, all microscopic values given in Table 2 results in the same value of macroscopic friction (µpw = µw ) except for the case of SP and µpw = 0.9, which is related to a macroscopic value of µw = 0.36. Considering Fig. 4b, Models M1/M3 are clearly expected to depict mass flow whereas model M8 is expected to depict funnel flow. In model M6 no clear funnel flow is expected, although it can be considered to be closer to this flow pattern. In all other cases the flow pattern could be one or the other. However, considering the values of MFI given in Fig. 4a, models M5/M7 would be likely to be mass flow and models M2/M4 would follow a funnel flow pattern. It should be noted that mixed flow can be perceived as a form of funnel flow in which the internal funnel widens to reach the walls, thus becoming mixed flow. This situation may occur when a funnel flow is expected but there is a relatively lager outlet in relation to the silo width and an expanding flow channel. Taking into account these considerations, the DEM simulations predicted mass flow for both M1 and M3, which agrees with the expectation. For models M6/M8, the present DEM results predicted mixed flow which also matches with the expectation. The results also suggest that non-spherical PP particles tend to produce flow pattern that is closer to the ones observed in experiments. When funnel flow is expected (M2, M4, M6, M8), MFI values for those models where PP were used (M4 and M8) are lower than those for models where SP were used (M2 and M6). When mass flow is more probable (M1, M3, M5, M7), MFI values for those models where PP were used 123 (M3 and M7) are higher than those for models when SP were used (M4 and M6). This observation highlights the key role of particle shape and particle interlocking in influencing the flow pattern. 3.2 Direct observation of the discharge process Figure 5 shows different snapshots of the discharge process for all the models considered for two discrete instants: T1 = 1 s and T2 = 3 s. Models M1, M3, M5 and M7 were likely to develop mass flow according to Fig. 4a,b, so the initial coloured layers would be expected to remain largely horizontal during most of the discharge process. This expected behaviour is clearly achieved in models M3 and M7 where non-spherical particles PP were used. However, in models M1 and M5 where single spheres were used, the initial horizontal layers tend to become somewhat irregular even at positions far from the hopper. Also, these layers become asymmetric and exhibit abrupt fluctuations. It seems that single spheres can develop into crystalline packing configurations with dislocation zones, leading to sudden local slip failures associated with the greater velocity fluctuations observed earlier. Since all particles within the silo are still in motion, this does not mean that flow pattern is not mass flow in these models, but it indicates that when non-spherical PP particles are used, the predicted flow pattern is much closer to the ones observed in experiments. Models M2, M4, M6 and M8 were likely to show funnel flow or mixed flow according to Fig. 4a,b. In an ideal funnel flow case, the coloured layers should clearly show the presence of an internal channel along the silo height. This situation is not observed in any of the models considered in the present paper. Discrete element modelling of grain flow in a planar silo: influence of simulation parameters Fig. 5 Different snapshots of the discharge process at two discrete times after discharge (T1 = 1 s and T2 = 3 s) for all possible combinations of the variables considered (Hopper angle of inclination (β), friction coefficient between particles and walls (µp-w ) and particle shape (SP–PP)) High Friction (µ p-w=0.9) (T1 - After 1s -) (T2 – After 3s -) PP β=60º SP PP β=45º SP Low Friction (µ p-w=0.13) (T1 - After 1s -) (T2 - After 3s -) 155 However, in cases M4 and M8, were PP were used, a channel reaching the silo walls at a transition point below the solid surface can be discerned from the snapshots. This could be considered as a mixed flow mode. A detailed observation of the entire discharge process shows that there is an area close to the transition junction where the material remains at rest until the surface reaches that point. In model M2, although funnel or mixed flow is expected, the flow pattern developed is not very different from the one shown in M1. The snapshots for model M6 in Fig. 5 could be considered similar to the ones for models M4 and M8. However, a detailed observation of the whole discharge process reveals that there is no material completely at rest anywhere in the silo. Besides, the curvature exhibited by the layers in 123 156 M6 is not as smooth as the ones in models M4/M8 due to the effect of the packing configurations. This again confirms that a flow pattern closer to the one expected is better achieved by using PP instead of SP particles. 3.3 Analysis of residence times Figure 6 shows the residence times associated with the particles located at three different levels (L1, L2, L3) at the beginning of the discharge process. The vertical coordinates above hopper outlet (z) of these levels are: z = 0.55 m, z = 0.40 m and z = 0.15 m respectively. The origin of times has been chosen to coincide with the instant when the first particle in these levels exits the silo. The earlier discussion relating to Figs. 4 and 5 have provided evidence that models M5/M7 should be considered as displaying mass flow. In that case, all particles at a certain height would exit the silo earlier than those particles placed above but later than those placed beneath, thus satisfying the “first in first out” rule. Particles at the bottommost level L3 satisfied this rule for both models. However it is evident that a significant part of the particles around the central axis at level L1 was discharged before the particles adjacent to the walls at the lower level L2. This observation provides the clearest evidence that models M5/M7 did not exhibit mass flow throughout the whole discharge process. Earlier deductions have concluded that M5/M7 exhibited mass flow during the initial discharge period. At some point during the discharge process, the predominantly mass flow mode changed over to an internal funnel flow mode. This switch from mass flow to funnel flow as the solid level drops to a small height above the hopper transition is a known phenomenon and has been reported before [23]. Models M6 and M8 are expected to show funnel or mixed flow according to earlier deductions. The residence time curves for these two models has a U shape representing an abrupt change in residence time for particles at some horizontal coordinates relating to the internal channel boundaries. Within the channel, the solids are all in motion obeying the “first in first out” rule whereas outwith the internal channel, “first in last out” rule is expected to apply. This situation is seen for the case of model M8 where non-spherical PP particles were used (Fig. 6). However, the residence time curves for model M6 show more complications suggesting that simple classical description of mass/funnel flow are not adequate to describe the flow pattern. The residence time curve for level L1 shows particles at mid-distance between the centre and the wall came out later than the ones near the wall. One possible explanation is the inclined solid surface causing those particles closer to the walls to roll down this surface reaching a central position and thus exiting the silo at an early time (See M6 at T2 in 123 C. González-Montellano et al. Fig. 5). This situation may again be due to the use of spheres instead of non-spherical particles. 3.4 Pressure distributions predicted through DEM In addition to the flow pattern, the predicted pressure distributions on the silo walls have been analysed. Figure 7 shows an example of the predicted normal pressure distribution on the lateral walls (W ) of the silo for models M7 and M8 at the end of the filling process and at the beginning of discharge. The values of the filling and discharge pressures shown in Fig. 7 came, respectively, from time-averaging the individual values of pressure within the last 0.05 s of the filling process and the first 0.05 s of the discharge process, At a qualitative level, the predicted pressure distribution in the vertical section is similar to those calculated using the Jansen analytical equation [6]. Additionally, noticeable increase in pressure is predicted just below the hopper transition during initial discharge, which is the expected behaviour in this zone. The DEM computations also predicted an overall increase in normal pressure in the vertical section when the discharge is initiated, an observation which matches numerous previous studies. In addition, although a rigorous quantitative analysis has not been presented, the authors have verified that the values obtained from DEM simulations are within the order of magnitude of the values given by Janssen’s equation. According to Fig. 7, there is not a clear influence of the flow pattern in the pressure distributions, showing only clear differences near the transition, where higher pressures have been predicted for model M7. However, this observation could be due to the much reduced scale used in the models. 4 Conclusions In this paper, DEM has been used to study the flow pattern and the wall pressure in a model planar silo. The main results are summarised below: 1. The velocity profile as well as the flow pattern developed in the silo under different conditions have been analysed and compared. 2. The DEM predicted flow patterns are generally in line with the expected flow pattern from Eurocode. The study has highlighted the fact that flow pattern evolves during the discharge and that the classic description of mass flow and funnel/mixed flow may not be adequate to describe some complex situations. 3. DEM models considered in this study did not predict a classic internal funnel flow when this pattern was expected (Model M8). The reason is that with the relatively large outlet in the model, the developing flow Discrete element modelling of grain flow in a planar silo: influence of simulation parameters Fig. 6 Residence times of particles located at levels L1, L2 at L3 for models M5, M6, M7 and M8. The vertical and horizontal axis represent, respectively, the residence time (s) and initial horizontal coordinate x (m) of each particle 157 High Friction (µp-w=0.9) Low Friction (µp-w=0.13) 8 8 L1 L1 L2 6 L2 6 L3 4 4 2 2 L3 SP 0 -0 .2 5 -0 .1 5 -0.0 5 0 .0 5 0 .15 0 .2 5 0 -0 .2 5 -0 .1 5 M5 -0 .0 5 0 .05 0 .1 5 0 .2 5 0 .1 5 0 .2 5 M6 8 8 6 L1 L1 L2 L2 6 L3 4 4 2 2 L3 PP 0 -0 .2 5 -0 .1 5 -0 .0 5 M7 0 .0 5 0 .1 5 0 .2 5 0 -0 .2 5 -0 .1 5 -0 .0 5 0 .05 M8 Fig. 7 Normal pressure distributions on walls W at the end of the filling and at the beginning of discharge processes, for models M7 and M8 123 158 channel quickly propagated into a wider flow channel that reached the walls. 4. DEM using sphere appears to be problematic in capturing realistic behaviour of granular solid. Particle shape and the resulting geometric interlocking is important for predicting a flow pattern closer to the experimental observations 5. The pressure distributions predicted by DEM were qualitatively in good agreement to those expected. Further work should be done to establish whether there is also a good quantitative agreement and a fuller description of the dynamic pressure during discharge. References 1. EN 1991–4.: Actions on structures. Silos and tanks (2006) 2. Ketterhagen, W.R., Curtis, J.S., Wassgren, C.R., Hancock, B.C.: Predicting the flow mode from hoppers using the discrete element method. Powder Technol. 195(1), 1–10 (2009) 3. Rotter, J.M.: Guide for the economic design of circular metal silos. Spon, London (2001) 4. Jenike, A.W.: Bulletin 108: Gravity flow of bulk solids. Salt Lake City, UT (1961) 5. Jenike, A.W.: Bulletin 123: Storage and flow of solids. Salt Lake City, UT (1964) 6. Janssen, H.A.: Versuch über Getreidedruck in Silozallen. Z Des Vereines Deutscher Ingenieure 39(35), 1045–1049 (1895) 7. Reimbert, A, Reimbert, M.: Silos théorie et pratique. Editions Eyrolles, France (1959) 8. Tejchman, J., Klisinski, M.: FE-studies on rapid flow of bulk solids in silos. Granul. Matter 3(4), 215–230 (2001) 9. Vidal, P., Gallego, E., Guaita, M., Ayuga, F.: Finite element analysis under different boundary conditions of the filling of cylindrical steel silos having an eccentric hopper. J. Constr. Steel Res. 64(4), 480–492 (2008) 123 C. González-Montellano et al. 10. Vidal, P., Guaita, M., Ayuga, F.: Discharge from cylindrical slender steel silos: finite element simulation and comparison with Eurocode 1. T. ASABE. 48(6), 2315–2321 (2005) 11. Chen, J.F., Yu, S.K., Ooi, J.Y., Rotter, J.M.: Finite element modelling of filling pressures in a full–scale silo. J. Eng. MechASCE 127(10), 1058–1066 (2001) 12. Ooi, J.Y., Pham, L., Rotter, J.M.: Sistematic and random features of measured pressures on full-scale silo walls. Eng. Struct. 12, 74–87 (1990) 13. Ramirez, A., Nielsen, J., Ayuga, F.: On the use of plate-type normal pressure cells in silos. Part 2: validation for pressure measurements. Comput. Electron. Agr. 71(1), 64–70 (2010) 14. Härtl, J., Ooi, J.: Experiments and simulations of direct shear tests: porosity, contact friction and bulk friction. Granul. Matter 10(4), 263–271 (2008) 15. Chung, Y.C., Ooi, J.Y.: Influence of discrete element model parameters on bulk behavior of a granular solid under confined compression. Particul. Sci. Technol. 26(1), 83–96 (2008) 16. Cundal, P., Strack, O.: A discrete numerical model for granular assemblies. Geotechnique 29(2), 47–65 (1979) 17. Tsuji, Y., Tanaka, T., Ishida, T.: Lagrangian numerical simulation of plug flow of cohesionless particles in a horizontal pipe. Powder Technol. 71(3), 239–250 (1992) 18. Anand, A., Curtis, J.S., Wassgren, C.R., Hancock, B.C., Ketterhagen, W.R.: Predicting discharge dynamics from a rectangular hopper using the discrete element method (DEM). Chem. Eng. Sci. 63(24), 5821–5830 (2008) 19. Foerster, S.F., Louge, M.Y., Chang, H., Allia, K.: Measurements of the collision properties of small spheres. Phys. Fluids 6(3), 1108 (1994) 20. Lorenz, A., Tuozzolo, C., Louge, M.: Measurements of impact properties of small, nearly spherical particles. Exp. Mech. 37(3), 292–298 (1997) 21. Johanson, J.R., Jenike, A.W.: Bulletin 116: Stress and velocity fields in gravity flow of bulk solids. Salt Lake City, UT (1962) 22. Härtl, J., Ooi, J., Rotter, J., Wojcik, M., Ding, S., Enstad, G.: The influence of a cone-in-cone insert on flow pattern and wall pressure in a full-scale silo. Chem. Eng. Res. Des. 86(4), 370–378 (2008) 23. Zhong, Z., Ooi, J.Y., Rotter, J.M.: The sensitivity of silo flow and wall pressures to filling method. Eng. Struct. 23, 756–767 (2001)