Analytic continuationfree Green’s function approach to
correlated electronic structure calculations
Abstract
We present a new charge selfconsistent scheme combining Density Functional and Dynamical Mean Field Theory, which uses Green’s function of multiple scatteringtype. In this implementation the manybody effects are incorporated into the KohnSham iterative scheme without the need for the numerically illposed analytic continuation of the Green’s function and of the selfenergy. This is achieved by producing the KohnSham Hamiltonian in the subspace of correlated partial waves and allows to formulate the Green’s function directly on the Matsubara axis. The spectral moments of the Matsubara Green’s function enable us to put together the real space charge density, therefore the charge selfconsistency can be achieved. Our results for the spectral functions (density of states) and equation of state curves for transition metal elements, Fe, Ni and FeAl compound agree very well with those of Hamiltonian based LDA+DMFT implementations. The current implementation improves on numerical accuracy, requires a minimal effort besides the multiple scattering formulation and can be generalized in several ways that are interesting for applications to real materials.
.pspdf.pdfps2pdf dEPSCrop dNOSAFER #1 \OutputFile
I Introduction
Density functional theory (DFT) Hohenberg and Kohn (1964) in conjunction with the KohnSham scheme Kohn and Sham (1965) and the local density approximation (LDA) Perdew and Wang (1992), or the generalized gradient approximation (GGA) Perdew et al. (1996), to the exchangecorrelation potential has shown great success in the computation of groundstate properties of real materials. However, the method cannot correctly describe materials where electronic correlations are important, such as the Mott insulators, transition metals and lanthanides. One successful approach to improve on the description of the electronic structure of strongly correlated materials is to merge DFT with Dynamical Mean Field Theory (DMFT) Metzner and Vollhardt (1989); Georges et al. (1996); Kotliar et al. (2006). Within DMFT the complicated manybody lattice problem is mapped selfconsistently into a single quantum impurity hybridized with an effective bath. Nowadays impurity problems are efficiently solved by various manybody techniques. Hence DMFT developed into a comprehensive, nonperturbative and thermodynamically consistent theoretical framework for the investigation of correlated electrons on the lattice. The combination of DMFT and DFT, referred to as LDA+DMFT and GGA+DMFT, respectively, has now become the stateoftheart method to study correlated materials Kotliar et al. (2006); Held (2007).
During the last decade, various LDA+DMFT implementations have been proposed. The early implementations employed a twostep procedure: in the first step the LDA problem was solved using an effective oneparticle KohnSham Hamiltonian and the singleparticle wavefunctions (KohnSham basis set) were integrated into the density functional variational approach. The corresponding Green’s function was then obtained using the spectral representation of the KohnSham Hamiltonian. In the second step the interaction problem was treated, i.e., the lowenergy effective Hamiltonian was formulated within a Wannierlike basis obtained through downfolding or, alternatively, by a suitable combination of KohnSham basis sets. This lowenergy Hamiltonian was solved using DMFT. Some of the initial LDA+DMFT implementations kept the effective KohnSham potential fixed, and considered only the selfconsistency of the local selfenergy. Therefore in these approaches the effect of the selfenergy on the electronic charge was neglected. Inserting the selfenergy back into the KohnSham iterative scheme allows one to converge towards selfconsistency in both the selfenergy and charge. Several fully selfconsistent Hamiltonian based implementations have been used within the framework of different basis sets, for example pseudopotential plane waves Leonov et al. (2008, 2010); Amadon (2012), linearized muffintin orbitals Savrasov and Kotliar (2004); Pourovskii et al. (2007); Grechnev et al. (2007); Di Marco et al. (2009); Grånäs et al. (2012) and augmented plane waves Haule et al. (2010). These procedures follow partly the spirit of the spectral density functional theory (SDFT) proposed by Savrasov and Kotliar Savrasov and Kotliar (2004), in which a selfconsistent solution of the Dyson equation is sought. This leads to a quasiparticle Schrödinger (or Dirac) equation with a nonhermitian part in the Hamiltonian.
An elegant way to avoid the difficulties involved in dealing with the nonhermitian Hamiltonian in the SDFT formulation of LDA+DMFT, is provided by the multiple scattering method based on Green’s functions. Green’s function methods have the attractive feature that they can be easily used to treat systems such as surfaces, defects and impurities Gunnarsson et al. (1983); Skriver and Rosengaard (1991). They can also be employed in connection with the coherent potential approximation (CPA) to study substitutional disorder Soven (1967). Common to many Green’s function methods is the problem that the electronic eigenvalue problem is formulated as an energydependent secular equation, from which it is difficult to extract the energy bands. Therefore, the charge density and the total energy, the relevant quantities for the DFT calculation, are obtained by integration of the Green’s function along contours in the complex energy plane Zeller et al. (1982). Some of the first charge selfconsistent implementations of LDA+DMFT with a Green’s function formulation of the KohnSham DFT were implemented within the exact muffintin orbitals (EMTO) method Chioncel et al. (2003) and the KorringaKohnRostoker (KKR) method Minár et al. (2005).
One of the major goals of any selfconsistent LDA+DMFT computation is to answer the question of how the effects of electronic correlation modify the equilibrium properties, like lattice parameters and bulk modulus, beyond the LDA. It is hence necessary to calculate accurate total energies within LDA+DMFT, from which the equilibrium quantities can be derived. Several of the ground state quantities and spectral properties have already been discussed Leonov et al. (2011, 2014) within the Hamiltonian framework. Despite the many successes of Green’s functionbased LDA+DMFT methods Minár (2011), several numerical difficulties still remain for total energy calculations. When the Green’s function based LDA+DMFT scheme is executed in practice, Padé approximants Vidberg and Serene (1977) (rational polynomials) are used to pass Green’s functions from the complex energy contour to the Matsubara frequencies, and to return with the selfenergy from the Matsubara frequencies back to the complex contour. Besides being sensitive to numerical noise Beach et al. (2000), Padé approximants may miss important features, that can only be captured by resummation of the continued fraction to infinite order Wall (1948); Baker (1975); Beach et al. (2000). In recent years some methods have been proposed in order to improve on the original Padé approximation technique Östlin et al. (2012); Schött et al. (2016); Nordström et al. (2016) to the analytic continuation of the Green’s function, but as of yet no fully satisfactory solution to this problem exists. Such numerical problems are presently a bottleneck for an accurate and stable selfconsistent Green’s function based LDA+DMFT method that can produce reliable total energies.
The success of LDA+DMFT consists in its ability to produce a selfconsistent, numerically manageable approximation for the spectral function and for lattice properties at equilibrium. It is desirable that LDA+DMFT developments be exact in principle, and that even approximate perturbative solvers should give good results, irrespective of whether a Hamiltonian or Green’s function method is used. For these reasons it is essential to pursue alternative methods that improve on the numerical accuracy. In general, for a Green’s function formulation of the LDA+DMFT the knowledge of the noninteracting Green’s function along the imaginary axis is required. Consequently, our primary objective of the present paper is to describe an approach which yields an accurate Green’s function in Matsubara frequencies which can be used in the DMFT part and, at the same time, in constructing the charge density.
Our novel method makes the analytic continuation during the selfconsistent KohnSham iterations unnecessary. The key observation that triggered this method development is that the charge density is the only ingredient needed to close the KohnSham selfconsistent loop. The charge density difference between correlated and noncorrelated calculations, evaluated on the imaginary Matsubaraaxis, is taken as the correction on the DFT level charge density. Quantities like eigenvalues, Green’s functions and selfenergies are only auxiliary quantities in this respect. In the method, MTO+DMFT, presented here the Green’s function in Matsubara frequencies is evaluated from the LMTO eigenstates, i.e., in the basis of linearized partial waves. The choice to take the character in the denomination MTO+DMFT is to remind of the fact that the Green’s functions in DFT are usually computed along a general complex contour mesh, i.e. , for a given muffintin potential. We implemented this scheme starting from our previous EMTO+DMFT method Chioncel et al. (2003), which has been successfully used to study correlated systems, such as transition metals and compounds Chioncel et al. (2003, 2007); Katsnelson et al. (2008), magnetic heterostructures Beiuşeanu et al. (2011) and transport properties through layered structures Chioncel et al. (2015). The use of a Green’s function method opens the possibility to study systems that deviate from perfect crystalline conditions, such as alloys and surfaces.
Ii Overview of the muffintin formalism
Muffintin based methods have in common that they partition space into spherical muffintins, centered around the ions in the lattice, and the interstitial, the area outside of the muffintins. Inside the muffintins the effective potential is assumed to be spherically symmetric, while it is taken to be constant in the interstitial. The KohnSham equations are solved separately within these regions, and the solution for the entire space is found by imposing boundary conditions between the muffintins and the interstitial. The algebraic formulation of the matching conditions takes the form of a secular equation, which is in general energydependent. Sec. II.1 describes this concept for the EMTO method. Sec. II.2 briefly reviews the concept of basis function linearization, which is important for the construction of the correlated orbitals in this work.
ii.1 Charge density and the complex contour Green’s function in the EMTO basis set
Within the muffintin formalism, the effective KohnSham potential ( denotes the spin) in the singleelectron KohnSham equations, labeled by state index ,
(1) 
is approximated by spherical muffintin wells centered at lattice sites . The exchangecorrelation part of will in the following always be approximated by the spinpolarized LDA, and we will from now on suppress the spin index for simplicity. For the EMTO basis set Andersen et al. (1994); Vitos et al. (2000); Vitos (2001, 2010), the oneelectron wavefunctions are expanded in exact muffintin orbitals ,
(2) 
where denotes the orbital and azimuthal quantum numbers respectively, and , where the vector notation for the index has been omitted. The superscript denotes the screening parameter. The orbitals are linear combinations of partial waves , which are normalized solutions of the radial Schrödinger eqution inside the muffintins with spherical potential ,
(3) 
and of the solutions in the interstitial region Vitos (2010). The angular momentum sum in Eq. (2) can in practice be truncated at , making the basis minimal. Since the orbitals are centered around the lattice sites , the basis is “local”, making it suitable as a basis for correlated orbitals within DMFT. The coefficients are determined from the condition that the expansion should fulfill Eq. (1) in all space, i.e. the orbitals should be everywhere continuous and have no derivative discontinuities (kinks) anywhere. In the EMTO formalism this leads to the kink cancellation equation:
(4) 
which is equivalent to the KKR tail cancellation equation Vitos (2010), written in a screened representation. The quantity defines the kink matrix for an arbitrary complex energy and has the form:
(5) 
denotes the EMTO logarithmic derivative function Vitos et al. (2000); Vitos (2001), and is the slope matrix Andersen and SahaDasgupta (2000). Note that Eq. (4) is an energydependent secular equation, which allows one to determine the eigenvalues . These are obtained using numerical search algorithms for the roots of the secular determinant along the real energy axis. To simplify the notation further, we suppress the index for the screening parameter .
For translation invariant systems, the index runs over the atoms in the primitive cell only, and the Fourier transformation of Eq. (5) produces a matrix equation in the reciprocal space:
(6) 
that is solved using Green’s function methods. Accordingly, the path operator is the unique solution of Eq. (6) (the inverse of the kink matrix ) that fulfills the combination of lattice symmetry and boundary conditions. The elements of the kink matrix are constructed from the Bloch wave vector () dependent slope matrix Vitos (2010). Since the energy derivative of the kink matrix, , gives the overlap matrix for the EMTO basis set Andersen and SahaDasgupta (2000), these are used to normalize the path operator and construct the matrix elements of the EMTO Green’s function Vitos et al. (2000); Vitos (2001)
(7)  
where accounts for the unphysical poles of Vitos (2001, 2010). The total number of states at the Fermi level is obtained as
(8) 
where the energy integral is carried out on a complex contour that cuts the real axis below the bottom of the valence band and at . The ksummation is performed over the Brillouin zone (BZ).
To close the KohnSham selfconsistency scheme requires the computation of the charge density. Within the EMTO method this is achieved through the real space path operator (corrected for unphysical poles similarly as in Eq. (7) Vitos (2010)) integrated over the same complex contour that is used to determine ,
(9)  
where is a real Gaunt number. Eq. (9) is valid within the muffintin spheres and for , and , where is a normalization function Vitos (2001, 2010). The specific set of real harmonics is denoted by .
ii.2 Charge density and the Matsubara Green’s function in the LMTO basis set
An alternative solution of Eq. (1) is obtained by the linearized muffintin orbitals (LMTO) Andersen (1975); Andersen et al. (1986) method. The same muffintin shape is used for the potentials as in the EMTO method, but with the additional approximation that the interstitial region is neglected, leading to the atomic sphere approximation (ASA). The LMTOs are constructed from the partial wave solutions inside the muffintin spheres, computed at an arbitrary energy (commonly chosen as the center of gravity of the occupied part of the band), and from the energy derivative of the partial wave, , viz.
(10) 
The omitted energy argument of the partial wave means that the function is evaluated at an energy . In Eq. (10), is defined as
(11) 
where is the KohnSham Hamiltonian in the socalled nearly orthogonal representation Andersen and Jepsen (1984); Andersen et al. (1986) viz.
(12) 
where are the LMTO structure constants, and the potential parameters and are computed from the partial waves according to the prescription given in Ref. Andersen et al., 1986. With the energy independent LMTO basis functions, Eq. (10), the lattice wave function (i.e. the linear muffintin wave function):
(13) 
follows the energyindependent eigenvalue problem:
(14) 
where the Hamiltonian eigenvalues provides the band structure, and the eigenvectors contain Bloch vector specific information.
ii.2.1 Moments from the LMTO eigenstates and complex contour
Once the LMTO Hamiltonian has been diagonalized, Eq. (14), the energy moments can be evaluated as
(15) 
where the and moments correspond to the orbitals occupation and oneelectron energies, respectively. Note that the moments computed with the help of Eq. (15), is along the real energy axis.
To make contact with DMFT we point out that the LMTO method has been already used to construct Green’s functions: either from the potential parameters directly or from the Lehmann (eigenvalue) representation Andersen et al. (1986); Skriver and Rosengaard (1991); Pourovskii et al. (2007):
(16) 
The energy moments can then be computed along a similar complex contour as in the EMTO method Skriver and Rosengaard (1991), using the site and orbital diagonal part of the Green’s function, viz.
(17) 
where we remind the reader of the definition . The eigenvalue summation done in Eq. (15), is now replaced with the complex contour integration Eq. (17). The knowledge of the moments and the partial waves allows the computation of the charge density Andersen et al. (1986), viz.
(18)  
and the DFT selfconsistency loop can be closed.
Note that one advantage of the LMTO Green’s function over a multiplescattering Green’s function is that its spectrum is discrete and upwards bound, i.e. it does not contain the freeelectron continuum Weinberger (1990).
ii.2.2 Moments from Matsubara LMTO Green’s function
Eq. (16) can be also defined for the Matsubara frequencies , where , and is the temperature. Pourovskii et al. Pourovskii et al. (2007), showed recently that the LMTO zeroth energy moments can be extracted also from the imaginary frequency domain by standard Matsubara summation Fetter and Walecka (1971), viz.
(19) 
with the resolved Green’s function given by the Lehmann representation
(20) 
The local Green’s function is computed as:
(21) 
The higher order moments can be calculated as products of the zeroth order moment , and ,
(22) 
The charge density can be computed again from Eq. (18). Note that a cutoff at a finite frequency will lead to inaccurate Matsubara sums Nicholson et al. (1994). This can be corrected to some extent by taking the analytic tail of the Green’s function into account Pourovskii et al. (2005, 2007).
ii.3 Incorporating the local manybody selfenergy
After the brief review of the energydependent and the energylinearized basis sets we proceed with discussing a combination of these methods which allows to include the local DMFT selfenergy in a charge selfconsistent way. The DMFT maps selfconsistently the manybody lattice problem to an impurity model, which can be solved by various manybody techniques and produces the impurity Green’s function and the local selfenergy Kotliar et al. (2006). The DMFT selfconsistency condition is obtained by imposing that the impurity Green’s function is the same as the lattice local Green’s function.
In the EMTO+DMFT method Chioncel et al. (2003), the selfconsistent procedure starts with a guess for the local selfenergy to be combined, through the Dyson equation, with the resolved LDA Green’s function, Eq. (7), which represents the “noninteracting” lattice Green’s function:
(23) 
The local Green’s function is extracted from Eq. (II.3) on the complex contour: . Its matrix elements are analytically continued to the Matsubara frequencies:
(24) 
In the next step one has to construct the bath Green’s function which specifies the impurity problem, which within EMTO+DMFT is computed from the analytically continued lattice local Green’s function and the selfenergy:
(25) 
The manybody problem is solved on the Matsubara axis, and the resulting selfenergy is then analytically continued to the semicircular contour:
(26) 
in order to close the LDA+DMFT loop. In Figure 1 we illustrate the contours used in the EMTO+DMFT calculations. Accordingly, the selfconsistency procedure requires two Padé analytic continuation Vidberg and Serene (1977); Chioncel et al. (2003) steps, that has to be controlled numerically.
In order to close the charge selfconsistent loop, the LDA+DMFT path operator is extracted from the interacting Green’s function (II.3), while the realspace charge density is computed according to Eq. (9) substituting the LDA path operator with the corresponding LDA+DMFT path operator. The new effective KohnSham potential is obtained by solving the Poisson equation, and the scheme is iterated until selfconsistency is achieved.
The LMTO method has previously been used as a choice for charge selfconsistent basis sets. In particular, Pourovskii et al. Pourovskii et al. (2007) implemented an LDA+DMFT scheme in the LMTOASA method. In the case of LMTOASA, the LDA level Green’s function is easily evaluated along the imaginary axis (Eq. (20)), and the selfenergy is embedded via the Dyson equation to obtain the LMTO LDA+DMFT level Green’s function. After performing the ksum, the bath Green’s function is given similarly as in Eq. (25), and is given as input to the DMFT impurity solver. In order to close the charge selfconsistent loop, the energy moments are computed as in Eqs. (19II.2.2), with the exception that the Green’s function in Eq. (19) is now on the LDA+DMFT level. The charge density is then computed from the energy moments as outlined in Eq. (18).
Iii Implementation of the new MTO+DMFT method
iii.1 Motivation
In this Section, we present a novel scheme that removes the need for the illposed analytic continuations Eqs. (24) and (26), during the selfconsistent loops. Two main ideas are used to achieve this: (i) the Green’s function can be well approximated by linearization of the muffintin orbitals, and (ii) the charge density can be calculated by Matsubara summation.
iii.1.1 Elimination of :
the benefit of a linearized basis set
A major difference between the Green’s function within EMTO, Eq. (II.3), and within LMTOASA, Eq. (16), is that the latter can be easily evaluated for any energy once the potential parameters are known. The EMTO Green’s function on the other hand requires the computation of the slope matrix and the solution of the radial Schrödinger equation at each energy point along the complex contour, and this is a numerically demanding task. The two Green’s functions should be equivalent up to the error in the linearization imposed on the kink cancellation condition Tank and Arcangeli (2000), reflecting the error of the linearization of the muffintin basis set.
Based on the formal equivalence of these methods, and the similar results for the corresponding quantities (Green’s functions and moments of these), we propose the following:

The EMTO Green’s function should be used for LDA calculations,

The LMTO Green’s function should be used for DMFT calculations.
This replaces the need of a Padé approximant with a linearization of the basis set, a more well controlled approximation.
To be specific, we outline the procedure: At each KohnSham iteration, the kink matrix in Eq. (5) is set up for the complex energies along the contour, and the EMTO Green’s function is used to solve the electronic structure problem as outlined in Sec. II.1. The partial waves are obtained by radially integrating the Schrödinger equation (3) for the linearization energy . From these partial waves, the LMTO potential parameters and can be obtained (see Ref. Skriver, 1984). The LMTO Hamiltonian (12) is constructed and diagonalized, providing eigenvalues and eigenvectors . In the next step, the noninteracting local LMTO Green’s function (21) is computed for the Matsubara frequencies . Correlation effects are generated by the interaction term, formally to be added to the noninteracting Hamiltonian . The explicit form of the four index Coulomb interaction matrix elements is discussed in Sec. IV. From the Green’s function formulated on the Matsubara axis the bath Green’s function (25) at the LMTO level is obtained, and passed into the DMFT manybody solver.
The error in the linearization can be assessed by comparing the density of states (DOS) arising from the EMTO and the LMTO Green’s functions, both at LDA level, see left panel of Fig. 2. The EMTO method was iterated selfconsistently for Ni (above left) and Fe (below left), using an basis set. The DOS was then evaluated from the imaginary part of Eq. (7) (black solid lines) and Eq. (16) (red dashed lines), along a horizontal contour close to the real energy axis. The curves are in good agreement with each other. The basis set linearization will introduce approximations, but these are easily controlled and can in principle be improved by including higher order MTOs (the NMTO method Andersen and SahaDasgupta (2000)).
iii.1.2 Elimination of :
Charge density difference
An essential step for the charge selfconsistency of LDA+DMFT with the EMTObasis set Chioncel et al. (2003) is the analytic continuation of the selfenergy back to the complex contour, which allows to update the path operator from which the real space charge density, Eq. (9), is computed. The correlation effects upon the real space charge density has been analyzed in the previous LDA+DMFT implementation for Fe, Ni and Cr Chioncel et al. (2003). In particular for Cr, LDA+DMFT charge density shows accumulation of electrons due to correlation effects inside the muffintin spheres and a depletion of density in the interstitial region. To capture these correlation induced corrections to the LDA charge density it seems natural for the current implementation to propose the following scheme:

The LDA charge density should be computed within EMTO, , on the complex contour,

The DMFT charge density correction, , should be computed within LMTO on the Matsubara axis.
To be specific, we outline the procedure: The LDA real space charge density is calculated from the complex contour, see Eq. (9). Once the LMTO Green’s function has been constructed on the Matsubara frequencies, the energy moments Eqs. (19)(II.2.2) are computed both on the LDA and the LDA+DMFT level. This allows to evaluate the charge density according to Eq. (18). The charge density difference is then simply defined as
(27) 
where the superscript of emphasize that this quantity is computed on the imaginary axis. The final LDA+DMFT real space charge density is obtained through
(28) 
and is used to close the selfconsistent cycle. Note that the charges computed along the Matsubara axis contain contributions from all orbitals, and not only from the correlated subset.
To assess the possible differences between the EMTO and LMTO charge density, at the LDA level, we plot in Fig. 2 (right column) the valence charge density for Ni/Fe in the upper/lower panel. The EMTO charge densities (black solid lines) were iterated to selfconsistency and evaluated according to Eq. (9). The LMTO charge (red dashed lines) was evaluated from the EMTO selfconsistent potentials by computing first the energy moments of the LMTO Green’s function Eq. (16) using the contour integration Skriver and Rosengaard (1991), and then applying Eq. (18). The charge densities are in a very good agreement.
iii.1.3 Total energy
Within the KohnSham scheme, the total energy functional can be expressed as
(29)  
where is the external ionic potential, is the exchangecorrelation energy and is the KohnSham singleparticle kinetic energy. The square brackets indicate that the energy components are functionals of the density . For the proposed new method, the charge density given as input is now computed on the LDA+DMFT level, Eq. (28), as outlined in the previous section. A slight change in the expression of the kinetic energy:
(30)  
has to reflect the change in the oneelectron energies caused by the presence of the real part of the selfenergy. Eq. (1) was used for the second equality, of the above equation. In order to account for this change in the oneelectron energies, the difference between the LDA and LDA+DMFT oneelectron energies is added to the expression for the kinetic energy. The total energy of a manybody system in the ground state includes also the GalitskiiMigdal contribution Fetter and Walecka (1971). This contribution is added into all LDA+DMFT computations. Other formulations such as the variational LuttingerWard functional may give improved energies Savrasov and Kotliar (2004); Kotliar et al. (2006); Haule and Birol (2015) but do not appear straightforward to implement in the present scheme. In the current implementation the GalitskiMigdal energy contribution is computed on the Matusbara axis in the LMTO formulation:
(31)  
where is on the LMTO LDA+DMFT level. The final expression for the LDA+DMFT total energy is
(32) 
The KohnSham oneelectron energies from the DFT (LDA) calculation already include some interaction effects through the Hartree and the exchangecorrelation potential terms. Including the interactions explicitly in the form of the Hubbard Hamiltonian, some interaction contributions would be counted twice. Consequently, some double counting correction has to be included. There is no universal solution to this problem, and most of the double counting schemes are empirical. In the present method we take over the schemes used in the previous implementation Chioncel et al. (2003), a detailed discussion is found in Ref. Petukhov et al., 2003.
iii.2 Flow Diagram of the selfconsistency calculation in Mto+dmft
The ideas presented in the previous section can be condensed in the following scheme that we call the MTO+DMFT method (see Fig. (3)):

The KohnSham iterations are initiated with a starting guess for the effective potential and the selfenergy .

The kinkcancellation equations are constructed for points along the complex contour, and the LDA level charge , Eq. (9), is obtained by integrating along the contour. At this stage, the LMTO potential parameters are also computed from the partial waves.

The Hamiltonian is constructed from the potential parameters from step (2) using Eq. (12), and the eigenvalue problem is solved.

The charges are obtained by Matsubara summation, and the difference according to Eq. (27) is evaluated.

The final LDA+DMFT charge density Eq. (28) is computed by adding from step (6) to the DFT charge density from step (2).

Return to step (2) until selfconsistency in both the charge and selfenergy is reached.
Once the selfconsistency has been reached, observables such as the total energy Eq. (32) and spectral functions can be evaluated. Note that the spectral functions are evaluated on a horizontal contour slightly shifted away from the real axis. To analyze the selfenergy along the real axis, a Padé approximant is can be used. Note however, that this does not affect the KohnSham loops, and has to be carried out only once at the end, after selfconsistency has been reached. In this case is also easy to identify spurious poles in the Padé approximant, as outlined in Ref. Östlin et al., 2012.
Iv Results
To assess the implementation electronic structure calculations have been performed according to the method proposed above. Transition metals and compounds in which the orbitals form the correlated basis set have been considered. For the DMFT impurity solver a fluctuation exchange (FLEX) Bickers and Scalapino (1989) type of approximation was used for the multiorbital case Lichtenstein and Katsnelson (1998); Katsnelson and Lichtenstein (1999); Pourovskii et al. (2005). In contrast to the original formulation of FLEX Bickers and Scalapino (1989), the spinpolarized matrix FLEX (SPTFLEX), used for the present calculations treats the particleparticle and the particlehole channel differently Lichtenstein and Katsnelson (1998); Katsnelson and Lichtenstein (1999); Pourovskii et al. (2005). While the particleparticle processes are important for the renormalization of the effective interaction, the particlehole channel describes the interaction of electrons with the spinfluctuations. In addition the advantage of such a computational scheme is that the electronelectron interaction term can be considered in a full spin and orbital rotationally invariant form, viz. . Here, annihilates/creates an electron with spin on the orbital at the lattice site . The Coulomb matrix elements are expressed in the usual way Imada et al. (1998) in terms of Slater integrals. Since specific correlation effects are already included in the exchangecorrelation functional, socalled “double counted” terms must be subtracted. To achieve this, we replace with Lichtenstein et al. (2001) in all equations of the DMFT procedure Kotliar et al. (2006). Physically, this is related to the fact that DMFT only adds dynamical correlations to the DFT result Petukhov et al. (2003).
iv.1 Transition metals: nickel and iron
Within the family of the late transition metals, nickel and iron are known to show in their band structures signatures of electronic correlation Lichtenstein et al. (2001). Nickel is wellknown for a “6eVsatellite” in its photoemission spectra Liebsch (1979), while a similar satellite in iron is debated Grechnev et al. (2007); SánchezBarriga et al. (2009).
For both Fe and Ni, a basis was used, and the and states were treated as valence. The core electron levels were recalculated at each KohnSham iteration (softcore approximation). The kink cancellation condition was set up for 16 energy points distributed around a semicircular contour with a diameter of 1 Ry, enclosing the valence band. The BZ integrations were carried out on an equidistant mesh with 285 kpoints (for Fe) and 240 kpoints (for Ni) in the irreducible BZ. For the exchangecorrelation potential the local spin density approximation with the PerdewWang parameterization Perdew and Wang (1992) was used. For the DMFT impurity calculations, the Matsubara frequencies were truncated after 2048 frequencies, and the temperature was set to K. The values for the Coulomb and the exchange parameters are discussed in connection with the presentation of the results in each case. The equations of state were obtained by fitting the energyversusvolume data to a BirchMurnaghan function Birch (1947). The densities of state were computed along a horizontal contour shifted a distance Ry away from the real energy axis. At the end of the selfconsistent calculations, to obtain the selfenergy on a real energy mesh, can be analytically continued into a horizontal contour by a Padé approximant constructed by the Thiele method Vidberg and Serene (1977).
Ni  LDA  eV  eV  Exp.  
This work  67.65  75.84  (0.12)  86.04  (0.27)  
FPLMTO (Ref. Di Marco et al., 2009)  67.88  76.20  (0.12)  89.48  (0.31)  73.79 
KKR (Ref. Di Marco et al., 2009)  66.86  76.28  (0.14)  85.53  (0.28)  
This work  259  162  (0.37)  99  (0.62)  
FPLMTO (Ref. Di Marco et al., 2009)  260  163  (0.37)  84  (0.68)  179 
KKR (Ref. Di Marco et al., 2009)  280  171  (0.39)  132  (0.53)  
Fe  LDA  eV  Exp.  
This work  70.09  86.21  (0.23)  
FPLMTO (Ref. Grånäs et al., 2012)  70.49  87.06  (0.24)  79.46  
This work  253  124  (0.51)  
FPLMTO (Ref. Grånäs et al., 2012)  234  90  (0.62)  163 
In the top left part of Fig. 4, the LDA and LDA+DMFT density of states for Ni is presented. The volume was set to the experimental value (73.79 a.u.). The new method compares well with previous DFT+DMFT studies employing the SPTFLEX impurity solver Chioncel et al. (2003); Minár et al. (2005); Grechnev et al. (2007), and captures the main correlation effects of Ni such as the satellite formation and band narrowing. Note that the correlation effects are stronger in the majority spin channel (more pronounced satellite, more narrow bandwidth) than in the minority spin channel, which is common for the late metals. For the case of eV (blue line), the position of the “6eV” satellite is at higher binding energy than in experiment. The value eV has previously given the correct satellite position when a quantum Monte Carlo impurity solver was used Lichtenstein et al. (2001), and the fact that the SPTFLEX solver overestimate the effect of correlation is thought to be due to the perturbative nature of the solver Katsnelson and Lichtenstein (2002). Recent spinpolarized positron annihilation experiments and LDA+DMFT calculations allowed to determine the value for the local electronelectron interaction strength in ferromagnetic nickel to the value of eV Ceeh et al. (2016). By decreasing the Coloumb parameter to eV (red line), the satellite is shifted to lower binding energy, in better agreement with experiment, as found previously Katsnelson and Lichtenstein (2002).
The top right part of Fig. 4 shows the equation of state of Ni as calculated within the new method, for various values of the Coulomb parameters and . The effect of correlation can be seen to increase the equilibrium volume from the value given by the LDA (corresponding to , black line). The equilibrium volumes are given in Table 1, together with the bulk moduli. As already mentioned in the discussion of the nickel DOS, the SPTFLEX solver overestimates the effect of correlation. This is seen for the equilibrium volume, where the commonly accepted value of eV (blue line) overestimates the equilibrium volume. eV (red line) gives a better agreement with the experimental volume. It should also be noted that the bulk modulus is softened as is increased, which corrects for the overestimation made by the LDA functional.
Fig. 4 shows the DOS (bottom left) and equation of state (bottom right) for bcc Fe, for the case of standard LDA () and for eV, eV. Similar values of and have previously been successfully used to describe the photoemission spectra and energetics of iron SánchezBarriga et al. (2009); Leonov et al. (2011); Grånäs et al. (2012). The effect of correlation is seen to broaden the peaks in the DOS, and create a satellite structure at eV binding energy, in agreement with previous SPTFLEX studies Grechnev et al. (2007); Grånäs et al. (2012). By including local correlation effects, the equilibrium volume is increased, similar as for Ni. This can be seen in the bottom right part of Fig. 4, where the equation of state is given. The effect of correlation also reduces the bulk modulus (see Table 1). The agreement between our results and the ones from the Ref. Grånäs et al., 2012 is very good, the slight differences are due to the spinorbit coupling explicitly present in Ref. Grånäs et al., 2012. On the other hand it is known that spinorbit effects are quite small for Fe Eriksson et al. (1990).
iv.2 Iron aluminium
The stoichiometric intermetallic compound FeAl has attracted the interest of the electronic structure community mainly due to its magnetic properties. While FeAl is paramagnetic in experiment, LSDA calculations within density functional theory predict an ordered ferromagnetic ground state with a magnetic moment of about . Mohn et al. Mohn et al. (2001) showed that including the effect of the local Coulomb interaction through the LDA+U method the nonmagnetic state can be stabilized for a narrow range of values. It was further argued that the reduction in the DOS at the Fermi level, caused by increasing values, will favor the nonmagnetic state through the Stoner criteria. Petukhov et al. Petukhov et al. (2003) pointed out the importance of dynamic effects by LDA+DMFT calculations of the spectral functions, showing that the nonmagnetic solution is stable within LDA+DMFT, and that the DOS is pinned to the Fermi level. Later on, Galler et al. Galler et al. (2015) confirmed this, while also computing susceptibilities for FeAl within LDA+DMFT using a continuoustime quantum Monte Carlo (CTQMC) impurity solver. None of the above previous LDA+DMFT studies presented total energies.
We have investigated the electronic structure of FeAl with our new method, in order to evaluate the density of states and the total energy for volumes around the experimental value. FeAl crystallizes in the B2 (CsCl) structure, i.e. a simple cubic lattice with Fe at position () and Al at (), where the experimental lattice constant is a.u. Mohn et al. (2001) (note that also the value a.u. is reported in the literature Sundararajan et al. (1995); Kulikov et al. (1999)). An basis was used, and a contour of diameter 1 Ry with 16 energy points was employed for the energy integrations. For the BZ integration 286 kpoints in the irreducible part was employed.
In the left part of Fig. 5 we present the nonmagnetic density of states for FeAl, computed assuming eV (black line) and eV (blue line). As pointed out in the previous studies Petukhov et al. (2003); Galler et al. (2015), the increasing of the Coulomb parameter, within LDA+DMFT has little effect on the density of states at the Fermi level, in contrast to LDA+U calculations Mohn et al. (2001), while it leads to a band narrowing. This is an indication that spinfluctuations, which are included on a perturbative level in the SPTFLEX solver, changes the simple picture of Stoner instability.
In the right panel of Fig. 5, our computed total energies for ferro and nonmagnetic FeAl are presented. In the case of LDA (, bottom right), the nonmagnetic total energy (black line) is never lower than the ferromagnetic total energy (red line), for all the studied volumes. In the lower volume range the ferromagnetic moment is lost, indicated by the coincidence of the two energy curves a.u.. The fact that a ferromagnetic ground state is favored in LDA is in agreement with previous DFT studies Kulikov et al. (1999). The equilibrium volumes for the respective curves are 73.95 a.u. ( a.u.) for the ferromagnetic curve, and 73.62 a.u. ( a.u.) for the nonmagnetic curve, and hence the ferro and nonmagnetic lattice constants differ by only. Previous DFT studies have found lattice constants of value a.u. (TBLMTO, nonlocal corrections to the LDA, Ref. Kulikov et al., 1999), a.u. (TBLMTO, BarthHedin parametrization of LDA, Ref. Sundararajan et al., 1995) and a.u. (fullpotential linearized augmented Slatertype orbital method using LDA, Ref. Watson and Weinert, 1998), using different basis sets and exchangecorrelation functionals. The previously reported lattice constants are all larger than the current results, but are consistent given the fact that different basis sets and exchangecorrelation functionals were used.
As local correlation effects are taken into account within LDA+DMFT ( eV, top right), the situation is reversed. In this case the ferromagnetic solution is always higher in energy compared to the nonmagnetic solution, indicating that the nonmagnetic solution is the ground state for the whole volume range. For volumes a.u., the magnetic moment is lost, and the two curves coincide. The equilibrium volumes for the respective curves are 80.99 a.u. ( a.u.) for the ferromagnetic curve, and 82.67 a.u. ( a.u.) for the nonmagnetic curve, which is in good agreement with experiment.
Associating the analysis of the DOS and equation of state, we see that LDA+DMFT is able to explain the experimentally observed fact that FeAl is in a nonmagnetic ground state, while at the same time providing an equilibrium lattice constant in better agreement with experiment than the LDA. By investigating the DOS, the Stoner criterion (an increased DOS at is leading to a magnetic instability) for ferromagnetism can be ruled out as an explanation for the magnetism in FeAl.
V Conclusion and Outlook
In this paper we have introduced a new computational scheme for LDA+DMFT calculations, using Green’s function methods. The new method is able to describe correlated systems such as transition metals and compounds, and shows results in very good agreement with previous LDA+DMFT implementations. At the heart of the current implementation is the formulation of the LDA Green’s function directly on the Matsubara axis, using the Lehmann representation in terms of the eigenvalues and eigenfunctions of the LMTO Hamiltonian. This simple procedure is essential for circumventing the analytical continuation of the Green’s function from the complex contour to the Matsubara frequencies (Sec. III.1.1). The real advantage of this construction appears in the computation of the charge density. Starting from the zeroth moment of the LMTO Green’s function, the extension to higher order moments becomes possible. From these moments the real space charge can be constructed. The difference between correlated and noncorrelated charge density allows for the selfconsistency and in the same time circumvent the second analytical continuation, that of the selfenergy from the Matsubara axis to the complex contour (Sec. III.1.2). The idea to consider charge density differences between LDA and LDA+DMFT might also prove useful for Hamiltonian based methods, since the operation of subtraction could help in reducing systematic errors coming from the numerically difficult Matsubara sums.
By sidestepping the illposed analytic continuation problems, a numerically stable implementation is possible, at the minor cost of performing basis set linearization for the calculations along the imaginary axis.
Numerical results are presented for Fe and Ni. A direct numerical comparison between the imaginary part of the EMTO and the LMTO Green’s functions along a horizontal contour close to the real axis is studied in Fig. 2. The agreement between the basis sets as well as for radially distributed real space charge are found to be excellent. The MTO+DMFT densities of states and total energy curves are then presented in Fig. 4, and are found to be in very good agreement with previous LDA+DMFT studies that were employing other basis sets. As a final example, the spectral functions and equations of state of the FeAl transition metal compound is studied. Similarly, an excellent agreement is found when comparing to previous LDA+DMFT methods Petukhov et al. (2003). For a Coulomb interaction strength of magnitude eV (on Fe in FeAl), the total energies for FeAl are seen to favor a nonmagnetic groundstate, in accordance with experiment.
As an outlook, we propose several possibilities to extend the current MTO+DMFT implementation. First, the downfolding of the linearized basis set can be included Pourovskii et al. (2007), in order to reduce the size of the minimal basis set even further. Second, the fullcharge density (FCD) technique Vitos et al. (1994) applied to the EMTO method has previously provided accurate total energies for lowsymmetry structures, while still keeping the efficiency of the spherical potential approximation (see Ref. Vitos (2010) The implementation of the FCD into the MTO+DMFT method would make it possible to study the energetics of lowsymmetry structures and anisotropic lattice distortions of correlated materials, which currently is work in progress. Finally, a major motivation is to enable a combination of the present method with the coherentpotential approximation Soven (1967), or with the typical medium theory for disorder Terletska et al. (2017). This would provide a method that could handle strong correlation and disorder in alloy systems, including the problem of Anderson localization Terletska et al. (2017).
In conclusion we have attempted to demonstrate by means of elementary examples that the current MTO+DMFT, in conjunction with the perturbative SPTFLEX solver, can successfully describe the electronic structure and energetics of transition metals and their compounds. Even though the SPTFLEX solver is numerically simple due to its algebraic structure, it is still sufficiently rigorous to deal with correlated electrons in condensed matter. A more sophisticated implementation using a variant of Continuous Time Quantum MonteCarlo, DMFT impurity solvers is in progress.
Acknowledgments
We greatly benefited from the discussions with D. Vollhardt and O. K. Andersen, whose advices are gratefully acknowledged. Financial support of the Deutsche Forschungsgemeinschaft through FOR 1346 is gratefully acknowledged. A.Ö. acknowledges helpful discussions with I. Di Marco. L.V. acknowledges financial support from the Swedish Research Council, the Swedish Foundation for Strategic Research, the Swedish Foundation for International Cooperation in Research and Higher Education, and the Hungarian Scientific Research Fund (OTKA 84078 and 109570). We acknowledge computational resources provided by the Swedish National Infrastructure for Computing (SNIC) at the National Supercomputer Centre (NSC) in Linköping.
References
 Hohenberg and Kohn (1964) P. Hohenberg and W. Kohn, Phys. Rev. 136, B864 (1964).
 Kohn and Sham (1965) W. Kohn and L. J. Sham, Phys. Rev. 140, A1133 (1965).
 Perdew and Wang (1992) J. P. Perdew and Y. Wang, Phys. Rev. B 45, 13244 (1992).
 Perdew et al. (1996) J. P. Perdew, K. Burke, and M. Ernzerhof, Phys. Rev. Lett. 77, 3865 (1996).
 Metzner and Vollhardt (1989) W. Metzner and D. Vollhardt, Phys. Rev. Lett. 62, 324 (1989).
 Georges et al. (1996) A. Georges, G. Kotliar, W. Krauth, and M. J. Rozenberg, Rev. Mod. Phys. 68, 13 (1996).
 Kotliar et al. (2006) G. Kotliar, S. Y. Savrasov, K. Haule, V. S. Oudovenko, O. Parcollet, and C. A. Marianetti, Rev. Mod. Phys. 78, 865 (2006).
 Held (2007) K. Held, Adv. Phys. 56, 829 (2007).
 Leonov et al. (2008) I. Leonov, N. Binggeli, D. Korotin, V. I. Anisimov, N. Stojić, and D. Vollhardt, Phys. Rev. Lett. 101, 096405 (2008).
 Leonov et al. (2010) I. Leonov, D. Korotin, N. Binggeli, V. I. Anisimov, and D. Vollhardt, Phys. Rev. B 81, 075109 (2010).
 Amadon (2012) B. Amadon, Journal of Physics: Condensed Matter 24, 075604 (2012).
 Savrasov and Kotliar (2004) S. Y. Savrasov and G. Kotliar, Phys. Rev. B 69, 245101 (2004).
 Pourovskii et al. (2007) L. V. Pourovskii, B. Amadon, S. Biermann, and A. Georges, Phys. Rev. B 76, 235101 (2007).
 Grechnev et al. (2007) A. Grechnev, I. Di Marco, M. I. Katsnelson, A. I. Lichtenstein, J. Wills, and O. Eriksson, Phys. Rev. B 76, 035107 (2007).
 Di Marco et al. (2009) I. Di Marco, J. Minár, S. Chadov, M. I. Katsnelson, H. Ebert, and A. I. Lichtenstein, Phys. Rev. B 79, 115111 (2009).
 Grånäs et al. (2012) O. Grånäs, I. D. Marco, P. Thunström, L. Nordström, O. Eriksson, T. Björkman, and J. Wills, Computational Materials Science 55, 295â302 (2012).
 Haule et al. (2010) K. Haule, C.H. Yee, and K. Kim, Phys. Rev. B 81, 195107 (2010).
 Gunnarsson et al. (1983) O. Gunnarsson, O. Jepsen, and O. K. Andersen, Phys. Rev. B 27, 7144 (1983).
 Skriver and Rosengaard (1991) H. L. Skriver and N. M. Rosengaard, Phys. Rev. B 43, 9538 (1991).
 Soven (1967) P. Soven, Phys. Rev. 156, 809 (1967).
 Zeller et al. (1982) R. Zeller, J. Deutz, and P. Dederichs, Solid State Communications 44, 993â997 (1982).
 Chioncel et al. (2003) L. Chioncel, L. Vitos, I. A. Abrikosov, J. Kollar, M. I. Katsnelson, and A. I. Lichtenstein, Phys. Rev. B 67, 235106 (2003).
 Minár et al. (2005) J. Minár, L. Chioncel, A. Perlov, H. Ebert, M. I. Katsnelson, and A. I. Lichtenstein, Phys. Rev. B 72, 045125 (2005).
 Leonov et al. (2011) I. Leonov, A. I. Poteryaev, V. I. Anisimov, and D. Vollhardt, Phys. Rev. Lett. 106, 106405 (2011).
 Leonov et al. (2014) I. Leonov, V. I. Anisimov, and D. Vollhardt, Phys. Rev. Lett. 112, 146401 (2014).
 Minár (2011) J. Minár, Journal of Physics: Condensed Matter 23, 253201 (2011).
 Vidberg and Serene (1977) H. J. Vidberg and J. W. Serene, J. Low Temp. Phys. 29, 179 (1977).
 Beach et al. (2000) K. S. D. Beach, R. J. Gooding, and F. Marsiglio, Phys. Rev. B 61, 5147 (2000).
 Wall (1948) H. S. Wall, Analytic Theory of Continued Fractions (Chelsea, New York, 1948).
 Baker (1975) G. A. Baker, Essentials of Padé Approximants (Academic Press, New York, 1975).
 Östlin et al. (2012) A. Östlin, L. Chioncel, and L. Vitos, Phys. Rev. B 86, 235107 (2012).
 Schött et al. (2016) J. Schött, I. L. M. Locht, E. Lundin, O. Grånäs, O. Eriksson, and I. Di Marco, Phys. Rev. B 93, 075104 (2016).
 Nordström et al. (2016) J. Nordström, J. Schött, I. L. Locht, and I. Di Marco, SoftwareX 5, 178â182 (2016).
 Chioncel et al. (2007) L. Chioncel, H. Allmaier, E. Arrigoni, A. Yamasaki, M. Daghofer, M. I. Katsnelson, and A. I. Lichtenstein, Phys. Rev. B 75, 140406 (2007).
 Katsnelson et al. (2008) M. I. Katsnelson, V. Y. Irkhin, L. Chioncel, A. I. Lichtenstein, and R. A. de Groot, Rev. Mod. Phys. 80, 315 (2008).
 Beiuşeanu et al. (2011) F. Beiuşeanu, C. Horea, E.V. Macocian, T. Jurcuţ, L. Vitos, and L. Chioncel, Phys. Rev. B 83, 125107 (2011).
 Chioncel et al. (2015) L. Chioncel, C. Morari, A. Östlin, W. H. Appelt, A. Droghetti, M. M. Radonjić, I. Rungger, L. Vitos, U. Eckern, and A. V. Postnikov, Phys. Rev. B 92, 054431 (2015).
 Andersen et al. (1994) O. K. Andersen, O. Jepsen, and G. Krier, Lectures on Methods of Electronic Structure Calculation (World Scientific, Singapore, 1994).
 Vitos et al. (2000) L. Vitos, H. L. Skriver, B. Johansson, and J. Kollár, Comp. Mat. Sci. 18, 24 (2000).
 Vitos (2001) L. Vitos, Phys. Rev. B 64, 014107 (2001).
 Vitos (2010) L. Vitos, Computational Quantum Mechanics for Materials Engineers (Springer, London, 2010).
 Andersen and SahaDasgupta (2000) O. K. Andersen and T. SahaDasgupta, Phys. Rev. B 62, R16219 (2000).
 Andersen (1975) O. K. Andersen, Phys. Rev. B 12, 3060 (1975).
 Andersen et al. (1986) O. K. Andersen, O. Jepsen, and M. Sob, Electronic Band Structure and Its Applications (Springer Verlag, Berlin, 1986).
 Andersen and Jepsen (1984) O. K. Andersen and O. Jepsen, Phys. Rev. Lett. 53, 2571 (1984).
 Weinberger (1990) P. Weinberger, Electron Scattering Theory for Ordered and Disordered Matter (Clarendon Press, Oxford, 1990).
 Fetter and Walecka (1971) A. L. Fetter and J. D. Walecka, Quantum Theory of ManyParticle Systems (McGrawHill, New York, 1971).
 Nicholson et al. (1994) D. M. C. Nicholson, G. M. Stocks, Y. Wang, W. A. Shelton, Z. Szotek, and W. M. Temmerman, Phys. Rev. B 50, 14686 (1994).
 Pourovskii et al. (2005) L. V. Pourovskii, M. I. Katsnelson, and A. I. Lichtenstein, Phys. Rev. B 72, 115106 (2005).
 Tank and Arcangeli (2000) R. Tank and C. Arcangeli, Phys. Stat. Sol. B 217, 89 (2000).
 Skriver (1984) H. L. Skriver, The LMTO Method (Springer, Berlin, 1984).
 Haule and Birol (2015) K. Haule and T. Birol, Phys. Rev. Lett. 115, 256402 (2015).
 Petukhov et al. (2003) A. G. Petukhov, I. I. Mazin, L. Chioncel, and A. I. Lichtenstein, Phys. Rev. B 67, 153106 (2003).
 Bickers and Scalapino (1989) N. E. Bickers and D. J. Scalapino, Ann. Phys. (N. Y.) 193, 206 (1989).
 Lichtenstein and Katsnelson (1998) A. I. Lichtenstein and M. I. Katsnelson, Phys. Rev. B 57, 6884 (1998).
 Katsnelson and Lichtenstein (1999) M. I. Katsnelson and A. I. Lichtenstein, J. Phys.: Condens. Matter 11, 1037 (1999).
 Imada et al. (1998) M. Imada, A. Fujimori, and Y. Tokura, Rev. Mod. Phys. 70, 1039 (1998).
 Lichtenstein et al. (2001) A. I. Lichtenstein, M. I. Katsnelson, and G. Kotliar, Phys. Rev. Lett. 87, 067205 (2001).
 Liebsch (1979) A. Liebsch, Phys. Rev. Lett. 43, 1431 (1979).
 SánchezBarriga et al. (2009) J. SánchezBarriga, J. Fink, V. Boni, I. Di Marco, J. Braun, J. Minár, A. Varykhalov, O. Rader, V. Bellini, F. Manghi, H. Ebert, M. I. Katsnelson, A. I. Lichtenstein, O. Eriksson, W. Eberhardt, and H. A. Dürr, Phys. Rev. Lett. 103, 267203 (2009).
 Birch (1947) F. Birch, Phys. Rev. 71, 809 (1947).
 Young (1991) D. A. Young, Phase Diagrams of the Elements (University of California Press, Berkeley, 1991).
 Katsnelson and Lichtenstein (2002) M. I. Katsnelson and A. I. Lichtenstein, Eur. Phys. J. B 30, 9 (2002).
 Ceeh et al. (2016) H. A. Ceeh, J. A. Weber, P. Böni, M. Leitner, D. Benea, L. Chioncel, H. Ebert, J. MinÃ¡r, D. Vollhardt, and C. Hugenschmidt, Sci. Rep. 6, 20898 (2016).
 Eriksson et al. (1990) O. Eriksson, B. Johansson, R. C. Albers, A. M. Boring, and M. S. S. Brooks, Phys. Rev. B 42, 2707 (1990).
 Mohn et al. (2001) P. Mohn, C. Persson, P. Blaha, K. Schwarz, P. Novák, and H. Eschrig, Phys. Rev. Lett. 87, 196401 (2001).
 Galler et al. (2015) A. Galler, C. Taranto, M. Wallerberger, M. Kaltak, G. Kresse, G. Sangiovanni, A. Toschi, and K. Held, Phys. Rev. B 92, 205132 (2015).
 Sundararajan et al. (1995) V. Sundararajan, B. R. Sahu, D. G. Kanhere, P. V. Panat, and G. P. Das, Journal of Physics: Condensed Matter 7, 6019 (1995).
 Kulikov et al. (1999) N. I. Kulikov, A. V. Postnikov, G. Borstel, and J. Braun, Phys. Rev. B 59, 6824 (1999).
 Watson and Weinert (1998) R. E. Watson and M. Weinert, Phys. Rev. B 58, 5981 (1998).
 Vitos et al. (1994) L. Vitos, J. Kollár, and H. L. Skriver, Phys. Rev. B 49, 16694 (1994).
 Terletska et al. (2017) H. Terletska, Y. Zhang, L. Chioncel, D. Vollhardt, and M. Jarrell, Phys. Rev. B 95, 134204 (2017).