3- Results and Discussion
\label{results-and-discussion}
In this section at first the collision model is validated by comparison
with the experimental results for a particle wall collision. Then, the
robustness of the method is studied for various cases of non-spherical
particles including in collisions with friction and without friction.
Finally, the model is tested for the simulation of particle-particle
cases.
3.1. Validation: Particle-Wall
Collision
\label{validation-particle-wall-collision}
In order to validate the method for the collision with a stationary
wall, the experimental results of Gondret et al.
[76] which were also implemented
later in [27,
of steel with \(D_{p}=3\ mm\), \(\rho_{p}=7800\ kg/m^{3}\), \(\nu=0.3\) and \(E=2.4\times 10^{11}\) is released for a free fall
in silicon oil RV10 with \(\mu=0.01\ Pa.s\) and \(\rho_{f}=935\ kg/m^{3}\). Here, \(d_{p}\), \(\rho_{p}\), \(\nu\), \(E\) are respectively solid particle’s diameter, density, Poisson’s
ratio, and Young modulus. \(\rho_{f}\) and \(\mu\) are the fluid’s
density and viscosity. The simulation is performed in non-dimensional
form by normalizing the distance with the particle diameter and the
non-dimensional parameters \(Re=164\), \(Fr=11.8\) and \(\frac{{\gamma=\rho}_{p}}{\rho_{f}=8.342}\). The simulation is made
in a \(10\times 10\times 10\ \)domain with different minimum grid
sizes to study the effect of grid refinement. The cut off value for the
overlapping function is set to 10% of the particle diameter.
Consequently, when a particle approaches the wall, the course of
deceleration is \(\sigma=0.1\).
The particle starts to fall and reaches the final steady velocity before
hitting the wall. It starts bouncing and losing the kinetic energy and
finally it stays motionless on the wall. The very low or zero spring
stiffness values related to the final stage of the collision means that
the particle may not have enough collision force to hold it over the
wall after ultimate reduction of spring coefficient. This happens when
the bouncing height is so low that the particle cannot escape the
overlapping area. Therefore, the spring coefficient is not updated to
higher \(k_{N,ij}\) values related to the next bounce as it still has
the low bounce spring coefficient value \(k_{N,ij,b}\). To avoid
penetration of the particle to the wall in these cases, an extra upward
force equal to negative of weight is applied.
The results are presented in Figure (1) and they show the overall
agreement of the method with the experiment. Moreover, the independency
of the results from the grid also is approved. There is a small
discrepancy among the current results and the experiment which is more
noticeable in the velocity figure after the second collision. The
maximum velocity just after the collision is smaller in the simulation
than the experiment. This is because in the current soft collision
model, the effect of viscous dissipation is applied through the
restitution coefficient. In reality there is also viscous dissipation
through the computations which becomes significant at low Stokes number
values. For this reason, Kempe and Fröhlich
[14] ignored the hydrodynamic drag on
the particle during collision course. The discrepancy it is not apparent
in the first collision because of the short duration of collision at
high impact velocity which results in a harder collision. However, it
grows in the subsequent collisions with diminution of the impact
velocity and lower Stokes numbers. Moreover, the accumulation of errors
by time also is the other reason for increase in the discrepancy.
3-2. Effect of Collision
Course
\label{effect-of-collision-course}
In the presented collision method, the collision course is the only
parameter that is not defined. We defined the collision course, \(\sigma\), as the distance between the point with \(\left.\ \text{dψ}_{i}\left(\mathbf{x},t\right)\right|_{\text{cut\ off}}=\varepsilon\) and the particle interface. This section studies the dependency of the
results to \(\sigma\). A problem like the one in the section 3.1 is
considered but with different collision courses (Figure 2). The results
are related to minimum grid size of \(\Delta x=0.075\) and \(CFL=0.25\).
For \(\sigma=0.3\), the error related to the simulation is not
significant as it is mostly limited to the instants of impact. After the
collision .i.e. when particle leaves the overlapping region, it reaches
the actual path and the error is partially compensated. This occurs in
spite of the fact the collision now has become too soft and the sudden
changes in velocity are smoothed. Moreover, after the second collision,
the particle does not escape the very wide overlapping region and the
error increases.
Variation of height with time for the first collision of the particle is
represented in Figure (2c) for the particle using different collision
courses. When \(\sigma=0.3\) the particle does not reach the wall
completely leading to emerge of the small error. For \(\sigma=0.05\),
the difference between the simulation and experiment becomes noteworthy
by showing a higher restitution coefficient for the collision. In this
case, the collision duration is so short that the simulation becomes
sensitive to the small errors. As it could be observed in Figure (2c),
even though the particle path is very close to the experiment, after the
collision it bounces with higher velocity.
3-3- Collision of an Elliptical Particle with a Frictionless
Wall
\label{collision-of-an-elliptical-particle-with-a-frictionless-wall}
When a particle is not spherical, the rotation becomes important even if
the walls are frictionless. This is because in general the collision
normal vectors do not point to the center of gravity. Here, we consider
the collision of a free-falling ellipse on a frictionless wall. The
ellipse with \(D_{p}=(2.0,0.5)\) and \(\frac{\rho_{p}}{\rho_{f}=10}\) is located at \((5.0,7.0)\) in a \(10\times 10\) channel with four no-slip walls. The particle is
released and settled under the gravity. At Froude number \(Fr=0.51\),
we study two Reynolds numbers of \(Re=20\) and \(Re=2000\). The
Reynolds and Froude numbers are measured with respect to the major
diameter of the particle just before the collision happens. The initial
angle of the particle is \(\theta=\frac{\pi}{4}\) with respect to the
horizontal wall. Due to high asymmetry of the problem, the particle
should start to rotate after the start of motion. However, to avoid more
complexities, the rotational motion of the particle is restricted until
the collision occurs. During the falling time, the particle is forced to
fall with a unit downward vertical velocity. The minimum grid size is \(x=0.0375\) and \(CFL=0.2\).
Figure (3) shows the sequence of the particle position in the channel.
Due to the tilted placement of the particle at the time of impact, the
particle reaches the wall on its sharp edge. Then it bounces slightly
and rotates towards the final position of the particle with minimum
gravitational potential on the flat side of the particle.
Figure (4) represents the time history of the particle kinematical
variables. It is clear from Figures (3) and (4) that at \(Re=20\) after the impact of the particle with wall and a small vertical bounce,
the particle settles down to the final vertical position very fast. This
is because of the high dissipation at low Re. The other hand, the
particle continues its motion in the horizontal direction by sliding
over the frictionless wall after gaining a horizontal momentum. This
momentum damps due to dissipation in the fluid. At \(Re=2000\), the
lower dissipation of energy leads a number of subsequent collisions
after the first collision as well as a horizontal motion in the other
direction.
3-4- Tangential Collision
Model
\label{tangential-collision-model}
To validate the tangential collision model, we study a case similar to
Biegert et al. [27]. A sphere with \(D_{p}=0.125\) and \(\frac{\rho_{p}}{\rho_{f}=2.5}\) is exposed to
a linear shear flow in a \(1\times 1\times 1\) domain. Two cases in
different fluids with \(\nu=0.1\) and \(\nu=0.02\) are considered.
The top wall is driven with \(u_{w}=1\) leading to \(\text{Re}_{w}=10\) and \(\text{Re}_{w}=50\). The static and
dynamic friction coefficient as \(\mu_{s}=0.8\) and \(\mu_{k}=0.15\) related to silicate materials are considered. The linear velocity
profile in the domain varies with height and is applied to all
boundaries except the outlet which satisfies the condition of zero
gradient for the variables. The particle is placed on the bottom wall
and is subject to gravity of \(g=9.81\). A minimum grid size of \(0.05\) is considered with \(CFL=0.5\).
The particle is initially at stationary and with the start of the
simulation the collision is detected. It is held for \(T=\frac{tu_{w}}{H}=0.01\ \) for the flow to develop slightly and
avoid the spikes in the forces at the start of the simulation. Then the
particle is released leading to acceleration and motion to reach the
steady condition. We compared our results with Biegert et al.
[27] in Figure (5). While the results
are very similar, there is a small discrepancy which is probably because
of using a different diffuse interface solver. Figure (5) also shows
that at similar to Biegert et al.
[27], at \(\text{Re}_{w}=10\) the
translational velocity (\(U_{p}\)) is more than the rolling velocity
(\(R_{p}\omega_{p}\)) which means that the particle is sliding on the
wall. However, at \(\text{Re}_{w}=50\), the particle does not slide
and a rolling condition happens because \(U_{p}=R_{p}\omega_{p}\).
Here, we did not distinct between \(R_{p}\omega_{p}\) and \(U_{p}\) as
they are very similar.
Figure (6) represents the streamlines and the vortex contours on the
symmetric vertical plane parallel to the flow. It shows that the
contours are very similar for the two Reynolds number values but at \(\text{Re}_{w}=50\) there is a small asymmetry due to the wake
effect. The streamlines also show that the particle is in rolling
condition at \(\text{Re}_{w}=50\) while it slips significantly at \(\text{Re}_{w}=10\).
3-5- Collision of Rectangular Particle with
Wall
\label{collision-of-rectangular-particle-with-wall}
To study the effect of tangential part of the collision in this section
we consider a rectangular block of copper with \(2\times 0.5\) dimensions standing vertically on a wall. We consider the bottom wall
with zero velocity and in two situations: one case without friction and
the other one with static and dynamic friction coefficients similar to
copper on steel friction i.e. \(\mu_{s}=0.53\) and \(\mu_{k}=0.36\).
A uniform unit velocity is applied for the input and the top wall in the \(x\) direction. The particle is located vertically at \((3,1)\) in a \(15\times 6\) channel. With \(\frac{\rho_{p}}{\rho_{f}=8.96}\), \(Re=\frac{U_{\infty}L}{\nu}=100\) and \(Fr=\frac{U_{\infty}^{2}}{\text{gL}}=0.102\) related to \(\nu=0.01\), \(g=9.81\) and \(L\) as the average dimension of the
rectangle. The rectangle is fixed at first until \(T=0.1\) to have the
flow developed over the particle and avoid a sudden spike in force at
the start of the simulation.
The sequence of particle motion is shown in Figure (7) for both friction
and frictionless collisions. Figure (8-9) display the quantitative
variation in particle’s motion for the two cases. With friction, the
particle shows a limited motion which leads to a final stop at about \(T\sim 2.2\). Without friction, it continues to move to reach the end
of the channel with a steady velocity of \(u_{\text{final}}\sim 0.56\).
Even though the contact is frictionless, the particle does not reach the
freestream velocity because of the boundary layer on the wall. Moreover,
the particle is losing momentum by time as it has less flow blockage in
comparison with the start of the motion.
The differences between the two cases is not limited to the
translational motions as the two particles have different rotational
motions. With friction, the particle pivots around the front collision
point and rotates until it contacts the wall on the longer side and
becomes stable. On the other hand, without friction, the particle is
free to rotate around any of the bottom edges. Therefore, when the
applied moment on the particle is positive due to more blockage of the
flow on the bottom, the particle rotates on the other direction leading
a positive angular velocity. However, for both cases the rotation is
restricted because after the particle reaches the long edge, the
stabilizing effect of the weight does not allow to particle to rotate
more.
The effect of grid refinement on the variation of angle is studied in
Figure (10) and it shows that with refinement of the grid, the results
follow the same overall trend. However, there is a small discrepancy for
the small fluctuations on the angle which belongs to the time that the
particle hits the wall. This time is the time that the collision forces
become dominant and the results are more sensitive to the initial
conditions and also the coarse grid.
Figures (11-12) display the contour plots of vorticity for the
frictionless and frictional collision at different time instants. The
evolution of vorticity from the tip of the rectangle is noticeable which
in the case of the frictionless collision is attached to the tip and
moves with the flow and the particle. For the frictional collision, it
detaches the tip due to fast downward motion while another vorticity is
generated on the other side of the particle.
3-6- Collision of Two
Particles
\label{collision-of-two-particles}
Now, we study a drafting, kissing, and tumbling case similar to
Glowinski et al. [77] who used a
repulsion force to apply the normal collision for the disk shape
particles. Two similar particles start to fall under the gravity of \(g=981\) in a vertical channel with size \(2\times 6\) containing a
fluid with \(\nu=0.01\). The particles with \(d_{p}=0.25\) and
density ratio of \(\frac{\rho_{p}}{\rho_{f}=1.5}\ \)are initially
located at \((1,4.5)\) and \((1,5)\). Having the same horizontal
position, motion of the leading particle generates a wake which results
in a drag reduction on the trailing particle. Therefore, the trailing
particle will have greater velocity and reach the leading particle.
After the collision, the instability in the system leads to deviation of
the particles from the centerline of the channel.
Here, a minimum grid size of \(x=0.0125\) with \(CFL=0.5\) is
considered. Figures (13-14) shows the results for a collision course
equal to 10% of the particle diameter. It is clear from the results
that having a diffuse interface and a soft sphere collision model with a
thick collision course results in an error in velocity during the
collision period. This error will be compensated for the vertical
velocity after the collision period. For the horizontal velocity, the
small errors in collision accumulates and leads to greater deviation
because of high sensitivity to the initial conditions. In overall, the
error in the position of the particles is not significant because the
change in the direction of velocity components results in compensation
of error in the position as the integral of the velocity.