Ab initio calculation of spin fluctuation spectra using time-dependent density functional perturbation theory, plane waves, and pseudopotentials

We present an implementation of time-dependent density functional perturbation theory for spin fluctuations, based on plane waves and pseudopotentials. We compute the dynamic spin susceptibility self-consistently by solving the time-dependent Sternheimer equation, within the adiabatic local density approximation to the exchange and correlation kernel. We demonstrate our implementation by calculating the spin susceptibility of representative elemental transition metals, namely bcc Fe, fcc Ni, and bcc Cr. The calculated magnon dispersion relations of Fe and Ni are in agreement with previous work. The calculated spin susceptibility of Cr exhibits a soft-paramagnon instability, indicating the tendency of the Cr spins to condense in an incommensurate spin density wave phase, in agreement with experiment.

Figure 1

(a) Calculated $\mathrm{Im}{\chi }^{+-}(\mathbf{q},\omega )$ of Fe along the $\mathrm{\Gamma }\mathrm{N}$ and $\mathrm{\Gamma }\mathrm{H}$ directions. The susceptibility is given in atomic units, $1/\left({\mathrm{Ry}\mathrm{bohr}}^3\right)$. (b) Calculated $\mathrm{Im}{\chi }_{\mathrm{KS}}^{+-}(\mathbf{q},\omega )$ along the $\mathrm{\Gamma }\mathrm{N}$ and $\mathrm{\Gamma }\mathrm{H}$ directions. We also show the magnon dispersion curves superimposed as dots. (c) Magnon bands of Fe along the $\mathrm{\Gamma }\mathrm{N}$ and $\mathrm{\Gamma }\mathrm{H}$ directions. We compare our results to the experimental data of Loong et al. [55], and the calculations of Rousseau et al. [18], and Buczek et al. [15]. (d) Calculated $\mathrm{Im}{\chi }^{+-}(\mathbf{q},\omega )$ of Fe, for selected wave vectors along the $\mathrm{\Gamma }\mathrm{N}$ direction. All peaks are normalized to their maximum for clarity, and the wave vectors are given as $\mathbf{q}=2\pi /a(1,1,0)q$.

Figure 2

(a) Calculated $\mathrm{Im}{\chi }^{+-}(\mathbf{q},\omega )$ of Ni along the $\mathrm{\Gamma }\mathrm{X}$ direction. (b) Calculated $\mathrm{Im}{\chi }_{\mathrm{KS}}^{+-}(\mathbf{q},\omega )$, with the magnon band superimposed as dots. (c) Calculated magnon dispersion curve of Ni along the $\mathrm{\Gamma }\mathrm{X}$ direction. Our results are compared with experimental data of Mook et al. [56, 57], and with the calculations of Rousseau et al. [18], Buczek et al. [15], using TD-DFPT, and Şaşıoğlu et al. [33] using many body perturbation theory (MBPT). (d) Calculated $\mathrm{Im}{\chi }^{+-}(\mathbf{q},\omega )$ of Ni for selected wave vectors along the $\mathrm{\Gamma }\mathrm{X}$ direction. The peaks are normalized to their maximum value, and the wave vector is given as $\mathbf{q}=2\pi /a(1,0,0)q$.

Figure 3

(a) Calculated $\mathrm{Im}\chi (\mathbf{q},\omega )$ of Cr along the $\mathrm{\Gamma }\mathrm{H}$ direction. (b) Calculated $\mathrm{Re}{\chi }_{\mathrm{KS}}(\mathbf{q},\omega )$ of Cr along the $\mathrm{\Gamma }\mathrm{H}$ direction. To highlight the effect of Fermi-surface nesting, we plot $\mathrm{Re}{\chi }_{\mathrm{KS}}(\mathbf{q},\omega )-0.2$ and only show the positive values.

Figure 4

(a) Calculated magnon energies of Ni together with full width at half maximum of corresponding magnon peaks. Wave vectors are given in units of $2\pi /a(1,0,0)q$. (b) Calculated magnon energies of Ni at $\mathbf{q}=(0,0,0)$ with respect to $\mathbf{k}$-points meshes used in the ground state calculations. The $\mathbf{k}$-points meshes are given in $N_k\times N_k\times N_k$. Linear response calculations are all performed with $50\times 50\times 50\mathbf{k}$-points mesh.