Local magnetic moments in iron and nickel at ambient and Earth’s core conditions

Some Bravais lattices have a particular geometry that can slow down the motion of Bloch electrons by pre-localization due to the band-structure properties. Another known source of electronic localization in solids is the Coulomb repulsion in partially filled d or f orbitals, which leads to the formation of local magnetic moments. The combination of these two effects is usually considered of little relevance to strongly correlated materials. Here we show that it represents, instead, the underlying physical mechanism in two of the most important ferromagnets: nickel and iron. In nickel, the van Hove singularity has an unexpected impact on the magnetism. As a result, the electron–electron scattering rate is linear in temperature, in violation of the conventional Landau theory of metals. This is true even at Earth’s core pressures, at which iron is instead a good Fermi liquid. The importance of nickel in models of geomagnetism may have therefore to be reconsidered.


SUPPLEMENTARY NOTE 1. THE COULOMB INTERACTION TENSOR
The general local Coulomb interaction Hamiltonian for a single five band impurity readŝ whered † i ,d i are the creation and annihilation operators for an electron in Wannier orbital φ i (r) and the Coulomb interaction tensor is defined the "Slater" way [1] U ijkl = dr dr φ * i (r)φ * j (r ) We have estimated the Coulomb interaction tensor U ijkl in a fully ab-initio way including the five 3d-and the 4sorbital in the target space, by means of constrained-random phase approximation (cRPA) [2]. To avoid confusion over the definitions and values of the ubiquitous average parameters U and J we discuss them here carefully. One possible definition of the average parameters (for a five band system) is the following Since this is in the spirit of Kanamori [3], we call these Kanamori averages and their values are reported in Supplementary Table 1. Another possible definition of average parameters is rooted in the Slater parametrization of the atomic Coulomb interaction [1,4] U ijkl = 4 α=0 a α (ik; jl)F α (4) and its density-density parts where a(ik; jl) is a shorthand for integrals of three spherical harmonics and can be expressed using Wigner 3j-symbols, and F α are the radial Slater integrals. Averaging over above expressions one obtains the "Slater" averages relating back directly to the Slater integrals. Here l = 2 is the angular momentum of the d shell. The two averages are exactly related in the spherical (atomic) case by which is still very well fulfilled for our cRPA matrix, see Supplementary To better understand how the choice of interaction affects the desctiption of magnetism, we rewrite the atomic Coulomb-interaction Hamiltonian in terms of the total particle number operatorN , total spinŜ and pseudo total angular momentumL. From the expression for density-densitŷ and for KanamoriĤ one can see, that density-density tends to maximizeŜ z , whereas Kanamori tends to maximizeŜ 2 . Thus densitydensity overestimates the magnetisation and Curie-temperature, which has been found in several model calculations [5,6].
As shown in Fig. 1 a and b of the manuscript, the former overestimates T C for both nickel and iron, as a consequence of the well-known absence of spin-flip processes [5,6]. For nickel, both Kanamori and full Coulomb reproduce the experimental result with remarkable accuracy. α-iron (bcc structure) however has lower coordination than nickel and it is closer to half-filling, i.e. charge fluctuations and correlation effects are stronger. For these reasons, sizeable non-local corrections to DMFT are expected in iron [7]. In fact, even with the most sophisticated form of the Coulomb tensor, we find T C about 30% larger than the experimental result [8][9][10].  The Hunds-rule coupling has a great influence on the critical temperatures in both iron and nickel. To explain the higher T C we obtain, as compared to earlier work in Refs. [11][12][13] we performed density-density calculations varying the Hunds coupling. To this end we spherically averaged the Coulomb matrix from cRPA to be able to parametrize it simply via Slater integrals. In this way a change of the Hunds coupling is readily achieved. We employ values ranging from J S = 1 14 (F 2 + F 4 ) = 0.8eV up to the cRPA estimate of 1.0eV for iron and 1.08eV for nickel. We see from Supplementary Table 2 that J has profound influence on the critical temperatures. A reduction of J S to 0.9eV, the value used in the previous study by one of us [11], reduces the critical temperature by a few hundred Kelvins and gets us close to the result published in aforementioned paper. The remaining difference can most probably be traced back to the different value of U and the different DFT implementations and ways to construct the local Hamiltonian.

SUPPLEMENTARY NOTE 3. DOUBLE COUNTING CORRECTION
The Coulomb matrix elements calculated within the cRPA contain naturally the cubic symmetry of the crystal. This means, that unlike in the usual spherically symmetric parametrization of the Coulomb matrix the diagonal entries U iiii are different for the t 2g and e g orbitals. This poses a challenge to the so-called double counting correction of the DFT+DMFT formalism, which is usually estimated as an orbitally averaged Coulomb interaction [14]. The most widely used approaches are the fully-localized-limit (FLL) and around mean-field (AMF) corrections, giving for the double-counting potential µ DC (using the Slater averages) where N is the total number of 3d electrons and n is the occupation per orbital and spin. Since these depend on the filling of the shell N d as well as the average Coulomb interaction U S different corrections are in order when the Coulomb matrix is split into two different blocks.
A generalization of the AMF and FLL can be obtained following Ref. [14]. One arrives at the following orbitally dependent AMF and FLL corrections with the average DFT occupancy n 0 = 1 2(2l+1) i,σ n i,σ . Using either of these criteria gives results for the critical temperatures as well as the orbital occupancies close to the selfconsistency requirement on the charges from Supplementary Eq. 15 below.
Another possible approach is to constraint the total charge on the impurity, which is based on the Friedel sum rule [15]. The Friedel sum rule gives at zero temperature, a relationship between the extra states induced below the Fermi level by a scattering center (an impurity) and the phase shift at the chemical potential. For the Anderson model the extra states induced are given by the occupation number of the impurity states, and the scattering potential is the hybridization that affects the conduction electrons. We thus require the following to be true in self-consistency The value of the imaginary time Green function at τ = β gives the orbital occupancy, thus the trace over them amounts to the total occupancy of the impurity. This works well in metallic systems [16], since in a metal the total particle number of the system and of the impurity are both very sensitive to small variations in the chemical potential µ and the double counting potential µ DC . Also the likeness to the Friedel sum rule, that applies to metals, indicates that such a constraint will work for metals only. Since this constraint affects the chemical potential, µ, as well, both µ and µ DC are adjusted self-consistently to fulfill it.
We have compiled orbital occupancies for the three different double counting corrections for one temperature close to T C for Fe and Ni in Supplementary Table 3.

SUPPLEMENTARY NOTE 4. NI UNDER PRESSURE
In contrast to iron (see e.g. Ref. [17]), the bulk elastic properties of nickel are already well described within DFT [18][19][20]. To account for the effects of pressure in the case of Ni, we thus proceed as follows. First, we recorded energy vs. volume curves for Ni using the VASP code, employing nine different functionals as shown in Supplementary Table  4 and Supplementary Fig. 1. We used a 20 × 20 × 20 point k-mesh centered at the Γ-point and a plane wave cutoff of 600eV. The resulting E(V ) curves are shown in Supplementary Fig. 1 a-c. Subsequently, these were fitted using the following equations of state (EOS): • Vinet [24]  In above equations E 0 is the energy minimum, V 0 the volume at the minimum, B 0 is the bulk modulus and B 0 is its derivative with respect to pressure. The fits to the data are very good in all cases, as shown in Supplementary Fig. 1 d. The resulting fit parameters are compiled in Supplementary Table 4. Subsequently the obtained parameters were used to produce the P (V ) curves for each equation of state: • Vinet One set of such P (V ) curves is shown in Supplementary Fig. 1 Supplementary Fig. 1 f. All other DFT functionals used here fall in between the PW91(FM) and the LSDA+U (FM) curves.

Birch-Murnaghan EOS are shown in
In the large amount of data we produced one can make out different combinations of functional and EOS that fit to the experimentally determined parameters of the textbook case fcc-Ni. It is, however, not our intention here to find the best fit to experiment, but to qualitatively estimate the effects of pressure on the band structure of Ni. Since the data show a relatively small spread we use the average of the curves shown in Supplementary Fig. 1 f as the aggregate "best guess" P (V ) curve. We choose the Birch-Murnaghan EOS, but as shown in Supplementary Fig. 1  Spectral information can be obtained from CT-QMC data via a transformation of the Green's function from imaginary time or Matsubara axis to real energy. Since this transformation is unstable and the QMC data is inherently noisy, one has to resolve to special techniques, such as maximum entropy [25] or stochastic analytical continuation [26,27].
Here we use the stochastic optimization method for analytical continuation developed by Mishchenko et al. [28][29][30] to obtain spectra from the imaginary time Green's function. Since the spectra are not the main focus of this work, we only discuss them here briefly. The results for Fe at T = 386K are shown in Supplementary Fig. 2 for the three parametrizations of the interaction used. The spectra generally agree well with photoemission experiments [31][32][33], showing a principal peak below the Fermi level, a secondary peak at about 2eV binding energy, as well as additional features at higher binding energies with a potential satellite at about 6-7eV. The spectra undergo an evolution as a function of the parametrization of the Coulomb interaction, showing more multiplet features when going form density-density (panels a,b in Supplementary  Fig. 2) towards the full inclusion of U ijkl (e,f in Supplementary Fig. 2). This is somewhat expected, due to the large effects of the Coulomb interaction on the electronic structure of iron in general, as discussed in the manuscript. For Ni the situation is complicated by the long history of the satellites in the photoemission spectrum. Our results for Ni (T = 386K) are shown in Supplementary Fig. 3 for the three parametrizations of the interaction used. Here, the spectra do not vary as strongly as in the case of iron. This is again not unexpected, since the different parametrizations do not have such a strong effect in Ni. Only small changes in the intensities and positions of features are observed.
Comparing to experiments we identify a large principal peak, a sholder next to it, as well as an additional peak at about 2eV binding energy [32,[34][35][36]. Apart from these, we see only one additional satellite around 6 − 8eV binding energy, which we identify as the "6eV valence band satellite". In Supplementary Fig. 4 we show an enlarged version of the satellite region. The whole spectrum below about 4eV binding energy shows stronger majority than minority character, in accordance with experiment [37,38]. We find the relative spin polarization (majority-spin minus minority-spin divided by their sum) to be between 40% (density-density) and 35% (Kanamori, Full) for the satellite depending on the parametrization of the Coulomb interaction. We do not see additional spectral features beyond the energy regions shown. Satellites reported at higher binding energies have been identified in Ref. [34] as potential artifacts of data subtraction procedures. We have also performed a preliminary comparison of our calculations with angle-resolved photoemission spectroscopy (ARPES) data, which is shown in Supplementary Fig. 5. Unlike the local spectral functions shown above, the data were obtained using a faster standard maximum-entropy approach. We compare our results with ARPES of Himpsel et al. [39] for Ni and of Schäfer et al. [40] for Fe. We have used the results of the density-density interaction, since the QMC data at low temperature have the least noise here. For Ni the agreement is quite satisfactory only close to the Fermi level along L-Γ, but worsens at larger binding energies as well as along Γ-X. From the analytically continued Greens function at the L-point we have estimated the exchange splitting to be about 0.25eV in good agreement with experiments that give values between 0.26 ± 0.05eV [41] and 0.31 ± 0.03eV [39]. Since a recent study [7], employing quasiparticle self-consistent GW (QSGW) to account for non-local correlations, found that using QSGW+DMFT leads to an improved description for the whole band structure, it appears that non-local effects play an important role in ferromagnetic Ni. For Fe the agreement between our data and ARPES is very good along the N-Γ-P lines of the band structure. Larger quantitative discrepancies appear along P-H-Γ, which is a trend also seen in the Gutzwiller-DFT of Schickling et al. [42]. Since we can reasonably describe the ARPES spectrum of the ferromagnetic phase of Fe, it appears that local correlations play an important role here. On the other hand it was shown in Ref. [7] that the ARPES spectrum of Fe is also very well described within QSGW alone, neglecting local correlations. In this case, only a simultaneous analysis of one-and two-particle quantities within the same theoretical scheme can clarify the relative role of local (DMFT) and non-local (QSGW) correlations.

SUPPLEMENTARY NOTE 6. CPA+DMFT
The coherent potential approximation (CPA) can be combined with DMFT to take into account substitutional disorder, see e.g. Ref. [43] and references therein. Take a binary alloy of two components A x B 1−x with concentrations x and 1 − x.
If the concentration of A is much larger than that of B, one assumes that the dispersion and lattice structure of A is preserved for the whole structure to a good approximation, and that the substituents B are randomly distributed over the lattice sites of A.  [44,45] to be 7.05Å/atom. The c/a ratio was fixed at 1.6. Since Mn has one 3d-electron fewer than Fe, and Ni has one 3d-electron more than Fe, we obtain shifts of CPA,Ni = −1.67eV for nickel, and CPA,Mn = +0.94eV for manganese. Within the DMFT the impurity problem is subsequently solved for each site A and B, and an average is performed over the Greens functions according to their concentrations This gives a disorder independent mean-field, that closes the DMFT self-sonsistency loop. We have benchmarked our implementation against Ref. [43] and reproduced some of the results reported there. Since the CPA can account only for substitutional disorder, we have performed additional molecular dynamics (MD) simulations combined with DMFT to estimate effects of thermal disorder. Following Pozzo et al. [46] we used a canonical ensemble at 6000K and a run time of 20 ps and a time step of 1 fs. Configurations were written every 50 fs, so in total 20.000 configurations were generated. These calculations were performed with an fcc Ni supercell at 330 GPa containing 27 atoms. Subsequently, for the last six MD configurations we performed an additional static calculation and a local orbital projection as for the bulk. These inputs were used to perform a DMFT calculation for all 27 Ni atoms in the unit cell. We used a temperature of about 6000K (β = 2eV −1 ) also in the DMFT. Since the crystal structure of the MD configurations is not cubic anymore we used the spherically averaged Coulomb interaction matrix for Ni, see Supplementary Note 1, with U S = 2.71eV and J S = 1.0eV. As expected from the structural changes in each configurations the spectrum, self-energies and fillings vary somewhat between atoms. The fillings of the atoms in each studied configuration vary between 8.8 and 8.9 electrons (bulk: 8.85), where the filling is equally distributed between the orbitals. The atom-averaged filling in each studied configuration gives 8.84-8.86, basically on top of the bulk result. Averaging the spectra and self-energies over all 27 impurity problems shows that the cell-averaged results are again very close to the bulk fcc Ni result. In Supplementary Fig. 6 we show the DMFT spectrum and electronic self-energy for the MD configuration at 19.80ps compared to the bulk result obtained using the same parameters. The data shown is a goods representative of all MD+DMFT results we obtained for this system.

SUPPLEMENTARY NOTE 8. NONINTERACTING SUSCEPTIBILITY OF IRON
In Supplementary Fig. 7 we show the noninteracting susceptibilty of iron. It has a Curie-Weiss behaviour coming from the e g orbitals, which have a van-Hove singularity above the Fermi level. Since the DOS is asymmetric around E F and there are minor features, there is also a tiny kink. The t 2g part has no sharp feature around E F and thus is of Pauli form. It is difficult to calculate the scattering rate from data on the imaginary axis. One has either to do an extrapolation of the self-energy or analytic continuation of the data to the real frequency axis. The extrapolation becomes doubtful for large temperatures, since the distance of the first Matsubara frequency iπ/β = iπk B T scales linearly with temperature. Analytic continuation is per construction ill-defined, as discussed above. One can, however, show, that there is a "first-Matsubara rule" [47], which states that for a canonical Fermi-liquid the scaling of the first Matsubara frequency in temperature has to be O(T ). In Supplementary Fig. 8 we show the result of application of this rule to the 3-band Bethe lattice (90% filled) with a bandwidth of 4.4eV. The result is a line through the origin and thus a Fermi-liquid. Nickel at ambient pressure instead begins to strongly deviate from linear behaviour at a few hundred Kelvins, whereas Nickel at 330GPa pressure begins to deviate at 600K from Fermi-liquid behaviour. This confirms our results for the non-Fermi-liquid scattering rates for nickel shown in the manuscript.