Coulomb Parameter U and Correlation Strength in LaFeAsO V. I. Anisimov,1 Dm. M. Korotin,1 S. V. Streltsov,1 A. V. Kozhevnikov,1, 2 J. Kuneš,3 A. O. Shorikov,1 and M. A. Korotin1 arXiv:0807.0547v1 [cond-mat.str-el] 3 Jul 2008 2 1 Institute of Metal Physics, Russian Academy of Sciences, 620041 Yekaterinburg GSP-170, Russia Joint Institute for Computational Sciences, Oak Ridge National Laboratory P.O. Box 2008 Oak Ridge, TN 37831-6173, USA 3 Theoretical Physics III, Center for Electronic Correlations and Magnetism, Institute of Physics, University of Augsburg, Augsburg 86135, Germany First principles constrained density functional theory scheme in Wannier functions formalism has been used to calculate Coulomb repulsion U and Hund’s exchange J parameters for iron 3d electrons in LaFeAsO. Results strongly depend on the basis set used in calculations: when O-2p, As-4p, and Fe-3d orbitals and corresponding bands are included, computation results in U =3÷4 eV, however, with the basis set restricted to Fe-3d orbitals and bands only, computation gives parameters corresponding to F 0 =0.8 eV, J=0.5 eV. LDA+DMFT (the Local Density Approximation combined with the Dynamical Mean-Field Theory) calculation with this parameters results in weakly correlated electronic structure that is in agreement with X-ray experimental spectra. PACS numbers: 74.25.Jb, 71.45.Gm Following the discovery of high-Tc superconductivity in iron oxypnictide LaO1−x Fx FeAs [1], a question of the influence of electronic correlation effects on the normal and superconducting properties of LaFeAsO has arisen. In striking similarity with high-Tc cuprates, undoped material LaFeAsO is not superconducting with antiferromagnetic commensurate spin density wave developing below 150 K [2]. Only when electrons (or holes) are added to the system via doping, antiferromagnetism is suppressed and superconductivity appears. As it is generally accepted that Coulomb correlations between copper 3d electrons are responsible for cuprates anomalous properties, it is tempting to suggest that the same is true for iron 3d electrons in LaFeAsO. Correlation strength in a system is determined by ratio of Coulomb interaction parameter U and band width W . If U/W is significantly less than 1 then the system is weakly correlated and results of the Density Functional Theory (DFT) calculations are reliable enough to explain its electronic and magnetic properties. However, if U value is comparable with W or even larger then the system is in intermediate or strongly correlated regime and Coulomb interactions must be explicitly treated in electronic structure calculations. For LaFeAsO the bands formed by Fe-3d states have width ≈4 eV (see shaded area in the lower panel of Fig. 1), so an estimation for Coulomb interaction parameter U should be compared with this value. In practical calculations, U is often considered as a free parameter to achieve the best agreement of calculated and measured properties of investigated system. Sometimes U value could be estimated from the experimental spectra. The first principles justified methods to determine Coulomb interaction parameter U value are constrained DFT scheme [3], where in DFT calculations the d-orbital occupancies are fixed to the certain values and U is numerically determined as a derivative of d- orbital energy over its occupancy and Random Phase Approximation (RPA) method, where screened Coulomb interaction between d-electrons is calculated via perturbation theory [4]. In Ref. 5 it was proposed to use in LaFeAsO U =4 eV obtained in RPA calculations for metallic iron [6]. This value for Coulomb parameter (with Hund’s exchange parameter J=0.7 eV) was used in Dynamical Mean-Field Theory (DMFT) [7] calculations for LaFeAsO [5, 8, 9]. Results of these works show iron 3d electrons being in intermediate or strongly correlated regime, as it is natural to be expected for Coulomb parameter value U =4 eV and Fe-3d band width ≈4 eV. The most direct way to estimate correlation effects strength in a system under consideration is to compare the experimental spectra with densities of states (DOS) obtained in DFT calculations. For strongly correlated materials additional features in the experimental photoemission and absorption spectra appear that are interpreted as lower and upper Hubbard bands absent in the DFT DOS. If no such additional features are observed and DOS obtained in DFT calculations satisfactorily describe the experimental spectra then the material is in weakly correlated regime. LaFeAsO was studied by soft X-ray absorption and emission spectroscopy [10], X-ray absorption spectroscopy (O K-edge) [11] and photoemission spectroscopy [12]. In all these works the conclusion was that DOS obtained in DFT calculations gave good agreement with the experimental spectra and the estimations for Coulomb parameter value are U <1 eV [11]. That contradiction with results of the DMFT calculations [5, 8, 9] shows that first principles calculation of Coulomb interaction parameter U value for LaFeAsO is needed to determine the correlation effects strength in this material. Results of such calculations by constrained DFT calculations are reported in the present work. We have obtained 2 10 8 6 4 2 0 -1 Density of states, (eV.atom) the value U <1 eV for Fe-3d band that agrees with the estimates from spectroscopy. Recently the RPA calculations for Coulomb interaction parameter U in LaFeAsO were reported, where U was estimated as 1.8÷2.7 eV [13]. It is important to note that Coulomb interaction parameter U value depends on the model where it will be used and, more precisely, on the choice of the orbital set that is taken explicitly into account in the model. For example, in constrained DFT calculations for high-Tc cuprates the resulting U value for Cu d-shell was found between 8 and 10 eV [14]. The U value in this range was used in cluster calculations where all Cu d-orbitals and p-orbitals of neighboring oxygens were taken into account and calculated spectra agree well with experimental data [15]. However, in one band model, where only x2 − y 2 orbital per cooper atom is explicitly included in the calculations, the U value giving good agreement with experimental data falls down to 2.5÷3.6 eV [16], that is 3-4 times smaller than constrained DFT value. The same situation occurs for titanium and vanadium oxides: the U value from constrained DFT calculations is ≈6 eV and cluster calculations where all d-orbitals and p-orbitals of neighboring oxygens were taken into account with U close to this value gave good agreement between calculated and experimental spectra [17]. However, in the model where only partially filled t2g orbitals are included, much smaller U value (corresponding to Slater integral F 0 =3.5 eV) gives the results in agreement with experimental data [18]. It is interesting that such a small U value can be obtained in constrained DFT calculations for titanates and vanadates where only t2g -orbital occupancies are fixed while all the other states (eg -orbitals of vanadium and p-orbitals of oxygens) allowed to relax in self-consistent iterations [18, 19]. So the calculation scheme used in constrained DFT (the set of the orbitals with fixed occupancies) should be consistent with basis set of the model where the calculated U value will be used. Another source of uncertainty in constrained DFT calculation scheme is a definition of atomic orbitals whose occupancies are fixed and energy calculated. In some DFT methods, like Linearized Muffin-Tin Orbitals (LMTO) [20], these orbitals could be identified with LMTO. However, in other DFT calculation schemes, where plane waves are used as a basis, like in pseudopotential method [21] one should use more general definition for localized atomic like orbitals such as Wannier functions (WFs) [22]. The practical way to calculate WFs for specific materials using projection of atomic orbitals on Bloch functions was developed in [23]. In Fig. 1 the total and partial DOS for LaFeAsO obtained in LMTO calculations are shown. Crystal field splitting for Fe-3d orbitals in this material is rather weak (∆cf =0.25 eV) and all five d orbitals of iron form common band in the energy region (−2, +2) eV relative to the Fermi level (see grey region on the bottom panel in Total 2 0 p-O 2 0 p-As 2 0 d-Fe -5 -4 -3 -2 -1 0 1 Energy, eV, EF=0 2 FIG. 1: Total and partial densities of states for LaFeAsO obtained in DFT calculation in frame of LMTO method. Fig. 1). There is a strong hybridization of iron t2g orbitals with p orbitals of arsenic atoms which form nearest neighbors tetrahedron around iron ion. This effect becomes apparent in the energy interval (−3, −2) eV (white region on the bottom panel in Fig. 1) where band formed by p orbitals of arsenic is situated. More week hybridization with oxygen p states reveals in (−5.5, −3) eV energy window (black region on the bottom panel in Fig. 1). We have calculated Coulomb interaction U and Hund’s exchange J parameters for WFs basis set via constrained DFT procedure with fixed occupancies for WFs of d symmetry. For this purpose we have used two calculation schemes based on the different DFT methods. One of them involves linearized muffin-tin orbitals produced by the TB-LMTO-ASA code [20]; corresponding WFs calculation procedure is described in details in Ref. 24. The second one is based on the plane waves obtained within the pseudopotential plane-wave method PWSCF, as implemented in the Quantum ESPRESSO package [21], and described in details in Ref. 25. The difference between the results of these two schemes could give an estimation for the error of U and J determination. The WFs are defined by the choice of Bloch functions Hilbert space and by a set of trial localized orbitals that will be projected on these Bloch functions. We performed calculations for two different choices of Bloch functions and atomic orbitals. One of them includes only bands predominantly formed by Fe-3d orbitals in the energy window (−2, +2) eV and equal number of Fe-3d orbitals to be projected on the Bloch functions for these bands. That choice corresponds to the model where only five dorbital per Fe site are included but all arsenic and oxygen 3 TABLE I: The constrained DFT calculated values of average Coulomb interaction Ū and Hund’s exchange J (eV) parameters for d-symmetry Wannier functions computed with two different sets of bands and orbitals. DFT method separate Fe-3d band full bands set TB-LMTO-ASA Ū =0.49, J=0.51 Ū =3.10, J=0.81 PWSCF Ū =0.59, J=0.53 Ū =4.00, J=1.02 FIG. 2: (Color online) Module square of dx2 −y 2 -like Wannier function computed for Fe-3d bands only (left panel) and for full set of O-2p, As-4p and Fe-3d bands (right panel). Big sphere in the center marks Fe ion position and four small spheres around it correspond to As neighbors. p-orbitals are omitted. Second choice includes all bands in energy window (−5.5, +2) eV that are formed by O2p, As-4p and Fe-3d states and correspondingly full set of O-2p, As-4p and Fe-3d atomic orbitals to be projected on Bloch functions for these bands. That would correspond to the extended model where in addition to d-orbitals all p-orbitals are included too. In both cases we obtained Hamiltonian in WF basis that reproduces exactly bands predominantly formed by Fe-3d states in the energy window (−2, +2) eV (Fig. 1), but in the second case in addition to that bands formed by p-orbitals in the energy window (−5.5, −2) eV will be reproduced too. However, WFs with d-orbital symmetry computed in those two cases have very different spatial distribution. In Fig. 2 the module square of dx2 −y2 -like WF is plotted. While for the case when full set of bands and atomic orbitals was used (right panel) WF is nearly pure atomic d-orbital (iron states contribute 99%), WF computed using Fe-3d bands only is much more extended in space (left panel). It has significant weight on neighboring As ions with only 67% contribution from central iron atom. The physical reason for such effect is p − d hybridization that is treated explicitly in the case where both pand d-orbitals are included. In the case where only Fe-3d bands are included in calculation p − d hybridization reveals itself in the shape of WF. Fe-3d bands in the energy window (−2, +2) eV correspond to antibonding combination of Fe-3d and As-4p states and that is clearly seen on the left panel of Fig. 2. The different spatial distribution for two WFs calculated with full and restricted orbital bases can be expected to lead to different effective Coulomb interaction for electrons occupying these states. The results of constrained DFT calculations of the average Coulomb interaction Ū and Hund’s exchange J parameters for electrons on WFs computed with two different set of bands and orbitals (and using two different DFT methods: LMTO and pseudopotential) are presented in Tab. I. One can see that very different Coulomb interaction strength is obtained for separate Fe-3d band and full bands set calculations. While the latter gives value 3÷4 eV, that agrees with previously used values [5, 8, 9], separate Fe-3d band calculation results in 0.5÷0.6 eV, that is much smaller but agrees with spectroscopy estimations [11]. The main reason for such a drastic difference between two calculations is very different spatial extension of the two WFs (see Fig. 2): nearly complete localization on central iron atom for “full bands set” WF (99%) and only 67% for “Fe-3d band set” WF. Another possible reason for strong reduction of the calculated Ū value in going from “full bands set” to “Fe-3d band set” WF is screening via p − d hybridization with As-4p band that is situated just below Fe-3d band (see Fig. 1). The effect of decreasing of the effective Ū value in several times going from full orbital model to restricted basis was found previously for high-Tc cuprates (U =8÷10 eV for full p − d-orbitals basis [14] and 2.5÷3.6 eV for one-band model [16]). In constrained DFT calculations one obtains an average Coulomb interaction Ū that can be estimated as Ū = F 0 − J/2. Hence, Slater integral F 0 can be calculated as F 0 = Ū +J/2 [3]. For “Fe-3d band set” WF that gives F 0 =0.8 eV at J=0.5 eV. With this set of parameters we performed the LDA+DMFT [26] calculations (for detailed description of the present computation scheme see Ref. 24). The DFT band structure was calculated within the TB-LMTO-ASA method [20]. Crystal structure parameters were taken from Ref. 1. The restricted basis set including only Fe-3d WFs was used in the LDA+DMFT calculations. The effective impurity model for the DMFT was solved by the QMC method in Hirsh-Fye algorithm [27]. Calculations were performed for the value of inverse temperature β=10 eV −1 . Inverse temperature interval 0 < τ < β was divided into 100 slices. 6 · 106 QMC sweeps were used in self-consistency loop within the LDA+DMFT scheme and 12 · 106 of QMC sweeps were used to calculate the spectral functions. The iron 3d orbitally resolved spectral functions obtained within DFT and LDA+DMFT calculations are presented in Fig. 3. The influence of correlation effects on the electronic structure of LaFeAsO is minimal: there are small changes of peak positions for 3z 2 − r2 , xy and x2 − y 2 orbitals (the shift toward the Fermi energy) and 4 0,5 LDA+DMFT calculation with these parameters results in weakly correlated nature of iron d bands in this compound. This conclusion is supported by several X-ray spectroscopic investigations of this material. Support by the Russian Foundation for Basic Research under Grant No. RFFI-07-02-00041, Civil Research and Development Foundation together with Russian Ministry of science and education through program Y4-P05-15, Russian president grant for young scientists MK1184.2007.2 and Dynasty Foundation is gratefully acknowledged. J.K. acknowledges the support of SFB 484 of the Deutsche Forschungsgemeinschaft. Spectral function, eV -1 xy 0 0,5 yz, zx 0 2 2 3z -r 0,5 0 2 x -y 0,5 0 -3 -2 0 -1 1 Energy, eV, EF=0 2 2 3 FIG. 3: Partial densities of states for different Fe-3d orbitals obtained within the DFT (filled areas) and LDA+DMFT orbitally resolved spectral functions for F 0 =0.8 eV, J=0.5 eV (bold lines). practically unchanged picture of spectral function distribution for yz, zx bands. There are no appearance of neither quasiparticle peak on the Fermi level nor Hubbard bands in the energy spectrum with such values of U and J. Hence LaFeAsO can be considered as weakly correlated material. This agrees with the results of soft X-ray absorption and emission spectroscopy study [10]. It was concluded there that LaFeAsO does not represent strongly correlated system since Fe L3 X-ray emission spectra do not show any features that would indicate the presence of the low Hubbard band or the quasiparticle peak that were predicted by the LDA+DMFT analysis [5, 8, 9] with the large U =4 eV. A comparison of the X-ray absorption spectra (O K-edge) with the LDA calculations gave an upper limit of the on-site Hubbard U ≈1 eV [11]. Photoemission spectroscopy study of LaFeAsO suggests [12] that the line shapes of Fe 2p core-level spectra correspond to an itinerant character of Fe 3d electrons. It was demonstrated there that the valence-band spectra are generally consistent with band-structure calculations except for the shifts of Fe 3d-derived peaks toward the Fermi level. In conclusion, we have calculated the values of U and J via constrained DFT procedure in the basis of WFs. For minimal model including only Fe-3d orbitals we have obtained Coulomb parameters F 0 =0.8 eV, J=0.5 eV. The [1] Y. Kamihara et al., J. Am. Chem. Soc. 130, 3296 (2008). [2] M. A. McGuire et al., arXiv: 0804.0796; C. de la Cruz et al., Nature 453, 899 (2008). [3] P. H. Dederichs et al., Phys. Rev. Lett. 53, 2512 (1984); O. Gunnarsson et al., Phys. Rev. B 39, 1708 (1989); V. I. Anisimov et al., ibid. 43, 7570 (1991). [4] I. V. Solovyev et al., Phys. Rev. B 71, 045103 (2005); F. Aryasetiawan et al., ibid. 74, 125106 (2006). [5] K. Haule et al., Phys. Rev. Lett. 100, 226402 (2008). [6] T. Miyake et al., Phys. Rev. B 77, 085122 (2008). [7] A. Georges et al., Rev. Mod. Phys. 68, 13 (1996). [8] L. Craco et al., arXiv: 0805.3636. [9] A. O. Shorikov et al., arXiv: 0804.3283. [10] E. Z. Kurmaev et al., arXiv: 0805.0668. [11] T. Kroll et al., arXiv: 0806.2625. [12] W. Malaeb et al., arXiv: 0806.3860. [13] K. Nakamura et al., arXiv: 0806.4750. [14] M. S. Hybertsen et al., Phys. Rev. B 39, 9028 (1989); M. S. Hybertsen et al., ibid. 41, 11068 (1990); A. K. McMahan et al., ibid. 42, 6268 (1990). [15] H. Eskes et al., Phys. Rev. B 44, 9656 (1991). [16] Th. Maier et al., Phys. Rev. Lett. 85, 1524 (2000); A. Macridin et al., Phys. Rev. B 71, 134527 (2005); W.G. Yin et al., J. Phys.: Conf. Series 108, 012032 (2008). [17] A. E. Bocquet et al., Phys. Rev. B 53, 1161 (1996). [18] I. A. Nekrasov et al., Phys. Rev. B 72, 155106 (2005). [19] I. Solovyev et al., Phys. Rev. B 53, 7158 (1996). [20] O. K. Andersen, Phys. Rev. B 12, 3060 (1975); O. Gunnarsson et al., ibid. 27, 7144 (1983). [21] S. Baroni et al., http://www.pwscf.org [22] G. H. Wannier, Phys. Rev. 52, 191 (1937). [23] N. Marzari et al., Phys. Rev. B 56, 12847 (1997); W. Ku et al., Phys. Rev. Lett. 89, 167204 (2002). [24] V. I. Anisimov et al., Phys. Rev. B 71, 125119 (2005). [25] D. M. Korotin et al., arXiv: 0801.3500. [26] V. I. Anisimov et al., J. Phys.: Condens. Matter 9, 7359 (1997); A. I. Lichtenstein et al., Phys. Rev. B 57, 6884 (1998); K. Held et al., Phys. Stat. Sol. (b) 243, 2599 (2006). [27] J. E. Hirsch et al., Phys. Rev. Lett. 56, 2521 (1986).