Band Gap in Magnetite above Verwey Temperature Induced by Symmetry Breaking.

Magnetite exhibits a famous phase transition, called Verwey transition, at the critical temperature TV of about 120 K. Although numerous efforts have been devoted to the understanding of this interesting transition, up to now, it is still under debate whether a charge ordering and a band gap exist in magnetite above TV. Here, we systematically investigate the charge ordering and the electronic properties of magnetite in its cubic phase using different methods based on density functional theory: DFT+U and hybrid functionals. Our results show that, upon releasing the symmetry constraint on the density but not on the geometry, charge disproportionation (Fe2+/Fe3+) is observed, resulting in a band gap of around 0.2 eV at the Fermi level. This implies that the Verwey transition is probably a semiconductor-to-semiconductor transition and that the conductivity mechanism above TV is small polaron hopping.

Introduction

Magnetite, the oldest magnetic material, shows a great prospect for applications in magnetic resonance imaging (MRI), drug delivering,[1] and spintronic devices[2] and thus has attracted a tremendous attention.[3] At room temperature, magnetite crystallizes in an inverse spinel structure with oxygen anions arranged in a slightly distorted face centered cubic lattice and iron atoms occupying tetrahedral and octahedral interstitial sites. It is assumed that both divalent and trivalent irons, in a ratio of Fe2+:Fe3+ = 1:2, exist in magnetite, with tetrahedral sites occupied by Fe3+ and octahedral sites occupied by an equal number of Fe2+ and Fe3+.[4] Below 858 K, magnetite is a ferrimagnet with the cations at octahedral sites coupling antiferromagnetically with the cations at tetrahedral sites.[5] The electrical conductivity of Fe3O4 is around 250 Ω–1 cm–1 at room temperature, like a poor metal.[6,7] On cooling from room temperature to 120 K, the electrical conductivity decreases gradually and suddenly drops by 2 orders of magnitude at around 120 K, when the crystal structure changes from cubic to monoclinic symmetry,[8,9] which is called Verwey transition.[6]

Attracted by the puzzling Verwey transition, many efforts have been devoted to the investigation of magnetite electronic properties. By photoemission spectroscopy (PES)[10−18] and scanning tunneling spectroscopy (STS),[19] researchers concluded that above TV, magnetite is half-metal because the valence band emission in PES[10−18] and the signal in STS[19] start from the vicinity of the Fermi level (EF) and a band gap of 0.14–0.3 eV exists at EF below TV. However, using high-resolution PES, Park et al. found that, on heating through TV, the band gap in magnetite did not collapse but was merely reduced by ∼50 meV.[20] This change in the gap perfectly accounts for the conductivity jump at TV. Later, Jordan also proposed that a ∼0.2 eV gap exists in magnetite above TV by an STS study.[21] Another STS investigation on magnetite nanocrystals shows that a band gap of 140–250 meV is present below TV and of 75 ± 10 meV above TV.[22] However, one may question that both PES and STS are surface sensitive techniques and the results may depend on the surface considered.[23] Schrupp et al.[14] probed deeper bulk layers by using higher energy photons and provided evidence for the existence of strongly bound small polarons, which were proposed to be responsible for the conductivity above TV.

In contrast to all these experimental observations, on the theoretical side, all DFT calculations for bulk magnetite in cubic phase lead to a coincident conclusion, i.e., magnetite above TV is half-metal ferrimagnet possessing 100% spin polarized charge carriers at EF.[10,24−29]

Another critical issue regarding magnetite is the charge ordering since Verwey originally proposed that the statistical distribution of Fe2+ and Fe3+ at octahedral sites, accounting for the high electronic conductivity, would lead to some type of order at lower temperature.[6] In recent years, both experiments[9,30−32] and DFT calculations[33−36] approved the charge ordering at octahedral sites in Fe3O4 below TV. Above TV, according to the Anderson’s condition,[37] the long-range charge ordering could be lost, whereas short-range charge ordering should be maintained. Recent X-ray diffuse scattering experiments conducted by Alexey et al. show that short-range charge ordering survives up to room temperature.[38] However, all the DFT calculations for magnetite above TV (i.e., in the cubic phase) result in an average valence state Fe2.5+ at octahedral sites,[10,24−29,39] except one unrestricted Hartree–Fock study, where charge disproportionation was proposed.[40]

Because the electronic properties of magnetite above TV are still under debate and an inconsistency between previous DFT calculations and diffuse scattering experiments on the charge ordering exists, it is a necessary and urgent issue to reinvestigate Fe3O4 above TV theoretically. Because of the strong correlation effects among Fe 3d electrons, standard DFT calculations fail to provide an accurate description of Fe3O4 electronic properties and predict magnetic moments for tetrahedral Fe ions (3.34–3.47 μB[26,27]) that are much smaller than the experimental value of about 4.2 μB.[41,42] In contrast, DFT+U calculations, including an additional orbital-dependent interaction for highly localized d orbitals, and hybrid functional calculations, including partial exact exchange, improve the magnetic moment to a value of 4.1 μB.[27] Thus, here we carefully investigate Fe3O4 above TV using both DFT+U and hybrid functional calculations. We impose the symmetry constraint (space group Fd3̅m) only on atomic positions but not on electron density. Our results, for the first time, predict a charge disproportionation (Fe2+/Fe3+) and a small indirect band gap (around 0.2 eV) in cubic phase magnetite, implying a small polaron hopping mechanism for the conductivity above TV, in agreement with high energy photons experiments.[14]

Methods

The spin polarized DFT+U calculations were performed using the plane-wave-based Quantum ESPRESSO package.[43] The projector augmented wave (PAW) potentials were adopted to describe the electron–ion interactions, with Fe (3s, 3p, 3d, 4s) and O (2s, 2p) treated as valence electrons. The exchange and correlation interaction was described by the Perdew–Burke–Ernzerhof (PBE) functional within the generalized gradient approximation (GGA). To test the U effect, different U values from 1 to 4 eV were used both for the geometry optimization and the electronic property calculations. Energy cutoffs of 64 Ryd and 782 Ryd (for kinetic energy and charge density expansion, respectively) were adopted for all calculations. The convergence criterion of 0.026 eV/Å for force was used during geometry optimization, and the convergence criterion for total energy was set at 10–6 Ryd for all the calculations. Bader charge analysis,[44−46] considering only the valence electrons, was used to estimate the charge of Fe ions in Fe3O4. A fine FFT grid (360 × 360 × 360) was carefully tested and used to accurately produce the charge.

Hybrid functional calculations were carried out using the CRYSTAL14 package based on DFT,[47,48] where the Kohn–Sham orbitals are expanded in Gaussian-type orbitals (the all-electron basis sets are O|8–411G* and Fe|8–6–411G*, according to the scheme previously used for Fe2O3[49]). The HSE06[50] and B3LYP[51,52] functionals were used for both geometry optimization and electronic properties calculations. The convergence criterion of 0.023 eV/Å for force was used during geometry optimization, and the convergence criterion for total energy was set at 10–7 Hartree for all the calculations. Mulliken population analysis was used to estimate the charge of Fe ions in Fe3O4.

For all the calculations, the primitive cell containing eight oxygen atoms and six iron atoms was adopted and the k points generated by the Monkhorst–Pack scheme were chosen to be 10 × 10 × 10. To check the polaron effect, the coordinates of all atoms in the primitive cell (with lattice parameters kept fixed at the optimized values) were relaxed without the symmetry constraint for both PBE+U and HSE methods. After that, the bond lengths of FeOct3+-O slightly decrease and the bond lengths of FeOct2+-O slightly increase compared with those in the symmetric structure, by only about ±0.03 Å for both PBE+U (U = 3.5 eV) and HSE calculations. The charge difference variation between FeOct3+ and FeOct2+ is negligible: 0.002 e for PBE+U and 0.067 e for HSE. The band gap increases to 0.450 eV for PBE+U (U = 3.5 eV) and 0.776 eV for HSE compared with 0.177 eV (U = 3.5 V) and 0.421 eV (HSE) for the symmetric geometry. Considering that above TV, Fe3O4 crystallizes in an exact inverse spinel structure with Fd3̅m symmetry, the slightly distorted lattice based on DFT calculations at 0 K may not suitable to model Fe3O4 above TV. This may be because that, at temperature above TV, electrons hop so frequently between Fe2+ and Fe3+ at octahedral sites (as discussed in the main text) that the lattice distortion cannot follow the electron transfer. Therefore, to partially consider the temperature effect and properly model magnetite in cubic phase above TV, the symmetric unit cell is used throughout the manuscript.

To check the charge ordering, single-point energy calculations were conducted on a conventional cell containing 32 oxygen atoms and 24 iron atoms using a uniform k points mesh of 6 × 6 × 6. We did not consider spin–orbit coupling in our calculations because spin–orbital splitting of the 3d band was found to be 2 orders of magnitude smaller than the crystal field splitting in previous calculations for cubic Fe3O4.[53]

Results and Discussion

The geometry of Fe3O4 was fully relaxed at the PBE+U and at the hybrid functional (HSE06 and B3LYP) level of theory, as shown in Figure S1 in the Supporting Information (SI). The Fd3̅m symmetry was kept during the geometry relaxation, otherwise it would have been lost, as already reported by Andrew et al.[36] As compared with the experimental value of 8.394 Å for the conventional cell,[42] PBE+U and B3LYP slightly overestimate the lattice parameter, i.e., 8.491 Å for U = 3.5 eV and 8.448 Å for B3LYP. In the case of PBE+U the overestimation is in accordance with a previous study.[27] In contrast, HSE06 gives a value (8.389 Å) that is much closer to the experimental one. The detailed structural information on Fe3O4 is listed in Table S1 and Table S2 in the SI.

According to recent X-ray diffuse scattering experiments,[38] the short-range (∼1.5 nm just above TV) charge ordering should be maintained above TV. Thus, for comparison, we performed self-consistent field (SCF) calculations for magnetite both with and without the symmetry constraint on the electron density and wave function. In contrast with standard PBE, the asymmetric ground state is found to have a much lower total energy than the symmetric state for HSE06, B3LYP, and PBE+U calculations. As shown in Table , the energy difference per Fe3O4 formula unit (ΔE) is around −0.15 eV for hybrid functional calculations (HSE06 and B3LYP) and is as large as −0.423 eV for PBE+U (U = 3.5 eV). Therefore, the ground state wave function should be asymmetric even though the geometry is symmetric. For the asymmetric electronic ground state, the four Fe ions at octahedral sites (FeOct) are divided into two groups, two FeOct3+ and two FeOct2+, according to the charge states listed in 1 and represented in Figure . To check whether the charge ordering exists, we conducted analogous calculations on the conventional cell. Different FeOct2+/FeOct3+ ions patterns can be obtained by both PBE+U and HSE calculations, as shown in Figure S2 in SI. Those only slightly differ in energy (less than ∼15 meV per formula unit). Electrons could hop between FeOct3+ and FeOct2+. Both experimental (selective excitation double mössbauer study)[54] and theoretical studies[40] show that the electron hopping time at the octahedral sites is less than 10–11 s at room temperature. Therefore, we may expect that, above TV, no charge ordering can be observed. Note that the charge difference between FeOct2+ and FeOct3+ is not 1.0 e but around 0.27–0.32 e, which is similar to the charger ordering in the low temperature phase as reported before.[9,30,42] The calculated magnetic moments of Fe ions at tetrahedral sites (FeTet) are −4.212 μB from HSE06, −4.139 μB from B3LYP, and −3.956μB from PBE+U (U = 3.5 eV), which are in good agreement with the experimental value of ∼4.2 μB.[41,42] In accordance with the charge disproportionation, the magnetic moments of FeOct ions are also divided into two groups. Experimentally, because of the frequent hopping of electrons between FeOct3+ and FeOct2+, an average magnetic moment at octahedral sites can be measured (3.97 μB),[42] which agrees well with our mean results (3.840 μB for U = 3.5 eV, 4.029 μB for HSE, and 3.977 μB for B3LYP). Note that previous HSE06 calculations by Junghyun et al. also led to a magnetic moment ordering at octahedral sites.[27] However, the authors did not give credit to this result and concluded that it was incorrect on the basis of some additional PBE+U calculations that did not show any charge disproportionation.

Table 1

Magnetic Moment and Net Charge of Fe Ions at Tetrahedral and Octahedral Sites Labeled As mTet, mOct, QTet, and QOct, Energy Difference per Fe3O4 Formula (ΔE) between Symmetric (Esy) and Asymmetric (Easy) Electronic Ground States (ΔE = Easy – Esy) and Band Gap Value (Eg) Derived from Different Methods: PBE, PBE+U, HSE06, and B3LYP

methods	mTet (μB)	mOct (μB)	QTet (e)	QOct (e)	ΔE (eV)	Eg (eV)
PBE	–3.584	3.626	+1.553	+1.540	0
PBE+U (U = 3.5 eV)	–3.956	4.058, 4.054	+1.655	+1.735, +1.732	–0.423	0.177
		3.624, 3.624		+1.464, +1.463
HSE06	–4.212	4.271, 3.789	+2.189	+2.225, +1.906	–0.157	0.421
		3.785, 4.270		+1.905, +2.224
B3LYP	–4.139	4.211, 4.197	+2.133	+2.173, +2.163	–0.146	0.340
		3.750, 3.750		+1.879, +1.879

Figure 1

Repeated primitive cell of Fe3O4 above TV with charge ordering pattern.

To investigate the U effect on the Fe3O4 system further, we calculated the magnetic moments, charge difference, and energy difference with different U values, as plotted in Figure . The magnetic moment of FeTet decreases from −3.988 μB to −3.584 μB, with the U value going from 4 to 0 eV. When the U value is smaller than 1 eV, the energy difference (ΔE) goes to zero and the charge disproportionation disappears with the four FeOct irons becoming identical, or all FeOct2.5+. This is because standard PBE and PBE with a very small U term fail to describe the strong correlation between d electrons that causes symmetry breaking of wave function, inducing the charge disproportionation. In addition, hybrid functional calculations with different percentage of exact exchange were also performed. Similar to the U effect, the magnetic moment decreases with the decreasing of the percent of exact exchange. The detail information is reported in Table S4 in the SI.

Figure 2

U value effect on magnetic moment mTet, charge difference between FeOct3+ and FeOct2+ (ΔQ), energy difference ΔE and band gap Eg as obtained by PBE+U calculations. The left Y axis indicates mTet, and the right Y axis indicates ΔQ, ΔE, and Eg. Note that in order to unify the unit, ΔQ is multiplied by 1 V.

To understand and rationalize the electronic properties of Fe3O4, projected density of states (PDOS) on 3d orbitals of Fe ions were calculated using PBE, PBE+U, and the hybrid functional HSE06, as shown in Figure . PDOS by standard PBE, based on a plane wave basis set (Figure a), agrees well with previous results based on LSDA[25] and PBE,[27] showing that some prominent states of minority FeOct 3d electrons cross EF, with a band gap (0.657 eV) in the majority spin-up channel. Note that PBE calculations with and without symmetry give the same PDOS due to the lack of charge disproportionation. When a U value of 3.5 eV is included in the SCF calculation with the symmetry constraint, the band gap in the spin-up channel increases to 1.539 eV (Figure b). This can be understood as follows: the U term corrects the strong Coulomb interaction between 3d electrons (self-interaction) in the PBE functional and thus shifts the filled 3d states into a deeper energy level, increasing the gap between the filled and empty 3d states. Still, the PDOS shows a half-metallic character, which agrees well with previous studies.[27] However, after removing the symmetry constraint on the electron density, a small band gap of 0.177 eV opens for the spin-down channel at EF (Figure c), in excellent agreement with the band gap obtained by PES (0.1 ± 0.03 eV)[20] and STS (0.2 eV)[21] experiments. The original states from FeOct2.5+ crossing EF split into two parts, i.e., states from FeOct2+ below EF and states from FeOct3+ above EF. PDOS with U = 4 eV exhibits a similar character (see Figure S3 in the SI) with a band gap of 0.474 eV, which is however much larger than the experimental value (0.1–0.2 eV).[20,21] As shown in Figure , the band gap decreases with the U values decreasing and goes to zero when U ≤ 3 eV. Considering that a big U (4 eV) gives a too-large band gap and small U (≤3 eV) results in a too-small magnetic moment for Fe ions, the proper U value should be set at around 3.5 eV. PDOSs with different U values are shown in Figure S3 in the SI.

Figure 3

Projected density of states (PDOS) of magnetite calculated with different methods. The left part (a–c) are from Quantum Espresso code, and the right part (d–f) are from CRYSTAL14 code. The U value used in (b) and (c) is 3.5 eV. Legend of colors is on the top of the panels.

The PBE calculation based on atomic orbitals basis sets (Figure d) provides analogous PDOS to that based on plane wave basis sets, indicating that the calculated electronic properties of Fe3O4 is not dependent on basis sets. Then hybrid functional (HSE06 and B3LYP) calculations, based on atomic orbitals basis sets, were conducted to be compared with the results from PBE+U. As shown in Figure e, symmetrical ground states from HSE06 exhibit a half-metal character with a much larger band gap for spin-up channel compared with that from PBE and PBE+U methods. This can be attributed to the large percent of exact exchange (25%) included in HSE06. Similar to PBE+U, the unsymmetrical ground state from HSE06 also shows a semiconductor character but with a much bigger band gap of 0.421 eV at EF (Figure f). The B3LYP calculation gives a band gap of 0.340 eV and PDOS very similar to that by HSE, as shown in Figure S4 in the SI. However, the band gap decreases quickly to zero when the percent of exact exchange included in HSE06 and B3LYP is reduced, as shown in Table S4 and Figure S5 in the SI.

Even though the band gaps given by U = 3.5 eV, HSE, and B3LYP are different, one prominent common character in the DOS (Figure c,f) is that a peak of spin-down states exists 0.5 eV below EF. This peak can also be found in many experimental PES.[13−17,20,55] Spin-resolved PES experiments[13,15] further confirm that the peak at −0.5 eV is from spin-down states, which agrees well with our results. As shown in Figure c, below EF, the spin-up states arise at around −0.5 eV, just following the peak of spin-down states. While in the DOS from HSE06 (Figure f), the spin-up states arise at about −1.9 eV, far away from the peak of spin-down states. In experiments, the intensity of spin-up states starts increase from around −1 eV,[13,15] and −1 eV is the watershed for the total spin polarization.[17] Thus, PBE+U may slightly underestimate the exchange energy between Fe d electrons and HSE06 may overestimate the exchange energy. Moreover, in experimental PES,[13−17,20,55] the intensity at EF is rather weak or zero, which is different from the DOS by the symmetric wave function that exhibits some prominent states at EF, as shown in Figure a,b,d,e. In contrast, the DOS by the asymmetric wave function exhibits a semiconducting character (Figure c,f), in very good agreement with the experiments. The semiconductor property can also reconcile the conflict between half-metallicity and the low spin polarization (40%–65%) at EF observed in experiments.[15,56,57] In additions, the conductivity of magnetite increases with the increasing of temperature above TV, which is not compatible with metallic properties. On the contrary, if magnetite is computed to be a semiconductor, the experimental conductivity behavior is readily understood in terms of small polaron hopping. To summarize, magnetite is found to be a semiconductor even above TV (i.e., in the cubic phase) with an indirect band gap, according to the DOS in Figure and also to the band structure shown in Figure S6 in the SI.

In Fe3O4, the spin-up bands and the spin-down bands are split by exchange energy and the five degenerate d-electron levels of the FeOct ions are further split by the crystal field into three t2g and two eg levels. In a previous study,[21] a simple schematic picture for d electrons in FeOct ions was proposed to explain the half-metallicity of Fe3O4, in which five d electrons of FeOct2+ and FeOct3+ occupy the spin-up states, and the sixth electron of FeOct2+ occupies the spin-down t2g band at EF, giving rise to conductivity. However, half-metal is not the ground state of Fe3O4, as discussed above, and the schematic picture must be modified by considering the symmetry breaking among FeOct ions. Here, the three-dimensional charge density at every peak from d electrons in the PDOS is calculated to investigate the arrangement of t2g and eg orbitals of Fe ions. As shown in Figure , the valence band maximum (VBM) (Figure a) and conduction band minimum (CBM) (Figure b) for spin-up electrons are contributed by eg electrons from FeOct2+ and t2g electrons from FeTet3+, respectively. The VBM (Figure c) and the CBM (Figure d) for spin-down electrons are contributed by t2g electrons from FeOct2+ and FeOct3+. Note that Figure d has a different orientation from Figure a–c in order to see the t2g character more clearly. States near EF are marked with different d orbitals in the PDOS in Figure c,f. Therefore, the energy barrier for electrons hopping from FeOct2+ to FeOct3+ is probably comparable with the band gap (because the electron must jump from the t2g of FeOct2+ to the t2g of FeOct3+) and the conductivity of Fe3O4 is mainly contributed by the hopping of t2g electrons between FeOct2+ and FeOct3+.

Figure 4

Charge density plots for the VBM and CBM of spin-up (a,b) and spin-down (c,d) electrons from HSE calculations. PBE+U (U = 3.5 eV) gives the same results. The isosurface level is 1.8 × 10–5 electron/bohr3.

The PDOS with peaks marked with t2g and eg for a broader energy spectrum are shown in Figure S7 in the SI. Accordingly, we propose a new schematic diagram of the d levels for the FeOct3+ and FeOct2+ cations, as shown in Figure . The spin-up d orbitals of FeOct are all occupied, leaving the spin-down orbitals unoccupied, except one t2g orbital of FeOct2+ occupied with the extra electron. The unoccupied d orbitals of FeOct3+ above EF are between the occupied and unoccupied t2g orbitals of FeOct2+, resulting in a band gap of around 0.2 eV at EF. Besides, for FeOct ions, one occupied eg orbital with even deeper energy than t2g orbitals involves mixing with O 2p states in a bonding mode and the other eg orbital with higher energy in an antibonding mode.

Figure 5

Schematic diagram of the d levels for the FeOct3+ and FeOct2+ cations in Fe3O4 obtained from PBE+U and HSE calculations.

Conclusions

In summary, high temperature phase of magnetite (cubic) was investigated using both DFT+U and hybrid functional calculations. Equivalent number of Fe2+ and Fe3+ exists at octahedral sites with a charge difference of about 0.3 e, instead of average Fe2.5+ charge states, as predicted by standard DFT calculations in this and in previous studies. The symmetry breaking among FeOct ions splits the t2g states around EF into two parts, t2g from FeOct2+ below EF and t2g from FeOct3+ above EF, and thus induces a small indirect band gap of about 0.2 eV. Electron hopping from FeOct2+ to FeOct3+, overcoming an energy barrier comparable to the band gap, gives rise to the conductivity above TV. Our results for the first time provide theoretical evidence for the existence of the charge disproportionation and of a small band gap in magnetite cubic phase, which is helpful to understand the experimentally observed conductivity behavior of magnetite above TV and rationalize it in terms of a small polaron hopping mechanism.