Experimental and theoretical examination of shock-compressed copper through the fcc to bcc to melt phase transitions

Copper is a widely used material; many industrial and technical applications exploit its high electrical and thermal conductivity. Within the high-pressure community, Cu is used as a pressure standard for x-ray diffraction (XRD) experiments. The copper volume, derived from shock wave studies, has been used to calibrate other pressure standards, including gold¹ and the well-known ruby fluorescence scale.² As a consequence, the behavior of copper at extreme pressures (⁠$P$⁠) and temperatures (⁠$T$⁠) has been widely studied both experimentally and theoretically for many years.

Different theoretical approaches have been used to explore the high $P$-$T$ equation-of-state (EOS) and phase stability of copper. A study by Neogi and Mitra,³ using density functional theory (DFT) and embedded atom model (EAM) molecular dynamics simulations, predicted a transition from the ambient fcc phase to a body-centered (bcc or bct) structure at pressures in excess of 80 GPa and temperatures above 1520 K. Recent ab initio studies by Smirnov,⁴ using the quasiharmonic approximation, and Baty et al.,⁵ using the inverse-Z method, predict a thermodynamic transition from fcc to bcc along the Hugoniot at 180 GPa at about 4000 K. These studies further confirm the bcc phase to be the high-temperature solid phase of Cu below the melt line at pressures above 80 GPa. Also, the classical molecular dynamics study of solidification of Cu liquid by Sadigh et al.⁶ found crystallization into the bcc phase at this pressure range.

Several experiments have investigated the evolution of crystal structure in copper along different thermodynamic compression paths. Room-temperature static-compression experiments by Dewaele et al. report that copper remains in the ambient fcc phase up to at least 153 GPa.⁷ More recent laser ramp-compression experiments by Fratanduono et al.,⁸ combined with in situ XRD, have shown the fcc structure to be stable to at least 1150 GPa (1.15 TPa). Along the higher temperatures accessed by shock compression, the fcc $→$ bcc transformation in Cu was first observed with in situ XRD measurements at $∼182$ GPa by Sharma et al.⁹ This apparently correlated with an earlier study on shock-compressed Cu by Hayes et al.,¹⁰ which reported a distinct kink at 185 GPa in the slope of the measured Poisson’s ratio vs shock pressure. The XRD study of Sharma et al.⁹ utilized diffraction line profile analysis to infer that the production of stacking faults was a precursor to the onset of the bcc transition. It was proposed that the observed fcc stability under ramp-compression to TPa pressures was due to a strain-rate effect, and that a critical loading rate may be required for the production of sufficient stacking faults to instigate transformation into the bcc phase.⁹ However, several theoretical calculations actually predict that the bcc phase is stable only at high temperatures,^3–6 consistent with the observation of the fcc phase stability at TPa pressures under low-temperature ramp compression.⁸ 

Here, we report on a combined experimental and theoretical approach to further explore the fcc $→$ bcc transition mechanism in Cu, and the P-T stability of the bcc phase. Cu foil samples were shock-compressed with a high-powered laser driver, and in situ XRD was employed to measure the evolution of micro-structural texture associated with the transformation to the high-pressure bcc phase. For metals, the transformation pathway (atomic mapping from the ambient to high-pressure crystal structure) associated with the martensitic fcc to bcc phase transformation has been studied for the last hundred years, and there are several reported mechanisms.^11–17 To date, no constraint of this transformation has been established during nanosecond shock compression. In this work, we find that the change in sample texture from the ambient fcc phase to the high-pressure bcc phase under shock compression is consistent with a martensitic transformation mechanism through the Pitsch distortion.^15,17 To determine the $P$-$T$ extent of the bcc phase, the fcc–bcc phase boundary was calculated with EAM molecular dynamics simulations by performing a thermodynamic integration (TI) using the Frenkel and Ladd technique at high pressures and high temperatures.¹⁸ These calculations confirm that the bcc phase is stabilized at high temperatures only. Based on the melting line obtained from two-phase (solid–liquid) coexistence simulations in the isobaric–isoenthalpic (NPH) ensemble following the procedure in Ref. 6, the bcc phase is predicted to persist to pressures greater than 500 GPa.

Laser-driven shock-compression experiments were conducted at the Dynamic Compression Sector (DCS) at Argonne National Laboratory.¹⁹ The experimental geometry is illustrated in Fig. 1. The target assembly consisted of a 37.5 $μ$m layer of polyimide, a 22.5 $μ$m thick polycrystalline copper foil (supplied by Goodfellow Inc.), and a 500-$μ$m thick single-crystal LiF [100] window with a 0.1 $μ$m Al coating to enhance the reflectivity of the Cu–LiF interface. All target components were held together with epoxy (⁠$∼1$ $μ$m thick). A 5 or 10 ns flat-top laser pulse from the 351-nm DCS drive laser of energies between 39 and 80 J was focused within a 500 $μ$m diameter focal spot to uniaxially compress the sample. This setup generated longitudinal stress states in Cu between 120 and 280 GPa.

The U17 undulator was used to generate x rays with a peak flux at $∼23.5$ keV (⁠$λ=0.5275$ Å) and a bandwidth of $∼3.5$%, as shown in Fig. 1. The scattered x rays were collected in transmission geometry on a Rayonix SX165 flat panel CCD detector with the x-ray beam aligned 52$°$ or 72$°$ from target normal. The diffraction measurements reported in this work were obtained using a $∼80$ ps duration x-ray pulse, with an on-sample focal spot size of $∼30$ $μ$m (vertical) $×$ 80 $μ$m (horizontal).¹⁹ The x-ray pulse was timed to probe the Cu sample during shock transit within the Cu layer. Therefore, each diffraction image samples both the compressed volume and the ambient volume (the region ahead of the shock front). A representative pre-shot diffraction image from an ambient-pressure Cu foil is shown in Fig. 1.

During each shock-compression experiment, a velocity interferometer system for any reflector (VISAR) was used to measure the particle velocity history (⁠$u$(t)) at the Cu/LiF interface. A 532-nm point-VISAR laser was focused to a $∼80$-$μ$m diameter spot on the Cu/LiF interface, and it was carefully aligned so that the time-history was recorded from the same location as the x ray and drive laser focal spots. To determine the correct interface velocity from the apparent velocity measured by the raw VISAR, a correction for the refractive index of the LiF window was applied.^35,36 A representative velocity profile from the setup is shown in Fig. 1 and several additional profiles are shown in Fig. S11 in the supplementary material.

Standard impedance matching techniques were used to determine the longitudinal stress in the Cu layer, from the measured Cu/LiF $u$(t), and linear $US$-$u$ (shock velocity–particle velocity) fits to previously published shock-compression data for LiF³⁶ [Eq. (1)] and Cu^20–32 [Eq. (2)],

We calculated the uncertainty in stress as the result of contributions from three major sources: (i) the distribution of stress states due to non-ideal temporal steadiness of the shock wave creates substantial error; the mean and standard distribution of velocity states above the initial shock were included in the impedance matching calculation (see Fig. S11 in the supplementary material); (ii) the accuracy to which fringe shifts can be measured in the point-VISAR system, which is taken here to be 0.024 km/s (2% of a fringe shift); and (iii) the uncertainty in the LiF and Cu Hugoniots. For some shots, an additional stress uncertainty was added to account for the expected spatial distribution of states associated with laser-driven compression over the 500 $μ$m diameter spot (see Fig. S12 in the supplementary material and associated discussion). Other contributors to the stress uncertainty, which are considered negligible, include the uncertainties in the refractive index of LiF, uncertainties in the timing of the x-ray probe with respect to the VISAR and uncertainties in the measurements of sample thickness.

Eight shock-compression experiments were performed which ranged in stress from 120 to 280 GPa. Sections of raw x-ray diffraction images from representative experiments are shown in Fig. 2(a), along with stress determinations based on VISAR analysis. Before structural analysis, azimuthal integration and several intensity corrections were performed on the raw XRD images using the HEXRD software package.³⁷ First, powder diffraction rings from known standards (CeO$2$ and Si) were used for detector calibration. The pink beam profile function detailed in Ref. 21 was used for this step. Figure 2(b) shows the intensity-corrected and azimuthally averaged diffraction profiles for shock-compressed Cu. The intensity corrections included: (i) subtracting the dark counts, (ii) subtracting the ambient Cu pattern and subtracting LiF single-crystal Laue peaks for all experiments, (iii) correcting for the solid angle subtended by each pixel on the detector, (iv) accounting for the polarization factor (for liquid diffraction profiles only), and (v) correcting for the attenuation due to varying path length of the diffracted x rays.

For these profiles, the XRD contribution of the unshocked Cu volume was subtracted from the shock-compressed XRD data, so only the compressed peaks are displayed. Here, line profiles are shown in momentum (Q) space, where $Q=[4πsin(θ)]/λ=2π/d$⁠, 2$θ$ is the diffraction angle, and $d$ is the atomic d-spacing. For each line profile, the shot number and shock stress are listed. Using the HEXRD software, as described in Ref. 33, we employ LeBail refinement to fit the diffraction data, and determine the lattice parameter and stacking fault probability in the fcc phase, during shock compression. For all shots, the fits to the measured diffraction profiles (including the contribution from the ambient uncompressed volume) are shown in Figs. S7 and S8 in the supplementary material. Based on the determined lattice parameters (see Table S1 in the supplementary material), we fit four fcc peaks for stresses below 185 GPa. At 185(20) GPa, we observed a mixed fcc + bcc phase. For pressures from 197(11) to 220(16) GPa, we only observed the bcc phase. At 280 GPa, only liquid diffraction is observed. This phase evolution is consistent with the published work of Sharma et al.⁹ The raw XRD data for all experiments shown are included in Figs. S2–S6 in the supplementary material. The experimentally determined stress-density data are plotted along with previous shock-compression studies^9,20–32 in Fig. 3.

Analyzing the azimuthal intensity distributions from two-dimensional XRD images can provide constraints on high-pressure deformation processes³⁸ and transformation pathways for any phase transitions, from ambient to high-pressure crystal structures, that occur.^39–42 Textured samples have a non-random distribution of crystallite orientations within the bulk sample, which results in non-constant azimuthal intensity distributions around individual Debye–Scherrer cones. Analysis of these intensity variations can provide useful information about the crystallite orientations within the ambient and shocked samples. The orientation relationships (ORs) between the ambient and high-pressure crystal structures is often identifiable when the ambient or shocked samples are single crystal or textured.^39–41 Although past experiments on copper used x-ray diffraction profiles to examine potential microscopic deformation mechanisms in the fcc phase,^9,43 to date, real-time in situ methods have not been used to evaluate previously proposed fcc to bcc orientation relations and transformation mechanisms during shock compression of copper.

To evaluate relationships between crystallite and lattice plane orientations in the ambient and shocked samples, we employed two separate techniques: (i) a ray-tracing analysis using the two-dimensional XRD images and (ii) a Rietveld model to estimate texture.

(i) The ray-tracing process consisted of two general steps. First, we identify the peak intensity locations (in polar coordinates) on each 2D x-ray diffraction image. Peak locations as a function of azimuthal (⁠$ϕ$⁠) angle, the angle around the Debye–Scherrer cones, were determined by integrating over a defined radial range using the cake function from Fit2D.⁴⁴ The function produces an intensity vs azimuthal angle lineout.

In the present analysis, the radial integration 2$θ$ ranges were chosen to isolate diffracted intensity from a single {hkl} family of lattice planes (ambient or shocked) while minimizing diffracted intensity from other lattice planes in the ambient or shocked material. Due to the broad x-ray spectra, however, it was impossible to completely avoid diffracted intensity from other lattice planes in some cases. The potential effects of ambient diffracted intensity on the azimuthal peak analysis of shocked diffraction peaks were explicitly examined and evaluated, as described in Sec. S4 in the supplementary material. A multi-peak fitting routine was performed using Igor Pro (WaveMetrics, Inc., Lake Oswego, OR, USA) to locate the azimuthal centers of the diffracted intensity. Radial positions (or 2$θ$ scattering angles) were determined using Dioptas.⁴⁵ 

Second, after establishing individual peak locations on the 2D XRD images, we determined the lattice plane normal directions consistent with the observed XRD peaks. The observed intensity on the 2D images fixes the scattered x-ray direction, $S^$⁠, based on the sample to detector distance and other detector parameters. The experimental geometry also establishes the incident x-ray direction, $S^0$⁠. Thus, using the peak wavelength, $λ$⁠, from the experimental flux spectrum (see Fig. 1), we can establish the reciprocal lattice vector, $Hhkl$⁠, necessary to satisfy the diffraction condition⁴⁶ given in Eq. (3),

Since a reciprocal lattice vector is normal to the corresponding lattice plane, Eq. (3) also defines the lattice plane normal directions. Generally, peak intensity locations from raw XRD images can be correlated to suggest mutual orientation of lattice planes. However, the analysis method used here provides a more rigorous determination of the angle between the lattice planes in ambient and shocked samples, removing possible ambiguities.

(ii) The Rietveld model was used to determine the orientation relationship between the fcc and bcc phases. Instead of applying the traditional axis-distribution functions employed in most Rietveld packages, we employed the discrete harmonics to explicitly enforce non-negativity in the calculated orientation distribution function (ODF).⁴⁷ The pink beam powder diffraction profile described in Ref. 48 was used to model diffracted peak shapes.

Figures 4 and 5 summarize the ray tracing analysis for representative experiments below and above the phase transformation stress, respectively. In each experiment, the starting material had significant texture, which is typical of rolled foils (see Figs. S1–S6 in the supplementary material). Figures 4 and 5 show the 2D $ϕ$ vs 2$θ$ scattering angle intensity patterns (known as cakes) from ambient and shock-compressed copper [(a) and (b), respectively]. In Figs. 4(b) and 5(b), 37$%$ and 19$%$⁠, respectively, of the ambient images were subtracted from the corresponding shocked images to highlight diffracted intensity from the shock-compressed material. The residual ambient fcc (111) peak intensity is due to the non-uniform background subtraction resulting from the depth-dependent sample microstructure. This issue and the potential effects on the azimuthal peak analysis are described in Sec. S4 in the supplementary material. By integrating over a small 2$θ$ range (as highlighted in the blue and red boxes), we obtain intensity vs $ϕ$ lineouts for specific hkl families of lattice planes, as shown in Figs. 4(c) and 5(c). The peak locations, determined from the 2$θ$ and $ϕ$ locations of the diffracted intensity spots shown in Figs. 4 and 5, were then used to determine the corresponding plane normal directions as described above.

For experiments below the phase transformation stress, e.g., Fig. 4, our analysis shows that the ambient and shocked fcc lattice planes remain approximately parallel, as indicated by the shaded regions in Fig. 4(c). For experiment S19C137 shown in Fig. 4, lattice planes in the {111} and {220} families of lattice planes are rotated by less than 4.5$°$ upon shock compression. Significant crystallite re-orientation seems highly unlikely based on our data because peak widths did not notably change and only slight grain re-orientation was observed across multiple hkl families of lattice planes. Slight re-orientation, indicated by azimuthal intensity peak splitting, was commonly observed for the (111) XRD peaks. This peak splitting is shown in Figs. 4(b) and 4(c). If more significant grain re-orientation did occur, as through twinning or some other mechanism, individual diffraction peaks would shift more drastically or disappear as the diffraction condition was broken. Diffraction peak shifts consistent with twinning of the fcc structure were not observed in this work.

The atomic mapping from an ambient fcc to high-pressure bcc crystal structure has several orientation relationships reported.^11–17 A review of the different proposed pathways is given in Ref. 17. To date, no constraint of this transformation has been established during ns shock-compression experiments. Our texture analysis provides the first such constraint.

The fcc $→$ bcc transition induced during the rapid quenching of high-temperature austenite (fcc) phase to the martensite (bcc) phase has been studied extensively in steels and other iron alloys. Many different ORs between the fcc and the bcc phase are reported in the literature, with the Kurdjumov–Sachs (KS) OR being the most commonly observed OR in fcc $→$ bcc transformation during rapid quenching experiments. However, no such OR measurement has been previously reported during ns-shock experiment.

There are a number of proposed pathways for the fcc $→$ bcc transformation: (1) brain distortion,¹¹ (2) the Kurdjumov–Sachs–Nishiyama shear model,⁴⁹ (3) the Bogers and Burgers model,⁵⁰ and (4) the one step Pitsch-distortion model.¹⁷ This section evaluates the consistency of our measurements with the predictions from these models.

For experiments above the phase transformation stress (Fig. 5), a notable correlation was observed; the ambient fcc (200) plane was oriented approximately parallel to the shock-formed bcc (110) plane. The gray regions in Fig. 5(c) highlight the general correlations between the diffracted intensity from the fcc (200) and bcc (110) lattice planes, a pattern observed across the three experiments performed above the phase transformation stress. To quantitatively determine the relative plane angles, we used the ray tracing approach described above. Due to the broad x-ray flux spectrum, diffracted intensity from bright fcc (111) peaks can overlap with fcc (200) peak intensity, artificially shifting the azimuthal peak locations. This has been accounted for as described in Sec. S4 in the supplementary material. For a quantitative comparison of the lattice plane angles, we analyze the peaks near 90$°$ and 270$°$—(200) peaks which are both more intense and minimally affected by ambient (111) intensity. By comparing these peaks for the three experiments, we determined an average lattice plane angle between the fcc (200) and bcc (110) planes of $0.88±0.52°$⁠. Because slight plane rotations are expected for some material deformation mechanisms,⁵² and well-defined orientation relations likely exist between the fcc and bcc phases,⁵³ slight misorientation between lattice planes between the different structures is not unexpected.

In the pathway analysis discussion, we use the standard notation for martensitic transformations, where the bcc phase is referred to as the $α$ phase, and the fcc phase is referred to as the $γ$ phase.¹⁷ From the XRD data in Fig. 5, we calculate the inverse pole figures (IPFs) for the ambient fcc and shock-compressed bcc phases. These IPF maps are shown in Fig. 6(c) and represent the distribution of crystallographic axis aligned with a particular sample axis. We have chosen the z or normal direction (ND) parallel to the x-ray axis, and x axis or the rolling direction (RD), which is perpendicular to ND. An OR of $[001]γ||[011]α$ and $[011]γ||[111]α$ is obtained, consistent with transformation into the high-pressure bcc phase along a pathway described by the Pitsch OR^15,17 [see Figs. 6(a) and 6(b)].

The Pitsch distortion is defined as a shear of the $(001)γ$ plane to the $(1¯10)α$ plane, which changes the angle between the diagonals from $90°$ to $70.5°$ while leaving one of the diagonal unchanged [see the $[001]γ$ projection in Fig. 6(a)]. The un-rotated and un-strained diagonal is the neutral $[110]γ||[111]α$⁠. This shear couple with a hydrostatic strain in the $[001]γ$ direction results in the final bcc unit cell. The principal strain values for this mechanism are $0%,−5.7%,$ and $15.5%$⁠. Notably, no atomic shuffle is required for this mechanism. Readers are referred to Ref. 17 for more details about this mechanism. This is consistent with the average lattice plane angle between the fcc (200) and bcc (110) planes of $0.88±0.52°$ reported earlier. It is important to note that due to the high temperature during shock experiments, the higher angle peaks have relatively low intensity, which can affect the quality of texture analysis. In addition, the difference in orientations between the different ORs are small (⁠$<6°$ to Kurdjumov–Sachs OR¹⁷) and multiple variants can be observed in a single experiment. Furthermore, there is a continuous pathway connecting the orientation relationship observed during Pitsch distortion to the Kurdjumov–Sachs (KS) and Nishiyama–Wasserman (NW) orientation relationship. However, there is additional evidence for the Pitsch distortion being active as described below.

Two additional experimental constraints provide evidence for the Pitsch-distortion mechanism. The angle between the close-packed plane of the two phases $(111)γ$ and $(110)α$ is predicted to be $5.25°$ using a hard sphere model of atoms.¹⁷ This is in good agreement with the average $5.73±3.24°$ angle between the $(111)γ$ and $(110)α$ lattice planes from these experiments (see Fig. 7). We considered the possibility that this angle could be a result of grain rotation due to plastic deformation. However, above the fcc $→$ bcc onset pressure, it is expected that the phase transformation is the dominant mechanism for relaxing the plastic strain.^53,54 Therefore, above the transition pressure, we do not expect dislocation mediated plastic flow to significantly contribute to the strain relaxation process. As a result, very little or no grain rotation is expected after the shock arrival and before the formation of the bcc phase. This has been observed in large scale molecular dynamics simulations of Ti and Zr^56–58 on ps timescales. Therefore, we conclude that the angle measured between $(111)γ$ and $(110)α$ is a result of the phase transition mechanism.

The Pitsch-distortion model for hard spheres also predicts a neutral or undistorted $[110]γ||[111]α$ direction [green atoms with black outline in Fig. 6(a)] such that $2aγ=3aα$⁠, where $a$ denotes the lattice parameters.¹⁷ In our experiments, this relationship was found to be accurate to within $∼2%$ for the shot where both the compressed $γ$ (fcc) and the $α$ (bcc) phases were simultaneously observed (185(20) GPa, shot ID s19C140). The small deviation from exact correspondence can be attributed to the slightly different pressure in the compressed $γ$ and the transformed $α$ phase or due to inaccuracies within the hard sphere model.

While our data are consistent with the Pitsch-distortion OR, we note that in shock studies on the $α−ω$ transformation in Ti,⁵⁷ there is evidence of different ORs being activated as a function of uniaxial shock loading along specific crystallographic axis. While no similar studies are reported for Cu, a similar OR-dependence on the shock loading direction may be possible. Further studies on highly oriented Cu samples are needed to clarify this point.

Experimentally, the bcc phase is first observed along the shock Hugoniot at 185 GPa (Fig. 3). This corresponds to an expected temperature of $∼3700$ K based on the Sesame equation-of-state (EOS) no. 3225 table⁵⁹ [Fig. 8(a)]. Along a cooler ramp-compression path, the fcc phase is observed through in situ XRD techniques to be stable to $>1150$ GPa⁸. While these datasets confirm that the bcc phase is stabilized at high temperatures, the extent of this $P$-$T$ stability remains unknown. Recently, the region of bcc $P$-$T$ stability has been explored theoretically.^3–6 The study of Sadigh et al.⁶ reports an fcc $→$ hcp $→$ bcc phase sequence along the Hugoniot,⁶ where stability of the hexagonal-close-packed (hcp) phase is predicted for very slow rates of compression. In the first-principles studies of Smirnov⁴ and Baty et al.,⁵ the high-temperature bcc phase is predicted to persist at all elevated pressures. In particular, Smirnov reported calculations for pressures as high as 1200 GPa [red dashed curve in Fig. 8(a)].

We use thermodynamic integration (TI) and melting line simulations to determine the extent of high P-T bcc phase stability. During the molecular dynamics, the forces are calculated with an Embedded Atom Model (EAM) interaction potential. As in the studies of Sadigh et al.⁶ and Neogi and Mitra,³ we use the EAM developed by Mishin et al.⁶⁰ The computational cost of EAM simulations allows us to explore various numerical convergence parameters and, in particular, the system size and the length of the simulations required to resolve the difference between the fcc and bcc free energies at high pressure and temperatures. The Helmholtz and Gibbs free energies were calculated for a set of pressure between 180 and 300 GPa and temperatures between 2000 and 5000 K. A system size of 32 000 atoms is used to study the bcc-fcc phase boundary. For additional details, see the supplementary material.

The EAM phase boundary between the fcc and bcc-structured copper is plotted as the bold blue curve in Fig. 8(a). We calculate the phase boundary with TI to 300 GPa and $∼5000$ K. Given the increasing slope of the fcc–bcc phase boundary between 200 and 300 GPa and the possibility that it might cross the melting line, we used a set of two-phase (solid–liquid) isobaric–isoenthalpic (NPH) coexistence simulations to establish the bcc-liquid and fcc-liquid phase boundaries for 300, 400, and 500 GPa. We used 216 000 and 160 000 atoms for the fcc-liquid and bcc-liquid simulations, respectively. For these three pressures, the fcc melting temperature was found to be lower than the bcc melting temperature (see Table I); hence, the bcc phase is predicted to persist to pressures greater than 500 GPa for EAM copper. Our boundary calculations are plotted with results from Baty et al. and Smirnov in Fig. 8(a).^4,5

By using the Rankine–Hugoniot equation, based on conservation of energy, $E=12P((ρ−ρ0)/ρρ0)$⁠,⁶² where P is the experimentally determined shock pressure, $ρ$ is the density from x-ray diffraction measurements, and $ρ0$ is the initial density (8.96 g/cm$3$⁠), we can compare our measured Hugoniot data with the calculated Hugoniot and melting curve and provide a useful description of the Cu phase diagram in pressure-energy space [Fig. 8(b)].

The computational modeling technique generally affects the determination of the phase boundary. EAM copper is a model potential. As such, whether it exhibits a phase diagram relevant to real copper is not a given. In fact Sadigh et al.⁶ show EAM copper favoring a hcp phase at high pressure.⁶ Such a phase is absent in real copper and in DFT-based calculations of copper. Neogi and Mitra³ quantified the crystal structures obtained upon shock compression also using the EAM model. Their simulations predicted a transformation from the fcc phase to an intermediate bct phase and to the bcc phase beginning at 80 GPa.³ We often use EAM and other model potentials to explore certain aspects of equilibrium and non-equilibrium systems. However, these models cannot be taken as predictive. This is in contrast to DFT-based calculations such as QMD. DFT is a first-principles (or ab initio) theory. It is an approximate formulation of quantum mechanics and only takes as input the type of atomic nuclei. For the construction of phase diagrams, one often starts with a quasi-harmonic DFT calculation. This is the approach taken by Smirnov.⁴ This consists of a set of 0-K DFT calculations of the total energy and of the phonon density of states at several densities. As such, this approach lacks any temperature-induced anharmonicity, which makes its use for the evaluation of phase transition lines a questionable endeavor at high temperatures. This temperature-induced anharmonicity could be captured by performing DFT-MD (aka QMD) calculations. However, as shown in the supplementary material, direct thermodynamic integration within DFT for the free-energy differences between the fcc and the bcc phases at these high-pressure and temperature conditions is computationally difficult. An attempt to take into account the effect of anharmonicity has been made by Baty et al.,⁵ who performed calculations in 512-atom supercells using the inverse Z method. This method relies on dynamics of solidification to be dominated by free-energy difference between the solid and the liquid phases. As a result, it is expected that at any given condition, the solid phase with the lowest free energy to be the phase that first solidifies. It is, however, well-known and is discussed at length by Sadigh et al.⁶ for solidification of Cu, that metastable bcc phase can be stabilized by kinetics of solidification due to lower solid–liquid interfacial free energy of this phase with the melt compared to that of the fcc phase with the melt. Hence, even this method suffers from uncontrolled errors. Nevertheless, we find near-quantitative agreement between the ab initio phase diagrams of Smirnov⁴ and Baty et al.⁵ Furthermore, the comparison of the results of the EAM copper calculations to the experimental results [Fig. 8(b)] shows that also EAM copper captures the energetics of actual copper for the fcc–bcc phase transition. Our experimental work does not support the existence of a bct³ or an hcp phase,⁶ and these may well be an artifact of the EAM potential.

The full azimuthal coverage and quality of the 2D x-ray diffraction images allows us to perform texture and structure analysis to better understand proposed deformation and phase transformation mechanisms. In copper at lower pressures (⁠$<60$ GPa), there is a change in the deformation mechanism from slip dominated to twinning dominated regimes.⁶³ However, this twinning dominated regime does not persist. We find little evidence of twinning based on examination of expected and observed angles between the lattice planes in ambient and shock-compressed copper. In particular, diffracted intensity at locations consistent with the ${111}⟨112¯⟩$ deformation twin, active in fcc metals,⁶⁴ was not observed. However, Schneider et al. suggested that twinning may require sufficient time that sub-microsecond experiments may not observe the twinning regime.⁶³ 

Sharma et al. proposed that the transition from fcc to bcc copper is facilitated through the generation of stacking faults,⁹ which was inferred by their observation of shifted (200) fcc XRD peaks. It was further proposed that the onset of the fcc $→$ bcc transition is strain-rate dependent, and a critical loading rate may be required for the production of stacking faults.⁹ This, it was argued, would explain the observed fcc stability at TPa conditions along a ramp-compression path.⁸ In our data, there was also a slight shift in the fcc (200) peak, consistent with stacking fault generation (see Fig. S9 in the supplementary material). We estimate a maximum stacking fault abundance, in the compressed fcc phase, of 6% (compared to 10% from Sharma et al.⁹). While stacking faults are present, there is no evidence that the transformation into the bcc phase requires stacking faults. Indeed, our computational work suggests that stacking fault generation is not required to mediate the transition to bcc copper; the bcc phase is the thermodynamically stable phase at high-$T$ and plastic deformation mechanisms may not be necessary for transformation. Additional high-pressure high-temperature experiments under different loading conditions could clarify this issue. We note also that the evidence for the Pitsch-distortion fcc $→$ bcc pathway (Fig. 6), which is not mediated by stacking faults, further suggests that their presence is not required for the transformation.

Additionally, our computational work in combination with the findings of Smirnov⁴ and Baty et al.⁵ suggest that the stability of the fcc phase in ramp compression experiments to much higher pressures than in shock experiments is rather due to thermodynamics and is not a kinetic effect. The theoretical calculations clearly show that the fcc phase is thermodynamically stable at low temperatures up to very high pressures. It is the higher temperatures accessed by the Hugoniot that allows for the bcc phase to become thermodynamically stable in the shock experiments.

We present a combined experimental and theoretical approach to determine the deformation processes associated with shock compression along the Cu Hugoniot and the stability of the high-temperature bcc phase.

The transformation mechanism between the fcc and bcc phases has been studied theoretically for the last hundred years with several candidate mechanisms proposed. Our study is the first experimental analysis of this transformation pathway, achieved using shock compression of Cu foils paired with synchrotron X-radiation. The diffraction texture analysis we present is consistent with this transformation occurring through the Pitsch-distortion mechanism. This is the likely transition mechanism for the other Group XI transition metals (Au and Ag), which are also observed to undergo an fcc–bcc transformation along the shock Hugoniot.^33,65,66 We use molecular dynamic calculations with an embedded atom model approach to determine the fcc–bcc phase boundary and the melt line in pressure–temperature energy space. Our calculations confirm that the bcc phase is stabilized by temperature and remains stable at high pressure to at least 500 GPa.

See the supplementary material for details of our thermodynamic integration and melting line simulations as well as detailed structural analysis of the fcc and bcc copper unit cells.

We thank Xiaoming Wang, Nicholas Sinclair, Adam Schuman, Paulo Rigg, Yoshi Toyoda, and the staff of the Dynamic Compression Sector for assistance with laser experiments and VISAR measurements. We are grateful for Carol A. Davis’ help with preparing the copper targets. M. Sims and J. Wicks were supported by the National Nuclear Security Administration for these experiments (Award No. DE-NA0003902). M. Sims was also supported by NSF Award No. 1952923 on this project. This publication is based upon work performed at the Dynamic Compression Sector, which is operated by Washington State University under the U.S. Department of Energy (DOE)/National Nuclear Security Administration Award No. DE-NA0003957. This research used resources of the Advanced Photon Source, a DOE Office of Science User Facility operated for the DOE Office of Science by Argonne National Laboratory under Contract No. DE-AC02-06CH11357. This work was performed under the auspices of the U.S. Department of Energy (DOE) by Lawrence Livermore National Laboratory under LLNL’s Laboratory Directed Research and Development (LDRD) Program under Grant Nos. 18-ERD-001, 18-ERD-012, and 21-ERD-032. Lawrence Livermore National Laboratory is operated by Lawrence Livermore National Security, LLC, for the U.S. Department of Energy, National Nuclear Security Administration under Contract No. DE-AC52-07NA27344.

This paper was prepared as an account of work sponsored by an agency of the United States Government. Neither the United States Government nor any agency thereof, nor any of their employees, makes any warranty, express or implied, or assumes any legal liability or responsibility for the accuracy, completeness, or usefulness of any information, apparatus, product, or process disclosed, or represents that is use would not infringe privately owned rights. Reference herein to any specific commercial product, process, or service by trade name, trademark, manufacturer, or otherwise does not necessarily constitute or imply its endorsement, recommendation, or favoring by the United States Government or any agency thereof. The views and opinions of the authors expressed herein do not necessarily state or reflect those of the United States Government or any agency thereof.

The authors have no conflicts to disclose.

Melissa Sims: Data curation (equal); Formal analysis (equal); Supervision (equal); Writing – original draft (equal). Richard Briggs: Formal analysis (equal); Writing – original draft (equal); Writing – review and editing (equal). Travis J. Volz: Formal analysis (equal); Writing – original draft (equal); Writing – review and editing (equal). Saransh Singh: Funding acquisition (equal); Writing – original draft (equal); Writing – review and editing (equal). Sebastien Hamel: Formal analysis (equal); Writing – original draft (equal); Writing – review and editing (equal). Amy L. Coleman: Formal analysis (equal). Federica Coppari: Formal analysis (equal); Writing – review and editing (equal). David J. Erskine: Formal analysis (equal). Martin G. Gorman: Formal analysis (equal). Babak Sadigh: Formal analysis (equal); Writing – original draft (equal); Writing – review and editing (equal). Jon Belof: Formal analysis (equal). Jon H. Eggert: Writing – original draft (equal); Writing – review and editing (equal). Raymond F. Smith: Formal analysis (equal); Funding acquisition (equal); Resources (equal); Writing – original draft (equal); Writing – review and editing (equal). June K. Wicks: Formal analysis (equal); Funding acquisition (equal); Writing – original draft (equal); Writing – review and editing (equal).

The data that support the findings of this study are available from the corresponding author upon reasonable request.