RKKY interactions of CeB6 based on effective Wannier model

We examine the RKKY interactions of CeB6 between multipole moments based on the effective Wannier model obtained from the bandstructure calculation including 14 Cef orbitals and 60 conduction orbitals of Ce-d, s and B-p, s. By using the f -c mixing matrix elements of the Wannier model together with the conduction band dispersion, the multipole couplings with the RKKY oscillation are obtained for the active moments in Γ8 subspace. Both of the Γ5g quadrupole Oxy and the Γ2u octupole Txyz couplings are largely enhanced with q = (π, π, π) which naturally explains the antiferro-quadrupolar phase of the phase II, and are also enhanced with q = (0, 0, 0) corresponding to the elastic softening of C44. Also the couplings of the Γ5u octupole T β z is quite large for q = (0, 0, π) which is related to the antiferro-octupolar ordering of a possible candidate for the phase IV of CexLa1−xB6.


Introduction
Electronic state of CeB 6 has been one of the central issues in the heavy Fermion systems, since it exhibits a rich phase diagram of the multipole orderings [1] due to the Γ 8 quartet ground state with the degrees of freedom of multipole moments as shown in Table I. Several multipole orderings have been observed in temperature and magnetic field (T, H) phase diagram such as the antifero-quadrupolar (AFQ) ordering of Γ 5g quadrupole moments (O yz , O zx , O xy ) with a critical transition temperature T Q = 3.2 K (phase II) and the antifero-magnetic ordering of Γ 4u magnetic multipoles (σ x , σ y , σ z ) with T N = 2.3 K (phase III). The antifero-octupolar (AFO) ordering of Γ 5u octupoles (T β x , T β y , T β z ) is also discussed as a possible candidate of the phase IV in the La-substitution system of Ce x La 1−x B 6 .
The 4 f electrons of CeB 6 are almost localized from several experiments. The Fermi-surface (FS) has been observed in the de Haas-van Alphen (dHvA) experiments [2], the angle resolved photoemission spectroscopy (ARPES) [3,4] and the high-resolution photoemission tomography [5], where an ellipsoidal FS centered at X point in the Brillouin zone (BZ) has been confirmed and is almost the same as that of LaB 6 with the 4 f 0 state. In such a localized f electron system, Ruderman-Kittel-Kasuya-Yosida (RKKY) interaction [6][7][8] plays an important role for the multipole ordering and must determine the ordering moments and wavevectors. The phonomenological RKKY models of CeB 6 [9,10] succeeded in reproducing the basic phase diagram in T -H plane and giving a great advance in the multipole physics in Γ 8 ground state system. However in these studies, only the nearest neighbor RKKY Hamiltonian with a symmetric coupling and/or asymmetric correction terms was used for explaining the experimental phase diagram. As shown in the original studies [6][7][8], the RKKY interaction has a decaying and oscillating function of 2k F R with the Fermi wavenumber k F of conduction (c) band and distance between multipole moments R, and thus these effects should be taken into account through the microscopic description of the f -c mixing and c band states. The explicit derivation of the RKKY multipole couplings based on the realistic bandstructure calculation still has been an important challenge for the microscopic IRR notation pseudospin rep. multipole , and g (u) means even (odd) time-reversal symmetry. Here all the multipole operators are normalized.
understanding of the multipole order.
In this study, we derive the RKKY interaction of CeB 6 based on the 74-orbital effective Wannier model derived from the bandstructure calculation directly [11]. By using the realistic c band dispersion together with the f -c mixing matrix elements from the Wannier model of CeB 6 , we calculate the RKKY couplings between the active multipole moments in Γ 8 subspace mediated by the realistic c band states explicitly. The obtained RKKY multipole interactions show that the 1st leading multipole mode is the q = (π, π, π)-AFQ ordering with Γ 5g quadrupole O xy together with Γ 2u octupole T xyz and the 2nd leading mode is the q = (0, 0, π)-AFO ordering with Γ 5u octupole T β z .

Model & Formulation
First we perform the bandstructure calculation of CeB 6 and LaB 6 by using the WIEN2k code [12], which is based on the density-functional theory (DFT) and includes the effect of the spin-orbit coupling (SOC) within the second variation approximation. The explicit bandstructures and FSs are shown in Fig.1 (a)-(d) of Ref. [11], where three FSs are obtained in CeB 6 , while an ellipsoidal FS centered at X point is obtained in LaB 6 , which well accounts for the experimental results [2][3][4][5].
Next we construct the 74-orbital effective Wannier model based on the maximally localized Wannier functions (MLWFs) method [13,14] from the DFT bandstructure of CeB 6 , where 14 f -states from Cef (7 orbital × 2 spin) and 60 c-states from Ce-d (5 orbital × 2 spin), Ce-s (1 orbital × 2 spin), B-p (6 site × 3 orbital × 2 spin) and B-s (6 site × 1 orbital × 2 spin) are fully included and the obtained bandstructure well reproduces the DFT band as shown in Fig.1 (e) and (f) of Ref. [11]. The obtained tight-binding (TB) Hamiltonian is given by the following form as, where f † km (c † k ) is a creation operator for a f (c) electron with wavevector k and 14 (60) spin-orbital states m ( ). Here 14 f states of m are represented by the CEF eigenstates as Γ 8 quartet and Γ 7 doublet with the total angular momentum J = 5/2, and Γ 6 , Γ 7 doublets and Γ 8 quartet with ] includes the f (c) energy levels, SOC couplings, CEF splittings and ff (c-c) hopping integrals, and V km is the f -c mixing element.
Here we consider the RKKY interaction between the multipole moments of Γ 8 quartet. For this purpose, we start from the localized f limit where a 4 f 1 state is realized in Γ 8 at each Ce site and the ff hopping and f -c mixing become zero. Hence the remained c electron Hamiltonian H c TB is diagonalized as follows, where a † ks is a creation operator for a electron with k and band-index s with the c band dispersion ε c ks and the eigenvector u c ks where a ks = u c ks c k . Figure 1 shows the c bandstructure [ Fig. 1(a)] and their FSs of the 21st band [ Fig. 1(b) & (c)] and the obtained c state well describes the experimental FSs [2][3][4][5] and is also almost same as the DFT-bandstructure of LaB 6 without f electron.
By using the second-order perturbation w. r. t. the f -c mixing V km in the third term of H TB , we obtain the multi-orbital Kondo lattice Hamiltonian which is given by, where f † im is a creation operator for a f electron with a Ce-atom i and 4-states in Γ 8 quartet |m = |1 ∼ |4 which are given with the J z -base of J = 5/2 |M explicitly as, where ε f Γ 8 is the bare f energy-level and ∆ε f Γ 8 is a energy shift due to the DFT potential which is of the order of a few eV. The Kondo coupling J k ,k imm can be written by the following simple form, where only f 0 -intermediate process is considered and the scattered c orbital energies are fixed to µ. The RKKY Hamiltonian can be obtained from the second-order perturbation w. r. t. the third term of H MKL together with the thermal average for the c states. The final form is given by, where K m 1 m 2 m 3 m 4 (R i j ) is the RKKY coupling between {m 1 , m 2 } at Ce-atom R i and {m 3 , m 4 } at R j and i j represents a summation for Ce-Ce vectors R i j = R i − R j and K m 1 m 2 m 3 m 4 (q) is given by, where f (x) is the Fermi distribution function f (x) = 1/(e (x−µ)/T + 1) and µ is a chemical potential.
Here v mm ks is a f -c mixing matrix between m and m via the c band state with k, s given by, which has all information about the f state scattering between {m, m } through the c state with k, s. The multipole interaction in q-space is explicitly given by the follwing form, where O Γ mm is the matrix element of the multipole operator and K loc m 1 m 2 m 3 m 4 = (1/N) q K m 1 m 2 m 3 m 4 (q). The mean-field multipole susceptibility χ O Γ (q) is written by, which is enhanced towards the multipole ordering instability for the ordering moment O Γ and wavevector q = Q, and diverges at a critical point of the multipole ordering transition temperature In the localized f limit, the Γ 8 ground state is degenerate and the single-site susceptibility for all multipole moments exhibits the Curie law, χ 0 O Γ (q) = 1/T , and then the transition temperature for a certain multipole ordering is determined by the condition . Therefore the sign and maximum value of K O Γ (q) plays an central role for the multipole ordering. Hereafter we set µ − ε f Γ 8 = 2 eV, and µ is determined so as to keep n tot = n c = 21 and T is set to T = 0.005 eV throughout the calculation.

Results
The RKKY couplings K O Γ (q) for several multipole moments as a function of the wavevector q along the high symmetry line in the BZ are plotted as shown in Fig. 2, where the positive (negative) coupling for a certain multipole with (O Γ , q) enhances (suppresses) the corresponding multipole fluctuation and its positive maximum value gives a leading multipole ordering mode. All leading modes are summarized in Table II of Ref. [11]. The couplings of the Γ 5g quadrupole O xy and Γ 2u ocutupole T xyz for q = (π, π, π) become largest among all moments and q, which corresponds to the AFQ ordering of CeB 6 as phase II. From the analysis of the real space couplings of O xy and T xyz , we have found that the main origin of this mode comes from the fact that the couplings with the 1st and 2nd neighbor Ce-Ce vectors exhibit an anti-ferro (AF) and ferro (F) interaction respectively, which indicates the realization of the RKKY oscillation as shown in Fig. 7 in Ref. [11] but the absolute value of the 2nd neighbor coupling is almost same or slightly larger than that of the 1st neighbor coupling. The 2nd neighbor F couplings of O xy and T xyz also increase the uniform mode with a substantial peak for q = (0, 0, 0) as shown in Fig. 2 corresponding to the elastic softening of C 44 [15].
The difference of the couplings between O xy and T xyz has often been discussed in the early studies [16,17], where the two couplings must have the same value within the 1st neighbor coupling due to the point group symmetry. In the present calculation, the same coupling value of the 1st neighbor O xy and T xyz has been obtained, while in the 2nd neighbor couplings, a slightly but finite difference between O xy and T xyz has been observed for the Ce-Ce vectors R = a(011), a(101) where a is the lattice constant, which enhances the peak of O xy over that of T xyz at q = (π, π, π). The anisotropy of the present interaction Hamiltonian including the origin of the above correction together with the so-called 'bond density' [16] will be presented in the subsequent paper [18].

√
3 [19]. The present development of the ζ z and O v X-point mode appears to be related to the recent inelastic neutron scattering experiments in Ce x La 1−x B 6 [20], where the intensity q = (π, 0, 0) is enhanced and becomes dominant mode for x < 0.8. As for the phase III of the AFM order, the Γ 4u magnetic multipole couplings of σ z and η z shall be dominant when the system enters into phase II. They does not become so large in the present paramagnetic system (phase I) and their maximum values are less than half of the 1st leading peak value of Γ 5g -(π, π, π). The realistic description of the successive transition from phase I (paramagnetic) to phase II (AFQ) and from phase II (AFQ) to phase III (AFM) is an important future problem.