X-ray absorption spectroscopy: state-of-the-art analysis

Natoli, C.R.; Benfatto, M.; Della Longa, S.; Hatada, K.

doi:10.1107/S0909049502017247

research papers

JOURNAL OF
SYNCHROTRON
RADIATION

ISSN: 1600-5775

Volume 10| Part 1| January 2003| Pages 26-42

doi:10.1107/S0909049502017247

X-ray absorption spectroscopy: state-of-the-art analysis

C. R. Natoli,^a ^* M. Benfatto,^a S. Della Longa ^b and K. Hatada ^a

^aLaboratori Nazionali di Frascati, INFN, CP13, 00044 Frascati, Italy, and ^bUniversita' dell'Aquila, Via Vetoio, loc. Coppito II, 67100 L'Aquila, Italy
^*Correspondence e-mail: natoli@lnf.infn.it

(Received 19 February 2002; accepted 29 August 2002)

State-of-the-art techniques for analysing X-ray absorption spectra are reviewed, with an eye to biological applications. Recent attempts to perform full spectral fitting of the XANES energy region and beyond for the purpose of structural analysis have met with encouraging success. The present paper analyses the theoretical motivations behind this success and indicates routes for future improvements. The theoretical background is not entirely new, although the point of view is, and some sections and appendices present material that the authors believe has never been published before. The aim of this paper is to provide a theoretical analysis that is as self-contained as possible.

Keywords: X-ray absorption spectra; XANES; biological applications.

1. Introduction

In the past 20 years, X-ray absorption spectroscopy (XAS) from inner-shell electrons has proved to be an invaluable tool in the study of the electronic and structural properties of condensed matter systems and, in particular, of the active sites of metalloproteins. In spectroscopic analyses, it has been common practice to treat separately the near-edge region (conventionally from below the edge up to ∼30–50 eV above it) and the high-energy part of the spectrum (above ∼30–50 eV), the so-called EXAFS (extended X-ray absorption fine structure) region. The motivation behind this separation is entirely empirical, in that the extraction of a structural `signal' from the absorption spectrum via subtraction of an approximate atomic background can be performed with a certain confidence and reliability only in the EXAFS region, whereas the low-energy part cannot be adequately background subtracted.

There exist several analysis packages based on spherical-wave multiple-scattering theory (MST) and the complex optical potential of the Hedin–Lundqvist type (Hedin & Lundqvist, 1969 , 1971 ). These packages are able to reproduce satisfactorily the EXAFS signal [for a review see Rehr & Albers (2000 )]. The same codes even offer the possibility of fitting the entire spectrum; however, it does not seem that this potential has been seriously exploited as yet. Some of the factors that have probably deterred the various practitioners of MST from extending the fitting procedure used in the EXAFS region to the edge region, via the calculation of the total cross section, are well known: the inadequacy of the muffin-tin approximation to the potential at low photoelectron energies, the lack of a satisfactory description of the screening and relaxation processes following the sudden creation of the core hole, and the need to include electronic correlations and, in particular, two-electron excitations (since they change the slope of the background absorption). Therefore, the analysis of the near-edge region has remained at a semi-quantitative level.

Nevertheless, quantitative analysis of the X-ray near-edge structure (XANES) in order to obtain structural and electronic information can be very relevant in many fields of scientific application, like extra-dilute systems, surface spectra, real-time measurements of dynamic systems, trace-element analysis, the local investigation of materials under extreme conditions and especially biological systems (enzymes), where the low S/N ratio and the weak scattering power of the light elements constituting the organic material limit the k-range of the available experimental data. In all these instances, the EXAFS region of the spectrum cannot be adequately exploited, since the usable data are very often below 200 eV from the absorption edge.

Very recently, a method for performing full spectral fitting of the XANES energy region and beyond (up to ∼200 eV above the edge) for the purpose of structural analysis has been proposed and applied to the K-edge spectra of a number of transition metal compounds, both organic and inorganic (Benfatto & Della Longa, 2001 ; Della Longa et al., 2001 ). The encouraging success of this attempt requires a reconsideration of the theory with a view to understanding why the method works and highlighting areas of possible improvement.

The present formulation reflects our long experience of the performance of multiple-scattering theory with a complex optical potential and might not duly acknowledge other groups' contributions to the subject. We apologize if this is the case.

2. Theoretical background

In this section we shall present the derivation of the photoemission cross section for the ejection of a photoelectron of final momentum $[{\bf k}]$ and kinetic energy k² along the direction $[\hat k]$ and illustrate the reduction of this many-body problem to an effective one-particle problem with complex energy-dependent optical potential. This process of reduction will help us to understand the validity of the necessary approximations to the optical potential and will give us guidance in choosing among various approximation schemes. The photoabsorption cross section is nothing more than the integration of the photoemission cross section over all the emission angles and all the final channels (elastic plus inelastic) with the same final energy. The reason we treat both cases together is duplex. Firstly, the mathematical formalism is the same, as is the reduction process to a one-particle problem; secondly, and more importantly, on purely physical grounds we can think of photoabsorption as a kind of photoemission that has the same electron source (the photoabsorber), in which the detector, instead of being outside the measured system, coincides with the source. The validity of this assumption will be apparent from the mathematical formalism in the following. By treating photoabsorption and photoemission together, we can judge the sensitivity to structural details of a particular potential in absorption by looking at its performance in photoelectron diffraction. In photoabsorption, because of the obvious impossibility of controlling the measured variables (except total energy), we are obliged to sum over all final states at a given energy.

2.1. Reduction of the many-body problem to an effective one-particle problem: the photoemission case

The photoemission cross section in the many-body case for the ejection of a photoelectron of final momentum $[\bf k]$ and kinetic energy k² along the direction $[\hat k]$ can be written in the dipole approximation as

$[\textstyle {\rm d}\sigma (\omega) / {\rm d}\hat k = 4\pi ^2 \alpha _0 \,\hbar \omega \left| {\left\langle {\Theta \Psi _{\bf k}^N } \right|{\boldvarepsilon} \cdot \sum\limits_{i = 1}^N {\bf r}_i \left| {\Psi_g^N } \,\right\rangle } \right|^2, \eqno(1)]$

where $[\Psi _{\bf k}^N]$ is the many-body final scattering state, normalized to one state per energy interval unit, for the N-electron system with one electron of momentum $[\bf k]$ travelling to infinity and $[\Psi _g^N]$ is its ground state; these states have respective energies E_k^N and E_g^N; $[\hbar \omega]$ is the incoming photon energy and $[\boldvarepsilon]$ is its polarization. Energy conservation imposes that $[\hbar \omega = E_k^N - E_g^N]$ . According to Breit & Bethe (1954 ), in order to satisfy the correct boundary conditions for the ejected photoelectron (no electron in a continuum state in the remote past), we must take the time-reversed scattering state by application of the time-reversal operator $[\Theta]$ . Throughout this paper we shall use atomic units of length and Rydberg units of energy. The quantity $[\alpha _0]$ is the fine-structure constant, equal to $[e^2 / \hbar c = 1 / 137]$ .

The photoabsorption cross section is obtained from the photoemission cross section by integrating over all the emission angles and all the final channels (elastic plus inelastic) with the same final energy E_k^N = k² + E_n^{N - 1}, where E_n^{N - 1} is the energy of the remaining N-1 electrons of the system. By definition, the elastic channel is the channel for which E_n^{N - 1} = E_g^{N - 1}, the ground state of the (N−1)-electron system. The states with energy E_n^{N - 1} can have one or more electrons in the continuum, as in the case of double photoionization (shake-off channels), or one excited electron in a valence excited state that remains in the system (shake-up channels).

Since we sum over all final states (f) we can eliminate the time-reversal operator and write the photoabsorption cross section in the usual way as

$[\sigma _{\rm abs} (\omega) = 4\pi ^2 \alpha _0 \,\hbar \omega \textstyle\sum\limits_f {} \left| {\left\langle \,{\Psi _f^N } \right|{\boldvarepsilon} \cdot \sum\limits_{i = 1}^N {\bf r}_i \,\left| {\Psi _g^N } \,\right\rangle } \right|^2 \delta (\hbar \omega - E_f^N + E_g^N).\eqno(2)]$

In the case of photoemission from a deep core state $[\varphi _{L_0 }^{\rm c}]$ of angular momentum L₀ = (l₀, m₀), we assume that, to a good approximation,

$[\eqalignno{ \Psi _g^N ({\bf r},{\bf r}_1,. \,.\,.\,,{\bf r}_{N - 1}) & = (N!)^{1/2} \, A\,\varphi _{L_0 }^{\rm c} ({\bf r})\textstyle\sum\limits_n {c_n } \Phi _n^{N - 1} ({\bf r}_1,.\,.\,.\,,{\bf r}_{N - 1}\,) \cr & = (N!)^{1/2} \, A\,\varphi _{L_0 }^{\rm c} ({\bf r})\,\Psi _g^{N - 1} ({\bf r}_1,. \,.\,.\,,{\bf r}_{N - 1}\,),& (3)}]$

where A is the usual anti-symmetrization operator [A = $[A(1 /N!)\textstyle\sum_P {(- 1)^P } P]$ , with A² = A], and $[\Phi _n^{N - 1} ({\bf r}_1,. \,.\,.\,,{\bf r}_{N - 1}\,)]$ are Slater determinants describing the configurations present in the ground state of the system. Normalization imposes $[\textstyle\sum_n {|c_n |^2 } = 1]$ if $[\langle {\varphi ^c} |{\varphi ^c } \rangle = 1]$ , and for simplicity we shall omit spin variables, since we are going to deal with non-magnetic systems.

In a similar way, we can write without loss of generality

$[\Psi _f^N ({\bf r},{\bf r}_1,. \,.\,.\,,{\bf r}_{N - 1}) = (N!)^{1/2}\, A \textstyle\sum\limits_\alpha {\varphi\, _\alpha ^f } ({\bf r})\,\tilde \Psi _n^{N - 1} ({\bf r}_1,. \,.\,.\,,{\bf r}_{N - 1}), \eqno(4)]$

where the functions $[\varphi \,_\alpha ^f]$ , ignoring exchange effects, can be thought of as describing the excited photoelectron while the $[\tilde \Psi _\alpha ^{N - 1}]$ states are eigenstates of the Hamiltonian H^{N - 1} describing the remaining (N−1)-electron system with eigenvalues $[E_\alpha ^{N - 1}]$ :

$[H^{N - 1} \,\tilde \Psi _\alpha ^{N - 1} = E_\alpha ^{N - 1} \,\tilde \Psi _\alpha ^{N - 1}. \eqno(5)]$

The tilde stands as a reminder that in the expansion (4) the relaxed states around the core hole are dominant. If needed, they can in turn be expanded in terms of Slater determinants that describe the intervening configurations in the final state. Borrowing the term from many-particle scattering theory, we can call the states $[\tilde \Psi _\alpha ^{N - 1}]$ final-state channels. Here and henceforth, the lower index f in the final state $[\Psi _f^N]$ can be replaced by $[\bf k]$ whenever we deal specifically with the scattering state $[\Psi _{\bf k}^N]$ .

The wave function $[\Psi _f^N]$ is an eigenstate of the total Hamiltonian H^N with eigenvalue $[E_f^N = E_g^N + \hbar \omega]$ , i.e.

$[H^N \Psi _{f}^{N} = {E}_{f}^{N} \Psi _{f}^{N}.\eqno(6)]$

Moreover,

$[H^N = - \nabla _{\bf r}^2 + \textstyle\sum\limits_i {V({\bf r},{\bf r}_i }) + H^{N - 1}, \eqno(7)]$

where $[V({\bf r},{\bf r}_i)]$ is the interaction potential of the excited photoelectron with the rest of the system.

By inserting (4) into (6), projecting onto the states $[\tilde \Psi _\alpha ^{N - 1}]$ and using (5), we obtain for the amplitude functions $[\varphi \,_\alpha ^f]$ the set of coupled equations

$[\textstyle(\nabla ^2 + k_\alpha ^2)\,\varphi\, _\alpha ^f ({\bf r}) = \sum\limits_\beta {\int {V_{\alpha \beta } } } ({\bf r},{\bf r}')\,\varphi\, _\beta ^f ({\bf r}')\,{\rm d}^3 r', \eqno(8)]$

where

$[k_\alpha ^2 = \hbar \omega - (E_g^{N - 1} - E_g^N) - (E_\alpha ^{N - 1} - E_g^{N - 1}) = \hbar \omega - I_{\rm c} - \Delta E_\alpha, \eqno(9)]$

I_c being the ionization potential for the core state and ΔE_α the excitation energy left behind in the (N − 1)-particle system. The non-local interchannel potentials $[V_{\alpha \beta } ({\bf r},{\bf r}')]$ are the matrix elements between states $[\tilde \Psi _\alpha ^{N - 1}]$ and $[\tilde \Psi _\beta ^{N - 1}]$ of the interaction potential $[V({\bf r},{\bf r}_i)]$ and include local terms coming from the Coulomb interaction as well as non-local exchange terms originating from the exchange interaction. The set of equations (8) is supplemented with the boundary conditions that are related to the behaviour of the photoelectron at infinity and to the state of the (N − 1)-electron system according to the partition of the total energy $[E_f^N = E_g^N + \hbar \omega]$ between the photoelectron and the system. Each different partition has a different set of boundary conditions, which lead to a different solution of (8). For example, if we are interested in a particular photoemission channel $[\beta]$ with kinetic energy $[k_\beta ^2]$ leaving behind the energy $[\Delta E_\beta]$ in the system, in the limit $[r \rightarrow \infty]$ we should impose the scattering conditions

$[\varphi _\alpha ({\bf r};{\bf k}_\beta) \simeq (k_\beta / \pi)^{1/2} \exp (i\,{\bf k}_\beta \cdot {\bf r}) \,\delta _{\alpha \beta } + f_\alpha (\hat r;{\bf k}_\beta)\exp(ik_\alpha r) / r, \eqno(10)]$

where we have made explicit the dependence of $[\varphi _\alpha]$ on $[k_\beta]$ as an argument rather than an upper index. Hence, as usual, δ_αβ is the Kronecker symbol, and f_α is the scattering amplitude. The factor $[(k_\beta /\pi)^{1/2}]$ normalizes the photoelectronic plane wave at the detector to one state per Rydberg. In the case of double electron transitions with two electrons in the continuum, (10) can easily be changed accordingly. For more details we refer the reader to Natoli et al. (1990 ), where a formal solution of (8) is given in the framework of the multichannel multiple-scattering theory. For the sake of the argument developed here, we need only the expression for the photoemission cross section, which is obtained from (1), (3) and (4) as

$[{\rm d}\sigma(\omega)/{\rm d}\hat k = 8\pi^2\alpha\hbar\omega\textstyle\sum\limits_{m_0}\left|\sum\limits_\alpha\left\langle S_\alpha^* \varphi _\alpha ^ - ({\bf r}\semi{\bf k}_\beta) \left| \boldvarepsilon \cdot {\bf r}\right| \varphi _{l_0 m_0 }^{\rm c} ({\bf r})\right\rangle\right|^2, \eqno(11)]$

which is valid if we orthogonalize the excited channels $[\varphi _\alpha ^ -]$ to all the one-particle states belonging to the configurations $[\Phi _n]$ that are present in the ground state $[\Psi _g^N]$ . Here we have introduced the overlap integrals $[S_\alpha = \langle \tilde \Psi _\alpha ^{N - 1} | \Psi _g^{N - 1} \rangle]$ of the `passive' electron and indicated by $[\varphi _\alpha ^ -]$ the time reversal of $[\varphi _\alpha]$ (in practice, the complex conjugate, if spin is neglected). Spin degeneracy has been taken into account by including an extra factor of 2.

The set of equations in (8) contains the complete description of all the outcomes of the photoemission process, be they of intrinsic (i.e. consequent to the relaxation of the system around the core hole) or extrinsic (excitations created by the photoelectron on its way out of the system) origin. The complete solution of (8) is out of the question; however, one can analyse the consequences of some solutions in particular cases. Since we are mainly interested in structural analysis, both in photoemission and in photoabsorption, we need only consider the completely relaxed or elastic channel (i.e. that for which $[\Delta E_\alpha = 0]$ ), since this is the only channel capable of carrying structural information. For simplicity we shall attribute to this channel here and in the following the index $[\alpha = 0]$ . This channel carries most of the weight since, as a typical value, $[| {S_0 } |^2 = | {\langle \tilde \Psi _0^{N - 1} | \Psi _g^{N - 1} \rangle } |^2 \sim 0.8]$ –0.9 (Mustre de Leon et al., 1991 ).

We can then think of solving the set of coupled Schrödinger equations (8) by eliminating all unwanted channels in favour of the elastic one. The result is a single equation for the channel function $[\varphi _0 ({\bf r})]$ with an effective complex energy-dependent non-local optical potential of the kind

$[[\,\nabla ^2 + k_0^2 - V_{\rm c} ({\bf r})]\,\varphi _0^{} ({\bf r}) = \textstyle\int {\Sigma ^{\rm opt} } ({\bf r},{\bf r}';\hbar \omega)\,\varphi _0 ({\bf r}')\,{\rm d}^3 r', \eqno(12)]$

where we have isolated the local Coulomb part of the potential (V_c) and indicated the energy dependence coming from the eliminated channels by the argument $[\hbar \omega]$ in $[\Sigma ^{\rm opt}]$ . Once this equation is solved, we can write each $[\varphi _\alpha ({\bf r})]$ in terms of $[\varphi _0 ({\bf r})]$ through a relation of the type $[\varphi _\alpha ({\bf r}) = \textstyle\int {A_\alpha } ({\bf r},{\bf r}';\hbar \omega)\,\varphi _0 ({\bf r}')\,{\rm d}^3 r']$ , which involves complicated inversions of the operators $[[{\nabla ^2 + k_\alpha ^2 - V_\alpha ({\bf r},{\bf r}')} ]]$ in (8). Again, a formal scheme of solution is given by Natoli et al. (1990). We can therefore write (11) as

$[\eqalignno{ {\rm d}\sigma (\omega) / {\rm d}\hat k_0 = {}&8\pi ^2 \alpha \hbar \omega \textstyle\sum\limits_{m_0 } \Big| \sum\limits_\alpha \Big\langle S_\alpha ^* \int A_\alpha ^ - ({\bf r},{\bf r}'\semi\hbar \omega) \,\varphi _0^ - ({\bf r}'\semi{\bf k}_0)\,{\rm d}^3 r' \cr& \times\big| \boldvarepsilon \cdot {\bf r} \big| \, \varphi _{l_0 m_0 }^{\rm c} ({\bf r}) \Big\rangle \Big| ^2, & (13) }]$

so that everything is expressed in terms of $[\varphi _0^ - ({\bf r})]$ . Notice that in the summation over α the most important term is $[S_0^* \varphi _0^ - ({\bf r})]$ , since by definition $[A_0^ - ({\bf r},{\bf r}';\hbar \omega) = \delta (r - r')]$ .

Our task is then to solve (12) with the asymptotic boundary conditions given by (10). This can be achieved in the framework of MST by transforming the integro-differential equation (12) into a Lippman–Schwinger equation with non-local potential, following the method illustrated by Natoli et al. (1986 , 1990). In these papers the derivation of the MS equations is given for local potentials, but it is clear that the same derivation is valid for non-local potentials as well. The essence of the method rests on the partition of the space into Voronoi polyhedra (equivalent to Wigner–Seitz cells for periodic systems) such that the diameters of the polyhedra are always smaller than the nearest distance between their centres. In each such polyhedra, a local solution $[\Phi _L ({\bf r})]$ of (12) is obtained, which behaves like $[J_L ({\bf r}) = j_l (k_0 r)\,Y_L (\hat r)]$ near the origin, where j_l (kr) is the usual spherical Bessel function of order l and $[Y_L ({\hat r})]$ are the spherical harmonics of type $[L \equiv (l,m)]$ (we shall use a real basis throughout). These functions are used to expand locally in the ith cell the overall solution of class C¹ (continuous together with first derivatives) in the whole space, satisfying the boundary condition (10). Thus we can write locally,

$[\Phi ^i ({\bf r}_i;{\bf k}_0) = \textstyle\sum\limits_L {A_L^i } ({\bf k}_0)\Phi _L^i ({\bf r}_i), \eqno(14)]$

provided the amplitudes $[A_L^i ({\bf k}_0)]$ satisfy the compatibility equations

$[\textstyle\sum\limits_{L'} {} C_{LL'}^i A_{L'}^i ({\bf k}_0) = A_L^o ({\bf k}_0) - \sum\limits_{jL'L''} {(1 - \delta _{ij}) } G_{LL'}^{ i j} S_{L'L''}^j A_{L''}^j ({\bf k}_0), \eqno(15)]$

where $[A_L^o ({\bf k}_0) = - 4\pi i^l \,Y_L ({\hat k}_0)\exp(i\,{\bf k}_0 \cdot {\bf R}_{io} ) (k_0/ \pi)^{1/2}]$ is the exciting amplitude originating from the plane wave in (10). A derivation is provided in Appendix A. Here, $[C_{LL'}^i]$ and $[S_{LL'}^i]$ are surface integrals over the surface S_i of the ith cell,

$[C_{LL'}^i = \textstyle\int\limits_{S_i } {\left[{\tilde H_L^ + ({\bf r}){\bf \nabla} \Phi _{L'}^i ({\bf r}) - \Phi _{L'}^i ({\bf r}){\bf \nabla} \tilde H_L^ + ({\bf r})} \right]} \cdot {\bf n}_i\,{\rm d}\sigma_i, \eqno(16)]$

$[S_{LL'}^i = \textstyle\int\limits_{S_i } {\left[{J_L ({\bf r}){\bf \nabla} \Phi _{L'}^i ({\bf r}) - \Phi _{L'}^i ({\bf r}){\bf \nabla} J_L ({\bf r})} \right]} \cdot {\bf n}_i\,{\rm d}\sigma_i, \eqno(17)]$

and $[G_{LL'}^{ i j}]$ are the KKR structure factors, which have the well known form

$[G_{LL'}^{ i j} = 4\pi \textstyle\sum\limits_{L''} i ^{l - l' + l''} C_{LL'}^{L''}\, \tilde H_{L''}^ + ({\bf R}_{ij}),\eqno(18)]$

where $[\tilde H_L^ + ({\bf r}) = - ik\,h_l^ + (k_0 r)\,Y_L (\hat r)]$ , h_l⁺ is the Hankel function with outgoing spherical-wave behaviour, $[{\bf R}_{ij} = {\bf R}_i - {\bf R}_j]$ is the vector connecting the origins of the two cells centred at $[{\bf R}_i]$ and $[{\bf R}_j]$ , and $[{\bf R}_{io} = {\bf R}_i - {\bf R}_o]$ connects site i with the origin of the coordinates o, which is assumed to coincide with the photoabsorber. The quantities $[C_{LL'}^{L''} = \textstyle\int {Y_L } (\hat r)\,Y_{L'} (\hat r)\,Y_{L''} (\hat r)\,{\rm d}\hat r]$ are known as Gaunt coefficients.

In order to be able to use the physical language of MST, we introduce the quantities $[B_L^i ({\bf k}_0) = \textstyle\sum_{L'} {S_{LL'}^i } A_{L'}^i ({\bf k}_0)]$ , and we define the suitably normalized local basis solutions $[\tilde \Phi _L^i = \textstyle\sum_{L'} {\Phi _{L'}^i (S^i })_{L'L}^{ - 1}]$ , so that $[\Phi ^i ({\bf r}_i;{\bf k}_0) = \textstyle\sum_L {\tilde \Phi _L^i } ({\bf r}_i)\,B_L^i ({\bf k}_0)]$ . The coefficients $[B_L^i ({\bf k}_0)]$ are easily seen to satisfy the MS equations

$[\textstyle\sum\limits_{L'} {(T^i\,)_{LL'}^{ - 1} } \,B_{L'}^i ({\bf k}_0) = A_L^o ({\bf k}_0) - \textstyle\sum\limits_{jL'} {}(1 - \delta _{ij})\, G_{LL'}^{i j}\, B_{L'}^j ({\bf k}_0), \eqno(19)]$

where we have introduced the quantity $[T_{LL'}^i = \textstyle\sum_{L''} {S_{LL''}^i (C^i)_{L''L'}^{ - 1} }]$ . In the case of cells of spherical shape it can be shown (Natoli et al., 1990) that $[T_{LL'}^i]$ is the scattering matrix of the non-spherical potential located inside the ith cell and $[B_L^i ({\bf k}_0)]$ is a total scattering amplitude impinging on that cell in response to a plane-wave excitation of momentum $[{\bf k}_0]$ . Indeed, (19) are the self-consistent equations for these amplitudes. We shall assume in the following that we are dealing with spherical cells, in order to make the discussion more physical. In the general case, (15) should be used.

For convenience, we define the cell T-matrix $[T_{\rm c} \equiv T_{LL'}^i \delta _{ij}]$ and the matrix $[G \equiv (1 - \delta _{ij}) G_{LL'}^{ij}]$ of the free-spherical-wave propagator. Then (19) are easily solved for the amplitudes $[B_L^i ({\bf k}_0)]$ to give

$[B_L^i ({\bf k}_0) = - 4\pi ({\bf k}_0 / \pi)^{1/2} \textstyle\sum\limits_{jL'} \tau _{LL'}^{ i j}\, i^{l'} \,Y_{L'} ({\hat k}_0)\exp(i\,{\bf k}_0 \cdot R_{jo}), \eqno(20)]$

where we have introduced the inverse τ of the MS matrix; $[\tau _{LL'}^{i j} = [(T_{\rm c}^{ - 1} + G)^{ - 1}] _{LL'}^{ij}]$ , which is known as the full-scattering path operator and gives the total amplitude of propagation from site i to site j, starting with angular momentum L and arriving with angular momentum [L'] .

We now have all the ingredients that are needed to calculate the photoemission cross section. Introducing the energy-dependent matrix element for the creation of the photoelectron,

$[M_{L_0L}^{\rm c}(\hbar\omega)= \textstyle\sum\limits_\alpha\left\langle S_\alpha^* \int A_\alpha^- ({\bf r},{\bf r}'\semi\hbar\omega)\,\tilde\Phi_L^- ({\bf r}')\,{\rm d}^3 r'|\boldvarepsilon \cdot {\bf r}| \varphi _{L_0 }^{\rm c} ({\bf r}) \right\rangle, \eqno(21)]$

we easily find that our quantity of interest takes the simple form

$[{\rm d}\sigma (\omega)/ {\rm d}\hat k_0 = 8\pi ^2 \alpha \hbar \omega \textstyle\sum\limits_{m_0 } {\left| {\textstyle\sum\limits_L {M_{L_0 L}^{\rm c} } (\hbar \omega)B_L^o ({\bf k}_0)} \right|} ^2. \eqno(22)]$

Remembering (20), we then see that the photoemission current along direction $[\hat k_0]$ is the modulus squared of the sum of all of the possible composite amplitudes that are obtained as products of an amplitude M_{L₀ L}^c for exciting a core electron at site o (the origin), times the amplitude of propagation $[\tau _{LL'}^{o j}]$ from this site to any other site (cell) j with final angular momentum [L'] , times the amplitude $[Y_{L'} (\hat k_0)]$ for emission along $[\hat k_0]$ , times the phase factor $[\exp(i\,{\bf k} \cdot R_{jo})]$ that takes into account the phase relation of the electronic wave between sites o and j. The photoemission current is therefore the result of a complicated interference process that requires the reliability of the optical potential.

Expression (22) has constituted the basis for the interpretation of photoelectron diffraction data at low photoelectron kinetic energies with some success (Gunnella et al., 1998 ). Computer programs in the MT approximation are also available (Gunnella et al., 2000 ). The emitted current depends on three variables: the photon energy and the two polar angles of the direction $[\hat k_0]$ . Therefore, as anticipated at the beginning of this section, by fixing the energy we can test the performance of a particular model optical potential in reproducing polar and azimuthal diffraction scans.

2.2. The photoabsorption cross section

Turning now to absorption, we need to integrate (11) over $[\hat k]$ and sum over all the final channels that have the same total energy $[E_{\rm tot} = \hbar \omega + E_g]$ . To this purpose it is expedient to make explicit the energy-conservation condition and to label the final state with f. By defining the Green's function matrix with outgoing wave boundary conditions,

$[G_{\alpha \alpha '}^ + ({\bf r},{\bf r}';E) = \textstyle\sum\limits_f \,\{\varphi \,_\alpha ^f ({\bf r})[\varphi \,_{\alpha '}^f ({\bf r}')]^* \}/(E - \Delta E_f - k_f^2 + i\eta ),]$

the photoabsorption cross section is easily seen to be

$[\eqalignno{ \sigma _{\rm abs} (\omega) = & - 8\pi \alpha \hbar \omega \, {\rm Im} {\textstyle\sum\limits_{m_0 }\Big[{ \int\!\!\int {\rm d}^3 r \,\varphi \,_{L_0 }^{\rm c} ({\bf r})\,{\boldvarepsilon}\cdot {\bf r} }} \cr & \,\times{\textstyle\sum\limits_{\alpha \alpha '} S_\alpha ^* S_{\alpha '} G_{\alpha \alpha '}^ + ({\bf r},{\bf r}'\semi\hbar \omega - I_{\rm c} ){\boldvarepsilon}} \cdot {\bf r}' \, \varphi \, _{L_0 }^{\rm c} ({\bf r}')\,{\rm d}^3 r' \Big], & (23)}]$

which is the same as the expression that we would have obtained if we had started from (1), (3) and (4). In Appendix C-4 of Natoli et al. (1990) it is shown that this matrix satisfies the set of coupled equations,

$[\eqalign{(\nabla ^2 + k_\alpha ^2) \,G_{\alpha \beta }^ + ({\bf r},{\bf r}'\semi E) - \textstyle\sum\limits_\gamma \int V_{\alpha \gamma } ({\bf r},{\bf r}'')\,G_{\gamma \beta }^ + ({\bf r}'',&{\bf r}'\semi E)\,{\rm d}^3 r'' \cr & = \delta _{\alpha \beta } \,\delta ({\bf r} - {\bf r}'), }]$

where we have written E for $[\hbar \omega - I_{\rm c}]$ . Following the same steps as in the photoemission case, the elimination of all channels in favour of the relaxed one ( $[\alpha = 0]$ ) leads to the following expression

$[\eqalignno{ \textstyle\sum\limits _{\alpha \alpha '} S_\alpha ^* \, &S_{\alpha '} \,G_{\alpha \alpha '}^ + ({\bf r},{\bf r}';E) =\cr & \textstyle\sum\limits_{\alpha \alpha '} {S_\alpha ^* \, S_{\alpha '} } \int\!\!\int A_\alpha ^* ({\bf r},{\bf x};E)\,G_{00}^ + ({\bf x},{\bf x}';E)\,A_{\alpha '} ({\bf x}',{\bf r}';E)\,{\rm d}^3 x\,{\rm d}^3 x', & (24)}]$

where the non-local operator $[A_\alpha ({\bf r},{\bf r}';E)]$ is the same as before and G₀₀ obeys an equation corresponding to (12) with the same optical potential:

$[\eqalignno{[\nabla ^2 + k_0^2 - V_c ({\bf r})]\,G_{00}^ + ({\bf r},{\bf r}';E) - {\textstyle\int} {\Sigma ^{\rm opt} } ({\bf r},{\bf x};E) \, G_{00}^ + &({\bf x},{\bf r}';E)\,{\rm d}^3 x \cr & = \delta ({\bf r} - {\bf r}'). & (25)}]$

The solution to this equation within MST, for $[{\bf r}]$ in cell i and $[{\bf r}']$ in cell j, can be written as (Faulkner & Stocks, 1980 )

$[\eqalignno{G_{00}^ + ({\bf r}, {\bf r}';E) = & \textstyle\sum\limits_{LL'} \tilde \Phi _L ({\bf r})\left({\tau _{LL'}^{i j} - \delta _{ij} T_{LL'}^i } \right)\tilde \Phi _{L'} ({\bf r}') \cr & + \delta _{ij} \textstyle\sum\limits_{LL'} \tilde \Phi _L ({\bf r}_ \lt)T_{LL'}^i \tilde \Psi _{L'}^ + ({\bf r}_ \gt), & (26)}]$

remembering that the functions $[\tilde \Phi _L ({\bf r})]$ do not carry the normalization to one state per Rydberg. The second term on the right-hand side is the singular part of Green's function, $[{\bf r}_ \lt \,({\bf r}_ \gt)]$ is the lesser (greater) of [r,r'] and $[\tilde \Psi _L^ + ({\bf r})]$ is the solution of (25) inside cell i that is irregular at the origin and matches smoothly to $[\tilde H_L^ + ({\bf r})]$ at the boundary. Making the reasonable assumption that the range of the functions $[A_\alpha ({\bf r},{\bf r}';E)]$ is of the order of the atomic dimensions (this assumption is obviously true for A₀) and using (24) and (26) in (23), we finally obtain for the photoabsorption cross section

$[\eqalignno{\sigma _{\rm abs} (\omega) = & - 8\pi \alpha \, \hbar \omega \cr & \times {\rm Im} \textstyle\sum\limits_{m_0 LL'} {M_{L_0 L}^{\rm c*} (\omega)\, ({\tau _{LL'}^{oo} - T_{LL'}^o } )} \,M_{L_0 L'}^{\rm c} (\omega) + \sigma _{\rm at} (\omega).& (27)}]$

This expression arises because of the localization of the core initial state. M_{L₀ L}^c is given by equation (21) and we define an `atomic' absorption given by

$[\eqalignno{\sigma _{\rm at} (\omega) = & - 8\pi \alpha \, \hbar \omega \cr & \times {\rm Im} \textstyle\sum\limits_{m_0 } {\int\!\!\int {{\rm d}^3 r} } \,\varphi _{L_0 }^{\rm c} ({\bf r})\, {\boldvarepsilon} \cdot {\bf r}\,M({\bf r},{\bf r}';\omega){\boldvarepsilon} \cdot {\bf r}'\varphi _{L_0 }^{\rm c} ({\bf r}')\, {\rm d}^3 r',& (28)}]$

where

$[\eqalignno{M ({\bf r},{\bf r}'\semi\omega) = {}&\textstyle\sum\limits_{\alpha \alpha '} S_\alpha ^* S_{\alpha '} \textstyle\sum\limits_{LL'} \textstyle\int\!\!\int A_\alpha ^* ({\bf r},{\bf x}\semi\omega)\tilde \Phi _L ({\bf x}_ \lt)\cr &\times T_{LL'}^o \tilde \Psi _{L'} ({\bf x}_ \gt)A_{\alpha '} ({\bf x}',{\bf r}'\semi\omega) \,{\rm d}^3 x\,{\rm d}^3 x'. & (29)}]$

Equation (27) is the expression that we wanted to arrive at. It is valid under quite general conditions, as is apparent from our derivation, and shows the natural partition into an `atomic' contribution and an MS contribution. By removing all cells of the cluster except that containing the photoabsorber, $[\tau _{LL'}^{oo}]$ reduces to $[T_{LL'}^o]$ so that the MS contribution is zero. Even though in a multiatom system we cannot define precisely an atomic entity, we shall continue to use the term `atom' for the central cell, because it is after all a reasonable approximation to the isolated atom. The fact that the second site index is equal to the first means that the electron source and the detector coincide, as already anticipated. Since the escape direction of the photoelectron has been integrated out, the only variable left is the energy and we can only test energy-dependent diffraction patterns. Therefore, these patterns might be distorted (compared with photoelectron diffraction patterns) by the energy dependence of the matrix elements that are necessary to create and detect the photoelectron at the absorption site.

However, we can hope to recover the structural information that we are interested in if we can neglect electronic correlations in the final state and if only one configuration (Slater determinant) is predominant in (3). Both of these conditions seem to be true, to a good approximation, for the K-edge spectra of metal ions in materials of biological interest. Indeed, the final photoelectronic states reached in these instances possess p symmetry around the photoabsorber and therefore are sufficiently delocalized not to suffer correlation effects with the electrons of the system. Moreover, very often the spin-restricted (unrestricted) Hartree–Fock approximation for the initial state provides a reasonably good description of the ground state for closed (open) shell configurations. Only in special cases do we have examples of configuration mixing in the ground state, as in the case of the Cu²⁺ ion in La₂CuO₄ (Wu et al., 1996 ) or valence fluctuating compounds. In such instances, a better strategy would be to solve the multichannel equations for the configurations that enter in the ground state and to eliminate all the other channels by expressing them as a function of the configurations of interest.

The predominance of one particular configuration implies that the spectral shape of the excited channels is rather featureless and similar to the ground state, since most of them consist of particle–hole or plasmon-like excitations, which do not drastically change the ground-state potential. Therefore, to a reasonable approximation, we can put $[A_\alpha ({\bf r} - {\bf x};\omega) = A_\alpha (\omega)\,\delta ({\bf r} - {\bf x})]$ . We can also write $[| {\textstyle\sum_\alpha {S_\alpha \,A_\alpha (\omega)} } |^2 =]$ $[| {S_0 (\omega)} |^2]$ , because $[| {S_0 } |^2]$ is the preponderant term of the series [by definition $[A_0 (\omega) = 1]$ ]. Thus, we obtain from (24)

$[\textstyle\sum\limits_{\alpha \alpha '} {S_\alpha ^* S_{\alpha '} \,G_{\alpha \alpha '}^ + ({\bf r},{\bf r}';E) = \left| {S_0 (\omega)} \right|} ^2 G_{00}^ + ({\bf r},{\bf r}';E).\eqno(30)]$

This equation, together with (23), tells us that the effect of the eliminated channels results in a shape function $[| S_0 (\omega) |^2]$ that modulates the absorption coefficient originating from the primary channel

$[\sigma _{\rm abs} (\omega) = \left| {S_0 (\omega)} \right|^2 \sigma _{\rm abs}^0 (\omega).\eqno(31)]$

The expression for $[\sigma _{\rm abs}^0 (\omega)]$ is easily obtained by observing that, in the same approximation, (29) gives

$[M({\bf r},{\bf r}';\omega) = \left| {S_0 (\omega)} \right|^2 \textstyle\sum\limits_{LL'} {\tilde \Phi _L ({\bf r}_ \lt)T_{LL'}^o \tilde \Psi _{L'} ({\bf r}_ \gt)} = \left| {S_0 (\omega)} \right|^2 \tilde M({\bf r},{\bf r}';\omega),\eqno(32)]$

whereas (21) becomes, dropping $[\hbar]$ ,

$[M_{L_0 L}^{\rm c} (\omega) = \langle {\tilde \Phi _L ({\bf r})} \right|{\boldvarepsilon} \cdot {\bf r}| {\varphi _{L_0 }^{\rm c} ({\bf r})} \rangle \textstyle\sum\limits_\alpha {S_\alpha } A_\alpha (\omega) = \tilde M_{L_0 L}^{\rm c} (\omega)\textstyle\sum\limits_\alpha {S_\alpha } A_\alpha (\omega),\eqno(33)]$

so that insertion of (33) into (27) and (28) provides the desired result, simply by replacing the transition matrix elements M with $[\tilde M]$ ,

$[\eqalignno{\sigma _{\rm abs}^0 (\omega) & = - 8\pi \alpha \hbar \omega \,{\rm Im} \textstyle\sum\limits_{m_0 LL'} {\tilde M_{L_0 L}^{\rm c*} (\omega)\left({\tau _{LL'}^{oo} - T_{LL'}^o } \right)} \tilde M_{L_0 L'}^{\rm c} (\omega) + \sigma _{\rm at}^0 (\omega) \cr & = \sigma _{\rm str}^0 (\omega) + \sigma _{\rm at}^0 (\omega). &(34)}]$

Even with these approximations, the first-principle calculation of $[| S_0 (\omega) |^2]$ is not an easy task. The most noticeable effect is provided by double electron transitions (that are included in our scheme), since they show up as a change of slope in the `atomic' background component of the cross section (Filipponi, 1995b ), which occurs at definite energies that depend on the photoabsorber. An analytical or numerical modelling of these effects would be highly desirable, although for our purposes it is sufficient to know that such `kinks' can be factorized into this shape function.

We are left with the problem of determining the nature of the optical potential $[\Sigma ^{\rm opt} ({\bf r},{\bf r}';E)]$ in (25). This potential contains the effect of both the intrinsic channels (excitations induced in the system by the sudden creation of the core hole) and the extrinsic channels (excitations created by the photoelectron on its way out of the system). Model systems to describe both types of processes, and their interference patterns have been studied by some authors (Fujikawa, 1999 ; Campbell et al., 2002 , and references therein). However, a practical scheme for realistic calculations has not yet been devised. The only optical potentials currently in use are those based on the Dyson self-energy of the photoelectron propagation in the system. This type of self-energy clearly accounts only for the extrinsic losses. Depending on the systems, various reasonable approximations have been devised.

It is well known, for example, that for metals we can obtain very good agreement with the observed absorption spectra by using a one-particle approach with an X-α potential and convoluting the calculated spectrum with a Lorentzian broadening function, which has an energy-dependent width that is related to the mean-free path of the photoelectron in the system by the relation (Müller & Wilkins, 1982 ; Müller et al., 1984 )

$[\Gamma (E) = [\hbar / \lambda (E)] (2E / m)^{1/2},]$

where $[\Gamma (E)]$ is the full width at half-maximum.

In the framework of the above multichannel approach, this finding can be rationalized by observing that in a metal the completely relaxed channel together with the plasmon excitation channels (whether intrinsic or extrinsic) almost completely exhaust the sum rule

$[\textstyle\sum\limits_\alpha {| {S_\alpha } |} ^2 = \textstyle\sum\limits_\alpha \left\langle {\Psi _g^{N - 1} } | {\tilde \Psi _\alpha ^{N - 1}} \right\rangle \left\langle {\tilde \Psi _\alpha ^{N - 1}} | {\Psi _g^{N - 1}} \right\rangle = \left\langle {\Psi _g^{N - 1}} | {\Psi _g^{N - 1}} \right\rangle = 1,]$

which holds because of the completeness of the intermediate relaxed states $[\tilde \Psi_\alpha ^{N - 1}]$ . Indeed, the intensity of the double-electron excitation channels is of the order of 10⁻²–10⁻³ times that of the main relaxed channel (Filipponi, 1995b). Therefore, an optical potential given by $[V_{X - \alpha } + i\Gamma (E)]$ is able to give a satisfactory picture of the absorption process for metals. The only discrepancy with experiment is that the calculated absorption maxima fall short of the observed maxima because of the energy independence of the X-α exchange.

A better approximation is provided by the Hedin–Lundqvist (H–L) potential (Hedin & Lundqvist, 1969, 1971), owing to its energy-dependent exchange and its imaginary part that is able to reproduce rather accurately the observed mean free path in metals (Penn, 1987 ). The H–L potential is the self-energy (based on the GW approximation and calculated at the local density of the system under study) of an electron that is propagating in a homogeneous interacting electron gas. Although initially devised to describe exchange and correlation corrections to the Coulomb potential due to the valence charge, Lee & Beni (1977 ) have extended its validity in the atomic-core region as well.

By neglecting the effect of the intrinsic processes, we can approximate $[\Sigma ^{\rm opt} ({\bf r},{\bf r}';E)]$ as the Dyson self-energy of the photoelectron in the final state. This approximation is consistent with the physical picture of the photoabsorption process, in which we add an electron to the ground state of the (Z + 1)-equivalent atom. From this point of view, $[G_{00} ({\bf r},{\bf r}';E)]$ describes the propagation amplitude of the excited photoelectron from point $[{\bf r}]$ to point $[{\bf r}']$ . $[G_{00} ({\bf r},{\bf r}';E)]$ is the probability amplitude that the added electron remains in the state in which it was added to the system. The imaginary part gives the total probability amplitude for scattering out of this initial state. In this scheme the Dyson self-energy $[\Sigma_{\rm D} ({\bf r},{\bf r}';E)]$ acts as a complex optical potential that describes the reduction of the wave-function amplitude of the elastic channel due to the transitions to all the other channels. The localization of the initial core state has the consequence that the optical paths of the photoelectron in the final state begin and end at the photoabsorbing site. We expect that neglecting the effects of the intrinsic processes on $[\Sigma^{\rm opt} ({\bf r},{\bf r}';E)]$ is a reasonably good approximation, since their main effect is already incorporated in the shape function $[| {S_0 (\omega)} |^2]$ .

We can therefore interpret the H–L potential as an effective optical potential that controls the propagation and damping of the excited photoelectron everywhere in the system. In this approach, this potential can be viewed as a local-density approximation to the self-energy of the photoelectron in real systems. Nowadays it is the potential that is most widely used in the calculation of the absorption and photoelectron diffraction cross sections of many systems, ranging from metals and semiconductors to ionic and covalent systems with varied success. In Appendix B we derive in some detail analytical expressions for the H–L potential, both for the sake of completeness and in order to make contact with other kinds of potentials, like the X-α and the Dirac–Hara potentials.

Finally, Fujikawa et al. (2000 ) have improved on the H–L potential by restricting its validity to the valence charge, as originally devised. They relied on the GW approximation for the photoelectron self-energy, $[\Sigma _{\rm GW} = GW]$ in the solid, where $[W = \varepsilon ^{ - 1} V]$ , V is the bare Coulomb interaction and $[\varepsilon = 1 - VP]$ is the dielectric response function of the system, and they split both the polarization propagator P = P ^v_o + P ^c and the one-electron Green's function G = G ^v + G ^c into core and valence parts. Since the core polarization was assumed to be much smaller than the valence polarization, Fujikawa et al. (2000) obtain an expansion in powers of P ^c for $[\Sigma _{\rm GW} = G\,^{\rm v} W\,^{\rm v} + V\,_{\rm ex}^{\rm c} + G\,^{\rm v} W\,^{\rm v} P\,^{\rm c} W\,^{\rm v} +. \,.\,.]$ . Here G ^v W ^v is the self-energy for the valence electrons, which, when calculated via the plasmon-pole approximation for the dielectric function, is equivalent to the H–L potential; V _ex^c = G ^c V is the bare Hartree–Fock exchange potential and G ^v W ^v P ^c W ^v is the screened polarization potential for the ion cores. Details of this latter potential are given in Appendix C. Preliminary calculations for photoelectrons with kinetic energy greater than about 100 eV show that this non-local potential gives the same scattering amplitude as the total H–L potential at large scattering angles and provides a better description of small-angle scattering (Fujikawa et al., 2000). However, more work is needed to establish its performance at low photoelectron energies.

Once the optical potential has been specified, we can proceed to the calculation of the key ingredients in (34), namely the cell functions $[\tilde \Phi _L ({\bf r})]$ and $[\tilde \Psi _L^ + ({\bf r})]$ in the absorbing sphere, which are needed to calculate the atomic absorption $[\sigma _{\rm at}^0 (\omega)]$ and the transition-matrix elements $[\tilde M]$ that create the photoelectron, and the scattering-path operator $[\tau _{LL'}^{i j} = [(T_{\rm c}^{ - 1} + G)^{ - 1}] _{LL'}^{i j}]$ . The inversion of the MS matrix (T_c^{- 1} + G) becomes time consuming at energies greater than about 50–150 eV, depending on the number of cells in the cluster. Fortunately, in most cases, above ∼50 eV we can invert the MS matrix by series expansion, i.e.

$[\eqalignno{\tau & = (T_{\rm c}^{ - 1} + G)^{ - 1} = (\,I + T_{\rm c} G)^{ - 1} T_{\rm c} \cr & = \textstyle\sum\limits_n {(- 1)^n } (T_{\rm c} G)^n T_{\rm c} = \textstyle\sum\limits_n {(- 1)^n } T_{\rm c} (GT_{\rm c})^n, & (35)}]$

so that the total absorption is seen to be made up of an atomic smooth contribution plus an infinite series of oscillatory signals (Benfatto et al., 1986 ; Tyson et al., 1992 ). This observation has been the basis for the development of the packages mentioned in §1, based on the empirical extraction of a structural EXAFS signal from the experimental data and the comparison with a theoretical signal.

A major ingredient in the extraction procedure is the definition of a background atomic absorption coefficient $[\mu _{\rm at} (\omega)]$ such that the structural signal is defined as

$[\chi (\omega) = [\mu (\omega) - \mu _{\rm at} (\omega)]/[\mu _{\rm at} (\omega)].]$

By identifying $[\mu _{\rm at} (\omega)]$ with $[\sigma _{\rm at} (\omega)]$ (ignoring the usual proportionality factor $[N\rho/A]$ , where N is the Avogadro number) we are led by (34) to the identification

$[\chi (\omega) = \sigma _{\rm str}^0 (\omega) / \sigma _{\rm at}^0 (\omega). \eqno(36)]$

In this formula, the shape function $[| {S_0 (\omega)} |^2]$ seems to have dropped out from the ratio. However, the function's disappearance is a consequence of the identification of $[\mu _{\rm at} (\omega)]$ with $[\sigma _{\rm at} (\omega)]$ , which is only approximate because the empirical definition of $[\mu _{\rm at} (\omega)]$ does not take into account the exact form of $[| {S_0 (\omega)} |^2]$ . The nearer $[\mu _{\rm at} (\omega)]$ is to $[\sigma _{\rm at} (\omega)]$ , the less dependent $[\chi (\omega)]$ becomes on the shape function. This is the reason why in many instances of EXAFS analysis $[| {S_0 (\omega)}|^2]$ seems to be constant and very near to 1. Also note that the atomic absorption $[\sigma _{\rm at}^0 (\omega)]$ does not factorize from the structural signal $[\sigma _{\rm str}^0 (\omega)]$ . Indeed, only if the optical potential is real can we show that

$[\sigma _{\rm at}^0 (\omega) = - 8\pi \alpha \hbar \omega \,{\rm Im} \textstyle\sum\limits_{m_0 LL'} {\tilde M_{L_0 L}^{\rm c*} (\omega)\,({T_{LL'}^o } )} \,\tilde M_{L_0 L'}^{\rm c} (\omega),]$

so that the structural signal is proportional to the atomic absorption. However, because of the presence of inelastic processes, the potential is complex, and the different behaviour of $[\sigma _{\rm at}^0 (\omega)]$ and $[\sigma _{\rm str}^0 (\omega)]$ should be taken into account in a refined treatment of the experimental data.

Therefore, in order to extract structural and electronic information on the system under study, we should not rely on an empirical definition of $[\mu _{\rm at} (\omega)]$ , which becomes increasingly difficult to define in the XANES part of the absorption spectrum. Rather, we should try to fit the entire experimental spectrum, especially since the relevant information is confined to the first 100–150 eV from the edge and is consequently well within the energy range where a complete inversion of the MS matrix is possible for a wide variety of systems, including biological material. Indeed from photoelectron diffraction studies we already know that the coherent diffraction process carrying the structural information is well described by an optical potential of the H–L type (Chen et al., 1998 ; Gunnella et al., 1998). Regarding the matrix elements for the creation of the photoelectron and its energy dependence, past and present experiences indicate that the smoothing action of the complex part of the potential helps us in reproducing the correct behaviour. Remaining discrepancies can be easily absorbed into a redefinition of this complex part. Even the effect of the shape function $[| {S_0 (\omega)} |^2]$ on the absorption spectrum, under the assumptions described above, can be mimicked via an additional damping. These statements will be substantiated by the applications of the method to particular systems in a following paper (Benfatto et al., 2003 ).

2.3. The mean free path

Expression (19) tells us that the structural part of the photoabsorption cross section $[\sigma _{\rm str} (\omega)]$ is proportional to the imaginary part of the product of the amplitude M for emitting the photoelectron, times the full scattering amplitude of propagation $[\tau]$ from the photoemitter site and back, excluding the contribution of the photoemitter itself, times another amplitude M for detecting the photoelectron. Since the optical potential is complex, we expect that the coherent part of the propagation (the part that only carries the structural information) will be attenuated by a damping factor related to the mean free path of the photoelectron in the system, whereas incoherent terms will appear that describe the effects of scattering out of the coherent channel due to the presence of inelastic processes. These incoherent terms are described by the complex part of the potential.

In this section we shall show explicitly that the mathematical structure of the theory contains the intuitive picture of the photoelectron damping in the expected way, by making contact with existing expressions for the mean free path $[\lambda]$ . At the same time, we shall find out how the incoherent terms due to the presence of inelastic processes find their place in the theory and, as a by-product, how we can mimic their effects when we ignore the details of their manifestation.

For the sake of illustration, we shall assume that the MS series for $[\tau]$ in (35) converges and that the optical potential is local and of the muffin-tin form. The general nth term of the series, dropping one factor T_c^o that factorizes into the amplitude for emitting and detecting the photoelectron, has the form

$[\left [{(T_{\rm c} G)^n } \right]_{LL'}^{ oo} = \textstyle\sum\limits_{ij\cdots k}\,\,\textstyle\sum\limits_{L_1 L_2 \cdots L_n } {t\,_l^o }\, G_{LL_1 }^{ oi} t\,_{l_1 }^i \, G_{L_1 L_2 }^{ij} \,t\,_{l_2 }^j \cdots G_{L_n L'}^{ko} \eqno(37),]$

where $[t\,_l^i = k^{-1}\exp(i\delta _l^i ) \sin \delta _l^i]$ is the atomic scattering matrix for the atom at site i and angular momentum l in terms of the corresponding phase shift $[\delta _l^i]$ ; $[G_{LL'}^{ij}]$ are the spherical wave propagators that are introduced in Appendix A. Without loss of generality, for simplicity we shall consider in (37) the case n = 3 with three sites o, i and j, so that j = k, L_n = L₂.

Since the potential is complex, the generic atomic phase shift $[\delta]$ is also complex, so we can write $[\delta = \delta _1 + i\delta _2]$ . Therefore,

$[\eqalignno{kt & = \exp(i\delta ) \sin \delta \cr&= \exp(-2\delta _2) \exp(i\delta _1 ) \sin \delta _1 + i\,[1 - \exp( - 2\delta _2 )]/2, &(38)}]$

which reduces to the known expression with real $[\delta]$ in the limit $[\delta _2 \rightarrow 0]$ . In an electron–atom scattering process with complex potential (Landau & Lifshitz, 1966 ), it is known that, within a factor k^-2, Im kt represents the total scattering cross section (elastic plus inelastic), whereas $[| {kt}|^2]$ gives the elastic cross section (without energy loss for the impinging electron). Based on the definition (25), the relation

$[{\rm Im} \,kt = \left| {kt} \right|^2 + \,[1 - \exp(- 4\delta _2)]/4]$

holds, so that $[[1 - \exp(- 4\delta _2)]/4]$ is the inelastic cross section. For convenience, we shall work henceforth in terms of the dimensionless quantity kt, since the propagator G in (37) is proportional to k, the photoelectron momentum.

From (38), it is clear that the coherent signal is obtained by choosing, for all the kt factors appearing in (37), the term $[\exp(- 2\delta _2) \exp(i\delta _1) \sin \delta _1]$ . All the other terms containing at least one factor $[i\,[1 - \exp(- 2\delta _2)]/2]$ describe inelastic processes with loss of coherence for the photoelectronic wave. Indeed, we have $[[1 - \exp(- 2\delta _2)]/2 \simeq [1 - \exp(- 4\delta _2)]/4 \simeq \delta _2]$ for sufficiently small $[\delta _2]$ , which is usually the case.

Concerning the coherent signal of interest, we see that it has the same form and value as the structural signal obtained for a real potential, except for a damping factor coming from each kt and each propagator G. To proceed further, we need to find an expression for this latter factor. Again, for illustrative purposes, we shall work in the so-called plane-wave approximation (PWA), but it can be shown that all the following considerations are valid even when spherical-wave corrections for the photoelectron propagating wave are taken into account. In the PWA we can write (Natoli et al., 1986)

$[G_{LL'}^{ij} \simeq \exp(i\,\kappa ^{\rm I} R_{ij}) \, Y_L (\hat R_{ij})\,Y_{L'} (\hat R_{ij})/ (\kappa ^{\rm I} R_{ij}),]$

where R_ij is the distance between sites i and j and $[\kappa ^{\rm I} = [{E - V_{\rm I} } ]^{1/2} = \kappa _1^{\rm I} + i\kappa _2^{\rm I}]$ is the photoelectron momentum in the interstitial region, in which the potential V_I is constant and can be complex.

Calling R_i the sphere radius of the potential at site i, it is easily seen that in the case of the three sites o, i and j the damping is represented by an exponential factor with exponent

$[- 2\delta _2^o - \kappa _2^{\rm I} R_{oi} - 2\delta _2^i - \kappa _2^{\rm I} R_{ij} - 2\delta _2^j - \kappa _2^{\rm I} R_{ijo}. \eqno(39)]$

We now need an expression for the phase shift $[\delta ^i]$ , which in the WKB approximation is given by (Hara, 1967 )

$[\delta ^i = \textstyle\int\limits_0^{R_i } [E - V_i (r)]^{1/2} {\rm d}r - \kappa ^{\rm I} R_i, \eqno(40)]$

where the potential V_i (r) inside the sphere R_i is assumed to be complex. Notice that in keeping with the plane-wave approximation for the propagator G we have dropped the centrifugal term in the WKB expression, so that in this case the phase shifts are independent of the angular momentum l. However, as already emphasized, the argument remains valid even when spherical-wave corrections are taken into account.

In the path going from site o, the photoabsorbing site, to site i and j and back to o, each atom is traversed for a length 2R, if R is the radius of the corresponding sphere. Writing R_path = R_oi + R_ij + R_jo for the total length of the path, we see that the damping factor is given by $[\exp (- \kappa _2 R_{\rm tot})]$ , where

$[\kappa _2 = R_{\rm path}^{-1}\,{\rm Im} \textstyle\int\limits_{\rm path} [E - V(r)]^{1/2} {\rm d}r, \eqno(41)]$

since in the interstitial region $[V(r) \equiv V_{\rm I}]$ and, for example, inside sphere j, $[V(r) \equiv V_j (r)]$ . Equation (41) is exactly the same as the expression that we would have obtained had we studied the propagation of an electronic wave in the potential V(r) by solving the problem in the WKB approximation. The mean free path is accordingly given by

$[\lambda = (2\kappa _2)^{-1}, \eqno(42)]$

which is consistent with the fact that $[\exp (- \kappa _2 R_{\rm tot}\,)]$ is an attenuation factor for an amplitude of propagation. A further simplification is achieved if we take into account that everywhere in the system $[E - V_1 (r) \gg V_2 (r)]$ , V₁ (r) and V₂ (r) being the real and imaginary parts of the potential. Then, expanding the square root in (41) as

$[[E - V(r)]^{1/2} \simeq [E - V_1 (r)]^{1/2} + (i /2) V_2 (r)\,[E - V_1 (r)]^{-1/2},]$

we obtain

$[\kappa _2 = 2R_{\rm path}^{-1} \textstyle\int\limits_{\rm path} V_2 (r) \,[E - V_1 (r)]^{-1/2} {\rm d} r \le (2 k R_{\rm path})^{-1} \int\limits_{\rm path} V_2 (r) \,{\rm d}r,]$

remembering that k = E^1/2. Therefore, from (42), we finally get an expression in atomic units,

$[\lambda \, ({\rm a.u.}) = k\,({\rm a.u.}^{ - 1}) / [\bar \Sigma _2 \,({\rm Ryd})], \eqno(43)]$

where $[\bar \Sigma _2 = R_{\rm path}^{-1} \textstyle\int_{\rm path} V_2 (r) \,{\rm d}r]$ , since V₂ (r) is actually a self-energy. Equation (43) is the same as that given by Penn (1987), who has, in conventional units,

$[\lambda = (\hbar ^2 k^2)/ (2m k {\bar \Sigma _2 }) = {E_k } / (k {\bar \Sigma _2 }),]$

taking into account that E_k = k² (Ryd). The only difference is that $[\bar \Sigma _2]$ is defined as a volume average of the imaginary part of the H–L self-energy (see Appendix B) instead of the line average found in this paper.

From the above discussion, we see that what in practice determines the photoelectron damping is a line average of the various MS paths. If the system is homogeneous enough, this average will not depend on the path, so in practice we can ignore the position dependence of the self-energy and replace it by a function of only the energy of the photoelectron: $[E = \hbar \omega - I_{\rm c}]$ . Since the cross section is analytical in the energy E, neglecting the cut at the Fermi energy, we can approximately take into account the effect of the absorptive part of the potential on the cross section $[\sigma _{\rm abs}^{\rm r} (E)]$ , calculated with the real part of the potential, by performing the following convolution

$[\eqalign{\sigma _{\rm abs}^{\rm c} (E) & = \pi^{-1} \textstyle\int\limits_{\varepsilon _{\rm F} }^\infty {\bar \Sigma _2 (E')} \,\sigma _{\rm abs}^{\rm r} (E') \, [(E - E')^2 + \bar \Sigma _2^2 (E')]^{-1} {\rm d} E' \cr & \simeq \sigma _{\rm abs}^{\rm r}\, [E - i\bar \Sigma _2 (E)],}]$

which holds if we are away from the edge, so that we can extend the integral to $[ - \infty]$ , and if $[\bar \Sigma _2 (E)]$ does not vary too rapidly with E. Explicit calculations via both methods substantiate this relation.

The effect of the core-hole lifetime can be easily incorporated into the theory by adding $[\Gamma _{\rm h}/ 2]$ to $[\bar \Sigma _2 (E)]$ , where $[\Gamma _{\rm h}]$ is the full width at half-maximum of the core hole and is related to its lifetime by $[\tau = \hbar /\Gamma _{\rm h}]$ . We assume here an exponential decay of the core hole, which is true in most of the cases that we are interested in. To estimate the effect of this added damping on the mean free path, we simply have to add $[\Gamma _{\rm h} /2]$ to $[\bar \Sigma _2 (E)]$ in (43), so that

$[\lambda _{\rm tot} = E_k / [k\,(\bar \Sigma _2 + \Gamma _{\rm h} /2)] \rightarrow \lambda _{\rm tot}^{ - 1} = \bar \Sigma _2 / k + \Gamma _{\rm h} /(2k) = \lambda _{\rm in}^{ - 1} + \lambda _{\rm h}^{ - 1}, \eqno(44)]$

where $[\lambda _{\rm in}]$ is the inelastic mean free path and $[\lambda _{\rm h}]$ is the mean free path corresponding to the finite lifetime of the core hole.

In the H–L approximation for the self-energy at energies lower than the plasmon energy, the mean free path would be infinite, if it were not limited by the core-hole lifetime, because of the plasmon-pole approximation that we have retained in the expression for the dielectric function of $[\Sigma _{\rm h}]$ in Appendix B. Damping originating from channels other than plasmon excitations (creation of electron–hole pairs) cannot be calculated by the H–L approximation. However, these damping terms have been evaluated (Hedin & Lundqvist, 1969) and found to be negligible compared with the core-hole lifetime and the experimental resolution. The damping terms can be introduced back into the equations using the formula given by Quinn & Ferrell (1958 ), as implemented in the FEFF code (Rehr & Albers, 2000). The photoelectron, on reaching the plasmon energy, interacts with the electron gas by creating plasmon excitations and suffering a reduction in the coherent wave-function amplitude. This fact would entail a rather sharp decrease in the photoelectron mean free path, which should show up as a localized feature in photoabsorption spectra. In reality, the setting in of the plasmon damping is not so sharp because of the quantum interference between intrinsic and extrinsic losses at the plasmon edge (Fujikawa, 1999; Campbell et al., 2002, and references therein). Therefore, we introduce into our fitting analysis an empirical model that takes into account both the smooth opening of the plasmon-loss channel and the effect of the electron–hole excitations.

2.4. The muffin-tin approximation

As a result of the reduction operated in §2.2, Green's function $[G({\bf r},{\bf r}';E)]$ obeys an effective one-particle Schrödinger equation, (25), which is better known as the Dyson equation. In (25) we have made explicit the appearance of the usual Coulomb or Hartree potential $[V_{\rm c} ({\bf r})]$ , which is given by

$[V_{\rm c} \left({{\bf r}} \right) = - \textstyle\sum\limits_k 2 Z_k \left| {{\bf r} - {\bf R}_k } \right|^{-1} + \, 2 \int {\rm d}^3 r' \rho ({\bf r}') \left| {\bf r} - {\bf r}' \right|^{-1}. \eqno(45)]$

Here $[{\bf R}_k]$ and Z_k indicate the position and the charge of the kth atomic nucleus and $[\rho ({\bf r})]$ is the charge density of the system under study, which, in an independent-particle scheme or in a local-density approach, is given by

$[\rho ({\bf r}) = \textstyle\sum\limits_i^{\rm occ} {\left| {\psi _i ({\bf r})} \right|} ^2.]$

In general, the solution of the Schrödinger–Dyson equation in three dimensions does not pose any particular problems in the framework of the MST, as illustrated in Appendix A. However, the generation of the local solutions $[\Phi _L ({\bf r})]$ and the calculation of the surface integrals $[C_{LL'}^i]$ and $[S_{LL'}^i]$ for all site and angular momentum indices are rather time consuming, especially if we have to generate many theoretical signals to compare with experimental data in a fitting procedure. A great simplification is achieved if we spherically average both the Coulomb potential and the self-energy $[\Sigma ^{\rm opt} ({\bf r},{\bf r}';E)]$ inside the atomic spheres. Both quantities depend on the electron density of the system under study, so, clearly, a self-consistent charge density would be highly desirable; however, this is often complicated to obtain in an amount of time that is appropriate to a fitting procedure. Therefore, usually the electron density is approximated by a superposition of spherically symmetric self-consistent atomic charge densities that are generated by currently available atomic programs. In order to make this approximation, we need to know how to expand a spherically symmetric function referred to one centre j around another centre o and take only the L = 0 component. This expansion is easily obtained based on the procedure proposed by Löwdin (1956 ) and later utilized by Mattheiss (1964 ). If $[\rho (r_j)]$ is the atomic radial charge density around centre j, normalized so that $[4\pi \textstyle\int_0^\infty {\rho (r_j) \, r_j^2 \,{\rm d}r_j } = Z_j]$ , then the component $[\rho _{L = 0} (r)]$ of the charge density at distance r from the centre o is given by

$[\rho _{L = 0} (r) = (2Rr)^{-1} \textstyle\int\limits_{|R - r|}^{R + r} {\rho (} r_j)\, r_j \,{\rm d}r_j, \eqno(46)]$

where R is the distance between the two centres. As a result, the total overlapped spherically averaged charge density around centre o is given by

$[\rho _{\rm tot}^o (r) = \rho ^o (r) + \textstyle\sum\limits_j (2R_j r)^{-1} \int\limits_{|R_j - r|}^{R_j + r} {\rho (} r_j) \, r_j \, {\rm d}r_j, \eqno(47)]$

where R_j is the distance of neighbour j from centre o and $[\rho ^o (r)]$ is the charge density of the atom at this centre. Similarly, because of the linearity of the Poisson equation, it is easy to see that the spherically averaged potential around centre o, V_c^o (r), generated by this overlapped charge density is given by

$[V_{\rm c}^o (r) = V^o (r) + \textstyle\sum\limits_j (2R_j r)^{-1} \int\limits_{|R_j - r|}^{R_j + r} {V^j (} r_j) \, r_j \, {\rm d}r_j, \eqno(48)]$

where V^j (r) is the atomic potential generated by the radial charge density $[\rho (r_j)]$ of the atom at site j:

$[V^j (r) = (2/r)\textstyle\int\limits_0^r {4\pi r_j^2 \rho (r_j })\, {\rm d}r_j + 2\int\limits_r^\infty {4\pi r_j \rho (r_j)} \, {\rm d}r_j. \eqno(49)]$

Now, the spirit of the muffin-tin (MT) approximation consists of partitioning the molecular space (i.e. the space occupied by the atomic cluster under consideration) into three regions. Region (I) is made up of atomic spheres around the physical atoms, and the radii of these spheres are determined as described below. Region (III) lies outside an outer sphere circumscribing the cluster, and the interstitial region (II) lies between the outer and the atomic spheres. Empty spheres can be added to (II) to minimize the interstitial volume and to better account for the local variation of the true potential. This latter is spherically averaged into the spheres and approximated in regions (II) and (III) to a constant, which is set equal to the volume average of the potential in region (II). The true long-range variation of the potential in region (III) (Coulomb tail, in the case of isolated molecules) or a different constant (if the cluster is embedded in a different medium) can be easily incorporated into the scheme without any problems, although these corrections are usually neglected because of the damping effects of the imaginary part of the optical potential when the radius of the outer sphere is greater than the photoelectron mean free path. The interstitial constant is easily calculated by expanding the potential around the centre of the outer sphere (and centre of the cluster). Assuming that o is such a centre and denoting $[\Delta \Omega]$ as the volume of the interstitial region, we can use (48) to obtain

$[\eqalignno{\bar V_{\rm c} & = \Delta \Omega^{-1} \textstyle\int\limits_{\Delta \Omega} {V_{\rm c} ({\bf r})} \,{\rm d}^3 {\bf r} \cr & = 4\pi\, \Delta \Omega^{-1} \left[ {\textstyle\int\limits_0^{R_o } {V_{\rm c}^o (r)\, r^2 \,{\rm d}r} - \textstyle\sum\limits_j {\int\limits_0^{R_j } {V^j (r_j)\, r_j^2 \,{\rm d}r_j }} } \right], &(50)}]$

where R_o is the radius of the outer sphere. Similarly for the constant interstitial charge we find

$[\eqalignno{\bar \rho _{\rm int} & = (\Delta \Omega)^{-1} \textstyle\int\limits_{\Delta \Omega } {\rho ({\bf r})} \,{\rm d}^3 {\bf r} \cr & = 4\pi \Delta \Omega^{-1} \textstyle\left[{\int\limits_0^{R_o } {\rho _{\rm tot}^o (r) \, r^2 \,{\rm d}r} - \textstyle\sum\limits_j {\int\limits_0^{R_j } {\rho ^j (r_j)\, r_j^2 \,{\rm d}r_j } } } \right]. & (51)}]$

Equations (47)–(51) define the Coulomb potential and the charge density throughout the molecular system in its muffin-tin form. Via the density, the exchange-correlation part of the electron self-energy can be calculated at each photoelectron energy, either in its local (H–L) or in the more elaborate non-local (GW) form.

In order to determine the MT radii of the atomic spheres, we can follow one of two prescriptions: one given by Norman (1974 ) and the other by Wille et al. (1986 ). In the Norman scheme, a Norman radius R_i^Nor is determined for each site i such that

$[\textstyle\int\limits_0^{R_i^{\rm Nor} } {4\pi \, r_i^2 \rho _{\rm tot}^i (r_i)}\, {\rm d}r = Z_i,]$

which is roughly proportional to the ionic radius of the atom in the periodic table.

Then, given two nearest-neighbour sites i and j, the touching-sphere MT radii are given by the following formula,

$[R_{i(j)} = [R_{i(j)}^{\rm Nor}\,R_{ij}]/ (R_i^{\rm Nor} + R_j^{\rm Nor}),]$

which implies that their ratio equals the ratio of the corresponding Norman radii. Empirically, this prescription works in the case of covalent bonds, so that both spheres should be full. In the case of several neighbours at slightly different distances, we can take average interatomic distances in place of R_ij. The use of averages might entail some overlap between the various MT spheres, but that overlap is not considered to be a problem. Actually, in practical calculations, we can allow an overlap of about 10–15% of the MT radii to empirically simulate the covalent bond and make up for the approximation of the potential in the interstitial region. In many cases this procedure seems to improve the agreement between theory and experiment, even though it is arbitrary on theoretical grounds. In the case of empty spheres or H atoms, ad hoc prescriptions can easily be devised.

On the other hand, Wille et al. (1986) suggest choosing the MT radii in such a way that the potential discontinuities at the boundaries of the MT spheres are minimized. This prescription seems to work for ionic compounds. In any case, even though the potential discontinuities can be reduced to an acceptable number, discontinuities remain between the boundaries of the atomic spheres and the interstitial region. These latter produce unphysical scattering, especially at photoelectron energies near the threshold, and can deform the shape of the calculated spectrum. Fortunately, by using the freedom of slightly varying the interstitial Coulomb potential and charge density and the overlap factor between MT spheres, we can nearly always reach a good agreement with experiments. However, the entire situation is highly unsatisfactory and calls for the elimination of the MT approximation.

In the generation of the charge density $[\rho ({\bf r})]$ , care must be taken in many cases to reproduce the relaxation of the atomic electrons following the creation of the core hole and the screening action of the valence electrons. This relaxation is a dynamical process that, for deep core holes, is simulated in a static way by considering for the photoabsorbing atom the SCF charge density obtained by promoting the core electron to the first non-occupied valence orbital. This simulation seems to work for many systems (metals, semiconductors etc.) for which electronic correlation in the final state can be neglected, including also the biological compounds that we are interested in.

Given the MT potential and exploiting the ensuing simplifications, it is then easy to generate the local regular and irregular solutions of the Schrödinger–Dyson equation and Green's function, as detailed in Appendices A and C. In this way, polarized and unpolarized spectra can be calculated to fit to the experimental data. In the above procedure it would be a step forward to use self-consistent charge densities (and therefore potentials) in the MT form. Indeed, starting from Green's function, a self-consistent loop could be initiated by first generating a new charge density based on the formula

$[\eqalign{\rho ^1 ({\bf r}) & = - (\pi)^{-1}{\rm Im} \textstyle\int\limits_{ - \infty }^{\varepsilon _{\rm F} } G^ + ({\bf r},{\bf r};E)\,{\rm d}E \cr & = (2\pi)^{-1} \left[ {\textstyle\int\limits_{- \infty} ^{\varepsilon _{\rm F} }} G^ - ({\bf r},{\bf r};E)\,{\rm d}E - \textstyle\int\limits_{\varepsilon _{\rm F }}^{ - \infty } G^ + ({\bf r},{\bf r};E)\,{\rm d}E \right] \cr & = (2\pi)^{-1} \textstyle\int\limits_L {G({\bf r},{\bf r};E)\,{\rm d}E},}]$

where L is any path in the complex energy plane containing the cut on the negative real axis due to the valence states and the poles due to the inner core states. By deforming the path so that it comprises two vertical lines at $[\varepsilon _{\rm F}]$ (the Fermi level) and $[\varepsilon _{\rm C}]$ (the energy value below any core state) and two horizontal lines at $[E + i\Gamma _{\rm u(l)}]$ in the upper and lower planes with an appropriate $[\Gamma]$ , the contour integral could easily be calculated. The Fermi level $[\varepsilon _{\rm F}]$ is easily determined by the condition

$[\textstyle\int\limits_V \rho ({\bf r})\,{\rm d}{\bf r} = \int\limits_V {\rm d}{\bf r}(2\pi)^{-1}\int\limits_L G({\bf r},{\bf r};E)\,{\rm d}E = N_{\rm el},]$

where V is the molecular volume (outside which the charge is negligible) and N_el is the total number of electrons. The new charge so obtained could serve to generate a new Coulomb and exchange-correlation potential (always of the MT form) that would generate a second Green's function, so as to repeat the loop until self-consistency was reached. This procedure would substantially improve the reliability of the fitting procedure at the cost of increasing the fitting time, although we could choose to repeat the loop only when really necessary. The self-consistent procedure has been implemented in the package FEFF8 [Ankudinov et al. (1998 ), to which we refer the reader for details] and seems to generate good final charge densities.

3. Future improvements and conclusions

The current status of techniques for calculating the inner-shell photoabsorption cross section for clusters, based on the real-space MS method and using the local density approximation for exchange and correlation, has been reviewed by many authors. A theoretical review mainly focused on EXAFS and to a lesser extent on XANES with emphasis on the path-by-path method has been given by Rehr & Albers (2000), whereas XANES analysis has been treated in the reviews by Kizler (1993 ), Binsted & Hasnain (1996 ), Ankudinov (1999 ) and, in the special case of narrow-band highly correlated materials, de Groot (1996 ). In particular, Binsted & Hasnain (1996) have advocated the fitting of the entire XAS spectrum rather than its components, EXAFS and XANES, and have reached the conclusion that this is a viable procedure. To our knowledge, they were the first to indicate such a procedure. These authors analyse four representative model compounds and by a comparison of matrix-inversion and finite-path-sum methods they reach the conclusion that the latter method is more promising for fitting the edge region. Our preference is for the matrix-inversion method, because the MS series might not converge near the edge, in which case the maximum order of scattering becomes a further parameter without physical significance, which could lead to false minima in the R factor. Rather, we should attempt to reduce the impact of the many theoretical approximations present nowadays in the calculation of XAS spectra. We summarize in the following the motivations behind this choice, having in mind not only biological applications but also the many X-ray spectroscopies that are used to obtain information on the many aspects of condensed matter physics.

Indeed, experience has shown that the sophistication of present-day spectroscopies has reached a point where a simple interpretation of the data is not sufficient to extract all the information that they conceal. Structural EXAFS analysis of absorption spectra is one such instance (Filipponi et al., 1995 ), and in the domain of surface physics photoelectron diffraction is another example (Chen et al., 1998). In such cases the possibility of calculating a good theoretical signal has been of paramount importance in establishing reliable tools for structural analysis, and the success of these techniques is due to the fact that the corresponding theories are asymptotic in the energy parameter. In fact, at high photoelectron energies simple approaches, like the description of the photoemission process in terms of one electron moving in an optical potential of the local density type and the muffin-tin approximation, are usually sufficient for a realistic description of systems of general type.

Unfortunately, while structural information dominates the high-energy spectral region, both structural and especially electronic information are present in the low-energy region. It would be extremely useful to have a reliable theory to help us extract the relevant information from the data, especially since the most interesting and promising spectroscopies for unraveling electronic structure (like X-ray magnetic and natural dichroism, resonant elastic and inelastic scattering, both magnetic and non magnetic) give most of their electronic information near the excitation edge. In this energy region, the muffin-tin approximation is no longer adequate for describing the geometrical details of the potential, since the excited electron is quite sensitive to them. Moreover, complicated many-body processes intervene to screen the core hole and to interfere with the direct excitation channel, as described in §1. Indeed, electron-correlation effects in the final state and the coupling of the excited photoelectron with the core hole concur to complicate the description of the near-edge region in many spectroscopies. Last but not least, the effects of atomic vibrations should be included in the theoretical simulations if we extend the fitting procedure beyond the near edge, since the damping effect of these vibrations interferes with the electronic damping due to the imaginary part of the optical potential.

Despite many efforts in the past to tackle these problems, they have not yet received a satisfactory solution. This is the reason why in many spectroscopies the near-edge region is not adequately exploited to extract the rich information, both structural and electronic, present in it. This latter information would be extremely valuable in the case of biological systems, concerning, for example, the spin and charge state of metal centres. Based on the above considerations, efforts should be concentrated in the following areas.

3.1. Inclusion of vibrational effects in XAS spectra

In the structural analysis of XAS spectra, much attention has been devoted to the inclusion of vibrational effects in the general terms of the MS series, but little effort has been put in to incorporating such effects directly into the MS matrix. The reason is both historical and methodological. On the one hand, interest was concentrated in analysing the spectral region where the MS series is converging and the damping effect of the atomic vibrations is sizable. On the other, the inclusion of these effects in the near-edge spectral region was theoretically more difficult and estimated to be of minor importance. From the point of view of fitting both the near-edge and the MS region of the spectrum using full matrix inversion, this inclusion becomes essential. Fortunately, related research areas, like low-energy electron diffraction (LEED), had to tackle the same problem, and various methods have been proposed to cope with it, up to the point that computer subroutines are available to calculate effective t-matrix elements that take into account in an approximate way the effect of atomic vibrations [see de Andres & King (2001 ) for details and references]. It would be very useful to incorporate such developments into the photoabsorption programs and continue the research in the field. Also, the availability of an experimentally determined vibrational spectrum would be highly desirable, since this would reduce the number of quantities to fit. Note that the experimental EXAFS data can support only a limited number of fitted parameters (Stern, 1993 ; Filipponi, 1995a ; Michalowicz & Vlaic, 1998 ). In this respect several attempts have been made to avoid the need to consider the Debye–Waller (DW) factors as parameters in the fitting procedure. These DW factors can be either calculated using semi-empirical and ab initio methods that consider a small cluster around the absorbing atom (Dimakis & Bunker, 1998 ; Poiarkova & Rehr, 1999 ) or derived from other types of experiments as in the work of Loeffen & Pettifer (1996 ). These authors performed an EXAFS calculation for the zinc tetraimidazole molecular cluster, which included thermal four-body correlations established from inelastic neutron scattering. In this case full vibrational information has permitted a quantitative comparison between theory and experiment, without recourse to fitting. They found that differences between theory and experiment persist, setting the limits on the systematic errors still present in the theory. Presumably the MT approximation (not good for short atomic bonds) and the over-damping of the imaginary part of the H–L potential (not suitable for such low-Z covalent systems) are the causes of the remaining discrepancy.

3.2. Elimination of the muffin-tin (MT) approximation in the solution of the Schrödinger–Dyson equation

This stage is preliminary to all spectroscopies, since a reliable solution of the one-electron Dyson equation moving in an effective optical potential is necessary for obtaining the related response functions. A method based on the finite-difference method (FDM) (and therefore free of the MT approximation) and a program for solving the associated equation are already available. This method has been developed by Joly (2001 ), and a collaboration with the Frascati group in the past few years has led to interesting applications, ranging from the electron population analysis of near-edge absorption spectra (Joly et al., 1999 ) to the interpretation of resonant elastic scattering in systems showing orbital ordering (Benfatto et al., 1999 ). However, such applications have also shown the method's practical limitations, which are due to the large memory and the long CPU time required by the computer code, especially for low-symmetry systems. In this scheme, fitting procedures to experimental data become prohibitive and many applications to interesting systems in biology, earth sciences, chemistry and physics are out of reach. An alternative method has been developed by Foulis et al. (1990 ), based on the non-MT (full potential) multiple-scattered-wave theory of Natoli et al. (1986). This method it is not completely free from approximations in treating the interstitial region. Nevertheless, this method, together with self-consistent-field electron densities, has demonstrated that quantitative accuracy for molecules and clusters of modest size can be achieved (Foulis et al., 1995 ). Probably, a combination of the two methods, as illustrated in Appendix A, is the way to eliminate the drawbacks of both. Moreover, in view of applications to magnetic dichroism, spin–orbit interaction and spin-polarized potentials can be included and indeed have been incorporated by Joly into his code (to be published).

3.3. Implementation of the multichannel MST (MMST) in the photoemission process

As already discussed, in the near-edge spectral region (either in absorption or in photoemission) intrinsic inelastic phenomena (shake-up, shake-off) very often modify the transition amplitude due to the primary excitation channel. The same phenomena (hole–particle excitation) also contribute to the dynamical screening of the core hole following the photoemission process and determine the type of effective optical potential acting in the final state. Likewise the coupling phenomena between the core hole and the excited electron (multiplet effects) need a multielectron description. The presence of all these effects has hindered up to now a satisfactory exploitation of the electronic and structural information contained in the near-edge spectra. The scheme usually employed in the literature to deal with the electron-correlation problems rests on the theory of atomic multiplets in a crystal field (CF). However, this approach is only static and relies on CF parameters. A way to deal with both the static and the dynamic correlation problem is provided by the multichannel multiple-scattering theory developed by Natoli et al. (1990). This theory is, in principle, `ab initio' and is able to incorporate into its description the local atomic multiplet aspect of the problem and the extended character of the excited photoelectron wave function. However, the theory has never been implemented in a computer code. If successful, this implementation should lead to a breakthrough in the analysis of near-edge spectroscopies concerning electronic properties.

3.4. Core hole dynamic screening and the problem of the optical potential in the final state

As we have shown, the calculation of the total absorption cross section (or the elastic channel in photoemission spectroscopy) can be reduced to a one-particle problem moving in an effective optical potential by elimination of the inelastic channels in the MMST. Since the problem is rather complicated, we usually resort to empirical ad hoc rules to construct such a potential. In particular, the most widely used rule is based on the work of Von Barth & Grossman (1982 ), which prescribes the use of a relaxed fully screened charge to calculate both the Coulomb and the exchange-correlation potential of the H–L type (the so-called Z + 1 screened approximation). However, in order to convince ourselves that the situation is more complicated, it suffices to consider some absorption spectra. The Z + 1 approximation may be necessary for the description of the K-edge absorption in diatomic molecules (N₂, O₂) in order to bind the $[\pi ^*]$ resonance in the final state, since this state lies in the continuum in the ground state of the molecule. However, this approximation is not necessary for the L_II, L_III circular absorption edges in the ferromagnetic phase of nickel. In this case the dichroic effect would disappear if we were to use a self-consistent spin-polarized relaxed potential, as has been checked by direct numerical calculations. Clearly the optical final-state potential is the result of a dynamic screening process that depends on the system and on the type of final state reached. In particular, in the case of nickel, the self-screening of the core hole by the excited photoelectron and the Pauli principle hinder the electronic relaxation and screening by the other electrons of the system, so that the optical potential closely resembles the initial-state unrelaxed and unscreened potential. In general, in the case of electronic correlations in the final state, the true amount of screening charge should be determined based on a dynamical theory, as discussed in §3.3 above. It is very important to have a theoretical model for the screening process since it affects the relative position of the pre-edge and rising-edge peaks in the absorption spectra.

In summary, even though decisive advances in the understanding of inner-shell X-ray absorption have been made over the past few decades, much remains to be done in order to obtain a versatile tool for extracting electronic and structural information from experimental spectra. However, we have also given theoretical arguments, supported by the applications made up to now (Benfatto et al., 2003), that even with the present approximations structural analysis by full-spectrum fitting is possible. Electronic analysis, like spin and charge state, electron population analysis etc., should probably await further improvements.

APPENDIX A

Derivation of the MS equations

We derive here the solution of the Dyson-like equation

$[[\nabla ^2 + k^2 - V_{\rm c} ({\bf r})]\,\varphi ({\bf r}) = \textstyle\int {\Sigma ^{\rm opt} } ({\bf r},{\bf r}';\hbar \omega)\,\varphi ({\bf r}')\,{\rm d}^3 r',\eqno(52)]$

subject to the asymptotic scattering boundary conditions

$[\varphi ({\bf r}) \simeq (k / \pi)^{1/2} \exp(i{\bf k} \cdot {\bf r}) + f(\hat r;{\bf k})\exp(ikr) / r \eqno(53)]$

in the framework of multiple-scattering theory and in the non muffin-tin case.

The main ingredients for this derivation are the single and two-centre expansions of the free Green's function with the following outgoing wave boundary conditions:

(a) Around one centre located at the origin,

$[G_0^ + ({\bf r} - {\bf r}';E) = - {{\exp (ik|{\bf r} - {\bf r}'|)}\over{4\pi |{\bf r} - {\bf r}'|}} = \textstyle\sum\limits_L {J_L ({\bf r}_ \lt)} \tilde H_L^+ ({\bf r}_ \gt), \eqno(54)]$

with $[\tilde H_L^+ ({\bf r}) = - ikH_L^ + ({\bf r})]$ and the same definitions as in §2.1.

(b) Around two centres located at $[{\bf R}_i]$ and $[{\bf R}_j]$ . Defining $[{\bf r}_i = {\bf r} - {\bf R}_i]$ and $[{\bf R}_{ij} = {\bf R}_i - {\bf R}_j]$ , we have, provided $[R_{ij} \, \gt \, r_i + r_j]$ ,

$[G_0^ + ({\bf r} - {\bf r}';E) = \textstyle\sum\limits_{LL'} {J_L ({\bf r}_i) }\, G_{LL'}^{i j} (E) \, J_{L'} ({\bf r}_j), \eqno(55)]$

where

$[G_{LL'}^{i j} = 4\pi \textstyle\sum\limits_{L''} i ^{l - l' + l''} C_{LL'}^{L''}\, \tilde H_{L''}^ + ({\bf R}_{ij}) \eqno(56)]$

are the real-space KKR structure factors already introduced by (18) in §2.1.

$[G_0^ + ({\bf r} - {\bf r}';E) = \textstyle\sum\limits_{LL'} {J_L ({\bf r}_i) } \,J_{LL'}^{i o} (E) \, \tilde H_{L'}^+ ({\bf r}_o),\eqno(57)]$

with

$[J_{LL'}^{i o} = 4\pi \textstyle\sum\limits_{L''} i ^{l - l' + l''} C_{LL'}^{L''}\, J_{L''}^{} ({\bf R}_{io}). \eqno(58)]$

We assume here that $[E \,\gt\, 0]$ in order to treat scattering states; however, the above relations still hold for $[E \,\lt \,0]$ , if we take the analytic continuation of the functions of energy appearing in the above equations and if we set the incoming plane wave to zero in (53). In this way, within the same formalism we can also treat bound states.

For simplicity we consider only short-range potentials (the extension to long-range potentials is, however, straightforward, as will be apparent from the following). There exists a volume $[\Omega _o]$ that encloses all the atoms of our system, and this volume is sufficiently large that in $[C\Omega _o]$ the solution of (52) is of the form of (53). We partition this `molecular' volume in N non-overlapping cells $[\Omega _j]$ with centres at $[{\bf R}_j]$ (these cells might even be empty, i.e. not enclosing a physical atom). The cells completely fill $[\Omega _o]$ (i.e. without interstices, so that $[\Omega _o = \textstyle\sum_{j = 1}^N {\Omega _j }]$ ) in such a way that:

(a) There is a finite, but small, neighbourhood around the origin of each cell that lies completely within the cell. If there is an atom inside the cell, its nucleus should coincide with the origin.

(b) The shortest inter-cell vector, joining the origins of nearest-neighbour cells, is larger than any intra-cell vector.

We then start from the obvious equality involving surface integrals

$[\eqalignno{\textstyle\sum\limits_{j = 1}^N & \,\,{\textstyle\int\limits_{S_j } {\left [{G_0^ + ({\bf r} - {\bf r}'){\bf \nabla} \varphi ({\bf r}') - \varphi ({\bf r}'){\bf \nabla} G_0^ + ({\bf r} - {\bf r}')} \right]} } \cdot {\bf n}_j \,{\rm d}\sigma _j \cr & - \textstyle\int\limits_{S_o } {\left [{G_0^ + ({\bf r} - {\bf r}'){\bf \nabla} \varphi ({\bf r}') - \varphi ({\bf r}'){\bf \nabla} G_0^ + ({\bf r} - {\bf r}')} \right]} \cdot {\bf n}_o \,{\rm d}\sigma _o = 0, & (59)}]$

in which the meanings of the symbols are obvious. This equality is valid for all $[{\bf r}]$ , provided that $[\varphi ({\bf r})]$ is continuous with its first derivatives, as required by the solution of the Schrödinger-like equation.

By taking $[{\bf r}]$ inside $[\Omega _i]$ such that $[r_i \rightarrow 0]$ , we use (54) in the surface integral over S_i, (55) in the integral over S_j and (57) in the integral over S_o. We then find

$[\eqalignno{\textstyle\sum\limits_L & J_L ({\bf r}_i) \Big\{\textstyle\int\limits_{S_i}\left[\tilde H_L^ + ({\bf r}_i){\bf \nabla} \varphi ({\bf r}_i) - \varphi ({\bf r}_i){\bf \nabla} \tilde H_L^ + ({\bf r}_i) \right] \cdot {\bf n}_i \,{\rm d}\sigma _i \cr & + \textstyle\sum\limits_{j \ne i} {\textstyle\sum\limits_{L'} {G_{LL'}^{i j} } } \int\limits_{S_j } {\left[{J_{L'} ({\bf r}_j){\bf \nabla} \varphi ({\bf r}_j) - \varphi ({\bf r}_j){\bf \nabla} J_{L'} ({\bf r}_j)} \right]} \cdot {\bf n}_j \,{\rm d}\sigma _j \cr & - {\textstyle\sum\limits_{L'} J_{LL'}^{i o} \int\limits_{S_o } \left [{\tilde H_{L'}^ + ({\bf r}_o){\bf \nabla} \varphi ({\bf r}_o) - \varphi ({\bf r}_o){\bf \nabla} \tilde H_{L'}^ + ({\bf r}_o)} \right] \cdot {\bf n} _o \,{\rm d}\sigma _o } \Big\} = 0, & (60)}]$

so that, because of the angular completeness of the set $[J_L ({\bf r})]$ , the expression inside the brackets {} should be zero for each L.

Now, inside each cell $[\Omega _j]$ , we introduce basis functions $[\Phi\, _L^j ({\bf r}_j)]$ , which are solutions of the Dyson-like equation (52) and behave at the origin like $[J_L ({\bf r}_j)]$ . These basis functions constitute a complete set so that the general solution $[\varphi ({\bf r})]$ can be locally expanded as $[\varphi ({\bf r}_j) = \textstyle\sum_L {A_L^j ({\bf k})\,\Phi\, _L^j ({\bf r}_j)}]$ . Likewise, in the outer domain $[C\Omega _o]$ , in order to impose the boundary conditions (53), we can take

$[\varphi ({\bf r}_o) = \textstyle\sum\limits_L {\left[ {{\underline A}_{\,L}^o ({\bf k})\,J_L ({\bf r}_o) + C_L^o ({\bf k})\,\tilde H_L^+ ({\bf r}_o)} \right]}, \eqno(61)]$

where $[{\underline A}_{\,L}^o ({\bf k}) = 4\pi i^l\, Y_L ({\hat k})\,(k / \pi)^{1/2}]$ , owing to the well known decomposition of a plane wave. Inserting these expressions into (60) using the relation

$[\textstyle\int\limits_{S_o } {\left [{\tilde H_{L'}^ + ({\bf r}_o){\bf \nabla} J_L ({\bf r}_o) - J_L ({\bf r}_o){ \nabla} \tilde H_{L'}^ + ({\bf r}_o)} \right]} \cdot {\bf n}_o \, {\rm d}\sigma _o = \delta _{LL'}]$

and the identity

$[\textstyle\sum\limits_{L'} {J_{LL'}^{i o} i^{l'} Y_{L'} (\hat k)} = i^l Y_L (\hat k)\exp(i\,{\bf k} \cdot {\bf R}_{io}),]$

which is obtained from (58) by observing that $[\textstyle\sum_{L'} {C_{LL'}^{L''} } \,Y_{L'} ({\hat k}) = Y_L ({\hat k})\,Y_{L''} ({\hat k})]$ , we derive (15) of §2.1,

$[\textstyle\sum\limits_{L'} {} C_{LL'}^i \,A_{L'}^i ({\bf k}) = A_L^o ({\bf k}) - \textstyle\sum\limits_{jL'L''} {(1 - \delta _{ij}) } \,G_{LL'}^{i j} \,S_{L'L''}^j \,A_{L''}^j ({\bf k}),\eqno(62)]$

where the definitions (16) and (17) give the surface integrals $[C_{LL'}^i]$ and $[S_{LL'}^i]$ .

In the case of long-range potentials (or if the short-range potential is substantially different from zero in the outer region) we have to supplement (62) with a similar equation obtained from (59) by taking $[{\bf r} \in C\Omega _o]$ to describe the effect of the scattering of the external potential. However, the solution inside the central cell is very often rather insensitive to the potential in the external region, provided the radius of the cluster is greater than the mean free path of the photoelectron at the energy considered. Clearly the above formalism can be easily extended to treat bound states.

The solution of (62) proceeds via the calculation of the surface integrals $[C_{LL'}^i]$ and $[S_{LL'}^i]$ and the structure factors $[G_{LL'}^{i j}]$ . These latter are routinely calculated and, in principle, there should be no problem in calculating the surface integrals $[C_{LL'}^i]$ and $[S_{LL'}^i]$ . However, their practical calculation for polyhedral cells of general shape, though feasible, might be rather complicated. For this reason, it is preferable to use cells of spherical shape, inscribed in the polyhedra, especially since the quantity $[T_{LL'}^i = \textstyle\sum_{L''} {S_{LL''}^i (C^i)_{L''L'}^{ - 1} }]$ , which was introduced in §3.1, can then be interpreted as the scattering amplitude of the cell. The filling of space in this case leaves the problem of the interstitial region (IR), which can be minimized by increasing the number of empty spheres. However, in order to treat the irregular shape of the remaining IR we are obliged to approximate the potential by a constant (e.g. its volume average). By absorbing this constant into the definition of energy the interstitial potential becomes zero, so that the surface integrals $[C_{LL'}^i]$ and $[S_{LL'}^i]$ and (62) remain unchanged by application of Green's theorem. A further simplification, though not a necessary one, is achieved by replacing the potential inside each sphere by its angular average, so that we have to deal with spherically symmetric potentials. In this case the basis functions $[\Phi _L^i ({\bf r}_j)]$ can be written as $[R_l^i (r_i)Y_L (\hat r_i)]$ , where R_lⁱ (r_i) is the solution of the radial Schrödinger equation inside sphere i. Consequently the surface integrals $[C_{LL'}^i]$ and $[S_{LL'}^i]$ reduce to

$[C_{LL'}^i = k\delta _{LL'} R_{\rm s}^2 \, W [- ih_l^ +, R_l^i] |_{r = R_{\rm s} }]$

and

$[S_{LL'}^i = \delta _{LL'} \, R_{\rm s}^2 \, W[\,j_l, R_l^i] |_{r = R_{\rm s} },]$

where we have introduced the Wronskian of two functions f(r) and g(r), $[W[\,f,g] = fg' - gf']$ , calculated at the sphere radius R_s. The atomic scattering amplitude $[T_{LL'}^i = \textstyle\sum\limits_{L''} {S_{LL''}^i (C^i)_{L''L'}^{ - 1} }]$ , also known as the atomic t-matrix, takes the form

$[T_{LL'}^i = \delta _{LL'} \,t\,_l^i = \delta _{LL'} k^{-1} W[\,j_l, R_l^i]/ W [- ih_l^ +, R_l^i], \eqno(63)]$

which is the expression that would be found by solving the scattering problem for a spherical wave of angular momentum l impinging on a spherical potential truncated at radius r = R_s. From scattering theory we also have $[t\,_l^i = k^{-1}\exp(i\delta _l ) \sin \delta _l]$ , where $[\delta _l]$ is the phase shift of the radial solution R_l (r) with respect to the free solution j_l (r) caused by the potential. If this latter is real, conservation of flux requires that $[k| {t_l } |^2 = {\rm Im} \,t_l]$ (optical theorem), so that $[\delta _l]$ is real; otherwise $[k| {t_l } |^2 \,\lt\, {\rm Im}\, t_l]$ , which implies that $[\delta _l]$ is complex and that the difference $[{\rm Im}\, t_l - k| {t_l } |^2]$ is related to the loss of flux due to the absorptive part of the potential, which becomes the source of the damping of the electronic wave (see §2.3).

APPENDIX B

The Hedin–Lundqvist potential

Much work has gone into approximating the exchange-correlation part of the optical potential $[\Sigma ^{\rm opt} ({\bf r},{\bf r}';E)]$ in a way suitable for numerical applications. Hedin & Lundqvist (1969, 1971), by incorporating the Sham–Kohn (Sham & Kohn, 1966 ) density-functional formalism for excited states within the single-plasmon-pole approximation of the electron gas dielectric function, have produced a useful local-density approximation to $[\Sigma]$ :

$[\Sigma {^{\rm opt} \left({{\bf r},{\bf r}';E} \right)} \simeq \Sigma_{\rm h} {\left[\,{p\left({{\bf r}} \right),E - V_{\rm c} \left({{\bf r}} \right);\rho \left({{\bf r}} \right)} \right]} \equiv V_{\rm exc} \left({{\bf r}} \right). \eqno(64)]$

Here $[\Sigma _{\rm h} (\,p,\omega; \rho)]$ is the self-energy of an electron in a homogenous interacting electron gas with momentum $[p = p({\bf r})]$ , energy $[\omega = E - V_{\rm c} ({\bf r})]$ and density $[\rho = \rho ({\bf r})]$ , the local density of the actual physical system.

Since $[E - V_{\rm c} ({\bf r}) \simeq p^2 ({\bf r})]$ , neglecting the small exchange and correlation correction, we can write for the local exchange and correlation potential (following Lee & Beni, 1977)

$[V_{\rm exc} \left({{\bf r}} \right) \simeq \Sigma_{\rm h} {\left[\,{p\left({{\bf r}} \right),p^2 \left({{\bf r}} \right) ; \rho \left({{\bf r}} \right)} \right]}. \eqno(65)]$

The local momentum $[p({\bf r})]$ is defined as

$[\eqalignno{p^2 ({\bf r}) & = k^2 + k_{\rm F}^2 ({\bf r}) - \Sigma_{\rm h} {(\,{p,p^2; \rho} ) }+ \Sigma { _{\rm h} } ({k_{\rm F}, k_{\rm F}^2; \rho} ) \cr & = k^2 + k_{\rm F}^2 ({{\bf r}} ) + \Delta \Sigma, & (66)}]$

where k² is the photoelectron kinetic energy measured from the Fermi level for an extended system or from the first ionization potential for finite systems (Hara, 1967) (atoms and molecules) and $[k_{\rm F}^2 ({{\bf r}} ) = [{3\pi ^2 \rho ({{\bf r}} )} ]^{2/3}]$ is the local Fermi energy. Usually, we can omit to perform the self-consistent procedure implicit in (66) for the determination of $[p({\bf r})]$ , except perhaps near the Fermi energy.

We follow Lee & Beni (1977) in extending the validity of (64) to the atomic core region, although the H–L approach was initially devised to describe exchange and correlation corrections to the Hartree potential due to the valence charge.

To calculate $[\Sigma _{\rm h} (\,p,\omega; \rho)]$ , we use equations (25.1), (25.14) and (25.15) of Hedin & Lundqvist (1969), which are based on the so-called GW approximation:

$[\eqalignno{\Sigma_{\rm h} (\,{p,\omega } ) = {}&i ({2 \pi } )^{-4} \textstyle\int \exp(i\omega '\delta) \,V(\,{{\bf p}'} ) \, \varepsilon (\,{{\bf p}',\omega '} )^{-1} \,\cr & \times G_0 \left(\,{{\bf p} + {\bf p}',\omega + \omega '} \right) \,{\rm d}^3 p'\,{\rm d}\omega '. & (67)}]$

Equation (67) corresponds to the self-energy of the test electron interacting with the charge fluctuations of the medium. Here $[G_0 (\,p,\omega) = [{\omega - p^2 + i\delta \sin(\omega - \varepsilon _{\rm F})} ]^{ - 1}]$ is the propagator for the test electron at energy $[\omega]$ measured from the Fermi energy $[\varepsilon _{\rm F}]$ ; $[\delta]$ is an infinitesimal positive quantity that is necessary for imposing the correct boundary conditions so that the integral converges when $[\omega \rightarrow \infty]$ .

The Fermi energy $[\varepsilon _{\rm F}]$ of the electron gas is given by $[\varepsilon _{\rm F} = k_{\rm F}^2]$ , where

$[k_{\rm F} = (3\pi ^2 \rho)^{1/3} = (\beta & r_{\rm s})^{ - 1}]$

is the Fermi momentum, which corresponds to the constant density $[\rho]$ of the homogeneous electron gas [eventually to be taken equal to the local density $[\rho ({\bf r})]$ of the inhomogeneous system]. The mean inter-particle distance r_s is given by

$[r_{\rm s} = \left[3 /(4\pi \rho) \right]^{1/3}]$

so that

$[\beta = \left[4 / (9\pi) \right]^{1/3}\,\,\simeq 0.52.]$

Moreover, the function $[\varepsilon (\,p,\omega)]$ is the frequency-dependent dielectric function of the electron gas in the plasmon-pole approximation (Lee & Beni, 1977):

$[\left [{\varepsilon (\,{p,\omega } )} \right]^{-1} = 1 + \omega _{\rm p}^2 \left [{\omega ^2 -\omega _1^2 (\,p )} \right]^{^{-1} },]$

where

$[\omega _{\rm p} = 4\left[(\,\beta r_{\rm s})/(3\pi)\right]^{1/2} \varepsilon _{\rm F} = 41.7 \left[ r_{\rm s}\, ({\rm a.u.} ) \right]^{-3/2} \,{\rm eV}]$

is the plasmon energy and

$[\omega _1^2 \left(\, p \right) = \omega _{\rm p}^2 + \varepsilon _{\rm F}^2 \left [(4/ 3)(\,p/k_{\rm F } )^2 + (\,p / k_{\rm F })^4 \right]]$

is the plasmon-pole dispersion relation. We have taken the coefficient 4/3 following Filipponi et al. (1988 ) rather than Lee & Beni (1977).

By performing the ω′ integral of (67), we obtain

$[\eqalign{\Sigma _{\rm h } (\,p,\omega ) = & - {\textstyle\int} {\rm d}^3 q \, \, (4\pi) \, f({\bf q} + {\bf p}) /\left\{ { (2\pi)^{3} q^{2}\, \varepsilon \left [{{\bf q},({\bf q} + {\bf p})^2 -\omega }\right] }\right\} \cr & + \omega _{\rm p}^2 \textstyle\int {\rm d}^3 q \, (4\pi ) / \left\{{ (2 \pi)^3 q^2 \, \left[\omega -\omega _1 (q)-({\bf q} + {\bf p})^2\right] } \right\} ,}]$

and hence

$[\eqalignno{{\rm Re}\Sigma_{\rm h} (\,p,\omega ) = {}&- \textstyle\int {\rm d}^3 q \, (4\pi) \, f({\bf q} + {\bf p}) / \left\{ {(2\pi)^3\,q^2\, \varepsilon \left [{{\bf q},({\bf q} + {\bf p})^2 -\omega } \right]} \right\} \cr & + \omega _{\rm p}^2 \textstyle\int {\rm d}^3 q \, (4\pi) / \left\{ { (2\pi)^3 q^2 \left[ {\omega -\omega _1 (q)-({\bf q} + {\bf p})^2 }\right]} \right\}, \cr&& (68)}]$

$[\eqalignno{{\rm Im} \Sigma_{\rm h} (\, p,\omega ) = & \,(\pi / 2) \, \omega _{\rm p} \textstyle\int \Big( {\rm d}^3 q \, (4\pi) \,\left[ (2\pi)^3 \, q^2\, \omega _1 (q) \right]^{-1} \cr & \times \left\{ f ({\bf q} + {\bf p}) \,\delta \,\left[ ({\bf q} + {\bf p})^2 - \omega _1(q) - \omega \right] \right. \cr & \,- \left. \left[1 - f({\bf q} + {\bf p}) \right]\, \delta\, \left[({\bf q} + {\bf p})^2 + \omega _1 (q) - \omega\right] \right\} \Big), &(69)}]$

where we should take the principal part of the integral in (68); $[f(\bf q)]$ is the Fermi distribution function. To proceed further we introduce the dimensionless quantities

$[X = {q \over {k_{\rm F}}},\quad E_r = {\omega \over {\varepsilon _{\rm F} }}, \quad X_k = {p \over {k_{\rm F} }} = \left [{1 + \left({{{k^2 } \over {k_{\rm F}^2 }}} \right) + \left({{{\Delta \Sigma } \over {k_{\rm F}^2 }}} \right)} \right]^{1/2},]$

so that

$[\omega _1 (q) /\varepsilon _{\rm F} = 2\left [{w_{\rm p}^2 + X^2/3 + X^4/4 } \right]^{1/2} \equiv 2 w\,\left(X \right)]$

with

$[w_{\rm p}^2 = 4 \beta r_{\rm s} / 3\pi.]$

Calculating the above integrals at $[\omega ^2 = p^2]$ (i.e. E_r = X_k²) we obtain for $[{\rm Re} \Sigma _{\rm h} (\,p,p^2)]$

$[\eqalignno{{1 \over {\varepsilon _{\rm F} }} {\rm Re}\Sigma _{\rm h} (\,p,p^2 ) = {}&- {{2\beta r_{\rm s} } \over \pi }\left [{1 + {{1 - X_k^2 } \over {2 X_k }} \,\ln {{X_k + 1} \over {X_k - 1}}} \right] + {{2\beta r_{\rm s} } \over \pi } \cr & \times\left\{ {{2\beta r_{\rm s} } \over {3\pi }}\,{1 \over {X_k }}\left [\int\limits_{X_k - 1}^{X_k + 1} F_{\rm sex} \left({X\semi X_k } \right)\, {\rm d}X\right.\right.\cr& \left.\left.\,- \int\limits_0^\infty F_{\rm ch} \left({X\semi X_k } \right)\,{\rm d}X \right] \right\}, &(70)}]$

where

$[\eqalign{F_{\rm sex} (X;X_k) = & \, [X \, w(X)]^{-1} \cr & \times\ln \left| { {{\left [{2w(X) + X_k^2 - 1} \right]\left [{2w(X) + ({X - X_k })^2 - X_k^2 } \right]} \over {\left [{2w(X) + X_k^2 + 1} \right]\left [{2w(X) - ({X - X_k })^2 + X_k^2 } \right]}} } \right|}]$

and

$[F_{\rm ch} (X;X_k ) = [X \, w(X)]^{-1} \ln {{\left [{2w(X ) + ({X - X_k })^2 - X_k^2 } \right]} \over {\left [{2w(X ) + ({X + X_k })^2 - X_k^2 } \right]}}.]$

Note that we have been considering excitations above the Fermi level, so that $[X_k \,\gt \,1]$ . We recognize in the first term of (70) the usual static and energy-dependent exchange obtained by Dirac (1930 ) and discussed by Hara (1967):

$[V_{\rm ex}^{\rm DH} = - {{2 k_{\rm F} } \over \pi }\left [{1 + {{1 - X_k^2 } \over {2 X_k }} \ln \left({{{X_k + 1} \over {X_k - 1}}} \right)} \right].]$

This comes from the constant part of the dielectric function in the first integral of (68). When multiplied by $[3\alpha /2]$ , where $[\alpha]$ is a parameter, we obtain the usual form of the X- $[\alpha]$ Slater static exchange potential (Slater, 1979 ; Schwarz, 1972 ),

$[V_{\rm ex}^{\rm S} = - 3 \alpha k_{\rm F} / \pi.]$

Usually a parametrized form is also assumed for V_ex^DH,

$[V_{\rm ex}^{{\rm DH}\alpha } = - {{3 \alpha k_{\rm F} } \over \pi }\left [{1 + {{1 - X_k^2 } \over {2 X_k }} \ln \left({{{X_k + 1} \over {X_k - 1}}} \right)} \right],]$

in such a way that V_ex^DH reduces to the X- $[\alpha]$ form for X_k = 1.

As already noted by Chou et al. (1985 ),

$[V_{{\rm ex}}^{{\rm DH}} \sim k_{\rm F}^3\, p^{-2} \simeq \rho k^{-2 }\quad(k \, \gg \,k_{\rm F })]$

for excitation energies much larger than the Fermi energy ( $[X_k \to \infty]$ ). Thus, for high-energy photoelectrons V_ex^DF becomes proportional to the density $[\rho]$ and decreases inversely with energy, instead of being proportional to k_F (or $[\rho ^{1/3}]$ ) and independent of energy, as for the X- $[\alpha]$ exchange V_ex^S. This behaviour is to be expected since, with increasing energy, the Pauli exclusion principle is less and less active, owing to the fact that at high energies the wave function of the test electron is already nearly orthogonal to that of the electrons in the system. Several authors have found that this energy dependence is crucial for a good calculation of EXAFS phase shifts (Chou et al., 1985, 1987 , and references therein).

The next contribution in (70) still comes from the first integral in (68) but contains the effect of the frequency dependence of the dielectric function. Therefore, this contribution is a kind of screened exchange. The contribution is positive and roughly constant in the interval $[1 \, \le \, X_k \, \le \, 2]$ and then decays rapidly like

$[V_{\rm sex} \mathop \sim\limits_{X_k \to \infty } \,16 \,\left[3\pi ^2 \, X_k^2 \, w(X_k) \right]^{-1} \simeq 16 \,(3\pi ^2 \, X_k^4 )^{-1}]$

since

$[\mathop {\lim }\limits_{X_k \to \infty} \,\textstyle\int\limits_{X_k - 1}^{X_k + 1} { F_{\rm sex} (X; X_k)} \, {\rm d}X = - 4 \,\left[w(X_k) \, X_k \right]^{-1} \sim - 4 X_k^{-3 }.]$

(This asymptotic relation actually holds already for $[X_k \,\ge \,3]$ ).

Finally, the last contribution in (70) comes from the last term in (68), which in turn arises from the poles of the inverse dielectric function. Therefore, this contribution is associated with a Coulomb correlation hole and stems from the fact that the test electron tends to keep away electrons of both spins because of the Coulomb repulsion. The contribution has a negative sign, like the exchange, but cannot be given an analytic form in terms of elementary functions. However, we can calculate the behaviour of the term at high energies,

$[\textstyle\int\limits_0^\infty {F_{\rm ch} } (X;X_k) \,{\rm d}X \mathop \sim\limits_{X_k \to \infty } - w_{\rm p}^{-1 } \left[{ \pi ^2 / 2 + 2 w_{\rm p} X_k^{-1} + O (X_k^{-2 })} \right],]$

so that the net contribution to V_exc is

$[V_{\rm ch} \mathop \sim\limits_{X_k \to \infty } - 4\, (3\pi ^2 \,X_k \,w_{\rm p})^{-1} \left[{ \pi ^2 / 2 + 2 w_{\rm p}X_k^{-1} + O (X_k^{-2})} \right].]$

This behaviour shows that the effect of the Coulomb correlation hole is more persistent at high energy than the D–H exchange. This effect causes the difference in the phase shifts that are calculated with the two potentials.

The numerical calculation of (70) does not pose particular problems, although it is a little time consuming. A good numerical approximation, valid in the entire energy range of the excited photoelectron, has been given by Mustre de Leon et al. (1991). Both options are incorporated into our computer code.

It is found that the imaginary part in (69) is amenable to quadrature. In fact, for $[X_k \,\ge\, 1]$ and $[\omega = p^2]$ the integral becomes

$[{1 \over {\varepsilon _{\rm F} }} {\rm Im} \Sigma_{\rm h} (\,p, p^2) = {{4(\beta r_{\rm s})^2 } \over {3 \pi }} {1 \over {X_k }} \int\limits_{0}^{X_{\rm M} } { {{{\rm d}X} \over {X \, w(X)}}} \,\Theta \left({X;X_1, X_2 } \right),\eqno(71)]$

where

$[X_{\rm M} = \left\{ { - (2/ 3) + 2\left [{(1/9) - w_{\rm p}^2 + ({X_k^2 - 1} )^2/4 } \right]^{1/2} } \right\}^{1/2}]$

is the maximum real root of the equation

$[w^2 (X) - ({X_k^2 - 1})^2/4 = 0.]$

Therefore, $[{\rm Im} \Sigma _{\rm h} (\,p,p^2)]$ is non-zero only if $[X_k^2 - 1 \,\gt\, 2w_{\rm p}]$ , i.e. for

$[k^2 \,\gt \,2 w_{\rm p} \varepsilon _{\rm F} = 4\,(\,\beta r_{\rm s} /2\pi)^{1/2} \varepsilon _{\rm F} = \omega _{\rm p}.]$

In other words there is no damping if the electron kinetic energy is not sufficient to excite the plasmon mode of the medium. Lower-energy excitations, like particle–hole excitations, are not possible because of the assumption of the plasmon-pole approximation for the dielectric function $[\varepsilon (\,p,\omega)]$ . Lower-energy excitations are negligible in a first approximation (see Hedin & Lundqvist, 1969, page 93). The function $[\Theta (X;X_1, X_2)]$ is defined to be 1 for $[X_1 \, \le \, X \,\le\, X_2]$ and zero outside. The interval [X₁, X₂] is any interval, lying inside [0,X_M], where $[w(X)\, \lt\, X(X_k - X/2)]$ or

$[f (X) \equiv X_k X^3 + (1/3 - X_k^2) \,X^2 + w_{\rm p}^2 \,\lt\, 0.]$

For $[X_k^2 - 1 \,\gt \,2w_{\rm p}]$ , this cubic polynomial has three real roots, of which only one, X₁, lies in the interval [0,X_M] and is such that $[f(X)\, \le\, 0]$ in [X₁, X_M]. Therefore, the integral becomes

$[\eqalignno{ \varepsilon_{\rm F}^{-1} \, {\rm Im} \Sigma_{\rm h} (\,p, p^2) & = 4 (\,\beta r_{\rm s })^2\, (3 \pi \, X_k)^{-1} \textstyle\int\limits_{ X_1 }^{X_{\rm M} } {\rm d}X \, \left[{X \, w(X)}\right]^{-1} \cr & = 3\pi\, (\,\beta r_{\rm s}/ 3\pi )^{3/2} X_k^{-1}\cr&\quad\times \ln \left [{ (X_{\rm M} /X_1)^2\, F(X_1) / F(X_{\rm M}) }\right], & (72)}]$

where

$[F(X) = X^2 + 6 w_{\rm p}^2 + w_{\rm p} w(X). \eqno(73)]$

The root X₁ is given by

$[X_1 = \left[{(3 X_k^2 - 1) / (3 X_1 )}\right] \left[{1 - \cos \,(\varphi / 3) + 3^{1/2} \sin \,(\varphi / 3)} \right]/3,]$

where

$[0\, \le \,\varphi = ar\cos \left\{ {1 - [27 \, w_{\rm p}^2 / (2\, X_k)] \left[{3 X_k / (3 X_k^2 - 1)} \right]^3 } \right\} \le \pi.]$

We could also exploit the facts that X_M is a solution of the equation w(X) = (X_k² - 1)/2 and X₁ is solution of w(X) = X(X_k - X/2) to eliminate w(X) in (73).

Equation (73) is the same as the expressions that were obtained by Penn (1987), apart from a volume average and a correction for some misprints, and by Mustre de Leon et al. (1991). In the spirit of the local-density approximation (66), the expressions (70) and (71) for the real and imaginary parts of $[\Sigma _{\rm h} (\,p,p^2; \rho)]$ become functions of the position $[{\bf r}]$ through the same dependence on the effective density $[\rho ({\bf r})]$ of the system under study.

APPENDIX C

Non-local core polarization potential

In order to calculate the screened polarization potential G^v W^v P^c W^v we need the explicit expression for the full RPA (random phase approximation) polarization propagator (Hedin & Lindquist, 1969),

$[P({\bf r},{\bf r}'\semi\omega) = - \textstyle\sum\limits_k^{\rm unocc} \textstyle\sum\limits_l^{\rm occ} \, 2\,(\varepsilon _k - \varepsilon _l\,) \left[ (\varepsilon _k - \varepsilon _l\,)^2 - \omega ^2 \right]^{-1} \,f_{kl} ({\bf r})\,f_{kl}^* ({\bf r}') ,]$

$[f_{kl} ({\bf r}) = \textstyle\int {\psi _k (x)\, \psi _l^* (x)\,{\rm d}\xi }. \eqno(74)]$

Here $[x = ({\bf r},\xi)]$ includes both space and spin variables, and $[\varepsilon _k]$ are the one-electron excitation energies of the system. The sum over k runs over unoccupied electron states, while l runs over the occupied core and valence electron states. By splitting the summation over l into core and valence contributions, P can be written as a sum of core and valence parts,

$[P = P\,^{\rm c} + P\,^{\rm v_o }, \eqno(75)]$

where v_o means occupied valence states. Similarly, we can split the summation over k in the expression for the one-electron Green's function,

$[G(x,x';\omega) = \textstyle\sum\limits_k^{\rm core + valence} \psi _k (x)\,\psi _k^* (x') /(\omega - \varepsilon _k ), \eqno(76)]$

to obtain

$[G = G\,^{\rm c} + G\,^{\rm v}. \eqno(77)]$

In this case, the index c runs again over the core states while v runs over both the occupied and the unoccupied valence states.

The question to consider is the effect of screening in the core-polarization term G ^v W ^v P ^c W ^v. For a free atom, the screening is small and we have to a good approximation G ^v VP ^c V, where V is the bare Coulomb potential. For long distances from the ion core, the core-polarization term further reduces to the well known local potential $[- \alpha e^2 r^{-4 }]$ , where $[\alpha]$ is the dipole polarizability and e is the electron charge (Hedin, 1965a ,b ). In a solid we expect, from simple physical considerations, that for long distances we should have a statically screened polarization potential, G ^v W ^v (0)P ^c W ^v (0). Since an r^{- 4} potential is already very weak, the additional screening should make it negligible outside the Wigner–Seitz cell of the ion under consideration. Inside the Wigner–Seitz cell, on the other hand, we do not expect much screening to take place because of the cost in kinetic energy that is needed to localize the screening charge.

The symbol G ^v VP ^c V stands for a convolution in energy space (Fujikawa & Hedin, 1989 ), which can be performed analytically, giving

$[[G\,^{\rm v} VP\,^{\rm c} V](x,x';\omega) = \sum\limits_k^{\rm unocc^{} } {\sum\limits_l^{\rm core^{} } {\sum\limits_{k'}^{\rm valence} {{{V_{kl} ({\bf r})\,\psi _{k'} (x)\,\psi _{k'}^* (x')\,V_{kl}^* ({\bf r}')} \over {\omega - \omega _{kl} - \varepsilon _{k'} }}} } }, \eqno(78)]$

where $[\omega _{kl} = \varepsilon _k - \varepsilon _l]$ and $[V_{kl} ({\bf r}) = \textstyle\int {V({\bf r} - {\bf r}')\,\psi _k^* (x')\,\psi _l (x')\,{\rm d}x'}]$ . The more tightly bound the core level l, the smaller its contribution to $[V_{kl} ({\bf r})]$ , because the overlap with the unoccupied function k is smaller. Thus the outermost core level will give the dominant contributions.

We replace $[\omega _{kl}]$ by a constant $[\Delta]$ , the average excitation energy. This approximation has been very successful in the free-atom case (Byron & Joachain, 1974 , 1977a ,b ). We define a function $[A({\bf r},{\bf r}')]$ ,

$[\eqalignno{A({\bf r},{\bf r}') & = \textstyle\sum\limits_k^{\rm unocc} {\textstyle\sum\limits_l^{\rm core} {V_{kl} ({\bf r})\, V_{kl}^* ({\bf r}')} } \cr & = \textstyle\int V({\bf r} - {\bf r}_1) \, V({\bf r}' - {\bf r}_2)\, [\delta (x_1 - x_2) - \rho (x_1, x_2)]\,\cr&\quad\times \rho ^{\rm c} (x_2, x_1)\, {\rm d}x_1 {\rm d}x_2, &(79)}]$

where the last equality follows by closure and $[\rho]$ and $[\rho ^{\rm c}]$ are the one-electron density matrices for all electrons and for core electrons, respectively. With $[\omega _{kl} = \Delta]$ we then have

$[\Sigma ^{\rm pol} (x,x';\omega) = \left [{G^{\rm v} VP^{\rm c} V} \right] (x,x';\omega) = A({\bf r},{\bf r}') \, G^{\rm v} (x,x';\omega - \Delta).\eqno(80)]$

With the use of closure we avoid the summation over the unoccupied states. Still, the density matrix $[\rho (x_1, x_2)]$ contains a sum over the occupied extended states. We will here take a simplified approach and represent the sum over Bloch functions in each atomic cell i by a sum over localized functions $[R_l^i (r)Y_{lm} (\hat r)]$ . This approach is well motivated for rare-gas solids but a more serious approximation for, say, metals. In future, a more accurate treatment should be considered.

From (80) and the expression for Green's function (26), we see that the range of $[\Sigma ^{\rm pol} (x,x';\omega)]$ is apparently not limited to atomic dimensions. However, there are damping effects present in Green's function, and the energy-independent term $[A({\bf r},{\bf r}')]$ defined by (79) decreases fairly rapidly with the distance $[|{\bf r} - {\bf r}'|]$ [see (101) below]. Therefore, in (26) we can neglect, to a good approximation, terms coming from spheres other than the sphere under consideration (labelled i). Since, for simplicity, we shall spherically average the optical potential inside each sphere, we can then write (see Appendix A and §2.2)

$[G^i ({\bf r},{\bf r}';E) = \textstyle\sum\limits_L {g_l^i (r,r';E)} \,Y_L (\hat r)\, Y_L (\hat r'),\eqno(81)]$

with

$[g_l^i (r,r';E) = \tilde R_l (r_ \lt) \, t_l^i \, \tilde R_l^ + (r_ \gt),\eqno(82)]$

where $[\tilde R_l (r) = R_l (r)/[R_{\rm s}^2 \, W(j_l, R_l)_{r = R_{\rm s} }]]$ is the regular solution of the radial Schrödinger equation with the optical potential matching smoothly to $[j_l (kr)\,(t_l^i)^{ - 1} + \tilde h_l^ + (kr)]$ at the sphere radius R_s; $[\tilde R_l^ + (r)]$ is the irregular solution matching smoothly to $[\tilde h_l^ + (kr) = - ikh_l^ + (kr)]$ , $[r_ \gt = \max (r,r')]$ , $[r_ \lt = \min (r,r')]$ and t_lⁱ is the atomic t-matrix for the atom at site i.

Therefore, in order to obtain an expression for the atomic optical potential at this site, we need to calculate the energy-independent term $[A({\bf r},{\bf r}')]$ . In each atomic region i we average the core charge density to obtain a spherically symmetric function

$[d\,_i^{\rm c} (r) = (4\pi)^{-1} \textstyle\int {{\rm d}\hat r\rho ^{\rm c} ({\bf r})}.\eqno(83)]$

Though both $[{\bf r}]$ and $[{\bf r}']$ are in the same atomic region, $[{\bf r}_1]$ and $[{\bf r}_2]$ are not necessarily in the same atomic region. From the first term in square brackets in (79) we have an intraatomic contribution to $[A({\bf r},{\bf r}')]$ :

$[I_i ({\bf r},{\bf r}') = \textstyle\sum\limits_L {\left[{4\pi /(2l + 1)} \right]^2 \, I_i (r,r')\,Y_L (\hat r)\,Y_L^* (\hat r')},\eqno(84)]$

where both $[{\bf r}]$ and $[{\bf r}']$ belong to atomic region i and $[{\bf r}_1 (= {\bf r}_2)]$ is also in region i. [I_i (r,r')_l] is

$[I_i (r,r'\,)_l = S_1^l + S_2^l + S_3^l,\eqno(85)]$

where

$[S_1^l = (rr'\,)^{ - l - 1} \textstyle \int\limits_0^{r_ \lt } {d\,_i^{\rm c} (r_1)\,r\,_1 ^{2l + 2} \,{\rm d}r_1 }, \eqno(86)]$

$[S_2^l = (r_ \lt ^l) / (r_ \gt ^{l + 1}) \textstyle\int\limits_{r_ \lt }^{r_ \gt } {d\,_i^{\rm c} (r_1)\, r_1 \,{\rm d}r_1 },\eqno(87)]$

$[S_3^l = (rr'\,)^l \textstyle\int\limits_{r_ \gt }^{R_i } {d\,_i^{\rm c} (r_1)\,r_1 ^{ - 2l}\, {\rm d}r_1 },\eqno(88)]$

and R_i is the radius of the atomic region i (muffin-tin radius).

In addition to this intraatomic contribution to the first term of (79) we have the interatomic term, where $[{\bf r}]$ and $[{\bf r}']$ are still in region i, whereas $[{\bf r}_1 (= {\bf r}_2)]$ is in region j $[(\,j \ne i)]$ . This interatomic term is given by

$[\eqalignno{I_i^j ({\bf r},{\bf r}') & = \textstyle\sum\limits_{LL'} \left[{4\pi/ (2l' + 1)} \right]^2 |F_{LL'} ({\bf R}_{ji})|^2 \,(rr')^l \,\langle {r^{2l'} } \rangle _j \, Y_L (\hat r)\,Y_L^* (\hat r') \cr & = \textstyle\sum\limits_L {F_L (r,r';{\bf R}_{ji})\,Y_L (\hat r)\,Y_L^* (\hat r')},& (89)}]$

where in terms of $[G(LL'|L'')\equiv \textstyle \int {Y_L (\hat r)\, Y_{L'} (\hat r)\, Y_{L''}^* (\hat r)\,{\rm d}\hat r}]$ (Gaunt's integral)

$[\eqalignno{F_{LL'} ({\bf R}_{ji}) = &{{4\pi (- 1)^{l'} (2l + 2l' - 1)!!} \over {(2l' - 1)!!(2l + 1)!!R_{ji}^{l + l' + 1} }} \cr & \times G(l + l',m' - m,L|L')\,Y_{l + l',m' - m} (\hat R_{ji}), & (90)}]$

$[F_L (r,r';{\bf R}_{ji}) = (rr'\,)^l \textstyle\sum\limits_{L'} {\left[{4\pi /(2l' + 1)} \right]^2\, |F_{LL'} ({\bf R}_{ji})|^2 \,\langle {r\,^{2l'} } \rangle _j },\eqno(91)]$

$[\langle {r\,^{2l} } \rangle _j = \textstyle\int\limits_0^{R_j } {d_j^{\rm c} (r)\,r\,^{2l + 2}\,{\rm d}r},\eqno(92)]$

$[{\bf R}_{ji} = {\bf R}_j - {\bf R}_i.]$

We note that $[\textstyle\sum_{m'} {|F_{LL'} ({\bf R}_{ji})|^2 }]$ depends on m. To avoid this difficulty, we use the spherically averaged value of $[\textstyle\sum_{m'} {|F_{LL'} ({\bf R}_{ji})|^2 }]$ . Finally, we have the averaged expression of F_L in (89) in terms of the Clebsch–Gordan coefficient $[\langle {{l0l'0}} | {{l + l'0}} \rangle]$ , which is denoted F_l⁰:

$[\eqalignno{F_l^0 (r,r'\semi R_{ji}) = {}&(rr')^l \sum\limits_{l'} \left[ {{4\pi (2l + 2l' - 1)!!} \over {(2l + 1)!!(2l' + 1)!!R_{ji}^{l + l' + 1} }} \right]^2 \cr &\times \langle r^{2l'} \rangle _j \,(2l' + 1)\, \langle l0l'0 |l + l'0 \rangle ^2 . & (93)}]$

The second term in square brackets in (79) is more difficult to calculate in general. We note that

$[\rho (x_1, x_2)\,\rho ^{\rm c} (x_2, x_1) = \textstyle\sum\limits_\alpha ^{\rm occ} {\textstyle\sum\limits_\beta ^{\rm core} {d_{\alpha \beta } (x_1)\, d_{\alpha \beta }^* (x_2)} }, \eqno(94)]$

where $[d_{\alpha \beta } (x) = \varphi _\alpha (x)\,\varphi _\beta ^* (x)]$ and we assume that the core orbital $[\alpha]$ is localized on site i. We take both $[{\bf r}_1]$ and $[{\bf r}_2]$ to be in the same region i. For the occupied states, we will not use Bloch functions but instead take a simplified approach and use localized functions $[R_l^i (r)\, Y_{lm} (\hat r)]$ . This approach is well motivated for rare-gas solids but a more serious approximation for, say, metals.

By spherically averaging $[d_{\alpha \beta } (x)]$ at each site i, we can obtain a simple representation for (94) as

$[\rho (x_1, x_2)\rho ^{\rm c} (x_2, x_1)\sim \textstyle\sum\limits_m^{\rm core} {d_m^i (r_1)\,d_m^i (r_2)} + \textstyle\sum\limits_m^{\rm core} {\textstyle\sum\limits_{n(\ne m)}^{\rm occ} {d_{mn}^i (r_1)\,d_{mn}^i (r_2)} } \eqno(95)]$

$[{\rm for}\quad{\bf r}_1, {\bf r}_2 \in i,]$

where m and n stand for one-electron atomic states in atom i that have the same angular quantum number L = (l,m), e.g. 2s, 3s and so on. The quantity d_mⁱ ( = d_mmⁱ) is the spherically averaged electron density of the mth atomic function at site i, and d_mnⁱ is the cross charge from the mth and nth atomic functions on site i written in terms of the radial parts of the atomic wave functions, d_mnⁱ = R_mⁱ (r)R_nⁱ (r)^*. R_mⁱ refers to a localized function, which for metals has a fractional occupation number. When we use this simple approximation in (79), the intraatomic contribution to A can be written as

$[J_i (r,r') = \textstyle\sum\limits_m {J_m^i (r)\,J_m^i (r')} + \textstyle\sum\limits_{m \ne n} {J_{mn}^i (r)\,J_{mn}^i (r')},\eqno(96)]$

with

$[J_m^i (r) = 4\pi \left[ r^{-1} \textstyle\int\limits_0^r {r_1 ^2 \, d_m^i (r_1)\, {\rm d}r_1 + \int\limits_r^\infty {r_1\, d_m^i (r_1)\,{\rm d}r_1 } } \right]. \eqno(97)]$

We see that J_i has only a spherically symmetric contribution.

The interatomic contribution J_i^j is not difficult to evaluate:

$[J_i^j (r,r') = \textstyle\sum\limits_m^{\rm core} {J_m^{ij} (r)\,J_m^{ij} (r')}, \eqno(98)]$

where, by taking r_i = r, J_m^ij (r) is given as

$[J_m^{ij} (r) = (4\pi)^{-1} \textstyle\int {{\rm d}\hat r_i \, n_m^j \, |{\bf r}_i - {\bf R}_{ji} |^{-1}} \quad (i \ne j). \eqno(99)]$

Here n_m^j is the number of electrons on site j and $[{\bf R}_{ji} = {\bf R}_j - {\bf R}_i]$ . Cross terms such as d_mn^j give no contribution to J_m^ij because of the orthogonality between the mth and nth shell functions. Therefore, $[A({\bf r},{\bf r}')]$ can be written as

$[A({\bf r},{\bf r}') = \textstyle\sum\limits_L {A(r,r')_l\,Y_L } (\hat r)\,Y_L^* (\hat r'),\eqno(100)]$

where A is a sum of one- and two-centre terms,

$[A(r,r')_l = A_i (r,r')_l + \textstyle\sum\limits_{j \ne i} {A_i^j (r,r')_l }, \eqno(101)]$

$[A_i (r,r')_0 = (4\pi)^2\, [I_i (r,r')_0 - 4\pi J_i (r,r')_0], \eqno(102)]$

$[A_i (r,r')_l = [4\pi/ (2l + 1)]^2 \, I_i (r,r')_l \quad (l \ge 1),\eqno(103)]$

$[A_i^j (r,r')_l = {{(4\pi)^2 } \over 3}\left({{{l + 1} \over {2l + 1}}} \right){{\langle {r^2 } \rangle _j } \over {R_{ji}^{2l + 4} }}(rr')^l + \cdots.\eqno(104)]$

We find a good convergence for the two-centre sum, the second term on the right-hand side in (101), when we include the surrounding atoms up to the third shell for the systems considered here.

The optical potential can be given by (105) after the spherical averaging of the potential over $[\hat r]$ and $[\hat r']$ in the same atomic region:

$[\Sigma ^{\rm pol} ({\bf r},{\bf r}';E_0 + \varepsilon _{\rm p}) = \textstyle\sum\limits_L {\Sigma _l^{\rm pol} (r,r';\varepsilon _{\rm p})\,Y_L } (\hat r')\,Y_L^* (\hat r'), \eqno(105)]$

where $[\Sigma _l^{\rm pol}]$ is expressed in terms of $[A_{l'}]$ , $[g_{l''}]$ and Clebsch–Gordan coefficients,

$[\Sigma _l^{\rm pol} (r,r';\varepsilon _p) = {{1}\over{4\pi}} \sum\limits_{l'l''} {{{(2l' + 1)(2l'' + 1)} \over {2l + 1}}\langle {{l'0l''0}} | {{l0}} \rangle ^2 A_{l'} (r,r')g_{l''}^i (r,r';\bar p)}. \eqno(106)]$

As previously demonstrated, it is enough to include only A₀ and A₁ in the expansion (106) (Fujikawa et al., 1993 ). We thus obtain an explicit expression for $[\Sigma _l^{\rm pol}]$ ,

$[\Sigma _l^{\rm pol} (r,r';\varepsilon _p) = (4\pi)^{-1} \left [{A_0 (r,r')\,g_l (r,r';\bar p) + A_1 (r,r')\, \tilde g_l (r,r';\bar p)} \right], \eqno(107)]$

where $[\tilde g_l]$ is defined by

$[\tilde g_l (r,r';\bar p) = {l \over {2l + 1}}g_{l - 1} (r,r';\bar p) + {{l + 1} \over {2l + 1}}g_{l + 1} (r,r';\bar p).\eqno(108)]$

Each g_l includes the radial solution for the potential $[\Sigma _l^{\rm pol} (\varepsilon _{\rm p})]$ to be determined. Since $[\Sigma _l^{\rm pol}]$ depends on g_l and $[g_{ \pm l}]$ , in principle we have to solve coupled self-consistent equations.

References

Andres, P. L. de & King, D. A. (2001). Comput. Phys. Commun. 138, 281–301. Web of Science CrossRef Google Scholar
Ankudinov, A. L. (1999). J. Synchrotron Rad. 6, 236–238. Web of Science CrossRef CAS IUCr Journals Google Scholar
Ankudinov, A. L., Ravel, B., Rehr, J. J. & Conradson, S. D. (1998). Phys. Rev. B, 58, 7565–7576. Web of Science CrossRef CAS Google Scholar
Benfatto, M. & Della Longa, S. (2001). J. Synchrotron Rad. 8, 1087–1094. Web of Science CrossRef CAS IUCr Journals Google Scholar
Benfatto, M., Della Longa, S. & Natoli, C. R. (2003). J. Synchrotron Rad. 10, 51–57. CrossRef CAS IUCr Journals Google Scholar
Benfatto, M., Joly, Y. & Natoli, C. R. (1999). Phys. Rev. Lett. 83, 636–639. Web of Science CrossRef CAS Google Scholar
Benfatto, M., Natoli, C. R., Bianconi, A., Garcia, J., Marcelli, A., Fanfoni, M. & Davoli, I. (1986). Phys. Rev. B, 34, 5774–5781. CrossRef CAS Web of Science Google Scholar
Binsted, N. & Hasnain, S. S. (1996). J. Synchrotron Rad. 3, 185–196. CrossRef CAS Web of Science IUCr Journals Google Scholar
Breit, G. & Bethe, H. (1954). Phys. Rev. 93, 888–890. CrossRef Web of Science Google Scholar
Byron, F. W. & Joachain, J. C. (1974). Phys. Rev. A, 9, 2559–2568. CrossRef CAS Web of Science Google Scholar
Byron, F. W. & Joachain, J. C. (1977a). Phys. Rev. A, 15, 128–146. CrossRef CAS Web of Science Google Scholar
Byron, F. W. & Joachain, J. C. (1977b). J. Phys. B, 10, 207–226. CrossRef CAS Web of Science Google Scholar
Campbell, L., Hedin, L., Rehr, J. J. & Bardyszewski, W. (2002). Phys. Rev. B, 65, 064107-1–064107-13. Google Scholar
Chen, Y., Garcia de Abajo, F. J., Chassé, A., Ynzunza, R. X., Kaduwela, A. P., Van Hove, M. & Fadley, C. S. (1998). Phys. Rev. B, 58, 13121–13131. Web of Science CrossRef CAS Google Scholar
Chou, S. H., Kutzler, F. W., Ellis, D. E., Shenoy, G. K., Morrison, T. I. & Montano, P. A. (1985). Phys. Rev. B, 31, 1069–1076. CrossRef CAS Web of Science Google Scholar
Chou, S. H., Rehr, J. J., Stern, E. A. & Davidson, E. R. (1987). Phys. Rev. B, 35, 2604–2614. CrossRef CAS Web of Science Google Scholar
Della Longa, S., Arcovito, A., Girasole, M., Hazeman, J. L. & Benfatto, M. (2001). Phys. Rev. Lett. 87, 155501-1–155501-4. Google Scholar
Dirac, P. A. M. (1930). Proc. Cambridge Philos. Soc. 26, 376–385. CrossRef CAS Google Scholar
Dimakis, N. & Bunker, G. (1998). Phys. Rev. B, 58, 2467–2475. Web of Science CrossRef CAS Google Scholar
Faulkner, J. S. & Stocks, G. M. (1980). Phys. Rev. B, 21, 3222–3241. CrossRef CAS Web of Science Google Scholar
Filipponi, A. (1995a). J. Phys. Condens. Matter, 7, 9343–9356. CrossRef CAS Web of Science Google Scholar
Filipponi, A. (1995b). Physica B, 208/209, 29–32. CrossRef Web of Science Google Scholar
Filipponi, A., Bernieri, E. & Mobilio, S. (1988). Phys. Rev. B, 38, 3298–3304. CrossRef CAS Web of Science Google Scholar
Filipponi, A., DiCicco, A. & Natoli, C. R. (1995). Phys. Rev. B, 52, 15122–15134. CrossRef CAS Web of Science Google Scholar
Foulis, D. L., Pettifer, R. F., Natoli, C. R. & Benfatto, M. (1990). Phys. Rev. A, 41, 6922–6927. CrossRef CAS PubMed Web of Science Google Scholar
Foulis, D. L., Pettifer, R. F. & Sherwood, P. (1995). Europhys. Lett. 29, 647–652. CrossRef CAS Web of Science Google Scholar
Fujikawa, T. (1999). J. Phys. Soc. Jpn, 68, 2444–2456. Web of Science CrossRef CAS Google Scholar
Fujikawa, T., Hatada, K. & Hedin, L. (2000). Phys. Rev. B, 62, 5387–5398. Web of Science CrossRef CAS Google Scholar
Fujikawa, T. & Hedin, L. (1989). Phys. Rev. B, 40, 11507–11518. CrossRef Web of Science Google Scholar
Fujikawa, T., Saito, A. & Hedin, L. (1993). Jpn. J. Appl. Phys. 32, 18–22. CrossRef CAS Web of Science Google Scholar
Groot, F. de (1996). X-ray and Inner-Shell Processes: Seventeenth International Conference, Hamburg, Germany, edited by R. L. Jonson, H. Schmidt-Böcking & B. F. Sonntag, p. 497. New York: AIP. Google Scholar
Gunnella, R., Bullock, E. L., Patthey, L., Natoli, C. R., Abukawa, T., Kono, S. & Johansson, L. S. (1998). Phys. Rev. B, 57, 14739–14748. Web of Science CrossRef CAS Google Scholar
Gunnella, R., Solal, F., Sébilleau, D. & Natoli, C. R. (2000). Comput. Phys. Commun. 132, 251–266. Web of Science CrossRef CAS Google Scholar
Hara, S. (1967). J. Phys. Soc. Jpn, 22, 710–718. CrossRef CAS Web of Science Google Scholar
Hedin, L. (1965b). Arkiv Fysik, 30, 231–258. CAS Google Scholar
Hedin, L. (1965a). Phys. Rev. 139, A796–A823. CrossRef Web of Science Google Scholar
Hedin, L. & Lundqvist, B. I. (1971). J. Phys. C, 4, 2064–2083. CrossRef Web of Science Google Scholar
Hedin, L. & Lundqvist, S. (1969). Solid State Physics, edited by F. Seitz, D. Turnbull & H. Ehrenreich, Vol. 23, pp. 1–181. New York: Academic Press. Google Scholar
Joly, Y. (2001). Phys. Rev. B, 63, 125120-1–125120-10. Google Scholar
Joly, Y., Cabaret, D., Renevier, H. & Natoli, C. R. (1999). Phys. Rev. Lett. 82, 2398–2401. Web of Science CrossRef CAS Google Scholar
Kizler, P. (1993). Phys. Lett. A, 172, 66–76. CrossRef Web of Science Google Scholar
Landau, L. D. & Lifshitz, E. M. (1966). Mécanique Quantique. Théorie non Relativiste, p. 632 and ff. Moscow: Éditions MIR. Google Scholar
Lee, P. A. & Beni, G. (1977). Phys. Rev. B, 15, 2862–2883. CrossRef CAS Web of Science Google Scholar
Loeffen, P. W. & Pettifer, R. F. (1996). Phys. Rev. Lett. 76, 636–639. CrossRef PubMed CAS Web of Science Google Scholar
Löwdin, P. (1956). Adv. Phys. 5, 3–172. Google Scholar
Mattheiss, L. F. (1964). Phys. Rev. 133, A1399–A1403. CrossRef Google Scholar
Michalowicz, A. & Vlaic, G. (1998). J. Synchrotron Rad. 5, 1317–1320. Web of Science CrossRef CAS IUCr Journals Google Scholar
Müller, J. E., Jepsen, O. & Wilkins, J. W. (1984). Solid State Commun. 42, 4331–4343. Google Scholar
Müller, J. E. & Wilkins, J. W. (1982). Phys. Rev. B, 29, 365–367. Google Scholar
Mustre de Leon, J., Rehr, J. J., Zabinsky, S. I. & Albers, R. C. (1991). Phys. Rev. B, 44, 4146–4156. CrossRef CAS Google Scholar
Natoli, C. R., Benfatto M., Brouder, C., Ruiz Lopez, M. F. & Foulis, D. L. (1990). Phys. Rev. B, 42, 1944–1968. Google Scholar
Natoli, C. R., Benfatto, M. & Doniach, S. (1986). Phys. Rev. A, 34, 4682–4694. CrossRef PubMed Web of Science Google Scholar
Norman, J. G. (1974). Mol. Phys. 81, 1191–1198. Google Scholar
Penn, D. R. (1987). Phys. Rev. B, 35, 482–486. CrossRef CAS Web of Science Google Scholar
Poiarkova, A. V. & Rehr, J. J. (1999). Phys. Rev. B, 59, 948–957. Web of Science CrossRef CAS Google Scholar
Quinn, J. J. & Ferrell, R. A. (1958). Phys. Rev. 112, 812–827. CrossRef Web of Science Google Scholar
Rehr, J. J. & Albers, R. C. (2000). Rev. Mod. Phys. 72, 621–654. Web of Science CrossRef CAS Google Scholar
Schwarz, K. (1972). Phys. Rev. B, 5, 2466–2468. CrossRef Web of Science Google Scholar
Sham, L. J. & Kohn, W. (1966). Phys. Rev. 145, 561–567. CrossRef CAS Web of Science Google Scholar
Slater, J. C. (1979). The Self-Consistent Field for Molecules and Solids; Quantum Theory of Molecules and Solids. New York: McGraw-Hill. Google Scholar
Stern, E. A. (1993). Phys. Rev B, 48, 9825–9827. CrossRef CAS Web of Science Google Scholar
Tyson, T. A., Hodgson, K. O., Natoli, C. R. & Benfatto, M. (1992). Phys. Rev. B, 46, 5997–6019. CrossRef CAS Web of Science Google Scholar
Von Barth, U. & Grossman, G. (1982). Phys. Rev. B, 25, 5150–5179. CrossRef CAS Web of Science Google Scholar
Wille, L. T., Durham, P. J. & Sterne, P. A. (1986). J. Phys. (Paris), 47, C8-43–45. CrossRef Google Scholar
Wu, Z., Benfatto, M. & Natoli, C. R. (1996). Phys. Rev. B, 54, 13409–13412. CrossRef CAS Web of Science 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 OF
SYNCHROTRON
RADIATION

ISSN: 1600-5775

Volume 10| Part 1| January 2003| Pages 26-42

doi:10.1107/S0909049502017247

Format		BIBTeX
		EndNote
		RefMan
		Refer
		Medline
		CIF
		SGML
		Plain Text
		Text

Format		BIBTeX
		EndNote
		RefMan
		Refer
		Medline
		CIF
		SGML
		Plain Text
		Text

Search IUCr Journals		doi		Advanced search
Author		volume	page

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

X-ray absorption spectroscopy: state-of-the-art analysis

1. Introduction

2. Theoretical background

2.1. Reduction of the many-body problem to an effective one-particle problem: the photoemission case

2.2. The photoabsorption cross section

2.3. The mean free path

2.4. The muffin-tin approximation

3. Future improvements and conclusions

3.1. Inclusion of vibrational effects in XAS spectra

3.2. Elimination of the muffin-tin (MT) approximation in the solution of the Schrödinger–Dyson equation

3.3. Implementation of the multichannel MST (MMST) in the photoemission process

3.4. Core hole dynamic screening and the problem of the optical potential in the final state

APPENDIX A

Derivation of the MS equations

APPENDIX B

The Hedin–Lundqvist potential

APPENDIX C

Non-local core polarization potential

References

research papers