# Initial Moving Mechanism of Densely-Packed Particles Driven by a Planar Shock Wave

Authors: Hua Lv; Zhongqi Wang; Yunming Zhang; Jianping Li
Journal: Shock and Vibration
 (2021)

---

## Abstract

The initial moving mechanism of densely packed particles driven by shock waves is unclear but vital for the next accurate calculation of the problem. Here, the initial motion details are investigated experimentally and numerically. We found that before particles show notable motion, shock waves complete reflection and transmission, and stress waves propagate downstream on particle skeleton. Due to the particle stress wave, particles successively accelerate and obtain an axial velocity of 6–8 m/s. Then, the blocked gas pushes the upstream particles integrally to move downstream, while the gas flow in the pores drags the downstream particles to separate dramatically and accelerate to the velocity of 60–70 m/s. This gas push-drag dual mechanism transforms densely packed particles into a dense gas-particle cloud, which behaves as the expansion phenomena of the dense particles.


---

## Body

## 1. Introduction

The movement of densely packed particles driven by shock waves has been widely used in many civil and military fields, such as Pressurized Fluidized Bed Combustion (PFBC), Fire Extinguishing Agent Dispersal, Needle-Free Drug Delivery, Supersonic Cold Spray, Fuel-Air Bomb, and Low Collateral Damage Ammunition. It is a typical compressible gas-particle two-phase flow problem that has essential academic and application values. The research studies of this problem were mainly focused on the overall parameters of the far-field particle cloud, such as the maximum accelerating velocity and scattering range. However, the internal concentration distribution, which has not been sufficiently investigated, largely restricts its accurate application [1]. The effectiveness of the internal concentration distribution for its application fundamentally depends on the temporal and spatial evolution of the particle concentration in the particle cloud evolving from the initial driving process. The driving process can usually be divided into three stages depending on the particle volume fraction ϕp: the granular flow (ϕp ≥ the particle volume fraction of packed particles ϕp,packed), the dense gas-particle flow (0.01 < ϕp < ϕp,packed), and the dilute gas-particle flow (ϕp < 0.01) [2]. However, the research studies about these two-phase flows were started from the simplest single-particle and dilute particles.In compressible dilute gas-particle flows, collisions between particles can be ignored, and particles are mainly affected by the flow field. Igra and Takayama [3] laid some small spheres on the shock tube floor to be driven by planar shock waves and photographed to investigate the drag force on it. Devals et al. [4] captured the trajectory of a single sphere hanged in a shock tube by a multiple exposure shadowgraph technique coupled with a high-speed drum camera. Jourdan et al. [5] similarly snapped the spheres trajectories by shadowgraph technique and gave a correlation for the nonstationary drag coefficient of a sphere. Parmar et al. [6] calculated the unsteady inviscid force in compressible flow. It is four times larger than that predicted using incompressible theory. Ling et al. [7, 8] emphasized that the effects of unsteady contributions get importance when the particle-to-gas density ratio is small. McGrath et al. [9] simulated the compressible particle-dispersed flow using an Euler-Euler method.In compressible dense gas-particle flows, the interparticle interactions and the interactions between particles and the reflected or diffracted waves from neighboring particles become important [10–13]. Li et al. [14] calculated the drag force in a dense gas-particle flow. Mehta et al. [15] presented a three-dimensional fully resolved numerical simulation of shock propagation through a fixed array of particles and investigated the unsteady drag force. Wagner et al. [16] designed a gravity-fed apparatus to produce a dense particle curtain (ϕp = 0.2) in a shock tube and observed the expansion of the particle curtain after the passage of the shock wave. Theofanous et al. [17] carried out a similar experiment, and the dense particle curtain (ϕp = 0.4) shows a similar expansion phenomenon. Wagner et al. pointed out that the expansion of the particle curtain may be associated with some gas property differences across the streamwise thickness of the particle curtain. They also addressed that the exact physical mechanism remained unclear. Ling et al. [18] simulated the interaction between the planar shock wave and the dense particle curtain by an Euler–Lagrangian method. In their simulation, the effects of gas compressibility and particle volume fraction on the drag force are all considered. Deng et al. [19] and Jiang et al. [20] simulated the shock-particle curtain interaction by direct numerical simulation and emphasized the importances of the viscosity and collision models.With regard to the densely packed particles, the research studies were mainly focused on the shock attenuation performance of the densely packed particle layer and on the final maximum velocity of particles [21–24]. However, the initial driving stage may show some special movement phenomena, for example, the particle jet in the explosion-driven particle dispersal [25]. These special phenomena will determine the subsequent temporal and spatial evolution of the particle concentration. Stewart et al. [26] used the discrete element method code LIGGGHTS to simulate the traversal of a planar shock into a granular bed sealed at the bottom of a tube and statistically investigated the traveling wave structure in the particle concentration downstream of the shock. Nevertheless, the force on the particles and the movement of them during the stage from initial densely packed particles to a dense gas-particle flow is still unclear. Therefore, the initial densely packed particles moving mechanism needs further investigation.When the initial particles are diluted, the particle acceleration is definitely due to the unsteady drag force from the gas phase. As the initial particles become dense and even packed, the particle collision, contact, and rearrangement are engaged in sequence. These engagements cause a problem which of the particle stress and the gas drag is the decisive reason for the initial motion of particles. In this paper, a densely packed particle wall is driven by a planer shock wave in a shock tube to capture the start-up characteristics of the densely packed particles. An Euler–Lagrangian coupling simulation of this shock-particle interaction is conducted to investigate the force on the particles and moving mechanism in the start-up process. The contribution of the particle stress on the initial motion of particles is calculated and compared with that of the gas drag.


## 2. Experiment of the Initial Movement of the Densely Packed Particles

### 2.1. Experiment Setup

The experiments were carried out in a horizontal multiphase shock tube cooperated with a high-speed schlieren imaging system and a synchronized pressure measurement system. A brief description of the experimental setup is presented here, and further details can be found in Lv et al. [27]. The driver and driven sections of the shock tube, which are separated by the diaphragm, are filled with high-pressure air and atmospheric air, respectively. When the pressures of the driver air P1 are 3, 4, and 5 MPa, the Mach numbers of the incident shock measured by the crossing time Ms are 2.18, 2.31, and 2.38, respectively.The front view schematic of the driven section of the shock tube is given in Figure1. The driven section made of C45E4 steel is totally 1640 mm long. It consists of four subsections, which are connected by flanges. The first subsection transforms the cross section from the circular tube which has the same cross section as the driver gas cylinder to the rounded-angle square tube as the same as the last 4 subsections. The locations of the pressure transducers PT1-PT6 with the type of KISTLER 211B also are listed in Figure 1. The third subsection of the driven section is designed as the test section. The area of the observation window in the test section for the high-speed schlieren is 130 mm × 40 mm. A densely packed sand wall carried by a fixed thin-walled particle carrier and wrapped by papers is located in the middle of the test section. The sands packed in the wall are quartz sands with the most frequent diameter Dm of 0.646 mm and a density of 2.64 g/cm3. The initial porosity of the densely packed particle wall is measured to be 0.561.Figure 1
Front view schematic of the driven section of the shock tube.


### 2.2. Experiment Results

#### 2.2.1. Shock Wave Structure

Figure2 shows the pressure histories of a typical experimental case. The pressure of the driver gas P4 is 4 MPa, the most frequent diameter of the sands Dm is 0.646 mm, and the length of the sand wall L is 10 mm. The trace of PT3 shows both the incident shock at t = 166.50 ms and reflected shock at t = 166.76 ms, while the trace of PT4 shows the transmitted shock at t = 166.70 ms. The shock waves at the same time are also seen in Figure 3, which are the typical images captured by the high-speed schlieren system.Figure 2
Typical pressure histories at the two side of the particle wall.Figure 3
Typical high-speed schlieren images (I, R, and T denote incident, reflected, and transmitted shock wave; UF and DF denote upstream and downstream front of particles; PC denotes fixed thin-walled carrier of packed particles).


#### 2.2.2. Particles Movement

Besides the shock wave structure, Figure3 presents the locations of two interfaces of the densely packed particles. At t = 166.56 ms, the incident shock wave impacts the upstream front of the sand wall. 0.29 microseconds after the onset of the interaction, i.e., at t = 166.85 ms, the downstream front of the densely packed particle wall begins to show noticeable downstream motion. At t = 167.65 ms, the upstream front begins to leave the sand wall carrier. At t = 167.85 ms and 168.65 ms, the upstream front is still in the region of the observation window.


#### 2.2.3. Velocities of the Upstream and Downstream Fronts of the Particles

Based on the arriving times of the upstream and downstream fronts of the dense particles at specific locations such as their initial positions, positions of PT4 and PT5, which are captured by the high-speed schlieren and synchronized pressure measurement systems, the positions of the two fronts are pointed and the trajectories of the two fronts are fitted in Figure4. After the onset of the interaction (at t = 166.56 ms), the downstream front of the densely packed particle wall firstly moves downstream, while the upstream front is still stationary. Based on the captured fronts’ positions, the trajectories of the fronts can be fitted. By taking a tangent line at a point on the trajectory, we can get the velocity of the front at the selected moment. For example, at t = 167.42 ms, the downstream front already gets the velocity of 62.96 m/s. As the velocity of the downstream front of the particles is greater than that of the upstream front, the densely packed particles are gradually pulled apart and transformed into a dense gas-particle cloud.Figure 4
Trajectories of the upstream and downstream fronts of the particles.


## 2.1. Experiment Setup

The experiments were carried out in a horizontal multiphase shock tube cooperated with a high-speed schlieren imaging system and a synchronized pressure measurement system. A brief description of the experimental setup is presented here, and further details can be found in Lv et al. [27]. The driver and driven sections of the shock tube, which are separated by the diaphragm, are filled with high-pressure air and atmospheric air, respectively. When the pressures of the driver air P1 are 3, 4, and 5 MPa, the Mach numbers of the incident shock measured by the crossing time Ms are 2.18, 2.31, and 2.38, respectively.The front view schematic of the driven section of the shock tube is given in Figure1. The driven section made of C45E4 steel is totally 1640 mm long. It consists of four subsections, which are connected by flanges. The first subsection transforms the cross section from the circular tube which has the same cross section as the driver gas cylinder to the rounded-angle square tube as the same as the last 4 subsections. The locations of the pressure transducers PT1-PT6 with the type of KISTLER 211B also are listed in Figure 1. The third subsection of the driven section is designed as the test section. The area of the observation window in the test section for the high-speed schlieren is 130 mm × 40 mm. A densely packed sand wall carried by a fixed thin-walled particle carrier and wrapped by papers is located in the middle of the test section. The sands packed in the wall are quartz sands with the most frequent diameter Dm of 0.646 mm and a density of 2.64 g/cm3. The initial porosity of the densely packed particle wall is measured to be 0.561.Figure 1
Front view schematic of the driven section of the shock tube.


## 2.2. Experiment Results

### 2.2.1. Shock Wave Structure

Figure2 shows the pressure histories of a typical experimental case. The pressure of the driver gas P4 is 4 MPa, the most frequent diameter of the sands Dm is 0.646 mm, and the length of the sand wall L is 10 mm. The trace of PT3 shows both the incident shock at t = 166.50 ms and reflected shock at t = 166.76 ms, while the trace of PT4 shows the transmitted shock at t = 166.70 ms. The shock waves at the same time are also seen in Figure 3, which are the typical images captured by the high-speed schlieren system.Figure 2
Typical pressure histories at the two side of the particle wall.Figure 3
Typical high-speed schlieren images (I, R, and T denote incident, reflected, and transmitted shock wave; UF and DF denote upstream and downstream front of particles; PC denotes fixed thin-walled carrier of packed particles).


### 2.2.2. Particles Movement

Besides the shock wave structure, Figure3 presents the locations of two interfaces of the densely packed particles. At t = 166.56 ms, the incident shock wave impacts the upstream front of the sand wall. 0.29 microseconds after the onset of the interaction, i.e., at t = 166.85 ms, the downstream front of the densely packed particle wall begins to show noticeable downstream motion. At t = 167.65 ms, the upstream front begins to leave the sand wall carrier. At t = 167.85 ms and 168.65 ms, the upstream front is still in the region of the observation window.


### 2.2.3. Velocities of the Upstream and Downstream Fronts of the Particles

Based on the arriving times of the upstream and downstream fronts of the dense particles at specific locations such as their initial positions, positions of PT4 and PT5, which are captured by the high-speed schlieren and synchronized pressure measurement systems, the positions of the two fronts are pointed and the trajectories of the two fronts are fitted in Figure4. After the onset of the interaction (at t = 166.56 ms), the downstream front of the densely packed particle wall firstly moves downstream, while the upstream front is still stationary. Based on the captured fronts’ positions, the trajectories of the fronts can be fitted. By taking a tangent line at a point on the trajectory, we can get the velocity of the front at the selected moment. For example, at t = 167.42 ms, the downstream front already gets the velocity of 62.96 m/s. As the velocity of the downstream front of the particles is greater than that of the upstream front, the densely packed particles are gradually pulled apart and transformed into a dense gas-particle cloud.Figure 4
Trajectories of the upstream and downstream fronts of the particles.


## 2.2.1. Shock Wave Structure

Figure2 shows the pressure histories of a typical experimental case. The pressure of the driver gas P4 is 4 MPa, the most frequent diameter of the sands Dm is 0.646 mm, and the length of the sand wall L is 10 mm. The trace of PT3 shows both the incident shock at t = 166.50 ms and reflected shock at t = 166.76 ms, while the trace of PT4 shows the transmitted shock at t = 166.70 ms. The shock waves at the same time are also seen in Figure 3, which are the typical images captured by the high-speed schlieren system.Figure 2
Typical pressure histories at the two side of the particle wall.Figure 3
Typical high-speed schlieren images (I, R, and T denote incident, reflected, and transmitted shock wave; UF and DF denote upstream and downstream front of particles; PC denotes fixed thin-walled carrier of packed particles).


## 2.2.2. Particles Movement

Besides the shock wave structure, Figure3 presents the locations of two interfaces of the densely packed particles. At t = 166.56 ms, the incident shock wave impacts the upstream front of the sand wall. 0.29 microseconds after the onset of the interaction, i.e., at t = 166.85 ms, the downstream front of the densely packed particle wall begins to show noticeable downstream motion. At t = 167.65 ms, the upstream front begins to leave the sand wall carrier. At t = 167.85 ms and 168.65 ms, the upstream front is still in the region of the observation window.


## 2.2.3. Velocities of the Upstream and Downstream Fronts of the Particles

Based on the arriving times of the upstream and downstream fronts of the dense particles at specific locations such as their initial positions, positions of PT4 and PT5, which are captured by the high-speed schlieren and synchronized pressure measurement systems, the positions of the two fronts are pointed and the trajectories of the two fronts are fitted in Figure4. After the onset of the interaction (at t = 166.56 ms), the downstream front of the densely packed particle wall firstly moves downstream, while the upstream front is still stationary. Based on the captured fronts’ positions, the trajectories of the fronts can be fitted. By taking a tangent line at a point on the trajectory, we can get the velocity of the front at the selected moment. For example, at t = 167.42 ms, the downstream front already gets the velocity of 62.96 m/s. As the velocity of the downstream front of the particles is greater than that of the upstream front, the densely packed particles are gradually pulled apart and transformed into a dense gas-particle cloud.Figure 4
Trajectories of the upstream and downstream fronts of the particles.


## 3. Simulation of the Initial Movement of the Densely Packed Particles

Based on the above experiment measurements, we can establish the following images. Immediately after an incident shock wave impacts against the densely packed particle wall, a reflected shock wave is reflected from the upstream front of the wall, and the sands will be compacted suddenly. Then, stress waves among the densely packed particles will arise and propagate through the skeleton formed by the compacted particles. Meanwhile, the shock wave will also enter into and pass through the pores between the particle skeleton. Thus, which is the primary mechanism that the densely packed particles start to move and scatter? We simulated the interaction between shock waves and densely packed particles by an Euler–Lagrangian coupling method in the explicit finite element solver AUTODYN 3D and then dug into the details of the initial movement.


### 3.1. Validation of the Euler–Lagrange Coupling in AUTODYN 3D

#### 3.1.1. Modeling of a Classic Single-Particle Accelerating Experiment

To validate the particle acceleration calculating ability of the Euler–Lagrange coupling in AUTODYN 3D, we chose a classic single-particle accelerating experiment in a shock tube. A similar validation has been conducted by Price et al. [28] in their shock-particles simulation work. Devals et al. [4] investigated the acceleration of a single sphere on a spider web driven by a planar shock wave in a horizontal shock tube. The tube has a 750 mm long high-pressure section and a 3000 mm long low-pressure section. The inner cross section of the tube is 80 mm × 80 mm. The downstream end of the low-pressure section is closed with a flange which reflects incident shock waves. A 2 mm diameter nylon sphere is hanged 30 mm from the downstream closing flange. The measured velocity of the single sphere in the unsteady flow agrees well with Igra and Takayama’s work, as shown in Figure 5. Igra and Takayama [3] just laid some small spheres on the shock tube floor to be tested.Figure 5
Particle velocities measured in literatures [3, 4] and simulated in this work.In shock tube experiments, the low-pressure section is usually designed long enough to gradually develop the generated shock wave into a plane shock wave, while in simulations, the generated shock wave is naturally planar. In addition, according to the basic theory of one-dimensional unsteady flow, the strength of the shock wave and the properties of the high-speed gas flow behind the shock wave depend only on the pressure and adiabatic properties of the gas of the high- and low-pressure sections rather than the length of the two sections. Thus, in this simulation, we designed the total length of the shock tube to be 1000 mm, as shown in Figure6. Moreover, in order to ignore the impact of the rarefactive waves reflected from the upstream end of the high-pressure section on the incident shock wave, a long high-pressure section is designed to avoid reflected rarefactive waves catching up with the incident shock wave. So, the high- and low-pressure sections are set to be 750 and 250 mm, respectively. Other gas properties of the high- and low-pressure sections are the same as in the experiments. Furthermore, to furtherly save computing time, we designed a slender computing domain of which the cross section is 6 mm × 6 mm. A nylon sphere that has a diameter Dp = 2 mm and density ρp = 1130 kg/m3 is positioned 30 mm from the downstream end, and the distance is the same as that in Devals et al.’s shock tube. The detailed material model of the nylon sphere is omitted here, which is similar to the sands model presented in Section 3.2.2. A gauge point is defined on the nylon sphere to record the particle motion information during the simulation.Figure 6
Computing domain and meshing withδ = 0.1Dp.


#### 3.1.2. Validation for the Computation of Shock-Driven Particle Acceleration

The 1000 mm × 6 mm × 6 mm computing domain is discretized with uniform cubic mesh. Different mesh sizesδ = 1Dp, 0.5Dp, 0.3Dp, 0.2Dp, and 0.1Dp are used to assess the calculating ability of the Euler–Lagrange coupling in AUTODYN 3D for shock-driven particle acceleration. For instance, Figure 6 shows the enlarged local details with a grid size of 0.1Dp. Figure 5 gives the simulated velocity-time curves of nylon particles under these five grid conditions. These five simulation results are also compared with experimental data of Devals et al. and Igra et al.In Figure5, particle motions under different meshing conditions all undergo two processes: the acceleration process under the incident shock and the deceleration process as a result of the reflected shock. The origin time is set at the moment when the incident shock wave reaches the particle. At t = 0 ms, the particle starts to accelerate. Due to the long high-pressure section and the short low-pressure section, the particle acceleration process is not affected by the rarefactive waves reflected from the upstream end of the high-pressure section. The particle keeps accelerating with constant acceleration until the reflected shock wave from the downstream end of the tube reaches the particle. While the particle accelerations in Devals et al. and Igra et al.‘s experiments finally show a decline because of the rarefactive waves. After the acceleration stage, the particle motion comes into the decelerating stage. As the grid size decreases from 1Dp to 0.1Dp, the calculated velocity curve of the particles gradually approaches the experimental data.Since the spacing of the computation grid is much smaller than the size of the particle, the high-speed postshock gas around the particle is resolved, and the particle moves under the hydrodynamic force exerted by the resolved surrounding gas flow. Figure5 shows that the Euler–Lagrange coupling in AUTODYN 3D gives consistent drag forces to the experimental ones when the grid size is smaller than 0.3Dp. On the other hand, as the grid condition of 0.3Dp under a large number of particles has reached the maximum calculating power of our computer, the grid of 0.3Dp is chosen to accomplish the calculation. Moreover, based on Figure 5, the more fine grid condition would give a bit higher particle velocity because of the drag force, which would enhance the difference between the two velocities because of particle stress and gas drag.In conclusion, under the above meshing conditions, this detailed simulation method is a particle-resolved direct numerical simulation. Besides, dealing with object collisions is the strength of AUTODYN. So, although this fine meshing condition restricts the application of Euler–Lagrange coupling in AUTODYN 3D to a large number of particles, this direct numerical simulation method is useful for understanding more details of shock-particle interaction when the number of particles is relatively small.


### 3.2. Numerical Simulations of the Interaction between a Shock Wave and Packed Particles

On account of the validation in Section3.1.2, we safely employed AUTODYN 3D to simulate the interaction between a planar shock wave and packed particles. This section describes the geometric and material model of the simulation. Then, the simulation results of shock wave strength and particle motion are presented.


#### 3.2.1. Geometric Model

Since the sidewall friction of the shock tube has no noticeable effect on particle acceleration [22], a slender computing domain is defined to save computing time. The length of the computing domain is the same as the length of the experimental shock tube, and the cross section of the domain is a small part of that of the real shock tube. The detailed dimension of the Euler domain is 1900 mm × 2.4 mm × 2.4 mm. On balance between the calculation accuracy discussed in Figure 5 and the computation time, the grid size is set to be 0.3Dp. The Euler mesh is 0.2 mm × 0.2 mm × 0.2 mm, as shown in Figure 7. The boundary condition of the downstream end of the Euler domain is “flow out.” Six fixed gauge points are defined within the Euler domain to monitor the gas pressure and velocity. They have the same x coordinate values as those of the real pressure transducers (PT1, …, PT6) mounted on the experimental shock tube, as described in Figure 1.Figure 7
Computational domain and grids.The sands used in the experiment were sieved by vibratory sieve shaker RETSCH AS200 and then analyzed by particle analyzer CAMSIZER P4. Based on the dynamic image analysis principle, CAMSIZER P4 measures the equivalent spherical diameter, which is the equivalent diameter of a sphere that has the same volume as the irregular particle. The measured particle size distribution of the sieved sands is of Gaussian distribution with a most frequent diameterDm of 0.646 mm and a variance of 0.038 mm, and the percentage of particles larger than 0.6 mm is 88.69%. For the above reasons, the spherical particle with a uniform diameter of 0.6 mm is used under the current simulation conditions.As AUTODYN 3D has a disadvantage on the complicated arrangement of large amounts of solid particles, the particles are initially arranged in the simple cubic packing model in this work. Figure8 illustrates the initial particle packing model used in the simulation. Nevertheless, the influence of this initial simple particle arrangement can be eliminated because of the following reason. When a planar shock wave impacts the particle wall, the particles inside the wall are further compacted and rearranged in both simulation and experiment. The porosity of the rearranged particles in the simulation will further reduce to approach the porosity of the real rearranged particles.Figure 8
Simple cubic packed particles.


#### 3.2.2. Material Model

The high-pressure and low-pressure gases in the shock tube obey the ideal gas state equation. The pressureP is related to the density and energy by(1)P=γ−1ρe,where γ is adiabatic index, ρ is the air density, and e is the specific internal energy.The adiabatic index of the air is a constant 1.4. According to the standard constants from the AUTODYN material library, the density of the atmospheric air is set to be 1.225 × 10−3 g/cm3, and the specific internal energy of that is 2.068 × 105 kJ/kg. The density of the high-pressure air is set to be 4.677 × 10−2 g/cm3, and the specific internal energy of that is 2.138 × 105 kJ/kg, and then the pressure in the driver chamber will be 4 MPa.The SAND state equation embedded in AUTODYN is used to describe the mechanical properties of sands under dynamic loadings. The SAND state equation is given by Laine and Sandvik [29] as a piecewise linear curve with ten points, namely, density ρ as a function of pressure P. And the bulk elastic modulus K in the model is calculated by(2)K=ρc2.The parameters of the SAND state equation are listed in Table1.Table 1
Parameters of the SAND state equation.
Pressure (kPa)Density (g/cm3)Sound speed (m/s)Density (g/cm3)01.6742.65E + 021.6744.58E + 031.7398.52E + 021.7451.50E + 041.8741.72E + 032.0862.92E + 041.9971.88E + 032.1475.92E + 042.1442.27E + 032.39.81E + 042.252.96E + 032.5721.79E + 052.383.11E + 032.5982.89E + 052.4854.60E + 032.6354.50E + 052.5854.63E + 032.6416.51E + 052.674.63E + 032.8


#### 3.2.3. Shock Wave Structure Compared with the Experiment Results

To compare with the experimental case in Section2.2 (P4 = 4 MPa, Dm = 0.646 mm, and L = 10 mm), the particle layers of simple cubic packing are set to be 16. Therefore, the length of the simple cubic packed particle wall is 9.6 mm. The pressure histories of PT3 and PT4 obtained from the simulations are compared with the experimental results in Figure 9. The time at which the incident shock wave reaches PT3 is set to be the origin time. Then, the reflected waves of the simulation and experiment return to PT3 almost at the same time. The transmitted shock waves of the simulation appear a little later than that of the experiment. The strengths of the incident shock waves of the simulation are equal to those of the experiment. The reflected and transmitted shock strengths of the simulation are almost equal to those of the experiment. In sum, the simulation which has high consistency with the experiment can be used to present much more detailed information about the shock-particle interaction.Figure 9
Shock wave structure comparison between simulation and experiment.


#### 3.2.4. Initial Movement Pattern of the Packed Particles

By means of the Euler–Lagrange coupling method, the initial movement of the packed particles driven by a shock wave is presented in Figure10. Under the impact of the shock wave, the particles are compressed and rearranged. When the compressed particles integrally move downstream, the particles of the most downstream layer (i.e., the 16th layer) firstly separate from the whole compressed particles and then the particles of the 15th layer separate. Thus, the particles of the whole compressed particles separate and accelerate towards downstream layer by layer. Why do the particles move in such a pattern? The forces exerting on the particles and the velocity of the particles need further exploration in the next section.Figure 10
Simulated initial movement pattern of packed particles.


## 3.1. Validation of the Euler–Lagrange Coupling in AUTODYN 3D

### 3.1.1. Modeling of a Classic Single-Particle Accelerating Experiment

To validate the particle acceleration calculating ability of the Euler–Lagrange coupling in AUTODYN 3D, we chose a classic single-particle accelerating experiment in a shock tube. A similar validation has been conducted by Price et al. [28] in their shock-particles simulation work. Devals et al. [4] investigated the acceleration of a single sphere on a spider web driven by a planar shock wave in a horizontal shock tube. The tube has a 750 mm long high-pressure section and a 3000 mm long low-pressure section. The inner cross section of the tube is 80 mm × 80 mm. The downstream end of the low-pressure section is closed with a flange which reflects incident shock waves. A 2 mm diameter nylon sphere is hanged 30 mm from the downstream closing flange. The measured velocity of the single sphere in the unsteady flow agrees well with Igra and Takayama’s work, as shown in Figure 5. Igra and Takayama [3] just laid some small spheres on the shock tube floor to be tested.Figure 5
Particle velocities measured in literatures [3, 4] and simulated in this work.In shock tube experiments, the low-pressure section is usually designed long enough to gradually develop the generated shock wave into a plane shock wave, while in simulations, the generated shock wave is naturally planar. In addition, according to the basic theory of one-dimensional unsteady flow, the strength of the shock wave and the properties of the high-speed gas flow behind the shock wave depend only on the pressure and adiabatic properties of the gas of the high- and low-pressure sections rather than the length of the two sections. Thus, in this simulation, we designed the total length of the shock tube to be 1000 mm, as shown in Figure6. Moreover, in order to ignore the impact of the rarefactive waves reflected from the upstream end of the high-pressure section on the incident shock wave, a long high-pressure section is designed to avoid reflected rarefactive waves catching up with the incident shock wave. So, the high- and low-pressure sections are set to be 750 and 250 mm, respectively. Other gas properties of the high- and low-pressure sections are the same as in the experiments. Furthermore, to furtherly save computing time, we designed a slender computing domain of which the cross section is 6 mm × 6 mm. A nylon sphere that has a diameter Dp = 2 mm and density ρp = 1130 kg/m3 is positioned 30 mm from the downstream end, and the distance is the same as that in Devals et al.’s shock tube. The detailed material model of the nylon sphere is omitted here, which is similar to the sands model presented in Section 3.2.2. A gauge point is defined on the nylon sphere to record the particle motion information during the simulation.Figure 6
Computing domain and meshing withδ = 0.1Dp.


### 3.1.2. Validation for the Computation of Shock-Driven Particle Acceleration

The 1000 mm × 6 mm × 6 mm computing domain is discretized with uniform cubic mesh. Different mesh sizesδ = 1Dp, 0.5Dp, 0.3Dp, 0.2Dp, and 0.1Dp are used to assess the calculating ability of the Euler–Lagrange coupling in AUTODYN 3D for shock-driven particle acceleration. For instance, Figure 6 shows the enlarged local details with a grid size of 0.1Dp. Figure 5 gives the simulated velocity-time curves of nylon particles under these five grid conditions. These five simulation results are also compared with experimental data of Devals et al. and Igra et al.In Figure5, particle motions under different meshing conditions all undergo two processes: the acceleration process under the incident shock and the deceleration process as a result of the reflected shock. The origin time is set at the moment when the incident shock wave reaches the particle. At t = 0 ms, the particle starts to accelerate. Due to the long high-pressure section and the short low-pressure section, the particle acceleration process is not affected by the rarefactive waves reflected from the upstream end of the high-pressure section. The particle keeps accelerating with constant acceleration until the reflected shock wave from the downstream end of the tube reaches the particle. While the particle accelerations in Devals et al. and Igra et al.‘s experiments finally show a decline because of the rarefactive waves. After the acceleration stage, the particle motion comes into the decelerating stage. As the grid size decreases from 1Dp to 0.1Dp, the calculated velocity curve of the particles gradually approaches the experimental data.Since the spacing of the computation grid is much smaller than the size of the particle, the high-speed postshock gas around the particle is resolved, and the particle moves under the hydrodynamic force exerted by the resolved surrounding gas flow. Figure5 shows that the Euler–Lagrange coupling in AUTODYN 3D gives consistent drag forces to the experimental ones when the grid size is smaller than 0.3Dp. On the other hand, as the grid condition of 0.3Dp under a large number of particles has reached the maximum calculating power of our computer, the grid of 0.3Dp is chosen to accomplish the calculation. Moreover, based on Figure 5, the more fine grid condition would give a bit higher particle velocity because of the drag force, which would enhance the difference between the two velocities because of particle stress and gas drag.In conclusion, under the above meshing conditions, this detailed simulation method is a particle-resolved direct numerical simulation. Besides, dealing with object collisions is the strength of AUTODYN. So, although this fine meshing condition restricts the application of Euler–Lagrange coupling in AUTODYN 3D to a large number of particles, this direct numerical simulation method is useful for understanding more details of shock-particle interaction when the number of particles is relatively small.


## 3.1.1. Modeling of a Classic Single-Particle Accelerating Experiment

To validate the particle acceleration calculating ability of the Euler–Lagrange coupling in AUTODYN 3D, we chose a classic single-particle accelerating experiment in a shock tube. A similar validation has been conducted by Price et al. [28] in their shock-particles simulation work. Devals et al. [4] investigated the acceleration of a single sphere on a spider web driven by a planar shock wave in a horizontal shock tube. The tube has a 750 mm long high-pressure section and a 3000 mm long low-pressure section. The inner cross section of the tube is 80 mm × 80 mm. The downstream end of the low-pressure section is closed with a flange which reflects incident shock waves. A 2 mm diameter nylon sphere is hanged 30 mm from the downstream closing flange. The measured velocity of the single sphere in the unsteady flow agrees well with Igra and Takayama’s work, as shown in Figure 5. Igra and Takayama [3] just laid some small spheres on the shock tube floor to be tested.Figure 5
Particle velocities measured in literatures [3, 4] and simulated in this work.In shock tube experiments, the low-pressure section is usually designed long enough to gradually develop the generated shock wave into a plane shock wave, while in simulations, the generated shock wave is naturally planar. In addition, according to the basic theory of one-dimensional unsteady flow, the strength of the shock wave and the properties of the high-speed gas flow behind the shock wave depend only on the pressure and adiabatic properties of the gas of the high- and low-pressure sections rather than the length of the two sections. Thus, in this simulation, we designed the total length of the shock tube to be 1000 mm, as shown in Figure6. Moreover, in order to ignore the impact of the rarefactive waves reflected from the upstream end of the high-pressure section on the incident shock wave, a long high-pressure section is designed to avoid reflected rarefactive waves catching up with the incident shock wave. So, the high- and low-pressure sections are set to be 750 and 250 mm, respectively. Other gas properties of the high- and low-pressure sections are the same as in the experiments. Furthermore, to furtherly save computing time, we designed a slender computing domain of which the cross section is 6 mm × 6 mm. A nylon sphere that has a diameter Dp = 2 mm and density ρp = 1130 kg/m3 is positioned 30 mm from the downstream end, and the distance is the same as that in Devals et al.’s shock tube. The detailed material model of the nylon sphere is omitted here, which is similar to the sands model presented in Section 3.2.2. A gauge point is defined on the nylon sphere to record the particle motion information during the simulation.Figure 6
Computing domain and meshing withδ = 0.1Dp.


## 3.1.2. Validation for the Computation of Shock-Driven Particle Acceleration

The 1000 mm × 6 mm × 6 mm computing domain is discretized with uniform cubic mesh. Different mesh sizesδ = 1Dp, 0.5Dp, 0.3Dp, 0.2Dp, and 0.1Dp are used to assess the calculating ability of the Euler–Lagrange coupling in AUTODYN 3D for shock-driven particle acceleration. For instance, Figure 6 shows the enlarged local details with a grid size of 0.1Dp. Figure 5 gives the simulated velocity-time curves of nylon particles under these five grid conditions. These five simulation results are also compared with experimental data of Devals et al. and Igra et al.In Figure5, particle motions under different meshing conditions all undergo two processes: the acceleration process under the incident shock and the deceleration process as a result of the reflected shock. The origin time is set at the moment when the incident shock wave reaches the particle. At t = 0 ms, the particle starts to accelerate. Due to the long high-pressure section and the short low-pressure section, the particle acceleration process is not affected by the rarefactive waves reflected from the upstream end of the high-pressure section. The particle keeps accelerating with constant acceleration until the reflected shock wave from the downstream end of the tube reaches the particle. While the particle accelerations in Devals et al. and Igra et al.‘s experiments finally show a decline because of the rarefactive waves. After the acceleration stage, the particle motion comes into the decelerating stage. As the grid size decreases from 1Dp to 0.1Dp, the calculated velocity curve of the particles gradually approaches the experimental data.Since the spacing of the computation grid is much smaller than the size of the particle, the high-speed postshock gas around the particle is resolved, and the particle moves under the hydrodynamic force exerted by the resolved surrounding gas flow. Figure5 shows that the Euler–Lagrange coupling in AUTODYN 3D gives consistent drag forces to the experimental ones when the grid size is smaller than 0.3Dp. On the other hand, as the grid condition of 0.3Dp under a large number of particles has reached the maximum calculating power of our computer, the grid of 0.3Dp is chosen to accomplish the calculation. Moreover, based on Figure 5, the more fine grid condition would give a bit higher particle velocity because of the drag force, which would enhance the difference between the two velocities because of particle stress and gas drag.In conclusion, under the above meshing conditions, this detailed simulation method is a particle-resolved direct numerical simulation. Besides, dealing with object collisions is the strength of AUTODYN. So, although this fine meshing condition restricts the application of Euler–Lagrange coupling in AUTODYN 3D to a large number of particles, this direct numerical simulation method is useful for understanding more details of shock-particle interaction when the number of particles is relatively small.


## 3.2. Numerical Simulations of the Interaction between a Shock Wave and Packed Particles

On account of the validation in Section3.1.2, we safely employed AUTODYN 3D to simulate the interaction between a planar shock wave and packed particles. This section describes the geometric and material model of the simulation. Then, the simulation results of shock wave strength and particle motion are presented.


### 3.2.1. Geometric Model

Since the sidewall friction of the shock tube has no noticeable effect on particle acceleration [22], a slender computing domain is defined to save computing time. The length of the computing domain is the same as the length of the experimental shock tube, and the cross section of the domain is a small part of that of the real shock tube. The detailed dimension of the Euler domain is 1900 mm × 2.4 mm × 2.4 mm. On balance between the calculation accuracy discussed in Figure 5 and the computation time, the grid size is set to be 0.3Dp. The Euler mesh is 0.2 mm × 0.2 mm × 0.2 mm, as shown in Figure 7. The boundary condition of the downstream end of the Euler domain is “flow out.” Six fixed gauge points are defined within the Euler domain to monitor the gas pressure and velocity. They have the same x coordinate values as those of the real pressure transducers (PT1, …, PT6) mounted on the experimental shock tube, as described in Figure 1.Figure 7
Computational domain and grids.The sands used in the experiment were sieved by vibratory sieve shaker RETSCH AS200 and then analyzed by particle analyzer CAMSIZER P4. Based on the dynamic image analysis principle, CAMSIZER P4 measures the equivalent spherical diameter, which is the equivalent diameter of a sphere that has the same volume as the irregular particle. The measured particle size distribution of the sieved sands is of Gaussian distribution with a most frequent diameterDm of 0.646 mm and a variance of 0.038 mm, and the percentage of particles larger than 0.6 mm is 88.69%. For the above reasons, the spherical particle with a uniform diameter of 0.6 mm is used under the current simulation conditions.As AUTODYN 3D has a disadvantage on the complicated arrangement of large amounts of solid particles, the particles are initially arranged in the simple cubic packing model in this work. Figure8 illustrates the initial particle packing model used in the simulation. Nevertheless, the influence of this initial simple particle arrangement can be eliminated because of the following reason. When a planar shock wave impacts the particle wall, the particles inside the wall are further compacted and rearranged in both simulation and experiment. The porosity of the rearranged particles in the simulation will further reduce to approach the porosity of the real rearranged particles.Figure 8
Simple cubic packed particles.


### 3.2.2. Material Model

The high-pressure and low-pressure gases in the shock tube obey the ideal gas state equation. The pressureP is related to the density and energy by(1)P=γ−1ρe,where γ is adiabatic index, ρ is the air density, and e is the specific internal energy.The adiabatic index of the air is a constant 1.4. According to the standard constants from the AUTODYN material library, the density of the atmospheric air is set to be 1.225 × 10−3 g/cm3, and the specific internal energy of that is 2.068 × 105 kJ/kg. The density of the high-pressure air is set to be 4.677 × 10−2 g/cm3, and the specific internal energy of that is 2.138 × 105 kJ/kg, and then the pressure in the driver chamber will be 4 MPa.The SAND state equation embedded in AUTODYN is used to describe the mechanical properties of sands under dynamic loadings. The SAND state equation is given by Laine and Sandvik [29] as a piecewise linear curve with ten points, namely, density ρ as a function of pressure P. And the bulk elastic modulus K in the model is calculated by(2)K=ρc2.The parameters of the SAND state equation are listed in Table1.Table 1
Parameters of the SAND state equation.
Pressure (kPa)Density (g/cm3)Sound speed (m/s)Density (g/cm3)01.6742.65E + 021.6744.58E + 031.7398.52E + 021.7451.50E + 041.8741.72E + 032.0862.92E + 041.9971.88E + 032.1475.92E + 042.1442.27E + 032.39.81E + 042.252.96E + 032.5721.79E + 052.383.11E + 032.5982.89E + 052.4854.60E + 032.6354.50E + 052.5854.63E + 032.6416.51E + 052.674.63E + 032.8


### 3.2.3. Shock Wave Structure Compared with the Experiment Results

To compare with the experimental case in Section2.2 (P4 = 4 MPa, Dm = 0.646 mm, and L = 10 mm), the particle layers of simple cubic packing are set to be 16. Therefore, the length of the simple cubic packed particle wall is 9.6 mm. The pressure histories of PT3 and PT4 obtained from the simulations are compared with the experimental results in Figure 9. The time at which the incident shock wave reaches PT3 is set to be the origin time. Then, the reflected waves of the simulation and experiment return to PT3 almost at the same time. The transmitted shock waves of the simulation appear a little later than that of the experiment. The strengths of the incident shock waves of the simulation are equal to those of the experiment. The reflected and transmitted shock strengths of the simulation are almost equal to those of the experiment. In sum, the simulation which has high consistency with the experiment can be used to present much more detailed information about the shock-particle interaction.Figure 9
Shock wave structure comparison between simulation and experiment.


### 3.2.4. Initial Movement Pattern of the Packed Particles

By means of the Euler–Lagrange coupling method, the initial movement of the packed particles driven by a shock wave is presented in Figure10. Under the impact of the shock wave, the particles are compressed and rearranged. When the compressed particles integrally move downstream, the particles of the most downstream layer (i.e., the 16th layer) firstly separate from the whole compressed particles and then the particles of the 15th layer separate. Thus, the particles of the whole compressed particles separate and accelerate towards downstream layer by layer. Why do the particles move in such a pattern? The forces exerting on the particles and the velocity of the particles need further exploration in the next section.Figure 10
Simulated initial movement pattern of packed particles.


## 3.2.1. Geometric Model

Since the sidewall friction of the shock tube has no noticeable effect on particle acceleration [22], a slender computing domain is defined to save computing time. The length of the computing domain is the same as the length of the experimental shock tube, and the cross section of the domain is a small part of that of the real shock tube. The detailed dimension of the Euler domain is 1900 mm × 2.4 mm × 2.4 mm. On balance between the calculation accuracy discussed in Figure 5 and the computation time, the grid size is set to be 0.3Dp. The Euler mesh is 0.2 mm × 0.2 mm × 0.2 mm, as shown in Figure 7. The boundary condition of the downstream end of the Euler domain is “flow out.” Six fixed gauge points are defined within the Euler domain to monitor the gas pressure and velocity. They have the same x coordinate values as those of the real pressure transducers (PT1, …, PT6) mounted on the experimental shock tube, as described in Figure 1.Figure 7
Computational domain and grids.The sands used in the experiment were sieved by vibratory sieve shaker RETSCH AS200 and then analyzed by particle analyzer CAMSIZER P4. Based on the dynamic image analysis principle, CAMSIZER P4 measures the equivalent spherical diameter, which is the equivalent diameter of a sphere that has the same volume as the irregular particle. The measured particle size distribution of the sieved sands is of Gaussian distribution with a most frequent diameterDm of 0.646 mm and a variance of 0.038 mm, and the percentage of particles larger than 0.6 mm is 88.69%. For the above reasons, the spherical particle with a uniform diameter of 0.6 mm is used under the current simulation conditions.As AUTODYN 3D has a disadvantage on the complicated arrangement of large amounts of solid particles, the particles are initially arranged in the simple cubic packing model in this work. Figure8 illustrates the initial particle packing model used in the simulation. Nevertheless, the influence of this initial simple particle arrangement can be eliminated because of the following reason. When a planar shock wave impacts the particle wall, the particles inside the wall are further compacted and rearranged in both simulation and experiment. The porosity of the rearranged particles in the simulation will further reduce to approach the porosity of the real rearranged particles.Figure 8
Simple cubic packed particles.


## 3.2.2. Material Model

The high-pressure and low-pressure gases in the shock tube obey the ideal gas state equation. The pressureP is related to the density and energy by(1)P=γ−1ρe,where γ is adiabatic index, ρ is the air density, and e is the specific internal energy.The adiabatic index of the air is a constant 1.4. According to the standard constants from the AUTODYN material library, the density of the atmospheric air is set to be 1.225 × 10−3 g/cm3, and the specific internal energy of that is 2.068 × 105 kJ/kg. The density of the high-pressure air is set to be 4.677 × 10−2 g/cm3, and the specific internal energy of that is 2.138 × 105 kJ/kg, and then the pressure in the driver chamber will be 4 MPa.The SAND state equation embedded in AUTODYN is used to describe the mechanical properties of sands under dynamic loadings. The SAND state equation is given by Laine and Sandvik [29] as a piecewise linear curve with ten points, namely, density ρ as a function of pressure P. And the bulk elastic modulus K in the model is calculated by(2)K=ρc2.The parameters of the SAND state equation are listed in Table1.Table 1
Parameters of the SAND state equation.
Pressure (kPa)Density (g/cm3)Sound speed (m/s)Density (g/cm3)01.6742.65E + 021.6744.58E + 031.7398.52E + 021.7451.50E + 041.8741.72E + 032.0862.92E + 041.9971.88E + 032.1475.92E + 042.1442.27E + 032.39.81E + 042.252.96E + 032.5721.79E + 052.383.11E + 032.5982.89E + 052.4854.60E + 032.6354.50E + 052.5854.63E + 032.6416.51E + 052.674.63E + 032.8


## 3.2.3. Shock Wave Structure Compared with the Experiment Results

To compare with the experimental case in Section2.2 (P4 = 4 MPa, Dm = 0.646 mm, and L = 10 mm), the particle layers of simple cubic packing are set to be 16. Therefore, the length of the simple cubic packed particle wall is 9.6 mm. The pressure histories of PT3 and PT4 obtained from the simulations are compared with the experimental results in Figure 9. The time at which the incident shock wave reaches PT3 is set to be the origin time. Then, the reflected waves of the simulation and experiment return to PT3 almost at the same time. The transmitted shock waves of the simulation appear a little later than that of the experiment. The strengths of the incident shock waves of the simulation are equal to those of the experiment. The reflected and transmitted shock strengths of the simulation are almost equal to those of the experiment. In sum, the simulation which has high consistency with the experiment can be used to present much more detailed information about the shock-particle interaction.Figure 9
Shock wave structure comparison between simulation and experiment.


## 3.2.4. Initial Movement Pattern of the Packed Particles

By means of the Euler–Lagrange coupling method, the initial movement of the packed particles driven by a shock wave is presented in Figure10. Under the impact of the shock wave, the particles are compressed and rearranged. When the compressed particles integrally move downstream, the particles of the most downstream layer (i.e., the 16th layer) firstly separate from the whole compressed particles and then the particles of the 15th layer separate. Thus, the particles of the whole compressed particles separate and accelerate towards downstream layer by layer. Why do the particles move in such a pattern? The forces exerting on the particles and the velocity of the particles need further exploration in the next section.Figure 10
Simulated initial movement pattern of packed particles.


## 4. Discussions of the Initial Motion Mechanism

Based on the above experimental motion measurement and simulated motion pattern, we noticed that packed particles do not show notable motion before the shock wave transmits the packed particles. Britan et al. [30] studied the motion of a single sphere in a shock tube employing the schlieren technique. They observed similar phenomena: the single sphere keeps stationary when the incident shock wave passes through it; then, the sphere begins to move after all the transient disturbances vanish; finally, the sphere reaches a relatively stable acceleration. So, in this section, we divided the initial motion of the packed particles into two stages to investigate the quasi-static stage and the acceleration stage.


### 4.1. The Quasi-Static Stage

Tens microseconds after the shock wave impacts on the packed particles, all the particles visually remain still, while shock waves and stress waves propagate.


#### 4.1.1. Shock Wave Propagation

In order to demonstrate the gas flow in the process of the interaction between the shock wave and the packed particles, the AUTODYN multimaterial Euler model is used here. The atmospheric air in the driven section of the shock tube is artificially simulated by two gases, which have two different colors but the same parameters as shown at the 0 ms of Figure11(a). In the following time, the red gas in front of the 1st particle layer is pushed into the pores between the particle skeleton and emerges from the last particle layer at 1.94 ms.Figure 11
Gas percolation and pressure distribution in the particle pores.
(a)(b)Figure11(b) shows the gas pressure distribution in the particle pores. Immediately after the incident shock wave impacts on the first particle layer, the reflected shock wave arises. At the same time, the shock wave also enters into the pores between the particles skeleton. Tens microseconds later, the shock wave is attenuated, and the transmitted shock wave recreates at the last particle layer. The strengths of the incident, reflected, and transmitted shock are recorded and compared with the experimental data measured by the pressure transducers in the same x positions in Figure 9.In this shock wave propagating process, the particles remain stationary. However, in the pores between the particles skeleton, the pressure gradient is established when the high-speed gas percolates through the packed particles, as shown in Figure11(b).


#### 4.1.2. Particle Stress Delivery

Due to the complexity of particle stress measurement, it is difficult to measure the particle stress during the interaction between the shock wave and the packed particles experimentally. In this section, the numerical method is used to calculate the particle stress. Figure12(a) displays the stress distribution among the particle layers, and Figure 12(b) gives the stress histories of every particle layer. # 1, # 2, …, and # 16 are 16 points in each particle layer which are used to monitor the particle velocity and stress. The points PT3 and PT4 have the same positions as the actual pressure transducers installed in the experimental shock tube. They are used to record the incident and transmitted shock waves at the two sides of the particle wall.Figure 12
Stress distribution on the packed particles and stress histories.
(a)(b)As the direct collision of the packed particles, the impact of the incident shock wave stimulates particle stress, which is delivered downstream on the particle skeleton. After the first stress wave goes downstream, the new stress wave arises at the upstream front and propagates downstream. The velocity of the first stress wave is about 9.6 mm/0.06 ms = 160 m/s.


#### 4.1.3. Particle Velocity during the Quasi-Static Stage

Although the particles do not show notable motion in this stage, the particle collision still causes particle acceleration because of the particle stress wave delivery. Figures13(a) and 13(b) are the stress and velocity histories of every particle layer in the time interval between 1.845 ms and 1.925 ms. At 1.845 ms, the incident shock wave impacts upon the upstream front of the densely packed particles, and at 1.925 ms, the first particle stress wave leaves the downstream front of the particles. During the interval, each particle layers successively accelerate and obtain an axial velocity of 6–8 m/s. As the velocity of the particles is too small to be recognized, this stage can accurately be called a quasi-static stage.Figure 13
Stress and velocity histories of every particle layers.
(a)(b)


### 4.2. The Acceleration Stage in the Subsequent Time

After the quasi-static stage, the particles start to accelerate significantly.


#### 4.2.1. Flow Field Pressure

Because of the push of the high-pressure gas blocked by the densely packed particles and the propagation of the particle stress wave, the particles are further compressed and rearranged, as shown in Figure14. The particles integrally move downstream; at the same time, the downstream particles begin to separate and accelerate. Figure 14 also presents the flow field pressure among the particle layers. As the high-pressure gas blocked by the densely packed particles percolates into the pores among the particle skeleton, the pressure gradient establishes in the pores. Driven by the pressure gradient in the pores, the particles of the most downstream layer (16th layer) are separated and accelerated and then the particles of the 15th layer. In the figure, the pressure gradient among particle layers exists all the time. Thus, the pressure gradient of the flow field among the particle layers is the direct reason for the notable acceleration of the downstream particles.Figure 14
Gas pressure gradient among the particle layers.


#### 4.2.2. Particle Velocities during the Acceleration Stage

Figure15 is the axial velocity histories of each particle layer (#1, #2, …, #16). The first red dotted line (at 1.845 ms) indicates the time that the incident shock wave impacts on the upstream front of the packed particles. The second red dotted line (at 1.925 ms) indicates the time that the first stress wave leaves the downstream front of the particles. Between two red dotted lines, each particle layer successively obtains an axial velocity of 6–8 m/s due to the stress wave among compressed particles. Then, each particle layer maintains the velocity for several microseconds. After that, each particle layer is significantly accelerated once again due to the gas pressure gradient. At 2.7 ms (0.855 ms later than the interaction beginning), each particle layer obtains an axial velocity of 60–70 m/s. Comparing with the experimental results in Figure 4, at t  = 167.42 ms (0.86 ms later than its interaction beginning), the velocity of the downstream front is 62.96 m/s.Figure 15
Particle velocity during the acceleration stage.
(a)(b)


## 4.1. The Quasi-Static Stage

Tens microseconds after the shock wave impacts on the packed particles, all the particles visually remain still, while shock waves and stress waves propagate.


### 4.1.1. Shock Wave Propagation

In order to demonstrate the gas flow in the process of the interaction between the shock wave and the packed particles, the AUTODYN multimaterial Euler model is used here. The atmospheric air in the driven section of the shock tube is artificially simulated by two gases, which have two different colors but the same parameters as shown at the 0 ms of Figure11(a). In the following time, the red gas in front of the 1st particle layer is pushed into the pores between the particle skeleton and emerges from the last particle layer at 1.94 ms.Figure 11
Gas percolation and pressure distribution in the particle pores.
(a)(b)Figure11(b) shows the gas pressure distribution in the particle pores. Immediately after the incident shock wave impacts on the first particle layer, the reflected shock wave arises. At the same time, the shock wave also enters into the pores between the particles skeleton. Tens microseconds later, the shock wave is attenuated, and the transmitted shock wave recreates at the last particle layer. The strengths of the incident, reflected, and transmitted shock are recorded and compared with the experimental data measured by the pressure transducers in the same x positions in Figure 9.In this shock wave propagating process, the particles remain stationary. However, in the pores between the particles skeleton, the pressure gradient is established when the high-speed gas percolates through the packed particles, as shown in Figure11(b).


### 4.1.2. Particle Stress Delivery

Due to the complexity of particle stress measurement, it is difficult to measure the particle stress during the interaction between the shock wave and the packed particles experimentally. In this section, the numerical method is used to calculate the particle stress. Figure12(a) displays the stress distribution among the particle layers, and Figure 12(b) gives the stress histories of every particle layer. # 1, # 2, …, and # 16 are 16 points in each particle layer which are used to monitor the particle velocity and stress. The points PT3 and PT4 have the same positions as the actual pressure transducers installed in the experimental shock tube. They are used to record the incident and transmitted shock waves at the two sides of the particle wall.Figure 12
Stress distribution on the packed particles and stress histories.
(a)(b)As the direct collision of the packed particles, the impact of the incident shock wave stimulates particle stress, which is delivered downstream on the particle skeleton. After the first stress wave goes downstream, the new stress wave arises at the upstream front and propagates downstream. The velocity of the first stress wave is about 9.6 mm/0.06 ms = 160 m/s.


### 4.1.3. Particle Velocity during the Quasi-Static Stage

Although the particles do not show notable motion in this stage, the particle collision still causes particle acceleration because of the particle stress wave delivery. Figures13(a) and 13(b) are the stress and velocity histories of every particle layer in the time interval between 1.845 ms and 1.925 ms. At 1.845 ms, the incident shock wave impacts upon the upstream front of the densely packed particles, and at 1.925 ms, the first particle stress wave leaves the downstream front of the particles. During the interval, each particle layers successively accelerate and obtain an axial velocity of 6–8 m/s. As the velocity of the particles is too small to be recognized, this stage can accurately be called a quasi-static stage.Figure 13
Stress and velocity histories of every particle layers.
(a)(b)


## 4.1.1. Shock Wave Propagation

In order to demonstrate the gas flow in the process of the interaction between the shock wave and the packed particles, the AUTODYN multimaterial Euler model is used here. The atmospheric air in the driven section of the shock tube is artificially simulated by two gases, which have two different colors but the same parameters as shown at the 0 ms of Figure11(a). In the following time, the red gas in front of the 1st particle layer is pushed into the pores between the particle skeleton and emerges from the last particle layer at 1.94 ms.Figure 11
Gas percolation and pressure distribution in the particle pores.
(a)(b)Figure11(b) shows the gas pressure distribution in the particle pores. Immediately after the incident shock wave impacts on the first particle layer, the reflected shock wave arises. At the same time, the shock wave also enters into the pores between the particles skeleton. Tens microseconds later, the shock wave is attenuated, and the transmitted shock wave recreates at the last particle layer. The strengths of the incident, reflected, and transmitted shock are recorded and compared with the experimental data measured by the pressure transducers in the same x positions in Figure 9.In this shock wave propagating process, the particles remain stationary. However, in the pores between the particles skeleton, the pressure gradient is established when the high-speed gas percolates through the packed particles, as shown in Figure11(b).


## 4.1.2. Particle Stress Delivery

Due to the complexity of particle stress measurement, it is difficult to measure the particle stress during the interaction between the shock wave and the packed particles experimentally. In this section, the numerical method is used to calculate the particle stress. Figure12(a) displays the stress distribution among the particle layers, and Figure 12(b) gives the stress histories of every particle layer. # 1, # 2, …, and # 16 are 16 points in each particle layer which are used to monitor the particle velocity and stress. The points PT3 and PT4 have the same positions as the actual pressure transducers installed in the experimental shock tube. They are used to record the incident and transmitted shock waves at the two sides of the particle wall.Figure 12
Stress distribution on the packed particles and stress histories.
(a)(b)As the direct collision of the packed particles, the impact of the incident shock wave stimulates particle stress, which is delivered downstream on the particle skeleton. After the first stress wave goes downstream, the new stress wave arises at the upstream front and propagates downstream. The velocity of the first stress wave is about 9.6 mm/0.06 ms = 160 m/s.


## 4.1.3. Particle Velocity during the Quasi-Static Stage

Although the particles do not show notable motion in this stage, the particle collision still causes particle acceleration because of the particle stress wave delivery. Figures13(a) and 13(b) are the stress and velocity histories of every particle layer in the time interval between 1.845 ms and 1.925 ms. At 1.845 ms, the incident shock wave impacts upon the upstream front of the densely packed particles, and at 1.925 ms, the first particle stress wave leaves the downstream front of the particles. During the interval, each particle layers successively accelerate and obtain an axial velocity of 6–8 m/s. As the velocity of the particles is too small to be recognized, this stage can accurately be called a quasi-static stage.Figure 13
Stress and velocity histories of every particle layers.
(a)(b)


## 4.2. The Acceleration Stage in the Subsequent Time

After the quasi-static stage, the particles start to accelerate significantly.


### 4.2.1. Flow Field Pressure

Because of the push of the high-pressure gas blocked by the densely packed particles and the propagation of the particle stress wave, the particles are further compressed and rearranged, as shown in Figure14. The particles integrally move downstream; at the same time, the downstream particles begin to separate and accelerate. Figure 14 also presents the flow field pressure among the particle layers. As the high-pressure gas blocked by the densely packed particles percolates into the pores among the particle skeleton, the pressure gradient establishes in the pores. Driven by the pressure gradient in the pores, the particles of the most downstream layer (16th layer) are separated and accelerated and then the particles of the 15th layer. In the figure, the pressure gradient among particle layers exists all the time. Thus, the pressure gradient of the flow field among the particle layers is the direct reason for the notable acceleration of the downstream particles.Figure 14
Gas pressure gradient among the particle layers.


### 4.2.2. Particle Velocities during the Acceleration Stage

Figure15 is the axial velocity histories of each particle layer (#1, #2, …, #16). The first red dotted line (at 1.845 ms) indicates the time that the incident shock wave impacts on the upstream front of the packed particles. The second red dotted line (at 1.925 ms) indicates the time that the first stress wave leaves the downstream front of the particles. Between two red dotted lines, each particle layer successively obtains an axial velocity of 6–8 m/s due to the stress wave among compressed particles. Then, each particle layer maintains the velocity for several microseconds. After that, each particle layer is significantly accelerated once again due to the gas pressure gradient. At 2.7 ms (0.855 ms later than the interaction beginning), each particle layer obtains an axial velocity of 60–70 m/s. Comparing with the experimental results in Figure 4, at t  = 167.42 ms (0.86 ms later than its interaction beginning), the velocity of the downstream front is 62.96 m/s.Figure 15
Particle velocity during the acceleration stage.
(a)(b)


## 4.2.1. Flow Field Pressure

Because of the push of the high-pressure gas blocked by the densely packed particles and the propagation of the particle stress wave, the particles are further compressed and rearranged, as shown in Figure14. The particles integrally move downstream; at the same time, the downstream particles begin to separate and accelerate. Figure 14 also presents the flow field pressure among the particle layers. As the high-pressure gas blocked by the densely packed particles percolates into the pores among the particle skeleton, the pressure gradient establishes in the pores. Driven by the pressure gradient in the pores, the particles of the most downstream layer (16th layer) are separated and accelerated and then the particles of the 15th layer. In the figure, the pressure gradient among particle layers exists all the time. Thus, the pressure gradient of the flow field among the particle layers is the direct reason for the notable acceleration of the downstream particles.Figure 14
Gas pressure gradient among the particle layers.


## 4.2.2. Particle Velocities during the Acceleration Stage

Figure15 is the axial velocity histories of each particle layer (#1, #2, …, #16). The first red dotted line (at 1.845 ms) indicates the time that the incident shock wave impacts on the upstream front of the packed particles. The second red dotted line (at 1.925 ms) indicates the time that the first stress wave leaves the downstream front of the particles. Between two red dotted lines, each particle layer successively obtains an axial velocity of 6–8 m/s due to the stress wave among compressed particles. Then, each particle layer maintains the velocity for several microseconds. After that, each particle layer is significantly accelerated once again due to the gas pressure gradient. At 2.7 ms (0.855 ms later than the interaction beginning), each particle layer obtains an axial velocity of 60–70 m/s. Comparing with the experimental results in Figure 4, at t  = 167.42 ms (0.86 ms later than its interaction beginning), the velocity of the downstream front is 62.96 m/s.Figure 15
Particle velocity during the acceleration stage.
(a)(b)


## 5. Conclusions

The initial movement of the densely packed particles driven by a planar shock wave has been investigated in a horizontal multiphase shock tube. To find out the main mechanism for the densely packed particles state to move and scatter, the interaction between the densely packed particles and the shock wave is simulated by an Euler–Lagrangian coupling method.The following basic conclusions can be drawn from the experiment and simulation. Firstly, due to the particle stress waves delivered downstream on the particle skeleton, each particle layer successively accelerates and obtains an axial velocity of the order of a few meters per second. Secondly, the gas flow in the pores among the particle skeleton drags the downstream particles to further accelerate to the velocity of the order of a few ten meters per second. The difference between the two velocities caused by the particle stress and gas drag proves that the gas drag is the main reason for the particle acceleration. Combined with the results of previous experiments, it can be inferred that the different porosities because of different particle sizes, Mach numbers, and wall thicknesses will not change the order of magnitude of the two velocities.The initial movement of the densely packed particles driven by the shock wave can be divided into two stages. In the quasi-static stage (tens microseconds), after the impact of the incident shock wave on the densely packed particles, the stress waves on the particle skeleton give the particles an unremarkable axial velocity. In the acceleration stage, because of the percolation of the blocked high-pressure gas into the pores among the particle skeleton, the gas flow in the pores drags the downstream particles to separate and accelerate dramatically. At the same time, the upstream particles integrally move downstream due to the push of the blocked high-pressure gas. Therefore, it is the gas push-drag dual mechanism that transforms the densely packed particles into a dense gas-particle flow, which behaves as the typical expansion phenomenon of the particle curtain or particle wall.



---

*Source: 8867615-2021-04-19.xml*