Open access peer-reviewed chapter

Discrete Element Modeling of a Projectile Impacting and Penetrating into Granular Systems

Written By

Waseem Ghazi Alshanti

Submitted: October 12th, 2017 Reviewed: February 16th, 2018 Published: November 5th, 2018

DOI: 10.5772/intechopen.75550

Chapter metrics overview

1,166 Chapter Downloads

View Full Metrics


From theoretical standpoint, it is difficult to analytically build a general theory and physical principles that critically describe the mechanical behaviour of granular systems. There are many substantial gaps in understanding the mechanical principles that govern these particulate systems. In this chapter, based on a two-dimensional soft particle discrete element method (DEM), a numerical approach is developed to investigate the vertical penetration of a non-rotating and rotating projectile into a granular system. The model outcomes reveal that there is a linear proportion between the projectile’s impact velocity and its penetration downward displacement. Moreover, depending on the rotation direction, there is a significant deviation of the x-coordinate of the final stopping point of a rotating projectile from that of its original impact point. For negative angular velocities, a deviation to the right occurs while a left deviation has been recorded for positive angular velocities.


  • discrete element method
  • granular bed
  • rotating projectile
  • penetration
  • vertical impact

1. Introduction

A granular material is any collection of many macroscopic discrete solids, whose typical size ranges from micrometres to centimetres such as sand, coal, sugar, and rice. It is obvious that these materials cannot be characterised as gas, fluid, or solid. But, they have hybrid state between these three phases depending on the average energy of the individual grains inside the granular system. Piles, snow avalanches, and planetary rings (interstellar dust) are, respectively, considered as solid, fluid, and gas phases for granular materials [1].

In 1895, the German engineer H. A. Janssen [2] made the first attempt to study the dynamical behaviour of the granular systems under specific conditions. He investigated the pressure propagation along the wall of a silo filled with granular solids. Nevertheless, the science of granular materials has a standing research history that has deep roots in time. A general theory that governs the mechanical behaviour of the granular materials has not yet been established due to the complex behaviours of these particulate systems. Most of the successful theories in physics are based on forming a linear differential equation that contains physical quantities from the system under study. But, the knowledge needed to learn about the dynamics of granular materials requires the construction of nonlinear differential equation that can be solved only in rare cases. Therefore, efforts have mainly been made over the past 50 years to build theories that are based on linear differential equations to study these particulate systems. The absence of comprehensive physical principles of granular materials makes the physics of granular materials still not completely understood.

It is well known that the experiment studies of granular materials are limited, and are not capable of capturing most of the internal features of the granular dynamical systems governed by cohesion, friction, and collisions of the particles. Thus, the exploit of such investigations in granular material studies become ineffective and unreliable in most cases. Accordingly, modellers start to employ new techniques and approaches in their investigations relevant to the mechanics of granular materials. Away from experimental and physical approaches, new techniques have been proposed to investigate the behaviour of granular systems. These new techniques can be classified into three main categories:

  • Theoretical-based techniques such as continuum mechanics methods and micromechanics analysis of particle collisions.

  • Statistical averaging-based techniques such as the dense gas kinetic theory methods.

  • Molecular dynamic simulation-based techniques such as the soft and hard particle discrete element methods (DEM), Monte Carlo techniques, and cellular automata techniques.

With the advance in computer technology, the molecular dynamic simulations or discrete approaches have recently been emerged. These methods have gained popularity and have been considered as powerful tools capable of handling the demands of the granular system research. In this chapter, based on the two-dimensional soft particle discrete element method (DEM), numerical simulations are carried out to study the dynamics of projectiles impacting and penetrating, vertically, beds of granular materials. First, we investigate the relationship between the penetration displacement of a projectile underneath a granular bed and its impact velocity. Then we investigate the effect of the rotation direction of a spinning projectile on its trajectory beneath the surface of a granular bed.


2. Development of discrete element model

Discrete element methods (DEMs) are numerical techniques capable of simulating the entire behaviour of systems of discrete interacting elements. Examples of discrete element simulations include hard and soft particle methods. The most common and flexible type of discrete element simulation methods is the soft particles models. They are able to model any kind of configuration, including static and dynamic situations [3]. Diluted and dense as well as long-lasting granular systems can be modelled through soft particle methods. The first soft particle contact model was carried out by Cundall and Strack [4]. The term soft particle refers to the fact that the particles are assumed to have collisions of (small) limited time and possibly experience elastic deformation. During the finite contact, the geometrical shape of the particle remains rigid, and the deformation is represented as small overlap at the surface which is taken into the account in the force model. The values of forces at the contact points are always varying as the particles are being deformed. At the same time, a particle may undergo multiple contacts with neighbouring particles and hence, numerous acting contact forces. Many types of forces can be simulated to act on a particle and various contact force models such as elastic and viscoelastic can be incorporated in this case. The general outline of any soft particle contact model is to determine, at one instant, each individual particle interchanges. Those interchanges are the surface (normal and frictional contact) forces with neighbouring particles or with the system boundaries. The normal and frictional interaction forces between particle-particle and particle-boundary are modelled by system of linear springs, shock absorbers (Dashpots), and sometime sliding elements are used. Due to the influence of these interaction forces and other forces like gravitational body force, the particles change their positions and velocities continuously throughout the simulation time. For each particle, the motion is governed by the principle of linear momentum and the principle of angular momentum. The resultant force acting on any individual particle can be calculated. Therefore, the acceleration of each individual particle at any simulation time step can be determined by using Newton’s second law. Finally, the new velocity and position of the particle can be obtained by integrating the differential equation of its motion over small simulation step time.

2.1. Mathematical formulation set-up

For a typical simulation, each granular particle within the granular system bears two types of forces: contact forces and gravitational body force. Any contact force between two particles is decomposed into normal and tangential components. The normal contact force is modelled by a damped linear spring, while the tangential contact force by a linear spring in series with a frictional sliding element. The formula that determines the contact force of particle i and particle j is



are unit vectors in the normal and shear directions of the contact plane,
are, respectively, the magnitudes of the normal contact force and shear contact force, namely,



are, respectively, the particle-particle normal and tangential spring coefficients,
are, respectively, the elastic contribution of the contact force between the particles i and j in the normal direction (
direction) and the friction coefficient of the granular particles,
are, respectively, the normal compression and the tangential displacement between the particles i and j over the time step
are the radii of the particles i and j . Under the contact forces and the gravitational body force, each particle has the following motion dynamic equations



, and
are, respectively, the mass, rotational moment of inertia, position, and rotational vectors of the centre of particle i ,
, are, respectively, contact force and moment acting on particle i due to particle j and external forces and moment acting on particle
is the number of particles within the granular bed. Hence,
we have a system of first-order ordinary differential equations as follows:


Therefore, by numerical integration of Newton’s equation of motion, the updated velocities and positions of all particles can be determined.

2.2. Numerical simulation structure

In general, numerical simulations are performed in two environments, namely, two dimensions and three dimensions. Most of the current numerical simulations are of two dimensions. They are popular due to the low computational cost by comparison with three-dimensional simulations. The behaviour of the entire granular system is also well observable. Furthermore, most of the granular dynamical problems require the employ of two-dimensional simulations rather than three-dimensional simulations. For instance, two-dimensional simulations work well for problems such as the motion of a projectile penetrating into granular beds. Furthermore, in three-dimensional simulations, more contact points and extra spatial degrees of freedom have to be considered for each individual particle inside the system, and thus more computational power and memory storage are required.

When developing a two-dimensional discrete element code to imitate a particular granular problem, several aspects have to be taken in consideration. For instance, penetration of vertical projectile into granular beds is a granular dynamical problem in which most of particles in the system experience considerable multiple, long duration contacts. In this case, sufficient number of particles, in atypical simulation, needs to be taken into account, and the collisions and particles’ overlaps should be easy to detect.

The current two-dimensional discrete element simulation code assumes the particles as a set of circular disks which are placed in the workplace at predetermined positions and velocities. The boundary conditions are assumed to be frictional. In the simulation workplace, any particle can be characterised through its radius R i , mass m i , initial moment I i , position vector of the particle centre r i , translational (linear) velocity r ̇ i , and the rotational (angular) velocity θ ̇ i , for i = 1 , 2 , 3 , , N p , where N p is the total number of particles in the simulation. During the simulation, three types of forces are considered to act on the particles. These forces are gravitational body force and particle-particle and particle-boundary contact forces. A proper simulation time step was chosen so that short contacts are not missed. The interstitial fluid, cohesive, electromagnetic forces were neglected. Therefore, this model simulates granular material as dry particles in a blanked space and gravity field.


3. Vertical penetration displacement of a projectile and its impact velocity

In order to describe a relationship between the vertical penetration displacement of a projectile and its impact velocity, a series of numerical simulations were conducted. The methodology was to vary the impact velocity of the projectile and keep the other simulation parameters alike. The projectile impact velocities range from 5 to 65 m/s. The particle bed was considered to be mono-sized with particle diameter equal to 1.5 mm . Figure 1 shows simulation sequential snapshots of the vertical penetration process. During the entire process, one can identify three distinct regimes, namely, impact, penetration, and collapse.

Figure 1.

Simulation snapshots show the penetration process.

Figure 2(a) shows the relationship between the impact velocities and the penetration distances underneath the granular bed’s surface for three different values of projectile’s diameter. The observations reveal that the projectile’s penetration downward displacement d pen increases linearly with the projectile’s impact velocity v imp .

Figure 2.

(a) Impact velocities versus penetration distances and (b) log-log plots of impact velocities versus penetration distances.

To extrapolate a power law that governs the relationship between the projectile impact velocity and its penetration displacement, we convert the obtained data to log-log plots, as shown in Figure 2(b). Obviously, the relationship between the penetration displacement and the impact velocity follows a power law of the form

d penetration v impact 1 2 . E10

The results from the simulation compare well with previous experimental results such as [5].


4. Vertical penetration of rotating projectiles into granular beds

To analyse the relationship between the angular velocity of a rotating projectile on its orthogonal penetration displacement, a set of numerical simulations was carried out with different aspects. The method is to equip the rotating projectile with various angular impact velocities and keep all other simulation parameters identical including its impact velocity (

). Figure 3 shows snapshots of simulation of a rotating projectile orthogonally penetrating mono-sized granular bed with

Figure 3.

Simulation sequential snapshots of normal penetration of a rotating projectile.

For series of numerical simulations, a rotating projectile with fixed impact vertical velocity

is given various angular velocities, namely,
. Then the projectile strikes the middle of the free surface of the mono-sized granular bed. Figure 4 shows the resultant trajectories of the rotating projectile underneath the granular bed’s surface for each different value of the angular velocities. In general, the trajectory profile of the projectile exhibits two different regimes. For negative angular velocities, namely,
the penetration trajectories of the projectile under the granular bed exhibit negative exponential like traces. Conversely, positive exponential-like traces have been recorded for positive angular velocities, namely,
The case when the projectile has no angular velocity, i.e.
can be considered as a turning value angular velocity between the two regimes. The resultant projectile’s trajectory, in this case, locates in the middle of the two regimes’ traces and the ultimate penetration point of the projectile is located, approximately, under its impact point. Moreover, it is found that, when the projectile comes to rest after achieving its maximum penetration depth, there is a considerable deviation for the x-coordinate of its ultimate stopping point from that of its original impact point. For negative angular velocities, namely,
a deviation to the right under the original impact point occurs while a left deviation happens for positive angular velocities, namely,

Figure 4.

Deviation of the projectile?s x-coordinate ultimate stopping point from that of its original impact point for various values of angular velocities. (a) Multisized particle bed. (b) Monosized particle bed.


5. Conclusion

This chapter focuses on the development of mathematical model based on the soft particle discrete element method for the study of two granular dynamical problems:

  1. The relation between the vertical penetration displacement of projectile and its impact velocity. It has been found that the scaling law of the penetration displacement of a projectile with its impact velocity follows a power law of the form

d penetration v impact 1 2 . E11
  1. We investigate the trajectories of a rotating projectile penetrating normally a granular system. Our numerical results show that the model is capable to simulate the normal penetration process for the various values of angular velocities. Moreover, it is found that the trajectory profile of the rotating projectile underneath the granular bed is affected the magnitude as well as the direction of its angular velocity. Depending on the rotation direction of the projectile, there is a relatively small change on the ultimate horizontal position of the projectile after penetration. For negative angular velocities, a right shift on the ultimate penetration point from the original projectile’s impact point is observed. On the other hand, the projectile’s eventual penetration point is located to the left of its original impact point for the positive angular velocities.


  1. 1. Goldhirsch I. Some open problems in granular media. Department of Fluid Mechanics and Heat Transfer. Faculty of Engineering. Tel-Aviv University. In: Subprogram on Uncertainty in Models of Granular Materials: Sources and Consequences; Kickoff Workshop; September 10-14, 2006
  2. 2. Janssen HA. Versuche über getreidedruck in silozellen (investigation on the pressure of grain in silos). Zeitschrift des Vereines Deutscher Ingenieure. 1895;(39):1045
  3. 3. Babic M, Shen HH, Shen HT. The stress tensor in granular shear flows of uniform, deformable disks at high solids concentrations. Journal of Fluid Mechanics. 1990;219:81-118
  4. 4. Cundall PA, Strack DL. A discrete numerical model for granular assemblies. Geotechnique. 1979;29:47-65
  5. 5. de Bruyn JR, Walsh AM. Penetration of spheres into loose granular media. Canadian Journal of Physics. 2004;82:439-446

Written By

Waseem Ghazi Alshanti

Submitted: October 12th, 2017 Reviewed: February 16th, 2018 Published: November 5th, 2018