research papers\(\def\hfill{\hskip 5em}\def\hfil{\hskip 3em}\def\eqno#1{\hfil {#1}}\)

Journal logoJOURNAL OF
SYNCHROTRON
RADIATION
ISSN: 1600-5775

Real-space multiple-scattering Hubbard model calculations for d- and f-state materials

CROSSMARK_Color_square_no_text.svg

aDepartment of Physics, University of Washington, Seattle, WA 98195, USA, and dTheoretical Division, Los Alamos National Laboratory, Los Alamos, NM 87545, USA
*Correspondence e-mail: vorwerk@physik.hu-berlin.de

Edited by S. M. Heald, Argonne National Laboratory, USA (Received 15 December 2014; accepted 19 May 2015; online 27 June 2015)

Calculations are presented of the electronic structure and X-ray spectra of materials with correlated d- and f-electron states based on the Hubbard model, a real-space multiple-scattering formalism and a rotationally invariant local density approximation. Values of the Hubbard parameter are calculated ab initio using the constrained random-phase approximation. The combination of the real-space Green's function with Hubbard model corrections provides an efficient approach to describe localized correlated electron states in these systems, and their effect on core-level X-ray spectra. Results are presented for the projected density of states and X-ray absorption spectra for transition metal- and lanthanide-oxides. Results are found to be in good agreement with experiment.

1. Introduction

Many-body perturbation theory and the GW approximation starting from density functional theory (DFT) calculations or the real-space multiple-scattering (RSMS) Green's function theory have become popular approaches to calculating the excited state electronic structure and spectra of condensed matter systems (Rehr & Albers, 2000[Rehr, J. J. & Albers, R. C. (2000). Rev. Mod. Phys. 72, 621-654.]; Onida et al., 2002[Onida, G., Reining, L. & Rubio, A. (2002). Rev. Mod. Phys. 74, 601-659.]) in recent decades. In particular, this approach has been widely successful in describing the electronic structure of weakly interacting sp bonded system. On the other hand, this formalism often falls short in the description of localized correlated 3d and 4f systems. In these cases the single-electron DFT or quasi-particle spectra are often in striking disagreement with spectroscopic data (Fabris et al., 2005[Fabris, S., de Gironcoli, S., Baroni, S., Vicario, G. & Balducci, G. (2005). Phys. Rev. B, 71, 041102.]). A straightforward approach to addressing these shortcomings is based on adding Hubbard model corrections, e.g. within the local density approximation (LDA) (Anisimov et al., 1991[Anisimov, V. I., Zaanen, J. & Andersen, O. K. (1991). Phys. Rev. B, 44, 943-954.], 1993[Anisimov, V. I., Solovyev, I. V., Korotin, M. A., Czyżyk, M. T. & Sawatzky, G. A. (1993). Phys. Rev. B, 48, 16929-16934.], 1997[Anisimov, V. I., Aryasetiawan, F. & Lichtenstein, A. I. (1997). J. Phys. Condens. Matter, 9, 767-808.]). This approach adds a spin- and orbital-dependent potential parameterized by two parameters, U and J, to the LDA-Hamiltonian to account for the screened Coulomb interaction between these localized states. Though highly simplified, the model accounts for the gap corrections observed in the density of states of these materials. Besides electronic structure, related spectra that include such Hubbard corrections are needed to interpret experimental data. To this end, Ahmed et al. (2012[Ahmed, T., Kas, J. J. & Rehr, J. J. (2012). Phys. Rev. B, 85, 165123.]) established a proof of principle for an LDA+U extension of the RSMS approach. This real-space Green's function method is now widely used to calculate X-ray spectra, for example in the FEFF codes (Rehr & Albers, 2000[Rehr, J. J. & Albers, R. C. (2000). Rev. Mod. Phys. 72, 621-654.]). The implementation of Hubbard corrections in this framework was, however, limited to d-state systems.

In this paper, we generalize the approach to both d- and f-electron systems, together with an automated way of calculating the Hubbard parameters U. The approach is implemented within the modular RSMS code FEFF9 (Rehr et al., 2010[Rehr, J. J., Kas, J. J., Vila, F. D., Prange, M. P. & Jorissen, K. (2010). Phys. Chem. Chem. Phys. 12, 5503-5513.]). The FEFF9 code is highly automated and includes many-body effects such as the electron self-energy and Debye–Waller factors. The efficiency of this code has been demonstrated by calculations for a wide range of materials. Moreover, the approach builds in all-electron full-relativistic initial states and self-consistent semi-relativistic final states, and does not depend on crystal or point-group symmetry. Consequently, the code is applicable to aperiodic materials such as nanostructures, biological materials or crystals with impurities throughout the periodic table. For these reasons, this method is both widely applicable and user-friendly compared with many other codes. A major approximation is the reliance of spherical scattering potentials, but, as seen in the examples presented here, is not a serious limitation.

The main product of the current investigation is a generalization of FEFF9 for d- and f-electron systems with new modules that calculate Hubbard model corrections in the electronic structure and X-ray absorption spectrum (XAS). In particular, the Hubbard model parameters are calculated within the constrained random phase approximation approach using a generalization of that introduced by Ahmed et al. (2012[Ahmed, T., Kas, J. J. & Rehr, J. J. (2012). Phys. Rev. B, 85, 165123.]). This procedure is in contrast to many works where the Hubbard parameters are either estimated or fitted to experiment. The new code is tested on both d- and f-state materials including transition-metal and lanthanide oxides. Compared with similar calculations without the Hubbard corrections, we find that this approach greatly improves the agreement with experiment from direct [X-ray photoelectron spectroscopy (XPS)] and inverse [bremsstrahlung isochromat spectroscopy (BIS)] photoemission spectroscopy for these materials. In the remainder of this paper, the theory and implementation are discussed in §2[link], results in §3[link], and conclusions in §4[link].

2. Theory and implementation

In this section, we briefly review the theoretical approach used to calculate the Hubbard potential and the implementation in the RSMS Green's function formalism.

2.1. RSMS Green's function approach

X-ray absorption spectroscopy measures the X-ray absorption coefficient [\mu(E)] which characterizes the exponential decay of the intensity of an X-ray beam passing through a material. In the RSMS approach, [\mu(E)] and other physical quantities are expressed in terms of the quasiparticle Green's function [G({\bf r},{\bf r}',E)]. This is in contrast to conventional approaches based on wavefunctions, and greatly simplifies the calculations. In this approach the absorption coefficient μ(E) is given by

[\mu(E) \propto -{{2} \over {\pi}}\,{\rm Im} \langle \varphi_{\rm{c}}({\bf r}) | \hat{\varepsilon}\cdot {\bf r}G({\bf r},{\bf r}',E+E_{\rm{c}})\hat{\varepsilon}\cdot{\bf r}' | \varphi_{\rm{c}}({\bf r}') \rangle,\eqno(1)]

which is formally equivalent to the Fermi golden rule. This can be seen formally by inserting the spectral decomposition of the Green's function,

[G({\bf r},{\bf r}',E)= \sum\limits_k \varphi_k({\bf r}) \, \varphi_k^*({\bf{r}}')\,\theta(\varepsilon_k-\varepsilon_F)/(E-\varepsilon_k + i \eta),]

into equation (1)[link] and taking the imaginary part. Here [|\varphi_{\rm{c}}({\bf r})\rangle] is the core state wavefunction, Ec the core level energy and [\hat{\varepsilon}] the photon polarization.

The spin and angular momentum projected density of states [\rho_{l\sigma}^{(n)}(E)] at atomic position n, which is of special interest in this paper, can be expressed in terms of the Green's function,

[\rho_{l\sigma}^{(n)}(E) = -{{1} \over {\pi}}\sum_{m}\textstyle\int\limits_{0}^{R_n} G_{Ln,L'n}^{\,\sigma,\sigma'}(r,r,E)\,r^2\,{\rm{d}}r,\eqno(2)]

with Rn the Norman radius about atomic position n. This can be derived, again using the spectral decomposition of G and recognizing that the imaginary part of [G({\bf r},{\bf r}',E)] is the density-matrix. [G_{Ln,L'n'}^{\,\sigma,\sigma'}] are the coefficients of the expansion of the Green's function about sites n and [n'] in terms of spherical harmonics,

[G({\bf r},{\bf r}',E) = \sum_{L,L',\sigma}Y_L(\hat{{\bf r}})\,G_{Ln,L'n'}^{\,\sigma, \sigma'}Y^{*}_{L}(\hat{{\bf r}}'),\eqno(3)]

where L = (l,m) are the orbital and azimuthal quantum numbers. Analogously we obtain the spin-dependent density matrix for site n as

[n^{\sigma\sigma'}_{Ln,L'n'}= -{{1} \over {\pi}}\int\limits^{E_f}{\rm{d}}E\int\limits_{\rm{cell}} {\rm Im}\, G^{\,\sigma\sigma'}_{Ln,L'n'}({\bf r},{\bf r},E)\,{\rm{d}}{\bf r}.\eqno(4)]

In the RSMS approach the Green's function G is naturally separated in the multiple-scattering series in contributions from the central (absorbing) atom G c and multiple-scattering (MS) contributions from the environment G sc, i.e. G = G c+G sc. The scattering part can be expressed as a sum over all MS paths that the photoelectron can take (Rehr et al., 2009[Rehr, J. J., Kas, J. J., Prange, M. P., Sorini, A. P., Takimoto, Y. & Vila, F. (2009). C. R. Phys. 10, 548-559.], 2010[Rehr, J. J., Kas, J. J., Vila, F. D., Prange, M. P. & Jorissen, K. (2010). Phys. Chem. Chem. Phys. 12, 5503-5513.]), i.e.

[G^{\,\rm{sc}}_{L,L'} = \left [\overline{G}^{\,0}T\overline{G}^{\,0}+\overline{G}^{\,0}T\overline{G}^{\,0}T\overline{G}^{\,0}+\ldots\right]_{L,L'},\eqno(5)]

where the successive terms represent the different orders of the scattering process. In this equation the scattering matrix T and the bare Green's function [\overline{G}^{\,0}] are represented in the angular-momentum and site basis set: [T_{ Ln,L'n'}] = [t_{Ln}\delta_{L,L'}\delta_{n,n'}] with the single-site scattering t matrix tLn and [\overline{G}^{\,0}] = [G^0_{Ln,L'n'}(1-\delta_{nn'})].

In the near-edge region the MS expansion is more efficiently carried out to full order [full multiple-scattering (FMS) expansion] by explicit matrix inversion,

[G^{\,\rm{sc}}_{L,L'}= \left[\left(1-\overline{G}^{\,0}T\right)^{-1}\,\overline{G}^{\,0}\right]_{L,L'}.\eqno(6)]

This can be seen by expanding the matrix inverse in powers of the scattering T-matrix. For more details on the RSMS Green's function, see Rehr et al. (2009[Rehr, J. J., Kas, J. J., Prange, M. P., Sorini, A. P., Takimoto, Y. & Vila, F. (2009). C. R. Phys. 10, 548-559.]).

2.2. Hubbard correction

The construction of the Hubbard model potential [V_{lm\sigma}^{\,U}(E)] in our approach directly follows the LDA+U formalism of Anisimov et al. (1991[Anisimov, V. I., Zaanen, J. & Andersen, O. K. (1991). Phys. Rev. B, 44, 943-954.], 1993[Anisimov, V. I., Solovyev, I. V., Korotin, M. A., Czyżyk, M. T. & Sawatzky, G. A. (1993). Phys. Rev. B, 48, 16929-16934.]). In this approach, a Hubbard correction E U, which takes the Coulomb interaction between localized correlated electrons into account, is added to the LDA energy functional E LDA[n] of the system:

[E\left[n^{\sigma}({\bf r}),{\bf n}^{\sigma}\right] = E^{\,\rm{LDA}}\left[n^{\sigma}({\bf r})\right]+E^{\,U}\left[{\bf n}^{\sigma}\right]-E_{\rm{dc}}\left[{\bf n}^{\sigma}\right],\eqno(7)]

with the charge density [n^{\sigma}({\bf r})], the density matrix [{\bf n}^{\sigma}] and the double-counting term Edc.

In order to obtain a rotationally invariant implementation, that takes the full density matrix into account, we perform a unitary transformation to diagonalize the density matrix. This yields a new basis [|l\alpha\rangle] for all states to which the Hubbard correction is applied. This approach is equivalent to a basis-independent LDA+U approach (Liechtenstein et al., 1995[Liechtenstein, A. I., Anisimov, V. I. & Zaanen, J. (1995). Phys. Rev. B, 52, R5467-R5470.]).

We write the total energy as

[\eqalignno{E = {}& E^{\,\rm{LDA}}+{{1}\over{2}} \sum\limits_{\alpha,\alpha',\sigma} U\left(n_{\alpha}^{-\sigma}-n^{0}\right) \left(n_{\alpha'}^{-\sigma}-n^{0}\right) \cr& +{{1}\over{2}} \sum\limits_{\alpha,\sigma'\,\neq\,\alpha,\sigma} (U-J)\left(n_{\alpha}^{\sigma}-n^{0}\right)\left(n_{\alpha'}^{\sigma}-n^{0}\right).&(8)}]

Here, the double-counting term is represented by the average occupation n0 = [({1}/{4l+2})\textstyle\sum_{\alpha\sigma}n^{\sigma}_{\alpha}]. From the energy functional the potential correction [V_{l\alpha\sigma}] to the LDA potential can be obtained, i.e.

[V_{l\alpha\sigma}^{\,U}= U\sum\limits_{\alpha'}\left(n_{\alpha'}^{-\sigma}-n^{0}\right) +(U-J)\sum\limits_{\alpha'\,\neq\,\alpha}\left(n_{\alpha'}^{-\sigma}-n^{0}\right).\eqno(9)]

2.3. Constrained random phase approximation (cRPA)

The cRPA calculation (Aryasetiawan et al., 2006[Aryasetiawan, F., Karlsson, K., Jepsen, O. & Schönberger, U. (2006). Phys. Rev. B, 74, 125106.]; Nilsson et al., 2013[Nilsson, F., Sakuma, R. & Aryasetiawan, F. (2013). Phys. Rev. B, 88, 125123.]) of the Hubbard parameter U starts by separating the non-interacting polarization operator into two contributions: P = Pl+Pr. Here Pl represents the contribution to the polarization due to the transitions between localized states to which the Hubbard corrections are applied, e.g. 3d or 4f states, and Pr is due to the remaining states. This separation is possible if the corresponding band is narrow, as is typical for strongly correlated materials. The effective Coulomb interaction [U({\bf r}, {\bf r}', \omega)] in the narrow band can then be calculated using the relation

[U({\bf r},{\bf r}',\omega) = \left [1-\upsilon P_{\rm{r}}({\bf r}, {\bf r}', \omega)\right]^{-1}\upsilon\eqno(10)]

where the polarization is given by

[\eqalignno{P_{\rm{r}} ={}& \sum_{i}^{\rm occ}\sum_{j}^{\rm unocc} \psi_i({\bf r})\,\psi_i^{*}({\bf r}')\,\psi_j^{*}({\bf r})\, \psi_j({\bf r}') \cr& \times \left [{{1} \over {\omega-\varepsilon_j+\varepsilon_i+i0^{+}}}+{{1} \over {\omega+\varepsilon_j-\varepsilon_i-i0^{+}}}\right].&(11)}]

Here, neither sum includes the narrow correlated band of states. The Hubbard parameter is then calculated by averaging the interaction over the atomic site n,

[U = \int_{0}^{R_n} {\rm{d}}{\bf r}\,{\rm{d}}{\bf r}' \,|\psi_{\rm{l}}({\bf r})|^2U({\bf r},{\bf r}',\omega = 0)|\psi_{\rm{l}}({\bf r}')|^2. \eqno(12)]

In the RSMS Green's function approach, we can express the zero frequency polarization operator in terms of the Green's function as

[\eqalignno{ P_{\rm{r}}(r,r^{\,\prime},\omega=0) = {}& -2{\rm Im} \int\limits_{-\infty}^{E_{\rm{f}}} {{{\rm{d}}\omega}\over{2\pi}} \cr& \times \Big[ \sum_{L\,\neq\,L_{\rm{H}}} G^{\,+}_{L,L'}(r,r^{\,\prime},\omega)\,G^{\,+}_{L,L'}(r,r^{\,\prime},\omega) \cr& +G^{\,+}_{L_{\rm{H}},L_{\rm{H}}}(r,r^{\,\prime},\omega)\,G^{\,+}_{L_{\rm{H}},L_{\rm{H}}}(r,r^{\,\prime},\omega) \cr& \times \Theta(r-R_{\rm{c}})\,\Theta(r^{\,\prime}-R_{\rm{c}}) \Big], &(13)}]

where LH = [|l\alpha\rangle] denotes the states the Hubbard correction is applied to, and [\Theta(r)] is a smooth function which goes to zero at r = Rc. This expression contains two approximations: we have assumed that it is sufficient to use a spherical average of the Coulomb interaction about a given site, since the correlated states are limited in spatial extent; we have also neglected off-diagonal angular-momentum elements of the Green's function.

2.4. Implementation

In our implementation we use a single-step LSDA calculation with the Barth–Hedin LSDA functional to obtain the occupation numbers [n_{lm}^{\sigma}]. For this calculation the spin polarization for the ith atom mi = [n_{i}^{\uparrow}-n_{i}^{\downarrow}] has to be known. We use the Hund's multiplicity rule for free atoms to determine the spin polarization. For a rotational invariant formalism, we diagonalize the occupation matrix [n_{l_{\rm{H}},m,m'}] [as in equation (4)[link]] and express it in the diagonal basis [|l_{\rm{H}}\alpha\rangle]. Here, lH is the orbital quantum number for which the Hubbard corrections is applied. Thus, we construct the unitary transformation matrix [\tau_{lm\alpha}^{\,\sigma}] from the eigenvalues of the occupation matrix [n_{mm'}]. The corrected self-consistent potential yields new solutions [R_{l\alpha}(E)] and [H_{l\alpha}(E)] of the Schroedinger equation within the muffin-tin spheres. Orbital dependent phase shifts [\delta_{l\alpha}^{\sigma}(E)] are found by matching these solutions with spherical Bessel functions at the sphere boundaries. From these phase shifts the scattering matrices are constructed as [t_{l\alpha}^{\,\sigma}] = [\exp(\delta_{l\alpha}^{\sigma})\sin(\delta_{l\alpha}^{\sigma})]. The scattering matrices are transformed back using the unitary transformation matrix to obtain tlm. Finally the Green's function is calculated using equation (5)[link], from which we can calculate the angular-momentum projected density of states, the X-ray absorption spectrum and other physical quantities. With the above implementation, Hubbard model calculations in FEFF9 are straightforward to make. They are governed by two keywords in the FEFF9 input file feff.inp. The Hubbard U parameter can be calculated automatically within the code using the cRPA approach (§2.3[link]), or a known value can be used. The FEFF9 Hubbard calculations typically take about 2–4 times longer than regular non-spin-polarized FEFF9 calculations due to the need to explicitly treat the spin variables and to calculate the DOS twice. However, as this new implementation is parallelized and scales near-linearly up to about 100 cores even on slow networks, the calculations are still very efficient and much faster than the original d-state version, as well as being more user-friendly. For the cRPA calculations, the Rcut parameter which describes the localization of the Hubbard correlated orbitals (§2.3[link]) needs to be chosen appropriately. Since the localization has been shown to depend only weakly on the crystal structure (Nilsson et al., 2013[Nilsson, F., Sakuma, R. & Aryasetiawan, F. (2013). Phys. Rev. B, 88, 125123.]), we expect elements with similar electronic structure to have the same Rcut. Differences in the Hubbard parameter U are due to different screening of the Coulomb interaction by the extended s- and p-states. This approximation allows us to compare the correlation in a group of elements, e.g. the lanthanides. We find that this approximation yields reasonable results over a wide group of elements and crystal structures.

3. Results

3.1. cRPA calculations

Calculations of the Hubbard parameter U were carried out using the procedure of §2.3[link]. The calculation of U depends on the radial cut-off Rc in equation (13)[link]. We choose a cut-off that yields a density of states (DOS) in good agreement with experiment, following the method used for d-state materials. This cut-off is then used to calculate U for a number of lanthanide oxides in the La2O3 structure (space group [P\bar3m1]) as well as other lanthanides in the close-packed hexagonal structure (space group P63mmc).

We find Rc = 1.49Rn which yields a Hubbard parameter of U = 6 eV for CeO2. This Hubbard parameter yields a DOS in good agreement with experimental results. It also agrees with Hubbard parameters used in LDA+U calculations (Fabris et al., 2005[Fabris, S., de Gironcoli, S., Baroni, S., Vicario, G. & Balducci, G. (2005). Phys. Rev. B, 71, 041102.]; Castleton et al., 2007[Castleton, C. W. M., Kullgren, J. & Hermansson, K. (2007). J. Chem. Phys. 127, 244704.], and reference therein) of U ≃ 5–7 eV. The DOS depends only weakly on the Hubbard parameters in this range (Table 1[link]).

Table 1
cRPA U-values compared with previously calculated values from Fabris et al. (2005[Fabris, S., de Gironcoli, S., Baroni, S., Vicario, G. & Balducci, G. (2005). Phys. Rev. B, 71, 041102.]) [1] and Castleton et al. (2007[Castleton, C. W. M., Kullgren, J. & Hermansson, K. (2007). J. Chem. Phys. 127, 244704.]) [2], and experimentally obtained values from Herbst & Wilkins (1987[Herbst, J. & Wilkins, J. (1987). High Energy Spectroscopy, edited by L. E. Karl A. Gschneidner Jr. and S. Hüfner, Vol. 10 of Handbook on the Physics and Chemistry of Rare Earths, pp. 321-360. Amsterdam: Elsevier.]) [3]

Material U (eV) Uref (eV) [1] [2] Uexp (eV) [3]
CeO2 6.0 5–7
Ce2O3 4.9
La2O3 4.1
αPr 6.7 4.8 5.4
αNd 7.3 4.8 6.3

We also investigate the dependence of the cRPA U-value on the unit cell volume for the example of CeO2. The crystal structure was not changed for the calculations. Results of these calculations are shown in Fig. 1[link]. The effect of the lattice constant can be seen in equation (11)[link]. With decreased lattice parameter, the d-states become less localized, increasing the transition matrix element in the polarization. On the other hand, the energy bands are pushed further apart by a decrease of the lattice constant, therefore diminishing the polarization (Tomczak et al., 2009[Tomczak, J. M., Miyake, T., Sakuma, R. & Aryasetiawan, F. (2009). Phys. Rev. B, 79, 235133.]). The total change of the Hubbard parameter is therefore material-dependent. A decrease of the Hubbard parameter with decreased lattice parameter has also been found for light actinides (Ahmed et al., 2014[Ahmed, T., Albers, R. C., Balatsky, A. V., Friedrich, C. & Zhu, J.-X. (2014). Phys. Rev. B, 89, 035104.]).

[Figure 1]
Figure 1
Dependence of the cRPA U on the unit cell volume. The line is inserted to guide the eye.

3.2. NiO

Hubbard corrections for transition metals have been applied with some success (Ahmed et al., 2012[Ahmed, T., Kas, J. J. & Rehr, J. J. (2012). Phys. Rev. B, 85, 165123.]; Jiang et al., 2010[Jiang, H., Gomez-Abal, R. I., Rinke, P. & Scheffler, M. (2010). Phys. Rev. B, 82, 045108.]) to correct the description of the localized d-states. As a first test of our generalized method, we calculate the DOS and XAS of NiO. We used a distorted crystal structure with a = b = 4.168 Å, c = 4.166 Å and [\alpha] = [\beta] = 90.055° and [\gamma] = 90.082°. The Hubbard parameters U = 8.0 eV and J = 0.9 eV were used. The Hubbard parameter U was calculated by the cRPA algorithm (cf. §3.1[link]).

Fig. 2[link] shows the DOS for NiO, calculated by FEFF9 and our FEFF9+U code. We also compare with experimental results (Zimmermann et al., 1999[Zimmermann, R., Steiner, P., Claessen, R., Reinert, F., Hüfner, S., Blaha, P. & Dufek, P. (1999). J. Phys. Condens. Matter, 11, 1657-1682.]) obtained from XPS–BIS measurements. It can be seen that the Hubbard correction opens a band gap, which does not occur in a FEFF9 calculation. Also the DOS in the valence band is broadened by the Hubbard correction, in good agreement with the experimental result. The Hubbard model still underestimates the width of the band gap. The cRPA U-value is a physically meaningful quantity, which can be considered as the on-site Coulomb interaction of the d-electrons in the static limit. In the ground state of transition metal oxides (e.g. NiO, MnO) with partially unoccupied d-orbitals, the Hubbard correction plays a major role in estimating the band gap. However, due to the presence of oxygen p-like states near the Fermi energy, and their non-trivial hybridization with localized d-states, one must consider the contribution from non-local dynamical self-energy corrections. Within the scope of FEFF, this is done using the GW self-energy formalism with a many-pole screening model (Kas et al., 2007[Kas, J. J., Sorini, A. P., Prange, M. P., Cambell, L. W., Soininen, J. A. & Rehr, J. J. (2007). Phys. Rev. B, 76, 195116.]). This frequency-dependent GW correction can account for the additional 1 eV band-gap correction on top of the Hubbard correction with cRPA U = 8.0 eV, as was earlier demonstrated in the O K-edge X-ray absorption near-edge structure (XANES) and X-ray emission spectroscopy (XES) for NiO and MnO (Ahmed et al., 2012[Ahmed, T., Kas, J. J. & Rehr, J. J. (2012). Phys. Rev. B, 85, 165123.]).

[Figure 2]
Figure 2
Total DOS of NiO calculated with the FEFF9 and FEFF9+U implementation compared with experimental results obtained from XPS–BIS measurements. Experimental results are extracted from Zimmermann et al. (1999[Zimmermann, R., Steiner, P., Claessen, R., Reinert, F., Hüfner, S., Blaha, P. & Dufek, P. (1999). J. Phys. Condens. Matter, 11, 1657-1682.]).

Our FEFF9+U implementation increases the agreement with experimental results (Kurmaev et al., 2008[Kurmaev, E. Z., Wilks, R. G., Moewes, A., Finkelstein, L. D., Shamin, S. N. & Kuneš, J. (2008). Phys. Rev. B, 77, 165127.]) for the O K-edge XAS of NiO (Fig. 3[link]). The Hubbard implementation yields a better description of the two peaks below 545 eV. These peaks are attributed to the O p-states, which are strongly hybridized with Ni d-states (Ahmed et al., 2012[Ahmed, T., Kas, J. J. & Rehr, J. J. (2012). Phys. Rev. B, 85, 165123.]). A better description of the correlated Ni d-states therefore yields the good agreement with experimental results.

[Figure 3]
Figure 3
NiO O K-edge XANES calculated with FEFF9 and FEFF9+U implementation compared with experimental results. Spectra are aligned at the main peak. Vertical lines are inserted at the experimental main peaks to guide the eye. Experimental results are extracted from Kurmaev et al. (2008[Kurmaev, E. Z., Wilks, R. G., Moewes, A., Finkelstein, L. D., Shamin, S. N. & Kuneš, J. (2008). Phys. Rev. B, 77, 165127.]) and all spectra are aligned at the second main peak.

3.3. CeO2

As a first example of an f material, we calculate the properties for CeO2 in the cubic fluorite lattice (space group [Fm\bar3m]) (Skorodumova et al., 2001[Skorodumova, N. V., Ahuja, R., Simak, S. I., Abrikosov, I. A., Johansson, B. & Lundqvist, B. I. (2001). Phys. Rev. B, 64, 115108.]) using the experimental lattice constant of 5.411 Å (Gschneidner, 2006[Gschneidner, K. A. (2006). Handbook on the Physics and Chemistry of Rare Earths, Vol. 36. San Diego: Elsevier.]). Converged RSMS results are obtained using a cluster of 137 atoms around the absorbing atom for the FMS calculation and a smaller cluster of 43 atoms for the calculation of the self-consistent muffin-tin potentials. For Hubbard calculations the calculated value U = 6.0 eV was used and J = 0.3 eV (Nilsson et al., 2013[Nilsson, F., Sakuma, R. & Aryasetiawan, F. (2013). Phys. Rev. B, 88, 125123.]). The oxygen K-edge XAS was calculated as well as angular-momentum projected density of states (LDOS).

Fig. 4[link] shows the total density of states of CeO2. The experimental spectrum (Wuilloud et al., 1984[Wuilloud, E., Delley, B., Schneider, W. D. & Baer, Y. (1984). Phys. Rev. Lett. 53, 202-205.]) is aligned at the upper valence band edge. The localized un­occupied f-states of Ce yields a sharp peak above the Fermi energy in the total density of CeO2. The valence band is mainly formed by the O p-states. In good agreement with experimental results (Wuilloud et al., 1984[Wuilloud, E., Delley, B., Schneider, W. D. & Baer, Y. (1984). Phys. Rev. Lett. 53, 202-205.]), the FEFF9+U calculation predicts a sharp peak with f-state character at the onset of the conduction band, as well as a less intense occupied f-state peak in the valence band. The f-state character of both peaks can also be seen in the angular-momentum projected DOS in Fig. 5[link]. The most important quantity to describe the influence of the Hubbard corrections is therefore the pf energy difference [\Delta{E}_{p-f}]. Experimentally, this gap is found to be [\Delta E_{p-f}^{\,\rm{exp}}] = 6.1 eV. Although the FEFF9 calculation does not yield a band gap between the valence states and the unoccupied f peak, the addition of Hubbard corrections increases the energy difference between the valence band and the unoccupied f peak. This difference [\Delta E_{p-f}] (Fig. 4[link]) is underestimated in calculations without the Hubbard correction by approximately 2 eV. The addition of Hubbard corrections increases the difference to [\Delta E_{p-f}^{\,U}] = 6.0 eV. Previous DFT calculations underestimate the pf band gap, even with the addition of Hubbard corrections, by more than 50% (Da Silva et al., 2007[Da Silva, J. L. F., Ganduglia-Pirovano, M. V., Sauer, J., Bayer, V. & Kresse, G. (2007). Phys. Rev. B, 75, 045121.]).

[Figure 4]
Figure 4
Total DOS of CeO2 compared with the experimental spectrum obtained from XPS–BIS measurements. The experimental spectrum is extracted from Wuilloud et al. (1984[Wuilloud, E., Delley, B., Schneider, W. D. & Baer, Y. (1984). Phys. Rev. Lett. 53, 202-205.]).
[Figure 5]
Figure 5
Total DOS for the Ce f-orbital states. Solid lines show spin-up DOS, dashed lines the spin-down DOS.

Fig. 5[link] shows the angular-momentum projected DOS for the Ce f-states. As expected (Anisimov et al., 1997[Anisimov, V. I., Aryasetiawan, F. & Lichtenstein, A. I. (1997). J. Phys. Condens. Matter, 9, 767-808.]), a gap opens up between the occupied and unoccupied states when the Hubbard corrections are applied. The spectrum varies only slowly with respect to the Hubbard parameter. Fig. 6[link] shows the dependence of the pf difference on the Hubbard parameter U. For all calculations J = 0.3 eV was used. As expected, the band gap is increasing linearly with the Hubbard parameter and the slope of this increase is smaller than 1, since the Ce f orbitals are partially filled.

[Figure 6]
Figure 6
Dependence of the pf energy difference [\Delta E_{p-f}] on the Hubbard parameter U. The line is inserted to guide the eye.

Fig. 7[link] shows the oxygen K-edge absorption spectrum for CeO2. The addition of a Hubbard potential shifts the first excitation to higher energies due to the increased gap between the conduction band and the unoccupied f peak. The Hubbard correction also introduces an absorption peak at 535 eV in good agreement with experimental results (Douillard et al., 1994[Douillard, L., Gautier, M., Thromat, N., Henriot, M., Guittet, M. J., Duraud, J. P. & Tourillon, G. (1994). Phys. Rev. B, 49, 16171-16180.]), that is missing in a FEFF9 calculation. The absorption at the main peak at 542 eV and at higher energies is not changed by the corrections.

[Figure 7]
Figure 7
CeO2 O K-edge XAS comparing the FEFF9+U with the FEFF9 calculation and experimental result extracted from Douillard et al. (1994[Douillard, L., Gautier, M., Thromat, N., Henriot, M., Guittet, M. J., Duraud, J. P. & Tourillon, G. (1994). Phys. Rev. B, 49, 16171-16180.]). A, B and C denote the main experimental peaks.

This is in good agreement with previous studies of CeO2 oxygen K-edge absorption (Douillard et al., 1994[Douillard, L., Gautier, M., Thromat, N., Henriot, M., Guittet, M. J., Duraud, J. P. & Tourillon, G. (1994). Phys. Rev. B, 49, 16171-16180.]). The peaks B at 538 eV and C at 542 eV represent transition of the O 1s states into O 2p states hybridized with Ce d- and p-orbitals, respectively. Peak A at 535 eV is interpreted as transitions into O 2p states hybridized with the Ce 4f states. Therefore the Hubbard corrected spectrum shows that these three distinct peaks and changes with respect to the FEFF9 spectrum occur mainly at the peak A.

3.4. Ce2O3

XAS and DOS were calculated for Ce2O3 in the hexagonal lattice (space group [P\bar3m1]) with lattice constants of a = 7.35 Bohr and c = 11.45 Bohr (Bärnighausen & Schiller, 1985[Bärnighausen, H. & Schiller, G. (1985). J. Less-Common Metals, 110, 385-390.]). The Hubbard parameter U = 4.9 eV was calculated by a cRPA calculation and we found the best agreement with experimental results for J = 1.0 eV. In the O K-edge XAS in Fig. 8[link], the FEFF9 calculation displays an additional peak at low energies, that does not occur in the experimental spectrum (Garvie & Buseck, 1999[Garvie, L. & Buseck, P. (1999). J. Phys. Chem. Solids, 60, 1943-1947.]). This additional peak is due to the underestimation of the Ce f-state energies in the FEFF9 calculation (Nishida et al., 2009[Nishida, I., Tatsumi, K. & Muto, S. (2009). Mater. Trans. 50, 952-958.]). Similar to the spectrum of CeO2 the spectrum shows the three peaks A at 531 eV, B at 534 eV and C at 536 eV, although the peak B has a greatly reduced oscillator strength compared with CeO2. The peak A at 535 eV once again represents transitions to the O p-states hybridized with Ce f-states and is only well described with the introduction of Hubbard corrections for the f-states. The Hubbard corrected spectrum is in good agreement with experimental results (Garvie & Buseck, 1999[Garvie, L. & Buseck, P. (1999). J. Phys. Chem. Solids, 60, 1943-1947.]).

[Figure 8]
Figure 8
Ce2O3 O K-edge XAS comparing the FEFF9+U with the FEFF9 calculation and experimental result extracted from Garvie & Buseck (1999[Garvie, L. & Buseck, P. (1999). J. Phys. Chem. Solids, 60, 1943-1947.]). Spectra are aligned at the main peak. A, B and C denote the main experimental peaks.

The correction of the Ce f-state energies can also be seen in the TDOS in Fig. 9[link]. FEFF9 and FEFF9+U calculations are compared with experimental results from XPS–BIS measurements (Allen, 1985[Allen, J. (1985). J. Magn. Magn. Mater. 47-48, 168-174.]). FEFF9 displays two distinct Ce f-state peaks around the Fermi level. The Hubbard corrections shift the unoccupied f-state peak to higher energies and the occupied peak to the O p-state valence band in good agreement with experimental results.

[Figure 9]
Figure 9
Total DOS of Ce2O3 compared with the experimental spectrum obtained from XPS–BIS measurements. The experimental spectrum is extracted from Allen (1985[Allen, J. (1985). J. Magn. Magn. Mater. 47-48, 168-174.]).

4. Conclusion

We have presented a generalized real-space multiple-scattering implementation of the Hubbard model for strongly correlated systems. This approach can treat localized states of d- or f-character. It can calculate the Hubbard U parameter automatically using a cRPA algorithm. Resulting values of U are in good agreement with other sources and with experiment. Our approach is included in the FEFF9 code for X-ray spectroscopy and yields density of states as well as many types of spectra. It is easy to use and efficient due to the RSMS formalism and efficient MPI parallelization. Results are presented for the f-state lanthanide CeO2 and Ce2O3, where we find good agreement with experiment.

Acknowledgements

This work was supported by DOE Grant DE-FG03-97ER45623 (JJR), Studienstiftung des deutschen Volkes (CV) and Deutscher Akademischer Austauschdienst (CV).

References

First citationAhmed, T., Albers, R. C., Balatsky, A. V., Friedrich, C. & Zhu, J.-X. (2014). Phys. Rev. B, 89, 035104.  CrossRef Google Scholar
First citationAhmed, T., Kas, J. J. & Rehr, J. J. (2012). Phys. Rev. B, 85, 165123.  CrossRef Google Scholar
First citationAllen, J. (1985). J. Magn. Magn. Mater. 4748, 168–174.  CrossRef CAS Google Scholar
First citationAnisimov, V. I., Aryasetiawan, F. & Lichtenstein, A. I. (1997). J. Phys. Condens. Matter, 9, 767–808.  CrossRef CAS Google Scholar
First citationAnisimov, V. I., Solovyev, I. V., Korotin, M. A., Czyżyk, M. T. & Sawatzky, G. A. (1993). Phys. Rev. B, 48, 16929–16934.  CrossRef CAS Google Scholar
First citationAnisimov, V. I., Zaanen, J. & Andersen, O. K. (1991). Phys. Rev. B, 44, 943–954.  CrossRef CAS Web of Science Google Scholar
First citationAryasetiawan, F., Karlsson, K., Jepsen, O. & Schönberger, U. (2006). Phys. Rev. B, 74, 125106.  CrossRef Google Scholar
First citationBärnighausen, H. & Schiller, G. (1985). J. Less-Common Metals, 110, 385–390.  Google Scholar
First citationCastleton, C. W. M., Kullgren, J. & Hermansson, K. (2007). J. Chem. Phys. 127, 244704.  CrossRef PubMed Google Scholar
First citationDa Silva, J. L. F., Ganduglia-Pirovano, M. V., Sauer, J., Bayer, V. & Kresse, G. (2007). Phys. Rev. B, 75, 045121.  Web of Science CrossRef Google Scholar
First citationDouillard, L., Gautier, M., Thromat, N., Henriot, M., Guittet, M. J., Duraud, J. P. & Tourillon, G. (1994). Phys. Rev. B, 49, 16171–16180.  CrossRef CAS Google Scholar
First citationFabris, S., de Gironcoli, S., Baroni, S., Vicario, G. & Balducci, G. (2005). Phys. Rev. B, 71, 041102.  CrossRef Google Scholar
First citationGarvie, L. & Buseck, P. (1999). J. Phys. Chem. Solids, 60, 1943–1947.  CrossRef CAS Google Scholar
First citationGschneidner, K. A. (2006). Handbook on the Physics and Chemistry of Rare Earths, Vol. 36. San Diego: Elsevier.  Google Scholar
First citationHerbst, J. & Wilkins, J. (1987). High Energy Spectroscopy, edited by L. E. Karl A. Gschneidner Jr. and S. Hüfner, Vol. 10 of Handbook on the Physics and Chemistry of Rare Earths, pp. 321–360. Amsterdam: Elsevier.  Google Scholar
First citationJiang, H., Gomez-Abal, R. I., Rinke, P. & Scheffler, M. (2010). Phys. Rev. B, 82, 045108.  CrossRef Google Scholar
First citationKas, J. J., Sorini, A. P., Prange, M. P., Cambell, L. W., Soininen, J. A. & Rehr, J. J. (2007). Phys. Rev. B, 76, 195116.  Web of Science CrossRef Google Scholar
First citationKurmaev, E. Z., Wilks, R. G., Moewes, A., Finkelstein, L. D., Shamin, S. N. & Kuneš, J. (2008). Phys. Rev. B, 77, 165127.  CrossRef Google Scholar
First citationLiechtenstein, A. I., Anisimov, V. I. & Zaanen, J. (1995). Phys. Rev. B, 52, R5467–R5470.  CrossRef CAS Google Scholar
First citationNilsson, F., Sakuma, R. & Aryasetiawan, F. (2013). Phys. Rev. B, 88, 125123.  CrossRef Google Scholar
First citationNishida, I., Tatsumi, K. & Muto, S. (2009). Mater. Trans. 50, 952–958.  CrossRef CAS Google Scholar
First citationOnida, G., Reining, L. & Rubio, A. (2002). Rev. Mod. Phys. 74, 601–659.  CrossRef CAS Google Scholar
First citationRehr, J. J. & Albers, R. C. (2000). Rev. Mod. Phys. 72, 621–654.  Web of Science CrossRef CAS Google Scholar
First citationRehr, J. J., Kas, J. J., Prange, M. P., Sorini, A. P., Takimoto, Y. & Vila, F. (2009). C. R. Phys. 10, 548–559.  Web of Science CrossRef CAS Google Scholar
First citationRehr, J. J., Kas, J. J., Vila, F. D., Prange, M. P. & Jorissen, K. (2010). Phys. Chem. Chem. Phys. 12, 5503–5513.  Web of Science CrossRef CAS PubMed Google Scholar
First citationSkorodumova, N. V., Ahuja, R., Simak, S. I., Abrikosov, I. A., Johansson, B. & Lundqvist, B. I. (2001). Phys. Rev. B, 64, 115108.  CrossRef Google Scholar
First citationTomczak, J. M., Miyake, T., Sakuma, R. & Aryasetiawan, F. (2009). Phys. Rev. B, 79, 235133.  CrossRef Google Scholar
First citationWuilloud, E., Delley, B., Schneider, W. D. & Baer, Y. (1984). Phys. Rev. Lett. 53, 202–205.  CrossRef CAS Google Scholar
First citationZimmermann, R., Steiner, P., Claessen, R., Reinert, F., Hüfner, S., Blaha, P. & Dufek, P. (1999). J. Phys. Condens. Matter, 11, 1657–1682.  CrossRef CAS Google Scholar

© International Union of Crystallography. Prior permission is not required to reproduce short quotations, tables and figures from this article, provided the original authors and source are cited. For more information, click here.

Journal logoJOURNAL OF
SYNCHROTRON
RADIATION
ISSN: 1600-5775
Follow J. Synchrotron Rad.
Sign up for e-alerts
Follow J. Synchrotron Rad. on Twitter
Follow us on facebook
Sign up for RSS feeds