Dislocation evolution and peak spall strengths in single crystal and nanocrystalline Cu

The response of metals under dynamic loading conditions (shock) is very complex and involves plastic deformation, damage creation and evolution, phase transformation, heat generation and transfer, etc.¹ The deformation response during shock compression is determined by the ability of the microstructure to distribute load as the material deforms plastically. This ability is related to the evolution of the defects that comprise interfaces, stacking faults, dislocations, twins, vacancies, and interstitials. The interaction and accumulation of these defects and their distribution during shock compression generate a heterogeneity in the microstructure that determines the failure behavior of the material under the action of the reflected tensile waves (spallation). The spall failure of metals initiates through nucleation of voids in the weak regions created in this heterogeneous microstructure during shock compression and the voids grow and coalesce to form a fracture surface.² A fundamental understanding of the links between the microstructure and the nucleation, evolution, and accumulation of defect structures that initiate failure is therefore critical to design impact-tolerant structures with simultaneously increased toughness and strength under dynamic loading conditions.³ 

Most of the experimental efforts towards the design of impact-tolerant structures are aimed at determining the spall strength (peak tensile pressure prior to failure) of the material for a variety of microstructures and loading conditions (peak pressure and strain rate). The experiments are carried out using plate impact experiments that generate strain rates of 10⁵s⁻¹–10⁷s⁻¹,^4–7 or using laser shocks which generate strain rates up to 10¹⁰s⁻¹.^8–14 These studies suggest that the strain rates generated during shock-wave propagation have a significant effect on the spall strength of the metal. For example, the spall strength for copper (single crystal and polycrystalline) is observed to increase slowly (0.7 GPa–3 GPa) with strain rate in the range 10³s⁻¹–10⁶s⁻¹, whereas it increases rapidly (3 GPa to 9 GPa) with strain rate in the range 10⁷s⁻¹–2.5 × 10⁷ s⁻¹.¹⁰ The microstructure of the metal is also observed to have an effect on the spall strength of the metal.⁶ In addition, nanocrystalline metals with ultra-fine grain sizes (d ≤ 30 nm) show increased strengths during deformation at high strain rates (10⁴s⁻¹)¹⁵ and render ultra-high strength values under shock loading conditions.^16,17 The nanocrystalline metals, due to their high strengths, show significant promise for applications to experience dynamic loading environments. The understanding of the microstructural response (nucleation, evolution, and accumulation of defects) that is responsible for the high strengths at the nanocrystalline grain sizes and the effect on the failure behavior and spall strengths, however, is still in its infancy. The generation of these defect distributions during shock compression, on one hand, can be beneficial in design of high strength materials if the related collective stress-strain response results in strain hardening of the material. On the other hand, these defect structure distributions can be detrimental if they result in weak regions to nucleate voids that grow and coalesce to form cracks and initiate failure. Thus, the nature of the evolution of defects (through their interaction with each other) and their collective response determines the strength and toughness characteristics for a given initial microstructure.

Continuum modeling of deformation and failure of ductile metals has been carried out using void-growth models.^18,19 An example is the Gurson model¹⁸ that is based on a pressure dependent yield criterion. These void growth models have been employed in several studies^20–25 that mostly address the mechanism of void growth and porosity to model failure under quasi-static loading conditions as well as dynamic loading conditions.^26–29 The continuum models for dynamic failure do not specifically incorporate the necessary details relevant to the creation of the atomic scale heterogeneity in the material microstructure through the evolution of defect structures (dislocations, stacking faults, and twins) that determine the nucleation sites for voids, their growth, as well as the coalescence to form cracks. This heterogeneity becomes more important at the nanoscale dimensions of the microstructure, i.e., in nanocrystalline metals. Here, plastic deformation mechanisms depend on the interplay between dislocation and grain boundary (GB) processes,³⁰ whereas plasticity in polycrystalline metals is dominated by the nucleation, multiplication, and propagation of dislocations.³¹ This understanding is critical for the design of continuum models aimed at modeling of dynamic failure as the existing void nucleation and growth models do not account for the grain size as well as the evolving complex defect structures or defect densities in the material.

Molecular dynamics (MD) simulations enable the investigation of the atomic level evolution of the microstructure and have therefore been used to investigate the failure behavior of single crystal Cu^32–34 and nanocrystalline Cu^35–37 at high strain rates as well as under shock loading conditions. Voids are observed to nucleate at intersections of stacking faults generated in the spall plane under the triaxial tensile stress conditions for the single crystal Cu system. Void growth is observed to proceed through heterogeneous nucleation of dislocations from the void surface in these simulations. The voids in nanocrystalline Cu are observed to nucleate at grain boundary junctions and grow and coalesce along the boundaries under the action of the reflected tensile waves. While these studies are aimed at the mechanisms of nucleation and growth of voids, the understanding of the effects of microstructure on the evolution of defect structures during shock loading and spall failure is still in its infancy. For example, the nanoscale dimensions of the grain size promote deformation twinning as the grain size decreases and therefore results in enhanced strength and ductility of the metal.^38–40 A high density of twins in the microstructure at the spall plane experiencing the reflected tensile waves is therefore also likely to affect the spall strength of the nanocrystalline metal.

This study aims to investigate the effect of microstructure on the evolution of defects during deformation and failure under shock loading conditions. The evolution of defect densities is first investigated for a single crystal Cu system to identify the characteristics that render ideal strengths for the metal under shock loading conditions. To investigate the effect of the microstructure, MD simulations are carried out for nanocrystalline Cu systems with an average grain size of 6 nm, 10 nm, 16 nm, 20 nm, and 30 nm and compared with that observed for single crystal Cu. The computational details are presented in Section II. The evolution of dislocation densities and the related spall strength of single crystal Cu is discussed in Section III. The effect of grain size and the evolution of dislocation densities and on the spall strength of nanocrystalline Cu is discussed in Section IV.

The initial single crystal Cu system is created with dimensions of 25 nm × 25 nm × 100 nm and ∼5.3 × 10⁶ atoms. The initial nanocrystalline Cu systems are created with five different average grain sizes of d = 6 nm, 10 nm, 16 nm, 20 nm, and 30 nm using the Voronoi construction method.⁴¹ The system size in the lateral direction (X and Y) is scaled to be at least twice the grain size and has a minimum value of 30 nm for the system with a grain size of 6 nm and corresponds to ∼7.7 × 10⁶ atoms. The largest system corresponds to a grain size of 30 nm and has dimensions of 60 nm × 60 nm × 100 nm and corresponds to ∼31.7 × 10⁶ atoms. The system size in the shock (Z) direction for all the nanocrystalline systems is also set to be 100 nm. Periodic boundary conditions are used in the X and Y directions and the Z (shock) direction is kept free to have a surface at the front and back end of the sample. The as-created systems are first equilibrated at 300 K for 50 ps under conditions of zero pressure and the density of the final configuration is calculated to be ∼99% of the bulk density of copper. The MD simulations are carried out using the Voter-Chen (VC)⁴² formulation of the embedded atom method (EAM) potential for copper. The VC potential is well suited to describe mechanical properties of nanocrystalline Cu as it provides a good description of the unstable and stable stacking fault energies, the surface and grain boundary energies,^43,44 as well as the spall strengths of single crystal and nanocrystalline Cu.^34–37 In addition, the dislocation nucleation energies predicted by the EAM potential show a very good agreement with the continuum models.⁴⁵ 

To achieve shock-induced deformation, the first 3 nm of the sample was chosen as the piston and the atoms in the piston were given an impact velocity $Up$⁠. The red colored atoms at the left end of the sample (as shown in Figure 1) correspond to the piston that is driven inward at the desired velocity as shown by the red arrow. All the MD simulations of impact reported here are performed using a constant piston velocity of 1000 m/s and a square pulse of 10 ps. Rear surface velocity profiles are computed by averaging the velocity of the green colored atoms in a width of 3 nm at the back end of the system as shown in Figure 1. The elements of the atomic-level stress tensor are calculated as

where α and β label the Cartesian components, Ω₀ is the atomic volume, F_ij is the force on atom i due to atom j, M_i is the mass of atom i, and v_i is the velocity of atom i. The pressure in the system is calculated as $P=−13(σxx+σyy+σzz)$⁠, where $σxx$⁠, $σyy$⁠, and $σzz$ are the stresses averaged over the entire system in the X, Y, and Z directions, respectively. The time step for all of the MD simulation runs was chosen to be 2 fs. No thermostat was applied during the MD simulations.

The characterization of defects (stacking faults, twins, and voids) is carried out using common neighbor analysis (CNA)⁴⁶ and the centro-symmetry parameter (CSP),⁴⁷ and the evolution of the dislocation lengths and densities are quantified using the dislocation extraction algorithm (DXA).^48,49 Based on the computed Burgers vectors, the identified dislocations are classified into perfect dislocations (⁠$a2⟨110⟩$⁠), Shockley partials (⁠$a6⟨112⟩$⁠), Frank partials (⁠$a3⟨111⟩$⁠), Stair-Rod (⁠$a6⟨011⟩$⁠), and Hirth (⁠$a3⟨001⟩$⁠) dislocations. The stacking faults and coherent twin boundaries in FCC crystals consist of atoms with HCP-like coordination, which can be detected using the CNA. However, to distinguish between these types of planar defects, an extended structure identification method is required. An algorithm implemented in the crystal analysis tool (CAT)⁵⁰ is therefore used to identify different planar defect types. Deformation twinning in FCC metals is attributed to glide of partial dislocations with the same Burgers vector on successive planes.⁵¹ Twinning dislocations in FCC metals therefore have a Burgers vector of $a6⟨112⟩$ and hence cannot be distinguished from Shockley partials just based on the Burgers vector alone. To resolve this ambiguity, the dislocation lines extracted by the DXA are mapped to the individual atoms that form their cores. The twinning dislocations are then discriminated from the Shockley partial dislocations as the ones that possess nearest neighbor atoms forming coherent twin planes. The total length (⁠$ld$⁠) of each type of dislocation is then computed as the sum over all corresponding dislocation segments for each snapshot generated during the simulation. The dislocation density (⁠$ρd$⁠) for each dislocation type is then computed by dividing by the volume of the system as $ρd=ldV$⁠.

The MD simulations are carried out by giving a constant inward velocity of $Up$ = 1 km/s to the atoms in the piston as shown in Figure 1 along the Z axis for 10 ps (square pulse). The inward impact results in a planar shock wave that travels with a velocity (U_s) towards the rear surface. To understand the temporal response of the metal to the shock wave propagation, the system is divided into sections along the Z axis and the values of pressure are averaged for all the atoms in each section. A contour plot of the evolution of pressure (P) as a function of time along the sample in the shock direction is shown in Figure 2(a). A positive value corresponds to compressive pressure and a negative value corresponds to tensile pressure. The peak shock pressure is calculated to be 48.9 GPa and the pressure in the material behind the shock front also remains constant. A tail of the pressure wave is observed at the end of the shock pulse (10 ps) that propagates to the rear end. The constant pressure shock wave reaches the rear surface at ∼19 ps and reflects back as a tensile wave. The reflected tensile wave interacts with the tail of the initial pressure wave as shown by the intersecting black arrows in Figure 2(a). The interaction results in an increase in the tensile pressure to a peak value at ∼30 ps that results in the nucleation of multiple voids to initiate spall failure. This peak value of the tensile pressure is defined as the spall strength of the material and is computed to be 12.2 GPa for single-crystal Cu for an impact velocity of 1 km/s and a square pulse of 10 ps. The propagation of the tensile pressure wave towards the front surface results in growth and coalescence of the voids to form the fracture surface. The corresponding rear surface velocity (in the shock direction) showing the various stages is plotted as a function of time in Figure 2(b).

The response of the metal can be classified into four stages: Stage I (SI) corresponds to shock compression of the metal for a pulse duration of 10 ps; stage II (SII) corresponds to the propagation of the shock wave from the end of the pulse till it hits the rear surface (∼19 ps); stage III (SIII) corresponds to the propagation of the reflected tensile wave till a peak tensile pressure is reached due to the interaction with the tail of the pressure wave (∼30 ps); and stage IV (SIV) corresponds to the onset of spall and the nucleation, growth, and coalescence of voids till the tensile wave reaches the front surface (∼ 40 ps). The four stages defined here are indicated in the pressure contour and the rear surface velocity profile in Figure 2. The snapshots generated at various times during the simulations are analyzed to identify various dislocations and their lengths are characterized using DXA. The coordinates of the dislocation lines are then mapped onto the coordinates of the atoms to identify the twinning partials. The snapshot of the system showing the evolution of microstructure of single crystal Cu at times of 10 ps (at the end of SI), 20 ps (beginning of SIII), 30 ps (at the end of SIII), and 40 ps (at the end of SIV) are shown in Figures 3(a), 3(b), 3(c), and 3(d), respectively. The color of the atoms corresponds to bulk FCC (light green), disordered (blue), stacking fault (green), twin fault (yellow), the twinning partial dislocation (light blue), and a surface (orange) as shown by voids.

The total dislocation density during the four stages of the simulation is plotted as a function of time in Figure 4 along with the dislocation densities of all types of dislocations. The dislocation density during stage I (shock compression) rapidly increases with time and is initially dominated by Shockley partials. The density of twinning partials increases at a slightly slower rate initially and increases more rapidly to almost equal the density of Shockley partials at the end of the shock pulse. Stage I also has a small contribution from Hirth dislocations. During stage II, a slight crossing over is observed in the density of twinning partials and that of Shockley partials. The twinning partials continue to increase at a slower rate till the shock wave hits the rear surface, whereas the Shockley partials begin to decrease slowly after the cross-over point. The density of Hirth dislocations is also observed to decrease during this stage. This decrease in density of the Hirth and Shockley dislocations is attributed to the propagation of the tail of the pressure wave at the front end of the sample towards the rear surface. As a result, the total density of dislocations as well as the density of twinning partials are observed to increase at a slower rate in stage II as compared with stage I and reaches a peak value when the shock wave reaches the rear surface. The reflection of compression wave as tensile wave (stage III) corresponds to a rapid decrease in all the densities of the dislocations and hence in the total dislocation density under the action of the reflected tensile wave and the tail of the pressure wave. The interaction of these two waves drives the material under tension and reaches a peak value for the tensile pressure that results in the nucleation of voids. The dislocation density reaches a minimum value with an equal density of Shockley partials and twinning partials at the peak value of the tensile pressure (spall strength). Stage IV corresponds to the nucleation, growth, and coalescence of voids and results in an increase in the dislocation density. This increase in the dislocation density is attributed to the nucleation of dislocations from the void surface during growth and coalescence of the voids^32,34 and is observed to be dominated by twinning partials and Shockley partials with a slightly greater density of twinning partials. The snapshot of the single-crystal system showing evolution of twins at times of 10 ps (at the end of SI), 20 ps (beginning of SIII), 30 ps (at the end of SIII), and 40 ps (at the end of SIV) along with the evolution of voids are shown in Figures 5(a), 5(b), 5(c), and 5(d), respectively. Only the atoms that correspond to the twin fault (yellow), the twin partial dislocation (light blue), and a surface (orange) are shown here.

Thus, these results demonstrate that the deformation response of single crystal Cu is dominated by the evolution of twinning partials and Shockley partials. The onset of spallation (nucleation of voids) in single crystal Cu occurs in a microstructure characterized by an equal density of Shockley and twinning partials. This equal density of twinning partials and Shockley partials renders high values of spall strength for the single crystal metal. It can be seen that the evolution contributions to dislocation density for single-crystal Cu has very little contributions from Hirth, perfect, and stair-rod dislocations during shock compression.

To study the effects of grain size of the nanocrystalline metal on deformation and failure response, the simulations were also performed for nanocrystalline Cu with an average grain size of 6 nm, 10 nm, 13 nm, 16 nm, 20 nm, and 30 nm. Contour plots of the evolution of pressure (P) as a function of time along the sample in the shock direction for the various grain sizes are shown in Figure 6 for a constant impact velocity (10 ps square pulse) of $Up$ = 1 km/s along the Z axis. The pressure profiles appear to be similar for all the grain sizes and the four stages are observed to occur at the same times as observed for the single crystal system. The values for the spall strength are computed using the pressure profiles and plotted in Figure 7 as a function of the average grain size in comparison with that for single crystal Cu (thick green line). It can be seen from the plot that the spall strength increases as the grain size decreases from 30 nm to 16 nm (Hall-Petch behavior) reaching a peak value of 11.92 GPa for an average grain size of 16 nm. A further decrease in the grain size transitions the system into the inverse Hall-Petch regime and hence the spall strength is observed to decrease for the 13 nm grain size system. A continued decrease in the grain size, however, results in an increase in the spall strength of the metal. This increase in spall strength is attributed to limited grain boundary sliding^9,16,17 and reaching a maximum of 11.82 GPa for an average grain size of 6 nm. These values for the 6 nm and the 16 nm system compare very well to the spall strength of 12.2 GPa for single-crystal Cu at an impact velocity of 1 km/s. Thus, it can be seen that the spall strength of the nanocrystalline metal reaches near-ideal strength values at grain sizes of 6 nm and 16 nm.

To understand the evolution of microstructure that renders the variations in the spall strength values (as plotted in Figure 7), the total dislocation density along with the densities of twining partials, Shockley partials and perfect dislocations is plotted as a function of time in Figure 8 for the various grain sizes considered. It can be seen from the plots that, in contrast to single crystal Cu, the nanocrystalline Cu systems have a substantial initial dislocation density due to dislocations distributed along the grain boundaries. These dislocations at the grain boundary are characterized to be perfect dislocations and Shockley partials. Stage I corresponds to an increase in the density of Shockley partials and twinning partials and hence an increase in the total dislocation density as observed for the single crystal system. The density of perfect dislocations, however, is observed to continuously decrease during the simulations. The evolution of dislocation density during stage II is observed to increase slowly as this stage is accompanied by the propagation of the tail of the pressure wave. The increase in dislocation density, however, is also observed to vary with the grain size of the metal. Depending on the evolution of the densities of the twinning partials, the peak value is either observed at the end of stage II (time when the shock wave hits the rear surface) or in the middle of stage II. The density of twinning partials is observed to be the highest for the case of the 16 nm followed by the 6 nm system at the end of stage II. The reflection of the compression wave as tensile waves in stage III results in a decrease in the density of all types of dislocations. The nucleation, growth, and coalescence of voids during stage IV is observed to occur along the grain boundaries and hence is accompanied by a decrease in the densities of all types of dislocations. A closer look at the densities of the twinning partials during stage IV for various grain sizes considered shows that only the 16 nm system shows a local minimum in the density of twinning partials at the onset of spallation as observed for the single crystal system. All the other grain sizes show a steady decrease in the density of twinning partials during stage IV. In addition, the dislocation density of the 16 nm system also has the highest value among all the grain sizes considered followed by the 6 nm system.

The variation of the spall strengths as a function of grain size as plotted in Figure 7 can be understood based on the evolution of the density of twinning partials for the various grain sizes. A comparison of the evolution of twinning partials during shock loading of nanocrystalline Cu using an impact velocity of 1 km/s for various grain sizes is shown in Figure 9 in comparison with that for single crystal Cu. It can be seen from the plots that while the density of twinning partials for nanocrystalline Cu with an average grain size of 6 nm and 16 nm are higher than that for single crystal Cu, the density has contributions from dislocations distributed at the grain boundary. The peak strengths are more likely to be determined by the evolution of twinning partials in the grain interior as for the case of single crystal Cu as compared with that at the grain boundary for the nanocrystalline metal. A closer look at the evolution of twinning partials for the 6 nm system shows that the density of twinning partials at the onset of spallation is very close to the initial density of twinning partials that are located at the grain boundaries. These dislocations in the initial as-prepared structure are likely part of growth twin boundaries or grain boundaries whose orientation relationship is close to the twin configuration and which relax to a structure with coherent twin facets. In contrast, the density of twinning partials for the 16 nm system at the onset of spallation is much larger than the initial density of twinning partials located at the grain boundary suggesting a higher density of deformation twins in the grain interior. The high density of twinning partials either at the grain boundary (as for the 6 nm system) or in the grain interior (as for the 16 nm system) affects the peak tensile pressure require to nucleate voids at the grain boundaries and hence render higher values for the spall strength. Illustrative snapshots of the failure of the nanocrystalline Cu system with an average grain size of 6 nm and 16 nm showing the growth and coalescence of voids along the grain boundaries at 40 ps (at the end of SIV) are shown in Figures 10(a) and 10(b), respectively. The corresponding microstructure of the 6 nm and the 16 nm system showing the distribution of twin faults (yellow atoms) and twinning partials (light blue atoms) at 40 ps (at the end of SIV) shown in Figures 11(a) and 11(b), respectively. The density of twin faults in the grain interior for the nanocrystalline systems, however, is anticipated to be less than the twinning partial densities in single crystal Cu and hence the spall strength is observed to be slightly less than that for single crystal Cu.

Thus, the spall strength of single crystal and nanocrystalline Cu is determined by the evolution of twinning partials and twin faults in the microstructure at the spall plane. The high density of twinning partials (twins) for the 16 nm system and the 6 nm system render near-ideal strengths as observed for single crystal Cu. Future efforts, therefore, need to be directed to increase the ability of the microstructure to twin in order to achieve ideal strengths for the metal. This ability to control the evolution of the density of twinning partials (twins) will provide pathways for the design of high strength and failure resistant metals. A similar trend in dislocation densities has been observed for Cu using another EAM potential⁴⁴ and the trend in the evolution of dislocation density is very similar to that reported here for the Voter-Chen EAM potential. A detailed comparison of the dislocation densities will be the focus of a future manuscript.

Large scale MD simulations were carried out to investigate the evolution of microstructure and defects in single crystal and nanocrystalline Cu under shock loading conditions. The results suggest that the ideal spall strength in the single crystal Cu system is attributed to the high density of twins in the material experiencing the peak tensile pressures. In contrast, the evolution of twinning partials is observed to vary with the grain size of the nanocrystalline metal. The evolution of twins for the nanocrystalline metal during stage II (propagation of the shock wave towards the rear end of the sample) and stage III (reflection of the waves as a tensile wave and the interaction with the tail of the pressure wave) is observed to vary with the grain size of the system. The density of the twinning partials is observed to be the highest for the 16 nm and the 6 nm system among all the grain sizes considered and hence have near-ideal spall strengths as observed for single crystal Cu. The higher strength for the 16 nm system is attributed to a higher density of twins in the grain interior as compared with the 6 nm system that is likely to have a larger fraction of coherent twin boundary segments due to the smaller grain size. The high density of twins at the grain boundary and/or grain interior for these two grain size systems affects the peak tensile pressures required to nucleate voids at the grain boundaries and hence render near-ideal spall strengths for nanocrystalline Cu.

This work acknowledges the support from NSF CMMI under Grant No. 1454547 and K.M. acknowledges the support in part from the Student Research Participation Program at the U.S. Army Research Laboratory (USARL) administered by the Oak Ridge Institute for Science and Education through an interagency agreement between the U.S. Department of Energy and USARL. Computational support through the Department of Defense (DOD), High Performance Computing Modernization Program (HPCMP), and the Booth Engineering Center for Advanced Technology (BECAT) is also acknowledged.