We present first-principles density functional theory calculations to study the *α*-*ω* phase transformation in Ti and Zr and its coupling to slip modes of the two phases. We first investigate the relative energetics of all possible slip systems in the *α* and *ω* phases to predict the dominant slip system that is activated during a plastic deformation under an arbitrary load. Using this and the crystallographic orientation relationships between *α* and *ω* phases, we construct low energy *α*/*ω* interfaces and study the energetics of the slip system at the interface between *α* and *ω* to compare to the slip systems in the bulk phases. We find that for a particular crystallographic orientation relationship, where $(basal)\alpha \u2225(prismatic-II)\omega $, and $[a]\alpha \u2225[c]\omega $, the slip at the interface is preferred compared to its bulk counterparts. This implies that the plastically deformed *α*/*ω* phase with this orientation relationship prefers to retain the interface (or coexisting phases) than transforming back to the pure phase after unloading. This is consistent with the observation that the *ω*-phase is retained in samples loaded in flyer plate experiments or under high-pressure torsion. Furthermore, calculation of the energy barrier for *α* to *ω* phase transformation as a function of glide at the *α*/*ω* interface shows significant coupling between the *α*-*ω* phase transformation and slip modes in Ti and Zr.

## I. INTRODUCTION

Titanium (Ti) has an exceptionally high strength to weight ratio, and it is a commonly used metal for many structural applications in aerospace, military, automotive, and several other industries. Zirconium (Zr), which belongs to the same group as Ti, has excellent resistance to corrosion and is extensively used in aggressive environments, such as those found in chemical and nuclear reactors.^{1–3} Consequently, the structural behavior of Ti and Zr as a function of temperature and pressure, as well as its mechanical properties under extreme conditions,^{4,5} has been the subject of much interest for several decades. At ambient pressure, the stable phase of Ti and Zr has a hexagonal close-packed (HCP) structure (*α*-phase); however, with pressure (in the range 2–10 GPa), it undergoes a structural phase transformation to a simple hexagonal (hex) structure (*ω*-phase).^{2,6,7} The *α*-*ω* phase transformation can occur under pressure generated either by static deformation,^{7–9} shock loading,^{10–13} or high-pressure torsion (HPT).^{14} The *α* phase is often not completely recovered after unloading at ambient temperature and pressure from the loaded and pressurized state, and the amount of the retained or metastable *ω* phase depends on the loading conditions. For example, samples obtained after shock loading retain nearly 80% *ω* phase,^{15} whereas under high pressure torsion, which involves shearing under high pressure, nearly 100%*ω* phase is retained.^{14}

Although the mechanical properties of the *α* phase for both Ti and Zr have been well studied, the *ω*-phase has received little study^{6,12,16–19} and only recently have atomistic simulations and density functional theory (DFT) been utilized to theoretically investigate the *ω*-phase.^{20,21} For example, to study the plastic behavior of the *ω*-phase, the dominant deformation or slip modes need to be characterized, and to understand the physics behind the retention and metastability of the *ω* phase on unloading, we need to understand how deformation modes couple with the energetics of the phase transformation process. Atomistic simulations have played an important role in understanding the properties of metals at small length scales. Molecular dynamics (MD) simulations are an efficient and useful technique to study such properties; however, we do not have very reliable interatomic potentials, particularly for Zr. On other hand, first-principles density functional theory (DFT) calculations, while computationally expensive, have been shown to be very accurate and reliable for studying mechanical properties, phase transformations, and energetics of defects and dislocations in metals and alloys.^{22–24} Moreover, information from first-principles calculations can be further incorporated into larger scale simulations or phenomenological models.^{19,25–30}

The energetics of the $\alpha -\omega $ phase transformation in Ti and Zr along different pathways has been studied using first-principles calculations, and it has been shown that the energy barrier for phase transformation varies for different pathways and decreases with pressure.^{31–34} Trinkle *et al.*, using a systematic symmetry analysis of transition pathways between *α* and *ω* phases for Ti, have shown that the TAO-1 pathway (with orientation relationship $(0001)\alpha \u2009\u2225(01\xaf11)\omega ;[112\xaf0]\alpha \u2009\u2225\u2009[011\xaf1]\omega )$ has the lowest energy barrier.^{31} Further, the pathways TAO-2 and Silcock with the same orientation relationship $((0001)\alpha \u2009\u2225\u2009(112\xaf0)\omega ;[112\xaf0]\alpha \u2009\u2225\u2009[0001]\omega )$ are energetically similar. Recently, several studies^{32,34} have examined the energetics of the phase transformation in Zr for TAO-1 and Silcock pathways and determined that the Silcock pathway has a larger energy barrier compared to TAO-1, similar to the findings in Ti.^{34}

We study in this work the Generalized Stacking Fault Energy (GSFE) surfaces for all possible slip systems in the *α* and *ω* phases of Ti and Zr using DFT calculations. Based on the relative energetics and ideal shear strength (ISS) obtained from the calculated GSFE curves, we show that the favored crystallographic slip modes in the *α*-phase are prismatic $\u3008a\u3009$, basal $\u3008101\xaf0\u3009$, basal $\u3008a\u3009$, and pyramidal-I $\u3008c+a\u3009$, whereas the favored crystallographic slip modes in the *ω*-phase are prismatic $\u3008c\u3009$, prismatic-II $\u3008101\xaf0\u3009$, and pyramidal-II $\u3008c+a\u3009$. These modes in *α* and *ω*-phases can accommodate any arbitrary plastic deformation state applied to the crystal. Comparison of the relative energetics of the slip systems in the two phases of Ti and Zr shows that slip in the *α* phase is much easier compared to slip in the *ω*-phase. This is consistent with experimental observations that the *ω*-phase is more brittle compared to *α*. We also investigate the GSFE for slip systems at the $\alpha /\omega $ interface and find that for a particular crystallographic orientation relationship, the slip at the interface is energetically favored compared to the corresponding slip systems in the bulk. The observation of an easy slip system at the interface provides an explanation of why the *ω*-phase is retained after unloading. We also performed calculation of energy barrier for *α* to *ω* phase transformation as a function of glide at $\alpha /\omega $ interface using the nudge elastic band (NEB) method and find large coupling between *α*-*ω* phase transformation and deformation modes in Ti and Zr. This change in the energy barrier with the deformation would change the phase transformation pressure for the plastically deformed samples compared to the one transformed with the hydrostatic pressure.

## II. METHOD

We performed first-principles DFT calculations using the generalized gradient approximation (GGA) for the exchange correlation functional with the Perdew-Becke-Erzenhof (PBE) parametrization^{35} as implemented in the VASP code.^{36,37} The interaction between the valence electrons and ionic cores is treated using projector augmented-wave (PAW) potentials.^{38,39} The number of valence electrons in the PAW potentials for Ti and Zr is four. In our DFT calculations, we used a plane wave energy cutoff of 400 eV and optimized the atomic structure until the force on each atom is smaller than 0.01 eV/$\xc5$. We used 19 × 19 × 11 and 19 × 19 × 25 Γ-centered Monkhorst-Pack^{40} k-point meshes to integrate the Brillouin Zone of the primitive unit cells of *α* and *ω* phases, respectively, to calculate the lattice and elastic constants for Ti and Zr. For calculation of the elastic constants, we used the finite difference distortion of the lattice as implemented in VASP (ISIF = 6) with a distortion step of 1.5% strain to calculate the elastic constants for both the systems.

The stable structure of Ti and Zr at ambient condition is hexagonal close packed (*α*-phase). However, using PAW PBE potential in DFT calculations, we find that energetically the *ω*-phase is more stable than the *α* phase for both Ti and Zr at ambient pressure. Therefore, to make the relative energies consistent with experiment for the *α* and *ω* phases in Ti and Zr (which are d-orbital electron metals), we use an on-site Hubbard parameter U.^{43,44} The use of the U parameter takes into account the electronic correlation in these metals. Figure 1 shows the relative energy of the *α* phase with respect to the *ω* phase as a function of pressure for different values of U. We find that the phase transformation pressure for *α* to *ω* increases with increase in U value (see Fig. 1). Based on the experimental phase transformation pressure for *α* to *ω* in Ti and Zr, which is in the range 2–10 GPa,^{2,6} we decided to use a U value of 2.2 eV and 1.2 eV for Ti and Zr, respectively. The calculated values of the lattice and elastic constants using DFT + U method for the *α* and *ω* phases of Ti and Zr are given in Table I. Comparison shows that the calculated values of the lattice constants *a* and *c* for both the *α* and *ω* phases of Ti and Zr are in very good agreement with experimentally measured values. Calculated elastic constants for the *α* phase are in good agreement with the experimental values. Experimental values of the elastic constants for the *ω* phase are, however, not available. All subsequent calculations were performed using DFT + U^{43,44} with the chosen values of U for Ti and Zr.

. | Ti . | Zr . | ||
---|---|---|---|---|

Parameter . | α
. | ω
. | α
. | ω
. |

a | 2.975 [2.951 (Ref. 41)] | 4.650 | 3.283 [3.233 (Ref. 1)] | 5.111 [5.036 (Ref. 9)] |

c | 4.729 [4.684 (Ref. 41)] | 2.865 | 5.228 [5.146 (Ref. 1)] | 3.190 [3.109 (Ref. 9)] |

C_{11} | 173.9 [160.0 (Ref. 42)] | 193.0 | 144.2 [144.0 (Ref. 42)] | 161.8 |

C_{33} | 195.4 [181.0 (Ref. 42)] | 245.3 | 166.2 [166.0 (Ref. 42)] | 197.3 |

C_{12} | 93.8 [90.0 (Ref. 42)] | 81.7 | 72.2 [74.0 (Ref. 42)] | 69.6 |

C_{13} | 71.0 [66.0 (Ref. 42)] | 52.0 | 62.8 [67.0 (Ref. 42)] | 46.3 |

C_{44} | 45.0 [46.5 (Ref. 42)] | 58.0 | 33.8 [33.4 (Ref. 42)] | 39.8 |

C_{66} | 40.0 [35.0 (Ref. 42)] | 55.7 | 36.0 [35.0 (Ref. 42)] | 46.1 |

Vol./atom | 18.12 | 17.88 | 24.395 | 24.053 |

E/atom | −5.829 | −5.825 | −7.471 | −7.464 |

. | Ti . | Zr . | ||
---|---|---|---|---|

Parameter . | α
. | ω
. | α
. | ω
. |

a | 2.975 [2.951 (Ref. 41)] | 4.650 | 3.283 [3.233 (Ref. 1)] | 5.111 [5.036 (Ref. 9)] |

c | 4.729 [4.684 (Ref. 41)] | 2.865 | 5.228 [5.146 (Ref. 1)] | 3.190 [3.109 (Ref. 9)] |

C_{11} | 173.9 [160.0 (Ref. 42)] | 193.0 | 144.2 [144.0 (Ref. 42)] | 161.8 |

C_{33} | 195.4 [181.0 (Ref. 42)] | 245.3 | 166.2 [166.0 (Ref. 42)] | 197.3 |

C_{12} | 93.8 [90.0 (Ref. 42)] | 81.7 | 72.2 [74.0 (Ref. 42)] | 69.6 |

C_{13} | 71.0 [66.0 (Ref. 42)] | 52.0 | 62.8 [67.0 (Ref. 42)] | 46.3 |

C_{44} | 45.0 [46.5 (Ref. 42)] | 58.0 | 33.8 [33.4 (Ref. 42)] | 39.8 |

C_{66} | 40.0 [35.0 (Ref. 42)] | 55.7 | 36.0 [35.0 (Ref. 42)] | 46.1 |

Vol./atom | 18.12 | 17.88 | 24.395 | 24.053 |

E/atom | −5.829 | −5.825 | −7.471 | −7.464 |

## III. RESULTS AND DISCUSSION

### A. Energy barrier for *α* to *ω* phase transformation in Ti and Zr

We investigated the energy barrier for the homogeneous phase transformation from the *α* phase to *ω* phase in both Ti and Zr. To calculate the energy barrier for *α* and *ω* phase transformation, we used fully relaxed structures of *α* and *ω* phases and performed a variable cell nudged elastic band (VCNEB) calculation using the DFT + U method starting from the linearly interpolated atomic configurations between the two phases. In VCNEB, both the atomic positions and the supercell dimensions are allowed to relax for each configuration, which allows pressure to be constant along the transition path. As discussed above, the two most common experimentally observed orientation relationships between the two phases are OR-I, $(0001)\alpha \u2225(112\xaf0)\omega ;[112\xaf0]\alpha \u2225[0001]\omega $ (i.e., $(Basal)\alpha \u2225(Prismatic-II)\omega $; and $[a]\alpha \u2225[c]\omega )$, and OR-II, $(0001)\alpha \u2225(01\xaf11)\omega ;[112\xaf0]\alpha \u2225[011\xaf1]\omega $ (i.e., $(Basal)\alpha \u2225(Pyramidal-I)\omega $; and $[a]\alpha \u2225[c+a]\omega )$.^{31–33} The OR-I and OR-II are also known as Silcock and TAO-1 orientation relationships. To calculate the energy barrier along the Silcock and TAO-1 pathways given by the two crystal orientation relationships as discussed in previously,^{31–33} we took two supercells consisting of 12 atoms and 6 atoms. For two supercells, there are many ways to shuffle the atoms, and in this work we used the same atomic shuffles as given by Ghosh *et al.*^{32} During the NEB calculation, we allow all the atomic degrees of freedom to relax along the pathway to calculate the energy barrier.

The calculated energy barrier for the homogeneous phase transformation along Silcock and TAO-1 pathways for Ti and Zr is given in Fig. 2 for two values of pressure (P = 0 and P = 10 GPa). The enthalpy at x = 0 corresponds to *α* phase and x = 1.0 corresponds to *ω* phase. The relative enthalpy shows that at P = 0, *α* phase is stable compared to *ω* phase, and with increasing pressure (at 10 GPa), *ω* phase becomes more stable compared to *α* phase, consistent with earlier studies.^{31–34} The energy barrier for the Silcock pathway decreases with pressure for both Ti and Zr, whereas the energy barrier for TAO-1 pathway does not show much sensitivity to pressure. In order to understand the coupling between plastic deformation modes in the *α* and *ω* phases with the phase transformation (i.e., energy barrier), in Secs. III B and III C we study the energetics of all possible deformation modes in the *α* and *ω* phases of Ti and Zr.

### B. Possible slip systems in the *α*-phase of Zr and Ti

To accommodate deformation in the $\u3008a\u3009$ and $\u3008c\u3009$ directions for the low-symmetry hcp crystal structure of the *α* phase, we have taken all possible slip planes.^{45} Their orientation and directions with respect to the hexagonal unit cell are shown in Fig. 3(b). These nine slip modes span five distinct planes, and four of these planes are associated with two independent slip modes. The number of slip systems per mode is determined by considering all crystal symmetries. The group of slip modes accommodating $\u3008a\u3009$ deformation are denoted as basal $\u3008a\u3009$, basal $\u3008101\xaf0\u3009$, prismatic $\u3008a\u3009$, prismatic-II $\u3008101\xaf0\u3009$, and pyramidal-I $\u3008a\u3009$, and the group accommodating $\u3008c\u3009$ deformation are referred to as prismatic $\u3008c\u3009$, prismatic-II $\u3008c\u3009$, pyramidal-I $\u3008c+a\u3009$, and pyramidal-II $\u3008c+a\u3009$.

The GSFE curve is the change in energy associated with the rigid displacement of one half of the crystal relative to the other half across a crystallographic plane along a crystallographic direction. A comparison of the GSFE curve for different slip systems can give information about dominant slip systems that will be activated under arbitrary applied load. Moreover, for some dislocations with planar dissociated states, it is possible to identify from the GSFE surface a set of partial dislocations that can produce the same displacement as a full dislocation. A local minimum, if exists, on the GSFE curve corresponds to the stable stacking fault energy (SFE) and maximum values correspond to unstable stacking fault energies.

In low symmetry crystal structures such as hcp and simple hexagonal, glide planes among the various slip modes are topologically different, some can be rumpled, and the atomic positions about the glide direction can be asymmetric.^{46} In such cases, additional mechanisms, like local atomic shuffling and glide non-parallel to the theoretical glide direction, are required to carry the shear. In the standard method of calculation of the GSFE, we divide the crystal into two halves and shift the upper half of the crystal with respect to the lower half of the crystal along the glide direction. At each displacement, we minimize the energy of the system by fixing all the positions of atoms in both the upper and lower crystals in the x and y directions and allow the positions in the z direction to relax.^{47–49} However, it has been shown that for low symmetric structures like hexagonal the standard approach sometimes does not give the local energy minima. Therefore in this work for the calculation of GSFE curves, we allowed additional relaxation of all atomic positions along the direction lying normal to the glide direction.^{46,50–52}

For the calculations of GSFE curves using the DFT + U method, we used slabs that are periodic in x and y directions. The periodic dimension along the z direction contains a 15 $\xc5$-thick vacuum layer. The values of *a* and *c* obtained from bulk *α*-phase of Ti and Zr (Table I) are used to construct the supercells for GSFE calculation. Since the dimensions of the unit cell depend on the slip plane, the appropriate number of atoms and the dimensions of the supercells is different for different slip planes. For each plane, the chosen supercell dimensions corresponded to the minimum number of layers for which convergence in energy of the system is attained.

Figures 3(c) and 3(d) show the calculated GSFE curves for all nine slip systems in the *α* phase of Ti and Zr using DFT + U method. We find that prismatic $\u3008a\u3009$, basal $\u3008101\xaf0\u3009$, basal $\u3008a\u3009$, pyramidal-I $\u3008c+a\u3009$, pyramidal-II $\u3008c+a\u3009$, and pyramidal-I $\u3008a\u3009$ are energetically favorable slip systems in the *α* phase for both the Ti and Zr, whereas prismatic-II $\u3008101\xaf0\u3009$, prismatic-II $\u3008c\u3009$, and prismatic $\u3008c\u3009$ slip systems are energetically very high and would be unfavorable. From these GSFE curves, we also estimate the ideal shear strength (ISS), which is a measure of the lattice resistance for gliding one part of the crystal with respect to other, by taking derivative of the GSFE curves with respect to displacement as $1b\u2202\Gamma \u2202u$ and choosing the value of lower of the two peaks from the derivative obtained from either shifting from left to right or in reverse from right to left on the GSFE curves,^{53} and list this in Table II. The knowledge of ideal shear strength (ISS) provides a measure for the energy barrier for shear and to identify the relative ease for the activation of different slip modes. The more likely modes on which dislocations would glide in an actual metal would tend to possess the lower energy barriers for shear or ideal shear stresses. This analysis, however, is a relative one in the sense that the actual barriers that would resist the motion of a dislocation along these pathways at deformation temperatures and pressures are typically orders of magnitude lower than the ideal shear stress at 0 K. To determine relative favorability of slip systems, we first rank these modes based on the ISS value from the lowest (most favorable) to the highest in Table II. We also provide the peak energy value from the GSFE curves corresponding to the lowest first peak. Based on a comparison of the ideal shear stress, the easier slip modes for accommodating $\u3008a\u3009$ deformation would be prismatic $\u3008a\u3009$ and basal $\u3008101\xaf0\u3009$ and for accommodating $\u3008c\u3009$ deformation those would be pyramidal-I $\u3008c+a\u3009$ and pyramidal-II $\u3008c+a\u3009$.

Index . | Slip . | $\Gamma usf(Ti)$ . | ISS (Ti) . | $\Gamma usf(Zr)$ . | ISS (Zr) . |
---|---|---|---|---|---|

1 | Prismatic $\u3008a\u3009$ | 282.9 | 3.86 | 249.1 | 3.02 |

2 | Basal $\u3008101\xaf0\u3009$ | 276.4 | 3.97 | 239.9 | 3.07 |

3 | Basal $\u3008a\u3009$ | 281.8 | 4.27 | 244.2 | 3.18 |

4 | Pyramidal-I $\u3008c+a\u3009$ | 539.0 | 5.46 | 512.7 | 4.60 |

5 | Pyramidal-II $\u3008c+a\u3009$ | 642.9 | 6.42 | 590.6 | 5.29 |

6 | Pyramidal-I $\u3008a\u3009$ | 729.8 | 7.82 | 620.5 | 5.97 |

7 | Prismatic $\u3008c\u3009$ | 1501.5 | 11.66 | 1334.7 | 8.81 |

8 | Prismatic-II $\u3008c\u3009$ | 1987.2 | 13.38 | 1688.3 | 9.61 |

9 | Prismatic-II $\u3008101\xaf0\u3009$ | 2840.4 | 18.46 | 2440.3 | 13.98 |

Index . | Slip . | $\Gamma usf(Ti)$ . | ISS (Ti) . | $\Gamma usf(Zr)$ . | ISS (Zr) . |
---|---|---|---|---|---|

1 | Prismatic $\u3008a\u3009$ | 282.9 | 3.86 | 249.1 | 3.02 |

2 | Basal $\u3008101\xaf0\u3009$ | 276.4 | 3.97 | 239.9 | 3.07 |

3 | Basal $\u3008a\u3009$ | 281.8 | 4.27 | 244.2 | 3.18 |

4 | Pyramidal-I $\u3008c+a\u3009$ | 539.0 | 5.46 | 512.7 | 4.60 |

5 | Pyramidal-II $\u3008c+a\u3009$ | 642.9 | 6.42 | 590.6 | 5.29 |

6 | Pyramidal-I $\u3008a\u3009$ | 729.8 | 7.82 | 620.5 | 5.97 |

7 | Prismatic $\u3008c\u3009$ | 1501.5 | 11.66 | 1334.7 | 8.81 |

8 | Prismatic-II $\u3008c\u3009$ | 1987.2 | 13.38 | 1688.3 | 9.61 |

9 | Prismatic-II $\u3008101\xaf0\u3009$ | 2840.4 | 18.46 | 2440.3 | 13.98 |

### C. Possible slip systems in the *ω*-phase of Zr and Ti

The crystallograhic slip planes in the *ω*-phase are the same as the *α*-phase. The main difference is that $c/a>1$ for the *α*-phase and $c/a<1$ for the *ω*-phase. We studied the GSFE curves for all nine slip systems that can accommodate deformation in the $\u3008a\u3009$ and $\u3008c\u3009$ directions in the *ω*-phase of Ti and Zr using the DFT + U method to identify the dominant slip modes. Similar to the *α*-phase, the five distinct crystallographic planes in *ω*-phase are shown in Fig. 4(a).

The calculated GSFE curves for all nine slip systems in the *ω* phase of Ti and Zr are shown in Figs. 4(b) and 4(c). In the *ω* phase of both Ti and Zr, we find that prismatic $\u3008c\u3009$, prismatic-II $\u3008c\u3009$, prismatic-II $\u3008101\xaf0\u3009$, pyramidal-II $\u3008c+a\u3009$, pyramidal-I $\u3008c+a\u3009$, and pyramidal-I $\u3008a\u3009$ are energetically favorable slip systems, whereas prismatic $\u3008a\u3009$, basal $\u3008101\xaf0\u3009$, and basal $\u3008a\u3009$ slip systems would be energetically unfavorable. Interestingly, the dominant slip modes in the *α* phase become unfavorable in the *ω* phase. From these GSFE curves in *ω* phase, we again estimated the ideal shear strength via $1b\u2202\gamma \u2202u$ taken from the lower of the two peaks obtained from the derivative of the GSFE curves^{53} (see Table III). To determine the favorability of slip modes, we first rank these modes based on the ISS value from the lowest (most favorable) to the highest in Table II. Based on a comparison of the ideal shear stress, the easier slip modes for accommodating $\u3008c\u3009$ deformation would be prismatic $\u3008c\u3009$ and prismatic-II $\u3008c\u3009$ and for accommodating $\u3008a\u3009$ deformation those would be prismatic-II $\u3008101\xaf0\u3009$, pyramidal-II $\u3008c+a\u3009$, pyramidal-I $\u3008c+a\u3009$, and pyramidal-I $\u3008a\u3009$.

Index . | Slip . | $\Gamma usf(Ti)$ . | ISS (Ti) . | $\Gamma usf(Zr)$ . | ISS (Zr) . |
---|---|---|---|---|---|

1 | Prismatic $\u3008c\u3009$ | 320.5 | 5.18 | 279.9 | 3.97 |

2 | Prismatic-II $\u3008c\u3009$ | 512.5 | 6.66 | 436.6 | 5.01 |

3 | Prismatic-II $\u3008101\xaf0\u3009$ | 796.3 | 7.29 | 710.1 | 6.19 |

4 | Pyramidal-II $\u3008c+a\u3009$ | 951.3 | 7.74 | 847.1 | 6.34 |

5 | Pyramidal-I $\u3008c+a\u3009$ | 901.1 | 7.88 | 827.3 | 6.54 |

6 | Pyramidal-I $\u3008a\u3009$ | 1060.6 | 10.16 | 803.9 | 7.57 |

7 | Basal $\u3008a\u3009$ | 2022.2 | 13.51 | 1703.5 | 10.21 |

8 | Basal $\u3008101\xaf0\u3009$ | 1968.2 | 14.63 | 1572.2 | 11.42 |

9 | Prismatic $\u3008a\u3009$ | 2700.8 | 19.52 | 2262.6 | 14.78 |

Index . | Slip . | $\Gamma usf(Ti)$ . | ISS (Ti) . | $\Gamma usf(Zr)$ . | ISS (Zr) . |
---|---|---|---|---|---|

1 | Prismatic $\u3008c\u3009$ | 320.5 | 5.18 | 279.9 | 3.97 |

2 | Prismatic-II $\u3008c\u3009$ | 512.5 | 6.66 | 436.6 | 5.01 |

3 | Prismatic-II $\u3008101\xaf0\u3009$ | 796.3 | 7.29 | 710.1 | 6.19 |

4 | Pyramidal-II $\u3008c+a\u3009$ | 951.3 | 7.74 | 847.1 | 6.34 |

5 | Pyramidal-I $\u3008c+a\u3009$ | 901.1 | 7.88 | 827.3 | 6.54 |

6 | Pyramidal-I $\u3008a\u3009$ | 1060.6 | 10.16 | 803.9 | 7.57 |

7 | Basal $\u3008a\u3009$ | 2022.2 | 13.51 | 1703.5 | 10.21 |

8 | Basal $\u3008101\xaf0\u3009$ | 1968.2 | 14.63 | 1572.2 | 11.42 |

9 | Prismatic $\u3008a\u3009$ | 2700.8 | 19.52 | 2262.6 | 14.78 |

## IV. POSSIBLE SLIP SYSTEMS IN THE *α*/*ω* INTERFACE OF Ti and Zr

To understand the coupling between slip modes and *α* to *ω* phase transformation in Ti and Zr, we took three bi-layer slabs as shown in Fig. 5 consisting of *α* and *ω* layers with three different crystallographic orientation relationships. The bilayer slabs are selected in such a way that the interface plane is parallel to the dominant slip planes of both *α* and *ω* phases. The dominant slip planes in the *α* and *ω* phases are prismatic; therefore, for the first bi-layer slab [Fig. 5(a)], we took the prismatic plane of the *α* layer parallel to the prismatic plane of the *ω* layer. In the second bilayer slab, we assumed that the basal plane of the *α* layer is parallel to the prismatic-II plane of the *ω* layer, which includes the second dominant slip planes of the *α* and *ω* phases and has the same crystallographic orientation relationship as the Silcock OR. In the third bilayer slab, we have the basal plane of *α* layer parallel to the pyramidal plane of the *ω* layer, which includes the second and third dominant slip planes of the *α* and *ω*, respectively, and it has the same crystallographic orientation as TAO-1. There is lattice mismatch for different orientation relationships between *α* and *ω* phases, and to make coherent bi-layer models for the three interfaces shown in Fig. 5, we first apply equal and opposite strains in both the *α* and *ω* layers. We then allow supercell for each bilayer model to fully relax (both atomic positions as well as supercell dimensions). The dimensions of the relaxed supercell models are used to calculate the in-plane strains for each layer compared to the bulk *α* and *ω* phases.

Next, we calculated the formation energy of these three interfaces using periodic supercells as shown in Fig. 5. Interface formation energy is given as

where $E\alpha /\omega super$ is the energy of supercell containing *n*_{1} atoms in *α* layer and *n*_{2} atoms in *ω* layer, $E\alpha bulk$ is the energy of *n*_{1} atoms in bulk *α*, $E\omega bulk$ is the energy of *n*_{2} atoms in bulk *ω*, $E\alpha elas$ is the elastic energy in *α* layer, $E\omega elas$ is the elastic energy in *ω* layer due to lattice mismatch and the formation of coherent interface, and A is the interface area.

The calculated formation energies for the three interfaces are shown in Table IV. The formation energy for interface 1 is $18.2\u2009mJ/m2$, interface 2 is $180.7\u2009mJ/m2$, and for interface 3 is relatively large $448.9\u2009mJ/m2$. Since the formation energy of interface 3 is energetically unfavorable compared to other two interfaces and does not have dominant slip systems parallel to the interface, next we only studied the GSFE curves for the interface 1 and interface 2.

Interface . | α-Strain (in %)
. | ω-Strain (in %)
. | No of atoms . | γ
. _{int} |
---|---|---|---|---|

1 | $\u03f5x=2.55;\u03f5y=0.7$ | $\u03f5x=\u22121.21$; $\u03f5y=\u22120.9$ | 72 | 18.2 |

2 | $\u03f5x=1.86;\u03f5y=2.18$ | $\u03f5x=\u22121.90$; $\u03f5y=\u22121.93$ | 120 | 180.7 |

3 | $\u03f5x=5.01;\u03f5y=4.8$ | $\u03f5x=\u22125.09$; $\u03f5y=\u22124.7$ | 96 | 448.9 |

Interface . | α-Strain (in %)
. | ω-Strain (in %)
. | No of atoms . | γ
. _{int} |
---|---|---|---|---|

1 | $\u03f5x=2.55;\u03f5y=0.7$ | $\u03f5x=\u22121.21$; $\u03f5y=\u22120.9$ | 72 | 18.2 |

2 | $\u03f5x=1.86;\u03f5y=2.18$ | $\u03f5x=\u22121.90$; $\u03f5y=\u22121.93$ | 120 | 180.7 |

3 | $\u03f5x=5.01;\u03f5y=4.8$ | $\u03f5x=\u22125.09$; $\u03f5y=\u22124.7$ | 96 | 448.9 |

To understand the slip activity of the interface plane compared to that of the bulk *α* and *ω* away from the interface in the bilayer system, we study the GSFE curves for the slip systems at the interface and compare them with the GSFE of slip systems away from the interface. As shown in Fig. 6(a), we calculated the GSFE for three slip planes parallel to the interface (i) in the *α* phase, (ii) at the $\alpha /\omega $ interface, and (iii) in the *ω* phase. For the interface 1, the *α* layer has crystallographic *a* direction, which is along the x axis of the supercell, and the *c* direction, which is along the y axis of the supercell, whereas the *ω* layer has crystallographic *c* direction which is along the x axis of the supercell and the *a* direction which is along the y axis of the supercell. Figures 6(b) and 6(c) show the GSFE curves for slip systems in Ti and Zr, respectively, for the interface 1. In the top panels, we show the GSFE curves for three planes with glide direction along the x axis of the supercell, and in the bottom panels we show the GSFE curves when the glide direction is along the y axis of the supercell. We find that it is relatively easy to glide in the *α* layer along the x axis compared to gliding at the interface in the *ω* layer. The GSFE curve for the interface plane lies between the GSFE curves for the *α* and *ω* layers. This is consistent with a direct comparison of bulk GSFE curves for *α* and *ω* phases, where prismatic $\u3008a\u3009$ of the *α* phase is the dominant slip mode in all slip modes of *α* and *ω* phases. Interestingly, all the results of the GSFE curves for Ti and Zr are qualitatively very similar but we present the results for both metals for completeness.

The GSFE curves along the y direction in supercell 1 show similar behavior as the GSFE along x direction except for the slip system at the interface for normalized displacement between 0.35 and 0.65. Analysis of the relaxed structure of these configurations shows some structural transformation in few layers near the interface in the *α* phase, which lowers the Interface GSFE curve compared to the GSFE curve for the glide in the *α* layer for both the Ti and Zr.

Similarly, we study the GSFE curves for three slip systems parallel to interface 2 (i.e., $(basal)\alpha \u2225\u2009(prismatic-II)\omega $). The calculated GSFE curves for slip systems in the *α* phase at the $\alpha /\omega $ interface and in the *ω* phase for interface 2 are shown in Fig. 7 for both Ti and Zr. In the supercell for interface 2, the crystallographic *a* direction of *α* and *c* directions of *ω* is along y axis, and crystallographic $[101\xaf0]$ direction of both *α* and *ω* is along x as shown in Fig. 7(a). We find that GSFE curves for slip system at interface for the interface 2 is energetically much favorable compared to GSFE in the *α* and *ω* phases. As shown in upper panels for the GSFE of Ti and Zr along the y direction, the GSFE curves for both the *α* and *ω* phases are higher in energy than the GSFE for interface. In addition to lower energy of the GSFE for the interface slip, the unstable stacking fault energy for partials vanishes, which implies that the full *a* or *c* Burgers vectors of *α* and *ω* phases, respectively, would split into two equal partials that can separate far apart from each other without any stacking fault between them.

Interestingly, the GSFE curve for the interface plane is lower than the GSFE curve for bulk *α* and *ω* phases for glide along both the x and y directions of the supercell and the fault energy vanishes for the partial Burgers vector at the interface. In the lower panels of Figs. 7(b) and 7(c), we show the GSFE for Ti and Zr for glide along x direction of the supercell. The periodic length in the *α* and *ω* crystal is different along the x direction of the supercell (i.e., $3\xd7[101\xaf0]\alpha =2\xd7[101\xaf0]\omega $); therefore, the normalized displacement on the horizontal axes are different for the *α* and *ω* GSFE curves. The GSFE curve for the interface is lower than the corresponding GSFE curves in the bulk and shows vanishing minimum at 0.17 normalized displacement.

Comparison of GSFE curves for slip planes at the interface 1 and interface 2 (Figs. 6 and 7) shows that the slip plane at interface 2 is the easiest compared to slip systems in both *α* and *ω* phases. The crystallographic orientation relationship for the interface 2 is the same as the Silcock OR, which has been observed during the phase transformation. If interface 2 slip is activated during loading of the mixed phase, the system would energetically prefer to maintain the mixed phase with $\alpha /\omega $ interface. Our finding of the dominant slip system in bi-layer at the interface 2 is consistent with the experimental observation and provides an explanation of why a certain fraction of the *ω* phase is retained after unloading the crystal.

## V. COUPLING BETWEEN PHASE TRANSFORMATION AND SLIP MODES IN THE *α* AND *ω* PHASES OF Ti AND Zr

To understand the coupling between phase transformation and slip mode, we study the change in the energy barrier for *α* to *ω* phase transformation as a function of glide at the interface. The preferable glide direction for the both interface 1 and interface 2 is along $\u3008a\u3009$ of the *α* phase (see Figs. 6 and 7). First, we calculated the energy barrier for moving the interface in the *α* phase by one periodic length along the z direction, which is $(3)a\alpha $ for the interface 1 and *c* for the interface 2. To calculate the energy barrier for moving the interface boundary, we perform the nudge elastic band calculations starting from initial and final configurations such that in the final configuration the interface boundary has moved by one periodic length in the *α* phase compared to the interface boundary in the initial configuration. The change in energy as function of moving the interface boundary by one periodic length along the z direction in the *α* phase for interface 1 and interface 2 is shown by ND = 0 in the Fig. 8 for Ti.

Next, we glide the *ω* part of the crystal along $\u3008a\u3009$ direction of the *α* crystal with respect to the *α* crystal at the interface for both initial and the final configurations by 0.25 of the normalized displacement. This corresponds to moving on the interface GSFE curve (shown by circle) in Figs. 6 and 7 to a point corresponding to a normalized displacement value of 0.25. At this value of the glide displacement, we again calculate the energy barrier for moving the interface boundary for one periodic length along the z direction in the *α* phase for interface 1 and interface 2 using the nudge elastic band method. The change in energy for moving the phase boundary at this value of glide for both the interfaces is shown by ND = 0.25 in Fig. 8. We find that energy barrier decreases for ND = 0.25 compared to ND = 0 for both the interfaces (I and II). We further glide the interface corresponding to the normalized displacement of 0.5 and repeat the nudge elastic band calculation to get the change in energy for moving the interface boundary in the *α* crystal, which is shown by ND = 0.50. We find that the energy barrier goes up compared to ND = 0.25. Comparison of the change in energy profile for three values of glide displacement shows that energy barrier first decreases as we glide from ND = 0 to ND = 0.25 and then increases as we glide ND = 0.25 to ND = 0.50 for both the interfaces. The energy barrier for *α* to *ω* phase transformation at interface 1 decreases by 243 mJ/m^{2} (28%) with the glide at interface [638 mJ/m^{2} for ND = 0.25 compared to 880 mJ/m^{2} for ND = 0 in Fig. 8(a)], whereas the energy barrier for *α* to *ω* phase transformation at interface 2 decreases by 177 mJ/m^{2} (60%) with the glide at interface [119 mJ/m^{2} for ND = 0.25 compared to 296 mJ/m^{2} for ND = 0 in Fig. 8(b)]. The change in energy barrier for *α* to *ω* phase transformation with glide shows that there is a significant amount of coupling between the slip mode and phase transformation.

As shown in Figs. 6 and 7, the energy barrier for the phase transformation is large compared to the coupling between deformation modes and phase transformation. To check if the glide of one part of crystal over the other can actually drive the *α* to *ω* phase transformation near the glide plane, we performed additional relaxation where we allowed atoms within one periodic length near the $\alpha /\omega $ phase boundary to relax. For example, the solid red curve with filled triangles in Fig. 9(a) is the GSFE curve for slip system in the *ω* layer of the interface 1. During the additional relaxation, we allowed the atoms between the slip plane in *ω* phase and the $\alpha /\omega $ phase boundary to fully relax during the glide to see if the *ω* phase atoms transform to the *α* phase atoms, which would move the phase boundary up in the *ω* phase. The relaxed energy of the system, shown by the dotted red curve with open triangles, shows that if glide is smaller than 0.3, the small dislocation created by glide in *ω* phase is relaxed by shearing the atoms in the relaxed layer. However, if the glide is larger than 0.3, then dislocation created by glide in the *ω* does not relax. However, during the additional relaxation, we do not see the phase transformation of the *ω* atoms to *α* phase. We found similar observation if we allowed the atoms between slip in the *α* phase and $\alpha /\omega $ phase boundary to relax during the glide at the interface as shown by dotted line with open circle in Fig. 9(a).

Interestingly, we found that the dislocation cross slips if we repeat the same relaxation for the second interface (basal $\u2225$ prismatic-II). The dotted red curve with open triangles in Fig. 9(b) is the energy of the system when we allowed the atoms between the slip plane in *ω* phase and the $\alpha /\omega $ phase boundary to fully relax during the glide in the *ω* phase. We find, when glide is small, that the additional relaxation of atoms between the slip plane in the *ω* phase and the $\alpha /\omega $ phase boundary allows cross slip of the dislocation glide plane from the *ω* phase to the $\alpha /\omega $ phase boundary. However, we do not find *ω* atoms transforming to *α* phase atoms, which is consistent with the observation that the energy barrier for the phase transformation is large compared to the coupling between deformation modes and the energy barrier between *α* to *ω* phases.

## VI. CONCLUSIONS

In summary, using first-principles density functional theory calculations, we have identified the preferred slip modes in *α* and *ω*-phases of Ti and Zr by calculating the GSFE curves for nine geometrically possible slip modes. From a knowledge of the dominant slip systems in the *α* and *ω* phases and crystallographic orientation relationships, we constructed low energy *α*/*ω* interfaces and studied the energetics of the slip systems at the interface compared to those in the bulk *α* and *ω* phases. We find that for interface 2 with $(basal)\alpha \u2225(prismatic-II)\omega $; and $[a]\alpha \u2225[c]\omega $ orientation relationship, which is one of the experimentally observed ORs, the slip system at this interface is preferred compared to the bulk counterparts. This implies that plastically the deformed *α*/*ω* coexisting phase with this orientation relationship would be preferred (or retained) compared to transforming back to the pure phase on unloading. This finding provides a possible explanation to the observation that the *ω* phase is retained under various loading conditions, such as high-pressure torsion (HPT). We also find evidence of significant coupling between the phase transformation and slip modes in the *α* and *ω* based on the change in the energy barrier associated with glide at the *α*/*ω* interface. Further analysis shows easy cross slip of the glide plane to the phase boundary in the course of glide in the *α* and *ω* phases.

## ACKNOWLEDGMENTS

We thank the DOE Advanced Science Computational (ASC) program at Los Alamos National Laboratory (LANL) for financial support and LANL Institutional Computing for computational resources.