In this example we show how to implement a simulation Discrete element simulation (DEM) Using the Lorenz-force contact model. We will also use dynamic load balancing to keep the simulation balanced during all simulation time. The simulation result is show in the figures and video below
A classical model for DEM simulations of spherical granular flows is the Silbert model, it includes a Herzian contact force and an elastic deformation of the grains. Each particles has radius \(R\), mass \(m\), polar momentum \(I\) and is represented by the location of its center of mass \(r_{i}\). When two particles i and j collide or are in contact, the elastic contact deformation is given by: \(\delta_{ij} = 2R-r_{ij}\)
where \(\vec{r}_{ij}=\vec{r}_i - \vec{r}_j\) is the distance vector connecting particle centers and \(r_{ij} = \|\vec{r}_{ij}\|_2\) its module. The normal and tangential components of the relative velocity at the point of contact is given by
whith \(\vec{n}_{ij}=\vec{r}_{ij}/r_{ij}\) is the normal unit vector in direction of the distance vector, \(\vec{\omega}_i\) is the angular velocity of a particle and \(\vec{v}_{ij}=\vec{v}_i-\vec{v}_j\) the relative velocity between the two particles. The evolution of the elastic tangential displacement \(\vec{u}_{t_{ij}}\) is integrated when two particles are in contact using.
Where $\delta t$ is the time step size. The deformation of the contacts points is stored for each particle and for each new contact point the elastic tangential displacement is initialized with \(\vec{u}_{t_{ij}} = 0\). Thus for each pair of particle interacting the normal and tangential forces become:
where \(k_{n,t}\) are the elastic constants in normal and tangential direction, respectively, and \(\gamma _{n,t}\) the corresponding viscoelastic constants. The effective collision mass is given by \(m_{\text{eff}}=\frac{m}{2}\). For each contact point in order to enforce Coulomb's law \(\|\vec{F}_{t_{ij}}\|_2 < \|
\mu\vec{F}_{n_{ij}}\|_2\) the tangential force is bounded by the normal component force. In particular the elastic tangential displacement \(\vec{u}_{t_{ij}}\) is adjusted with \(\vec{F}_{t_{ij}} \leftarrow
\vec{F}_{t_{ij}}\frac{\|\mu\vec{F}_{n_{ij}}\|_2}{\|\vec{F}_{t_{ij}}\|_2}
\)
This adjustment induce a truncation of the elastic displacement. The Coulomb condition is equivalent to the case where two spheres slip against each other without inducing additional deformations. Thus the deformation is truncated using:
Considering that each particle \(i\) interact with all the particles $j$ is in touch with , the total resultant force on particle $i$ is then computed by summing the contributions of all pair particles $(i,j)$. Considering that the grains are also under the effect of the gravitational field we obtain that the total force is given by
where \(\vec{g}\) is the acceleration due to gravity. Because particles has also rotational degree of freedoms, the total torque on particle $i$ is calculated using
\(\vec{r}_i\) and angular velocities \(\vec{\omega}_i\) for each particle \(i\) at time step \(n+1\), We integrate in time the equations using leap frog scheme with time step given by $\delta t = 10^{-6}\,$s
where \(\vec{r}_i^{n},\vec{v}_i^{n},\vec{\omega}_i^{n}\) denotes respectively the position, the speed and the rotational speed of the particle $i$ at time step \(n\), and \(\delta t\) the time step size.
We simulate an avalanche down an inclined plane
Constants
This simulation require several constants to work, here we define mass, radius, delta time elastic constants and stop time, simulation domain and ghost.