Penetration Depth of Free Falling Intruder into a Particles Bed in Fluid-Immersed Two-Dimension Spherical Particle System

In a system of spherical particles which is entirely immersed in fluid, an intruder is let free fall into the particles bed consisted of 220 particles with height about 13 particles diameter and width about 15 particles diameter, both are in average. The system is placed in a container with height of 40 cm and width of 15 cm. Bed particles have the same diameter of 1 cm, while the intruder has 4 cm. Water is the fluid with density of 1 g/cm3 and viscosity of 8.90×10−4 Pa·s, which is assumed to constant during the simulation. Intruder initial height is always the same and bed particles are let to be relaxed for about 2 s from its initial random configuration of 20×11 in height and width of bed particle diameter, before the intruder is dropped into the system. Bed particles has always the same density about 2 g/cm3, while intruder density ρ int is varied from 2 to 4.5 g/cm3 with increment of 0.5 g/cm3. It is observed that higher ρ int gives higher penetration depth after the density of 2.5 g/cm3. Density of 4.5 g/cm3 and beyond will give similar final result, since the intruder already reached bottom of the container. An empirical model for penetration depth as function of ρ int is proposed.


Introduction
Penetration of a larger object into a bed of grains is still interesting to investigage. It ranges from modeling of crater forming of due to meteor impact on planet surface [1], investigation of depthdependent resistance of penetration process [2], study of confinement influence [3], observing unified force of granular impact cratering [4], formulating granular bed viscosity [5], until seperation of two grains by one grain above [6]. In this work a whole system of spherical particles in two-dimension is immersed in fluid, so that both bed particles and intruder will be affected by drag force due to fluid viscosity and also by buoyant force due to fluid density. The intruder is released from a certain height z 0 and gives impact to the bed. How far it can penetrate the bed is then reported.

Simulation
A spherical particle or grain i will have density of ρ i , diameter of D i , and mass of where i = 1, 2, .., N, with index N is for the intruder, while 1 .. N -1 are for the remaining particles. Due to its surrounding fluid environment with viscosity η f and density ρ f particle i will have viscous force (or drag)  (2) and buoyant force with g  is earth gravity. Since the system is simulated on earth surface it will have also gravitational which differs from Eq. (3) in the sign and also the density. Position of particle i relative to particle j is distance between two particles and unit vector where for relative velocity and its unit vector can also be obtained using similar way. Using Eq. (6) and diameter of two particles, D i and D j , an overlap can be calculated with Using Eqs. (8) and (10) normal force on particle i due to particle j can be defined as [7] Then total force acted on particle i will be is the acceleration. At every time t Eq. (14) can be written in its full form as that shows a coupled differential equation, which is difficult to solve it analytically. Using forward difference method velocity of particle i at time t + Δt can be obtained and also the position with initial condition   must be known for all particles. And for interaction between particles and the container walls, it is formulated similar to Eq. (12), where a flat wall can approximated with a spherical grain with large diameter, so that only one grain diameter is used and unit vector of relative positition will become the normal vector of the wall. The simulation will be performed from t = t beg until t = t end time step Δt.
During the simulation three parameters are observed, which are vertical position of interuder z int , maximum vertical position of bed particle z max , and average vertical position of bed partices where vertical position of the intruder (with index N) is excluded.

Results and discussion
Values of parameters used in the simulation is shown in folowing Table 1, where all are represented in SI units. Simulation begins at t = t beg with creating spherical particles in N y × N z grid with slightly pertubated randomly and let the system condensed until t = t int . After that time an intruder is let free falling from height z 0 , hit the bed, and begin to penetrate it, as shown in Fig. 1.  Since particles bed are randomly generated before it condensed for each density ρ int there are five times repetition of simulation, which give average value of z avg , z max , and z int as shown in Fig. 4.  Fig. 3 (top) is one of the results contributing to average value in Fig. 4 (left), while Fig. 3 (bottom) is related to Fig. 4 (right). We can also calculate the difference between final and condensed configuration, as given in Fig. 5 with standard deviation from 1.2×10 -4 until 1.7×10 -2 for the right figure. It can seen from Fig. 5 (left) that Δz max and Δz avg has order of 10 -2 , while from Fig. 5 (right) that Δz int has order of 10 -1 . It can understand since Δz max and Δz avg have not change much, but Δz int is very dependent on initial configuration, that always be different. Monotonic increasing of Δz avg as increasing of ρ int indicates an interesting result since it is already average values.

Configuration in
Polynomial regression is used in fitting Δz int as function of ρ int and following relation 1 int