Surface deformation under overlapping impacts of solid particles

To improve the accuracy of erosion prediction, the effect of subsequent particles impacting the same area while the first single particle rebounds from the substrate must be considered. This issue has rarely been considered in studies pertaining to erosion damage. In the present study, the ABAQUS software is used to investigate the erosion crater morphology and stress distribution on a target material subjected to overlapping impacts of spherical particles. Subsequently, the validated model is applied to investigate the effect of the overlapping impacts of particles on the target. Accordingly, the correlation between erosion severity and the impact locations of the two incident particles is quantified. The results show that the horizontal distance between two solid particle impact locations can significantly affect the erosion magnitude and pattern. The interactions of the resulting craters diminish when the horizontal distance exceeds 0.6 times the particle diameter. When the horizontal distance is approximately 0.06 times the particle diameter, the energy loss originating from collisions reaches the maximum, which modifies the crater morphology. The present study is expected to provide in-depth insights into erosion mechanisms and erosion modeling.


Introduction
Studies have shown that erosion is one of the most significant causes of material damage and equipment failure [1][2][3]. Erosion is a typical form of damage to oil and gas pipelines [4][5][6]. Solid particles transported by fluid impact the pipe wall, which may cause local damages such as small cracks, plastic deformation, and material shedding to the inner wall of the pipe. The continuous impact of particles will cause wall damage to accumulate gradually until the material fails, thereby resulting in significant economic losses and safety issues. Previous investigations indicated that erosion is a complex problem that is affected by particle velocity, angular velocity, impact angle, particle shape, and particle size. Theories based on experiments and simulations have been proposed to model and predict material damage [7]. To accurately predict and manage erosive events in industrial processes, the erosion process must be investigated comprehensively, including the examination of the impact process and its relevant affecting factors, such as the particle diameter, impact kinetics, and energy conversion. The objectives are to determine the relationship between the aforementioned factors and the impact damage as well as to further generalize a user-friendly model to predict future events. In this regard, most of the currently available theoretical models focus on the particle shape, impact velocity, and impact angle [8][9][10][11], including the micro-cutting model [12], deformation wear theory [13,14], plastic deformation accumulation theory [15][16][17], and adiabatic shear theory [18]; these models reveal the effect of abrasive particles on erosion and wear mechanisms. Constant related to velocity of particle before collision K 0 Empirical constant related to material properties m Thermal softening of material m p Quality of abrasive particle N c Effective number of turns, which is equal to the total number of turns minus two n 1 Empirical constant related to material properties n 2 Empirical constant related to impact velocity P Pressure on surface of material Q p Formation of permanent indentation, which requires a certain amount of energy R Notch radius of smallest cross-section of pattern T m Melting temperature T r Reference temperature T * Normalized temperature t Duration of impact process u Horizontal distance between the two impact events v p Average impact velocity of particles v 0 Impact velocity v' Velocity of particle after impact W D Wear due to repeated deformation W t Initial kinetic energy of particles W u Increased internal energy of target material α Impact angle of particles ε Equivalent plastic strain rate ε B Deformation wear factor ε el Elastic stress ε f Fracture strain ε pl Equivalent plastic strain ε 0 Reference strain rate η Energy conversion rate σ Static yield stress σ c Total elastic stress σ u Elastic stress σ * Stress triaxiality Researchers have experimentally investigated different impact particles and target materials. Hutchings and Winter [19] investigated the impact of a single square particle on a low-carbon steel specimen.
Different impact parameters such as the impact angle and particle azimuth were investigated, and it was discovered that the impact process comprised two stages: forward spin impact and backspin impact. Sahoo et al. [20] performed sandblasting experiments to statistically evaluate the behavior of particles impacting the target material, as well as investigated the effects of different parameters, including microstructural degradation, impact velocity, impingement angle, and erodent size on the material. Tirupataiah et al. [21] incorporated the effects of hardness, elastic modulus, and impact velocity into a theoretical model and established a functional correlation for the crater recovery degree. Papini and Dhar [22,23] used a pneumatic launcher to investigate the erosion mechanism of angular particles based on an impact velocity of less than 100 m/s. Despite the achievements of the study, it was limited to the "rigid-plastic" framework.
Generally, experiments provide limited information, such as the impact velocity and crater depth, which cannot provide adequate information for analyzing dynamic erosion. Meanwhile, studies have indicated that the finite element method [24][25][26][27][28][29] can provide detailed information, including stress, interaction force, and crater morphology of particles and target materials [30][31][32][33]. In addition, it can be used to effectively investigate the effect of the abovementioned parameters on the target material. Many scholars [34][35][36][37][38] applied the finite element method to establish an erosion calculation model for angular particles that impact viscoplastic materials. The results obtained were verified via collision tests using an air gun. The effects of the impact angle and azimuth angle of the particles on the deformation mechanism of the target material were reflected in their numerical simulation and experiment. Previously, material removal [39,40] and accumulation deformation [41][42][43][44][45] mechanisms caused by abrasive particles have been investigated comprehensively. Most previous investigations focused on non-overlapping particle impacts on the target material, i.e., erosion modeling based on singleparticle impact [46][47][48]. Figure 1(a) shows a schematic illustration of a completely non-overlapping impact with an extremely large horizontal distance, which hindered the crater interaction and its effect on the erosion profile from being considered.
The ultimate goal of this study is to quantify the neighboring impact of particles on erosion, which can be incorporated to the erosion model based on the single-particle impact. Therefore, we assumed that the particle impact velocity can be quantified accurately. Under this condition, the general form of a single-particle erosion model is typically expressed as Eq. (1): where ER is the ratio of the mass loss of the target material to the mass of the erosion particles; F s is a coefficient related to the shape of the particles; BH is Brinell hardness of the material; K 0 and n 1 are the empirical constants related to the material properties; and v p is the average impact velocity of the particles. Moreover, n 1 is an empirical constant related to the material properties; n 2 is an empirical constant related to the impact velocity; and F(α) is a function of the impact angle. As shown, the interaction between adjacent craters is rarely considered. Therefore, in this study, a numerical model was established to investigate the overlapping impact of spherical particles, as shown in Fig. 1(b). The complex particle impact scenario was simplified to the impact of two neighboring particles. The expansion and   www.Springer.com/journal/40544 | Friction addition of an impacting particle couplet were assumed. Therefore, it is critical to first investigate erosion damage by two overlapping impact particles comprehensively to understand the trend and results quantitatively. Different impact parameters, including the impact velocity, impact angle, and horizontal distance between two impacts, were considered to investigate the resulting erosion profile under multiparticle impacts. The present study is expected to provide a scientific basis for interpreting material erosion originating from multi-particle impacts as well as for elucidating the erosion mechanism.

Experimental device
The experimental setup used in this study was designed and manufactured based on the slingshot ejection device proposed by Dhar et al. [49]. Figure 2 shows the configuration of the test setup. During the experiment, the ejection rod was manually loaded by pressing the launching arm down until it was securely locked in place to the release mechanism (rubber arm stop). The ejection rod was connected to two springs in series and was limited by an adjustable plate. Moreover, the particles were fixed by a bracket (particle holder). The target fixing plate was fixed on the frame of the ejection rod, and the target frame (target holder) was adjustable to change the incident angle.
The ejector rod pressed the spring down to convert the elastic potential energy of the spring into the kinetic energy of the particles. It is noteworthy that different ejection speeds can be obtained by changing where K is the stiffness coefficient representing the attribute of the spring, N/mm; G denotes the rigidity modulus of the wire, N/mm 2 ; d w is the wire diameter, mm; D m is the middle diameter, which is equal to the outer diameter minus the wire diameter, mm; and Nc is the effective number of turns, which is equal to the total number of turns minus 2. The specific parameters of the spring are shown in Table 1.

Experimental process
In the experiments, particles were synthesized using high-speed steel (W18Cr4V), and a target material composed of 6061 aluminum alloy (AA6061) were used. The target material was selected to achieve the excellent plastic deformation performance. Therefore, the hardness of the target material should be lower than that of the particles. In this case, the particles can generate pits on the target material after the impact. Prior to the experiment, the target material was polished to obtain a smooth surface. The friction coefficient between the two materials was 0.1. The simulation results were consistent with the experimental results. Table 2 presents the main parameters of the particles and target materials. Confocal microscopy (MicroPhase, PhaseView, France) was performed in the present study to scan the target material. Table 3 lists the experimental data. Meanwhile, Fig. 3 illustrates the morphology of the crater when the target was impacted by particles of 5 mm in diameter, where the impact angle and velocity were 45° and 20 m/s, respectively. The red   and yellow areas indicate the non-impacted surface and the crater that is deeper than the non-impacted area, respectively.

Numerical simulation
Particle impact is a complex phenomenon. Accordingly, several theories have been proposed to model contact and damage mechanics [50,51]. The Johnson-Cook (JC) viscoplastic model is used to model the flow stress behavior of typical metals subjected to large strains [52,53]. To reduce the computational cost, a two-dimensional (2D) model was applied to the calculations involved. In fact, the same calculations were performed by Takaffoli and Papini [54]. The modelling details are provided in the following sections.

Constitutive models
The strain of the target material changes significantly during dynamic erosion. To characterize materials at higher strain rates, the JC constitutive model, which is suitable for ductile metal materials [55][56][57][58], is used to characterize the nonlinear dynamic mechanical behavior of the materials. The present study focuses on ductile material erosion behavior. Therefore, the material constants for the JC constitutive model of AA6061 were obtained.

JC hardening criterion
Johnson and Cook [52] combined viscoplastic and continuum damage mechanics, and proposed an empirical constitutive model under a large deformation, a high strain rate, and thermal insulation [59,60]. In this model, the yield flow stress of the material can be expressed as a function of the strain, strain rate, and temperature, as Eqs. (3) and (4): where σ and ε pl are the static yield stress and equivalent plastic strain, respectively. Moreover, A, B, and n denote the material parameters at the transition temperature; C and m are the strain rate and thermal softening of the material, respectively; ε is the equivalent plastic strain rate; ε 0 is the reference strain rate. Finally, T, T * , T m , and T r denote the operating, normalized, melting, and reference temperatures, respectively. The defined parameters were acquired via tensile and split Hopkinson bar tests on AA6061. The parameters of the JC hardening criterion are listed in Table 4. Table 4 Parameters of JC hardening criterion.

JC damage criterion
The JC damage model is based on the cavitation growth and damage evolution principles of the material. In this model, the effects of stress triaxiality, strain rate, and temperature softening are considered. The failure criterion and calculation equations for the JC model are expressed as Eqs. (5)- (7): where ε f and σ * denote the fracture strain and stress triaxiality, respectively; D 1 , D 2 , D 3 , D 4 , and D 5 are the damage constants; a and R denote the radius and notch radius of the smallest cross-section of the pattern, respectively; d 0 and d 1 are the initial and final diameters after deformation, respectively. It is noteworthy that constants D 1 -D 5 can be obtained using Eq. (5). Moreover, Eqs. (6) and (7) show the expressions for obtaining σ * and ε f , respectively. Table 5 lists the values of these parameters. Table 5 Parameters of JC damage criterion [61].

Boundary conditions
In all the simulations, it was assumed that the particles were rigid bodies. In other words, the particle strain was set to zero during the impact. The initial velocity was applied to the particles through a predefined field, and the input velocity was decomposed along the X-and Y-directions. The movements of both the left and right sides were restricted along the X-direction, and the bottom plane was fixed. Figure 4 shows the imposed boundary conditions.

Mesh independence test
In numerical simulations, the use of a finer grid yields more accurate results, although at the expense of higher computational costs. In this study, a systematic grid independence test was performed to obtain a reasonable mesh size to ensure accuracy while maintaining a reasonable computational cost. In this regard, the von Mises stress and maximum depth were considered for selecting the optimized mesh. As shown in Table 6, five mesh sizes were considered in the grid-independence test. It is noteworthy that quadrilateral structural meshes were applied in all the simulations. Figure 5(a) shows significant deviations in the von Mises stress contours between the results obtained from Test Nos. 1 and 2, indicating that the meshes were coarse. Table 6 shows that for Test Nos. 1-3, the von Mises stress near the deformation ranged from 110 to 160 MPa after the second impact, whereas the von Mises stresses in Test Nos. 4 and 5 were 172.8 and 173.2 MPa, respectively, and their distributions were similar. As shown in Fig. 5(b), for Test Nos. 1-3, the depths of the target material ranged from 30 to 65 μm after the   Figure 6 shows the mesh generated for the 2D target material. A quadrilateral structural mesh was applied in the simulations. It is noteworthy that CPS4R is a four-node bilinear plane stress quadrilateral element with reduced integration. Therefore, only one integral point exists in the center of the unit and hourglass control. Consequently, the element can resist deformation without involving the stiffness. In this type of unit, the analysis accuracy is affected less when the mesh is deformed, and the displacement can be solved accurately. To improve the calculation accuracy and conserve computational   www.Springer.com/journal/40544 | Friction resources, only finer meshes were applied in the contact area, whereas coarse meshes were applied in other areas. The total number of elements in the refined and fine meshes in the contact region were 62,436 and 55,611, respectively.

Experimental verification
In  differences between the simulation and experiment results of the crater were 5.0 and 35.1 μm, respectively, which were only 8% and 3% greater compared with the experimental data. Therefore, the simulated value was consistent with the experimental value. It is noteworthy that the roughness data shown in Fig. 7 are expressed in units of micrometers, which are three orders of magnitude smaller compared with the units of the particle; therefore, they are negligible. Moreover, the effect of roughness is not considered in erosion models. Therefore, scholars typically characterize the effect of roughness using the friction coefficient [62,63]. Figures 8 and 9 show the penetration depths and lengths under different incident velocities, respectively. The experiments were repeated five times for each group. The results show that the test values of the five groups (involving different speeds) agreed well with the numerical values obtained from the ABAQUS software. When the speed was 22 m/s, the maximum difference between the penetration depths was approximately 8 μm, indicating a 6% error compared with the experiment. Additionally, when the speed was 16 m/s, the maximum difference between the pitting lengths was approximately 49.0 μm, indicating a 4% error compared with the experiment.
A 2D model is a simplified model. However, as verified using data, it can capture the prediction of the maximum erosion depth and profile within the   cross section. Owing to its relative accuracy, ability to represent the general trend, and reduced computational cost, it was used to present the results of the current study, similar to other previous studies [28,40]. Figure 10 shows a comparison between the experimental results and simulation data for the overlapping impacts of spherical particles. In this experiment, the particle size and horizontal distance of the impact were 5 and 1 mm, respectively. Moreover, the impact velocity and angle were 20 m/s and 85°, respectively. After the two impact events, the depths of the crater yielded by the first and second impact events were 87.6 and 97.0 μm, respectively; moreover, the depths obtained from the simulation results were 82 and 91 μm, respectively. Therefore, the deviation between the two results remained within 6 μm. In summary, the 2D model can be used to effectively calculate material damage during impact processes. Figure 11 shows the impact process of the two particles. It was assumed that the material and size of the particles were the same in the two collisions. The red and blue balls represent the particles of the first and second impact events, respectively, with different positions of impact on the target material. D and u denote the particle diameter and horizontal distance between the two impact events, respectively. α and  v 0 denote the impact angle and impact velocity, respectively.

Impact process
The entire impact process can be categorized into the following four stages: www.Springer.com/journal/40544 | Friction 1) The first particle impacts the surface of the target material.
2) The compression and deformation of the material surface result in Crater 1. Meanwhile, the first particle departs from the material surface.
3) The second particle impacts the surface of the target material.
4) The compression and deformation of the material surface result in Crater 2, and the second particle departs from the material surface.
When the particles successively impacted the material surface, the compression and deformation of the material surface reached a maximum when the relative impact velocity between the particles and target material was 0. When the particles collided with the target, the target material deformed elastically. Moreover, an impulse was applied such that the particles departed from the target material. It is noteworthy that when the surface of the target material deformed plastically, a crater was formed.
A study [39] was conducted to determine a time interval that allowed the stress waves in the material to be negligible after each impact. In that study, a time interval of 200 μs was selected, which is consistent with the current assumption that the second impact is due when the first impact is fully completed.

Numerical simulation analysis 4.1 Effect of impact distance between two impact events
In the study pertaining to the horizontal distance between the two impact events, the particle diameter, particle impact velocity, and incident angle were 5 mm, 20 m/s, and 45°, respectively. Figure 12 shows the effect of different incident positions on the crater depth. The horizontal distances of the incident positions were 0D, 0.06D, 0.1D, 0.14D, 0.2D, 0.4D, 0.6D. It is noteworthy that the second impact location was on the left side of the first impact. In Fig. 12, O and O' represent the centers of the first and second particles,  respectively. The dotted yellow line demarcates the edges of the two craters. As shown in Fig. 12, the horizontal distance between the two impact events on the target material can be classified into four ranges. Figure 13 shows that for Range A (u ≤ 0.1D), Crater 2 created by the second impact almost completely covered Crater 1. It was observed that the impact point of Particle 2 was located in Crater 1 generated by the first impact, whereas the maximum von Mises stress was located in Crater 2. In this case, the second impact imposed a greater effect on the target material. For Range B (0.1D < u < 0.2D), the impact point was located outside Crater 1. The two craters affected large areas, and the area of Crater 2 was larger than that of Crater 1. As the horizontal distance between the two impact events increased, the length of Crater 2 caused by the second impact decreased gradually. In Range C (u = 0.2D), the two craters resulted in overlapping areas. However, the lengths of the two crater edges were the same. In Range D (0.2D < u < 0.6D), the overlapping area of the two craters disappeared gradually, whereas the affected area remained. As the horizontal distance between the two impact events increased, the von Mises stress distribution of the two craters remained constant. Moreover, the effect of the impact on the area between the two craters was negligible. Figure 14 shows the depth difference between two craters during the impacts of a particle couplet (H 1 − H 2 , where H 1 is the penetration depth of the first single particle, and H 2 is the penetration depth of the second particle that can impact the neighboring areas). As the horizontal distance between the two impact events changed from 0.06D to 0.6D, the depth difference between the two impact events varied significantly with respect to the impact distance. When the horizontal distance between the two impact events was 0.06D (the second impact of the particle was located in Crater 1), the maximum depth difference was 57.5 μm, i.e., H 2 was 57.5 μm less than H 1 , which indicates a significant decrease in material damage caused by the second impact on the eroded material. Additionally, as the horizontal distance between the two impact events increased from 0.1D to 0.6D (the second impact of the particle was located outside Crater 1), the depth difference between the two craters decreased from 6.15 to 0.11 μm, i.e., crater-crater interactions diminished. At a distance of 0.6D, the depth of the second crater became the same as that of the first crater, which implied a complete transition from overlapping to non-overlapping impact effects.
As indicated by the residual stress distribution of the target material, the length of the affected area after the second impact was 0.6 times the particle diameter D. In other words, when the horizontal distance between two collisions exceeded 0.6D, the interaction between the two collisions was negligible. Figure 15 illustrates the affected area in this case. In this section, energy analysis is presented to elucidate the impact process. The energy transformation is vital for describing the particle impact process. The energy balance [31,61] for the entire model can be mathematically expressed as where E U is the dissipated portion of the internal energy; E K is the kinetic energy; E FD is the energy dissipated by friction (dissipative energy refers to the energy lost during energy transfer); E W is the work performed by external forces on the body; E QB is the energy dissipated by the damping effect of solid medium infinite elements, i.e., E QB = 0; E I is the remaining energy, which is known as the internal energy; E V is the energy dissipated by viscous effects, i.e., E V = 0; E S is the applied elastic strain energy; E P is the energy dissipated by plasticity; E C is the energy dissipated by time-dependent deformation (creep, swelling, and viscoelasticity), i.e., E C = 0. If damage occurs in the material, not all the applied elastic strain is recoverable. At any given time, the stress σ c can be expressed in terms of the "undamaged" stress σ u and the continuum damage parameter d as Eq. (11): The damage parameter d is initiated from the value of 0 in the undamaged material, and increases to a maximum value of 1 in the fully damaged material. It is assumed that during unloading, the damage parameter remains fixed at the value achieved at time t. Therefore, Eq. (12) is obtained: where ε el denotes the elastic stress, V is the volume of the sample. Moreover, the energy dissipated through the damage is expressed as Eq. (13) where d t denotes the damage parameters at time t. Therefore, the energy conversion equation for this case is as Eqs. (14) and (15): where E total is the total energy. Equation (14) shows that the final energy loss is caused by the joint action of friction dissipation energy and damage dissipation energy. Figure 16 and Table 7 illustrate the energy change. During the entire impact process, the initial kinetic energy was 124.97 mJ. The target material was subjected to the first impact of the particles at 0-0.1 ms, which resulted in plastic deformation and damage. Moreover, the second impact of the particles occurred at 0.3-0.4 ms, which further intensified the deformation and damage of the material. After the particles impacted the target, more than 80% of the energy was converted into internal energy, including a slight damage dissipation energy, and the energy dissipated by the damage was 0.06 mJ (see Fig. 16), which accounted for only 0.07% of the energy conversion. The remaining energy dissipated in the form of friction. Therefore, the change in the target morphology is mainly caused by plastic deformation, and the effect of material damage on the target morphology is negligible. Figure 17 shows that as the horizontal distance between the two impact events increases, the internal energy of the target material initially decreased and then increased. Equation (16) and Fig. 18 show that during the impact process, owing to the energy dissipation and kinetic energy conversion, the energy conversion rate remained in the range of 48.91%-53.35%. Moreover, the kinetic energy of the particles was primarily stored in the target material in the form of internal energy.
where η, W u , and W t denote the energy conversion rate, increased internal energy of the target material, and initial kinetic energy of the particles, respectively.  Figure 18 shows that when the horizontal distance between the two impact events was approximately 0.1D, the energy loss was the highest. Subsequently, the energy loss decreased as the horizontal distance increased. This is mainly caused by the change in the contact area between the particles and target. According to the formula P = F/A 0 (P is the pressure on the surface of the material; F is the applied force, where the impact velocity of the particle is constant, i.e., F is a constant value; A 0 is the contact area between the particle and material), as the contact area decreases, the resistance of the material to the particle increases gradually. This results in an increased energy loss in the material. Figure 19 presents the contact area between the particles and target.
When the horizontal distance of the second impact was within 0.06D (Range A), the contact area was large, and the energy loss was small. As the horizontal distance increased, the second impact occurred near the lip of Crater 1 (between A and B). Meanwhile, more kinetic energy dissipated from the material surface containing the two particles. It is noteworthy that the friction dissipation energy and damage dissipation energy remained at 11.74-12.64 mJ and 0.030-0.052 mJ, which are negligible. As the horizontal distance between the two impact events increased, the final impact of the particles occurred in the range between C and D, and the change in energy    | https://mc03.manuscriptcentral.com/friction conversion reduced. When the horizontal distance exceeded 0.6D, the impact was similar to that of a single particle impacting the target.

Effect of impact velocity
To investigate the effect of impact velocity, the particle diameter, horizontal distance between the two impact events, and impact angle were set to 5 mm, 0.2D, and 45°, respectively. Moreover, the impact velocities of the two impact events were the same and were set to 14,16,18,20, and 22 m/s. Figure 20 shows that as the speed increased, the depth of the crater and the difference between the depths of the two craters increased. This is because the increase in the particle velocity increased the particle forces in the normal direction, which subsequently increased the crater depth. During the second impact, as the deformation of Crater 2 increased, the material at the bottom of Crater 1 was continuously squeezed, causing it to move upward and hence resulting in a decrease in the depth of Crater 1.
The higher stress region indicates a larger deformation of the target, which resulted in a more significant erosion. Figures 20 and 21 show that as the velocity increased, the maximum stress remained in Crater 2. Moreover, as the impact velocity increased from 14 to 22 m/s, the maximum depths of the two craters increased from 38.2 to 94.5 μm. Furthermore, Fig. 20 Depths of two craters vs. particle velocity based on the particle size of 5 mm, impact angle of 45°, and horizontal distance of 0.2D. the sizes of the two craters increased proportionally, and the contours of the craters remained unchanged. It is noteworthy that the change in the impact velocity did not significantly affect the deformation mechanism of erosion. Figure 22 shows that as the velocity increased, the damage dissipation energy and damage to the target material increased. However, apart from the impact energy introduced to the target material that was converted into internal energy, only an insignificant portion of the impact energy dissipated in the form of damage dissipation energy. This is because the particle was spherical, and the initial velocity was low. Moreover, the damage caused by the particle impact to the target material was limited. The material mainly underwent plastic deformation instead of material removal. Therefore, impact damage to the target material accumulated gradually and can create a macro erosion morphology, e.g., in an actual oil and gas pipeline. It is noteworthy that three factors affect material damage caused by erosion as follows: number of particles, impact type (e.g., continuous), and duration.

Effect of impact angle
To investigate the effect of impact angle, the particle diameter, horizontal distance between the two impact events, and impact velocity were set to 5 mm, 0.2D, and 20 m/s, respectively. It is noteworthy that the characteristics of the overlapped impact were conspicuous and can be distinguished from completely overlapped and non-overlapped impact events. Moreover, the impact angles of the two impact events were set to 25°, 35°, 45°, 55°, 65°, 75°, and 85°. Figures 23  and 24 show that as the impact angle increased, the depth of the first and subsequent impact increased. When the impact angle exceeded 75°, the change in the impact depth was slight. This is because as the impact angle increased, the speed along the vertical direction increased, and the impact intensity on the target material increased as well. Because the inertial force changes the direction of the particles instead of their magnitude, the inertial force along the vertical , is correlated to the sine value. As the impact angle increased from 25° to 85°, the depth difference between the two craters increased from 0.11 to 14.40 μm. Moreover, the depth difference increased significantly from 45° to 55°. When the angle exceeded 55° and was affected by work hardening, the depth difference increased gradually. Figures 25 and 26 show that when the impact angle α ≤ 45°, three stress concentration areas appeared at the bottom and both lips of Crater 2. Moreover, when 45° ≤ α ≤ 85°, the two stress concentration regions were located near the bottom and left lip of the crater. No stress concentration area was observed on the lip between the two craters. This is because the large    plastic deformation of the target resulted in more plastic deformation after the formation of Crater 1. Therefore, surface work hardening prevented the further deformation of the target material. Figure 27 shows that as the impact angle increased, the damage dissipation energy increased. When the angle increased from 75° to 85°, compared with the impact angle increment from 25° to 35°, the change in energy dissipation reduced from 0.01 to 0.003 mJ. When the impact angle exceeded 75°, the change in the angle did not significantly affect the crater depth. Meanwhile, the crater morphology on the surface of the target material stabilized.

Modified Bitter model
To reveal the effects of adjacent craters, the internal energy variation of the target material was redefined in the modified model.

Modified model
The Bitter model is expressed mathematically as shown in Eqs. (17) and (18). where Q p is the formation of the permanent indentation, which requires a certain amount of energy; m p is the quality of the abrasive particle; α is the impact angle of the particle; K B is a constant related to the velocity of the particle before collision; W D is the wear due to repeated deformation; and ε B is the deformation wear factor. The effects of the impact velocity, impact angle, and horizontal distance between two craters on the internal energy were revised using the First Optimization (1stOpt) fitting software. Hence, the modified energy (Q pz ) was replaced by the Bitter model.

Model verification
As shown in Fig. 28, the calculation results of the fitting formula are similar to the calculation data of the model. In addition, four evaluation indicators, i.e., RMSE (root mean square error), SSE (sum of squares due to error), R (correlation coefficient), and R 2 (square of correlation coefficient), were used to assess the quality of the fitting results. The values of RMSE and SSE were approximately 0, whereas the values of R and R 2 were approximately 1. Hence, this formula can characterize the simulation results accurately.

Conclusions
In the present study, considering the JC hardening criterion and JC damage criteria, a numerical model for simulating overlapping collisions of spherical | https://mc03.manuscriptcentral.com/friction particles on AA6061 aluminum alloy was established. Moreover, the accuracy of the model was verified by comparing the results obtained with the experimental results. The following conclusions were obtained from this study: 1) The impact location of the neighboring particles, which is defined as the horizontal distance between the solid particle impacting location, significantly affects the erosion magnitude and pattern. As the horizontal distance between the two impact events increased from 0 to 0.6 times the particle diameter, the interaction between the formed craters reduced gradually. When the horizontal distance exceeded 0.6 times the particle diameter, the interaction between two craters disappeared.
2) As the horizontal distance between the two impact events changed from 0.06D to 0.6D, the depth difference between the two impact events changed significantly with respect to the impact distance. When the horizontal distance between the two impacts was 0.06D, the maximum depth difference was 57.5 μm. The impact on existing craters can significantly decelerate erosion. As the horizontal distance between the two impact events increased from 0.1D to 0.6D, the depth difference between two craters decreased from 6.15 to 0.11 μm, which implies diminished cratercrater interactions.
3) As the horizontal distance between impact events increased, the energy conversion rate during the impact initially decreased rapidly. Subsequently, it increased gradually after it reached a critical point, which was 0.1 times the particle diameter. 4) As the impact velocity increased from 14 to 22 m/s, the depth of the crater increased from 38.2 to 94.5 μm, whereas the profile of the two craters remained unchanged.
5) The removal of the target material by the impact particles was limited by plastic deformation erosion. Moreover, the change in the target morphology was mainly caused by plastic deformation. The continuous impact and long-term accumulation of numerous particles contributed significantly to pipeline erosion.
6) The second impact was hindered by work hardening as the depth increased. This resulted in an area without stress concentration on the lip between the craters, which decelerated the erosion of the material.