A Stillinger-Weber Potential for InGaN

Reducing defects in InGaN films deposited on GaN substrates has been critical to fill the “green” gap for solidstate lighting applications. To enable researchers to use molecular dynamics vapor deposition simulations to explores ways to reduce defects in InGaN films, we have developed and characterized a Stillinger-Weber potential for InGaN. We show that this potential reproduces the experimental atomic volume, cohesive energy, and bulk modulus of the equilibrium wurtzite / zinc-blende phases of both InN and GaN. Most importantly, the potential captures the stability of the correct phase of InGaN compounds against a variety of other elemental, alloy, and compound configurations. This is validated by the potential’s ability to predict crystalline growth of stoichiometric wurtzite and zinc-blende InxGa1-xN compounds during vapor deposition simulations where adatoms are randomly injected to the growth surface.


Introduction
Solid-state lighting (SSL) has the potential to broadly replace conventional light sources with tremendous energy advantages (Schubert, Kim, Luo, & Xi, 2006;Shur, & Zukauskas, 2005;Bergh, Craford, Duggal, & Haitz, 2001;Tsao, 2004;Steigerwald et al., 2002).However, the "green gap" (Krames et al, 2007;Wu, 2009) in SSL efficiency remains a major hurdle in this technological transformation.In principle, red, blue, and green lights are needed in order to create the white light suitable for human eyes.While high efficiencies have been achieved for blue and red light, the emission efficiency for green-to-yellow light is significantly lower.The devices used to emit green light are primarily based on InGaN films grown on GaN substrates.The poor performance of the devices has been attributed to various defects created due to the large lattice mismatch between InGaN and GaN (Liliental-Weber et al., 2001;Song, 2005;Lu, Li, Li, & Zhang, 2004;Pereira, 2006;Liu et al., 2006;Jahnen et al., 1998;Liu et al., 2006;Iida et al., 2013;Srinivasan et al., 2003).Molecular dynamics (MD) simulation of growth of InGaN films on GaN substrates provide an effective theoretical means to explore defect reduction strategies.This requires an In-Ga-N potential capable of crystalline growth simulations.
Stillinger-Weber (SW) potentials (Stillinger, & Weber, 1985) use an energy penalty for non-tetrahedral bond angles to ensure the lowest energy for tetrahedral structures, such as a diamond-cubic, zinc-blende, or wurtzite crystals.As a result, SW potentials normally predict crystalline growth of tetrahedral phases.SW potentials have been widely used for semiconductors including multi-element Ga-In-As-Se-Te (Grein et al., 1997), Al-Ga-In-P-As-Sb (Ichimura, 1996), Cd-Te-Zn-As-Si (Zhang et al., 2011), In-Ga-N (Lei, Chen, Jiang, & Nouet, 2009), and Cd-Zn-Hg-Se-S-Te (Zhou et al., 2013) systems.However, the literature In-Ga-N (Lei, Chen, Jiang, & Nouet, 2009) SW potential does not provide a complete set of parameters (in particular, only parameters between different species are given but parameters between the same species are missing).The purpose of the present work is to provide a new complete set of SW parameters for In-Ga-N system so that the community can begin to use MD simulations to explore improved InGaN/GaN films for solid-state lighting applications.

Stillinger-Weber Interatomic Potential
In the SW potential (Stillinger, & Weber, 1985), the total energy of a system composed of N atoms is written as where i 1 , i 2 , …, i N is a list of neighbors of atom i, θ jik is the bond angle formed by atoms j and k at the site of atom i, ( ) V r and ( ) V r are, respectively, pairwise repulsive and attractive functions, u IJ (r ij ) is another pair function, and subscripts ij and IJ indicate, respectively, the pair of atoms and the species of the pair of atoms.Here the three pair functions V r , and u IJ (r ij ) are expressed respectively as ( ) ( ) ) where ε, σ, a , λ, γ, A and B are seven pair dependent parameters, and p and q are two additional parameters that are usually set to 4 and 0 respectively.Note that a •σ represents the cutoff distance because all three pair functions V R (r), V A (r) and u IJ (r) vanish at r = a •σ.

Parameterization
As described above, SW potentials are fully defined by seven parameters ε, σ, a, A, B, λ, and γ for each atomic pair (Stillinger, & Weber, 1985).Not all these parameters need to be fitted.For instance, A is an over-specified parameter that can be merged into ε so that we fix A = 7.9170 for pairs InIn and GaGa and A = 7.0496 for the remaining pairs NN, InN, GaN, and InGa.Note that these different A values are purely due to historical reasons and they do not impact results.Furthermore, we assume λ = 32.5 and γ = 1.2 since these values have been found to give good results (Grein et al., 1997;Ichimura, 1996;Wang, & Stroud, 1990;Kioseoglou et al., 2003).We then fit ε, σ, a, B of each of the six pairs IJ (I,J = In, Ga, N) to target properties of elements (In, Ga, and N), alloys (InGa) and compounds (InN, and GaN).
Discussion of the parameterization methods have been documented in detail previously (Zhou et al., 2013) and hence will not be repeated.A complete set of the parameters we obtained is listed in Table 1.An electronic copy of the potential in a format compatible with the MD code LAMMPS (LAMMPS, 2017;Plimpton, 1995) is also provided (http://www.ccsenet.org/journal/index.php/jmsr/article/view/70840/38621).Characteristics of the potential, including its comparison with the target properties, will be provided in the next section.

Relative Energies of Different Phases
To examine the lowest energy phase prescribed by the potential, relaxed cohesive energies (per atom unit) of various elemental (In, Ga, N) and binary (In-Ga, In-N, Ga-N) phases are calculated using molecular statics simulations.Here elemental phases include face-centered-cubic (fcc), hexagonal-close-packed (hcp), bodycentered-cubic (bcc), simple-cubic (sc), and diamond-cubic (dc) crystals, and binary phases include wurtzite (wz), zinc-blende (zb), B1 (NaCl type), and B2 (CsCl type) crystals.Our calculated results are summarized in Table 2, where the lowest energy phases are marked in bold.It can be seen that fcc or hcp elemental crystals, and wz or zb binary crystals, have lower energies as compared with the other phases included in Table II.Hence, our potential captures correctly the observed wz and zb phases for the compounds InN and GaN.The equilibrium phases of elemental In, Ga, and N are I4/mmm, Cmca, and N 2 respectively (Donnay, & Ondik, 1973).SW potentials cannot capture these phases; however, the purpose of the potential is to study the growth of InGaN films under stoichiometric conditions where elemental phases never occur.As a result, particular phases of pure In, Ga, N are not important as long as the energies of the lowest-energy model-phases reproduce the energies of the experimentally observed phases.This ensures that the relatively stabilities of compounds with respect to elements are correctly captured and elements indeed do not occur during simulations.In this sense, our potential captures the stability of the compounds very well.

Agreement with Target Energetics
To provide users with an idea on the thermodynamic transferability, relaxed cohesive energies of the lowest-energy phases for elements and compounds are compared with the relevant target values in Table 3, where target energies of In, Ga, InN, and GaN are experimental values (Barin, 1993), and target energy for solid N is from our density function theory (DFT) calculations on an fcc lattice.It can be seen that our potential captures exactly the energies of all the characteristic phases.Note that although experimental and model phases are different for elemental In, Ga, and N, the exact reproduction of the experimental energies of the lowest-energy phases ensure that these model elemental phases do not occur during growth simulations just as during real deposition where the experimental elemental phases do not occur.

Lattice and Elastic Constants for zb and wz Phases of InN and GaN
The discussions above suggest that the potential is likely to correctly maintain the zb and wz phases of both InN and GaN during near-equilibrium MD simulations.To give users an idea as how well the InN and GaN compounds are simulated, relaxed lattice and elastic constants of zb and wz phases of InN and GaN compounds are calculated, and the results are summarized in Table 4.For the wz phases, the data is compared with available experiments (Donnay, & Ondik, 1973;Simmons, Wang, 1971;Derby, 2007), or DFT (Derby, 2007).The zb phases are not compared because like wz, zb also has tetrahedral bonding so that the zb phases are not expected to be far off once the wz phases are captured.It can be seen that when compared with the experimental lattice constants of the wurtzite phases, our potential predicts slightly smaller lattice constant a, and slightly larger lattice constant c, but the predicted atomic volume (scales with a 2 •c) reproduce exactly the corresponding experimental values.Our potential also reproduces exactly the experimental or DFT bulk moduli for wurtzite InN and GaN.Other elastic constants of the wurtzite InN and GaN predicted by our potential are also close to the corresponding experimental or DFT values.a: experimental data cited in (Donnay, & Ondik, 1973).
Compared with more complex potentials that can capture properties of a variety of phases (van Duin, Dasgupta, Lorant, & Goddard 2001), SW potentials can only ensure the correct lowest energy compound phases and their basic properties such as atomic volume, cohesive energy, and bulk modulus (Zhou et al., 2013).SW potential also cannot distinguish energies of zb and wz structures as can be seen in Table 2.This is a fair, perhaps even a good, approximation since first principle calculations indicate that the energy difference between wz and zb phases is small (~10 meV/atom) for both InN and GaN (Yeh, Lu, Froyen, & Zunger, 1992).Despite limitations of SW potentials, our InGaN SW potential is sufficiently good to explore defect formation during growth simulations as will be demonstrated in the next section.

Vapor Deposition Simulation Tests
The results shown in Table 2 only indicate that energies of wz and zb structures of the InN and GaN compounds are lower than those of the other phases listed in the table (or mixture of other phases at the same composition).The results do not prove that wz and zb structures of the InN and GaN compounds are more stable than any other configurations (which include an infinite number of amorphous structures).One effective method to demonstrate this stability is to perform MD simulations of vapor deposition.The advantages of vapor deposition simulations are that they sample a variety of configurations stochastically formed on the growth surface.If the potential predicts abnormally low energies for any wrong configuration, the simulations are likely to trigger an amorphous growth regardless temperature and growth rate.Likewise, if the potential correctly captures the lowest energy phase, the simulations will predict crystalline growth of the lowest energy phase near the melting temperature of the selected phase and at a growth rate feasible with the MD simulations (Zhou et al., 2012).Following the same approach used previously (Zhou et al., 2013;Ward et al., 2012a;Ward et al., 2012b;Ward, Zhou, Wong, & Doty, 2013), we performed four representative MD simulation tests to grow compound films in the [0001] direction of wurtzite substrates.In the first case, an InN film is grown on an InN substrate at a temperature of T = 1200 K, a growth rate of R = 0.088 nm/ns, and an incident adatom energy of 0.1 eV.In the second case, a GaN film is grown on GaN substrate at a temperature of T = 1600 K, a growth rate of R = 0.121 nm/ns, and an incident adatom energy of 0.1 eV.In the third case, an alloyed In 0.26 Ga 0.74 N layer is grown on a GaN substrate at a temperature of T = 1600 K, a growth rate of R = 0.110 nm/ns, and an incident adatom energy of E i = 0.1 eV.In the fourth case, a GaN/InN bilayer is grown on a GaN substrate at a temperature of T = 1800 K, a growth rate of R = 0.116 nm/ns, and an incident adatom energy of E i = 0.1 eV.Here the choice of different temperatures is to simply test different growth conditions.The atomic configurations obtained at the end of simulation (time t = 15.92 ns for the first case, 11.96 ns for the second case, 15.92 ns for the third case, and 31.84ns for the fourth case) are show respectively in Figs. 1 -4.Figs. 1 and 2 show that for homogenous InN-on-InN and GaN-on-GaN growth, large portion of the films are crystalline with the same orientation as the substrate, albeit defects exist under the extremely high deposition rates that have to be applied in the simulations.In particular, a majority of atoms do have tetrahedral bonds with neighboring atoms.Fig. 3 indicates that the growth of alloyed In 0.26 Ga 0.74 N compound is also crystalline.Although in the front view, the deposited film appears to be different from the substrate, this is caused by the formation of defects (e.g., misfit dislocations) that changed the alignment of atoms in the direction of the projection of the front view.In fact, the top view shows that the film is composed of domains of hexagonal regions with high crystal quality as indicated by the red circles, separated by defective boundaries.Similar results can be seen from Fig. 4 for the GaN/InN bilayer case.
The crystalline growth observed in Figs. 1 -4 provide strong validation that our potential captures the stability of the wz phase of In x Ga 1-x N compounds.Additional growth simulations using zb substrates also validated that the potential captures the stability of the zb phase of In x Ga 1-x N compounds.Hence, the potential is suitable for studying defect formation during growth of InGaN films on GaN substrates.

Conclusions
We have developed a Stillinger-Weber potential for In-Ga-N ternary systems.We demonstrate that this potential captures the lowest energies for the wz and zb structures of the InN and GaN compounds, prescribes lattice constants, cohesive energies, and elastic constants of the wz phases of the InN and GaN compounds that are in good agreement with the corresponding experimental values.Molecular dynamics vapor deposition simulations indicate crystalline growth of In x Ga 1-x N films, validating the potential's capability for exploring defect reduction strategies for In x Ga 1-x N-on-GaN solid-state lighting devices.

Figure 1 .
Figure 1.Atomic scale structure of the InN on InN obtained from molecular dynamics simulations.The pink shaded area in the front view highlights the initial InN substrate

Table 1
. SW potential parameters for In-Ga-N (energy in unit eV and length in unit Å).

Table 2 .
Relaxed energies (eV/atom) of selected structures.Note that bold numbers highlight the lowest energies

Table 3 .
Predicted and target cohesive energies (eV/atom) of the lowest-energy phases of elements and compounds * : from our DFT calculations.