A Comparative Ab Initio Thermodynamic Study of Oxygen Vacancies

Article pubs.acs.org/JPCC

A Comparative Ab Initio Thermodynamic Study of Oxygen Vacancies in ZnO and SrTiO3: Emphasis on Phonon Contribution Denis Gryaznov,*,†,‡ Evgeny Blokhin,† Alexandre Sorokine,‡ Eugene A. Kotomin,†,‡ Robert A. Evarestov,§ Annette Bussmann-Holder,† and Joachim Maier† †

Max Planck Institute for Solid State Research, Stuttgart, Germany Institute for Solid State Physics, University of Latvia, Riga, Latvia § Department of Quantum Chemistry, St. Petersburg State University, Peterhof, Russia ‡

ABSTRACT: Using a hybrid Hartree−Fock (HF)-DFT method combined with LCAO basis set and periodic supercell approach, the atomic, electronic structure and phonon properties of oxygen vacancies in ZnO and SrTiO3 were calculated and compared. The important role of a ghost basis function centered at the vacant site and defect spin state for SrTiO3 is discussed. It is shown that the use of hybrid functionals is vital for correct reproduction of defects basic properties. The Gibbs free energy of formation of oxygen vacancies and their considerable temperature dependence has been compared for the two oxides. These calculations were based on the polarizability model for the soft mode temperature behavior in SrTiO3. The supercell size effects in the Gibbs free energy of formation of oxygen vacancies in the two oxides are discussed. The major factors for the quite different behavior of the two oxides and the degree of electron delocalization nearby the oxygen vacancy have been identified.

More than 30 first-principles calculations dealt so far with the atomic and electronic structure of the oxygen vacancy in SrTiO3 (see ref 15 and references therein) and even more in ZnO (see, for example, an excellent review paper by Janotti and Van de Walle7 and references therein). However, the common problem of most of these calculations is use of the standard DFT exchange-correlation functionals (LDA/GGA) or parameter-dependent functionals like hybrid Heyd−Scuseria−Ernzerhof (HSE) or DFT+U. The LDA/GGA approach strongly underestimates the band gap in both ZnO and STO. As a result of such calculations, the energy level of the vacancy in STO (which is a shallow donor) falls into the conduction band instead of being in the band gap. On the other hand, different fractions of Hartree−Fock (HF) exact exchange were used for ZnO within the HSE functional (see below). Below we present a comparison between our results and those obtained by means of HSE functional implemented in plane wave codes, if not otherwise stated. Lastly, the DFT+U functional in STO suggested several different scenarios for the position of the defect level, depending on the on-site Coulomb interaction (Hubbard U-parameter).16,17 Its choice affects not only the band gap value but also the energetic preference of the oxygen vacancy in different charge states. Another common feature of the standard DFT calculations combined, as a rule, with the plane wave (PW) basis set is negligible electron density concentration in the oxygen vacancy, in contrast to the

1. INTRODUCTION Zinc oxide (ZnO) and strontium titanate (SrTiO3, hereafter STO) belong to the most prominent functional materials with numerous applications.1,2 Both have similar band gaps around 3.4 eV but reveal different atomic structure: wurtzite in ZnO3 and ABO3-type cubic perovskite in STO.4 The oxygen vacancy is a very common defect in these and related materials, controlling mass transport and other device properties, e.g. permeation membranes and solid oxide fuel cell cathodes.5−8 The two oxides are very different with respect to their dielectric properties. SrTiO3 is an incipient ferroelectric material with extremely high dielectric constant (the values of the order 104 were reported at low temperatures9) related to the soft phonon mode behavior with temperature. At low temperatures controlled by quantum statistics the ferroelectric instability is suppressed, and the compound has been defined as a quantum paraelectric.10,11 The very high dielectric constant of the bulk STO decreases, however, significantly in thin films (from several thousands to several hundreds12,13). Moreover, STO thin films14 demonstrate a change of the dielectric constant by a factor of 2, i.e. from 1000 to 450 (at 77 K) also for different concentrations of oxygen vacancies. On the other hand, the dielectric constant of ZnO is orders of magnitude smaller. We, thus, expect to see an influence of the materials dielectric properties on the phonon contribution to the vacancy formation energy in the two oxides. In addition, the same defect, i.e. the neutral oxygen vacancy, has quite different properties in the above-mentioned materials which, as shown below, have an additional effect on the oxygen vacancy formation energy. © XXXX American Chemical Society

Received: January 18, 2013 Revised: June 5, 2013

A

dx.doi.org/10.1021/jp400609e | J. Phys. Chem. C XXXX, XXX, XXX−XXX

The Journal of Physical Chemistry C

Article

eV per primitive unit cell. The Monkhorst−Pack scheme29 for 8 × 8 × 8 and 6 × 6 × 4 k-point mesh in the Brillouin zone (BZ) was applied for STO and ZnO, respectively, together with the tolerances 8, 8, 8, 8, and 16 for the Coulomb and exchange integrals calculations. Furthermore, the SCF convergence threshold on the total electronic energy is 10−7 au for ZnO and 10−10 au for STO for the lattice structure optimization and 10−10 au for the phonon frequency calculations in both oxides. High-frequency dielectric constants were estimated using the coupled perturbed HF (Kohn−Sham) method as implemented in the CRYSTAL code.30 The threshold on the energy first derivative change was 10−4. It allows us to determine with high accuracy the longitudinal optical (LO) phonon frequencies and static dielectric constants. We use the supercell approach31 to consider several defect concentrations of the neutral oxygen vacancies VO. The defect concentrations of 6.25% and 12.5% were modeled in STO with 80 atom and 40 atom supercells (see ref 15 for details). The analogous calculations for ZnO may require considering asymmetric extension of the unit cell due to its hexagonal structure. Thus, the supercells 2 × 2 × 2 and 3 × 3 × 2 were used for the VO calculations suggesting, consequently, the concentrations 6.25% and 2.78%. Only the Γ-point was taken into account in the phonon calculations of defective systems. The full structure optimization of the defective supercells was performed prior to the phonon frequency calculations. In the calculations of VO in both STO (oxygen site symmetry D4h) and ZnO (oxygen site symmetry C3v) we use neutral supercells and the so-called ghost BS as implemented in the CRYSTAL code.25 It assumes that the BS of a missing O atom remains in the vacancy with a zero core charge. This means that the electron density of two electrons stemming from a missing oxygen atom is divided between the oxygen site and neighboring cations. Thus, use of a ghost BS provides a reasonable initial guess for the electronic density distribution. B. Defect Formation Energy. The standard Gibbs free energy of oxygen vacancy formation, ΔG0F, as a function of temperature at constant pressure was calculated using the following two equations (oxygen-rich conditions):

considerable electron density in the vacancy obtained using more flexible calculations based on the atomic basis set.15,18 This effect can substantially affect the optimized atomic and electronic structure of the point defects. Most of the calculations so far dealt with the phonon properties of the bulk defect-free STO and ZnO. The phonon properties of defective STO (containing Fe impurities and oxygen vacancies) were discussed for the first time in ref 19, and, to our knowledge, no such calculations were performed for ZnO. It has been shown by us20,21 and other authors22 that the phonon properties are sensitive to the choice of the functional. The phonon frequencies in the solid, when properly calculated, allow us to obtain the correct temperature dependence of oxygen vacancy formation energies (thus, the Gibbs formation energies). In the present study, we aimed at a comparison of phonon properties including infrared-active modes for bulk defect-free STO and infrared and Raman-active modes for bulk ZnO as well as calculation of the Gibbs formation energies of oxygen vacancies within the same computational scheme without imposing any adjustable parameters. This paper is organized in the following way. In section 2, the method and computation parameters are discussed, together with the approach to calculate the Gibbs formation energy of oxygen vacancy. The main results in section 3 for the lattice structure and the electronic properties and phonon frequencies in STO and ZnO are subdivided into two subsections dealing with pure materials and defects therein. Lastly, the conclusions are given in section 4.

2. COMPUTATIONAL DETAILS A. Method and Computational Parameters. In the present calculations, we use a basis set (BS) of the linear combination of atomic orbitals (LCAO), Perdew−Burke− Ernzerhof (PBE0)24 hybrid exchange-correlation functional as implemented in the CRYSTAL09 computer code.25 The standard PBE0 hybrid functional supposes 25% of the exact exchange.24 It has been shown by us19,20,26 and other authors27 that PBE0 correctly reproduces the basic properties and defect behavior in different oxides. In our calculations on STO the BS for Sr, Ti and O atoms were taken from refs 19 and 28. The BS optimization for STO has already been discussed in our previous paper,20 whereas in the present study the all electron BS by Jaffe et al.23 was chosen for ZnO and reoptimized for the PBE0 functional. The exponents of Gaussian type orbitals smaller than 0.7 Bohr−2 were reoptimized. Thus, we are concerned only with the virtual orbitals and their contraction in comparison to original Jaffe’s BS (Table 1). Note that the BS optimization in ZnO led to a considerable energy gain of 0.29

VO VO VO ΔG F0(T ) = [Etot + Evib (T ) − TSvib (T ) + pV VO] p p p − [Etot + Evib (T ) − TSvib (T ) + pV p] 1 + μO0 (T ) 2 2

and O2 μO0 (T ) = Etot + EO2 − TS O2 + kT 2

where superscripts p and VO indicate perfect and defective (one VO per supercell) crystals, Etot total electron energies, Evib the phonon contribution to the internal energy including vibrations at T = 0 K, Svib the corresponding entropy, V the supercell volume, T the temperature, p the standard pressure and k Boltzmann’s constant. All the energies are given per supercell. For the calculation of the standard chemical potential μ0O2(T) of O2 molecule, we used the free energy of a gas-phase O2 molecule calculated using the LCAO method within an ideal gas model.19,33 The thermal energy EO2 and entropy SO2 include the translational, rotational and vibrational contributions to the chemical potential of O2. These were calculated using the CRYSTAL09 code (Table 2). The reference energy of the O atom in O2 molecule in eq 1 is traditionally presented as a sum

Table 1. The Exponents of Gaussian Type Orbitals (Bohr−2) after Optimization of Jaffes’s BS23 for PBE0 Functional in ZnO type of orbital sp sp d sp sp

Jaffe BS Zn 0.62679 0.15033 0.51592 O 0.536420 0.23973

(1)

present study 0.62668 0.14392 0.50107 0.40904 0.13710 B



Article

3. RESULTS AND DISCUSSION A. Bulk Properties. Atomic Structure and Electronic Properties. The calculated lattice parameters and bulk moduli agree very well with the experimental data and previous hybrid functional calculations for both oxides (Table 3). The Mulliken atomic charges evidence a considerable covalency of the Ti−O chemical bonding in STO (2.33 e (Ti), −1.46 e (O), 2.03 e (Sr)) and Zn−O bonding in ZnO (±0.95 e for Zn and O, respectively). As mentioned above, the underestimated band gaps in STO and ZnO within the standard LDA, GGA-type functionals were a subject of many discussions in the literature. The band gap ΔEg is significantly underestimated in both plane wave and LCAO calculations: 0.73 eV for ZnO45 and 2.10 (1.8) eV for direct (indirect) band gap in STO.20,47 In order to solve this problem, one should go beyond the standard DFT functionals, e.g. by using the hybrid functional, the DFT+U approach7,48 or self-interaction correction (SIC)-LDA technique.49 It is also well established for both STO and ZnO that the DFT+U approach is unable to reproduce correctly ΔEg by a simple variation of the U-parameter for the d electrons of the transition metal, i.e. the value of ΔEg remains underestimated by almost 1.5 eV in comparison with the experimental value for ZnO and by almost 1.0 eV for STO.17,45 A very sophisticated DFT+U study on ZnO was recently performed by Boonchun and Lambrecht50 where they separately adjusted the Uparameters for s, p, d electrons of Zn and s, p electrons of O, in order to obtain the correct ΔEg from their accurate GW calculation. On the other hand, Oba et al.45 have shown that the fraction of HF exact exchange must be increased to 0.375 for the HSE functional, in order to reproduce ΔEg (3.42 eV in Table 3) in ZnO, if the screening parameter is set to 0.2. A very similar value for the screening parameter and 0.25 for the HF exact exchange gave a very accurate band gap for STO22 (3.47 eV for the direct ΔEg in Table 3). A very elegant way to correct “the band gap problem” is suggested in the (SIC)-LDA

Table 2. The Properties of a Free O2 Molecule Calculated with PBE0 and HSE0 Exchange-Correlation Functionals and LCAOa equilibrium distance, Å binding energy, eV rotational temp, K vibrational temp, K a

PBE0

HSE32

expt33

1.20 5.30 (5.19) 2.11 2478.60

1.19 5.25 (5.14)

1.21 5.12 2.07 2230

2475.70

Numbers in parentheses are corrected for zero-point vibrations.

of the isolated O atom energy and half the binding energy of the O2 molecule.34 Unfortunately, the accuracy of the calculated properties of O2 for the defect formation energies is often overlooked in the literature. As is known, the binding energy of O2 is usually significantly overestimated by 1−1.5 eV using the standard DFT functional,35−37 which consequently could result in a defect formation energy error as large as 0.7 eV. In contrast, the LCAO method and the hybrid PBE0 functional used here allow us to reproduce very accurately the vibrational properties of O2 as well as the equilibrium bond length and the binding energy. As seen from Table 2, the properties of O2 as calculated in the present study are in good agreement with the experimental findings and previous HSE LCAO calculations. Such high accuracy LCAO calculations allow us to avoid any experimental data for defect formation energies. We have estimated that the VO formation energy is reduced by ∼2.0 eV within the temperature from zero to 1660 K, if the vibrations in the solid phase are ignored. The standard expressions for Evib and Svib are discussed in many reviews (e.g., refs 38 and 39). They require the knowledge of the phonon frequencies for perfect and defective crystals, which were calculated here by the direct method39 within the harmonic approximation and for optimized equilibrium crystal structures.

Table 3. Bulk Properties of Defect-Free ZnO (Wurtzite) and SrTiO3 (Cubic Perovskite)a ZnO PBE0 LCAO a, Å c, Å u, frac zMull(A), e zMull(B), e zMull(O), e zBorn(A), e zBorn(B), e zBorn(O1), e zBorn(O2), e ε0xx ε0zz ε∞ xx ε∞ zz ΔEg, eV B, GPa

3.26 5.22 0.38 0.95

HSEb PW 3.25 5.25

SrTiO3 expt

PBE0 LCAO

e,f

3.91

3.25 5.21e,f 0.38e,f

−0.95 2.00/2.04

2.06/2.11

2.16/2.22h

−2.00/−2.04

−2.06/−2.11

−2.16/−2.22h

6.67 7.31 3.25 3.28 3.56 153

7.72 8.37 3.66 3.73 2.5 (3.42l) 144

7.80i 8.75i 3.70i 3.75i 3.44i 143e, 183f

2.03 2.33 1.46 2.19 7.02 −5.60 −1.80 3753

HSEc (LDAd) PW 3.90

3.91

(2.54) (7.12) (−5.66) (−2.00) 3333j at T = 50 K 5.82k 5.18d

4.64 4.17 [3.87m] 195

expt g

3.47 [3.07m] 192

3.75 [3.25m]n 179o

a a, c, u: lattice parameters. zMull: Mulliken atomic charge. zBorn: Born effective atomic charge (two components for ZnO are separated by /). ε0: static dielectric constant. ε∞: high frequency dielectric constant. ΔEg: the band gap. B: bulk modulus. A, B, O denote A-cation, B-cation and O ion, respectively. Note that the ion site symmetry defines whether several inequivalent Born effective charges exist, i.e. there are two inequivalent directions for oxygen ions in the cubic perovskite structure. bRef 40. cRef 22. dRef 46. eRef 3. fRef 41. gRef 4. hRef 42. iRef 9. jRef 13. kRef 43. lHF exact exchange 0.375 and screening parameter 0.2 from ref 45. mIndirect band gap. nRef 20 and references therein. oRef 44.

C



Article

Table 4. Phonon Frequencies (cm−1) Calculated at the Γ-Point for ZnO and SrTiO3 ZnO

a

SrTiO3 HSEa (LDAb) PW

exptb

symmetry

PBE0 LCAO

e2 (TO1) b1 (silent) a1 (TO2) e1 (TO3) e2 (TO4) b1 (silent)

100 285 407 432 458 562

100 268 400 412 443 570

98 259 378 412 439 552

a1 (LO1) e1 (LO2)

607 619

(564) (556)

574 593

symmetry

PBE0 LCAO

HSEc (LDAd) PW

expte

Transverse (TO) t1u (TO1) t1u (TO2) t2u (silent) t1u (TO3)

25 180 268 547

74i 162 250 533

15 at T = 15 K 175, 170 265 545, 547

Longitudinal (LO) t1u (LO1) t1u (LO2) t1u (LO3)

180 480 824

(158) (454) (829)

171 474 795

Ref 40. bRef 51 measured at 10 K. cRef 22. dRef 46. eRefs 52−54.

scheme.49 This includes the so-called weight factors which turn on the self-interaction correction for the chosen electrons. In our calculations with PBE0 (Table 3) ΔEg = 3.56 and 4.17 eV (indirect gap 3.87 eV) for ZnO and STO, respectively. The calculated value ΔEg for ZnO is only by 0.12 eV larger than the experimental one. A slightly larger difference was observed for STO, viz., 0.4 (0.6) eV for the direct (indirect) gap. Interestingly, the accurate GW calculations for STO overestimate ΔEg in STO,47 suggesting 4.38 (3.82) eV for its direct (indirect) gap. Phonon and Dielectric Properties. As mentioned above, the STO cubic phase represents an example of the ideal perovskitetype cubic ABO3 structure whereas ZnO crystallizes in wurtzite structure. To analyze the symmetry of phonon states, the method of induced representations of space groups was used.20,21 In Table 4 the phonon symmetry and frequencies at the Γpoint of the BZ (except for acoustic modes) are given for the two oxides and compared to the experimental and theoretical data. As is known,20 there are four t1u modes and one t2u (silent) mode in STO. According to its cubic symmetry, three t1u modes are infrared active and one t1u is acoustic. A peculiar feature of STO is the existence of the two low frequency (soft) modes at the Γ- and R-points which are responsible for a trend toward the ferroelectric transition (never achieved due to quantum effects) and the anti-ferrodistortive phase transition observed at 105 K,2,55 respectively. Both modes show strong temperature dependencies which cannot be reproduced in the harmonic approximation. Phonon frequencies calculated with PBE0 at the Γ-point are all real and agree well with the experimental data. This is in contrast to most previous calculations where imaginary soft TO1 frequency was obtained employing either the standard DFT (PBE) or HSE functional.20,22,46,56 This mode is also responsible for the high dielectric constant in STO and known to follow the Cochran law53 ω2 ∝ (T − To) with the extrapolated transition temperature To ≈ 31 K. It is worth mentioning that this frequency is sensitive to the eigenvalue level shift technique commonly used to speed up the convergence in the calculations. Its use increases the soft mode frequency from 25 to 72 cm−1.20 In the present study we do not employ this technique. In a crystal with hexagonal structure like ZnO the cation and anion occupy 2b (1/3, 2/3, z) Wyckoff position. The following set of optical phonon modes at the Γ-point is 2e2 + 2b1 + 2a1 + 2e1 (Table 4). Two phonon modes b1 are silent whereas one pair e1 + a1 is acoustic. The symmetry of the crystal requires

that two phonon modes e1 are Raman active and phonon modes e2 and a1 are both Raman and infrared active. As a general rule, the calculated phonon frequencies are slightly higher than the measured ones.51 The highest discrepancy (the maximum relative error is 7.1%) is observed for low frequency modes (except for the phonon mode e1). A comparison to the results of HSE functional for ZnO (Table 4) shows a good agreement between the two types of calculations. The HSE results are slightly better than the present ones when compared to the measured frequencies except for the phonon mode b1 with highest frequency. We calculated also the longitudinal optical (LO) phonon frequencies for STO and ZnO. The ability to calculate properly the LO−TO splitting of the phonon modes is primarily related to the correct reproduction of the high-frequency dielectric constant ε∞ and the band gap.42 As discussed in ref 42 on ZnO, the underestimated band gap may be a reason of the large discrepancy between the calculated and measured ε∞ and, consequently, the LO−TO splitting. Note that the equilibrium unit cell volume and the Born effective charges can also influence the LO−TO splitting. Our calculated values of ε∞ in ZnO (Table 3) are smaller than the experimental ones by 0.86 and 0.75 for xx- and zz-components, respectively. Similar deviations for ε∞ (1.25 and 0.44, respectively) were obtained in ref 42 using the standard DFT functional. The LO frequencies as calculated with PBE0 (Table 4) reveal a different sequence order than those with the standard DFT functional,51 i.e. the phonon mode a1 is lower than the b1 one. The largest difference between the experiment and the calculated LO frequency with PBE0 is 33 cm−1 for the a1 phonon mode. A similar error (37 cm−1) was observed in the calculation with the standard DFT functional51 for the e1 phonon mode. Both components of the absolute values of the Born effective charges for Zn as well as O atoms are very similar to each other (this was found not only in the present study with PBE0 functional but also in the experiments42 and HSE calculations39). Note that STO, in contrast to ZnO, has a very large value of the static dielectric constant ε09,12−14 though high-frequency dielectric constants ε∞ for these crystals are similar (Table 3). Our calculated value of ε0 = 3753 in STO is close to the experimental one of 333313 at 50 K. The value of ε0 is mainly determined by (and sensitive to) the soft-mode (TO1). B. Oxygen Vacancy in ZnO and SrTiO3. Ground State Properties. The basic properties of a neutral VO were calculated for the defect concentration of 6.25% in STO (Table 5) and 2.78% in ZnO (Table 6). These results include D



Article

Table 5. Structure Relaxation around the Oxygen Vacancy and the Gibbs Free Energy of Formation ΔG0F at 0 K in SrTiO3 for the Defect Concentration of 6.25%a zMull (e)/μ (μB)

distances (Å), displacements in parentheses

ΔG0F (eV)

dV−V (Å)

magnetic configuration

VO−2 × Ti

Ti−4 × Opl

Ti−Oax

VO

2 × Ti

11.08

FM

2.04 (+0.08)

1.85 (−0.11)

0.50/0.4

2.22/0.2

6.26 (6.24)

NM

1.89 (−0.07)

2 × 1.94 (−0.02) 2 × 1.95 (−0.01) 1.94 (−0.02)

2.04 (+0.08)

0.70

2.22

6.59

a

dV−V: the distance between the oxygen vacancies. zMull: Mulliken atomic charge/spin magnetic moment for Ti and vacancy electron occupation due to the ghost basis functions (see the text for details). The values of ΔG0F in parentheses include the zero-point vibrations in the solid phase and oxygen molecule. Positive and negative signs for the displacements denote outward and inward displacement of the ions with respect to the oxygen vacancy site. NM and FM stand for the nonmagnetic and ferromagnetic cases.

moment of VO is 0.4 μB, twice as large as that of each of the two Ti ions (0.2 μB). The FM case is energetically more favorable than the NM and AFM ones, by ∼0.32 and 0.27 eV per supercell, respectively. Consequently, the calculation for the FM state and the ghost basis functions is energetically more preferable than without the ghost, by almost 0.30 eV per supercell. This is why we used the triplet spin state and ghost basis functions in our further study when calculating the phonon contribution to the Gibbs free energy of formation of VO in STO. Note that our results for the magnetic moments of Ti ions and the energy differences between the magnetic states above are similar to GGA+U calculations by Hou and Terakura.17 However, in addition to their findings, our calculations show the defect symmetry reduction from D4h to D2h (due to the Jahn−Teller effect) with almost singly ionized oxygen vacancy (one electron in the vacancy). The four planar O atoms around each of the two Ti centers near VO, being equivalent in the NM solution, split into the two Ti−O pairs in the FM solution and are displaced differently around the vacancy (Table 5). We compared the defect formation energies at 0 K using the 40 atom and 80 atom supercells. For the 40 atom supercell, the calculated ΔG0F = 6.50 eV in STO at 0 K, which is larger by 0.4 eV than ∼6.1 eV extracted from the experimental temperature dependence extrapolated toward low temperatures.19 In contrast, for the 80 atom supercell, ΔG0F = 6.26 eV, which is very close to the experimental value. This demonstrates that the employment of the 80 atom supercell is necessary to obtain the convergence for the defect formation energy at 0 K caused by a strong electron density delocalization around VO in STO. Interestingly, the contribution of zero-point vibrations to the formation energy in STO is very small, 0.02 eV (Table 5). The analysis of band structure for defective STO reveals (Figure 2a) that due to the Jahn−Teller effect the defect energy splits into two states, each being occupied by one electron (both with parallel spins). As Figure 2a shows, one defect level at the Γ-point lies ∼0.5 eV, whereas another one lies ∼1.0 eV, below the conduction band bottom, suggesting a strong electron density delocalization. A qualitatively similar pattern at the Γ-point has been discussed by Hou and Terakura.17 Note, however, that our band structure calculations show that the upper defect energy band has a large dispersion and falls into the bottom of the conduction band at other symmetry points of the BZ, which could be a consequence of the finite supercell size. Each oxygen atom in a perfect wurtzite structure of ZnO is surrounded by three Zn ions at the same distance, one Zn ion at slightly larger distance and 6 O ions at almost twice larger distances. The Zn and O ions are moved by 0.16 Å and 0.04− 0.07 Å inward the vacancy direction in the defective ZnO,

Table 6. The Same as Table 5 for ZnO for the Defect Concentration of 2.78% distances (Å), displacements in parentheses

zMull (e)

dV−V (Å)

VO−Zn

VO−O

VO

Zn

ΔG0F (eV)

9.78 10.44

1.83 (−0.16) 1.84 (−0.16)

3.15 (−0.07) 3.18 (−0.04)

0.83

0.91

4.20 (4.09)

among other properties for STO the distances between VO and two nearest Ti ions (VO − 2 × Ti). The five remaining O ions around these Ti ions are classified in Table 5 and Figure 1a as

Figure 1. The atomic structure of STO (a) and ZnO (b) with VO. The arrows show the relaxation pattern around VO.

four planar ions (Ti − 4 × Opl) and one axial ion (Ti − Oax). As pointed out above, VO in SrTiO3 has two electrons due to the missing oxygen ion being partly localized on nearest neighbor Ti ions. In such a case the open 3d electron shells of Ti ions suggest a spin-restricted solution for STO with VO. We considered a triplet (ferromagnetic, FM) state with parallel spin alignment on the two Ti ions, a singlet (antiferromagnetic, AFM) state with antiparallel spin alignment and a nonmagnetic (NM) solution. Interestingly, there is a considerable electronic charge concentrated within the vacancy due to the use of the ghost basis functions (Table 5). Moreover, the magnetic E



Article

0.11 eV) than in STO. In contrast to STO, ZnO demonstrates (i) the defect level lying very close to the top of valence band (Figure 2b), (ii) weak electron density delocalization and (iii) the valence band is mainly formed by the states of Zn whereas the O states form a separate band below the Zn states (Figure 2b) (even though the presence of VO leads to their enhanced hybridization). C. Phonon Contribution to the Formation Energy. In the following we consider the phonon contribution to the free defect formation energy in the two oxides with very different static dielectric constant ε0 and degree of the electron density delocalization around VO. Obviously, periodically distributed oxygen vacancies at concentrations as large as 12.5% in STO interact with each other, e.g. the difference in the formation energies for the defect concentrations of 12.5% and 6.25% is quite pronounced, namely, ∼0.18 eV at 0 K. The calculation of the Gibbs free energy of defect formation requires integration over the phonon spectrum of perfect and defective solids. As discussed above, bulk SrTiO3 is an incipient ferroelectric material and is characterized by the presence of two soft modes with the temperature dependence frequently discussed.52 The phonon mode frequencies above 400 cm−1 decrease very smoothly with increasing temperature, whereas the two soft modes mentioned above show strong temperature dependencies. Since this effect is beyond the harmonic approximation used here, the temperature dependence of the soft modes has been calculated within the polarizability model.59,60 Its results have in turn been used as an input for the evaluation of the defect formation energies. In Figure 3 the temperature

Figure 2. The band structure and density of states (DOS) with respect to the Fermi level EF of defective STO (a) and ZnO (b). The concentration of VO was 2.78% and 6.25% for ZnO and STO, respectively. The projected DOS of VO is multiplied by a factor of 50. Solid and dashed lines denote α and β electrons, respectively, whereas EF denotes the Fermi energy.

Figure 3. The calculated high temperature dependence of soft phonon modes at both Γ- and R-points of the BZ (filled squares and circles) and their comparison to the experimental data. Filled squares correspond to calculated values by the polarizability model at the Γpoint, filled circles calculated values at the R-point, open squares experimentally determined (ref 53) soft mode at the Γ-point according to Cochran law up to 300 K, open circles experimentally determined soft mode at the R-point at high temperatures extrapolated to higher temperatures using eq 1 of ref 61, triangles experimentally determined soft mode at the Γ-point at high temperatures.52.

respectively. As in the case of STO, the oxygen vacancy in ZnO shows a considerable Mulliken charge of 0.83 e confirming the important role of the ghost basis functions in the studies of oxygen vacancies. Our calculated value of ΔG0F (4.20 eV) obtained with PBE0 at T = 0 K is close to that obtained with HSE functional by Oba et al.45,57(∼4.1 eV), (SIC)-LDA study of Koerner and Elsaesser49 (4.49 eV) and the screened exchange functional by Clark et al.58 (∼4.2 eV). The improved DFT+U study of Boonchum and Lambrecht50 suggested a value which is by almost 1.0 eV smaller, i.e. 3.32 eV. We cannot exclude here the impact of the O2 binding energy obtained within the standard DFT. We found that the effect of zeropoint vibrations for ΔG0F is stronger in ZnO (ΔG0F changes by

dependence of soft phonon modes is shown. Obviously, the polarizability model results are in very good agreement with the experimental data. In contrast, the mean-field Cochran law prediction53 shows a considerable deviation above 300 K from them.52,60,61 Nevertheless, we plotted the Cochran law in Figure 3 for a comparison to the results of the polarizability model. Since our calculations of ΔG0F address high temperF



Article

atures, the deviations from the mean field behavior are important and have to be taken into account (Figure 3). The soft mode at the R-point saturates faster with temperature than that at the Γ-point, which is also in agreement with the experimental observations.61 Thus, ΔG0F as calculated using eq 1 for STO with the defect concentration of 12.5% and the temperature dependent soft modes taken from the polarizability model decreases by ∼1 eV within a broad temperature range, from 300 to 1660 K (Figure 4a). However, there is a

at high temperatures62 is also shown. The agreement between the calculated ΔGF0 for the 80 atom supercell and the experiment is quite good. Note that the experimental curves refer to electron-free vacancies. Consideration of the small energetic effect of the electron trapping (0.1−0.3 eV62,64) would slightly improve the agreement. ΔG0F is also plotted in Figure 4b for the two concentrations of vacancies in ZnO with and without the phonon contribution in the solid phase. In the case of phonon contribution in the solid phase included ΔG0F decreased by ∼1.3 eV within the same temperature range as used for STO with almost no difference between the two concentrations (see also refs 45 and 63). Thus, the oxygen vacancies in ZnO do not interact strongly (due to much smaller dielectric constant of ZnO and weak electron density delocalization around VO as compared to STO). The temperature dependence of ΔG0F is mainly caused by the variation of the oxygen chemical potential.

4. CONCLUSIONS By means of careful first-principles calculations we have shown that the hybrid PBE0 exchange-correlation functional (with the standard HF exact exchange of 0.25) and the LCAO basis set can correctly reproduce the electronic as well as phonon properties of both ZnO and SrTiO3 without use of any additional parameters. The calculated Raman- and infraredactive phonon frequencies in defect-free hexagonal ZnO and infrared-active phonon frequencies in a cubic SrTiO3 are in very good agreement with experimental results. In defective SrTiO3 the oxygen vacancy induces the Jahn−Teller lattice distortion and gives rise to the triplet spin ground state. In particular, we have shown that the additional “ghost” basis set at the vacancy site considerably affects the defect formation energies and leads in both oxides to almost singly ionized oxygen vacancies (one electron trapped in/around the vacancy). The corresponding defect formation energies at 0 K (6.26 eV in SrTiO3 and 4.20 eV in ZnO) are well converged with respect to the supercell size and consistent with the experiments. The temperature dependences of the Gibbs free energies of formation of oxygen vacancy in ZnO and SrTiO3 have been calculated by taking into account the phonon contributions. The temperature dependence of the soft modes as calculated for a defect-free SrTiO3 within the polarizability model is of key importance for the correct calculation of defect formation energy. SrTiO 3 evidences supercell finite-size effects on the Gibbs free energy of defect formation, in contrast to ZnO, which is caused by a considerable delocalization of electronic density from oxygen vacancy in SrTiO3 as well as inefficient number of phonon frequencies in calculating the Gibbs free energy of formation using standard 40 atom supercell.

Figure 4. The calculated Gibbs free energy (per vacancy) of formation of oxygen vacancy ΔG0F with and without phonon contribution in the solid phase as a function of temperature in STO (a) and ZnO (b). The defect concentrations were 2.78 (6.25)% and 12.5% (6.25%) for ZnO and STO, respectively. The distances between vacancies in STO are 7.83 and 11.08 Å for 12.5% and 6.25%, respectively. The distances between vacancies in ZnO are 6.52 (10.44) Å and 9.78 (10.44) Å for 6.25% and 2.78%, respectively. Note that the experimental curves refer to electron-free vacancies. Consideration of the small energetic effect of the electron trapping (0.1−0.3 eV62,64) would slightly improve the agreement.

■

clear difference between this dependence and that neglecting the phonon contribution in the solid phase. The latter case was calculated using only the chemical potential of the oxygen molecule. The energy difference between these two curves increases considerably with temperature for the 40 atom supercells but becomes much smaller for the larger, the 80 atom supercell (dash−dotted curve in Figure 4a). These results clearly evidence supercell finite-size effects, i.e. slow convergence of the Gibbs free energy of defect formation ΔG0F in perovskite materials as STO with respect to the supercell size (related to the defect concentration). In Figure 4a ΔGF0 extracted from the experimental conductivity measurements

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. Mailing address: Max Planck Institute for Solid State Research, Heisenbergstr. 1, Stuttgart, 70569, Germany. Phone: +49-711-6891771. Notes

The authors declare no competing financial interest.

■

ACKNOWLEDGMENTS This study was partly supported by ERAF 2010/0272/2DP/ 2.1.1.1.0/10/APIA/VIAA/088 project, St. Petersburg State University Grant 12.0.108.2010 and Juelich supercomputer G



Article

(21) Gryaznov, D.; Evarestov, R. A.; Maier, J. Hybrid DensityFunctional Calculations of Phonons in LaCoO3. Phys. Rev. B 2010, 82, 224301. (22) Wahl, R.; Vogtenhuber, D.; Kresse, G. SrTiO3 and BaTiO3 Revisited Using the Projector Augmented Wave Method: Performance of Hybrid and Semilocal Functionals. Phys. Rev. B 2008, 78, 104116. (23) Jaffe, J. E.; Hess, A. C. Hartree-Fock Study of Phase Changes in ZnO at High Pressure. Phys. Rev. B 1993, 48 (11), 7903. (24) Perdew, J. P.; Ernzerhof, M.; Burke, K. Rationale for Mixing Exact Exchange with Density Functional Approximations. J. Chem. Phys. 1996, 105, 9982−9985. (25) Dovesi, R.; Saunders, V. R.; Roetti, C.; Orlando, R.; ZicovichWilson, C. M.; Pascale, F.; Civalleri, B.; Doll, K.; Harrison, N. M.; Bush, I. J.; er al. CRYSTAL09 User’s Manual; University of Torino: 2009. (26) Evarestov, R. A.; Bandura, A. First-Principles Calculations on the Four Phases of BaTiO3. J. Comput. Chem. 2012, 33, 1123−1130. (27) Alkauskas, A.; Broqvist, P.; Pasquarello, A. Defect Energy Levels in Density Functional Calculations: Alignment and Band Gap Problem. Phys. Rev. Lett. 2008, 101, 046405. (28) Catti, M.; Vallerio, G.; Dovesi, R. Theoretical Study of Electronic, Magnetic, and Structural Properties of a-Fe2O3 (Hematite). Phys. Rev. B 1995, 51, 7441. (29) Monkhorst, J.; Pack, J. D. Special Points for Brillonin-Zone Integrations. Phys. Rev. B 1976, 13, 5188−5192. (30) Ferrero, M.; Rérat, M.; Orlando, R.; Dovesi, R. The Calculation of Static Polarizabilities of 1−3D Periodic Compounds. The Implementation in the CRYSTAL Code. J. Comput. Chem. 2008, 29, 1450−1459. (31) Evarestov, R. A. Quantum Chemistry of Solids; Springer: Heidelberg, 2013. (32) Heyd, J. Screened Coulomb Hybrid Density Functionals; PhD Thesis, Rice University, 2004. (33) Bokstein, B.,; Mendelev, M. I.; Srolovitz, D. J. Thermodynamics and Kinetics in Materials Science: A Short Course; Oxford Univ. Press: New York, 2005. (34) Kotomin, E. A.; Mastrikov, Yu. A.; Heifets, E.; Maier, J. Adsorption of Atomic and Molecular Oxygen on the LaMnO3(001) Surface: Ab Initio Supercell Calculations and Thermodynamics. Phys. Chem. Chem. Phys. 2008, 10, 4644−4649. (35) Bocharov, D.; Gryaznov, D.; Zhukovskii, Yu. F.; Kotomin, E. A. Ab Initio Modeling of Oxygen Impurity Atom Incorporation into Uranium Mononitride Surface and Sub-surface Vacancies. J. Nucl. Mater. 2011, 416, 200−204. (36) Lee, Y.-L.; Kleis, J.; Rossmeisl, J.; Morgan, D. Ab Initio Energetics of LaBO3(001) (B = Mn, Fe, Co, and Ni) for Solid Oxide Fuel Cells. Phys. Rev. B 2009, 80, 224101. (37) Dorado, B.; Freyss, M.; Martin, G. GGA+U Study of the Incorporation of Iodine in Uranium Dioxide. Eur. Phys. J. B 2009, 69, 203−209. (38) Stoffel, R. P.; Wessel, C.; Lumey, M.; Dronkowski, R. Ab Initio Thermochemistry of Solid-State Materials. Angew. Chem., Int. Ed. 2010, 49, 5242−5266. (39) Grimvall, G. Thermophysical Properties of Materials; NorthHolland: Amsterdam, 1986. (40) Wróbel, J.; Kurzydłowski, K.; Hummer, K.; Kresse, G.; Piechota, J. Calculations of ZnO Properties Using the Heyd-Scuseria-Ernzerhof Screened Hybrid Density Functional. Phys. Rev. B 2009, 80, 155124. (41) Desgreniers, S. High-Density Phases of ZnO: Structural and Compressive Parameters. Phys. Rev. B 1998, 58, 14102. (42) Serrano, J.; Romero, A. H.; Manjón, F. J.; Lauck, R.; Cardona, M.; Rubio, A. Pressure Dependence of the Lattice Dynamics of ZnO: An Ab Initio Approach. Phys. Rev. B 2004, 69, 094306. (43) Dore, P.; Paolone, A.; Trippetti, R. Mid-Infrared and NearInfrared Properties of SrTiO3 from Transmission and Reflection Measurements. J. Appl. Phys. 1996, 80 (9), 5270−5274. (44) Ferroelectrics and Related Substances; Hellwege, K. H., Hellwege, A. M., Eds.; Landolt-Börnstein, New Series, Group III; SpringerVerlag: Berlin, 1969; Vol. 3.

center (project HSS13.) The authors are grateful to Manuel Cardona, Mike Finnis, Juris Purans, Vladimir Trepakov, Rotraut Merkle and Eugene Heifets for numerous stimulating discussions.

■

REFERENCES

(1) Schmidt-Mende, L.; MacManus-Driscoll, J. L. ZnO-Nanosructures, Defects, and Devices. Mater. Today 2010, 10 (5), 40−48. (2) Lines, M. E.; Glass, A. M. Principles and Applications of Ferroelectrics and Related Materials; Oxford Univ. Press: New York, 1977. (3) Karzel, H.; Potzel, W.; Köfferlein, M.; Schiessl, W.; Steiner, M.; Hiller, U.; Kalvius, G. M.; Mitchell, D. W.; Das, T. P.; Blaha, P.; et al. Lattice Dynamics and Hyperfine Interactions in ZnO and ZnSe at High External Pressures. Phys. Rev. B 1996, 53, 11425−11438. (4) Cao, L.; Sozontov, E.; Zegenhagen, J. Cubic to Tetragonal Phase Transition of SrTiO3 Under Epitaxial Stress: An X-Ray Backscattering Study. Phys. Status Solidi A 2000, 181, 387−404. (5) Merkle, R.; Maier, J. How Is Oxygen Incorporated into Oxides? A Comprehensive Kineitc Study of a Simple Solid-State Reaction with SrTiO3 as a Model Material. Angew. Chem., Int. Ed. 2008, 47, 3874− 3894. (6) Sunarso, J.; Baumann, S.; Serra, J. M.; Meulenberg, W. A.; Liu, S.; Lin, Y. S.; Diniz da Costa, J. C. Mixed Ionic-electronic Conducting (MIEC) Ceramic-Based Membranes for Oxygen Separation. J. Membr. Sci. 2008, 320, 13−41. (7) Jannoti, A.; Van de Walle, C. G. Fundamentals of Zinc Oxide as a Semiconductor. Rep. Prog. Phys. 2009, 72, 126501. (8) Jung, E. S.; Kim, H. S.; Kong, B. H.; Cho, H. K.; Park, N. K.; Lee, H. S. The Effect of Oxygen Content on the Electric Charecteristics of ZnO. Phys. Staus. Solidi B 2007, 244 (5), 1553−1559. (9) Springer Materials, The Landolt-Börnstein Database (2012). (10) Müller, K. A.; Burkhard, H. SrTiO3: An Intrinsic Quantum Paraelectric below 4 K. Phys. Rev. B 1979, 19, 3593−3602. (11) Burkhard, H.; Müller, K. A. Evidenz für die Quantenmechanische Unterdrückung der Ferroelektrischen Phasenumwandlung in SrTi03 aus Suszeptibilitäts- und Elektrostriktionsmessungen. Helv. Phys. Acta 1976, 49, 725. (12) Li, H.-C.; Si, W.; West, A. D.; Xi, X. X. Thickness Dependence of Dielectric Loss in SrTiO3 Thin Films. Appl. Phys. Lett. 1998, 73 (4), 464−466. (13) Sirenko, A. A.; Bernhard, C.; Golnik, A.; Clark, A. M.; Hao, J.; Si, W.; Xi, X. X. Soft-Mode Hardening in SrTiO3 Thin Films. Nature 2000, 404, 373−376. (14) Liu, X. Z.; Li, Y. R. Dielectric Properties of Multilayered SrTiO3 Thin Films with Graded Oxygen Vacancy Concentration. Appl. Phys. A: Mater. Sci. Process. 2006, 83, 67−72. (15) Zhukovskii, Yu. F.; Kotomin, E. A.; Evarestov, R. A.; Ellis, D. E. Periodic Models in Quantum Chemical Simulations of F Centers in Crystalline Metal Oxides. Int. J. Quantum Chem. 2007, 107, 2956− 2985. (16) Cuong, D. D.; Lee, B.; Choi, K.; Ahn, H.; Han, S.; Lee, J. Oxygen Vacancy Clustering and Electron Localization in OxygenDeficient SrTiO3: LDA+U Study. Phys. Rev. Lett. 2007, 98, 115503. (17) Hou, Z.; Terakura, K. Defect States Induced by Oxygen Vacancies in Cubic SrTiO3: First-Principles Calculations. J. Phys. Soc. Jpn. 2010, 79 (11), 114704. (18) Alexandrov, V. E.; Kotomin, E. A.; Maier, J.; Evarestov, R. A. First-Principles Study of Bulk and Surface Oxygen Vacancies in SrTiO3 Crystal. Eur. Phys. J. B 2009, 72, 53−57. (19) Evarestov, R. A.; Blokhin, E.; Gryaznov, D.; Kotomin, E. A.; Merkle, R.; Maier, J. Jahn-Teller Effect in the Phonon Properties of Defective SrTiO3 from First Principes. Phys. Rev. B 2012, 85, 174303. (20) Evarestov, R. A.; Blokhin, E.; Gryaznov, D.; Kotomin, E. A.; Maier, J. Phonon Calculations in Cubic and Tetragonal Phases of SrTiO3. Phys. Rev. B 2011, 83, 134108. H



Article

(45) Oba, F.; Choi, M.; Togo, A.; Seko, A.; Tanaka, I. Native Defects in Oxide Semiconductors: A Density Functional Approach. J. Phys.: Condens. Matter 2010, 22, 384211. (46) Zhong, W.; King-Smith, R. D; Vanderbilt, D. Giant LO-TO Splitting in Perovskite Ferroelectrics. Phys. Rev. Lett. 1994, 72 (22), 3618−3621. (47) Benrekia, A. R.; Benkhettou, N.; Nassour, A.; Driz, M.; Sahnoun, M.; Lebégue, S. Structural, Electronic and Optical Properties of Cubic SrTiO3 and KTaO3: Ab Initio and GW Calculations. Physica B 2012, 407, 2632−2636. (48) Lany, S.; Zunger, A. Accurate Prediction of Defect Properties in Density Functional Supercell Calculations. Modell. Simul. Mater. Sci. Eng. 2009, 17, 084002. (49) Körner, W.; Elsässer, C. First-Principles Density Functional Study of Dopant Elements at Grain Boundaries in ZnO. Phys. Rev. B 2010, 81, 085324. (50) Boonchun, A.; Lambrecht, W. R. L. Critical Evaluation of the LDA+U Approach for Band Gap corrections in Point Defect Calculations: The Oxygen Vacancy in ZnO Case Study. Phys. Status Solidi B 2011, 248 (5), 1043−1051. (51) Serrano, J.; Manjón, F. J.; Romero, A. H.; Ivanov, A.; Cardona, M.; Lauck, R.; Bosak, A.; Krisch, M. Phonon Dispersion Relations of Zinc Oxide: Inelastic Neutron Scattering and Ab Initio Calculations. Phys. Rev. B 2010, 81, 174304. (52) Servoin, J. L.; Luspin, Y.; Gervais, F. Infrared Dispersion in SrTiO3 at High Temperatures. Phys. Rev. B 1980, 22, 5501−5506. (53) Petzelt, J.; Ostapchuk, T.; Gregora, I.; Rychetský, I.; HoffmannEifert, S.; Pronin, A. V.; Yuzyuk, Y.; Gorshunov, B. P.; Kamba, S.; Bovtun, V.; et al. Dielectric, Infrared, and Raman Response of Undoped SrTiO3 Ceramics: Evidence of Polar Grain Boundaries. Phys. Rev. B 2001, 64, 184111. (54) Stirling, W. G. Neutron Inelastic Scattering Study of the Lattice Dynamics of Strontium Titanate: Harmonic Models. J. Phys. C 1972, 5, 2711−2730. (55) Yuzyuk, Yu. I. Raman Scattering Spectra of Ceramics, Films, and Superlattices of Ferrelecrtic Perovskites: A Review. Phys. Solid State 2012, 54 (5), 1026−1059. (56) Lebedev, A. I. Ab Initio Calculations of Phonon Spectra in ATiO3 Perovskite Crystals (A = Ca, Sr, Ba, Ra, Cd, Zn, Mg, Ge, Sn, Pb). Phys. Solid State 2009, 51 (2), 362−372. (57) Oba, F.; Togo, A.; Tanaka, I.; Paier, J.; Kresse, G. Defect Energetics in ZnO: A Hybrid Hartree-Fock Density Functional Study. Phys. Rev. B 2008, 77, 245202. (58) Clark, S. J.; Robertson, J.; Lany, S.; Zunger, A. Intrinsic Defects in ZnO Calculated by Screened Exchange and Hybrid Density Functionals. Phys. Rev. B 2010, 81, 115311. (59) Bussmann-Holder, A. The Polarizability Model for Ferroelectricity in Perovskite oxides. J. Phys.: Condens. Matter. 2012, 24, 273202. (60) Bussmann-Holder, A.; Bilz, H.; Benedek, G. Applications of the Polarizability Model to Various Displacive-type Ferrelectric Systems. Phys. Rev. B 1989, 39, 9214−9223. (61) Shirane, G.; Yamada, Y. Lattice-Dynamical Study of the 110 K Phase Transition in SrTiO3. Phys. Rev. 1969, 177, 858−863. (62) Moos, R.; Haerdtl, K. H. Defect Chemistry of Donor-Doped and Undoped Strontium Titanate Ceramics between 1000° and 1400 °C. J. Am. Ceram. Soc. 1997, 80, 2549−2562. (63) Vidya, R.; Ravindran, P.; Fjellvåg, H.; Svensson, B. G.; Monakhov, E.; Ganchenkova, M.; Nieminen, R. M. Energetics of Intrinsic Defects and Their Complexes in ZnO Investigated by Density Functional Theory. Phys. Rev. B 2011, 83, 045206. (64) Maier, J. Physical Chemistry of Ionic Materials: Ions and Electrons in Solids; Wiley: Chichester, England, 2004.

I


A Comparative Ab Initio Thermodynamic Study of Oxygen Vacancies

Recommend Documents