Dynamic evolution of microstructure during laser shock loading and spall failure of single crystal Al at the atomic scales

The capability to predict the response of metals subjected to extreme loading environments of high strain rates relies on a fundamental understanding of the microstructural evolution in the metal at the atomic scales. This comprises the investigation of the links between loading conditions and the dynamic evolution of defects (dislocations) of the metal. Such an understanding is particularly challenging for conditions of shock loading that result in significantly high pressures generated in the metal in very short times. The capability to vary the pressures generated in the metal and time scales for the microstructure to respond using shock loading can result in microstructures with variations in defect distributions and hence the dynamic response of the metal. The understanding of the spatial and temporal evolution of defects and their interactions in such environments is therefore a critical challenge, as their interactions can either render strengthening mechanisms or render defect nucleation sites to nucleate damage and initiate failure.

The failure strength under shock loading (spall strength) is typically used as a metric to characterize the dynamic response of the metal and is defined as “the peak tensile pressure the material can withstand prior to nucleation of voids and initiate failure.” The spall strength of metals is typically investigated by impacting a flyer plate on to a metal target and renders strain rates ranging from 10⁵ to 10⁷ s⁻¹.^1–4 More recently, shock loading using short pulse lasers has allowed pushing the experimental limits to investigate the dynamic response of materials and spall failure at strain rates greater than 10⁷ s⁻¹.^5–12 The current advancements of these experiments allow one to investigate the shock compression and spall strengths of single crystal and polycrystalline samples with dimensions as low as several hundred nanometers that render strain rates of up to 10¹⁰ s⁻¹.^13–15 The experimental investigations of the shock response are largely focused on the characterization of the shock pressures, shock velocities, and the spall strength of metals using rear surface velocity profiles. While these studies suggest a clear strain rate dependence and a shock pressure dependence of the spall strength of the metal, it should be noted that the loading conditions under plate impact and laser shock result in large variations in the microstructural response. For example, laser interaction with metals can result in ablation, spallation (front end), and melting that can result in shock wave structures,^13–16 propagation, and hence the microstructural response which are significantly different from that observed during plate impact experiments. The reader is directed to Ref. 16 for a fundamental understanding of the time and length scales of the processes during laser shock loading using experimental and computational methodologies.

As a result, the evolution of microstructure (defect density) becomes a critical challenge that needs to be investigated to determine the dependence on spall response. The short length and the time scales of the processes, however, limit the experimental capabilities to unravel the links between defect evolution and spall strengths. Classical molecular dynamics (MD) simulations can provide critical insights toward the deformation and failure mechanisms.^17–26 These MD studies typically use a piston that is driven into the metal at a constant velocity and renders a constant pressure shock wave in the metal and deformation at strain rates of 10⁹ s⁻¹ and higher. While the MD simulations reported in these studies can characterize defect nucleation and evolution, the strain rates of the MD simulations and sample dimensions correspond to experimental capabilities achieved using laser shock loading. As a result, a direct comparison of MD simulations for the length and time scales of the phenomena of spallation requires the capability to model the phenomena of laser material interactions (ablation, spallation, and melting) to have an accurate prediction of the shock wave generation and spall failure behavior.

The capability to model the microstructural evolution during laser shock loading and spall failure requires an accurate incorporation of the mechanisms of laser energy deposition, heat generation, and dissipation in the metal microstructure to have an accurate description of the shock response. The conduction band electrons absorb the laser energy deposited in the metals, and the electron-phonon coupling transfers this energy to lattice vibrations and results in heating of the material that flows from the surface (interacting with the laser) into the bulk of the metal. A critical challenge in the capability to model shock propagation and spallation behavior during laser shock loading, therefore, is the ability to retain the underlying mechanisms that govern the laser energy deposition. The incorporation of accurate physics for these mechanisms requires a bottom-up approach that allows the atoms to absorb the deposited energy, heat generation and transfer, melt nucleation and propagation, and related defect nucleation and evolution.

MD simulations combined with a two-temperature model (TTM) create a hybrid atomistic-continuum approach that is able to describe laser energy deposition, heat conduction, melting, ablation, phase transformation, and microstructural evolution at the atomic scales.^27–32 This paper demonstrates the capability of the MD-TTM applicability of the hybrid MD-TTM method to investigate the spall failure of single crystal Al samples during laser shock loading and unravel the links between laser fluence as well as loading orientation and melting, defect evolution, and spall failure behavior at the atomic scales. The MD-TTM simulations are carried out for a range of laser fluence and loading orientations of single crystal Al to investigate the variations in the shock wave structures, defect nucleation and evolution, as well as spall failure behavior. As will be discussed later, the laser induced shock wave structures and microstructural evolution is observed to be different from that reported for the case of plate impact loading (piston shock). The “computational details” are discussed in Sec. II and the atomic scale mechanisms of “laser shock induced spall failure” are discussed in Sec. III. The “effect of laser fluence” on the melting behavior, defect evolution, and spall strengths is discussed in Sec. IV and the “effect of loading orientation” is discussed in Sec. V. A comparison of the predicted shock response with that observed in experiments is summarized in Sec. VI.

The capability of the MD-TTM method to model microstructural changes during short pulse irradiation has been successfully demonstrated to investigate shock pressures,^30,32–34 ablation/melting of metal surfaces,^35–39 and nanostructuring of thin films during laser processing.^40–43 The MD-TTM simulations investigate the laser-induced shock compression and spallation failure behavior of metals at the atomic scales and the role of laser parameters (fluence) and loading orientation on the defect generation and mechanisms of spall failure. The typical laser shock setup that results in spall failure is shown in Fig. 1. Here, laser (for a given fluence and duration) interacts with the metal at the front end (left) of the sample and results in the mechanical removal of the material at the front end through either ablation or spallation as well as melting of the material. While the mechanisms of ablation and spallation of the metal at the front end have been studied before, this paper investigates the defect evolution during laser shock compression, the wave reflections, and interactions to initiate spall failure. The thermal response generates a shock wave that propagates toward the rear end (right) of the metal. The deformation of the metal in the current simulations corresponds to a uniaxial strain, i.e., strain only in the shock direction (Z) and zero strain in the lateral directions (X and Y) for both tension and compression, i.e., $ε X =0, ε Y =0, and ε Z ≠0$⁠. The state of such a tensile strain corresponds to a state of triaxial tensile stress, i.e., $σ X ≠0, σ Y ≠0, and σ Z ≠0$⁠.

The MD-TTM simulations are carried out for a 25 nm × 25 nm × 500 nm [001] single crystal Al system (∼19 × 10⁶ atoms), wherein the embedded atom method (EAM) interatomic potential⁴⁴ is used to define the Al interactions and the parameters for TTM were taken from Ref. 32. The single crystal system is subjected to laser irradiation with fluences of 0.5 kJ/m², 1 kJ/m², 6 kJ/m², and 13 kJ/m² and a duration (pulse) of τ = 150 fs (full-width half maximum). The MD-TTM simulations are also carried out for the [110] and [111] orientations to compare the shock pressures, wave structures, as well as spall strength for the same laser fluence (13 kJ/m²) and duration (150 fs). The temporal response of the pressures and temperatures in the metal is analyzed by dividing the system into several sections along the laser loading direction and averaging the values for pressure and temperature in each section. The microstructural characterization uses “centro-symmetry parameter” (CSP)⁴⁵ and “common neighbor analysis” (CNA)⁴⁶ algorithms to analyze defects and damage in the microstructures generated. In addition, densities of Shockley partials (1/6<112>), Perfect dislocations (1/2<110>), Frank partials (1/3<111>), Stair-rod partials (1/6<110>), and Hirth locks (1/3<001>) are quantified at various times during the simulation using the “dislocation extraction algorithm” (DXA)^47,48 implemented in the “crystal analysis tool” (CAT).⁴⁹ 

The interaction of the metal with lasers for a fluence of 13 kJ/m² and a duration of τ = 150 fs results in ablation and melting at the front end and generates a compressive pressure wave that travels toward the rear end of the metal. The predicted shock wave structures (variation of pressure along the loading direction) in the metal are shown in Fig. 2 for various times. As compared to the response under piston shock simulations,^17,18,23,25 the shock wave structure due to laser interactions does not reveal any constant pressure zones in the metal behind the shock front. The peak value of compressive pressure generated in the liquid phase is calculated to be ∼50 GPa and is observed to decay as the compressive wave travels through the metal toward the rear surface. The structure of the shock wave is observed to have a two-wave structure: a plastic wave traveling with a speed of 6260 ± 250 m/s and an elastic wave with a speed of 7530 ± 200 m/s. The shock wave structure (as observed for times of 40 ps to 50 ps) shows a change in the shape of the compressive pressure wave. A wide bump as compared to a sharp rise in pressure is observed at these times and is attributed to an artificial phase transformation of the FCC microstructure to the BCC phase predicted by the potential. The equation of state for Al as predicted by the EAM potential for Al⁴⁴ indicates a crossover of the curves for the FCC and BCC phases at a pressure about 30 GPa. As will be discussed later, this phase transformation does not affect the spall failure of the metal as the tail of the compressive wave transforms the BCC phase back to the FCC phase.

The evolution of pressures generated due to laser interaction in the metal as a function of time is shown by plotting the average pressures along the length of the sample in Fig. 3(a). The compressive shock wave [red regions in Fig. 3(a)] that travels toward the rear surface reflects and interacts with the plastic wave (red arrows) as well as the tail of the pressure wave (black arrows). These wave interactions generate tri-axial tensile stresses in various regions in the system. Voids are observed to nucleate in these regions under triaxial tensile stresses and grow and coalesce to initiate spall failure. The “maximum tensile pressure” generated in the metal due to the interactions of the waves that result in void nucleation is computed as the “spall strength” of the metal. The maximum tensile pressure is calculated to be 6.7 GPa for the spall region closer to the rear surface (attributed to the interaction shown by red arrows), whereas the interaction shown by black arrows generates tensile pressures of 6.2 GPa for the spall region in the interior.

The evolution of temperatures with time generated due to laser interaction in the metal is shown by plotting the temperatures along the sample length in Fig. 3(b). The temperatures generated due to laser energy deposition result in rapid melting at the front end experiencing the laser energy deposition that propagates toward the rear surface. The kinetics of the melting can be investigated by the propagation of the solid-liquid interface as indicated by the dashed-dotted line in the plots in Fig. 3. The melt region is characterized using the average values for the centrosymmetry parameter and the potential energy for each bin along the length of the sample. The atoms are characterized as liquid for CSP values that are greater than 0.05 and a potential energy value that is greater than −3.0 eV/atom. These values are determined based on the pressure dependence of the phase stability of the liquid phase. The velocity of the rear surface showing the spall signal is shown in Fig. 3(c) and can be used to compute the strain rates generated as well as the spall strength. The rear surface velocities render a value of 7.1 GPa for the spall strength and compares well with the value of 6.7 GPa computed from the pressure profiles above. The spall strength values computed here compare well with the experimentally reported values of 5.8 GPa for the same parameters of laser shock loading and sample dimensions (depth).¹³ 

Four stages of shock loading and spall failure are identified using the pressure profiles and the rear surface velocity profile. Stage “T1” corresponds to laser energy deposition resulting in melting and generation of the shock compression wave of the metal that travels toward the rear surface. Stage “T2” starts when the elastic wave reaches the rear surface and results in the acceleration of the rear surface till it reaches a peak value. Stage “T3” corresponds to the deceleration of the rear surface and interactions of the waves (as discussed before) that result in maximum tensile pressures and nucleation of voids. Stage “T4” corresponds to the growth/coalescence of voids in the metal as the triaxial tensile wave travels toward the front surface. These four stages can be used to compare the evolution of microstructure and defect densities in the metal for the various simulations considered here.

Illustrative snapshots that show the microstructure predicted by the MD-TTM simulations of laser shock loading are shown in Fig. 4. Only part of the system (Z coordinate > 1500 Å) is shown here, and the front end of the metal that undergoes laser interaction induced melting is neglected. The liquid atoms are colored as red, FCC stacked atoms as green, disordered atoms as blue, and surface atoms as orange. The dark blue bands at a depth of ∼3000 Å to 3500 Å at a time of 50 ps as shown in Fig. 4(a) represent the phase transformed regions (an artifact of the potential) in the metal and are re-transformed by the tail (release) of the pressure wave to the FCC phase at a time of ∼80 ps. The first voids nucleate at a time of ∼90 ps [attributed to interactions near the rear surface as shown by red arrows in Fig. 3(a)] and more voids continue to nucleate in the metal as the tensile pressure wave moves across the sample toward the front end of the metal. The solid-liquid interface is observed to travel slowly toward the rear surface throughout the simulation.

The densities of the dislocations in the metal are plotted in Fig. 5 at various times. The deformation of the metal during stages T1 and T2 is observed to be dominated by Shockley partials. A peak value of dislocation density is observed till the plastic wave arrives at the rear surface [at about 70 ps as shown by the red arrow in Fig. 3(c)] during T3; after which, the acceleration of the rear surface results in a considerable annihilation of the density of Shockleys during T2. The deceleration of the rear surface in T3 results in a further reduction in the density of Shockleys till it reaches a minimum value. Voids nucleate during T3 and growth of voids results in an increase in the density of Shockley partials as the tensile pressure reaches a peak value at the end of T3. Shockley partials continue to increase during T4 till the tensile pressures exist in the crystalline phase of the metal. As a result, a decrease in the density of Shockley partials is observed as the tensile wave travels toward the front surface where crystalline regions do not exist. The contributions of other dislocations to plasticity are observed to be considerably low as compared to Shockley partials. Illustrative snapshots (at times corresponding to in Fig. 4) showing the distribution of Shockley partials, twin faults, stacking faults, and voids are shown in Fig. 6.

To investigate the role of laser fluence, MD-TTM simulations are carried out using laser fluence values of 0.5 kJ/m², 1 kJ/m², and 6 kJ/m² in addition to 13 kJ/m² and a duration of τ = 150 fs. The pressures predicted in the metal along the sample length are compared in Fig. 7, where the dashed-dotted lines indicate the time evolution of the solid-liquid interface. The values for shock pressures, shock velocities, strain rates, and spall strengths are tabulated in Table I for the various loading laser fluences. The location of the spall plane from the rear surface (L_spall), the location of the solid-liquid interface (L_melt), and the temperature in the spall region at the maximum tensile pressure are tabulated in Table I. It can be seen that high laser fluence results in increased shock pressures and a larger depth of penetration of the liquid phase; a clear dependence of the spall strength on the fluence is not observed for the systems studied here. While no spall failure is observed for a fluence of 0.5 kJ/m² for the pulse duration chosen, the spall strengths for the other fluences are computed to be 6.8 GPa and 6.35 GPa for loading conditions of laser fluences of 1.0 kJ/m² and 6 kJ/m², respectively. The strain rates are computed to be 3.2 × 10⁹ s⁻¹ and 3.9 × 10⁹ s⁻¹ for 1 kJ/m² and 6 kJ/m² laser fluences, respectively. The tabulated values also suggest that a direct correlation does not exist between the spall strength and the strain rates or with the temperature generated in the spall region (as the lowest spall strength does not correspond to the highest temperatures generated at the spall plane). As a result, this study investigates additional factors (initial microstructure as well as its evolution) that can explain the dependence of the spall strength observed here.

The evolution of dislocation densities is plotted in Fig. 8 for the various laser fluences to investigate a correlation between the observed variations in the values of spall strengths. No plasticity is observed at a fluence of 0.5 kJ/m², and hence, the evolution of dislocation densities is not discussed for this fluence. The evolution of dislocation densities for the 6 kJ/m² and 13 kJ/m² fluences is observed to be very similar with a larger density of dislocations observed for the higher fluence. For the case of 1 kJ/m², no plasticity is observed during stages T1 and T2 and dislocations are observed to nucleate at the end of stage T3 wherein voids nucleate under the action of triaxial tensile stresses. The dislocation density shows a sharp increase during stage T4 attributed to void growth/coalescence to initiate spall failure. The calculated densities of various types of dislocations in the metal at the peak values of tensile pressure (onset of spallation) for the various laser fluences are tabulated in Table II. As will be discussed later, a direct correlation is observed between the spall strength and the dislocation density in the metal at the onset of spallation. Stacking fault intersections that result in the creation of stair-rod partials become potential void nucleation sites in the single crystal systems.¹⁷ As a result, spall strength values are highest when the microstructure shows the lowest values for densities of stair-rod partials in the metal. Higher density of stair-rod partials increases the sites for nucleation of voids and lowers the spall strength values. The results obtained for the 1 kJ/m² fluence validate this behavior wherein shock compression does not result in any plasticity in the metal and results in stair-rod partials with a significantly low density (and hence void nucleation sites) under tensile pressures rendering the highest spall strength values.

Illustrative snapshots of the systems showing the microstructure at a time of 120 ps for the various loading fluences are shown in Fig. 9. It can be seen that the laser interaction with the front end of the metal results in melting and spallation at the front end for laser fluences of 0.5 kJ/m² and 1 kJ/m², but only results in spallation closer to the rear surface for the 1 kJ/m² fluence. For the case of 6 kJ/m² and 13 kJ/m², ablation (explosive-like material removal) is observed at the front end of the metal and spallation at the rear end of the metal. The character of the front-end ablation and spallation processes observed here correspond well with those reported previously ⁵⁰ using MD-TTM simulations. In addition, the solid liquid interface is observed to proceed to larger depths in the metal for higher laser fluences. No spall (back-end) is observed for a fluence of 0.5 kJ/m² attributed to the interaction of the waves in the sample. A very thin spall region (region subjected to nucleation of voids) is observed for the 1 kJ/m² fluence, suggesting that the tensile pressures generated are not enough to nucleate voids in Al. The width of the spall region is observed to increase with laser fluence, and the location of the spall region is closer to the rear surface for larger values of laser fluence.

To investigate the effect of loading orientation, MD-TTM simulations are also carried out for laser shock loading and spall failure for the [110] and [111] orientations and compared with that observed for the [001] orientation. The same laser fluence (13 kJ/m²) and duration (τ = 150 fs) are used for this comparison. The predicted shock wave structures for the [110] and [111] orientations at various times are plotted in Fig. 10. The initial stages (20 ps) of laser material interaction are observed to result in a similar response (shock wave structure and melting) for the three orientations. The peak values of compressive pressure generated are calculated to be ∼50 GPa at a time of 5 ps for both the orientations and is observed to decay as the compressive wave travels through the metal toward the rear surface. One interesting aspect for the two orientations is that the shock wave structure transitions to a two-wave elastic-plastic structure at a time of 50 ps, which is shown in Fig. 11. Such a two-wave structure is not observed for the [001] orientation. The elastic (plastic) wave velocity is computed to be 8100 ± 100 (4300 ± 300) m/s and 9200 ± 300 (5700 ± 300) m/s for the [110] and [111] orientation, respectively. The two-wave shock structure and the corresponding distribution of the dislocations in the microstructure for the two orientations are shown in Fig. 11. These variations observed in shock wave structures for the three orientations are consistent with those computed using MD simulations obtained using a piston to shock the metal.²⁵ 

The pressure evolution in the system with time for the three orientations is shown in Fig. 12. The spall strengths are computed to be 6.7 GPa, 6.9 GPa, and 6.4 GPa for the [001], [110], and [111] orientations, respectively. The values for shock pressures, shock velocities, strain rates, and spall strengths are tabulated in Table III for the three orientations. In addition, the location of the spall plane from the rear surface (L_spall) is different for the three orientations and is attributed to the variations in the defect evolution behavior. The difference in the distance of the spall plane from the liquid phase as well as the heat propagation rate leads to a variation in the temperature at the onset of spall. Thus, a correlation between the strain rate of loading or pressure and the spall strength is not observed here, and a critical parameter is the temporal evolution of microstructure (defects) that determines the spall strength of the metal. A comparison of the variations in the temporal evolution of dislocation densities for the three orientations is shown in Fig. 13. It can be observed that the evolution of dislocation densities is very different for all the three orientations. A peak in the dislocation density is observed in T1 for the [110] orientation, whereas the peak value for the [111] orientation is observed in T3 in contrast to the peak in T2 for the [001] orientation. In addition, a significantly high density of stair-rod partials is observed for the [111] orientation in T3 before the peak tensile pressure is reached. The calculated densities of various types of dislocations in the metal at the peak values of tensile pressure (onset of spallation) for the three orientations are tabulated in Table IV. The spall strength values are correlated to the density of stair-rod partials at peak tensile pressures. For varying microstructures and deformation mechanisms under consideration and the same loading conditions, the highest spall strength values are observed for the lowest values for densities of stair-rod partials. Higher density of stair-rod partials increases the sites for nucleation of voids and lowers the spall strength values.

Illustrative snapshots showing the comparison of the spall failure microstructures of the three orientations are shown in Fig. 14. The voids are observed to be closest to the rear surface for the [001] orientation and farthest for the [110] orientations wherein void nucleation is observed to occur in the liquid phase for the same system dimensions and loading laser fluence and pulse duration. Such variations in spall failure phenomena for the three orientations suggest that a microstructural dependence of the spall strength needs to be investigated to enable predictive capabilities of damage tolerance of metallic materials.

The current results cannot be validated using experiments as no experimental data are found on ultra-short laser shock loading induced spall behavior of sc-Al at the length scales of the simulations. The comparison is therefore made with studies for experimental studies on polycrystalline Al thin film microstructures with a film thickness of 500 nm and the same loading conditions (a pulse of 150 fs and a fluence of 13 kJ/m²).¹³ The spall strength values predicted for all the three orientations considered here are slightly higher than the experimental values as is expected for single crystal systems. A more direct comparison with experiments will require MD-TTM simulations of the spall failure of polycrystalline Al microstructures. The predicted dimensions of the melting front compare very well with that observed in the experiment (∼200 nm of material undergoes phase transformation); however, the final front-surface damage is not investigated in this study. The predicted microstructural evolution, as well as detailed pressure and temperature evolution of Al systems, reported here provides new insights that can complement the experimental studies of laser shock loading and spall failure.

Large scale hybrid MD-TTM simulations investigate the microstructural evolution during laser shock compression of single crystal Al at the loading conditions and sample dimensions used experimentally. Laser interactions with the metal result in melting of the metal, and the predicted shock wave structures are observed to be significantly different from that predicted using MD simulations of shock using piston impact. The MD-TTM simulations are carried out to investigate the role of laser fluence and the role of loading orientation on the spall failure behavior of single crystal Al systems. A decrease in the laser fluence results in reduced depths of melting and lower peak pressures of the shock wave. A reduction in the laser fluence results in reductions in shock wave velocities and a reduced width of the spall region wherein voids are nucleated. The response of the shock wave propagation for the three loading orientations indicates variations in shock wave structures as well as defect evolution behavior in the metal. These variations in defect densities determine the spall strength values. The observed spall failure for the various laser fluences and the three loading orientations suggests that a correlation between shock pressure or loading strain rates and spall strengths does not exist and clearly demonstrates the need for a correlation between the microstructure (evolution of defects) and the spall strength values. This paper demonstrates that a correlation exists between the density of dislocations (stair-rod partials) in the microstructure and the spall strength of the metal. The spall strength values are highest for the lowest densities of stair-rod partials that serve as nucleation sites for voids. The higher the density of stair-rods formed in the spall region, the lower the spall strength of the metal. The predicted response of the metal compares very well with the experimentally observed melting behavior and slightly overpredicts the spall strengths. A more direct comparison of the spall strength will require the investigation of the response of polycrystalline Al microstructures under laser shock loading conditions in the future.

S. Galitskiy and A. M. Dongare acknowledge funding support from the US Army Research Office under Contract/Grant No. W911NF-14-1-0257. Any opinions, findings, and conclusions or recommendations expressed in this publication are those of the authors and do not necessarily reflect the views of the US Army Research Office or of the US Government. The US Government is authorized to reproduce and distribute reprints for Government purposes notwithstanding any copyright notation herein. D.S. Ivanov acknowledges the financial support due to grants DFG IV 122/4-1 and RFMEFI57816X0197, and computing time by Lichtenberg Computer Center, Darmstadt.