The Bethe-Slater curve revisited; new insights from electronic structure theory

The Bethe-Slater (BS) curve describes the relation between the exchange coupling and interatomic distance. Based on a simple argument of orbital overlaps, it successfully predicts the transition from antiferromagnetism to ferromagnetism, when traversing the 3d series. In a previous article [Phys. Rev. Lett. 116, 217202 (2016)] we reported that the dominant nearestneighbour (NN) interaction for 3d metals in the bcc structure indeed follows the BS curve, but the trends through the series showed a richer underlying physics than was initially assumed. The orbital decomposition of the inter-site exchange couplings revealed that various orbitals contribute to the exchange interactions in a highly non-trivial and sometimes competitive way. In this communication we perform a deeper analysis by comparing 3d metals in the bcc and fcc structures. We find that there is no coupling between the E g orbitals of one atom and T 2g orbitals of its NNs, for both cubic phases. We demonstrate that these couplings are forbidden by symmetry and formulate a general rule allowing to predict when a similar situation is going to happen. In γ-Fe, as in α-Fe, we find a strong competition in the symmetry-resolved orbital contributions and analyse the differences between the high-spin and low-spin solutions.


I. VOLUME DEPENDENCE OF THE RKKY OSCILLATIONS IN AFM CR.
Here we present the results obtained for AFM Cr by varying the volume of the unit cell. We have considered two lattice parameters: 5.46 a.u. and 5.67 a.u. The former value corresponds to the experimental one and the latter was used in our previous work 1 due to convergence difficulties in RS-LMTO-ASA. Calculated long-ranged exchange interactions along the NN direction are shown in Fig. S1. As one can see, the RKKY oscillations are damped at FIG. S1. Panel "a": Orbitally-resolved JijR 3 ij along the direction of the NN in AFM Cr for two different choices of the lattice parameter. Panels "b" and "c": The corresponding band dispersions along Γ − R direction in the BZ. Dashed circle underlines the feature which gives rise to the long-ranged magnetic couplings. Note that in this case the BZ is simple cubic, due to an AFM order.
higher volume of the unit cell, but remain well pronounced at lower ones. Quite remarkably, the period of the oscillations looks the same for both volumes. Such a different behaviour is a result of the significant changes in the band structure. It was already shown in e.g. Ref. 2 that the formation of the AFM order in Cr is accompanied by opening of the band gap. In our calculations this gap appears for both values of the lattice parameters. However, an inspection of Fig. S1(a,b) reveals that at a lat =5.46 a.u. the valence band touches the E F , which, gives rise to RKKY oscillatory exchange, according to Eq. (1) in the main text. On the other hand, at higher volumes the band dispersion is reduced and the gap completely opens. Another manifestation of this drastic change is the variation in the sublattice magnetization. For a lat =5.46 a.u. the magnetic moment is about 0.65 µ B , whereas for a lat =5.67 a.u. it surges to 1.53 µ B , which is in line with an increase of the NN J ij value at larger volumes.   dxy dyz dxz d x 2 −y 2 d 3z 2 −r 2 dxy 0.051 0.000 0.000 0.000 0.000 dyz 0.000 0.017 0.000 0.000 0.000 dxz 0.000 0.000 0.017 0.000 0.000 d x 2 −y 2 0.000 0.000 0.000 0.083 0.000 d 3z 2 −r 2 0.000 0.000 0.000 0.000 -0.078 TABLE S4. Orbital-decomposed next NN Jij (in mRy) in bcc Fe, corresponding to the bond vector Rij=(0,0,1)a. dxy dyz dxz d x 2 −y 2 d 3z 2 −r 2 dxy 0.016 0.000 0.000 0.000 0.000 dyz 0.000 0.332 0.000 0.000 0.000 dxz 0.000 0.000 0.332 0.000 0.000 d x 2 −y 2 0.000 0.000 0.000 0.018 0.000 d 3z 2 −r 2 0.000 0.000 0.000 0.000 -0.074

III. SYMMETRY ANALYSIS
Eq. (3) in the main text defines the formula used to calculate the exchange constant, J ij , where one has to take the trace of the expression∆ iĜ ↑ ij∆ jĜ ↓ ji . The∆ i and∆ j matrices are diagonal matrices in orbital space. Hence, it is enough to focus on the Green's functionĜ σ ij , which is a matrix over the d-sector with orbital indices m 1 and m 2 that run over the following orbitals: d xy , d yz , d xz , d x 2 −y 2 and d 3z 2 −r 2 .
We first note that the sites i and j define a "bond" in the system that can be characterized by a vector n ij pointing from site i to j. Then J ij has to belong to the invariant IR of the point group (PG) that consists of the subset of the crystallographic PG, in this case O h , that preserves the bond direction n ij . There are four types of symmetric bond directions: 1) parallel to a cubic axis, n ij 00n , with PG C 4v , 2) parallel to body diagonal, n ij nnn , with PG C 3v , 3) parallel to a side diagonal, n ij nn0 , with PG C 2v , and finally 4) lying in a mirror plane, n ij nnk or n ij nk0 with PG C s (n and k are integers).
In Table S5 the IR's ofĜ σ ij are given as obtained through subductions from the IR of O h for the different possible PG. Since orbitals that belong to the same type of IR generally mix, we can deduce that only in the case of C 4v symmetry the mixed term J Eg−T2g ij in Eq. (6) in the main text is zero. For lower symmetry directions it will be non-zero since there are contributions from the same IR for both the subduction of T 2g and E g . For C 3v and C 4v our finding is consistent with the structure of Tables S1, S2, S3 and S4.
For C 2v , the relevant group of e.g. the NN interactions in fcc lattice, we note that the mixing occurs only between one orbital from E g and one from T 2g IR. If we specifically consider the direction [110], the two orbitals that belong to A 1 representation are d 3z 2 −r 2 and d xy , respectively. Then we can deduce that the corresponding Green's function matrix takes the generic form: with the presence of M A element originating from this mixing, while the off-diagonal T C element is due to the fact that the linear combination d yz + d xz belongs to A 2 representation and d yz − d xz belongs to B 1 . This is consistent with the actual form of the Green's function matrices obtained in the DFT calculations.