DEM Fundamentals

Introduction

The discrete-element method was first described by \cite{cund79} as an application to granular assemblies. DEM's have been used to investigate the properties of granular shear zones \cite[]{morg99a,morg99b}, the propogation of blind thrusting has been tackled in 2D by \cite{finc03} and in 3D by \cite{stra02}. A numerical simulation of seismic waves using a discrete particle scheme was performed by \cite{toom00} with the aim of developing a tool to investigate wave propogation in multiphase media and they have also proven particularly useful in studies of granular materials where the bulk properties emerge as a result of individual contact interactions and their histories at a finite scale. Recent work has produced FEM-DEM coupled models to study the interaction of particles suspended in a gas or liquid flow \cite[]{limt03}.

Elements can represent either discrete objects, such as individual grains of sand, or bulk materials. In the experiments presented in this study, each element represents some mass of rock. The particle sizes define the resolution of structures within the model. They will allow a full exploration of parameter space appropriate to geological and sandbox problems and have the potential to self organise structures.

Generally, contact forces are determined by some function of the separation of the discs in the normal and tangential directions. The collision of two balls produces non-zero forces only when particles overlap, where as an electrostatic repulsion would work over a larger range.

The forces exerted on individual particles are summed in a Cartesian frame of reference and integrated, as outlined in Fig.~\ref{fig:general_scheme}.

Disc contact variables

The key particle variables required for determining the repulsive force due to particle overlap are shown in the figure below.
Regular hexagonal packing

Contact force variables

Efficient contact seach scheme

Computationally, the slowest process in running a DEM is contact detection. If we were to compare every particle with every other particle, the speed of the program would decrease approximately with the square of the number of particles. To reduce this, we can linearly sort the particles into a cell structure (Fig. \ref{fig:Cell_structure}) by their $x$ and $y$ coordinates (Fig. \ref{fig:Sort_particle}) \cite[]{comsimliq}, storing this information in a Linked List. Then we only need compare particles in adjacent cells. Forces are equal and opposite, so we need only compare Cell~A with Cell~B and not the reverse (Fig. \ref{fig:eq_and_op}).

Random clustered disc

Box sorting for efficient contact seaching

Random clustered disc Random clustered disc

Cells need only be compared once as forces are equal and opposite

The size of the cell to attain the best performance is a balancing act between small cells which contain few particles but a short particle travel distance before the list needs to be updated and a large cell which contains more particles but does not need to be updated as frequently. This scheme works well for a packing involving discs of similar sizes. The size of the cell must be larger than the largest particle so a constant cell size may be inconvenient for particles with a large range of sizes, a dynamic cell size can then be useful.

Boundary conditions

A number of DEM experiments have been investigated as part of this study. This section describes the terminology used in describing the experimental setups.

Rigid walls interact with particles using the same repulsive force law as the particle-particle interactions (\ref{eqn:force}), where the magnitude of the force is proportional to the overlap with the boundary (\ref{eqn:wall_overlap}), the boundary is undeformable. A \emph{smooth} wall imparts no tangential force at its contacts. A \emph{rough} surface allows a tangential component which resists motion parallel to the surface through a frictional force law, in this study a Coulomb failure criteria has been applied.

A no-slip boundary is formed where the particles in contact with that surface are fixed in space, relative to some frame of reference. Thus it can form a constant velocity condition along a surface or provide a high resistance to motion. This is characteristically different to a Coulomb failure condition, which does allow slip on the boundary.

For some experiments, it is neccessary to completely remove particles from the system. This may be due to erosion of a surface or the flow of particles out of the system. So that many arrays need not be re-allocated, the particles are stored outside the model space with their radii set to zero so that they cannot interact and their force components perminantly set to zero, so they cannot move. It is useful to store them as a function of their particle number, so that they are never stored on top of each other. Similarly, this technique can be used as a reservoir of particles to be added to the system. This type of boundary will be refered to as an \emph{open} boundary, other boundaries we have described are therefore \emph{closed}.

Particle packing structures

The inital packing arrangement used in the experiments is generated in two stages. Firstly, a loose distribution of non-dimensional particles is produced; some typical particle arrangements are described below. Secondly, the packing is allowed to settle under gravity given a particular parameterization, specific to the experiment to be performed. This outputs a settled particle arrangement containing information on particle position, force balance, velocities and contact geometry. This is the inital particle data read into each experiment.

Regular hexagonal packing

Regular hexagonal packing

This regular packing structure Fig.~\ref{fig:pack1}, using balls of uniform size, was employed to study accretionary wedges and fold-and-thrust belts \cite[]{burb02}. This packing has pre-defined planes of weakness to shearing and artificial strength along long linear chains, which explains the emergence of pop-ups at a long distance from the toe of the wedge. This packing is also the most compact, and thus any disruption would weaken it. It was argued that the anisotropy in the layering would better simulate sedimentary packings, which are often layered, and facilitate strain weakening. We believe that the more fundamental behaviour is of a homogeneous material with no predefined planes of weakness or fault geometries should first be investigated first, before any anisotropy is introduced, to understand which properties result from the boundary conditions and which from the packing arrangement. Subsequent experiments can then relate the structural evolution to this reference. Additionally, it would be useful if the degree of strain weakening along fault planes could be definable through some inital parameter; this cannot be achieved if strain weakening is employed as a function of the packing.

Random single disc

Random single disc

To address some of the issues outlined above, a random medium was generated as outlined by Mora. This material consists of balls with a random size, between limits, positioned randomly (Fig. \ref{fig:pack2}). We found that this material immediately improved critical wedge-type growth over the regular hexagonal packing. A critical wedge was seen to form at the moving backstop with oscillations in the localisation of deformation dependent upon the basal friction. However the angle of the wedge was low due to balls frequently rolling and avalanching down the surface due to over-rotation. Various methods have been identified to overcome this problem, which mainly focus on increasing particle angularity. This is done through more complex particle geometries including ovular and polygonal particles, these more complex shapes increase the time for contact detection which rapidly becomes more complicated for non-circular particles.

Random clusters

Random clustered disc This study applied an alternative approach, proposed by \cite{jens99}, to tackling the problem of over-rotation. Particles consisting of clusters of 3 random sized discs are used to increase angularity (Fig. \ref{fig:pack3}). These particles will remain intact for all time by redefining the contact force for clustered particles such that,

This scheme is a simple extension of the standard contact force law, the expense is in increasing the degrees of freedom for each particle by the number of balls in each particle. \cite{jens99} used three regular sized particles to form an equilateral triangle; in our study the isotropic nature has been increased by using particles made from three random sized balls, between limits, linked together in a triangular formation. This increased the material shear strength and reduced the particle rotation.