Atomistic assessment of structural evolution for magnesium during hypervelocity nanoprojectile penetration

In the present study, investigating the effect of ballistic penetration of spherical projectiles on a monocrystalline magnesium specimen is performed using embedded atom method (EAM) potential in molecular dynamics (MD) simulation. The dynamic investigations of structural evolution based on common neighbor analyses and Wigner–Seitz defect analysis are carried out for the varying depths of penetration and velocities of the projectile (v = 2 km s−1, 6 km s−1, and 10 km s−1). It is found that the extent of amorphization in the specimen is more in the case of higher depth and lower projectile velocity. Voronoi cluster analyses are also done to identify cluster distribution and their transformation during ballistic penetration, which is accompanied by atomic strain and displacement vector evaluation to give light to the effect of shear strain and displacement of atoms respectively. According to Voronoi cluster analysis, Voronoi polyhedra having < 0, 4, 4, 6 > and < 0, 6, 0, 8 > clusters exhibit a higher population during hypervelocity projectile penetration. The findings have potential applications in hypervelocity applications such as defense and space technologies.


Introduction
The ballistic or hypervelocity impact is defined by relative velocities of the projectile and target material which determines the loading level of varied substrates and their physical response to the impact. The phenomenon is called ballistic when the relative impact velocity between the projectile and target material is 2 km s −1 or above [1]; the mechanical properties and failure mechanisms of the substrate depend on the loading rate, their response, the generation of kinetic energy, and entropy change [2].
The hypervelocity impact of nanosized specimens on metals has gained much popularity in recent years due to its enormous applications in the design of armor [3,4], modification of any surface [5,6], impact with micrometeorites [7,8], ballistic penetration [9,10], and the development of cluster technology [11]. Experimental and theoretical studies revealed the dynamic mechanism of a microjet from a grooved metal produced from the hypervelocity impact [12]. The evolution of plasma and radiofrequency due to hypervelocity impacts of meteoroids and dust was studied along with the exploration of parameters in the creation of dust and its dynamics [13]. Experimental studies were conducted to study the hypervelocity impact on biopolymerbound soil composites (BSC) at the National Aeronautics and Space Administration's (NASA) White Sands Testing Facility to verify numerical simulations of shock physics code developed by Sandia National Laboratories [14].
Though the uses of ballistic penetration are enormous, experimental studies are yet insufficient to analyze the collision of clusters as they offer limited information and narrow scope to understand the hypervelocity collision dynamics. To address the issue, atomistic simulations are employed as their length and timescales efficiently capture the movement of ultrafast particles [15]. Several atomistic simulations are done on hypervelocity impact to study the ballistic behavior of multilayer graphene polymer composites (MGPCs) using molecular dynamics (MD) simulations [16]. Fragmentation characteristics were studied at macro-and micro-scales at high and low shocks using MD and smoothed particle hydrodynamics (SPH) [17]. The investigation of the effect of interlayer spacings and the degree of graphene oxide nanochannel oxidation was performed using MD simulations [18]. MD simulation studies were done on single-crystal aluminum to study the plasma phase transition due to the hypervelocity impact [19]. Atomistic simulations were also done on graphene to figure out the mechanical and chemical properties using the reactive force field by the bombardment of nanoprojectiles [20]. Therefore, MD simulations are reliable and widely popular for the study of the transient ballistic phenomenon at varied hypervelocities.
There have been several nanoindentation processes devoid of the MD process which involves the use of an indenter to penetrate to a certain depth of the specimen, generally less than the specimen radius [21][22][23][24]. The strengthening mechanism has been studied for coated graphene in metal nanocomposites by nanoindentation by MD [25]. MD simulation was performed on Ni using the many-body tight-binding potential method to perform nanoindentation for evaluating Young's modulus and hardness values [26]. Indentation of the SiO 2 /Si bilayer composite was studied to determine the stress relaxation behavior using MD [27]; mechanical properties of single-crystal diamonds were also determined using MD during nanoindentation [28]. Techniques similar to nanoindentation, including cold-spraying, have been performed using MD simulation on the single-crystal copper to study the effect of crystal orientation [29].
In the present study, the ballistic penetration concept is an extrapolation of the indentation process where the hypervelocity penetrator substitutes the indenter affecting the surface and also penetrates the bulk of the specimen. This modified ballistic penetration process gives substantial information regarding the material characteristics during the penetration effect which lacks any literature, especially for hexagonal close-packed (HCP) materials such as magnesium (Mg). The primary difference between the indentation and penetration process is the domain of study of the specimen. The former is more concerned about the surface characteristics whereas the ballistic process indicates the bulk resistance to deformation. Furthermore, the deformation strain rate is much lower for indentation as compared to a hypervelocity ballistic process. Thus, the hypervelocity penetration process gives information about the change in material structural distribution due to the rapid deformation process having a higher penetrator speed as compared to slower deformation processes like nanoindentation.
Currently, the popularity of Mg-based materials is increasing due to their future potential for lightweight and cost-saving applications [29][30][31][32]. The demand for non-exhausting, energy-efficient, and advanced resources is rising which can be potentially fulfilled by Mg [33]. Mg and its alloys have a potential wide scope of application in biomedical, aerospace, energy sectors, etc. in the modern era [33][34][35][36]. In the present study, the hypervelocity penetration of monocrystalline Mg is done at different projectile speeds. The purpose of the simulations is to understand the behavior of magnesium specimens to hypervelocity impact by analyzing the rate of amorphization, change in polyhedral structures, and its atomic displacement and strain with the change in projectile speed. Four different techniques are used for the analysis of the specimen during the ballistic process: (a) common neighbor analysis [37] is performed to understand the amorphization, and distribution, (b) atomic displacement analysis [38] to analyze the distribution pattern of displaced atoms, (c) shear strain analysis [27] to evaluate the shear distribution in the specimen during deformation, and (d) Voronoi analysis [39] to find the local atomic cluster population distribution.

Simulation details
The hypervelocity ballistic penetration of pure monocrystalline Mg is done for velocities (v) ranging from 2 to 10 km s −1 . Applications for high-velocity impacts fall into two categories: low-velocity range applications (< 3 km s −1 ) and high-velocity range applications (> 8 km s −1 ). Military armor, impact erosion by combustor products in gas turbines for low-velocity impacts, and spacecraft colliding with meteoroids for highvelocity impacts are some of these uses [16]. Hence, different speeds across the spectrum are chosen: (a) v = 2 km s −1 for a low-velocity range, (b) v = 10 km s −1 for a high-velocity range, and (c) an intermediate speed of v = 6 km s −1 .
An embedded atom method (EAM) potential is used for simulation [40] whose potential energy (U) is given in Eq. 1.
where r i,j is the distance of separation between atoms i and j, N is the total number of atoms in the system, φ is the pairwise potential, and Φ is the embedding energy function. The spherically averaged atomic electron density, i , is given by Eq. 2.
where the is the density function. For the hexagonal closed-pack (HCP) system, the solid-liquid interface ( ) free energy (SLI) is calculated by expanding it to second-order harmonics as shown in Eq. 3 [41].
where 0 , 20 , 40 , 60 , and 66 are material-dependent expansion coefficients that describe how the depends on the inclination. Furthermore, the discrepancy (r) is quantified between an empirical relation and MD simulation data as shown by Eq. 4 [40].
where the γ rel i and γ MD i are the free energies from empirical relation and MD simulation respectively and n is the total number of the material tested. Thus, from the deviation values, the best model was compared to that predicted by the Turnbull, Laird, and Ewing relations, and it was found that the Ewing relation best describes the simulation results [40].
The specimen of pure Mg was prepared in Atomsk software having dimensions 187.91 Å × 187.56 Å × 186.66 Å containing 288,864 atoms [42]. The specimen was initially relaxed using an isothermal-isobaric (NPT) ensemble and the simulation is done at a low temperature (10 K). Shrink-wrapped with movable boundary is taken in the direction of penetration (Z-axis) which is acting on the basal plane and the remaining two axes are taken periodically. The projectile used for ballistic penetration is spherical with a radius of 25 Å and the force exerted by the projectile is given by Eq. 5 [43].
where k is the force constant, r is the distance from the center of the indenter to the atom, R is the radius of the indenter, and F(r) is the repulsive force exerted by the projectile on the atoms to be penetrated. The force here, unlike in the nanoindentation process, will be used both for the surface and bulk of the specimen until the projectile exits the specimen from the rear face. The timestep of 0.001 ps is taken and minimization is done in the specimen by the conjugate gradient (CG) method. To carry out the simulation, NVT (canonical) ensemble is used to maintain the temperature by the Nose-Hoover thermostat algorithm [44,45]. To stabilize the thermal fluctuation, rescaling is done at every 5 ps. To keep the Mg substrate fixed, a vertical region of width 15 Å is fixed perpendicular to the y-axis as shown in Fig. 1. The entire simulation is performed using Large-scale Atomic/Molecular Massively Parallel Simulator (LAMMPS) software [46]. Furthermore, to visually analyze the specimen Open Visualization Tool (OVITO) is used where displacement analyses [38], shear strain analyses [27], common neighbor analyses [37], and Voronoi analyses [39] are done.

Results and discussion
The depth and atomic displacement are presented in Fig. 2 to evaluate the underlying physics of the possible structural evolution of the Mg specimen during ballistic penetration, The simulations were performed for 2 km s −1 , 6 km s −1 , and 10 km s −1 where the hypervelocity penetration has caused a trace of an amorphized layer of atoms along the path of the projectile. The thickness of the amorphized layer at the sides which is at the trailing edge is represented by d1, and the thickness at the leading edge of the projectile is given by d2. The d1 and d2 signify the amorphization thickness which indicates the extent of turbulent structural change due to the penetration process. Thus, these thickness values indicate the extent of material damage in the specimen and correspondingly give a quantitative parameter during comparison with different projectile speeds and depths of penetration.
As the depth of penetration increases, the number of vacancies and interstices increases. This is because the projectile displaces the atoms in its path outwards as shown by red-colored vectors in Fig. 2b creating interstices while and Wigner-Seitz defects analysis [47]. The path taken by the projectile tends to taper off due to the displaced atoms moving back in direction of the parent site. The magnitude of displacement vectors at the leading end of the projectile is more than at the trailing end as shown in Fig. 2b which indicates that the net displacement of atoms at the leading side is more than that of the trailing side. Thus, due to penetration, there is an increase in the number of vacancies and interstices while at the same time the recovery takes place with the decrease in the point defects as displaced atoms tend to move back to their parent sites. The effect of the creation of point defects when the projectile is in the specimen is much faster than the recovery effect. There is a positive net generation of vacancies and interstices as the projectile goes deeper into the specimen. The recovery effect cannot be ignored due to which the d1 value increases over time. The recovery effect is insignificant for d2 due to the sheer lack of time available on the leading end of the projectile as compared to d1 indicating d1 > d2. The recovery process is solely present as the projectile leaves the Mg specimen from the rear end with a continuous decrease in vacancies and interstices with time.
The diameter of the path traversed reduces due to the migration of interstices to fill up the vacancies previously created.
However, the d1 and d2 values are seen to decrease with the increase in projectile velocity (v) as shown in Table 1.
The thickness of amorphization decreases as the time required for transformation decreases with increasing v, and hence, at any depth, we can observe less percentage of atoms transform from their crystalline to amorphous state. It is also seen that the d1 and d2 values would consistently increase with the increase in depth of the projectile at any given velocity as shown in Table 1.
As stated earlier, to analyze the displacement of atoms due to penetration, displacement analysis is done with a depth of penetration of 50 Å at different speeds as shown in Fig. 3. It can be seen that as v increases from 2 to 6 km s −1 , the magnitude of atomic displacement sharply decreases and this further decreases at v = 10 km s −1 . This is because with the increase in v there is lesser time available for atoms near the projectile to impart displacement to their neighboring atoms. Furthermore, the top view at v = 2 km s −1 shows a symmetrical pattern that tends to follow a six-fold symmetry on the XY plane along the directions [0 1 0], [1 1 0], and [1 10] as shown in Fig. 3d.
With the increase in v, there is a decrease in the shear strain effect as can be seen from the top view of the Mg specimen in Fig. 4. It can be seen that the influence of shear strain in the case of v = 2 km s −1 is more than that of    Fig. 4a. The average atomic displacement (S) and shear strain (γ) values are illustrated as shown in Fig. 5. It can be seen that with the increase in penetration depth, S and γ values increase for all penetration speeds. However, the rate of increase varies for the three given velocities where it can be observed that S and γ values for v = 2 km s −1 are much higher as compared to the other two velocities. This is evident due to the difference in time of contact with the change in v. The lesser the v, the more will be the time available for atoms to displace, and hence, the S value will be significantly higher at lower speeds. Moreover, the time of contact will determine the shear force generated in order to displace the neighboring particles. Thus, the available time for strain generation will be much lesser for higher speeds. There is a lack of significant difference in values at v = 6 km s −1 and v = 10 km s −1 . For higher velocities, the values will decrease slightly with an increase in the

Cluster analysis
Voronoi analysis helps us to understand the local atomic arrangements at every step of the simulation [39,48]. Voronoi polyhedrons (VPs) are the smallest polyhedrons in which their planar faces will bisect at right angles to the lines connecting particles such as atoms, ions, or molecules [49]. A VP is identified by the number of faces (n i ≥ 0) where i is the number of edges of the VP around the central atom. The coordination number is given by Σn i and the VP can be expressed by the indices < n 3 , n 4 , n 5 , n 6 > . This can be explained by an example < 0, 3, 6, 5 > which indicates a coordination number of 14 (Z14) having 3 polyhedral faces with 4 polyhedral edges, 6 polyhedral faces with 5 edges, and 5 polyhedral faces with 6 edges. To identify the varied clusters and their transformation with different speeds of the projectile, Voronoi analysis is done in OVITO [38].
The Voronoi analysis is studied for projectile velocities of 2 km s −1 , 6 km s −1 , and 10 km s −1 at depths of 50 Å, 85 Å, and 120 Å to identify the compositions of VPs and their relative change with respect to reference configuration as shown in Fig. 7. Due to ballistic penetration, there is seen to be a significant increase in < 0, 4, 4, 6 > and < 0, 6, 0, 8 > clusters while there is seen to be a reduction in the rest of the VPs. The VPs that have majorly reduced in numbers are < 0, 3, 6, 5 > , < 0, 3, 6, 6 > , < 0, 2, 8, 4 > , and < 0, 4, 4, 8 > . The amount of change in VP composition depends on the speed of the projectile where there is seen to be a significant change for lower speeds like 2 km s −1 while the least change for 10 km s −1 . The speed of the projectile dictates the time available for necessary atomic displacement and atomic strain as seen in Figs. 3 and 4 respectively and thus, the required displacement is less for the complete transformation of the VPs.
The specimen in Fig. 8 is divided into five planes A, B, C, D, E, and F which are located at depths of 20 Å, 50 Å, 80 Å, 110 Å, 140 Å, and 170 Å respectively whose top views are shown in Fig. 9 to visualize the VP distribution. It is evident from Fig. 8 that an increasing depth of the planes reduces < 0, 6, 0, 8 > clusters while < 0, 4, 4, 6 > clusters increase as shown in Fig. 8d and then start to decrease in Fig. 9e and f. < 0, 6, 0, 8 > and < 0, 4, 4, 6 > follow similar to a six-fold symmetry in the XY plane along the directions [0 1 0], [1 1 0], and [1 10] up to 90 Å depth where the former (shown in green) plays the major role at 20 Å depth. < 0, 4, 4, 6 > has the highest population (shown in yellow in Fig. 8c and d) following a distribution similar to six-fold symmetry at a depth of 80 Å and 110 Å. The degree of symmetry is lost at greater depths. < 0, 3, 6, 4 > and < 0, 4, 4, 7 > as shown in black and red, respectively, tend to follow six-fold symmetrical lines on the XY plane along the directions [0 1 0], [1 1 0], and [1 10] as shown in Fig. 9c and d. Hence, it can be concluded that the symmetrical patterns exist along the path traversed by the projectile and also some depth ahead of the projectile which is due to the effect of the shock wave influencing the atomic displacement. This symmetry would decrease as we move to greater depths where the influence of the projectile is insignificant. There is a lack of symmetry and fewer penetration effects at depths above 140 Å where < 0, 3, 6, 5 > (shown in light blue) is found to be significantly present. < 0, 3, 6, 5 > would greatly reduce in population with an increasing penetration depth of projectile due to its transformation to < 0, 4, 4, 6 > .

Conclusion
In the present MD study, hypervelocity penetration simulations are used to investigate the atomistic structural parameters of pure monocrystalline Mg. The hypervelocity penetration is done at different velocities to assess the possible structure-property correlation-based performance of Mg. The effect of velocity plays a crucial role in the time available for amorphization, the effect of atomic displacement, shear strain during penetration, and finally on the percentage transformation of VPs from their initial configuration. The findings are summarized as shown below.
1. The amorphization thickness decreases with an increase in the projectile velocity (v), while it increases with greater depth of projectile penetration. d1 is greater than d2 and is independent of v and the depth of penetration. 2. Point defects increase with the increase in the penetration depth but decrease after the projectile leaves the specimen from the other end of the penetration side. This creation of vacancies and interstices is accompanied by a reduction in point defects due to the recovery effect. The annihilation of point defects during recovery is much slower than the creation of the same during penetration. 3. The higher the v value, the less the magnitude and area of atomic displacement. The atomic displacement area tends to take a six-fold symmetrical shape at v = 2 km s −1 , while the area is seen to be circular at higher v values. 4. An increase in v causes a decrease in magnitude and area of shear strain. Similar to atomic displacement, the shear strain area takes a six-fold symmetrical shape at v = 2 km s −1 , while the area is circular at v = 6 km s −1 and 10 km s −1 . Fig. 10 Pathways of VPs during ballistic penetration at 10 K for pure Mg 5. Mixed and ICO-like VPs are reduced and crystal-like polyhedrons are increased during the penetration process. An increase in < 0, 4, 4, 6 > , and < 0, 6, 0, 8 > creates a symmetrical distribution. < 0, 4, 4, 6 > is seen to be the most active polyhedron with the highest number of transformations and population in the specimen.
The present study provides an insight into the atomistic structural behavior of Mg under hypervelocity penetration which can provide the base work for future ballistic applications in space technology and defense. This study can also be extrapolated to examine the property of Mg alloys and composites to prepare more lightweight and energy-efficient materials for high-velocity applications in the near future.