Subscriber access provided by Macquarie University
C: Surfaces, Interfaces, Porous Materials, and Catalysis 2
3
Temperature Dependence of Self-Diffusion in CrO from First Principles Bharat K Medasani, Maria L. Sushko, Kevin M. Rosso, Daniel K. Schreiber, and Stephen M. Bruemmer J. Phys. Chem. C, Just Accepted Manuscript • DOI: 10.1021/acs.jpcc.9b03218 • Publication Date (Web): 26 Aug 2019 Downloaded from pubs.acs.org on August 29, 2019
Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.
is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.
Page 1 of 38 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
The Journal of Physical Chemistry
Temperature Dependence of Self-diffusion in Cr2O3 from First Principles Bharat K. Medasani,∗,†,‡ Maria L. Sushko,∗,‡ Kevin M. Rosso,‡ Daniel K. Schreiber,¶ and Stephen M. Bruemmer¶ †Delaware Energy Institute, University of Delaware, Newark, DE 19702, USA ‡Physical and Computational Sciences Directorate, Pacific Northwest National Laboratory, Richland, WA 99354, USA ¶Energy and Environment Directorate, Pacific Northwest National Laboratory, Richland, WA 99354, USA E-mail:
[email protected],
[email protected];
[email protected] Abstract Understanding and predicting the dominant diffusion processes in Cr2 O3 is essential to its optimization for anti-corrosion coatings, spintronics, and other applications. Despite significant theoretical effort in modeling defect mediated diffusion in Cr2 O3 the correlation with experimentally measured diffusivities remains poor partly due to the insufficient accuracy of the theoretical approaches. Here an attempt to resolve these discrepancies is made through high accuracy density functional theory simulations coupled with grand canonical formalism of defect thermochemistry. In this approach, point defect formation energies were computed using hybrid exchange correlation functional. This level of theory proved to be essential for achieving the agreement with experimental self-diffusion coefficients. The analysis of the resulting self-diffusion coefficients indicate that chromium has higher mobility at low temperatures and high oxygen
1
ACS Paragon Plus Environment
The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
partial pressures, in particular at standard temperature and pressure conditions. At high vacuum, high temperature conditions, oxygen diffusion becomes dominant. At Cr/Cr2 O3 interfaces O vacancies were found to be more mobile than Cr vacancies at all temperatures. Cr diffuses preferentially along the c-axis at low temperatures but switches to basal plane at higher temperatures. O diffusion is primarily bound to basal plane at all temperatures.
Introduction Cr2 O3 finds applications in both functional and structural applications such as catalysis, 1,2 electro-optics, 3,4 magnetoelectronic, 5,6 and corrosion protection. 7 It serves as a protective oxide coating on widely used high-temperature Ni and Fe alloys. 8–10 Point defects play an important role in the effectiveness of Cr2 O3 for corrosion protection and other applications and significant experimental effort has been devoted to understanding the nature of point defects and their diffusion mechanisms. However, there is a large spread in the reported selfdiffusion coefficients and defect formation energies from various experimental studies. 11–20 Complementing the experimental investigations, diffusion mediated by both intrinsic and extrinsic point defects in Cr2 O3 have been extensively studied theoretically by many groups, including ours. 21–28 Due to the advances in density functional theory (DFT) accuracy in calculating many properties of materials, recent studies employed DFT to evaluate formation energies, migration barriers and charge localization of point defects. 22–24,27,28 The 0K energies obtained from DFT were combined with Einstein-Smoluchowski random walk thermal diffusion formalism 29 to predict the self-diffusion coefficients in Cr2 O3 . However, the agreement between calculated and experimental self-diffusion coefficients had, so far, been poor, which prompted the use of modified diffusion equations. 22,27 This leaves the question of the root cause of the discrepancy between theory and experiment unresolved. A possible source of discrepancy is the insufficient accuracy of the DFT in modeling transition metal oxides, which requires the use of empirical corrections to the exchange2
ACS Paragon Plus Environment
Page 2 of 38
Page 3 of 38 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
The Journal of Physical Chemistry
correlation functionals to induce electron localization. For example, to correct for the selfinteraction error and induce electron localization in modeling Cr2 O3 highly correlated Cr3d (and sometimes O-2p) electrons the Hubbard on-site correction is often employed. 30–32 However, the resulting defect formation energies vary significantly 20,22–25,33 depending on the choice of Hubbard U parameters and the type of the semi-local generalized gradient approximation (GGA) functionals. Adding to the confusion, corrections to minimize the electrostatic errors inherent in the periodic supercell formulation of defects and the errors due to the under-prediction of bandgap by semi-local functionals were either not applied or inconsistently applied in many of these studies. Rak and Brenner found that the barrier energies for defect migration were not sensitive to the variations in U value. In contrast, they showed that varying U value results in a variation of defect formation energies. 27 The likely reason for such variations is in electronic dissimilarities in the local environments in defected and bulk systems. For example, the electronic density surrounding a point defect such as a vacancy is much lower compared to the electronic density in the corresponding bulk system and the relative effect of U is different. Selecting the U parameter in the Hubbard correction scheme is neither trivial nor consistent. Often, U selection is based on fitting a bulk property to the experimental value. The chosen optimal U value is highly dependent on the selected bulk property. Another approach which avoids the empirical route in the selection of U is based on the linear response theory (LRT). 34 LRT based U values often tend to be smaller than the U values obtained from regressing bulk properties. To overcome the uncertainties associated with the empirical U value of Hubbard correction scheme, we employed hybrid functional to accurately evaluate defect electronic levels and defect formation energies. Here we report self-diffusion coefficients computed using Einstein-Smoluchowski formalism and high accuracy point defect energies. The resulting diffusion coefficients have a reasonable agreement with the experimental data. This work can be considered as the cul-
3
ACS Paragon Plus Environment
The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
mination of our previous studies on the atomistic diffusion mechanisms in Cr2 O3 mediated by vacancies and interstitials. 23,28 For barrier energies, the transition state and ground defect state exhibit similar local electronic environment, suggesting an insensitivity to the U value used in this prior work. Hence, by utilizing the transition barrier energies obtained in those studies in addition to the revised more accurate defect energies and defect thermodynamics, we evaluated Fermi level variations, stoichiometry deviations and self diffusion coefficients in undoped bulk Cr2 O3 as a function of temperature and oxygen partial pressure. The computed self diffusion coefficients were analyzed using Brouwer diagrams to reveal the dominant defects responsible for diffusion.
Methods Density Functional Theory Approach Density functional theory (DFT) calculations for defect formation energies were performed using the Vienna ab initio simulation package (VASP). 35–37 We utilized HSE 38 hybrid functional with 25% mixing of Hartree - Fock exchange at short range and 100% PBE correlation. Projector augmented wave (PAW) 39,40 optimized for PBE 41 functional was used. Wave functions preoptimized from GGA+U calculations with PBEsol 42 functional was used as input to the HSE calculations. The HSE screening parameter of 0.6, which is different from the commonly used values of 0.2 and 0.3 corresponding to HSE06 and HSE03, respectively, was chosen after fitting the computed bandgaps with the experimental bandgap. This point is further discussed in the Results section. The reciprocal space of the primitive cell of Cr2 O3 (R3c space group) with 10 atoms was sampled with 5×5×5 Γ-centered k -point grid. The primitive cell was fully relaxed (both size and shape were relaxed) until the forces converged to 0.01 eV/˚ A. The atomic positions in the 2×2×1 defect supercells were relaxed at constant volume and fixed cell shape until the individual forces on each atom were minimized to 0.03 eV/˚ A. For supercell calculations, a cut-off value of 400 eV was used for the plane-wave basis 4
ACS Paragon Plus Environment
Page 4 of 38
Page 5 of 38 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
The Journal of Physical Chemistry
set and 2×2×1 Γ-centered k -point grid was used to sample the reciprocal space. Spin polarization with anti-ferromagnetic (AFM) ordering of Cr spins was used. Gaussian method with a width of 0.01 eV was used for electronic smearing.
Diffusion Coefficient Calculation In materials which support charged defects, the diffusion coefficient of an element s can be defined as Ds =
XX
d(Xsq )
(1)
q
X
where d is the diffusion coefficient pertaining to a defect Xsq , with X denoting the defect type, and q representing the defect charge. X can be a regular point defect, such as a vacancy or an interstitial, or it can be a complex defect, such as the Frenkel defect. Based on Einstein’s random walk theory, d can be computed as 1 X mp lp2 Γp . d= c 2 p
(2)
Here c is the concentration of defect, Xsq , and p represents one of the migration pathways. mp , lp , and Γp indicate the multiplicity, length, and jump frequency associated with path p, respectively. The defect concentration, c, is a function of the defect formation energy Ef , c = c0 e−βEf ,
(3)
were c0 is the number of lattice sites per cell volume. According to Vineyard’s transition state theory (TST), 43 the jump frequency is defined as a function of the attempt frequency νp and the migration barrier energy Em,p as Γp = νp e−βEm,p .
5
ACS Paragon Plus Environment
(4)
The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 6 of 38
The attempt frequency is evaluated under harmonic approximation using the phonon modes of defect ground and transition states as Q3N −3
i ν = Q3N −4 i
νi νj0
,
(5)
where νi and νj0 represent the real phonon frequencies of the ground and transition states of the defect respectively and N is the number of atoms in the bulk supercell.
Finite Temperature Defect Concentrations The effect of pressure and temperature is accounted for in the diffusion coefficients through the chemical potential term in the expression for defect formation energies. Assuming dilute concentrations for defects, we utilized constrained grand canonical formalism to compute the elemental chemical potentials at finite temperatures. This formalism has been successfully applied to predict defect concentrations in ZrO2 under dilute conditions. 44 In this formalism, bulk Cr2 O3 is assumed to be in contact with an infinite O2 reservoir at temperature T and partial pressure pO2 . The oxygen chemical potential then is defined as µO (pO2 , T ) = µ0O + ∆µ0O (1atm, 298K) + ∆µO (pO2 , T ).
(6)
In the above equation, µ0O is the 0 K DFT computed O chemical potential, ∆µ0O (1 atm, 298 K) is the correction to obtain the experimental chemical potential of O at standard temperature and pressure (STP: 273 K and 1 atm pressure) conditions, and ∆µO (pO2 , T ) is given by
∆µO (pO2 , T ) = kT ln(pO2 ).
6
ACS Paragon Plus Environment
(7)
Page 7 of 38 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
The Journal of Physical Chemistry
The latter two quantities are obtained from NIST-JANAF tables. 45 Cr chemical potential is obtained from µO as µCr (pO2 , T ) = ∆HCr2 O3 − µO (pO2 , T ),
(8)
where ∆HCr2 O3 is the formation enthalpy of Cr2 O3 . Since the overall system is charge neutral, we impose the charge neutrality condition after accounting for the charged defects and free carriers such as electrons and holes using the expression q
X
q cXs + p − n = 0,
(9)
Xsq
where p and n represent the number densities of holes and electrons, respectively. Given that the formation energy of a charged defect can be calculated as 23,28 q
q
Xs bulk EfXs = Etot − Etot +
X
q
Xs ns0 (Xsq )µs0 + qEF + Ecorr ,
(10)
s0
we can rewrite the concentration of defect Xsq as q
q
q Xs (p
s −Ef cXs (pO2 , T, EF ) = cX 0 e
O2 ,T,EF )/kT
.
(11)
fF D (E, T, EF )g(E)dE,
(12)
The number density of free electrons, n, is given by Z
∞
n(T, EF ) = CBM
where g(E) is the density of states, fF D is the Fermi-Dirac distribution function given by
fF D (E, T, EF ) =
7
1 1+
e(E−EF )/kT
ACS Paragon Plus Environment
.
(13)
The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 8 of 38
Similarly the number of free holes, p, is given by Z
V BM
(1 − fF D (E, T, EF ))g(E)dE.
p(T, EF ) =
(14)
−∞
Substituting Eqns. 11, 12, and 14 into Eqn. 9, and by making T and pO2 as free variables, we can solve for EF in terms of T and pO2 self-consistently. The resulting defect concentrations are given by q
q
q Xs (p
s −Ef cXs (pO2 , T ) = cX 0 e
O2 ,T )/kT
.
(15)
Results and Discussion HSE Screening Parameter To eliminate the uncertainty surrounding the defect formation energies computed with GGA+U formalism, we utilized HSE hybrid functional. The bandgap of Cr2 O3 predictions by the two widely used flavors of HSE, of HSE06 and HSE03, 38 hybrid functional, of 4.27 eV and 3.85 eV, respectively, were much higher than the experimental values of 3.2-3.4 eV. In order to lower the calculated bandgap closer to the experimental value of 3.4 eV, we modified the parameters of the HSE functional. We retained the 25% mixing of Hartree-Fock exchange at short range and increased the screening parameter till the computed bandgap matched with the experimental value. For bandgap fitting, we used a fixed cell size that was optimized with PBEsol+U approach in our previous studies. 23,28 The bandgaps obtained with different screening parameters are listed in Table 1. The results indicate that calculations with the screening parameter of 0.6 provide the best match between the computed and experimental bandgap. The primitive cell was then optimized with the chosen HSE parameters. The resulting cell dimensions, magnetic moment and bandgap are given in Table 2. Specifically, both a and c are smaller compared to those obtained using GGA+U and closer to the experimental values. Similarly, the computed magnetic moment of Cr has a better
8
ACS Paragon Plus Environment
Page 9 of 38 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
The Journal of Physical Chemistry
match with the experimental value. After relaxation, the bandgap increased negligibly from 3.377 eV to 3.384 eV. The resulting electronic density of states (DOS) of bulk Cr2 O3 is given in Figure 1. Cr-O phase diagram was evaluated by optimizing the Cr, O2 , and CrO2 phases with HSE parameters identical to those used in the relaxation of Cr2 O3 cell. The resulting phase diagram (given in SI Figure S1) yielded the formation energies of Cr2 O3 and CrO2 at −2.54 eV/atom and −2.17 eV/atom, respectively. The computed formation energy of Cr2 O3 compares well with the corresponding experimental value of −2.35 eV/atom.
Defects Defect calculations were performed using a 2×2×1 supercell containing vacancies and interstitials and generated from the optimized primitive cell. In addition to vacancies and interstitials, Cr vacancy triple defect (or split-vacancy) identified in our earlier work, 23 was also considered. For completeness, electronic structure and formation energies of complex defects such as Frenkel and Schottky defects were evaluated. When electronic structure is evaluated, these defect complexes were modeled after Lebreau et al., 22 where the individual vacancy and interstitial defects comprising the defect complex were in close proximity. However, these defect complexes were not accounted for in the diffusion studies because of their low mobility due to bigger sizes. The ionic positions in the defect supercells were optimized under fixed volume and shape conditions. To gain an understanding of how HSE affected defect properties, we evaluated electronic density of states (DOS), distribution of excess electrons and holes for charged defects, and the formation energies of the optimized defects.
Electronic Structure of Defects Previous studies of defects in Cr2 O3 mainly employed GGA+U method. Given that U and the width of the bandgap often have monotonic relationship U value is often fitted to recover 9
ACS Paragon Plus Environment
The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
experimental bandgap. The choice of the U-value not only affects the bandgap but also the position of the VBM and to a certain extent CBM with respect to the core bands. To get an understanding of the influence of U on the defect induced electronic states, we focus on the defect levels of the neutral Cr and O vacancies. In particular, in the studies of Lebreau et al., 22 and Rak and Brenner, 27 the U value of 5.0 eV within Dudarev formalism was used for Cr-3d electrons. According to the DOS of the vacancies reported in these two studies, neutral Cr vacancy exhibits two defects levels, one that is slightly above the VBM (0.2 or 0.3 eV, respectively) , and another that is unoccupied and at the bottom of the CBM. Differing from these studies, Gray et al. 24 used {U, J} = {4.5, 1} eV within the Liechtenstein formalism and Medasani et al. 23 used U=3.7 eV within Dudarev formalism. Both studies indicate that VCr introduces multiple defect levels. Gray et al.’s results point to four distinct defect levels (with two of them very close to each other and slightly below the CBM) in the bandgap, 24 and Medasani et al. reported 5 distinct defect levels. 23 If defect levels that are energetically closer are grouped together, then Gray et al and Medasani et al’s data predict the existence of three defect level zones: a) just above the VBM, b) 1 eV above the VBM, and c) just below the CBM. The HSE results confirm this pattern (Figure 2). This indicates that Uef f of ~3.7 eV is able to describe defect levels of VCr more accurately at the GGA+U level of theory. For VO0 , HSE predicts two defect levels that are spin degenerate at 0.8 and 2.6 eV above VBM (Figure 3). Similar electronic structure for VO0 was obtained in the study of Gray et al., 24 where the spin degenerate levels are at 1.3 and 2.6 eV and in the study of Carey et al., 46 who used UCr−3d = 5 eV and UO−2p = 5.5 eV, where the corresponding defect levels are at 1 and 2 eV above the VBM. Medasani et al. instead observed no degeneracy, but the resulting four levels are still within the bandgap. In contrast, Lebreau et al. 22 observed no defect levels within the bandgap for VO . These results indicate that the energetic position of defect electronic states with respect to the VBM are highly dependent on the choice of U parameters and on the +U formalism selected. Overall, comparison with the HSE results shows that although some GGA+U data qualitatively match the HSE results the difference
10
ACS Paragon Plus Environment
Page 10 of 38
Page 11 of 38 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
The Journal of Physical Chemistry
is often quite noticeable especially for neutral vacancies. Interestingly, for charged defects −3 such as VCr , VO2 , for interstitials, and for Cr Frenkel defect, eigen levels obtained from our
GGA+U calculations are in good agreement with those obtained from HSE calculations (see Figures 2-6). It is noteworthy that unlike the study by Gray et al 24 of native vacancies in Cr2 O3 that assumed fixed defect levels in the bandgap independent of the defect charge state, our earlier GGA+U and the present HSE studies on vacancies and interstitials indicate that defect levels in the bandgap are not fixed but vary depending on the defect charge state. Electronic structure of defect complexes with individual defects in close proximity is also evaluated. For Cr Frenkel defect, HSE predicts multiple defect eigen states near the VBM and one near the CBM within the bandgap as shown in Figure. 7. These eigen states are in good agreement with the GGA+U predicted defect eigen states reported by Medasani et al. 23 But with Uef f = 5 eV, GGA+U predicts only one Cr Frenkel defect eigen state within the bandgap. For Schottky defect complex, GGA+U predicts only three eigen states (one near VBM and two in the center) within bandgap. But for the same defect, HSE predicts multiple defect states within the bandgap as shown in Figure 8.
Charge distribution The distribution of the charge of excess hole(s) or electron(s) for the various charge states of the native vacancies and self-interstitials is shown in SI Fig. S2 - S5. These plots indicate that both holes and electrons are delocalized. This finding on the delocalization of the holes is consistent with the earlier studies. 23,33,47 However, comparison of charge distributions of Cri1 and Cri2 holes shows that delocalization is more pronounced when GGA+U is used. 0-K Defect Formation Energies After structure optimization, the spurious electrostatic interactions inherent in the periodic boundary formalism of charged defects were corrected using the anisotropic FNV (Freysoldt, Neugebauer, and Vande Walle) method 48,49 as implemented in PyCDT 50 (see Table S1 in
11
ACS Paragon Plus Environment
The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
SI for the resulting corrections). The charge transition levels of the defects with respect to the valence band maximum are presented in Table S2 in the SI. When compared to the corrected transition levels reported in our previous works, the corresponding transition levels obtained with HSE are overall at a higher Fermi level. Some of the transition levels such as Cri (3/2) transition, which is at 2.6 eV above the VBM, differ by nearly 1 eV (1.64 eV in Ref 28 ). Further, HSE predicts that only +3 and +2 charge states are stable for Cri within the bandgap. On the other hand, the Oi (0/-2) transition level at 3.07 eV differs only by 0.14 eV from the corresponding transition level computed using GGA+U coupled with bandgap correction. This shows that while bandgap correction overall improves the position of the transition levels computed with GGA+U, the improved agreement between the HSE computed and badgap corrected GGA+U formation energies is not systematic. The position of some transition levels such as Oi (0/-2) and VO (1/0) are almost the same, while the position of other transition levels such as Cri (3/2) and VO (2/1) exhibit differences that are greater than 0.4 eV when computed at the HSE and GGA+U levels of theory. The defect transition levels for Cr2 O3 in Figure S6 in SI show that whereas VO acts as a deep donor, Cri acts as a shallow donor. The Oi and VCr defects (regular and vacancy triple defect) behave as deep acceptors. Qualitatively, GGA+U computed transition levels of corresponding defects in our prior works 23,28 exhibit similar behavior. Comparing the defect states in Cr2 O3 with those of other corundum oxides such as Al2 O3 and α-Fe2 O3 can help in understanding the defect behavior in corundum oxides, a technologically important class of materials. The bandgap of Cr2 O3 at 3.2 eV lies in between those of Al2 O3 and Fe2 O3 at 8.8 eV and 2.2 eV respectively. Choi et al. 51 evaluated the defect thermodynamics in Al2 O3 using HSE. To reproduce the high bandgap of Al2 O3 , Choi et. al. used a mixing parameter of 0.32 which differs from the value of 0.25 used in this study. Choi et. al. identified that VAl and Oi were deep acceptors, Ali was a deep donor and VO could act as either a donor or an acceptor depending on the Fermi level. The point defects in Fe2 O3 were evaluated by Lee and Han 52 at GGA+U level of theory. Lee and Han used a Uef f of 4.3 eV within
12
ACS Paragon Plus Environment
Page 12 of 38
Page 13 of 38 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
The Journal of Physical Chemistry
Dudarev formalism to obtain bandgap that matches with the experimental results. Their results indicate that both VO and F ei are deep donors, but F ei is relatively a shallow donor when compared to VO . Their results also indicate that Oi has no transition states within Fe2 O3 bandgap and VF e is a shallow acceptor. The results indicate that the behavior of defects in both Al2 O3 and Fe2 O3 is markedly different from that of Cr2 O3 . The dependence of 0K defect formation energies on the position of the Fermi level for both O-rich and Cr-rich conditions shows that not only transition levels, but also the formation energies differ between HSE and GGA+U even after corrections (Figure S6 in SI). At Cr-rich phase boundary, HSE predicts lower formation energies than GGA+U for the dominant Cri and VO defects. At O2 -rich phase boundary HSE predicted formation energies are lower for neutral Cr vacancies and higher for Oi . In contrast, for charged Cr vacancies, the HSE formation energies are in good agreement with those calculated using the GGA+U. Representative defect formation energies for O-rich conditions are given in Table 3 keeping the Fermi level at the middle of the bandgap. For defect complexes comprising of individual defects in close proximity, when compared to GGA+U, HSE predicts higher formation energies. HSE predicted formation energy of Cr Frenkel defect predicted at 4.56 eV is 2 eV higher than the corresponding GGA+U prediction at 2.6 eV. Similarly, the formation energy of Schottky defect predicted with HSE at 9.24 eV is 5.24 eV higher than that of the GGA+U prediction. These results validate our choice of using HSE to obtain more accurate formation energies needed for the calculation of self-diffusion coefficients. Defect energetics of the Schottky and Frenkel defect complexes comprising of individual defects under isolated conditions were also evaluated. The Frenkel defect energy was calculated as Etot,vac,i + Etot,inter,i − 2Etot,bulk , where Etot,vac , Etot,inter , and Etot,bulk refer to the total energies of the vacancy, interstitial and bulk supercells, and i refers to one of the Cr and O species. The resulting formation energies of O and Cr Frenkel defects are 11.92 eV and 6.71 eV respectively. The formation energy of Schottky defect complex comprsing of two Cr vacancies and 3 O vacancies under non-interacting conditions is evaluated as 2 ∗ Etot,V 3− + Etot,V 2+ − 5 ∗ (N − 1)/N Etot,bulk . The Cr
13
ACS Paragon Plus Environment
O
The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
resulting formation energy of Schottky defect is 18.12 eV.
Self Diffusion Coefficients We evaluated the self-diffusion of Cr and O species mediated by vacancies and interstitials along the diffusion pathways identified in earlier works. 22,23,28 Activation energies for defect migration are sum of transition state barrier energies computed at the the GGA+U level of theory, as reported in our previous works, 23,28 and the defect formation energies computed at the HSE level of theory. We expect this approach to yield reliable activation energies at lower computational cost due to insensitivity of GGA+U defect migration barrier energies to the variation in U value. 27 The details of the transition states and the barrier energies of the diffusion pathways identified in our previous works 23,28 are reproduced in the Supplementary Information for reader’s convenience. The resulting self-diffusion coefficients are plotted as a function of temperature in Figure 9 at four different conditions: a) standard conditions (1 atm), b) moderate vacuum (1 × 10−8 atm), c) high vacuum (1 × 10−14 atm), and d) metal/metal oxide interface conditions. Electronic effects such as band bending and the finite electric field arising at the metal-metal oxide and vacuum-metal oxide interfaces are ignored in this study. The plots indicate that at O2 partial pressure (pO2 ) of 1 atm, Cr has higher mobility than O. The difference in the mobilities is higher at lower temperatures. As pO2 is reduced to 1 × 10−8 atm, O become the dominant diffusing species at temperatures higher than 1100 K. At lower temperatures (< 1100 K), Cr still has a higher mobility, but the magnitude of mobility difference is reduced. As pO2 decreases further to 1 × 10−14 atm, O mobility becomes higher than that of Cr at even lower temperatures (800 K). At Cr/Cr2 O3 interface, relative O mobility is significantly higher than that of Cr at all temperatures studied. Sabioni reported the experimental diffusion coefficients of Cr and O as (2.2 − 7.2)×10−18 cm s−2 and (3.2 − 8.4)×10−18 cm s−2 measured at 1300° C and 1100° C respectively under high vacuum conditions. 12,13 The experimental diffusion coefficients reported by Sabioni et 14
ACS Paragon Plus Environment
Page 14 of 38
Page 15 of 38 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
The Journal of Physical Chemistry
al. were nearly the same for Cr and O, and relatively independent of pO2 . Sabioni proposed the suprising independence of the diffusion coefficients from pO2 could be due to the presence of impurities in the Cr2 O3 samples at around 1 ppm level. The resulting high concentration of extrinsic defects may have affected the diffusion coefficients in the Cr2 O3 samples. While the diffusion coefficients reported by Sabioni could not be considered representative of intrinsic Cr2 O3 , they could serve as higher limits on one of the diffusion coefficients by promoting the compensating defects of either Cr or O species. At 1373 K and at 1 × 10−14 atm, which are representative of experimental conditions used to measure O diffusion, the computed diffusion coefficient of O is found to be at 1 × 10−18 cm s−2 . This indicates that O diffusion coefficient was not affected by the extrinsic defects. At 1573 K, the maximal value of Cr diffusion coefficient observed by us is at 1 × 10−18 cm s−2 (at standard conditions), which matches with the experimental value of Cr diffusion. 12 To identify the nature of the predominant defects facilitating the diffusion of each species we plotted the defect concentrations in Brouwer diagrams at different temperatures (Fig. 10). At low temperatures (< 800 K), Cr vacancies (including Cr triple defects denoted as VCr−T D ) have higher concentrations at all pO2 . At 1200 K, the concentration of VO is higher at pO2 below 1 × 10−8 atm. Above 1 × 10−8 atm, Cr vacancies become the dominant defects. Finite temperature defect thermodynamics calculations by Lee and Han for Fe2 O3 52 indicate that at high temperatures (1373 K), F ei and VO are the dominant defects at low and high pO2 ’s respectively. However at low pO2 , the difference in the concentrations of F ei and VO in Fe2 O3 is negligible. For Cr2 O3 also at high temperature and low pO2 conditions, the concentration of Cri defects is closer to that of VO . The qualitative behavior of the Brower plots at high temperatures looks similar for both Cr2 O3 and Fe2 O3 . This could be ascribed to the higher concentrations of free carriers at high temperatures resulting in similar defect concentration profiles due to the charge balance constraint. Due to the high concentration of Cr vacancies at pressures close to 1 atm and temperatures above 1200 K, the self diffusion of Cr is equal to the one measured by Sabioni, 12 suggesting that self-compensating Cr vacancies could be
15
ACS Paragon Plus Environment
The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
responsible for the high Cr diffusion coefficient (when compared to the intrinsic Cr diffusion coefficients at corresponding pO2 ). These results indicate that vacancies are the mediators of reported self-diffusion coefficients in Cr2 O3 . If the impurity concentration in Cr2 O3 samples could be minimized to very low levels, Fig. 10c indicates that at high temperatures (> 1200 K) and low pO2 (< 1 × 10−12 atm), VO and Cri could be responsible for O and Cr diffusion respectively. However the relative mobility of O is very high in such cases, indicating again that vacancies dominate the diffusion process in Cr2 O3 . The Brouwer diagrams also reveal that it would be misleading to use commonly reported 0-K defect formation energies at the phase boundary edges to predict dominant defects responsible for experimentally measured diffusion coefficients. For example at Cr-rich boundary, the 0-K defect energy diagram (Figure S6 in SI) indicates that Cri could be the dominant defect if the Fermi level, EF , lies between 1.5 to 1.9 eV. However, the allowed Fermi level range is only 1.97 − 2.04 eV (Figure 11) and neutral vacancies are the dominant defects. Fermi level variation is plotted as a function of pO2 for different temperatures in Figure 12a. The plots indicate that Cr2 O3 is intrinsically p-type at low temperatures and high pO2 , where Fermi level can become as low as 1.2 eV. Under such conditions, the dominant defects are VCr (Figure S6 in SI). At low temperatures and low pO2 , the Fermi level is around 1.75 eV making Cr2 O3 slightly n-type intrinsically. Stoichiometry deviation plotted as a function of pO2 in Figure 12(b) shows that the deviation is negligible in intrinsic Cr2 O3 . Figure 13 shows the anisotropic ratio of diffusion for Cr and O at various operational conditions. The results indicate that Cr2 O3 basal plane is the preferential diffusion pathway for O under all conditions. For Cr, however, the preferred diffusion orientation is temperature dependent. At low temperatures, Cr diffuses mainly along the c-axis. The mode of diffusion switches to basal plane pathway at around 1200 K under high pO2 . In ultra high vacuum conditions of pO2 ≤1 × 10−14 atm the primary mode of Cr diffusion switches back to c-axis at around 1600 K. Under such conditions, Cr diffusion coefficients along basal plane and
16
ACS Paragon Plus Environment
Page 16 of 38
Page 17 of 38 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
The Journal of Physical Chemistry
c-axis become nearly the same.
Summary Towards a comprehensive understanding of the diffusion processes in Cr2 O3 , we computed the self diffusion coefficients and Brouwer diagrams in Cr2 O3 . Diffusion coefficients were computed using the DFT at two accuracy levels and the Einstein random walk formalism. The defect formation energies were reevaluated with HSE hybrid functional for high accuracy, and were corrected for spurious electrostatic interaction errors arising in the periodic supercell method. The resulting self diffusion coefficients have a good agreement with the experimental data. Our results reveal that vacancies are primarily responsible for both Cr and O diffusion in intrinsic Cr2 O3 . At high temperatures and low oxygen partial pressures, O has higher mobility. Cr has higher mobility at lower temperatures at moderate to high oxygen partial pressure. The preferred path for O diffusion is along the basal plane. Cr diffusion takes place along the c-axis at low temperatures and switches to basal plane at temperatures above 1200 K. The revealed mechanism for self-diffusion in Cr2 O3 in a wide range of temperature and pressure conditions opens new avenues for the design of Cr2 O3 materials.
Acknowledgement This work was supported by the U.S. Department of Energy (DOE), Office of Science, Basic Energy Sciences, Materials Sciences and Engineering Division. Simulations were performed using PNNL Institutional Computing facility. PNNL is a multi-program national laboratory operated by Battelle for the U.S. DOE under Contract DEAC05-76RL01830. The data and the associated scripts used to compute the diffusion coefficients are available from https: //github.com/mbkumar/Cr2O3_diffusion.
17
ACS Paragon Plus Environment
The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Supporting Information Available Cr-O phase diagram, defect charge distributions, defect formation energies at 0 K, all computed with HSE functional. Transition state structures, and diffusion barrier energies computed at GGA+U level reproduced from our previous works. 23,28
References (1) Bates, M. K.; Jia, Q.; Ramaswamy, N.; Allen, R. J.; Mukerjee, S. Composite Ni/NiOCr2 O3 Catalyst for Alkaline Hydrogen Evolution Reaction. J. Phys. Chem. C 2015, 119, 5467–5477. (2) Xie, Z.; Fan, J.; Cheng, Y.; Jin, L.; Hu, G.; Lu, J.; Luo, M. Cr2 O3 Catalysts for Fluorination of 2-Chloro-3,3,3-Trifluoropropene to 2,3,3,3-Tetrafluoropropene. Ind. & Eng. Chem. Res. 2013, 52, 3295–3299. (3) Arca, E.; Kehoe, A. B.; Veal, T. D.; Shmeliov, A.; Scanlon, D. O.; Downing, C.; Daly, D.; Mullarkey, D.; Shvets, I. V.; Nicolosi, V.; Watson, G. W. Valence Band Modification of Cr2 O3 by Ni-Doping: Creating a High Figure of Merit p-type TCO. J. Mater. Chem. C 2017, 5, 12610–12618. (4) Farrell, L.; Fleischer, K.; Caffrey, D.; Mullarkey, D.; Norton, E.; Shvets, I. V. Conducting Mechanism in the Epitaxial p-type Transparent Conducting Oxide Cr2 O3 :Mg. Phys. Rev. B: Condens. Matter Mater. Phys. 2015, 91, 125202. (5) Mu, S.; Wysocki, A. L.; Belashchenko, K. D. Effect of Substitutional Doping on the N´eel Temperature of Cr2 O3 . Phys. Rev. B: Condens. Matter Mater. Phys. 2013, 87, 054435. (6) Pati, S. P.; Shimomura, N.; Nozaki, T.; Shibata, T.; Sahashi, M. N´eel Temperature of
18
ACS Paragon Plus Environment
Page 18 of 38
Page 19 of 38 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
The Journal of Physical Chemistry
Cr2 O3 in Cr2 O3 /Co Exchange-Coupled System: Effect of Buffer Layer. J. Appl. Phys. 2015, 117, 17D137. (7) Zurek, J.; Young, D.; Essuman, E.; H¨ansel, M.; Penkalla, H.; Niewolak, L.; Quadakkers, W. Growth and Adherence of Chromia Based Surface Scales on Ni-base Alloys in High- and Low-pO2 Gases. Mater. Sci. Eng.: A 2008, 477, 259 – 270, 3rd International Conference on Spray Deposition and Melt Atomization (SDMA 2006) and the 6th International Conference on Spray Forming (ICSF VI). (8) Schreiber, D.; Olszta, M.; Saxey, D.; Kruska, K.; Moore, K.; Lozano-Perez, S.; Bruemmer, S. Examinations of Oxidation and Sulfidation of Grain Boundaries in Alloy 600 Exposed to Simulated Pressurized Water Reactor Primary Water. Microsc. Microanal. 2013, 19, 676–687. (9) Schreiber, D.; Olszta, M.; Bruemmer, S. Grain Boundary Depletion and Migration During Selective Oxidation of Cr in a Ni-5Cr Binary Alloy Exposed to High-Temperature Hydrogenated Water. Scr. Mater. 2014, 89, 41 – 44. (10) Brozek, V.; Pokorny, P.; Bouska, P.; Stoulil, J.; Mastny, L. Corrosion Properties of Chromia based Eco-friendly Coatings on Mild Steel. Metalurgija 2016, 55, 675–678. (11) Sabioni, A. C. S.; Huntz, A. M.; Philibert, J.; Lesage, B.; Monty, C. Relation Between the Oxidation Growth Rate of Chromia Scales and Self-Diffusion in Cr2 O3 . J. Mater. Sci. 1992, 27, 4782–4790. (12) Sabioni, A.; Lesage, B.; Huntz, A.; Pivin, J.; Monty, C. Self-Diffusion in Cr2 O3 I. Chromium Diffusion in Single Crystals. Philos. Mag. A 1992, 66, 333–350. (13) Sabioni, A.; Huntz, A.; Millot, F.; Monty, C. Self-Diffusion in Cr2 O3 II. Oxygen Diffusion in Single Crystals. Philos. Mag. A 1992, 66, 351–360.
19
ACS Paragon Plus Environment
The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
(14) Sabioni, A.; Huntz, A.; Millot, F.; Monty, C. Self-diffusion in Cr2 O3 III. Chromium and Oxygen Grain-Boundary Diffusion in Polycrystals. Philos. Mag. A 1992, 66, 361–374. (15) Latu-Romain, L.; Parsa, Y.; Mathieu, S.; Vilasi, M.; Galerie, A.; Wouters, Y. Towards the Growth of Stoichiometric Chromia on Pure Chromium by the Control of Temperature and Oxygen Partial Pressure. Corros. Sci. 2017, 126, 238 – 246. (16) Huntz, A. M.; Tsai, S. C. Diffusion in Oxide Scales: Application to Cr2 O3 Scales. J. Mater. Sci. Lett. 1994, 13, 821–825. (17) Kofstad, P.; Lillerud, K. Chromium Transport Through Cr2 O3 Scales I. On Lattice Diffusion of Chromium. Oxid. Met. 1982, 17, 177–194. (18) Tsai, S.; Huntz, A.; Dolin, C. Growth Mechanism of Cr2 O3 Scales: Oxygen and Chromium Diffusion, Oxidation Kinetics and Effect of Yttrium. Mater. Sci. Eng. A 1996, 212, 6–13. (19) Hoshino, K.; Peterson, N. Cation Self-Diffusion in Cr2 O3 . J. Am. Ceram. Soc. 1983, 66, c202–c203. (20) Schmucker, E.; Petitjean, C.; Martinelli, L.; Panteix, P.-J.; Lagha, S. B.; Vilasi, M. Oxidation of Ni-Cr Alloy at Intermediate Oxygen Pressures. I. Diffusion Mechanisms Through the Oxide Layer. Corros. Sci. 2016, 111, 467–473. (21) Catlow, C. R. A.; Corish, J.; Hennessy, J.; Mackrodt, W. C. Atomistic Simulation of Defect Structures and Ion Transport in α-Fe2 O3 and α-Cr2 O3 . J. Am. Ceram. Soc. 1988, 71, 42–49. (22) Lebreau, F.; Islam, M. M.; Diawara, B.; Marcus, P. Structural, Magnetic, Electronic, Defect, and Diffusion Properties of Cr2 O3 : A DFT+U Study. J. Phys. Chem. C 2014, 118, 18133–18145.
20
ACS Paragon Plus Environment
Page 20 of 38
Page 21 of 38 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
The Journal of Physical Chemistry
(23) Medasani, B.; Sushko, M. L.; Rosso, K. M.; Schreiber, D. K.; Bruemmer, S. M. Vacancies and Vacancy-Mediated Self Diffusion in Cr2 O3 : A First-Principles Study. J. Phys. Chem. C 2017, 121, 1817–1831. (24) Gray, C.; Lei, Y.; Wang, G. Charged Vacancy Diffusion in Chromium Oxide Crystal: DFT and DFT+U Predictions. J. Appl. Phys. 2016, 120, 215101. (25) Cao, P.; Wells, D.; Short, M. P. Anisotropic Ion Diffusion in α-Cr2 O3 : An Atomistic Simulation Study. Phys. Chem. Chem. Phys. 2017, 19, 13658. (26) Vaari, J. Molecular Dynamics Simulations of Vacancy Diffusion in Chromium(III) Oxide, Hematite, Magnetite and Chromite. Solid State Ion. 2015, 270, 10 – 17. (27) Rak, Z.; Brenner, D. W. First-Principles Investigation of Diffusion and Defect Properties of Fe and Ni in Cr2 O3 . J. Appl. Phys. 2018, 123, 155105. (28) Medasani, B.; Sushko, M. L.; Rosso, K. M.; Schreiber, D. K.; Bruemmer, S. M. FirstPrinciples Investigation of Native Interstitial Diffusion in Cr2 O3 . J Phys. Chem. C 2018, 122, 12984–12993. (29) Mehrer, H. Diffusion in Solids; Springer-Verlag: Berlin, 2007. (30) Liechtenstein, A. I.; Anisimov, V. I.; Zaanen, J. Density-Functional Theory and Strong Interactions: Orbital Ordering in Mott-Hubbard Insulators. Phys. Rev. B: Condens. Matter Mater. Phys. 1995, 52, R5467–R5470. (31) Dudarev, S. L.; Botton, G. A.; Savrasov, S. Y.; Humphreys, C. J.; Sutton, A. P. Electron-Energy-Loss Spectra and the Structural Stability of Nickel Oxide:
An
LSDA+U Study. Phys. Rev. B: Condens. Matter Mater. Phys. 1998, 57, 1505–1509. (32) Rohrbach, A.; Hafner, J.; Kresse, G. Ab initio Study of the (0001) Surfaces of Hematite and Chromia: Influence of Strong Electronic Correlations. Phys. Rev. B: Condens. Matter Mater. Phys. 2004, 70, 125426. 21
ACS Paragon Plus Environment
The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
(33) Kehoe, A. B.; Arca, E.; Scanlon, D. O.; Shvets, I. V.; Watson, G. W. Assessing the Potential of Mg-Doped Cr2 O3 as a Novel p-type Transparent Conducting Oxide. J. Phys. Condens. Matter 2016, 28, 125501. (34) Cococcioni, M.; de Gironcoli, S. Linear Response Approach to the Calculation of the Effective Interaction Parameters in the LDA + U Method. Phys. Rev. B: Condens. Matter Mater. Phys. 2005, 71, 035105. (35) Kresse, G.; Hafner, J. Ab initio Molecular Dynamics for Liquid Metals. Phys. Rev. B: Condens. Matter Mater. Phys. 1993, 47, 558–561. (36) Kresse, G.; Hafner, J. Ab initio Molecular-Dynamics Simulation of the Liquid-Metal– Amorphous-Semiconductor Transition in Germanium. Phys. Rev. B: Condens. Matter Mater. Phys. 1994, 49, 14251–14269. (37) Kresse, G.; Furthm¨ uller, J. Efficient Iterative Schemes for Ab Initio Total-Energy Calculations Using a Plane-Wave Basis Set. Phys. Rev. B: Condens. Matter Mater. Phys. 1996, 54, 11169–11186. (38) Krukau, A. V.; Vydrov, O. A.; Izmaylov, A. F.; Scuseria, G. E. Influence of the Exchange Screening Parameter on the Performance of Screened Hybrid Functionals. J. Chem. Phys. 2006, 125 . (39) Bl¨ochl, P. E. Projector Augmented-Wave Method. Phys. Rev. B: Condens. Matter Mater. Phys. 1994, 50, 17953–17979. (40) Kresse, G.; Joubert, D. From Ultrasoft Pseudopotentials to the Projector AugmentedWave Method. Phys. Rev. B: Condens. Matter Mater. Phys. 1999, 59, 1758–1775. (41) Perdew, J. P.; Burke, K.; Ernzerhof, M. Generalized Gradient Approximation Made Simple. Phys. Rev. Lett. 1996, 77, 3865–3868.
22
ACS Paragon Plus Environment
Page 22 of 38
Page 23 of 38 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
The Journal of Physical Chemistry
(42) Perdew, J. P.; Ruzsinszky, A.; Csonka, G. I.; Vydrov, O. A.; Scuseria, G. E.; Constantin, L. A.; Zhou, X.; Burke, K. Restoring the Density-Gradient Expansion for Exchange in Solids and Surfaces. Phys. Rev. Lett. 2008, 100, 136406. (43) Vineyard, G. H. Frequency Factors and Isotope Effects in Solid State Rate Processes. J. Phys. Chem. Solids 1957, 3, 121 – 127. (44) Youssef, M.; Yildiz, B. Intrinsic Point-Defect Equilibria in Tetragonal ZrO2 : Density Functional Theory Analysis with Finite-Temperature Effects. Phys. Rev. B: Condens. Matter Mater. Phys. 2012, 86, 144109. (45) Chase, M. W. J. NIST-JANAF Thermochemical Tables, 4th Edition; American Institute of Physics: New York, 1998. (46) Carey, J. J.; Legesse, M.; Nolan, M. Low Valence Cation Doping of Bulk Cr2 O3 : Charge Compensation and Oxygen Vacancy Formation. J. Phys. Chem. C 2016, 120, 19160– 19174. (47) Lany, S.; Zunger, A. Assessment of Correction Methods for the Band-Gap Problem and for Finite-Size Effects in Supercell Defect Calculations: Case Studies for ZnO and GaAs. Phys. Rev. B: Condens. Matter Mater. Phys. 2008, 78, 235104. (48) Freysoldt, C.; Neugebauer, J.; Van de Walle, C. G. Fully Ab Initio Finite-Size Corrections for Charged-Defect Supercell Calculations. Phys. Rev. Lett. 2009, 102, 016402. (49) Kumagai, Y.; Oba, F. Electrostatics-Based Finite-Size Corrections for First-Principles Point Defect Calculations. Phys. Rev. B: Condens. Matter Mater. Phys. 2014, 89, 195205. (50) Broberg, D.; Medasani, B.; Zimmerman, N. E. R.; Yue, G; Canning, A.; Haranczyk, M.; Asta, M.; Hautier, G. PyCDT: A Python toolkit for modeling point defects in semiconductors and insulators. Comp. Phys. Comm. 2018, 226, 165 – 179. 23
ACS Paragon Plus Environment
The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
(51) Choi, M.; Janotti, A.; Van de Walle, C. G. Native point defects and dangling bonds in α-Al2 O3 . J. Appl. Phys. 2013, 113, 044501. (52) Lee, J.; Han, S. Thermodynamics of native point defects in α-Fe2 O3 : an ab initio study. Phys. Chem. Chem. Phys. 2013, 15, 18906.
24
ACS Paragon Plus Environment
Page 24 of 38
Page 25 of 38 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
The Journal of Physical Chemistry
Table 1: Bandgap of Cr2 O3 obtained with different values of HSE functional screening parameter. µ
Bandgap (eV)
0.2
4.27
0.3
3.85
0.4
3.64
0.5
3.59
0.6
3.38
Table 2: Bulk properties of Cr2 O3 obtained with HSE functional. Property This Work
Literature Theory (GGA+U) Experiment
a
4.925
c
13.54
13.68, 23 13.85 22
13.566
µ
2.81
2.91, 23 3.1 22
2.8
Bandgap
3.38
2.95, 23 2.8 22
3.4
25
ACS Paragon Plus Environment
The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Table 3: Defect formation energies (in eV) of Cr2 O3 with Fermi level (EF ) at the middle of the bandgap under O-rich conditions. Defect VO
VCr
Oi
Cri Cr Frenkel Schottky
q HSE GGA+U 23,28 0 5.55 6.01 1 5.62 6.34 2 6.16 7.47 0 0.90 1.93 -1 0.93 1.69 -2 1.30 1.37 -3 2.00 1.66 0 2.94 2.18 -1 4.77 3.69 -2 5.76 4.38 0 11.65 10.13 1 9.58 8.82 2 7.85 7.91 3 6.92 7.97 0 4.56 2.6 0 9.24 3.97
26
ACS Paragon Plus Environment
Page 26 of 38
Page 27 of 38
Density of states
total
Cr
Density of states Density of states
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
The Journal of Physical Chemistry
−2
O d p s
−1
0
1 2 Energies (eV)
3
4
5
Figure 1: Electronic DOS of bulk Cr2 O3 computed with HSE functional with 0.6 screening parameter.
27
ACS Paragon Plus Environment
The Journal of Physical Chemistry
DOS
0 VCr
DOS
−1 VCr
DOS
−2 VCr
−3 VCr
d p s
DOS
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 28 of 38
−2
−1
0
1 2 Energies (eV)
3
4
5
Figure 2: Orbital projected DOS averaged over all sites for Cr vacancies in different charge states. Dashed vertical line marks the position of the Fermi level. The zero level has been fixed at the VBM.
28
ACS Paragon Plus Environment
Page 29 of 38
DOS
VO0
DOS
VO1
VO2 d p s
DOS
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
The Journal of Physical Chemistry
−2
−1
0
1 2 Energies (eV)
3
4
5
Figure 3: Orbital projected DOS averaged over all sites for O vacancies in different charge states. Dashed vertical line marks the position of the Fermi level. The zero level has been fixed at the VBM.
29
ACS Paragon Plus Environment
The Journal of Physical Chemistry
DOS
0 VCr − split
DOS
−1 VCr − split
DOS
−2 VCr − split
−3 VCr − split
d p s
DOS
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 30 of 38
−2
−1
0
1 2 Energies (eV)
3
4
5
Figure 4: Orbital projected DOS averaged over all sites for Cr vacancy triple defects in different charge states. Dashed vertical line marks the position of the Fermi level. The zero level has been fixed at the VBM.
30
ACS Paragon Plus Environment
Page 31 of 38
DOS
Cri0
DOS
Cri1
DOS
Cri2
Cri3
d p s
DOS
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
The Journal of Physical Chemistry
−2
−1
0
1 2 Energies (eV)
3
4
5
Figure 5: Orbital projected DOS averaged over all sites for Cr interstitials in different charge states. Dashed vertical line marks the position of the Fermi level. The zero level has been fixed at the VBM.
31
ACS Paragon Plus Environment
The Journal of Physical Chemistry
DOS
Oi0
DOS
Oi−1
DOS
Oi−2 d p s
−2
−1
0
1 2 Energies (eV)
3
4
5
Figure 6: Orbital projected DOS averaged over all sites for O dumbbell interstitials in different charge states. Dashed vertical line marks the position of the Fermi level. The zero level has been fixed at the VBM.
Cr Frenkel d p s
DOS
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 32 of 38
−2
−1
0
1 2 Energies (eV)
3
4
5
Figure 7: Orbital projected DOS averaged over all sites for Cr Frenkel defect. Dashed vertical line marks the position of the Fermi level. The zero level has been fixed at the VBM.
32
ACS Paragon Plus Environment
Page 33 of 38
DOS
Cr Frenkel d p s
−2
−1
0
1 2 Energies (eV)
3
4
5
Figure 8: Orbital projected DOS averaged over all sites for Schottky defect complex with 2 Cr and 3 O ions removed. Dashed vertical line marks the position of the Fermi level. The zero level has been fixed at the VBM.
pO2 1.0e+00 atm
0
log D (cm2.s)
log D (cm2.s)
Cr O
−10
−20 −30
−20 −30 −40
−40 −50
pO2 1.0e-08 atm
0
Cr O
−10
6
7
8
9
104/T
10
(K−1)
11
12
13
−50
14
6
7
8
9
12
13
14
Cr/Cr2O3 In erface
0
Cr O
−10
Cr O
log D (cm2.s)
−10
−20 −30 −40 −50
11
(b)
pO2 1.0e-14 atm
0
10
104/T (K−1)
(a)
log D (cm2.s)
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
The Journal of Physical Chemistry
−20 −30 −40
6
7
8
9
10
104/T (K−1)
11
12
13
−50
14
6
7
(c)
8
9
10
104/T (K−1)
11
12
13
14
(d)
Figure 9: Self diffusion coefficients of Cr and O at vacuum/oxide interface at O2 partial pressures of a) 1 atm, b) 1 × 10−8 atm, c) 1 × 10−14 atm, d) and at Cr/oxide interace.
33
ACS Paragon Plus Environment
The Journal of Physical Chemistry
300K
800K
10−26
10−42 10−50 10−58 10−66
10−12
VCr VO VCr − TD Oi Cri
Defect concentration
10−34
Defect concentration
n h
10−15
10−21
10−82
10−30 10−12
10−10
10−8
pO2 (atm)
10−4
10−6
n h
10−24 10−27
10−14
VCr VO VCr − TD Oi Cri
10−18
10−74
10−14
10−12
(a)
10−10
10−8
pO2 (atm)
10−6
10−4
(b) 1200K 10−9
Defect concentration
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 34 of 38
10−12 10−15 10−18 10−21
VCr VO VCr − TD Oi Cri
n h
10−24 10−14
10−12
10−10
10−8
pO2 (atm)
10−6
10−4
(c)
Figure 10: Brouwer diagrams depicting the defect concentrations vs. O2 partial pressure in Cr2 O3 for a) 300 K, b) 800 K, and c) 1200 K.
34
ACS Paragon Plus Environment
Page 35 of 38
2.05 2.04 2.03 2.02
EF (eV)
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
The Journal of Physical Chemistry
2.01 2.00 1.99 1.98 1.97
5
10
15
20
104/T
(K−1)
25
30
35
Figure 11: Fermi level in intrinsic Cr2 O3 at Cr/Cr2 O3 interface (Cr-rich condition) as a function of temperature.
35
ACS Paragon Plus Environment
The Journal of Physical Chemistry
1e−10 300 K Δ00 K 1200 K
0.5 1.7
0.0
ΔX of Cr
1.6
EF (eV)
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Page 36 of 38
1.5 1.4 1.3
300 K 400 K 600 K 800 K 1000 K 1200 K
10−15
−0.5
−1.0
−1.5 10−13
10−11
10−9
pO2 (atm)
10−7
10−5
10−15
10−3
(a)
10−13
10−11
10−9
pO2 (atm)
10−7
10−5
10−3
(b)
Figure 12: (a) Fermi level in intrinsic Cr2 O3 , and (b) stoichiometry deviation expressed in terms of relative change in Cr composition form stoichiometry as functions of temperature and O2 partial pressure.
36
ACS Paragon Plus Environment
Page 37 of 38
3.5 3.0
pO2 1.0e+00 a m
3.0
2.0 1.5
2.5 2.0 1.5
1.0
1.0
0.5
0.5 6
Cr c-axis Cr ab-plane O c-axis O ab-plane
3.5
2.5
0.0
pO2 1.0e-08 atm
4.0
Cr c-axis Cr ab-plane O c-axis O ab-plane
Anisotropic Ratio
4.0
Aniso ropic Ra io
8
104/T
10
(K−1)
12
0.0
14
6
8
(a)
10
104/T (K−1)
12
14
(b) pO2 1.0e-14 atm
4.0
Cr c-axis Cr ab-plane O c-axis O ab-plane
3.5 3.0
Anisotropic Ratio
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
The Journal of Physical Chemistry
2.5 2.0 1.5 1.0 0.5 0.0
6
8
10
104/T (K−1)
12
14
(c)
Figure 13: Anisotropic ratio of diffusion of Cr and O at vacuum/oxide interface at O2 partial pressures of a) 1 atm, b) 1 × 10−8 atm, c) 1 × 10−14 atm.
37
ACS Paragon Plus Environment
The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60
Graphical TOC Entry
38
ACS Paragon Plus Environment
Page 38 of 38