Object Impacted and Transported by Dry Granular Flow: 3D MPM-SDEM

. Granular geophysical flows pose a significant threat to human lives and infrastructure in mountainous regions globally. Property and victims impacted by these flows often get carried away and buried. Quickly locating the potential location of victims is critical for post-disaster recovery and rescue. However, the trajectory of object impacted by granular flows is difficult to predict because of the scale disparity between the sizes of granular particles and the object being transported. To deal with this multi-scale problem, a continuum-discrete solver, specifically the material-discrete element method (MPM-SDEM), was developed to simulate the trajectory of a cube (the object) impacted by a dry granular flow. Simulated results are used to understand the trajectory of objects under different flow and impact dynamics.


Introduction
After a granular geophysical flow occurs, first responders typically have a 72-hour window to recover survivors. Disaster recovery is challenging because the affected area is several hectares in area [1] and debris volumes can reach up to millions of cubic meters [2]. Evidently, it is extremely difficult to find the final location of the transported victims due to the scale of an event and the complex dynamics of granular flows [3].
The trajectory of objects impacted by a debris flow is governed by the driving forces from impact and the resisting forces from the object and its boundaries. Previous studies have examined the impact force required to mobilise [4][5] and transport an object [6]. But the force on an object could be influenced its submerged state in flow [7]. A governing equation for boulder transportation in fully submerged, partially submerged, and subaerial (not in contact with flow) states has been proposed [8]. But this analytical model could only deal with one-dimensional movement, while the movement of an object impacted by the flow is a 3D problem. Besides, only a one-dimensional model could not deal with diverse moving modes (sliding, rolling and saltation) of an object. To solve the problem about moving modes of the object impacted in flow, the onedimensional model has been improved [6]. The movement of an object impacted by flow is the interaction between the object and flow, but these models only solve the impact force on an object, without considering the reaction force on the flow from this object. Besides, they focus on the water flow instead of frictional geophysical flows, which are fundamentally different.
The trajectory of an object impacted by a flow could be abstracted as the process of interactions between a * Clarence Edward Choi: cechoi@hku.hk continuum material and a discrete body [9]. The force on a moving object varies spatiotemporally, so it is impossible to solve sliding, rotation, and saltation accurately at each time section based on analytical solutions. Many researchers have studied the interaction between a rigid body and flow based on numerical methods. Liu [10] has proposed a unified MPM framework to solve the problem about force interaction. But the method is limited by the size of grid, so it is difficult to balance the size of the impacted object and computational efficiency. Then an improved 2D coupling algorithm, named MPM-SDEM [11], was proposed to solve this problem. For an object impacted by a flow, the interaction is 3D problem, because the flow exerts external forces on the object from all sides of an object.
To solve the 3D multi-scale granular-solid interaction, we adopt the hybrid continuum-discrete solver, where the material point method is coupled with the discrete element method. Furthermore, we did nontrivial extensions on the discrete solver using threedimensional spheropolyhedron element (SDEM) to accurately capture the trajectory of non-spherical objects. The shape of object has an influence on the force on the object, so it could be represented by only a sphere. Here this method could utilize realistic particle shape [12].
In this extended abstract, details of a new continuum-discrete method, 3D MPM-SDEM, will be presented to simulate an object impacted and transported by the flow, in both partially and fully submerged states, to the trajectory of objects under different flow and impact dynamics.
The test object is considered as a rigid cube for simplicity. Meanwhile, to simplify the numerical model and reduce the computational cost, dry granular is adopted as the flow material, modelled as it is wellbenchmarked in the literature.

Modelling granular flow
MPM is the combination of the Lagrangian and Eulerian methods. It solves the movement of a continuum body by discretizing it into a set of particles [13], so, as a continuum body, the deformation and displacement of granular flow could be simulated based on MPM. The first step is discretization, meaning that all mass of the flow is distributed on the particles. After discretization, the position of a material point at each time step is updated according to the explicit method [5]. The Drucker-Prager yield criterion [14] (D-P model) is used to describe the failure of initial slope, and the μ(I) rheology [15] is used to describe movement of flow material.

Modelling impacted object
We first introduce the procedures -shrinkage and expansion -for generating a spheropolyhedron in the discrete solver, as shown in Fig. 1. In the shrinkage step, the original polyhedron's geometry is shrunk volumetrically by a distance of r. The polyhedron will then be expanded through Minkowski addition [12] with the sphere in r radius, which is in the expansion step. In this way, a contact layer -a volume that detects and calculates contact forces -is coated outside the polyhedron.
A SDEM particle has three elements for contact handling namely, vertexes {V}, edges {E}, and faces {f}. Among the three types of elements from different particles, there are three possible contact modes (shown in Fig. 2) between two SDEM particles (i.e., particle and particle ): (i) vertex-face contact, (ii) vertex-edge contact, (iii) edge-edge contact. Then the contact force ij which particle exerts on the particle can be expressed as follows: where ( i , j ) is the force produced by contact between a vertex of particle and a face of particle ; ( i , j ) is a force produced by contact between a vertex of particle and an edge of particle ; and ( i , j ) is the force produced by contact between an edge of particle and an edge of particle . To solve the contact force c between particle and particle in Fig.3, the overlapping depth δ between the vertices of particle and edges of particle is calculated [11]. For example, δ between a vertex of particle and a face of particle can be expressed as: ( i , j ) = i + j − ( i , j ) , (2) where ( i , j ) is the Euclidean distance between the vertex of particle and the closest surface of particle .
i and j are the radii of spheres used to shrink and expand particles and , respectively. After ( i , j ) is solved, the contact force, ( i , j ), is calculated based on Hertzian mechanics [16].

Fig. 4
A two-dimensional schematic diagram of contact between an MP and expanded boundary of a SDEM particle (object) In this section, the interaction between MPM and SDEM is achieved through the contact between material point and SDEM particle. And in the coupling, the SDEM particle is the impacted object here.
A material point (MP), when in contact with an object, is assigned a radius. The coupling contact force can be conveniently calculated based on the overlapping distance between the MP radius and the contact layer on the object. The MP radius, , can be solved as: where m is the mass of a MP, and ρ is the density of a MP. Fig. 4. shows the coupling contact force Fc calculated based on the overlapping distance δ between an MP and an object. For the MP that contacts the object, an extra normal force term − c is added to its momentum conservation equation at each simulation time step. In turn, the contact force Fc is added to the resultant force on the object. It is worthwhile to note that an object is most likely in contact with multiple material points. Hence, the total coupling force should be summed as: where n is the number of MPs that contact the object, and i is the -th MP. ci ( i , ) is the contact force https://doi.org/10.1051/e3sconf/202341502010 , 02010 (2023) E3S Web of Conferences 415 DFHM8 between -th MP and one face of the object, and ci ( i , ) is the contact force between -th MP and one edge of object.

Simulation setup and preliminary results
A simulation of a channelized dry granular flow impacting a cube in partially submerged state is simulated (Fig. 5). The channel includes an inclined part and horizontal part. θ is the inclination. 0 is the horizontal length of the inclined channel and 1 is the length of the horizontal channel. The material properties of dry granular was adopted from [9]. In this simulation, the dry granular begins to deform and flow because of gravity force (simulated by MPM). After accelerating on the inclined part, the dry granular flow impacts the cube (object simulated by SDEM). With different velocity and depth of the dry granular flow, the interaction between the cube and the granular flow, and the whole trajectory progress of the cube could be studied. In this simulation, the length of the cube is 0.04m, the horizontal length of the inclined part, l0, is 1m, and the length of the horizontal part, l1, is 2m. The aerial view of the simulation results is shown in Fig.6. The blue area is the granular flow, and the red square is the location of cube.
(a) (b) Fig.6 The movement of object impacted by granular flow Fig.7 The displacement and velocity of object centre Fig.6 shows the object began to move under the impact of the granular flow. Fig.7 shows the change of location and velocity and object centre during the simulation process. The horizontal axis is the simulation time, whose unit is s. The left vertical axis is the location of the object centre, and the right vertical axis is the velocity of the object centre. At about 0.17s, the object began to move under impact from the granular flow, then its velocity increased to the maximum value, 1.08m/s. Because the difference between cube velocity and flow velocity decreased, the impact force on cube shrank as well. With the basal frictional resistance, the cube began to decelerate and stop.
In this unsubmerged simulation case, the square object is placed at the central line of the flow path, so the net force on it is along the flow path, making it move only along the flow direction. For further study, different parameters should be considered to study the relationship between trajectory of object and the different flow conditions.