DEM - gnomeCreative/HYBIRD GitHub Wiki
Outline of the DEM model
What follows is adapted from Leonardi (2015). DEM is a method that has become over the last years the most common choice for the representation of granular systems (Bicanic, 2007). Every grain $\textrm{p}$ is represented as a Lagrangian point $x_\textrm{p}$ and is idealized as a particle with mass $m_\textrm{p}$ and moment of inertia $J_\textrm{p}$. In a three-dimensional space, the state of a particle is defined by six degrees of freedom, position $x_\textrm{p}$ and orientation $\theta_\textrm{p}$. The evolution of these is governed by Newton's equations of motion,
$m_\textrm{p}\frac{d^2 x_\textrm{p}}{dt^2} = F_\textrm{p}$,
$J_\textrm{p}\frac{d^2 \theta_\textrm{p}}{dt^2} = M_\textrm{p} - \frac{d \theta_\textrm{p}}{dt} \times \left(J_\textrm{p}\frac{d \theta_\textrm{p}}{dt} \right) $
where the force $F_\textrm{p}$ and the moment $M_\textrm{p}$ are the results of the interactions acting on the particle. These are the fluid-particle coupling interactions $F_\textrm{hydro}$ and $M_\textrm{hydro}$ coming from LES, the collisions between particles, and external force fields such as gravity. In general, the interaction resultants are a function not only of position and orientation of the particles, but also of their linear velocity $u_\textrm{p}$ and angular velocity $\omega_\textrm{p}$:
$F_\textrm{p}=F_\textrm{p}\left(x_\textrm{p},u_\textrm{p},\theta_\textrm{p},\omega_\textrm{p}\right)$,
$M_\textrm{p}=M_\textrm{p}\left(x_\textrm{p},u_\textrm{p},\theta_\textrm{p},\omega_\textrm{p}\right)$.
Representation of the DEM particles and of the contact law.
The simplest form of the DEM considers particles to be spheres of radius $r_\textrm{p}$. In this case, checking if two particles $\textrm{p1}$ and $\textrm{p2}$ are in contact is very simple, the overlap between them being
$\xi=r_\textrm{p1}+r_\textrm{p2}-d_{\textrm{p1},\textrm{p2}}$
where $d_{\textrm{p1},\textrm{p2}}=x_\textrm{p2}-x_\textrm{p1}$ is the vector distance between the two centers, which also defines the collision reference system (see figure above), whose first component is
$n=\frac{d_{\textrm{p1},\textrm{p2}}}{d_{\textrm{p1},\textrm{p2}}}$.
The second and third components, $t$ and $b$, are defined according to the component collision velocity vector $u_\textrm{coll}$. The component of $u_\textrm{coll}$ parallel to $n$, the normal collision velocity, can be computed as a function of the translational velocities of the particles $u_\textrm{p1}$ and $u_\textrm{p2}$
$u_\textrm{n,coll}=((u_\textrm{p2}-u_\textrm{p1})\cdot n)n$,
while the component aligned to $t$, the tangential velocity, is computed as a function of the rotational velocities of the particles $\omega_\textrm{p1}$ and $\omega_\textrm{p2}$ too, as
$u_\textrm{t,coll}=u_\textrm{p2}-u_\textrm{p1}-u_\textrm{n,coll}-r_\textrm{p1}\omega_\textrm{p1}\times n-r_\textrm{p2}\omega_\textrm{p2}\times n$.
From the definition of the relative tangential velocity, the remaining two vectors of the collision reference system are defined as
$t=u_\textrm{t,coll}/|u_\textrm{t,coll}|$,
$b=n\times t$.
A 2-sphere and a 3-sphere cluster.
The overlap represents the elastic deformation of the particles. If $\xi$ is positive, the two particles exchange a repulsive force $F_\textrm{coll}$, whose expression as a function of $\xi$ is defined according to the contact constitutive model. This representation of debris material as spherical particles is motivated by computational convenience, but is only acceptable as a first approximation of their actual behavior. One of the possibility to overcome this simplification is to represent grains as clusters of spheres, rigidly connected (see figure above). This preserves the simplicity of the contact identification function, but allows for a more realistic representation.
The collision force $F_\textrm{coll}$ is computed by calculation of the moduli of its normal component $F_\textrm{n,coll}$ and tangential component $F_\textrm{t,coll}$, which are both dependent on the employed constitutive model. The final system of interactions during a collision is
$F_\textrm{coll,p1}=-F_\textrm{n,coll}\cdot n+F_\textrm{t,coll}\cdot t$;
$F_\textrm{coll,p2}=F_\textrm{n,coll}\cdot n-F_\textrm{t,coll}\cdot t$;
$M_\textrm{coll,p1}=r_\textrm{p1} F_\textrm{t,coll}\cdot b$;
$M_\textrm{coll,p2}=r_\textrm{p2} F_\textrm{t,coll}\cdot b$.
If the particles are spherical, only the tangential component generates a moment. If the particles are non-spherical, Eq.~\ref{eq:allInteractions} must be adapted to take into account the secondary moment generated by the normal component.
The contact models are described in the following pages: