Collision models

There are several approaches to treat collisions in particle simulations. However, there is the assumption that in finite-time particle-particle and particle-wall collisions in a continuous setting do not take place. This results out of repulsive and lubrications forces preventing collisions in terms of sufficiently viscous fluids (Glowinski et al., 2001; Glowinski, 2003). Nevertheless, collisions have to be modeled for numerical simulations to avoid decreasing the time step size, as well as the element size, if particles come too close. This would lead to a massive increase in computational costs as soon as particle collisions occur in the simulation.

For common collision models, one can distinguish repulsive force, lubrication force and models based on momentum conservation. The first two add a repulsive force to the right-hand side of the NEE (2.7). The latter is based on the equations of motion and a balance of total momentum at the point of contact. Here, the models which are described can be extended to particles with arbitrary shape quite easy. However, the extension shifts the computational effort to determining the distance between both particles efficiently, which is more expensive than for, e.g. spherical particles.

Repulsive force

The following collision model was developed by Glowinski et al. (2001) and is based on repulsive forces between colliding particles. It is assumed, that collisions have a smooth nature, meaning their particle velocities coincide at the point of contact.

Additionally, it has to be ensured that the particles do not overlap each other at any time in the simulation. For circular particles the repulsive force~_Frep has to fulfill the following properties

~_Frepparallel to−−→_X

iX_j, (3.1a)

||~_Frep||=0 ifd_ij ≥^Ri+R_j+ρ, (3.1b)

||~_Frep||= ^cij

e ifd_ij = R_i+R_j. (3.1c) In (3.1), ~_Xi is the center of mass particlei, R_i is the radius and d_ij = ||−−→_X

iX_j||^{is the} distance between both center of masses. c_ij is a scaling factor ande a small positive number depending on grid size. Additionally to the properties of (3.1), the model has to satisfy the behavior of Figure 3.1. This means a decreasing behavior between maximum and minimum value of repulsive force. In Figure 3.1,ρis the range where the repulsive force model is acting and is chosen to be∆h, the minimum mesh size for the spatial discretization of the velocity in the original literature. Further information on how to choose this parameter can be found in Glowinski et al. (2001).

||~_Frep||

d_ij R_i+R_j+ρ R_i+Rj

Figure 3.1: Force depending on distance in the model of Glowinski et al. (2001).

The aforementioned method does not allow the particles to overlap each other. How-ever, overlapping can occur in numerical calculations, if the time step size is not adapted when particles come close to each other. As a result, Wan and Turek (2006) proposed a modified definition of the short-range repulsive force model where particles

can come arbitrary close and are even allowed to overlap. A schematic representation of the repulsive force with overlapping region can be seen in Figure 3.2.

~_Frep =0 ifd_ij ≥^Ri+R_j+ρ (3.2a)

~_Frep = ¹

e_p(~_Xi−~_Xj)(R_i+R_j+ρ−^dij)² = ifd_ij ≥ ^Ri+R_j (3.2b)

~_Frep = ¹

e⁰_p(~_Xi−~_Xj)(R_i+R_j−^dij)ifd_ij ≤ ^Ri+R_j (3.2c) In (3.2)ρis again the range of the repulsive force and can be chosen to be a factor of 0.5 to 2.5 of the meshsize∆h. Further,e_pande⁰_pare small stiffness parameters which are chosen to bee_p ≈(∆h)²and e⁰_p ≈^∆hif the ratio of both densities is around 1 and the fluid is sufficiently viscous.

||~_Frep||

d_ij R_i+R_j+ρ

R_i+Rj 0

ρ²

ep||~_Xi−~_Xj||

R_i+R_j

e⁰_p ||~_Xi−~_Xj||

overlapping region

Figure 3.2: Force depending on distance in the model of Wan and Turek (2006), here overlap is possible.

The model can also be used for particle-wall collisions with some small modifications.

In (3.3)d_i⁰ =|~_Xi−~_X⁰_i|is the distance between the center of mass and the center of an imaginary particle at the wall. The stiffness parameters for particle-wall collisions are chosen to be half of the particle-particle parameters, eg. e_w =e_p/2 ande⁰_w =e⁰_p/2.

~_Frep =0 ifd⁰_ij >2R_i+ρ (3.3a)

~_Frep= ¹

e_w(~_Xi−~_Xi⁰)(2R_i+ρ−^d0i)² = if 2R_i ≥^d0i ≥^Ri+R_j+ρ (3.3b)

~_Frep= ¹

e⁰_p(~_Xi−~_Xi⁰)(2R_i−^d0_i)ifd_i ≤^2Ri (3.3c)

Lubrication force

If particles come very close to each other, a Poiseuille-type flow develops between them, leading to high local stress and pressure values. Here, the magnitude of the lubrication force strongly depends on the distance between the particles.

A very popular example of lubrication force models is the one of Maury (1997). Let

~C_iand_C~_j be the material points on the surface of two particles which are closest to each other. Further, ˙~_{Ciand ˙C}~_jbeing the velocity at those points and~n_iand~_{tibeing the} normal and tangential vector on the particle surface∂Ω_i, resp. ∂Ω_j. All quantities are visualized in Figure 3.3. The lubrication force can be estimated by

~_Flub = [−^κn(d_ij)~n_c⊗~n_c−^κt(d_ij)~_tc⊗~_tc]·(_C~^˙_i−~_C^˙_j), (3.4) where

κ_n(d) = µ_n1

d andκ_t(d) = µ_t ln(^d0

d ). (3.5)

In (3.5),µ_n andµ_t depend on geometrical aspects like curveture of the particle surface and on the viscosity of the fluid. Further,d0will be chosen according to the character-istic size of the particles. In his publication Maury (1997) assumesκ_n(d)andκ_t(d)to be given funnctions. It is also important that, like in the repulsive force models, a force will only act ifdis greater than a given valued0.

d_ij C_j C_i

Ω_j

∂Ω_j Ω_i

∂Ω_i

Figure 3.3: Particle quantities used by Maury (1997).

Conservation of momentum

Using the basics of rigid particle dynamics from Section 2.4, collision models based on the conservation of momentum can be constructed. Those models do not depend on an a-priori stiffness parameter. Ardekani and Rangel (2008) published a method based on the conservation of momentum for circular shaped particles with a Lagrange-Multiplier method. They consider only centric collisions, where both collision forces

act in the direction of the center of masses. This results in no torque acting on the rigid particle.

For this method, the velocities at the points of contact are split into normal and tangential parts of the particle surfaces_U~ =U_n~n_c+U_t~_tt. For tangential components of the velocitiesU_t1andU_t2momentum is conserved as well as the angular momentum for each particle separately. For both normal components the momentum balance renders down to

U_n1⁺ = ^e(U_n2⁻ −^U_n1⁻)M2

M1+M2 + ^M1Un1+M2U_n2

M1+M2 (3.6a)

U_n2⁺ = ^e(U_n2⁻ −^U_n1⁻)M₁

M1+M2 + ^M1Un1+M2U_n2

M1+M2 (3.6b)

Here, the superscripts−^/+denote quantities right before and after the collision and Mis the mass of the rigid particle. eis the coefficient of restitution.

For the numerical treatment of collisions Ardekani and Rangel (2008) simply calculate the contact force between both particles based on the rigidity force. If its magnitude is negative the collision process is ended since particles can not apply a tensile force to each other. The process is triggered if the distance between both particles is equal to

∆r_min. Here,∆r_min denotes the roughness height of the particle. If no collision occurs and the particles are not in contact, the sum of all forces over each particle, excluding gravitational and hydrodynamical forces, renders to be 0.