2p x-ray absorption spectroscopy of 3d transition metal systems

This review provides an overview of the different methods and computer codes that are used to interpret 2p x-ray absorption spectra of 3d transition metal ions. We first introduce the basic parameters and give an overview of the methods used. We start with the semi-empirical multiplet codes and compare the different codes that are available. A special chapter is devoted to the user friendly interfaces that have been written on the basis of these codes. Next we discuss the first principle codes based on band structure, including a chapter on Density Func- tional theory based approaches. We also give an overview of the first-principle multiplet codes that start from a cluster calculation and we discuss the wavefunction based methods, including multi-reference methods. We end the review with a discussion of the link between theory and experiment and discuss the open issues in the spectral analysis.


Introduction
This review provides an overview of the various routes to calculate the 2p x-ray absorption (XAS) spectral shape of 3d transition metal ions in solids and in molecular complexes. We focus on a general description of the various methods. For readers interested in the detailed mathematical and numerical methods used we refer to the original references.
The excitation of the 2p core electron is also known under the name L 2,3 edge. The core-hole makes XAS an element-specific technique as well as a local probe that derives detailed information without the necessity of long-range order. The 2p XAS spectra of 3d transition metals are positioned in the soft x-ray range from 400 eV to 1000 eV. This energy range allows in-situ measurements, though they are more challenging than for the hard x-rays.
If a transition metal system is exposed to X-rays, fine structures appear at energies related to the core binding energies. The lowest possible transition relates to an energy from the ground state to the lowest empty state, or more precisely to the lowest electron-hole exciton. The L 2,3 edge relates to a 2p core-level electron being excited to an empty state where the edge is dominated by transitions to the empty 3d states. The electric field of the x-ray photon interacts with the core electron and their overlap determines the chance for an excitation. This is captured by the Fermi golden rule that is based on the interaction of electromagnetic radiation with matter as described in relativistic quantum mechanics. The Golden Rule states that the x-ray absorption intensity (I XAS ) between a system in its initial state Φ i and a final state Φ f is proportional to: The delta function takes care of the energy conservation and a transition takes place if the energy of the final state equals the energy of the initial state plus the X-ray energy (ℏω). The intensity is given by the matrix element of the dipole operator T 1 between the initial and final states. The transition operator T 1 describes one-photon transitions such as X-ray absorption. For 2p XAS, the quadrupole transitions can be neglected and T 1 can be written as ε⋅r, where ε is the polarization vector.
The parity selection rule states that the final state Φ f must have different parity as the initial state, implying that from a 2p core state only s and dsymmetry final states can be reached. The initial state (Φ i ) and final state (Φ f ) wave functions are not precisely known and in practical calculations one must make approximations to calculate the x-ray absorption cross-section.
In this review we will limit ourselves to 2p x-ray absorption spectra of 3d transition metal ions. The methods discussed here are also applicable to other systems such as 4d, 5d, 4f and 5f systems, and to experimental methods such as X-ray Photoemission spectroscopy (XPS), Resonant Inelastic X-ray Scattering (RIXS) and X-ray Magnetic Circular Dichroism (XMCD) and related techniques.

Overview of the methods
We start with an introduction of the multiplet models that are central for the interpretation of 2p XAS of 3d transition metal ions. Section 2 introduces the different multiplet models, section 3 discusses the multiplet codes and section 4 the user-friendly graphical interfaces. First principle (or ab initio) multiplet calculations computes the 2p XAS spectra based on first principle electronic structure calculations. We note that the semi-empirical multiplet codes just mentioned can be turned into first principle codes by calculating all necessary parameters accordingly. For example, in the case of an octahedral transition metal ion, the ligand field can be calculated from first principles and combined with the first principle atomic parameters. This type of first principle approach contains several choices and approximations, for example treating the 3d3d interactions spherically and adding an effective field to treat non-spherical perturbations. But as we will see, also most first principle multiplet calculations have to make choices and approximations in their treatment of the 2p XAS spectra. In other words, in practical applications, there is not always a large difference between the semi-empirical and first principle multiplet approaches.
We divide the first principle approaches into four groups: • Section 5 on band structure multiplet codes that start from a band structure calculation in reciprocal space and via a projection to localized states perform a multiplet calculation on a small cluster. • Section 6 on band structure calculations based on the Bethe-Salpeter equation (BSE) and Time-Dependent DFT (TDDFT) approaches, including also the multi-channel multiple scattering method. • Section 7 on first principle cluster multiplet codes for solids that use molecular orbitals in a real space cluster embedded by point charges. • Section 8 on multireference (MR) methods for molecules that use restricted calculations to be able to perform first principle calculations.

History
In 1945 Rule gave an analysis of the 3d XAS (M 4,5 edge) of samarium. He attributed the pre-edge structure to localised 3d to 4f transitions, followed by an edge [1]. In 1966, Williams realized that the whole edge could be attributed to 3d-to-4f transitions, that is as transitions from a 4f N ground state to a 3d 9 4f N+1 final state [2]. This idea of "atomic transitions in solids" was theoretically developed for rare earth 4d XAS [3]. The usage of the "atomic transitions in solids" for 3d transition metal edges started from the experimental 3p XAS data from DESY [4] and especially from data collected on 3d metal halides [5]. Nakai wrote "the detailed structures in the low energy structures may be due to the 3p 5 3d N+1 multiplet, including the crystalline field splitting". This is essentially, the crystal field multiplet model, which then was used theoretically to explain the 3p XAS spectra of transition metal systems [6]. The charge transfer model was then used to explain the screening effects in x-ray photoemission [7]. The combination of the crystal field model and screening effects evolved into the charge transfer multiplet model, which has been developed by Thole and coworkers [Thole et al., 1988;de Groot et al., 1990; [8,9]; Okada et al. 1992]. Cluster multiplet calculations were also performed by the Groningen group [10], Hiroshima group [11] and from the Saclay group [12].

The physical origin of the multiplet effects
In the final state of 2p XAS there is a strong interaction between the 2p core-hole and 3d electrons. In a 3d transition metal system, for example NiO, the final state will have an incompletely filled 3d-band. The initial state of the nickel sites in NiO can be approximated with a 3d 8 configuration, yielding a 2p 5 3d 9 final state. The 2p-hole and the 3dhole have radial wave functions that overlap significantly. This wave function overlap is an atomic effect. It creates final states that are found after the vector coupling of the 2p and 3d wave functions. In order to understand the Coulomb interaction between the core hole and valence electrons one can make a multipole expansion of the interaction. It is found that the monopole part of the core-hole potential (Q) is largely screened to a value of 6− 8 eV, while the higher order multipoles of the interactions are largely unscreened in the solid state. The eigenstates generated due to these higher-order multipole interactions are known as the atomic multiplet effects and they are of the same order of magnitude in atoms as in solids and molecules.
In the case of the Ni 2+ 2p 5 3d 9 final-state configuration, the 2p spinorbit coupling is approximately 17 eV, yielding an energy splitting of ~20 eV between the L 3 and L 2 edges. For multiplet effects to have a significant effect on the mixing of the L 3 and L 2 edges, the multiplet splitting due to Coulomb interactions need to be of the same order of magnitude as the spin-orbit coupling splitting. We can quantify the energy scale of the multiplet splitting by looking at the parameters arising from the multipole expansion of the Coulomb interaction. These parameters are known as the Slater-Condon parameters or integrals. In order to find substantial mixing between the spin-orbit split L 3 and L 2 edges the value of the Slater-Condon parameters must be of the same order of magnitude as the core-hole spin-orbit coupling separating the two edges. The quadrupole part of the Coulomb interaction, parameterized by the F 2 Slater integral is approximately 8 eV for Ni 2+ , large enough to mix the L 3 and L 2 edges and certainly large enough to modify the energies of the 3d 8 configuration. Whether a multiplet effect will actually be visible in x-ray absorption further depends on the lifetime broadening due to the core-hole decay and the amount of charge fluctuations of the valence shell. The core-hole lifetime leads to a line-width which approximately has a full width at half maximum value of 0.4 for the L 3 edge, implying that both the 2p core spin-orbit coupling and the 2p3d multiplets are larger than the lifetime broadening. Strong charge fluctuations can lead to an overlap of many different local configurations and a superposition of many different multiplet structures. This can lead to broad structures where the underlaying multiplets are not trivially visible. The interaction between the 2p and 3d orbitals is larger than between the 2p and 4d or 5d orbitals. At the same time the 2p core-hole lifetime of 4d and 5d transition metals is shorter and the amount of charge fluctuations often larger than for 3d transition metals. This leads to broader peaks, thereby hiding the multiplet structures.

The crystal field multiplet model
An effective method to analyze the L 2,3 edges of 3d metal systems is based on the crystal field multiplet model, where in this review we will use the terms ligand field (in coordination complexes) and crystal field (in solids) as equivalent. The crystal field multiplet model describes the atomic interactions, where the surroundings are treated as a perturbation using an effective electric field. The crystal field multiplet model is justified because the 2p3d transition creates 2p 5 3d N+1 self-screened excitonic states.
The crystal field multiplet model is an effective model Hamiltonian for the description of all charge conserving excitations of ionic transition metal systems. In a transition metal ion described with a 3d N configuration, the most important interaction is the 3d3d electronelectron interaction. The 3d3d interactions determine the ground state symmetry, which is best described using a basis of L and S quantum numbers, expressed in short as 2S+1 L. For example, the ionic 3d 8 configuration has L = 3 and S = 1, i.e. a 3 F ground state, as determined by the Hund's rules, i.e. (1) maximum S and (2) maximum L. The other four term symbols are labeled 1 S, 3 P, 1 D and 1 G. Their excitation energies are respectively 1 D at 1.7 eV, 3 P at 2.0 eV, 1 G at 2.6 eV and 1 S at 6.4 eV. The 1 S state is positioned at a higher energy position because it relates to two holes located in the same orbital. The 3d spin-orbit coupling is small with a value of less than 100 meV, but this energy is large enough to be effective for temperatures at and below room temperature. Following the third Hund's rule, maximum J (if the shell is more than half-filled), the 3d spin-orbit coupling leads to a 3 F 4 ground state with total-angular momentum J = 4 for a Ni 2+ ion.
In a molecular complex or a solid, the atomic states are modified by the point group of the transition metal ion. In other words, based on the self-screened excitons, the transition metal ion is approximated as an isolated 3d N ion surrounded by a distribution of charges that mimic the system. This crystal field multiplet model is able to explain a large range of experiments, including x-ray absorption, x-ray emission and optical transitions.
We will use the 3d 8 configuration as an example to show the effects of crystal field theory. The above mentioned five atomic states of 3d 8 in spherical symmetry split into eleven crystal field states that are further split by the small 3d spin-orbit coupling. The Tanabe-Sugano diagram captures the changes in energy of the electronic states as a function of the crystal field strength. Consequently, the 3 F ground state splits into 3 T 1 , 3 T 2 and 3 A 2 states. The 3 A 2 state identified by all t 2g -states plus the spin-up e g -states being occupied, is the ground state. Fig. 1 shows the Tanabe-Sugano diagram for a 3d 8 configuration. The 3 F atomic ground state is split into 3 A 2 , 3 T 1 and 3 T 2 states and in octahedral symmetry 3 A 2 is the ground state. The 2p XAS spectrum is calculated as all transitions from the 3 A 2 state in Fig. 1 at a crystal field of 1 eV, with the resulting spectrum shown in Fig. 2. The interactions that create this spectrum are (a) the 3d-3d electron interactions, (b) the 3d spin-orbit coupling, (c) the octahedral crystal field, and for the final state also the 2p3d electron-electron interactions and the 2p spin-orbit coupling.

The charge transfer multiplet model
The approximation of a transition metal ion as 3d N neglects all other electrons and interatomic interactions. The crystal field multiplet model can be improved by adding more configurations. This includes ligandmetal charge transfer, where an electron is transferred from the ligand to the metal indicated with 3d N+1 L, where L denotes a hole on the ligands. The charge transfer energy is indicated with Δ and relates to a band, where in the most drastic approximation this band is simplified to a single state. Other charge transfer channels include the two-metal interaction 3d N 3d N ↔ 3d N+1 3d N− 1 related to the Hubbard U parameter. In the charge transfer multiplet model, one effectively combines the ground state configuration with other low-lying configurations.
Self-screened excitons: From the charge transfer multiplet model one can explain why the crystal field multiplet model is an effective model Hamiltonian for 2p XAS. Fig. 3 shows the 3d 8 and 3d 9 L configuration of a Ni 2+ system, split by the charge transfer parameter Δ (middle). In 2p XAS the final states are respectively 2p 5 3d 9 and 2p 5 3d 10 L, split by an energy Δ F =Δ+U-Q, where Q is the core-hole potential. Since U is roughly equal to Q, the final state charge transfer energy Δ F is equivalent to Δ, implying that the bonding combination of 3d 8 and 3d 9 L has similar charge transfer effects for both the initial and final state. The consequence is that the great majority of the intensity originates from the bonding initial state to the bonding final state. In other words, charge transfer satellites are weak and the transitions in 2p XAS are self-screened excitons. This reasoning is valid for all neutral spectroscopic transitions, but not for ionizing experiments such as 2p XPS. In 2p XPS the electron is emitted from the local atom, which has the consequence that the ordering of states drastically changes due to the core-hole potential. This implies that 2p XPS is dominated by charge transfer satellites. The advantage of the charge transfer multiplet theory is that one can also describe both XAS and XPS from a unified description [13].
The energy configurations of 3d 8 +3d 9 L is indicated in Fig. 4. On the right the limit where the charge transfer energy is set to 10.0 eV. The crystal field multiplets of 3d 8 are found between 0 and 4 eV and the multiplets of 3d 9 L are found between 12 and 14 eV. Decreasing Δ the 3d 9 L states move towards lower energy. Due to their interaction with the 3d 8 states, they effectively compress the 3d 8 multiplet structure, which is known as the nephelauxetic effect. At negative values of Δ, the 3d 9 L states start to dominate the ground state. An ion with large charge transfer effects is Cu 3+ and its ground state can be described as a linear combination of 3d 8 and 3d 9 L with a Δ value of approximately -1.0 eV [14]. In this case, it is important to perform the full charge transfer model as discussed for rare-earth nickelates [Bisogni et al., 2016].

Important interactions in semi-empirical calculations
Below we will describe a series of semi-empirical and first principle calculations and in order to compare their treatment of the 2p XAS problem, we first separate the interactions into different categories, based on the local crystal field model Hamiltonian as a starting point. The interactions in the charge transfer multiplet model can be separated into (I) local electron-electron interactions, (II) spin-orbit interactions, Fig. 1. The Tanabe-Sugano diagram for a 3d 8 ground state, as a function of the octahedral crystal field parameter 10Dq. The atomic 3 F ground state (red) is split by 3d spin-orbit coupling into 3 F 4 , 3 F 3 and 3 F 2 . The octahedral crystal field splits the 3 F ground state into 3 A 2 , 3 T 2 and 3 T 1 states. The excited state atomic multiplets 1 D, 3 P and 1 G are given in blue. The 1 S state has a high energy. The spectrum is calculated with atomic parameters and with a crystal field parameter 10Dq = 1 eV.  (III) crystal field and molecular field interactions, (IV) charge transfer interactions and (V) the broadening.
(I) Local electron-electron interactions: These include (a) the core-hole binding energy and (b) the higher order terms of the interactions between the 3d valence electrons and (c) the multiplet effects, in other words the higher order terms of the interaction of the 2p corehole with the 3d valence electrons. The absolute XAS excitation energy is usually not calculated with an accuracy better than 1.0 eV. This implies that the theoretical spectral shape is shifted to align with the experiment, which in turn can suffer from a miscalibration of the x-ray energy. The core-hole valence-hole exchange interaction is in multiplet theory described by the G 1 pd Slater integral. G 1 pd is 4.6 eV for 3d 8 Ni 2+ and 1.5 eV for 4d 8 Pd 2+ , in other words the 2p3d Slater integrals are three times larger than the 2p4d Slater integrals. The higher order term of the core-hole valence-hole exchange interaction is described with the G 3 pd Slater integral. G 3 pd is 2.6 eV for 3d 8 Ni 2+ and 0.88 eV for 4d 8 Pd 2+ . The core-hole (2p) -valence-hole (3d) multipole interaction is described with the F 2 pd Slater integral. F 2 pd is 6.2 eV for 3d 8 Ni 2+ and 1.8 eV for 4d 8 Pd 2+ . It can be concluded that 4d systems have small multiplet effects and large 2p spin-orbit coupling, which brings the L 2,3 edge of 4d systems to the limit where multiplet effects can in first approximation be neglected, in contrast to the 3d systems.
(II) Spin-orbit interactions: The spin-orbit interactions include the 3d valence spin-orbit coupling and the 2p core spin-orbit coupling. The 2p core spin-orbit coupling creates two separated features in the x-ray absorption spectrum, historically running under the names L 3 and L 2 edges for 2p core-holes. Without the multiplet effects and the inclusion of the ground state spin-orbit coupling, the L 3 and L 2 spectra correspond to a 2p 3/2 core-hole and a 2p 1/2 core-hole respectively, which would imply that they have an integral ratio of 2:1. The 2p core-hole coupling is 11.5 eV for 3d 8 Ni 2+ and 107 eV for 4d 8 Pd 2+ .
(III) Field interactions: The crystal field and the molecular exchange field are effective field interactions in the local model Hamiltonians that are supposed to capture the electron-electron interactions with the neighboring atoms in molecules and solids. In the case of highsymmetry (octahedral, tetrahedral) only one effective parameter is needed for the 3d shell, simplifying empirical approaches. Lowsymmetry systems needs a number of empirical parameters, in which case a combination with a first principle calculational approach is suggested for a reliable result.
(IV) Charge transfer interactions: Charge transfer interactions describe the interaction of the 3d N ionic model Hamiltonian with other configurations, in order to describe the ground state electronic structure more completely and also to capture the core-hole screening processes more efficiently. As such, charge transfer interactions bear much resemblance to configuration interaction. Charge transfer interactions can also be calculated with densities-of-states based on band structure calculations and in this way the translation symmetry (dispersion) can be included into local cluster calculations.
(V) Broadening parameters: Traditionally 2p XAS is broadened with a Lorentzian function in order to capture the lifetime of the 2p core excitons, plus a Gaussian function to describe the experimental broadening. The lifetime broadening of 2p core excitons in the L 3 edge is relatively constant and usually approximated with a value of 0.4 eV fullwidth half-maximum, where the L 2 edge states have an additional broadening due to the super-Coster-Kronig Auger decay. If the 2p excitons interact strongly with continuum states, a Fano line shape is used to capture this. One open question regarding 2p XAS is if the intrinsic (Lorentzian) broadening used is only due to lifetime effects, given the fact that calculations of the core-hole lifetime typically yield a broadening of 0.2 eV. More likely the intrinsic broadening contains a significant effect from vibrations and as such captures the combination of lifetime and vibrational broadening. This vibrational broadening concerns the final state excitonic vibrational broadening and not the effects of vibrations on the ground state.

The THOLE multiplet code
The THOLE multiplet program is a suite of self-contained codes that have been adapted and modified by Theo Thole and his many coworkers. The first usage of the THOLE multiplet code was in the calculation of the 3d XAS spectra of the 4f rare earths [15]. A multiplet calculation from 4f 6 to 3d 9 4f 7 has 48048 fin. l states, implying that the matrix diagonalizations become very large. These atomic calculations were performed using the Cowan code [16], which computes the reduced matrix elements of the atom in spherical symmetry. This needs only a few effective parameters: (a) the Slater integrals to describe the core-valence and valence-valence electron-electron interactions in the atomic potential and (b) the core and valence spin-orbit interaction [17,18]. The THOLE multiplet code is comprised of several self-contained codes that are run in sequence. The RCN program calculates the initial and final-state wave functions in intermediate coupling using the atomic Hartree-Fock method with relativistic corrections, together with the Coulomb and exchange integrals (F k and G k ) and spin-orbit parameters. After empirical scaling of the Slater integrals, Cowan's RCG program calculates the electric dipole and quadrupole transition matrix elements from the initial state to the final-state levels of the specified configurations [17]. Since the Cowan code uses spherical wave functions and Wigner-Racah tensor operators, which specify the magnetic quantum number, the Wigner-Eckart theorem can be used to obtain the spectra for linear and circular polarization in the presence of a Zeeman field [19].
The next program of the THOLE code is based on Butler's point group program [20]. Starting from the atomic reduced matrix elements, all the required reduced matrix elements can be obtained in any of the 32 different point groups, using coefficients calculated from group theory. The advantage of the Butler method, with respect to older methods [21] is that the coupling coefficients (isoscalar factors) are fully consistent over all point groups. The program can calculate the transition probabilities between any two configurations, such as the 2p XAS spectra of 3d transition metal ions. This version of the THOLE codes effectively treats crystal field theory, and over the last 35 years many calculations have been published using this approach [15,22,23,8,9]. The program BANDER, developed by Thole in collaboration with Kotani and coworkers, extends the crystal field multiplet calculations with the hybridization between different configurations. The number of configurations is only limited by the memory size and computing speed. This enables the study of the interplay between atomic multiplet structure and solid-state hybridization using the Anderson impurity model. The BANDER program exists in a number of varieties based on either exact diagonalization or the Lanczos method. The Lanczos method speeds up a large matrix calculation, by effectively describing the full spectrum by a specific limited set of levels. This speeds up the calculation, but if the set is chosen too small it has the disadvantage that specific final states are lost, which can lead to inaccuracies for second-order processes such as RIXS.
Temperature dependent spectra are obtained taking a Boltzmann distribution over the calculated energy levels of the initial state. Furthermore, using the initial-state wave functions the expectation values of any ground state operator can be calculated, such as the spin and orbital moments and magnetic dipole term. The x-ray polarization vector and magnetic field can be chosen along arbitrary directions with respect to the crystalline field orientation. The Cowan code also efficiently calculates the transitions to continuum states, which enables to compute x-ray photoemission, resonant photoemission, Auger spectroscopy, constant initial state spectroscopy, and Bremsstrahlung isochromat spectroscopy (BIS). This provides a natural way to obtain the Fano line shape of the 2p XAS in 3d transition metals by calculating the coherent resonance between the 2p XAS mediated resonant photoemission 2p 6 3d N → 2p 5 3d N+1 → 2p 6 3d N− 1 ε process and the direct photoemission 2p 6 3d N → 2p 6 3d N− 1 ε, where ε represents a continuum state.

The TANAKA multiplet code
The TANAKA multiplet code was developed by Arata Tanaka in the early 90's. [11]. The program calculates transition metal clusters in the Anderson impurity model, including clusters with arbitrary number of orbitals and sites. It was one of the first full multiplet programs which extensively used the Lanczos and the method based on Krylov projection to calculate many-body Green's function of the final and intermediate states. The program can calculate all spectroscopies that involve a core-hole creation or decay, including XAS and RIXS.

The XCLAIM multiplet code
The XCLAIM code makes use of an uncoupled product basis of atomic or real orbitals. In fact, the use of a product basis allows one to go one step further and make a code where the starting point is the Hamiltonian and not the basis. The code only needs to know how different "shells" are coupled via the one-and two-particle interactions in the Hamiltonian without any knowledge of atomic orbitals. A shell is defined as a set of states that are grouped together. The shell can contain atomic orbitals of a particular angular momentum, a subset of orbitals, spins, etc. The total system is then a combination of those shells. Physically, these shells can be on the same site or different sites, although the code only needs to know the coupling within and between the shells via the matrices for the one-and two-particle interactions. Such a second-quantization based approach was initially incorporated in conventional programming languages, such as Fortran (as used in the underlying code for XCLAIM) and C, and later extended to script-based languages, such as Mathematica. In addition, to treating atomic spectroscopy, this approach allows one to build arbitrary clusters with the major restriction being the size of the many-body Hamiltonian. The Hamiltonian in matrix form is solved using Lanczos and tridiagonalization methods. The code itself is run by an external script that allows one to perform different types of calculations without changing the code. Finally, the code also calculates expectation values of importance for comparing X-ray spectroscopy with the sum rules for integrated intensities [Thole & van der Laan, 1988; [24,25].
One of the first advantages of this approach was the ability to study X-ray spectroscopy on systems containing more than one transitionmetal ion by making use of the flexibility in constructing systems and the possibility of using a reduced basis. This allowed the study of nonlocal screening effects in 2p core-level photoemission absorption site [26], the effects of doping [27], and the influence of exchange interactions between neighboring sites [28]. The graphical user interface for XCLAIM, discussed in section 4.2, provides a user-friendly way to create the input and script file for the underlying Fortran code and subsequently broadens and displays the output.

The QUANTY multiplet code
This section reviews the empirical multiplet calculations using QUANTY [29][30][31]. Parameter free calculations using QUANTY are described in section 4.1.
Philosophy: QUANTY implements a script language that can solve general problems in quantum chemistry and physics, focusing on spectroscopy and dynamics of correlated electron systems. The language allows one to define operators in second quantization and calculate Eigenstates or (generalized n-point) Green's functions for these operators. In principle one is completely free to define any operator and calculate any response function, given the memory and time constraints set. The idea behind the script language is that the user can focus on the physics whereas the numerical considerations are dealt with internally with only one parameter to control convergence, namely the required accuracy of the calculation. The program gives the user a great deal of freedom to choose the physical model which comes at the price of the need to do some programming.
For the case of L 2,3 edge core-level x-ray absorption spectroscopy on transition metal compounds there are several scripts available that only need minor modification to fit the situation at hand. Depending on the size of the basis set used one can change the level of the theoretical model used. We will describe here the specific case of 3d transition metal L 2,3 edges, but modification to other edges can be made straightforwardly. If one includes the 2p core orbitals as well as the 3d valence orbitals (i.e. a basis of 16 orbitals in total), and includes the effects of the solid as an effective potential on the 3d orbitals, one works on a level of crystal field multiplet theory. Extending the basis to include also 10 ligand orbitals and explicit covalence allows one to do charge transfer multiplet calculations. In cases where there is an interest to study π backbonding one can add a second ligand shell to account for the hybridization of the correlated 3d shell with both the valance and conduction orbitals of the solid. For the simulation of metals it can be necessary to include a full band, which can be realized by including not one, but several ligand orbitals, which represent a discretized version of the continuum states. Interactions between transition metal 3d shells can be included by creating a double cluster [32]. The only limit on the calculation is memory and computation time. The size of the Hilbert space is given by the binomial coefficient of the number of electrons present in the problem and the number of orbitals included. Five electrons in a 3d-shell lead to a Hilbert space of 252 states in crystal-field multiplet theory, 15,504 states in charge-transfer multiplet theory and 847,660,528 for a double cluster with two ligand shells. Due to several optimizations and the fact that QUANTY does not store operators as matrices nor wave-functions as vectors, but stores them as true operators and functions that allow one to calculate expectation values and new wave-functions after an operator acted on a function, the memory needed is often much smaller. These optimizations can be made even more effective by reducing the required accuracy to less than standard machine precision.
In order to stay in line with previous model calculations, the configuration included in the calculation can be restricted. For a ligand field theory calculation one can limit the occupations to 3d N and 3d N+1 L. This is very useful to reproduce earlier results as well as to decrease computation time. For accurate calculations it is advisable to include at least configurations with two or three ligand holes. For larger basis sets this allows one to implement several forms of restricted active space calculations and compare the accuracy of each of these approximations.
Parameters: For empirical calculations, the parameters needed are the size of the Coulomb interaction, spin-orbit coupling strength, possible exchange or magnetic fields and, depending on the level of theory, either the crystal-fields as potentials acting on the 3d-shell (crystal field multiplet theory) or the hopping strengths between the 3d and ligand orbitals, as well as potentials acting on the 3d-shell and ligand shell (charge transfer multiplet theory). It is common to take the multipole part of the Coulomb interaction from scaled atomic calculations, which works quite well as the multipole part of the Coulomb interaction is hardly screened in a solid. The monopole part is strongly screened and thus is treated as an empirical parameter.
Implementation: The start of any QUANTY calculation is the definition of the basis orbitals used. Although the program in principle does not need this, it gives a physical meaning to the different Fermion indices or quantum numbers. For a ligand field calculation one could define that index 0-5 refers to the 2p shell, index 6-15 to the 3d shell and index 16-25 to the ligand shell, that the even (odd) indices refer to states with spin down (up) and that the atomic shells are given on a basis of spherical harmonics with the angular momentum quantized in the z direction ordered from -l to l. The input in QUANTY would be: With this basis set one can define several operators. For the Coulomb interaction within the 3d shell one can use one of the many standard operators defined: Besides predefined operators for standard models one can create operators using creation and annihilation strings as well as add or multiply different operators.
Once the required operators are defined one can calculate the lowest N eigenstates of these operators (Hamiltonian) using the function Eigensystem. Note that operators are unaware of the number of electrons in the system. The most straightforward call to the function Eigensystem will thus calculate the lowest N eigenstates within a grand canonical potential, i.e. it determines the occupation with the lowest energy at the same time as it determines the eigenstates. For most models considered the chemical potential is not set and one needs to fix the occupation. This can be done with the use of restrictions. There are two different ways to use restrictions possible. The first way sets a restriction on the starting point of the iterative procedure that evolves to the ground-state of the operator whose eigensystem is calculated. As the iterative procedure conserves the particle number (symmetries included in the Hamiltonian) one can use the starting position to restrict the possible outcomes. The second way to use restriction is used throughout the calculation and allows one to restrict the number of configurations included or do a restricted active space calculation. Depending on the size of the Hilbert space and problem at hand different numerical routines will be used to find the eigensystem. The methods implemented include several dense methods, restarted Lanczos routines as well as a Block version of the Lanczos routine. During the calculation the Hilbert space will be enlarged iteratively until the required accuracy is reached. The aim of the code is to fully automatically choose the best numerical method and needed Hilbert space to deliver the result with the required accuracy. Although this is implemented and functional in several cases, there still is ongoing work in progress to automate the choice of best numerical method and gain full error control. The automatic mode will be useful for most users, at the same time expert users do have full control on which method is used with the use of options.
Once the ground-state wave-function is calculated one can calculate the spectra: which calculates the function where |Ψ 〉 is a starting eigenstate (often the ground state), T is a transition operator (for example a 2p3d dipole excitation for L 2,3 XAS), H is the Hamiltonian operator, ϖ is the energy relative to the energy of state |Ψ 〉, and Γ represents the core-hole lifetime. The variable G contains the spectrum object and can be modified printed or saved to disk. Additional options can be set defining for example the numerical mesh used for the energy grid. The algorithm calculates T|Ψ〉 and uses this state as a starting vector for the calculation of a Krylov basis of H. Within this Krylov basis the inverse of the operator is given as a continued fraction.
A full manual of all functions included and standard operators defined as well as the different options each function takes can be found at www. quanty.org. Output: QUANTY is a script language based on the programming language LUA. Variables containing a spectrum can be saved to file in Ascii format listing the real and imaginary part of the spectrum. The imaginary part is the absorption spectrum, whereas the real part enters in the calculation of resonant x-ray diffraction. One can also use the scripting properties of QUANTY / LUA to define additional functions which in combination with for example GNUPLOT can create immediate plots as output. The latter has the advantage of being very transparent and guarantees complete reproducibility, especially if the experimental data is read by the script and plotted with the same file. The most considerable advantage of having a script language is that one can write small scripts that automate part of the calculation or loop through parameter space.

Additional multiplet codes
In addition to the THOLE, TANAKA, XCLAIM and QUANTY multiplet codes a number of additional codes have been written. Crocombette wrote a cluster-based charge transfer multiplet code, using the Slater integrals and spin-orbit coupling plus the Slater-Koster hopping, Hubbard U, charge transfer as parameters. No symmetry restrictions have been used. The CROCOMBETTE multiplet code has been used to calculate the 2p XAS of 3d transition metal oxides such as TiO 2 and Fe 2 O 3 , with equivalent results to the THOLE calculations [12].
Stepanow et al. have written a charge transfer multiplet code that is based on Cowan, but that is not limited by symmetry. The STEPANOW multiplet code is a single-atom symmetry-free multiplet code that is mainly used to calculate XMCD spectra and spin-and orbital expectation values related to 2p XAS spectra of 3d systems and 3d XAS spectra of rare earths [33]. Multiplet codes that are based on ab-initio methods are discussed in section 5.
Krüger wrote a charge transfer multiplet code applied to 2p and 3 s XPS spectra [34,35] and a general "symmetry-free" CF multiplet code for XAS [36].

The CRISPY interface
CRISPY is a graphical user interface for the simulation of core-level spectra using the semi-empirical multiplet approaches implemented in QUANTY, shown in   the PYTHON ecosystem. The user interface is developed using the PyQt5 library which provides a comprehensive set of python bindings for the Qt5 framework. Spectra manipulations are done using a dedicated plotting widget from the SILX library. Computationally demanding tasks, such as interactive broadening of the spectra, are done using the FFT module from the NUMPY library.
Using CRISPY it is possible to perform XAS, XMCD and RIXS calculations for transition metals and lanthanides considering various level of approximations for the semi-empirical Hamiltonian and different site symmetries of the absorbing atom. In the case of non-centrosymmetric site symmetries, the inclusion in the Hamiltonian of elaborate pd mixing terms is also possible. In addition, for the simulation of L 2,3 edges in transition metals energy depended broadenings and angular dependence can be considered. More information about CRISPY, including detailed installation instructions and tutorials on how to use it, can be found at http://esrf.fr/computing/scientific/crispy. To make the interface more accessible for non-expert users, easy to use package installers for Windows and macOS operating systems are also provided.

The XCLAIM interface
XCLAIM [37] is a graphical interface written in PYTHON that interfaces with a compiled Fortran code as described in section 3.3 and calculates spectra using a multiplet model that includes a general crystal-field based on Wybourne parameters. The XCLAIM multiplet code can be downloaded at: https://subversion.xray.aps.anl.gov/xclaim/xclaim. html. It also possible to include charge-transfer within TMO 6 and TMO 4 clusters. The program can calculate XAS, RIXS, XPS, valence photoemission and inverse photoemission. Fig. 6 shows the main window with the controls for setting the model parameters. In addition to calculating spectroscopy, it is also possible to get information from the ground and excited states by calculating expectation values of quantum operators, like the spin and orbital momenta, and occupation of the d-orbitals that give information on crystal field excitations [38].
The expectation values are also useful to estimate the errors made when using the sum rules. Several factors can limit the application of the sum rules to experimental spectra: the magnetic dipolar term [25], the mixing of spectral weight between the L 2 and L 3 edges occurring in transition metals [12,39] or the fact that the isotropic intensity is approximated as the average of left and right circularly polarized absorption [40]. For cases that are not implemented in the graphical interface, it is possible to use the PYTHON modules of XCLAIM directly from the PYTHON interpreter, together with the Numpy, Scipy and Matplotlib libraries for plotting. The code saves the Hamiltonian matrices and transition operators and includes sample python scripts to calculate XAS and RIXS that can be edited and modified by the user. This is useful to repeat the spectra calculations with arbitrary polarizations, or change the population distribution of the initial state, for looking at temperature dependence of the spectra or to use the initial intermediate states that are reached in pump-probe experiments.

The CTM4XAS interface
The THOLE multiplet program has been used in the CTM4XAS interface. CTM4XAS can calculate the XAS, XPS, XES and RIXS spectra of transition metal systems and rare earths [41]. The XAS calculations include the 1 s, 2p and 3p core excitations of the 3d, 4d and 5d transition metal ions. In addition the 3d and 4d XAS spectra of uranium and of the rare earths can be calculated. The XPS calculations include all possible core excitations, where the interface limits the number of configurations to two. XPS calculations with more configurations must be calculated with the THOLE program. CTM4XAS also can be used to calculate the XES spectra and the matrix elements for 1 s2p, 1 s3p and 2p3d RIXS calculations, that can be performed with CTM4RIXS.
The CTM4XAS6 software includes an option to switch from the THOLE to the QUANTY program [29]. CTM4XAS6 is programmed by Mario Delgado. It can be used to calculate the 2p and 3p XAS spectra of 3d transition metal ions. The THOLE and QUANTY codes can both be used, and as far as has been tested provide the same result. The code can also be used to generate the LUA scripts for QUANTY calculations, that then can be further modified as required. CTM4XAS can be found at: http://www.anorg.chem.uu.nl/ CTM4XAS/ An additional program from Mario Delgado is CTM4DOC [42]. This program analyses the ground state that has been used to calculate the XAS spectra. The ground state can be analyzed in terms of (a) an orbital description, (b) an atomic term symbol description or (c) a crystal field term symbol description, where the respective components that make up the total ground state are given.
In addition the Differential-Orbital-Covalency (DOC) is given for charge transfer calculations and charge-transfer Tanabe-Sugano diagrams can be generated as a function of all parameters used. DOC calculations represent a method to project the ground state from multiplet simulations fitting experimental data into a molecular orbital picture, which can be compared to DFT calculations to assess theoretical models. This approach is particularly useful to better understand molecular systems in chemical catalysis and bioinorganic chemistry.

The MISSING interface
MISSING (Multiplet Inner-Shell Spectroscopy Interactive GUI) is a userfriendly interface to the THOLE program. The MISSING stays close to the original THOLE program and can be used as alternative to CTM4XAS. MISSING is written by Riccardo Gusmeroli and can be found at http://www .esrf.eu/computing/scientific/MISSING/

First principle band structure multiplet calculations
Two similar band structure multiplet codes are being developed, respectively QUANTY written by Maurits Haverkort (sections 5.1 and 5.2) and the LDA + DMFT program by Atsushi Hariki (section 5.3).

The QUANTY first principle multiplet code
In section 2.4 we described how to use QUANTY to calculate core-level spectroscopy using empirical models. Models including crystal field multiplet theory, charge transfer multiplet theory, multiple site clusters or Anderson impurity models have been used for many years with success and QUANTY allows for an implementation of each of these. For high symmetry crystals, i.e. those where the transition metal ion is on a site with cubic point-group symmetry, the number of parameters needed in these calculations is moderate, especially if the Coulomb interaction is taken from atomic Hartree-Fock or Density Functional Theory calculations. For low symmetry systems determining all crystal field or ligand field parameters becomes cumbersome.
The bonding within a crystal is a result of a gain of kinetic energy and a loss of potential energy with respect to an assembly of neutral atoms. In order to capture these energies one needs a minimal basis that includes covalence, i.e. explicit hopping between orbitals to represent the gain in kinetic energy. Within crystal field theory the kinetic energy gain is represented by an effective potential. Although efficient for the description of low energy states, we have not been able to device a reliable ab initio method creating crystal-field Hamiltonians. On a level of charge transfer multiplet theory where the gain in kinetic energy is explicitly included by hopping to the ligand orbitals leading to covalence we do find vary satisfactory ab initio descriptions of many materials.
Following methods devised by Gunnarson et al. [43,44] we implemented a method using Wannier functions to determine the crystal field parameters, i.e. the orbital dependent onsite energy, the hopping parameters and the orbital dependent Coulomb interaction [29]. The basic procedure is to start from a density-functional theory or Hartree-Fock calculation of the full periodic solid. We then reduce the basis set by picking a small set of Wannier orbitals for bands around the Fermi energy. For ligand field theory these include the transition metal 3d derived bands as well as the ligand (oxygen) 2p derived bands. Note that these orbitals are close to, but not exactly equal to the atomic orbitals. In most materials, independent of the amount of covalence, the transition metal 3d Wannier orbitals can be made closely resembling the atomic orbitals. The ligand Wannier orbitals must contain at least tails of transition metal 4 s character representing the transition metal 4 s, ligand 2p bonding states. Note that the choice here is not to get maximally localized Wannier functions, but to find a set of Wannier functions where the transition metal 3d states are close to the atomic orbitals and the ligand orbitals are chosen such to represent the occupied bands exactly.
Once the band-structure is calculated and a Wannier representation of the bands is found that describes the band eigen-states well and contains the transition metal 3d Wannier orbitals that are atomic like, one can create the ligand field Hamiltonian. The transition from a tight binding Hamiltonian to a ligand field cluster representation is done using a block tri-diagonalization procedure starting from the five transition metal 3d orbitals. This results in a unitary transformation of the ligand orbitals from Ligand atom centered Wannier states to transition metal centered ligand orbitals. This procedure not only leads to a substantial decrease in computational cost, it also leads to an ab initio Hamiltonian that is directly comparable to the model Hamiltonians described above. Any band-structure code that can create a Wannier representation of a subset of the bands can be used with this procedure, including the Stuttgart LMTO code [45], Wien2k [46] and FPLO [47]. Interfaces to other codes are planned.
Besides the calculation of the one-particle interactions, i.e. the hopping and onsite energy, one needs to include the Coulomb interaction in the calculation. This can be done by evaluating the Coulomb integrals between the different Wannier orbitals in the basis. As the transition metal Wannier orbitals are near atomic-like one can express within a very good approximation these integrals in terms of Slater integrals. As our basis sets are small, the direct calculation of the Coulomb interaction is expected to result in much larger values than needed in the empirical methods discussed in section 2.4. In order to remedy this one needs to calculate the screened interaction between the low energy effective orbitals due to all states not included in the basis. Such calculations are done using for example the random phase approximation [48].
We used a different method based on the fact that screening in a solid is very effective for spherical interactions, but that multipole interactions are much less screened. For the Coulomb interaction, this means that the monopole part of the Coulomb interaction, i.e. that part that depends on the total occupation of the d shell, but not in which orbital or spin state the atom is, needs to be screened. The multipole part, i.e. the part that describes the difference in the Coulomb interaction between a 3d xy orbital and a 3d x2-y2 orbital or a 3d z2 orbital is much less screened. The latter determines the multiplet splitting in the spectra and is for many cases much more decisive in predicting the experiment than the monopole part of the Coulomb interaction. Using the ab initio values for the multipole part of the Coulomb interaction (F 2 and F 4 ), but empirical values for the monopole part of the Coulomb interaction we find good agreement between theory and experiment, with one side condition. As the size of the monopole part of the Coulomb repulsion shifts the onsite energy of the transition metal d orbitals with respect to the ligand orbitals, one needs to treat also the energy difference between these orbitals as a parameter (the charge transfer energy, Δ). In practice we set the average energy of the d-shell (ε d ) and the average energy of the ligand shell (ε p ) to zero and add the charge transfer energy (Δ) to the many-body Hamiltonian. The resulting two parameters (much like those in an LDA + U or DMFT calculation) are well determined using corelevel photo-electron spectroscopy.
The calculation of spectra follows from here on the same procedure as is needed for the calculation of spectra from the empirical models. We generally advise that empirical model calculations start with an ab initio determination of the parameters and if deviations between theory and experiment are found parameter optimization can be made. This can be done, for example, by choosing different ground-state symmetries in rare earth compounds. These deviations between theory and experiment either show the numerical accuracy of the theory or the fact that some parts in the model used are not sufficient to represent the material. Fig. 7 illustrates several examples of ab-initio QUANTY calculations compared to experiment. For the L 2,3 XAS of high symmetry oxides SrTiO 3 , MnO, and NiO shown in Fig. 8(a-c), respectively, the agreement between the ab initio theory and experiment is excellent. In the lower symmetry case of TiOCl shown in Fig. 8(e), the agreement is still very good and the detailed polarization dependence of the XAS arising from the low symmetry is captured as well. Fig. 8(d) shows a corresponding 2p XPS spectrum.

The RSPt interface to QUANTY
QUANTY can also be used in conjunction with the full-potential linear muffin-tin orbital (FP-LMTO) code RSPt (http://fplmto-rspt.org). Here the impurity problem is formulated in the language of dynamical mean field theory (DMFT), which operates with several sets of effective ligand orbitals [49].
The first step is to construct a set of localized wavefunctions describing the 3d impurity states. There is no unique way of obtaining these states and a possible strategy is to perform Wannierization of the band structure by selecting the bands belonging to a certain energy range and possessing particular atomic character [29]. RSPt adopts another procedure and defines more general projection schemes, which are independent of the band structure details and do not require any a priori knowledge about the system. There are two projection schemes [50] muffin-tin-heads and Löwdin-orthogonalized LMTO's. The former projection takes the density within the muffin-tin sphere of a certain l-character and neglects the interstitial part. The latter one mixes the orbitals centered on different sites, making it reminiscent to Wannier functions. Once the projection is defined, one can calculate the local Hamiltonian. The projected Hamiltonian will contain all information about the crystal field and spin-orbit coupling directly from the DMFT calculation.
An impurity problem with several sets of bath orbitals is then efficiently solved using QUANTY or an alternative impurity solver, which also uses Krylov projection scheme [https://github.com/JohanSchott/impu rityModel]. It is worth mentioning that for a long time the number of bath states (N b ) was restricted to just a few orbitals, but recent development allowed to consider even hundreds of bath states. A relatively moderate number of N b is sufficient in order to converge the 2p XAS spectra in transition metal oxides. However, large N b would be vital e.g. for organometallic compounds where the hybridisation is spread among many different orbitals. As mentioned above, the developed approach takes advantage of a versatile projection scheme, which does not involve the fitting of the band structure, which can be cumbersome for large systems. Thanks to this feature, one can study, for example, the impact of the actual presence of the core-hole on the Slater integrals and other calculated parameters. The results for monoxides suggest that the transition metal 3d-orbitals become slightly less hybridizing and less crystal-field-split [49].

The LDA + DMFT first principle multiplet code
The LDA + DMFT first-principle multiplet code for calculating corelevel spectra builds on the Anderson impurity model that incorporates the valence states described by the local-density approximation (LDA) + dynamical mean-field theory (DMFT) [51][52][53][54]. LDA + DMFT is an ab-initio many-body method and has been successful for describing electronic, magnetic, and structural properties in materials with strongly correlated electrons, such as 3d transition-metal oxides [55][56][57]. Our computational scheme for calculating core-level spectra consists in the post-processing of the LDA + DMFT self-consistent calculation. For 2p excitations in 3d transition-metal oxides, the Anderson impurity model with the DMFT hybridization function V(ε) is extended to include explicitly the 2p core orbitals and their interaction with transition-metal 3d orbitals. Our method can be viewed as an extension of the cluster model (composed of a transition-metal ion and nearest-neighboring ligands [58]) to include the effect of the long-distant atoms, while retaining the single-impurity description. This can be carried out by replacing discrete ligand states of the cluster model by a continuous spectrum which represents atoms (host) surrounding a transition-metal ion. While this description is exact for a noninteracting host, the host in 3d transition-metal compounds is strongly correlated due to the electron-electron interaction between localized 3d electrons. The LDA + DMFT scheme provides an optimized effective description of the host in that case [51], [Ghiasi et al., 2019], where the interaction in the host is represented by the local self-energy Σ(ε) calculated by the Anderson impurity model and determined self-consistently in the valence states of the entire lattice. The advantages of the LDA + DMFT approach over the cluster model are following: (1) accurate description of the charge transfer effect from not only the neighboring ligand states but also from/to the correlated 3d bands. The latter charge-transfer channel, traditional called nonlocal charge transfer [26], is missing in the cluster model description.
(2) elimination of the ambiguities in the choice of the parameters (e. g., amplitude of the metal-ligand hybridization and crystal-field splitting), which are replaced by the (almost) parameter-free LDA + DMFT calculation.
Our code implements a configuration-interaction (CI)-based exact diagonalization method. It allows us to include the full multiplet interaction and a large number of the bath levels (about 25 bath levels per spin and orbital) for describing the band features in the DMFT hybridization function V(ε) in 2p XAS, 2p XPS and 2p-3d RIXS calculations for 3d transition-metal compounds. A careful consideration on the CI basis expansion may be necessary in materials with a large charge fluctuation due to e.g., a strong covalent bonding with ligand states or with a mixed valence. Practically, the inputs to the code are (a) discretized bath levels, (b) a double-counting correction (which relates to the charge transfer energy), and Coulomb interaction parameters for (c) valence-valence and (d) valence-core channels. The inputs (a)-(c) are fixed in the DMFT self-consistent calculation, independent of the core-level spectral calculations. The monopole part of the valence-core interaction needs to be fixed by fitting to the experimental core-level spectrum. In standard LDA + DMFT implementations, Wannier functions constructed from the LDA bands (usually including transition-metal 3d bands and oxygen 2p bands for 3d transition metal oxides) are used as a basis set. As the Wannier functions of the 3d orbitals are atomic-like and well localized, the Slater integrals estimated by the atomic Hartree-Fock calculation could be used for describing the multipole part of the valence-core interaction. A more consistent procedure with the Wannier functions for evaluating the Slater integrals is also proposed [29]. Fig. 8 shows 2p XPS and 2p XAS, respectively, calculated by the LDA + DMFT approach for NiO. In the XPS spectrum, besides the global structure, such as the spin-orbit splitting and the charge transfer satellite, the double-peak structure in the 2p 3/2 main line, reported in hard xray photoemission experiments [59], is well reproduced. The comparison to the spectrum calculated by the cluster model (dashed line) emphasizes the importance of non-local charge transfer effect in the formation of the double peaks. The nonlocal charge transfer is a common contribution factor to the 2p XPS spectra of transition metal oxides and must be taken into account when interpreting the spectra [51]. In the 2p XAS spectra, we find a reasonable agreement to the experimental data in both LDA + DMFT and the cluster model. This is expected since the charge transfer effect in the final states is weaker in 2p XAS compared to that in 2p XPS, especially in Mott insulators with a large charge gap [53]. However, the nonlocal charge transfer effect in 2p XAS becomes important in high-valence transition metal oxides with negative (or very small) charge-transfer energy. The LDA + DMFT approach reproduces 2p XAS fine features due to the nonlocal charge transfer in LaCuO 3 and LuNiO 3 [52].

Density Functional Theory based approaches
A key aspect of the various methods discussed thus far is the onerous scaling of computation time as the Hilbert space (or active space) is increased. Therefore, it has been important to find effective ways to reduce this space, such as the band down-folding procedure used in the DFT-based multiplet calculations. However, this approach remains difficult to realize when a system has many overlapping bands, such as for metals, or when both localized and delocalized excitations are important. In this chapter we highlight several approaches that are able to treat extended systems, such as the pure-element 2p XAS of 3d transition metals [60,61]. In this class of approach, the many thousands of electron states preclude the use of the many-body techniques highlighted in earlier chapters. Instead various simplifications to the interaction must be employed.
We note that in case of closed shell systems, one can describe the 2p XAS spectra accurately with DFT codes directly. For example using the XSPECTRA routine within QuantumEspresso, an accurate description can be found for Cu metal and Cu 2 O that both have their 3d-band filled, implying that the large 2p3d and 3d3d interactions are not necessary to include for the 2p XAS calculations [62]. This approach based on a DFT calculation of the ground state, followed by a projected augmented wave method is also accurate for the calculation of all condensed matter systems without an open shell 3d, 4f or 5f shell and in addition for the oxygen K edges of transition metal oxides [63][64][65].

Linear response time-dependent DFT
The extension of the linear response TDDFT formalism to treat core electron excitations has been developed first at non-relativistic level [66] and then at relativistic levels, which can be scalar relativistic or spin-orbit [67]. A reliable description of the 2p XAS spectrum requires the inclusion of the configuration mixing among different excitation channels due to the degeneracy of the 2p core-hole [68]. Furthermore the inclusion of spin-orbit coupling is of mandatory importance to properly describe both the L 2 and L 3 edges. In the scalar-relativistic scheme all the relativistic contributions are kept except the spin-orbit coupling. For closed-shell electronic structures, singlet-singlet and singlet-triplet excitations can be calculated, while for open shells a spin polarized scheme is possible but no more a spin eigenfunction. The space spanned by the solutions of the eigenvalue equation of the TDDFT Casida formulation [69] corresponds to the one-hole-one-particle excited configurations, so it is possible to approximate this space by selecting a subset of configurations, keeping only those necessary for an accurate description of the phenomenon. To describe core electron excitations, the indexes which span the occupied spinors can be limited to run only over the core shell under study. Therefore, the original Casida Ω matrix equation is solved for a (much smaller) submatrix of the original matrix; the Davidson algorithm can be used because core excitations correspond to the lowest states. TDDFT calculations of 2p XAS spectra have been performed for titanium compounds [66,67,70] and the series of the closed-shell ionic compounds VOCl 3 , CrO 2 Cl 2 and MnO 3 Cl [71]. Open-shell systems have been calculated for Mn [72] and Fe [73]. A systematic study on the bulk metals from Ca to Ni has been published [Bunȃu et al. 2012]. Nevertheless, due to the wrong asymptotic behavior of current exchange-correlation kernels for solid systems, strongly bound excitons, which occur e.g. in the 2p XAS spectra of transition metal oxides are not properly described [74]. This fact hampers the applicability of current TDDFT approaches to core spectroscopy in solids.
The Ti 2p XAS spectrum of TiCl 4 molecule is dominated by transitions from the Ti 2p core orbitals to the valence virtual orbitals with mainly metal 3d character. In the T d symmetry of TiCl 4 these orbitals are split in the 3e and 10t 2 manifolds by the crystal field. The TDDFT calculations are able to describe the strong mixing of configurations responsible for the intensity redistribution between the two Ti 2p valence transitions, with the t 2 transition 20 times more intense than the e transition [66]. The major limitation of this level of theory is the Fig. 9. Ti 2p XAS spectrum of TiCl 4 from TDDFT calculations with the LB94 exchange correlation potentials. Upper panel: scalar relativistic TDDFT results; lower panel spin-orbit TDDFT results [67]. Experimental spectrum from [66].
inability to treat the spin-orbit effects which prevents to distinguish between the L 3 and L 2 edges. The inclusion of the spin-orbit coupling in the calculation allows to describe both thresholds, shown in the lower panel of Fig. 9 [67].

Bethe-salpeter equation
The Bethe-Salpeter equation (BSE) approach of many-body perturbation theory is a general theory that describes the bound states of a twoparticle system [75]. BSE treats radiation-induced excitations using an effective interaction, where the computation of the absorption spectrum is performed in three steps: First, the electronic structure is calculated within DFT. While DFT codes are routinely employed to determine the electronic structure around the band gap, for 2p XAS, the electronic structure and wavefunctions of the deep-lying 2p states have to be calculated as well. On top of these DFT calculations, quasiparticle energies and wavefunctions can be obtained to account for electron-electron correlation. In practice, though, the electronic structure is corrected by many-body perturbation theory calculations in the GW scheme typically in the conduction region, while no correction is applied to the 2p states. As a third step, the absorption spectrum is described by the electron-hole Green's function. In the BSE formalism, this Green's function is determined by a Dyson-like equation, which determines the interacting electron-hole Green's function to the independent-particle ones through the four-point interaction kernel. In common approximations to the interaction kernel, this equation can be mapped to an effective two-particle problem of an excited electron and a core-hole. In this effective problem, the electron and the core-hole interact through the exchange interaction and the direct interaction, which represents the attraction of the oppositely charged particles.
Inclusion of these Coulomb interactions accounts for multiplet effects, rigorously in the 3d • case, but only approximately for open shells. BSE spectra tend to agree well with measurements of early 3d transition metal oxides, but show discrepancies for oxides toward the middle of the 3d series. The BSE for both core-and valence-level excitations is implemented within the OCEAN code [76,77]. Due to the localized nature of the core-hole, computational effort in evaluating the screening response may be limited to the neighborhood of the absorbing site by working in real-space [78]. Finally, rather than solving explicitly for the final states by full diagonalization, the spectrum is generated by Lanczos iteration. An all-electron version of the BSE formalism [79,80] is implemented within the all-electron full-potential code package EXCITING [81], which has recently been applied to calculations of core-spectroscopy [82]. Kohn-Sham energies and wavefunctions are obtained from all-electron full-potential DFT calculations employing the linearized plane-wave basis set. The non-local screening is fully calculated using currently a plane-wave expansion for the non-local operators. The effective BSE eigenvalue problem is solved by explicit diagonalization, which not only yields the absorption spectrum but also the excitonic eigenstates. As an all-electron code, EXCITING calculations for both core and valence spectra can be obtained on an equal footing. Both EXCITING and OCEAN are open-source and freely available, see www. exciting-code.org and www.ocean-code.com, respectively. Bethe-Salpeter calculations for core-level excitations can be remarkably efficient. Unlike core-hole constrained DFT approaches used to calculate absorption spectra for K-edges, the use of ground-state Green's functions within the BSE obviates the need for an explicit core-hole or supercell.
In Fig. 10 [84], due to the localized nature of both the initial 2p and final 3d states, the exchange interaction in the BSE Hamiltonian requires a careful treatment of local field effects. Within EXCITING these effects, which originate from the local variations of the Coulomb potential, are described by the |G + q| max parameter. The essential role of local field effect for an accurate description of the Ca 2p XAS can be seen in Fig. 10, where the calculated spectrum is displayed for increasing values of the |G + q| max parameter. Without local field effects (|G + q| max = 0), the spectrum is broad and featureless, in disagreement with the experimental one. At intermediate values for |G + q| max , the correct spectral shape is obtained, as the screened Coulomb interaction between electron and core-hole is captured correctly. The exchange interaction on the other hand, requires a more careful description. As the exchange interaction mixes the transition of the sub-edges, the relative peak heights are only captured at high values for |G + q| max . When the exchange interaction is converged,  Our study shows that the BSE formalism is capable of an accurate first principle description of the localized excitations in the 2p XAS spectrum of CaO, where it reproduces both the spectral shape as well as the off-statistical branching ratio. For these localized excitations, an accurate treatment of the local variations of the potentials is key, as the exchange interaction significantly influences the relative peak intensities.
As a further example, we present the iron 2p XAS of FeSe and FeTe in Fig. 11. Disagreement remains on the strength and nature of electronic correlations in the superconducting iron-chalcogenides and ironpnictides. Comparisons of spectra calculated for different ground-state electronic structures can provide additional constraints on theory. To test the appropriateness of band theory, we calculated the BSE spectra based on ground-state DFT electronic structures within the local density approximation to the exchange-correlation functional. DFT Kohn-Sham eigenvalues were corrected with a band-dependent, k-point averaged complex-valued GW self-energy prior to solving the BSE. Agreement with experiment would likely improve if full, state-by-state self-energy corrections were applied, though remaining discrepancies could indicate the importance of additional, neglected local correlation effects.
Because BSE and various other calculations build on an underlying DFT calculation, certain aspects of chemical realism can be included in exchange for limitations to the complexity of correlation effects accessible in other treatments. Ligand-field or crystal-field splittings emerge in BSE spectra by the separation of peaks, whereas the underpinning source for this splitting is encrypted in the band structure and crystal field that is not so obvious for spatially extended Bloch states. Spin-orbit and atomic multipolar Coulomb interactions are naturally included in a BSE treatment or a multi-channel multiple-scattering treatment. First principles DFT calculations can provide parameters suitable for building in Franck-Condon or Jahn-Teller effects at post-processing stages of a BSE calculation [85,86]. Additional post-processing effects can include incorporation of the cumulant representation of the Green's function [87], which is another means by which charge-transfer satellites are included in an approximate fashion, but without adjustable parameters. Finally, subtle effects can arise such as the polarization dependence of the position of a multiplet peak in rutile TiO 2 [88], which BSE calculations seem to report correctly, and other subtle changes that take place when transition between two phases of SrTiO 3 [89].

Multi-channel multiple scattering
Multiple scattering theory is widely used for x-ray absorption fine structure calculations in the single-particle approximation. Multichannel multiple scattering (MCMS) was devised as a generalization of the multiple scattering method to correlated many-electron wave functions [90]. In multichannel scattering, the N-electron wave function on a given site is written in close-coupling expansion, i.e. including configuration interaction with a continuum wave. For a x-ray absorption final state the wave function on the absorber atom is expanded as an eigenstate of the N-1-electron ion. The photo-electron wave components undergo multiple scattering with the environment and thereby probe the electronic structure of the entire system. At the absorber site, the components are coupled by scattering in and out the different channels corresponding to different local many-electron states such as multiplet levels. MCMS thus combines the single-particle electronic structure of the extended systems with local electron correlation effects. [91] have implemented MCMS with a particle-hole wave function and the variational R-matrix method for the computation of the multichannel T-matrices [92]. In this particle-hole scheme, the multiplet coupling between the excited electron and the core-hole is well accounted for, but configuration mixing with other electrons is disregarded. As a consequence, the MCMS method, similar to the Bethe-Salpeter equation approach [93], yields best result for systems with an empty shell ground state, e.g. compounds with nominal (3d • ) metal ions. MCMS has been successfully applied to 2p XAS spectra of various early transition metal compounds, including free molecules, metals and the ionic minerals CaO, CaF 2 , Ca(CO) 3 , TiO 2 , SrTiO 3 and V 2 O 5 [91,[94][95][96][97]. Fig. 12 shows the Ti 2p XAS spectra of titanium dioxide in its main polymorphs rutile (a) and anatase (b) [95]. While the independent particle approximation (IPA, Fig. 1a,b) completely fails, the MCMS calculation (c,d) reproduces well all spectral features, including the asymmetric doublet (D-E) whose origin had been debated for two decades. The MCMS results show that the D-E splitting is a non-local effect related to the arrangement of the TiO 6 octahedra in the crystal rather than to their (weak) distortion as suggested previously.

The DFT-CI multiplet code
The DFT-CI multiplet code has been applied to various transition metal oxides with different numbers of 3d-electrons, spin-states, and symmetries [98][99][100][101]. This method is a wavefunction-based approach on top of the configuration interaction theory in quantum chemistry. A DFT-CI calculation is performed in real-space using a cluster model composed of a single transition metal ion and neighboring ligand ions. The cluster is embedded into the array of point charges placed at external atomic sites of clusters so as to mimic the solid-state effects. For a given cluster model, relativistic DFT calculations are performed. All the relativistic effects including the 2p spin-orbit coupling are taken into account by solving the Dirac equation, where each molecular orbital is a four-component spinor. In the present DFT-CI code, a numerical basis set is used to expand the molecular orbitals, where the radial distributions of the basis atomic spinors are computed numerically. The radial distribution functions for Ti-2p 1/2 and Ti-2p 3/2 for SrTiO 3 and the molecular orbital energies of valence O-2p and Ti-3d levels are shown in Fig. 13, where P 2pj (r) and Q 2pj (r) (j = 1/2, 3/2) represent the radial distribution of upper and lower two components of the Ti-2p atomic spinors, respectively. The difference of spatial distribution of 2p 1/2 and 2p 3/2 core orbital results in the difference of Coulomb and exchange interactions between 2p and 3d orbitals, and hence has an influence on spectral shapes of 2p XAS. Once the molecular orbitals are obtained, the one-electron integrals and two-electron integrals are evaluated over all possible combinations of molecular orbitals, which are required to construct the many-electron Hamiltonian. The latter corresponds to the Slater integral in the THOLE multiplet calculation. The values of the two-electron integrals over the molecular orbitals are usually smaller than the Slater integrals for atomic orbitals, because metal 3d orbitals are spatially spread due to the formation of covalent bonds to ligand orbitals [Ikeno et al., 2011].
A many-electron wavefunction in the configuration interaction is expressed as a linear combination of Slater determinants composed of the relativistic orbitals. Due to the limitation of computational power, we usually restrict the number of Slater determinant. The condition for 2p XAS corresponds to the electronic transition from the 3d N configuration to the 2p 5 3d N+1 configuration. The charge transfer effects can be taken into account by including the additional electronic configurations. However, the computational cost will be drastically increased as the number of Slater determinants increases exponentially as increasing the number of orbitals considered in the configuration interaction. The many-electron Hamiltonian matrix is fully diagonalized to obtain all possible initial and final states under the restricted electronic configurations. The 2p XAS spectra are obtained in accordance with the Fermi's golden rule, where the transition matrix elements are computed using the wavefunctions for the initial and final states. Alternatively, transition metal 2p XAS spectra are evaluated solely from the initial state as the imaginary part of the Green's function formalism. The latter scheme is suitable for handling systems with large Hamiltonian matrices. Several iterative algorithms, such as the generalized Davidson, Lanczos, and reduced-shifted conjugate-gradient method with seed switching, have been implemented in the DFT-CI method [102,103].
This configuration interaction method is known to overestimate the absolute transition energies systematically. This can be ascribed to the incompleteness of our basis orbitals and the limitation on the number of Slater determinants. In the present DFT-CI method, the additional molecular orbital calculation is made for the Slater's transition state, where a partial core-hole (0.5 in most cases) is introduced. The transition energy is corrected by the considering the energy difference between orbitals for the Slater's transition state as a reference. Fig. 14 shows the Ti 2p XAS of SrTiO 3 obtained by the DFT-CI calculations without charge transfer (b) and with charge transfer (c), which are compared with the experimental spectrum (a) [Uehara. et al., 1997]. Widely spreading multiplet levels appear when the charge transfer is  included especially at the L 2 edge. As a result of this, the L 2 peaks are much more broadened.

The MULTIX multiplet code
The MULTIX multiplet method [104] is a self-contained, accessible ligand-field/crystal-field approach to calculate XAS, RIXS and inelastic neutron spectra (INS), as well as ground state properties including anisotropic magnetic properties of the ion. The MULTIX website http://multiplets.web.psi.ch/ provides downloads and documentation for MULTIX and auxiliary codes for converting crystal structures to ligand cluster files for the application of the crystal-field point-charge model. MULTIX works in the active space of the core-hole shell and the local orbital valence shell as a self-contained first principles code. This aspect calculates the key matrix elements for electron-electron interaction, relativistic effects including spin-orbit coupling and single-electron orbitals with minimal input. The atom label and the specification of the core-hole shell and the partially occupied valence shell is enough input to define the multi-determinant active space based on Dirac relativistic atomic orbitals from the built in density functional atomic differential equation solver. This alone produces rough XAS spectra with reasonable excitation energy and core spin-orbit splitting. States with the proper orientation are conveniently generated by specifying a quasi-crystal field as a point charge model of the molecular or crystal environment. The charges are considered as semi-empirical parameters to define the effective ligand field matching a corresponding experimental spectrum. An accurate parametrization of the electric crystal field should not be expected from a point charge model with conventional valence charges. And even an accurate crystal field would miss the effects of hybridization on the single-electron resonances, which crucially enter the ligand field-model as single-electron levels. A small number of neighbor shells with their charge parameters is sufficient to model all possible splittings of the single-electron levels for the localized open valence shell. Such semi-empirical parameters can be found by optimizing fits for one or several experimental spectra. The charge parameters relevant to the local orbital levels can be alternatively inferred from band structure calculations or from cluster calculations of electronic structure.
The connection of the Hamiltonian with the Cartesian coordinates of the crystal environment is useful to specify polarization directions. MULTIX handles the full Hamiltonian in the active space without recourse to symmetry, so any low symmetry environment can be treated with equal ease. The code automatically switches into Lanczos type spectral calculations, when this method is expected to be faster. The Lanczos method also allows to do the calculations requiring significantly less computer memory. The finite active space of the open localized shells leads to a multiplet line spectrum. For a useful comparison with measured spectra it is unavoidable to include the effects of an intractable number of further states in a semi-empirical fashion. The coupling into continuum states limiting the core-hole lifetime is modeled by Lorentzian core-hole lifetime broadening. The coupling with states responsible for screening is modeled semi-empirically by scaling the first principles parameters for electron-electron interaction, S-O splitting and by introducing a threshold correction for the X-ray absorption edge. As emphasized above, the coupling of the local valence orbitals with the bands of ligand states in a solid tends to become intractable in a seriously first principles manner. Thus, MULTIX allows to adjust the open-shell atom single-electron resonances by considering point charges as semiempirical parameters.
The case of SrTiO 3 may serve again for illustration of semi-empirical spectrum fitting with the MULTIX approach. Because of the cubic environment of the Ti atom, it is sufficient to scale the crystal field splitting as part of the spectral fitting. The point charge model of the crystal field was truncated after the second neighbour shell, to show a second subset of symmetry equivalent atoms for illustration. In fact, the same cubic splitting can be introduced to the Hamiltonian by using only the first neighbor shell, with a different crystal field scaling parameter value. The threshold correction is another single-electron parameter, in this case less than 2% correction for the value from the approximate first principles part of the calculation. In the case of spin-orbit coupling, a renormalization by less than 10 % is introduced by its semi-empirical parameter, which is quite a typical value. The scaler Coulomb parameter provides a scaling to the first principles electron-electron interaction matrix elements. A renormalization by less than 20 % is needed in the case of electron-electron interaction. This reduction of the electronelectron interaction by increased screening in a compound as compared to the free atom has been termed the nephelauxetic effect in an inorganic chemistry context [105]. For the 2p 6 3d • ground-state, electron-electron interaction only acts as core-valence interaction for the excited state in the active space. This electron-electron scaling factor is in a value range fitting many compounds. A detailed study of core-hole decay is not usually done for core-hole spectroscopy interpretation and is beyond the scope of the MULTIX approach. A Lorentzian core-hole broadening is specified as a semi-empirical parameter. The excitation energy dependence of the core-hole is an important effect for experimental comparisons. Fig. 15 shows the resulting fit of the experimental spectrum by Uehara (red) with the MULTIX model (green). The model is complemented by a Gaussian resolution profile of 0.25 eV width for the X-ray spectroscope and a background model based on convolution of the theoretical spectrum with a Heaviside function above an excitation gap of 4.5 eV. The case of SrTiO 3 illustrates that the MULTIX crystal/ligand field multiplet model can represent peak positions often with very good accuracy and intensities with fair accuracy in a complex XAS spectrum with modest semi-empirical adjustments. The small remaining discrepancies in intensities highlight mainly shortcomings of the core-hole ramp model. An early realization that atomic-like multiplets show up prominently in some types of core-hole spectroscopies came with the investigation of Moser et al. [106].
Model calculations predating the MULTIX development focus on the interplay of the localized (multiplet) levels with the continuum of band states based on the Anderson impurity model. As summary is in [107]. Multifaceted use of MULTIX is illustrated in a number of articles: An article predating the MULTIX method publication compares the excitations of the closed shell Ti 4+ ion as calculated with TDDFT to the MULTIX multiplet calculation. For this case of a non-degenerate ground state, scalar-relativistic TDDFT generates the LS coupled multiplet with qualitative correctness and in similarly good agreement with experiment as the LS limit for MULTIX [108]. A study on iron pnictides shows among other findings, that the 2p3d XAS and also RIXS is represented well by the crystal/ligand field multiplet [109]. Combined XAS and RIXS studies  [110]. Magnetic ground-states were studied by circularly polarized XAS [111] and with XMCD on M 4,5 edges of rare earth elements [112,113]. In the adsorption case, it was not possible to match experimental findings with point charge parameters on atomic sites. It was necessary to involve charges, with proper site symmetry, in the graphene p-cloud above the nuclear plane to semi-empirically model the ligand-field interaction with the f-electron levels. MULTIX has also been applied to inelastic neutron scattering on a rare earth element [114].

The LFDFT-ADF multiplet code
The Ligand Feld Density Functional Theory (LFDFT) methodology has been recently implemented in the Amsterdam Density Functional (ADF) program packages for the calculation of multiplet energies and ligand field effects of open-shell electron configurations. ADF2019 is available at htttp://www.scm.com. It is a DFT-based model that aimed at modeling the electronic structure and XAS spectra of metal compounds, namely transition metal L 2,3 -edge XAS [115] and lanthanide M 4,5 -edge XAS [116] and XMCD. [117] LFDFT-ADF utilizes a conceptual DFT procedure by using the average of configuration concept to resolve the multiplet structure, i.e. a non-standard procedure that is possible in ADF by controlling the electron occupation of the reference Kohn-Sham orbitals. That is, the Kohn-Sham orbitals calculated with predominant 2p and 3d characters of the transition metal ion are populated with fractional electron charge. This insures a totally symmetric electron density, under which the effective ligand field Hamiltonian is invariant and can be operated. The effective ligand field Hamiltonian is used to treat explicitly the configuration interaction within the restricted subspace of the reference Kohn-Sham orbitals. The success of LFDFT-ADF resides in the proper construction of the reference Kohn-Sham orbitals, which must contain predominant transition metal 2p and 3d atomic orbitals. In LFDFT-ADF, the multiplet structure is non-empirically calculated, and reasonable agreement with experiments are obtained for transition metal 2p XAS model [115]. Fig. 16 shows selective Slater-Condon parameters of series of divalent transition metal ions within configuration 2p 5 3d N+1 , compared with available reference data taken from the literature [23]. The represented parameters include the Slater-Condon integrals F 2 pd , G 1 pd and G 3 pd and spin-orbit coupling constant, as obtained from LFDFT-ADF calculations using a GGA functional [115].
The LFDFT-ADF methodology also enables an accurate description of ligand field effects by taking into consideration a small cluster of molecules or crystal structures into the calculation. The small cluster is composed of the transition metal ion center together with its coordination sphere, which is often composed of electronegative ligating elements [115]. Fig. 13 shows the change of the radial functions of the Kohn-Sham orbitals with predominant 3d characters of Mn 2+ ion from the free ion to the molecular (MnX 6 ) 4− cluster with octahedral symmetry, with X = F-, Cl-and Br-, obtained from LFDFT-ADF calculations using a GGA functional. Note that the 3d orbitals of Mn 2+ ion split into three lower and two upper energy molecular orbitals, usually designated as t 2g and e g irreducible representations of the O h point group. Therefore, Fig. 13 shows averaged radial functions that are take into consideration the five components of the Kohn-Sham orbitals with predominant Mn 3d characters. That is in LFDFT-ADF, there is no m l dependency of the manifold of the atomic basis functions, i.e. the electrons are supposed to move in a central field.    [126]. The active space is divided into three subspaces: RAS1 including Fe 2p orbitals, RAS2 including valence σ3d bonding, n3d non-bonding, and σ*3d antibonding combinations, RAS3 containing π*(CO)+Fe 3d orbitals. Only one core-hole is allowed in RAS1 and only one electron can be excited to RAS3, whereas a full configuration interaction calculation is done within the RAS2 subspace.
It is seen in Fig. 17 that the Kohn-Sham orbitals with predominant Mn 2+ 3d characters are associated with radial functions of different shapes, which depends on the bonding regime and covalency between the Mn 2+ ion and ligands orbitals. More particularly, the variation of the shape of the radial functions illustrates the phenomenology of the nephelauxetic effect, which describes the expansion of the electroncloud toward the position of the ligands. Naturally, the parameters representing the Slater-Condon integrals and spin-orbit coupling constants are strongly influenced by the ligands, so that they are reduced if compared to the free ion (see Fig. 13 inset). Hence in LFDFT-ADF, scaling and reduction factors are not required. One can obtain theoretical ligand field parameters, that are reasonably good for the emulation of the multiplet energies and ligand field effects associated with the initial and final states of the x-ray absorption process in transition metal compounds, i.e. 3d N and 2p 5 3d N+1 .

Multireference methods based on restricted active space
Another group of ab initio methods, multi-reference wavefunction methods, has been developed with a focus on accurate descriptions of valence excited states and can be considered as a "golden standard" for this type of calculations [118]. However, this approach has been systematically applied to 2p XAS spectra of transition metals only relatively recently [119][120][121]; for recent reviews see [122,123]. In the multi-reference ansatz [124], the reference wavefunction of the many-body state is represented as a linear combination of, e.g. Slater determinants, which are generated in configuration interaction fashion. In contrast to the ordinary configuration interaction procedure, the energy of the state is minimized variationally not only with respect to the coefficients but also with respect to the set of molecular orbitals.
The practical application of MR calculations is demanding if not impossible for medium-sized systems. A solution is provided by the concept of the active space. At the most basic level it corresponds to the subset of molecular orbitals for which different levels of configuration interaction expansions are performed. In complete active space (CAS) calculations, electron configurations, arising from all possible excitations of electrons in the active orbitals, are included in the configuration interaction calculation. More suitable for core-level states is the restricted active space (RAS) approach, which allows restrictions in the number of excitations in the selected sub-sets of orbitals, e.g., the core orbitals. This makes it possible to target in particular highly excited states without having to calculate all lower-lying states. The first multireference wavefunction calculations of transition metal 2p XAS were therefore performed using the RAS formalism in the MOLCAS quantum chemistry program suite [119][120][121]125].
To study dipole allowed L 2,3 edge transitions one needs to include 2p and 3d orbitals into the active space as illustrated in Fig. 18. The active space is further subdivided into RAS1, RAS2, and RAS3 subspaces. In absorption spectroscopy, one is interested in singly excited core-states and thus the number of holes is limited to one in the RAS1 space comprising 2p orbitals. For highly correlated 3d electrons included in RAS2 it is usually necessary to perform CAS-type calculations. Additionally, one can include ligand orbitals, e.g., in RAS3 space as in Fig. 18, with the possibility to limit the number of electrons excited into this subspace [121,[126][127][128]. As the energy difference between core-hole and valence states is large, the number of configurations can be further reduced by invoking the core-valence separation technique [Delcey et al., 2019]. Active space methods imply a local correlation treatment whereby the active orbitals usually belong to a small fragment of the system. However, the choice of the active space can become quite involved. To reproduce 2p3d transitions the active space should comprise of 2p and all orbitals containing notable 3d character. For the simple cases of cationic [M(H 2 O) 6 ] n+ systems this can be three t 2g non-bonding and two e g anti-bonding orbitals (denoted in octahedral symmetry) [129,130]. In case of stronger covalent bonding, ligand orbitals increasingly mix with the metal ones. Filled orbitals can be added to describe e.g., shake-up and charge-transfer transitions as in [FeCl 6 ] 3− [121] and Mn n+ complexes [131], while empty ligand-dominated orbitals can be included to describe effects of backdonation and explicit charge transfer excitations as in Fe(CO) 5 [126,127] and [Fe(CN) 6 ] 3-/4- [121,132]. To investigate dynamics pathways after ligand-to-metal charge-transfer excitations in aqueous ferricyanide by determining the spectral signature of transient chemical species, additional t 1u ligand orbitals can be added to the active space [133].
Valence 3d→3d excited states can be included as well. This allows addressing RIXS spectra which, being integrated along the outgoing photon energy, can also be used to model XAS in the PFY mode [134,135]. In principle, any other occupied or unoccupied ligand or metal orbital can be added, including Rydberg orbitals, thus switching certain effects on and off. For instance, the 3 s orbital has been added to describe 3 s→2p XAS-PFY spectrum in addition to 3d→2p XAS-PFY in aqueous Fe (II) (Fig. 19) [134].
By virtue of the variational principle, in the RAS self-consistent field (SCF) calculations, optimal MOs and optimal configuration interaction coefficients are obtained at the same time [124]. With this, the one-electron MO basis consistently adapts to the description of the electron correlation effects. As orbitals can be variationally optimized in both initial and core-excited states it is possible to correlate spectra to electron and spin density changes in the x-ray process [136]. Note that this optimization procedure searches for orbitals that minimize the energy of the targeted states, not necessarily those that contribute most strongly to the X-ray spectra. Including the right number of final states can therefore be crucial to reproduce the XAS spectra [137].  6 ] 2+ [134].
The sophisticated treatment of electron correlation comes with a major drawback, namely the complexity of simultaneous optimization of a large number of MO and configuration interaction coefficients. The number of MO coefficients increases with the number of atoms and the number of basis functions per atom. RAS calculations are therefore ideally suited to describe medium-sized gas phase or solvated metal complexes. Periodic systems are not well represented, which is one of the major limitations of the method. The number of configuration interaction coefficients mainly depends on the number of orbitals in the active space. If more than one metal is present, the active space and number of electronic states rapidly grow, eventually hitting an "exponential wall" of computational cost. To cope with this situation, a Frenkel exciton model has been applied to x-ray spectra of the hemin dimer [Preusse et al., 2016] being a model for a general molecular aggregate. This extends the applicability of RAS to a weakly interacting multi-center system with reduced computational costs.
In the RASSCF method, all MOs outside the active space are treated at the Hartree-Fock level with the respective moderate costs. To account for dynamic correlation effects more systematically, second order Restricted Active Space Perturbation Theory (RASPT2) can be applied on top of RASSCF [118]. This is very costly for large molecules. The results are also less stable because of the so-called intruder states, although that can be partly mitigated by complex level shifts [137]. The effect of PT2 depends on the system, with relatively small effects for aqueous metal ions [119,131,138] but large effects for systems, like metal hexacyanides, where the active space is not big enough to include all the important ligand-centered orbitals [137,139].
Spin-orbit coupling in these one-component schemes is treated using the numerically efficient Atomic Mean-Field Integral (AMFI) approach, where the scalar relativistic RASSCF states of different multiplicity are coupled in an LS-coupling manner. This coupling simplifies a chemically intuitive analysis of the SOC mixing patterns in the core-excited states [121]; Bokarev et al., 2015; [134]. As the RAS orbitals are optimized without SOC, constraining the 2p 3/2 and 2p 1/2 spinors to have equal spatial distribution, the L 3 /L 2 ratio and energy splitting are approximate, although there is often reasonably good agreement (see Fig. 15).
An active space with five orbitals of metal 3d character, in addition to the 2p orbitals, conceptually corresponds to a ligand-field calculation, while including also ligand-dominated orbitals corresponds to a chargetransfer calculation. It should be noted that the ligand field effect is introduced not as an external electrostatic potential. Instead it is obtained self-consistently within first principles variational procedure allowing for an accurate description of chemical bonding. Moreover, charge-transfer effects are treated on the same footing as the local transitions. RAS explicitly calculates ligand-field strengths, correlation and covalency corrections, as well as relative energies and coupling matrix elements for metal and charge-transfer configurations. Spin-orbit coupling is also included ab initio, but is treated as a perturbation and not as part of the optimization of the MOs.
As a first principles approach, the metal atom of interest is embedded explicitly in its chemical environment. Extended regions of the surrounding can be described as polarizable continuum or using frozen density embedding techniques. There are no restrictions with regard to symmetry. Multiconfigurational methods have typically been used to describe bond dissociation, distorted geometries, near-degenerate ground states, and electronically excited states. The RAS approach is therefore ideally suited to probe excited-state chemical dynamics with xrays, such as bond dissociation or photochemical paths involving conical intersections. A prominent example is the use of RAS method to discriminate between pathways for energy dissipation during photodissociation of iron-pentacarbonyl in ethanol solution [127,140]. However, symmetry can be exploited to reduce computational cost or to target core-and valence-excited states of certain symmetries, as done for aqueous Ni(II), iron-porphyrin model complexes, and studies of both ground state properties and photo-induced dynamics of iron hexacyanides [119,128,121,133].
In summary, the active space approach is exact in principle and systematically improvable by increasing the quality of the basis set and the number of orbitals in the active space. It contains no unknown density functional, no empirical parameters, and no double counting of correlation effects as in DFT-CI. However, due to the complexity and cost of the optimization procedure, RAS calculations cannot be done in a "black box" manner and many user choices, especially the active space design, can affect the final results [119,121,131,141]. The current restrictions on system size and complexity call for the development of more efficient algorithms. Recent progress in multiconfigurational methods, such as the density matrix renormalization group, could lead to improvements in the applicability of RAS to larger and more complex systems.

Methods using the ORCA code (ROCIS and MRCI)
The group in Mülheim is pursuing an approach from the starting point of a restricted-open-shell configuration interaction with singles (ROCIS) [142,143], based on the ORCA platform of first principle calculations [144]. ROCIS starts from a single high-spin open shell determinant and allows all single excitations (at the orbital level) into all unoccupied orbitals. The space of determinants spanned by ROCIS is larger than the particle/hole space that is treated by TD-DFT, which is necessary to describe the multiplet structure. In principle charge transfer is included in such approach non parametrically. In order to include dynamic correlation, this is parameterized in ROCIS/DFT with three empirical parameters for the entire periodic table. This parametrization allows the treatment of large systems.
The ROCIS/DFT method has been combined with pair natural orbitals (PNO), known as the PNO-ROCIS method, which is much faster than original so-called canonical ROCIS method. This allows the calculation of 2p XAS spectra of large metal-organic frameworks, without the need for any parameter tuning [145,146]. Looking at the 2p XAS spectrum for Cr(III)(acac) shows that ROCIS and PNO-ROCIS calculations are essentially identical, but the differences with experiment are still significant [145].
Maganas and coworkers also investigated wavefunction-based multireference calculations to simulate the metal 2p XAS spectral shape. The spectra of high-spin Fe 2+ and Fe 3+ chlorides were calculated with two MR methods, respectively MR configuration interaction (MRCI) and MR Equation of Motion Coupled Cluster (MREOM-CC). The MRCI method shows very good agreement with experiment when the ligand-based orbitals are included. The MREOM-CC calculations yield excellent agreement with small active spaces, suggesting that they efficiently capture the major effects that determine the spectral shape [147].

Experimental aspects of 2p XAS
In this section some of the experimental aspects of 2p XAS spectra are briefly discussed. In particular we discuss some experimental aspects that have to be considered for comparison with accurate calculations. 2p XAS spectra of 3d transition metal systems are mainly measured with synchrotron radiation sources. Since 1987, the experimental resolution of 2p XAS spectra is better than the intrinsic spectral broadening of approximately 0.4 eV full-width half-maximum. In fact, the published 2p XAS spectral shapes have hardly improved after the introduction of high energy resolution grating monochromators [148].
Energy calibration: Synchrotron measurements are usually calibrated to reference spectra to within 0.1 eV accuracy, though not all publications followed such procedure or sometimes used different calibration spectra. Energy calibration is not paid much attention to in the theoretical papers on 2p XAS spectra, as most, if not all, theoretical procedures are not able to do an energy position determination that is more accurate than ~0.3 eV.
Absolute intensity: Synchrotron experiments are usually measured in yield mode, either electron yield or fluorescence yield. This implies that there is no absolute intensity scale and the spectra are published normalised to one. There are very few experiments that provide the absolute cross-section, for example with transmission experiments. The measurements known are in good agreement with the cross-section calculated for atoms [17].
Probing depth: The large intensity of the 2p XAS spectra also implies that the penetration depth is only of the order of 20-30 nm at the peak positions. Electron yield has a measurement depth of approximately 5 nm due to the electron escape length [149]. Fluorescence yield has a measurement depth of several 100 nm [150].
Spectral distortions: Transmission measurements can be affected by the pinhole effect, i.e. the thickness variations in the very thin sample. Fluorescence yield suffers from saturation and self-absorption effects in non-dilute samples. In addition, fluorescence yield is affected by an energy dependence in the integrated x-ray emission intensity, which implies that the 2p XAS spectra appear distorted [151]. Total fluorescence yield is affected by fluorescence from the non-absorbing elements that can cause negative spectral features [152].
EELS: Spectra that are identical to 2p XAS spectra can be measured with electron energy loss (EELS) in an electron microscope. Under the conditions of (i) a primary electron energy larger than ~10 KeV and (ii) small electron scattering angles, EELS measures a spectral shape that is effectively the same as the 2p XAS spectra [153]. EELS has the option to manipulate the electric dipole selection rule and measure at scattering angles where non-dipole transitions become large.

Concluding remarks and open issues
Before discussing the differences between the methods, we would like to distinguish two goals of the multiplet calculations.
1 Theory to derive useful information from experimental 2p XAS spectral shape. 2 Theory to understand all details of the 2p XAS spectral shape.

Using the 2p XAS spectral shape
The level of theory needed to derive useful information form the 2p XAS spectral shape depends largely on the goal of the experiment.
In case of operando studies where a condition is slowly changed, or in case of time-resolved studies, it is often easier (and more accurate) to use a completely experimental approach. Using experimental references, ideally measured at the same beamtime, one can accurately fit the spectral changes. Theory is not as accurate and a basic understanding of the spectral shape is often enough to derive useful information.
In many cases it is important though to provide an interpretation of the spectral shape. For this purpose the user-friendly interfaces have been developed. All interface codes and also all other semi-empirical multiplet codes use essentially model Hamiltonians based on the crystal field multiplet or the crystal field multiplet approaches as indicated in Fig. 20. An important aspect of these calculations is that they are fast and that they are user-friendly. We again note here that the use of the crystal field and the charge transfer multiplet model Hamiltonians is made possible due to the self-screened 2p3d excitation. In contrast, the accurate description of 2p XPS needs more elaborate descriptions of the screening process.
The crystal field multiplet code is most efficient in high-symmetry ionic systems, where the self-screened 2p3d transitions allows a local model Hamiltonian as provided by the crystal field multiplet code. In ionic cubic systems there is essentially only one empirical parameter, the octahedral crystal field 10Dq. Everything else is provided by atomic theory. In more covalent systems one can often remain in crystal field theory, by using the screening of the electron-electron interactions (Slater integrals) by the nephelauxetic effect, equivalent with the use of crystal field theory in optical spectroscopy.
In case of high-valent and/or covalent systems the ionic description is not accurate anymore for the ground state, implying that the charge transfer multiplet model is needed. This is in particular important for molecular systems with π-backbonding due to the inverted screening processes [154]. The charge transfer multiplet model needs more empirical parameters as it combines symmetry effects (in the crystal field) and screening effects (in the charge transfer). Usually only two configurations are needed separated by the charge transfer energy Δ, assumed to be effectively reduced by ~1.0 eV in the final state. In addition the hopping parameters are needed, though they can be estimated from DFT.
In case of low-symmetry covalent systems, it becomes unpractical to estimate all parameters empirically. It then is an option to map first principle calculated parameters into the charge transfer multiplet model. This can be done in many ways. One can map the charge transfer multiplet calculations onto the first principle calculated states and as such determine the charge transfer parameters [155]. A more general method, as indicated in Fig. 21, is to determine general rules to map first principle results onto the charge transfer multiplet model, as was discussed in section 5.  Graphical summary of the calculation of the 2p XAS spectral shape from first principles. The calculation of the complete system involves all states, including translation symmetry in solids. The 2p XAS process is largely a localised phenomenon taking place on one atom.

Understanding 2p XAS spectral shapes
An important aspect for understanding the 2p XAS spectral shape is that the 2p core-hole localizes the final state, implying the need for a local real-space model for an accurate description of the spectral shape. This "revenge of real space" is a major drawback for the accurate description of 2p XAS spectral shapes for solids, where the first principle codes are performed in reciprocal space. This problem does not occur in the description of the 2p XAS spectral shapes for molecular systems, as they usually perform real space first principle calculations. Therefore, we divide the systems into molecules and solids, where we assume that the ground state of the molecules is described by real space quantum chemistry codes and the ground state of solids in reciprocal space.

Understanding the 2p XAS spectral shape with a minimal model Hamiltonian
The minimal Hamiltonian is essentially defined by the self-screened character of the 2p3d excitation, where the final state lifetime broadening of ~200 meV (hwhm) sets the required accuracy regarding energy splitting effects. An example of an accurate simulation of a 2p XAS spectrum of a 3d transition metal ion is the crystal field multiplet fit of Cr 3+ impurities in Al 2 O 3 [156]. The reason for the high accuracy is that the empirical parameters were optimized for the 2p3d RIXS experiments, allowing the separate optimization of the Racah's parameters B and C, which are an alternative representation of the F 2 and F 4 Slater integrals. The model Hamiltonian that was used in these simulations is the crystal field multiplet Hamiltonian. The ground state is described with three parameters 10Dq, Racah B and Racah C. The 3d spin-orbit coupling and the trigonal site distortion are not effective for a 3d 3 ground state with a 4 A 2 ground state. The calculation is limited to the 2p3d dipole transition. Transitions to 4 s and higher states of s and d character are neglected and also quadrupole transitions are neglected. The final state is described with the same three parameters (10Dq, B and C) plus the F 2 , G 1 and G 3 multiplet effects that have been reduced by the same (average) amount as Racah B and C.

Understanding the 2p XAS spectral shape from first principles
Molecules are described in real space using quantum chemical methods. The most accurate methods describe the complete system. The main complication is the calculation of the final state including the 2p core-hole. This is possible with (almost) the same level of accurate theory but this requires the correct description of the screening of the core-hole excitation. One range of methods is based on restricting the active space of the calculation to a sub-set of all states. The spectral shape is for ~90 % determined by strong local interactions and it turns out to be complicated to restrict an ab-initio calculation to the correct active space in order to correctly keep all these strong local interactions. Alternative approaches are based on a two-step approach. The molecule is calculated ab-initio and then all parameters are projected onto the crystal field multiplet or charge transfer multiplet model.
Solids are described in reciprocal space, which implies that the calculation of the 2p XAS spectral shape faces the additional complication that the reciprocal space results much be projected upon a real space description. This "real space revenge" adds an additional step in the process from the ground state calculation to the calculation of the 2p XAS spectral shape.

Open issues and future developments
2p XAS spectra have a lifetime broadening of 200 meV, while the experimental resolution of XAS (and also of EELS) has been better than this lifetime broadening since 30 years. This implies that 2p XAS experimental spectra have not changed over the last 30 years. In this time period, there are a number of theoretical challenges that have not yet been solved.
(1) Excitation energy: No method can calculate the excitation energy within 100 meV accuracy. The relative excitation energies can be calculated more accurately. (2) Cross-sections: The absolute values of the XAS cross-section are calculated accurately from atomic models. The effects of chemical bonding on the cross-sections is likely small, but little is known quantitatively. (3) Broadening: Most methods assume a constant lifetime broadening over the L 3 edge and another larger lifetime broadening of the L 2 edge. The 2p core-hole decays for 99 % via Auger decay and by adding all Auger channels this lifetime can be calculated. The L 2 edge has addition (super Coster Kronig) Auger channels creating the larger lifetime broadening. (4) Autoionization: It is assumed that autoionization is small for 2p XAS but likely it is not completely absent. The M 4 edges in rare earths have significant effects that are described by the Fano broadening parameters. Minor effects should be present in some 2p XAS spectra. (5) Vibrations: Vibrations are not explicitly included in the calculation of 2p XAS spectral shapes. Vibrations are effectively assumed to cause a symmetric broadening of the peaks and as such they are "hidden" in the lifetime broadening and experimental resolution broadening parameters. Data on systems with empty 3d-band, for example SrTiO 3 , seems to indicate the presence of vibrational broadening. (6) Translational symmetry: The effect of band structure on 2p XAS seems limited. The 2p core-hole is such a localized transition that dispersion effects only indirectly enter the simulation of the spectral shape, for example via the screening (charge transfer) channels. Systems where larger dissipative effects are likely, for example pure transition metals, lack spectral details to accurately check theoretical predictions. (7) Dichroism: In this review we did not discuss polarization dependence, angular dependence and magnetic circular dichroism (MCD). It turns out that MCD spectral shapes contain more information than the 2p XAS itself, which might need more detailed theoretical simulations, in addition to the description of the MCD itself. (8) Automation: Other recent developments include the possibility of fitting multiplet simulations to experimental data in order to a) explore the solution space of radial parameters and of Lorentzian/Gaussian broadenings to evaluate uncertainties and reduce user-bias; and b) to evaluate compositions from materials where several metallic sites contribute to the total experimental spectra. These developments include the possibility of directly optimizing underlying radial 3d and radial 2p functions (with respect to experiment) from which, Slater Integrals and crystalfield parameters are calculated. In other words, the experiment serves as the reference (equivalent to the SCF of theoretical first principle methods) to which the radial functions need to be optimized in order to match multiplet simulations to experimental data.
In this review we have focused on 2p XAS spectral shapes. Related spectroscopies do see much progress in the obtainable resolution. This is particularly the case for 2p3d RIXS experiments, where the final states are optical electron-hole excitations that have broadenings that vary from below 1 meV to ~50 meV. Many new features are visible in 2p3d RIXS spectra, particularly at low-loss, which need theoretical description additional to the description of 2p XAS, including phonon and magnon dispersions. In addition the electronic dd-excitations are measured with improved resolution, including their (mixed) angular, polarization and momentum dependences. Apart from the new effects visible in RIXS, the improved resolution from ~20 meV to ~100 meV needs more accurate theoretical calculations.

Declaration of Competing Interest
The authors report no declarations of interest.