Thermodynamics of native defects in In₂O₃ crystals using a first-principles method

Abstract

The stability of the intrinsic point defects in bixbyite In₂O₃, including oxygen vacancies, oxygen interstitials, indium vacancies and indium interstitials, under a range of temperatures, oxygen partial pressures and stoichiometries has been studied by computational methods. The calculated results indicate that the a-position is not a suitable site for interstitials and both the b-position and d-position are favorable for indium vacancies with similar defect formation energies. Both donors (oxygen vacancies and indium interstitials in the c-position) and acceptors (indium vacancies in the b/d-position and oxygen interstitials in the c-position) are predicted to have shallow defect transition levels. Then defect formation energies of all possible charged states are used in thermodynamic calculations to predict the influence of temperature and oxygen partial pressure and varied Fermi level in a limited area on the relative stabilities of the point defects. The combined formation energies of point defect complexes, including Frenkel pairs, anti-Frenkel pairs and Schottky pairs, are calculated to predict relative stability in the paper.

---

1. Introduction

Indium oxide (In₂O₃) is a wide band gap transparent conducting oxide (TCO) which is attracting a great deal of attention because of its high electrical conductivity, high optical transparency in the visible region, and high reflectivity in the infrared.^1,2 Bixbyite In₂O₃ is currently used in a wide range of optoelectronic applications, light-emitting diodes (OLEDs)³ including flat-panel displays, and solar cells. In addition, In₂O₃ is a suitable material for gas sensing applications especially of ozone and oxygen.⁴

The electron concentration and conductivity in n-type In₂O₃ have been studied experimentally since the 1970s.^5–7 It has been well concluded that the intrinsic defect contributes to the pronounced non-stoichiometry which results in n-type semiconductivity.⁶ Experimental techniques, such as thermogravimetric and conductivity measurements, have long been used to determine the deviation from non-stoichiometry in In₂O₃ as a function of temperature and oxygen partial pressure. From the experimental point of view, both the nonstoichiometry and the n-type conductivity have been assigned to doubly charged oxygen vacancies.⁵

Further research about anion vacancies in In₂O₃ was discussed.⁸ Lany and Zunger presented that conductivity should originate from a photo-induced process: the oxygen vacancy acts as color center binding two electrons.⁹ Gan and Lu et al. further demonstrated that oxygen vacancies promoting photoelectrochemical performance of In₂O₃ nanocubes.¹⁰ In addition, it is confirmed that oxygen vacancies can quantitatively account for the rather moderate conductivity in thin films surface.¹¹ From above, it is mostly implied that the native donor oxygen vacancies is the origin of the conductivity.

On the contrary, it is concluded that the origin of the native donor in undoped In₂O₃ is interstitial indium, but with the coexistence of V_O. Interstitial indium can generate a shallow donor level, close to the conduction band, and an even shallower donor level is formed when it associates with an oxygen vacancy.¹² V_O cannot generate free electrons by itself, but can facilitate emergence of the donor In_i.¹²

In the case of In₂O₃, several DFT studies have been applied to examine defect structures and stabilities. These studies have focused on, for example, Frenkel or oxygen vacancies and indium interstitials defect complexes,^12,13 extrinsic point defects^14,15 and the electronic structures of intrinsic point defects.¹⁶ Recent calculations of intrinsic defect formation energies (DFE) in In₂O₃ (ref. 17) find that oxygen vacancies are the predominant point defect under indium-rich conditions and oxygen dumbbell interstitials are the predominant point defect under oxygen-rich conditions. However, these studies have been restricted to special conditions of oxygen rich or meal rich limit. In addition, none of these studies address the relative stabilities of intrinsic point defect under varied experiment conditions, such as temperature and partial oxygen pressure.

In present paper, quantitative predictions of the stabilities of charged intrinsic point defects in bixbyite In₂O₃ have been performed using a combination of ab initio and thermodynamic calculation. We have performed the calculations on the basis of density of theory (DFT) using the GGA + U method to achieve a better description of the structure of the material. In particular, DFT calculations are used to study the structure and active sites about both perfect and defective crystal. This information is then used to determine DFEs. Importantly, a quantitative link is made to temperature and oxygen partial pressure and Fermi level, which are the key parameters for controlling the type and concentration of dominant defects in In₂O₃. The resulting self-consistent set of DFEs are not only necessary for predicting bulk point defect chemistry, but also are crucial input parameters for equilibrium, space-charge segregation models that can predict the equilibration of point with higher dimensional defects such as surfaces and grain boundaries.

2. Computational methods

2.1. Electronic structure computational

Density functional theory (DFT) calculations are all performed using Vienna ab initio simulation package (VASP).^18,19 Electron-ion interactions are treated using projector augmented wave (PAW) method.^20,21 The generalized gradient approximation (GGA) in the parameterization by Perdew, Burke and Ernzerhof (PBE)²² was selected to express the exchange–correlation potential. The electron wavefunctions are described using a plane-wave basis set with a cut-off energy of 500 eV.¹⁷ It is proved to be sufficient for defect calculations. The valence electrons for indium are generated in 4d¹⁰5s²5p¹ configuration and those for oxygen are in the 2s²2p⁴ configuration. Since the ‘standard’ GGA functional underestimates the binding energy of In 4d electrons, which leads to a considerable overestimation of the O 2p⁴–In 4d¹⁰ hybridization. Therefore we adopted the GGA + U method in the formulation by Dudarev et al.²³ with U − J = 7 eV.^17,24 This method has been successfully used in a recent study of the valence band structure of In₂O₃.²⁴ The Brillouin-zone integrations are carried out using the Monkhorst–Pack scheme with a k-mesh of 3 × 3 × 3 points, which proved to be sufficient for the defect calculation.

The defect calculations are carried out with fixed cell parameters as obtained for the ideal structure while the maximum force is converged to less than 10 meV Å⁻¹. In₂O₃ (bixbyite) has the symmetry point groups Ia3̄, space number 206 and lattice parameter a = 10.117 Å. The calculations are carried out of a supercell with 80 atoms. Fig. 1 shows the possible sites for interstitials and vacancies. Here, the calculated band gap of In₂O₃ is 1.79 eV, which is close to the experiment values.³¹

Fig. 1 Arrangement of atoms in In₂O₃. b and d site represent possible vacancy site for Indium atoms. Oxygen vacancy may locate on e site. a and c site are the interstitial site for indium or oxygen atoms.

2.2. Defect formation energies

The Gibbs free energy for point defect formation in In₂O₃ is given by a function of defect species α, charge state q, oxygen partial pressure P and temperature T.^25,26

ΔG_f(α, q, T, P) ≅ E^total(α, q) − E^total(perfect) + Δn_iμ_i(T, P) + q(E_VBM + E_F + ΔV) − TΔS(T, V₀) (1)

where ΔG_f(α, q, T, P) is the defect formation energy. E^total(α, q) is the relaxed total energy of supercell including the defect α of charge state q obtained from the DFT calculations, and E^total(perfect) is the relaxed total energy of the corresponding ideal supercell crystal also obtained from DFT calculations. Δn_i denotes the change of the number of atoms of type i from the ideal to defective crystal. For example, Δn_i = Δn_O = −1 for an oxygen interstitial, and Δn_i = Δn_In = 1 for an indium vacancy. μ_i(T, P) is the chemical potential of the defect atom i as a function of temperature and oxygen partial pressure, which is obtained from combination of DFE and thermodynamic calculations. E_VBM is the maximum of valence-band, and ΔV is the correction based on potential alignment correction method between ideal cell and defected crystal.³² E_F is the Fermi level which varied from 0 eV to band gap. V₀ is crystal volume which is constant in calculation. Lastly, ΔS is the vibrational entropy which is calculated from VASP.

2.3. Thermodynamic component

Chemical potential μ_i(T, P) is determined by the system environment. Oxygen chemical potential^27,28 is calculated as follows:

In this field, and are the chemical potential of In₂O₃ and In, which are calculated with VASP code, respectively. While is the Gibbs formation energy per mole in the standard state obtained from standard thermodynamic data.²⁹ Δμ⁰_O(T) is the difference chemical potential of oxygen between any temperature and standard temperature gained from thermodynamic data.³⁰ Combining eqn (1) and (2), we can get defect formation energies under different temperature and oxygen partial pressure. Table 1 shows a comparison of oxygen chemical potential, which calculated from eqn (2) and the defect formation energies for oxygen vacancies of +1 and +2 charges under different temperatures and oxygen partial pressures. This comparison indicates that temperature and oxygen partial pressure all play an important role, which should not be ignored in defect formation calculations.

Table 1 Calculated oxygen chemical potential from eqn (2) and defect formation energy from eqn (1) for vacancies with +1 and +2 charges at various temperatures and oxygen partial pressures

T (K)	q	EF (eV)	log(P/P^0)	μ (O)	qEF (eV)	Sum of eqn (1) (eV)
300	2	0	−2	−5.14	0.00	0.45
1400	2	0	−2	−6.70	0.00	−0.8
300	2	0	−10	−5.37	0.00	0.21
1400	2	0	−10	−7.81	0.00	−1.92
300	2	1.4	−2	−5.14	2.80	3.25
1400	2	1.4	−2	−6.70	2.80	2
300	2	1.4	−10	−5.37	2.80	3.01
1400	2	1.4	−10	−7.81	2.80	0.88
300	2	2.8	−2	−5.14	5.60	6.05
1400	2	2.8	−2	−6.70	5.60	4.8
300	2	2.8	−10	−5.37	5.60	5.81
1400	2	2.8	−10	−7.81	5.60	3.68
300	1	0	−2	−5.14	0.00	2.11
1400	1	0	−2	−6.70	0.00	0.86
300	1	0	−10	−5.37	0.00	1.87
1400	1	0	−10	−7.81	0.00	−0.26
300	1	1.4	−2	−5.14	1.40	3.51
1400	1	1.4	−2	−6.70	1.40	2.26
300	1	1.4	−10	−5.37	1.40	3.27
1400	1	1.4	−10	−7.81	1.40	1.14
300	1	2.8	−2	−5.14	2.80	4.91
1400	1	2.8	−2	−6.70	2.80	3.66
300	1	2.8	−10	−5.37	2.80	4.67
1400	1	2.8	−10	−7.81	2.80	2.54

---

3. Results and discussion

3.1. Charge transfers and band gap corrections

Fig. 2 shows the calculation results of the defect transition level for oxygen vacancy, indium interstitial on c and a sites (In_i–c, In_i–a), oxygen interstitial on c and a sites (O_i–c, O_i–a), indium vacancy on b and d sites (V_In-b, V_In-d). The results indicate that both donors and acceptors have shallow defect transition levels except for high ionized defects, such as and , which manifests that both donors and acceptors are prefer to remain in low ionized or neutral charge state. These shallow levels further indicate that the thermodynamic transitions are easily activated.

Fig. 2 Calculated defect transition levels ε (defect^q1/q2) for the major defects. The circles below and above the transition levels denote the occupancy of the corresponding defect charge state. The experiment band gap value is 2.8 eV.³¹

Usually, all the calculated transition levels that are scaled by the ratio of the experimental to the calculated data of the band gap using the scissor operator. This correction is only applied on conduction states and the DFEs of donors increased by an energy equal of the difference in one particle energies times the occupation of the corresponding defect levels. The choice to only correct the conduction band, however, is usually not justified, and furthermore the change in energy does not include defect level relaxation and the double counting correction term using this scheme.¹⁷ So the transition levels correction is omitted in this work.

3.2. Sites for vacancies and interstitials

Within the bixbyite lattice there are four distinct defect sites: a, c, b and d sites. The a-position and c-position are usually referred to as the structural vacancy sites, while b-position and d-position are acted as the structural interstitial sites (Fig. 1). The DFEs of interstitial in a-position is always much higher than that in c-position for both indium and oxygen interstitial. So the c-position is most favorable for indium (In_i–c) and oxygen interstitial (O_i–c) and a-position is not a suitable interstitial site although the temperature and oxygen partial pressure varied. However, the b-position and d-position are both favorable for indium vacancy as the defect formation energy are much closed to each other. The DFEs of indium vacancy in d-position is a little lower than in b-position when Fermi level varied from 0 eV to 0.85 eV and a reverse phenomenon is found when Fermi level position locates in the 0.85–2.8 eV region.

3.3. Defect formation energy

The calculation of predicting defect formation energies may be used to obtain defect concentrations and well predict doping behavior. Therefore, the defect formation energy of different charged defects are calculated using our approach as a function of system conditions. Fig. 3a–c shows the calculated DFEs of V_In-b, V_In-d, In_i–b, In_i–d, V_O, O_i–c and O_i–a in pure In₂O₃ as a function of Fermi level E_F under reducing conditions (log(P_O2/P⁰) = −10) at three different temperatures. The factual Fermi level of In₂O₃ is varied by dopants or impurities, in addition to the system temperature (T) and oxygen partial pressure (P_O2). It is therefore necessary to determine the DFEs over a physically reasonable range of Fermi energies to increase the correlation of the calculated values. So the Fermi level of all possible values is considered here.

Fig. 3 Calculated DFE of point defects (V_O, In_i, O_i and V_In) is a function of Fermi level, oxygen partial pressure and temperature (a–f). (a–c) shows formation energies at different temperatures (T = 300 K in (a), T = 800 K in (b) and T = 1400 K in (c) in the reduced state (log(P_O2/P⁰) = −10). (d–f) shows formation energy at different temperatures T = 300 K in (d),T = 800 K in (e) and T = 1400 K in (f) in the near-neutral state (log(P_O2/P⁰) = −2).

When the Fermi level varies from 0 to 2.8 eV at room temperature and low oxygen partial pressure, the most stable intrinsic point defects are shown to be , and , respectively (see Fig. 3a). In contrast, when the Fermi level is approximately equal to E_g/2 at room temperature, the stable charge states for these same point defects are , , and . The results also show that, as the Fermi level decreases, the thermodynamically preferred charge states of the indium interstitials change from neutral to +1, +2 and +3. Intriguingly, the p-type defect and are predicted to be the most stable under conditions of Fermi level E_F = 1.3–2.8 eV, T = 300 K and P_O2 = 10⁻¹⁰ atm. While n-type behavior always be find when experimental measurements are normally carried on In₂O₃ (ref. 5–7).

Over almost the full range of Fermi levels considered here, oxygen vacancies of various charge states are predicted to be the most stable intrinsic defect in In₂O₃ at high temperature, except at very high Fermi levels.

The most stable intrinsic defect as a function of Fermi level, oxygen partial pressure and temperature are determined in the three-dimensional defect formation energy diagrams which are divided into three regions (Fig. 4). The graph shows the most stable point defect is oxygen vacancies when the external environment is located in the region 1. When E_F is in the range of 0–1.66 eV, 1.66–2 eV and 2–2.8 eV, the oxygen vacancies of +2, +1 and 0 charge states are predicted to be the most stable intrinsic defect in the region 1, respectively. The most stable point defect has a transfer to oxygen interstitials and indium vacancies when the experimental environment fit to the region 2 and region 3, respectively. The most stable point defects of oxygen interstitials change from −1 to −2 charge state in region 2 when E_F changes from 1.26 to 2.8 eV and the most stable point defects of indium vacancies maintain −3 charge state in whole region 3. And, intriguingly, when E_F is less than 1.05 eV, the most stable intrinsic defect is always despite the temperature and oxygen partial pressure varying in Fig. 4.

Fig. 4 Three-dimensional defect formation scheme as a function of temperature, Fermi level and oxygen partial pressure. The space in graph is divided into three regions, region 1, region 2 and region 3.

Fig. 5 shows two-dimensional defect formation scheme as a function of temperature and oxygen partial pressure with the E_F = 1.4 eV, 2.4 eV respectively. Fig. 5a (E_F = 1.4 eV) show that the ordering of most stable defects are oxygen interstitials → oxygen vacancies with the temperature increases and partial pressure decreases (going from bottom right to top left in each figure). In addition, Fig. 5b (E_F = 2.4 eV) show the ordering of most stable defects are indium vacancies → oxygen interstitials → oxygen vacancies. The predicted preferred defects with different Fermi levels are disparate, although it has the same temperature and oxygen partial pressure. For example, and are predicted to be the most stable intrinsic defects in the low temperature region with the E_F = 1.4, 2.4 eV respectively, whereas, and are predicted to be the most stable intrinsic defects in the region of high temperature and reduced conditions.

Fig. 5 Two-dimensional defect formation scheme as a function of temperature and oxygen partial pressure calculated at two different Fermi level (E_F = 1.4 eV in (a) and E_F = 2.4 eV in (b), respectively).

The most stable charge states and corresponding formation energies of distinct intrinsic point defects in In₂O₃ at three typical conditions are listed in Table 2. Since it is unconcern the equilibrium concentration of defects or impurities can be achieved for any experimental conditions which have an influence on Fermi energy. Therefore, it is significant to computer all possible charge states under different external experiment and explore the relative stabilities of point defects in indium oxide.

Table 2 Calculated defect formation energies for their most stable charge states of defects under three typical conditions: standard condition (T = 300 K, E_F = 0.4 eV, P_O2 = 1 atm), reduced condition (T = 1200 K, E_F = 1.4 eV, P_O2 = 10⁻⁶ atm) and oxidized condition (T = 700 K, E_F = 2.4 eV, P_O2 = 10² atm)

	Standard condition		Reduced condition		Oxidized condition
	Charge	DFE (eV)	Charge	DFE (eV)	Charge	DFE (eV)
VIn-b	−1	4.87	−3	4.59	−3	−0.01
VIn-d	−1	4.81	−3	4.66	−3	0.06
VO	+2	1.31	+2	1.77	0	3.97
Ini–c	+3	3.95	+3	4.32	+1	8.39
Ini–a	+3	4.56	+3	4.96	+1	8.88
Oi–c	0	3.31	−2	3.16	−2	0.18
Oi–a	0	3.77	−2	3.96	−2	0.98

---

In the preceding discussion, the relative stabilities of individual defects are considered and discussed. However, point defects do not occur alone but as part of charge-neutral complexes. The DFEs of Frenkel, anti-Frenkel and Schottky defect complexes are therefore also calculated using the following equations:

ΔG^F_f = E^total(V_In, −q) − E^total(perfect) + E^total(In_i, q) − E^total(perfect) − T(ΔS_VIn(T, V) + ΔS_Ini(T, V)) (3)

ΔG^AF_f = E^total(V_O, q) − E^total(perfect) + E^total(O_i,−q) − E^total(perfect) − T(ΔS_VO(T, V) + ΔS_Oi(T, V)) (4)

ΔG^S_f = 2(E^total(V_In,−3q) − E^total(perfect)) + 3(E^total(V_O, 2q) − E^total(perfect)) + 2μ_In + 3μ_O − T(2S_VIn(T, V) + 3S_VO(T, V)) (5)

The definitions of E^total(α, q) and E^total(perfect) are the same as in eqn (1).

These defect complexes are treated as a combination of charged or neutral defects and the total charge states of the defect complex is neutral. For instance, four different types of Frenkel pair modes are considered: , , and . Similar combination is considered for anti-Frenkel and Schottky pairs. A nice feature of Analysis of the defect complexes is that, since the defect complexes are charge neutral, the calculations doesn't consider the Fermi level, as the formation energy form a transverse line in Fig. 6.

Fig. 6 Calculated DFE of defect complexes at 1400 K: (a) Frenkel pairs and (b) anti-Frenkel pairs. The DFEs of the complexes are compared with the DFEs of individual point defects at 1400 K and oxygen pressures of (log(P_O2/P⁰) = −10). The formation energy of combination of charged indium vacancy and interstitial in b-position or in d-position is much closed. For example, means that the DFE of Frenkel is a little lower than that of the Frenkel .

The two more stable types Frenkel pairs including V_In-b + In_i–c and V_In-d + In_i–c compared with other Frenkel pairs are shown in the Fig. 6a. The most stable Frenkel pair is and combination is also stable since its DFE is only a little higher than the former's about 0.1 eV.

From the Fig. 6b, it can be seen that all the DFEs of anti-Frenkel pairs are lower than that of Frenkel pairs. The combination of V_O + O_i–c is more stable than V_O + O_i–a for any charged state, which is well agreement with the fore-discussion that a-position is not a suitable interstitial site.

From Fig. 6a and b, it can be concluded that charged defect pairs have lower DFEs. For instance, the combination () has a lower DFE. The anti-Frenkel pair is more stable than the Frenkel pair. In addition, the DFEs of Schottky pairs is much higher than that of Frenkel and anti-Frenkel pairs. Therefore anti-Frenkel pairs are the most stable and Schottky pairs are lowest thermodynamically stable. However, Aron believes that the DFEs of the anti-Frenkel pair is moderately higher than that of the Schottky defect energy, nevertheless, it is still lower than that of the Frenkel pair.¹³ The reason to the argument is that their calculation formula of the DFEs has some difference from ours.

4. Conclusions

In summary, the native point defects have been studied at varied temperature and confined oxygen partial pressure in In₂O₃ using ab initio methods based on the GGA + U formalism. As to the In₂O₃ structural, we find the c-position is the favorable interstitial site and b/d-position has same role for indium vacancies with similar DEFs. It indicates that oxygen vacancies is the most stable point defect in n-type behavior and indium vacancies or oxygen interstitials is the most stable defect in p-type behavior. In addition, as to defect complexes, we conclude that anti-Frenkel pairs is most thermodynamic stable than Frenkel pairs and Schottky pairs and charged defect complexes are more likely to occur, such as . Of considering key parameters to calculate defect formation energy, it is therefore generally applicable for the systematic evaluation of defect formation in electronic ceramics under the full range of temperature and atmospheric conditions, such as are important for catalysis (low temperature, high pressure) and gas sensing (high temperature, low pressure) applications.

Acknowledgements

This work is supported by the Hujiang Foundation of China (A14001).

References

1. R. G. Gorden, MRS Bull., 2000, 25, 52–57.
2. H. Hartnagel, A. L. Dawar, A. K. Jain and C. Jagadish, Semiconducting transparent thin films, Institute of Physics Pub., England, 1995, pp. 110–126.
3. L. S. Hung and C. H. Chen, Mater. Sci. Eng., R, 2002, 39, 143–222.
4. V. Golovanov, M. A. Mäki-Jaskari, T. T. Rantala, G. Korotcenkov, V. Brinzari, A. Cornet and J. Morante, Sens. Actuators, B, 2005, 106, 563–571.
5. J. H. W. De Wit, G. Van Unen and M. Lahey, J. Phys. Chem. Solids, 1977, 38, 819–824.
6. J. H. W. D. Wit, J. Solid State Chem., 1975, 13, 192–200.
7. J. H. W. D. Wit, J. Solid State Chem., 1973, 8, 142–149.
8. I. Tanaka, K. Tatsumi, M. Nakano and H. Adachi, J. Am. Chem. Soc., 2002, 85, 68–74.
9. S. Lany and A. Zunger, Phys. Rev. Lett., 2007, 98, 045501.
10. J. Gan, X. Lu, J. Wu, S. Xie, T. Zhai, M. Yu, Z. Zhang, Y. Mao, S. C. I. Wang, Y. Shen and Y. Tong, Sci. Rep., 2013, 3, 1021.
11. S. Lany, A. Zakutayev, T. O. Mason, J. F. Wager, K. R. Poeppelmeier, J. D. Perkins, J. J. Berry, D. S. Ginley and A. Zunger, Phys. Rev. Lett., 2012, 108, 016802.
12. T. Tomita, K. Yamashita, Y. Hayafuji and H. Adachi, Appl. Phys. Lett., 2005, 87, 051911.
13. A. Walsh, C. R. A. Catlow, A. A. Sokol and S. M. Woodley, Chem. Mater., 2009, 21, 4962–4969.
14. O. Warschkow, D. E. Ellis, G. B. González and T. O. Mason, J. Am. Chem. Soc., 2003, 86, 1707–1711.
15. T. O. Mason, G. B. Gonzalez, D. R. Kammler, N. Mansourian-Hadavi and B. J. Ingram, Thin Solid Films, 2002, 411, 106–114.
16. O. Warschkow, D. E. Ellis, G. B. González and T. O. Mason, J. Am. Chem. Soc., 2003, 86, 1700–1706.
17. P. Ágoston, P. Erhart, A. Klein and K. Albe, J. Phys.: Condens. Matter, 2009, 21, 455801.
18. G. Kresse and J. Furthmüller, Phys. Rev. B: Condens. Matter Mater. Phys., 1996, 54, 11169.
19. G. Kresse and J. Furthmüller, Comput. Mater. Sci., 1996, 6, 15–50.
20. P. E. Blöchl, Phys. Rev. B: Condens. Matter Mater. Phys., 1994, 50, 17953.
21. G. Kresse and D. Joubert, Phys. Rev. B: Condens. Matter Mater. Phys., 1999, 59, 1758–1775.
22. J. P. Perdew and K. Burke, Phys. Rev. Lett., 1996, 77, 3865–3868.
23. S. L. Dudarev, G. A. Botton, S. Y. Savrasov, C. J. Humphreys and A. P. Sutton, Phys. Rev. B: Condens. Matter Mater. Phys., 1998, 57, 1505–1509.
24. P. Erhart and A. Klein, Phys. Rev. B: Condens. Matter Mater. Phys., 2007, 75, 153205.
25. S. B. Zhang and J. E. Northrup, Phys. Rev. Lett., 1991, 67, 2339–2342.
26. G. X. Qian, R. M. Martin and D. J. Chadi, Phys. Rev. B: Condens. Matter Mater. Phys., 1988, 38, 7649–7663.
27. I. G. Batyrev, A. Alavi and M. W. Finnis, Phys. Rev. B: Condens. Matter Mater. Phys., 2000, 62, 4698–4706.
28. M. W. Finnis, A. Y. Lozovoi and A. Alavi, The oxidation of NiAl: What can we learn from ab initio calculations?, Annu. Rev. Mater. Res., 2005, 35, 167.
29. W. M. Haynes and D. R. Lide, Handbook of chemistry and physics, C.R.C. Press, New York, 2010.
30. N. I. ST. NIST Chemistry WebBook. NIST, 2005.
31. A. Klein, Appl. Surf. Sci., 2000, 77, 2009–2011.
32. J. Lee and S. Han, Phys. Chem. Chem. Phys., 2013, 15, 18906–18914.

---