DMD - formalism

The central problem in applying the MF-RPA is the calculation of the dynamical susceptibility $\chi_{\delta\beta}^{s''s'}({\mathbf Q},\omega)$ from equation (73). The standard procedure is to substitute equation (74) into (73), which is then solved for each desired value of $\omega$ and ${\bf Q}$ by a matrix inversion. In order to avoid a numerical divergence, it is necessary to add to $\hbar\omega$ a small imaginary constant $\hbar\omega\to\hbar\omega +i\Upsilon$ and insert this into equation (74) leading to a susceptibility which is equivalent to equation (67). This method is inefficient and time demanding, however , because a $3N_b\times 3N_b$ matrix has to be inverted for each $(Q,\omega)$ in the calculation.

In order to minimise the computational effort an algorithm was developed [35], which requires only the solution of a single generalised eigenvalue problem at each scattering vector $\mathbf Q$. This dynamical matrix diagonalisation (DMD) resembles the standard approach to lattice dynamics. This approach is very fast and will allow for more complex systems to be investigated¹⁶.

In the following we describe the DMD for a single excitation $\epsilon_- \rightarrow \epsilon_+$ of each subsystem $s$, i.e. we assume that each subsystem is a two level system with a single transition only. Other transitions (terms in equation (74)) can be considered in the DMD formalism by assigning to each of these transitions an additional value of the index $s$ and increasing the total number of subsystems ($N_b$) correspondingly. This procedure is different from adding other terms (which are present in equation (74)) to the right hand side of equation (76). However, both procedures lead to the same results, this is shown in [21].

For readability it is convenient to adopt the following matrix notation: a $m \times m$ matrix is indicated by a bar on top of the symbol, e.g. $\overline{\chi}^s$ refers to the matrix $\chi_{\alpha\beta}^s$ with $\alpha,\beta=1,\dots,m$. A $N_b \times N_b$ matrix is denoted by a bar below the symbol. Making use of these two conventions the dynamical susceptibility $\chi_{\alpha\beta}^{ss'}({\mathbf Q},\omega)$ can be written as $\overline{\underline{\chi}}({\mathbf Q},\omega)$.

Considering only a single excitation $\epsilon_- \rightarrow \epsilon_+$ in the subsystem susceptibility, equation (74) can be rewritten as

$$\overline{\chi}^{s}(\omega)= \frac{\overline{M}^s}{\Delta^s-\hbar \omega}$$ (73)

with $\Delta^s\equiv\epsilon^s_+-\epsilon^s_-$ and the transition element matrix

$$M^s_{\alpha\beta}=\langle-\vert\hat \mathcal I^s_{\alpha}-... ...angle\hat \mathcal I^s_{\beta}\rangle\vert-\rangle (p_--p_+)$$ (74)

Note for experts on programming single ion modules: that any external single ion module has to provide the matrix $M^s_{\alpha\beta}$ for every transition $\vert-\rangle \rightarrow \vert+\rangle$which is to be taken into consideration in the calculation. If the energy of this transition is zero, i.e. $\Delta^s=0$ (diffuse scattering), the expression (77) would be zero because $(p_--p_+)$ vanishes. In this case the single ion module should calculate $(p_+/kT)$ instead of $(p_--p_+)$.

The $m \times m$ matrices $\overline{M}^s$ may be diagonalised giving eigenvalues which are all zero except for one real eigenvalue $\gamma^s$ (which has the same sign as $\Delta^s$):

$$\gamma^s=\rm Tr \{ \overline{M}^s \}$$ (75)

Now, the MF-RPA problem (72) may be simplified by using the unitary transformation ${\overline{\mathcal U}^s}$ ( ${\overline{\mathcal U}^{s\dag }\overline{\mathcal U}^s}=\overline{1}$), which diagonalises $\overline{M}^s$. Note that the first column of this matrix $\overline{\mathcal U}^s$ (the eigenvector with the eigenvalue $\gamma^s$) is simply

$${\mathcal U^s_{\alpha1}}=\sqrt{(p_--p_+)/\gamma^s}\langle-... ...hat \mathcal I^s_{\alpha}\rangle_{\mathbf % H,T}\vert+\rangle$$ (76)

This property is useful as most of the equations below require only knowledge of this first column. Following the procedure outlined in [35] one may transform the subsystem interaction ${\mathcal J}_{\alpha\beta}^{ss'}({\mathbf Q})$

$$\overline{\mathcal L}^{ss'}(\mathbf Q) \equiv \overline{\... ...ine{\mathcal J}^{ss'}({\mathbf Q}) \overline{\mathcal U}^{s'}$$ (77)

Now the $N_b \times N_b$ Hermitian dynamical matrix may be defined as

$$A^{ss'}({\mathbf Q})=\Lambda^{ss'}\Delta^s-\sqrt{\gamma^s} ... ... L^{ss'}({\mathbf Q}) \left (\sqrt{\gamma^{s'}} \right )^\ast$$ (78)

Figure 13: Illustration of the DMD procedure (see sections 9.3 and 9.4 of text).

The energies of the system may be calculated by solving the following generalised eigenvalue problem:

$$\underline{A}({\mathbf Q})\mathbf t=\hbar\omega \underline{\Lambda} \mathbf t$$ (79)

where the matrix $\underline{\Lambda}$ is defined as

$$\Lambda^{ss'}=\delta_{ss'}{\rm sgn}(\Delta^s)$$ (80)

The solution of the generalised eigenvalue problem (82) yields the eigenvectors $\underline{\mathcal T}=\left( \mathbf t^1, \mathbf t^2, \dots,\mathbf t^r,\dots\right)$ and eigenvalues $\hbar \omega^r$. These may be written as the eigenvalue matrix $\Omega^{rr'}=\delta_{rr'}\hbar \omega^r$, and correspond to the excitation energies of the system at the wavevector ${\bf Q}$ for which $\underline{A}({\mathbf Q})$ was calculated. The solution of the eigenvalue problem (82) corresponds to the diagonalisation of the dynamical matrix in the case of phonons and therefore this method for calculating magnetic excitations is called dynamical matrix diagonalisation (DMD). Note, that the matrix $\underline{A}({\mathbf Q})$ does not change when the number of dimensions $m$ of the subsystem susceptibility is increased (for example to include quadrupolar degrees of freedom), unless there is an interaction coupling these degrees of freedom between different ions (this can be seen by the definition (80) of the matrix $\underline{\mathcal L}$ having in mind the first column of the transformation matrices $\overline{\mathcal U}^s$, see equation (79)).

The eigenvector matrix $\underline{\mathcal T}$ provides a unitary transformation, which may be used to obtain the dynamical susceptibility. If the eigenvectors are normalised as $\underline{\mathcal T}^\dag\underline{A}\underline{\mathcal T}=\underline{1}$, then equation (73) may be transformed using $\overline{\mathcal U}^s$ and $\underline{\mathcal T}$(see [21]):

$$\chi^{ss'}_{\alpha\beta}(\mathbf Q,\omega)=\left (\sqrt{\g... ...}(\mathbf Q)} \mathcal U^{s'\dag }_{1\beta}\sqrt{\gamma^{s'}}$$ (81)

By its definition the generalised susceptibility gives information about the correlated movement of the operators $\langle \hat \mathcal I_{\alpha}^s \rangle$ for a specific excitation and contains the relative phases and amplitudes of the different operators. The procedure for the calculation of excitation energies $\hbar \omega^r$ and physical observables (such as correlation functions and spectra) outlined above is very fast, because it involves only a single diagonalisation (determination of the matrix $\underline{\mathcal T}$) for every scattering vector of interest. The dynamical susceptibility does not need to be calculated for each energy transfer by inverting equation (72) saving much computation time. Therefore the module McDisp of the McPhase package uses this method by default. We want to emphasize, that the procedure outlined in this section is general and allows to treat any number and combination of multipolar interactions just by letting the index $\alpha$ in $\hat \mathcal I_{\alpha}^n$ take values between 1 and the number of multipolar operators considered ( $\alpha=1,\dots,m$).

Figure 13 illustrates the DMD algorithm. We have described the first three parts shown in the figure, obtaining the eigenvectors $\mathbf t^{r}$ and eigenvalues $\omega^r$. The following parts are described in the next section where the dynamical susceptibility $\overline{\underline{\chi}}(\mathbf Q,\omega)$ is used to calculate a general susceptibility $\overline{\underline{X}}(\mathbf Q,\omega)$ corresponding to an arbitrary observable, from which physical properties may be calculated.