Melt-growth dynamics in CdTe crystals

We use a new, quantum-mechanics-based bond-order potential (BOP) to reveal melt-growth dynamics and fine-scale defect formation mechanisms in CdTe crystals. Previous molecular dynamics simulations of semiconductors have shown qualitatively incorrect behavior due to the lack of an interatomic potential capable of predicting both crystalline growth and property trends of many transitional structures encountered during the melt $\rightarrow$ crystal transformation. Here we demonstrate successful molecular dynamics simulations of melt-growth in CdTe using a BOP that significantly improves over other potentials on property trends of different phases. Our simulations result in a detailed understanding of defect formation during the melt-growth process. Equally important, we show that the new BOP enables defect formation mechanisms to be studied at a scale level comparable to empirical molecular dynamics simulation methods with a fidelity level approaching quantum-mechanical methods

We use a new, quantum-mechanics-based bond-order potential (BOP) to reveal melt-growth dynamics and fine-scale defect formation mechanisms in CdTe crystals. Previous molecular dynamics simulations of semiconductors have shown qualitatively incorrect behavior due to the lack of an interatomic potential capable of predicting both crystalline growth and property trends of many transitional structures encountered during the melt → crystal transformation. Here we demonstrate successful molecular dynamics simulations of melt-growth in CdTe using a BOP that significantly improves over other potentials on property trends of different phases. Our simulations result in a detailed understanding of defect formation during the melt-growth process. Equally important, we show that the new BOP enables defect formation mechanisms to be studied at a scale level comparable to empirical molecular dynamics simulation methods with a fidelity level approaching quantum-mechanical methods.
CdTe-based materials have been instrumental for technological breakthroughs in the semiconductor industry. The largest science and technology impact has been their widespread use in solar cells, radiation detectors, and medical imaging devices. While CdTe solar cells currently have the lowest cost compared to any other photovoltaic technologies [1], the material is far from optimum as the record energy-conversion efficiencies achieved today are only 16%, significantly below the 29% theoretical value [2,3]. CdTe-based Cd 1−x Zn x Te (CZT) alloys are currently the leading semiconductors for γ-ray detection, but their application is limited by low manufacturing yield (and, therefore, high cost) of detector-grade materials. The under-achievement of CdTe solar cells and CZT detectors can both be attributed to charge-trapping defects formed during CdTe growth [4][5][6].
Direct molecular dynamics (MD) simulations provide a fundamental understanding of the CdTe growth dynamics and defect formation. However, such simulations are extremely challenging because they sample a large number of metastable configurations not known a priori. If the interatomic potential used in a simulation over-stabilizes a metastable configuration, that configuration will likely persist, leading to an amorphous growth. Previous MD simulations of semiconductor growth have shown qualitatively incorrect behavior due to the lack of an interatomic potential capable of predicting both crystalline growth and property trends of many transitional structures encountered during the growth. In addition, past MD simulations of semiconductor crystalline growth were limited to vapor deposition [7][8][9][10][11][12][13], while cases for melt-growth (which are likely to be more challenging because the crystallization occurs from a condensed liquid phase rather than the atom-by-atom assembly during vapor deposition) have yet to be demonstrated.
A vast majority of successful MD simulations of semiconductor vapor deposition [7][8][9][10][11] used Stillinger-Weber (SW) [14] potentials. We established previously [15] that while SW potentials can easily ensure the lowest energy for the equilibrium tetrahedral semiconductor crystal and its crystalline growth during vapor deposition, they cannot satisfactorily capture the property trends of other configurations. Hence, they cannot accurately reveal defect formation. Tersoff potentials [16], on the other hand, capture property trends more accurately. However, there is no obvious way to ensure the lowest energy for the tetrahedral structure. Because the tetrahedral structure must have a lower energy than any other structures, it is unclear which phases should be included in the potential parameterization. As a result, crystalline growth is difficult to achieve with Tersoff potentials unless the potential parameterization is done iteratively with crystalline growth used as a criterion to actively select and modify the phases to be included in the fitting. Not surprisingly, we found that many literature Tersoff potentials [17][18][19] predict amorphous growth during vapor deposition simulations.
The objective of this work is to fill the missing research areas identified above by not only demonstrating the crystalline growth of semiconductor compounds from melt, but also advancing beyond Tersoff potentials on capturing properties of different phases using a new bond-order potential (BOP) [20][21][22][23]. This work will begin to allow detailed investigations of defect formation mechanisms during melt-or vapor-phase growth. While we focus on CdTe compounds, the methods are applicable to a broad range of semiconductors.
Unlike the Stillinger-Weber [14] and Tersoff/Brenner [16,24] potentials commonly used for semiconductors, BOP is analytically derived from quantum mechanical theories [20][21][22][23]. In particular, It incorporates primary (σ) and secondary (π) bonding and the valencedependence of the heteroatom interactions, with the functional forms of the potential derived from tightbinding theory under the condition that the first two levels of the expanded Green's function are retained. Details of the parameterization of the CdTe BOP are discussed elsewhere [25]. This parameterization considers properties of a variety of elemental and compound configurations (with coordination from 1 to 12) including clusters, bulk lattices, defects, and surfaces, in addition to the crystalline growth of vapor deposition.
There are currently two other CdTe interatomic potentials available in the literature, one [26] uses the Stillinger-Weber format [14], and the other [17] uses a Rockett modification [27] of the Tersoff format (TR). To evaluate different potentials, cohesive energies of a variety of Cd, Te, and CdTe phases calculated using various models are compared with the corresponding values obtained from our high-level density functional theory (DFT) calculations in Fig. 1. Here, various lattices are abbreviated as diamond-cubic (dc), simplecubic (sc), body-centered-cubic (bcc), face-centeredcubic (fcc), hexagonal-close-packed (hcp), graphite (gra), graphene (grap), γ-Se (A8), zinc-blende (zb), wurtzite (wz), NaCl (B1), and CsCl (B2). For clarity, these structures are sorted to give monotonic DFT energy trends. In Fig. 1, the unfilled stars show the experimental cohesive energies [28] of the equilibrium phases, while the straight lines connecting the neighboring data points merely guide the eye. Because the DFT calculations typically give accurate energy trends but not the absolute energies, the cohesive energies obtained from DFT calculations are scaled to match the experimental values for the equilibrium phases. Fig. 1 indicates that the cohesive energies calculated from BOP (solid lines) slightly oscillate (i.e., are non-monotonic) around the DFT benchmarks (thick light lines); however, these variations are quite minor. In fact, the BOP energy trends are considerably closer to those predicted by DFT than the corresponding results of the SW and TR parameterizations. Most importantly, BOP correctly specifies the lowest energies for the equilibrium phases of both elements and the compound, namely, the hcp Cd, the A8 Te, and the zb CdTe, and the calculated cohesive energies of the lowest energy phases also match the corresponding experimental values. In sharp contrast, the lowest energy phases are calculated to be dc Cd, dc Te, and zb CdTe by the SW parameterization and dc Cd, bcc Te, and B2 CdTe by the TR parameterization, with SW having the only correct result of zb CdTe. These results indicate that the TR parameterization cannot be used to study any of the equilibrium Cd, Te, and CdTe phases as the structures will not even be stable in MD simulations. While the SW parameterization can be used in some sort of MD simulation to study the equilibrium CdTe phase, caution should be taken in explaining the results involving defects as the potential is not transferrable to Cd and Te (and hence the defective) regimes. As a result, our new CdTe BOP approach significantly improves over other widely-used potentials on energy trends of different configurations leading to better description of defects.
The geometry of the melt-growth MD simulations is shown in Fig. 2. An initial zb CdTe crystal containing  (044) layers (about 30Å) in the z-direction was first created using the equilibrium lattice constant. Periodic boundary conditions were used in all three coordinate directions so that the system can be viewed as infinitely large. The two ends in the x-direction containing 2 (400) planes (about 1.8Å) were maintained at 0 K temperature (i.e., atom positions were fixed) so that these two regions acted as seeding crystals. The two regions containing 24 (400) planes (about 41Å) adjacent to the fixed ends were controlled at a low temperature of T low = 1000 K. A region containing 128 (400) planes (about 220Å) in the middle of the sample was controlled at a high temperature of T high . The remainder of the sample was left free. The middle portion of the sample was then melted by first giving a random displacement to all the atoms not in the fixed regions, and then annealing the system using a MD simulation with T high > 2200 K. The melt-growth of the crystal was then simulated in a second MD run at a different desired T high temperature where the T low and T high regions grew/shrank at a nominal growth rate of R = 0.2Å/ps, as shown in Fig. 2 (the accelerated growth rate is required for MD simulations).
Simulations were performed at two different high temperatures of T high = 2200 K and T high = 1800 K. The projected x-y configurations are shown in Fig. 3 as a function of time. The temperature profiles along the xdirection are superimposed on the atomic configurations. Fig. 3(a) shows the configuration prior to the growth where the middle section of the sample was melted at a temperature of T high = 2200 K. Fig. 3(b) shows that at t = 0.3 ns, the crystal/melt interface moved to a location corresponding to a temperature T i between 1300 K and 1500 K in both the T high = 2200 K and T high = 1800 K growth conditions, in good agreement with the melting temperature of CdTe, T m = 1365 K. In contrast and temperature profiles (along the x-direction). (a) prior to growth; (b) 0.3 ns after growth started using two growth conditions of T high = 2200 K and T high = 1800 K; and (c) 0.9 ns after growth started using T high = 2200 K and T high = 1800 K.
to conventional SW or TR-based potentials, our BOPbased MD method allows a physically-correct crystalline growth from the melt while at the same time predicts accurate energy trends of various configurations. Fig. 3(c) indicates that the temperature profiles dropped to 1000 K across the sample length at t = 0.9 ns for both simulation conditions. Interestingly, the sample became entirely crystalline at T high = 1800 K but the middle portion remained amorphous at T high = 2200 K. This amorphous zone did not change with further simulation at 1000 K, but was found to crystallize when the temperature of the middle portion was raised to 1500 K and then slowly cooled to 1000 K. Fig. 3, therefore, reveals how defects are trapped at a high growth rate.
To further explore defects in the grown crystals, the crystallinity parameter developed previously [12] was calculated along the x-direction for the two configurations shown in Fig. 3(c), and the results are shown in Figs. 4(a) and 4(b) respectively. Fig. 4(a) shows that the crystallinity for the sample obtained from the T high = 2200 K condition is uniform along the sample length except between 170Å and 280Å where the crystallinity exhibits a drop, in agreement with the trapped amorphous phase in this region. Surprisingly, the sample obtained at T high = 1800 K also exhibits a crystallinity drop near the center of the sample even though no amorphous zone is identified. To understand this, the low crystallinity region, which is framed in Fig. 4(b), is magnified and examined in Fig. 5(a). Fig. 5(a) indicates a planar defect where the projected Cd → Te stacking, shown as an arrow in the negative x-direction, is rotated by about 109.47 o clockwise. Such a rotation can originate from a [211]/6 slip on the (111) plane. The defect, therefore, corresponds to a stacking fault, and causes the drop in crystallinity observed in Fig. 4(b). Similar stacking faults have been observed in experiments [29,30]. A time-resolved analysis surprisingly indicates that the growth front is not always perpendicular to the growth direction; rather, it forms local trailing {111} facets, suggesting a higher stability and a slower solidification on these planes. In particular, the formation of such facets was often accompanied by the nucleation of stacking faults, Fig. 5(b). This discovery is further supported by a previous MD study on gold melt-growth, where {111} planes were found to cause stacking faults and slow solidification kinetics [31]. While mechanical twinning is known to be caused by stresses, our simulations provide a mechanistic explanation for the formation of a stacking fault: when growth occurs on a hexagonal {111} plane, say, plane A, atoms can occupy two different sites B and C. One of these sites corresponds to the lattice sites and the other corresponds to stacking fault sites. The energy difference between these two sites is small, and hence, it is not surprising that the condensation on a {111} plane may nucleate a stacking fault defect. Since the (100) and (110) planes have only one low energy (lattice) site, such stacking faults do not form if growth occurs strictly on (100) or (110) planes. These results suggest that using {100} or {110} growth planes can reduce stacking faults if the growth technique permits a sufficiently high temperature gradient and a sufficiently low growth rate to prevent the formation of local {111} interfaces. However, if the growth technique does not permit a sufficiently high temperature gradient and a sufficiently low growth rate, the {111} growth planes might be beneficial despite its high propensity for stacking faults because this interface is more likely to remain flat to minimize dendritic growth.
In conclusion, we demonstrate that (a) new BOP-based MD simulations can accurately predict melt-growth of semiconductors; (b) strictly derived from a quantum-  Fig. 3(c). (a) T high = 2200 K and (b) T high = 1800 K. The framed region will be further analyzed in Fig. 5. mechanical formalism, BOP enables defect formation mechanisms to be studied at a scale comparable to empirical MD methods and a fidelity approaching quantummechanical methods; and (c) amorphous defects can be trapped, and stacking faults can nucleate on {111} facets surprisingly formed during the melt-growth of semiconductors in non <111> directions. This work is supported by the NNSA/DOE Office of Nonproliferation Research and Development, Proliferation Detection Program, Advanced Materials Portfolio. Sandia National Laboratories is a multi-program laboratory managed and operated by Sandia Corporation, a wholly owned subsidiary of Lockheed Martin Corporation, for the U.S. Department of Energy's National Nuclear Security Administration under contract DE-AC04-94AL85000.