Formalism

As noted in equation 4, in the case of intermediate coupling, we calculate the single-ion Hamiltonian considering the electron-electron (Coulomb) interaction and the spin-orbit (SO), which are respectively the first¹⁴ and second lines of the equation, in addition to the crystal field (CF) interaction which cfield treats.

The term intermediate coupling is somewhat a misnomer, because it does not refer to a specific coupling scheme to determine the total angular momentum $J$, such as Russel-Saunders coupling ($LS$-coupling), or $jj$-coupling. Rather, it refers the case where all three of the above single-ion interactions are treated on an equal footing. This is to distinguish between the practice of using only the lowest energy manifold of levels with the same $J$ in either of the other two coupling schemes to calculate the low temperature (low energy) physical properties¹⁵.

Finaly in some cases, intermediate coupling, also refers to the narrower case where although the Coulomb and SO terms are included, their strengths are fixed by parameters determined from experimental (optical) spectra, or ab-initio (Hartree-Fock) calculations, rather than being freely varied.

In the case where the Coulomb interaction dominates, and the spin-orbit interaction is sufficiently stronger than the CF, then only the lowest energy multiplet - that is, the degenerate energy levels produced by the Coulomb and SO interactions - need to be taken into account. This multiplet may be determined from Hund's Rules, and hence we need only calculate the CF interaction (as cfield does). This case is termed Russel-Saunders or $LS$-coupling. Thus Russel-Saunders coupling corresponds to the limit where both the Coulomb and spin-orbit interaction are large (but where $\mathcal{H}_{\mathrm{Coulomb}}\gg \mathcal{H}_{\mathrm{SO}}\gg \mathcal{H}_{\mathrm{CEF}}$), i.e. the crystal field interaction is treated as a perturbation. Only the lowest energy spin-orbit multiplet is considered, and its degeneracy is split by the CF. These levels are labeled by the total angular momentum quantum numbers $L, S, J$ and $J_z$¹⁶. This is the case treated in the modules so1ion/cfield, and is sometimes called the weak field case.

The other limit corresponds to $jj$-coupling, where the spin-orbit interaction dominates $\mathcal{H}_{\mathrm{SO}}\gg \mathcal{H}_{\mathrm{Coulomb}}$, and there is no term or multiplet structure, in this case only the single electron total angular momentum quantum number $j$ serves to label the states. In between these two limits is thus, intermediate coupling.

In addition to the limits represented by these coupling schemes, the crystal field literature also introduces several other limiting cases, depending on the strength of the crystal field. The weak field limit corresponds to $LS$-coupling. The intermediate field is the case where $\mathcal{H}_{\mathrm{Coulomb}}\gg\mathcal{H}_{\mathrm{CEF}}\gg \mathcal{H}_{\mathrm{SO}}$ and this case may be calculated using the icf1ion module, which ignores $\mathcal{H}_{\mathrm{Coulomb}}$, and treats only $\mathcal{H}_{\mathrm{CEF}}$ and $\mathcal{H}_{\mathrm{SO}}$ within the lowest term $^{2S+1} L$ determined by Hund's rules.

The strong field limit is where $\mathcal{H}_{\mathrm{CF}}$ dominates even over the Coulomb interaction. In this case, the crystal field is considered to act first on single-electron states, splitting the degenerate orbital levels $m_l=-l,...,l$, as dictated by the point symmetry of the ion (e.g. for a transition metal ion in cubic symmetry this is usually into a doublet $\vert m_l=0,2\rangle$, labelled by the group theoretical irreducible representation of the cubic group, $E_g$; and a triplet $\vert m_l=-2,-1,1\rangle$ labelled $T_{2g}$). Electrons may then fill these single electron crystal field states, and their degeneracy may be lifted by the Coulomb and spin-orbit interactions. This case is not explicitly considered in the module, but one may calculate the resulting energy levels by increasing the crystal field above the Coulomb interaction. ic1ion however, uses the total orbital and spin angular momenta as labels for the states, rather than single-electron angular momenta or irreducible representations, so it is difficult to intepret ic1ion output in terms of this limit.

Finally, ic1ion may be used in the general case, where all the three main single ion interactions are considered with variable strengths $\mathcal{H}_{\mathrm{Coulomb}}\sim \mathcal{H}_{\mathrm{SO}}\sim \mathcal{H}_{\mathrm{CEF}}$, but may also specifically refer to the situation in which the strengths of the Coulomb and spin-orbit interaction is fixed by parameters determined either from experimental (optical) spectra or ab-initio (Hartree-Fock) calculations. If the CF interaction is relatively weak, there will still be a multiplet structure but the eigenstates of the system is a mixture of the $LS$-basis states (due to the SO interaction), with only $J$ being a good quantum number which can be used to distinguish the states. In the case of a large CF interaction, however, J-mixing can occur, which results in the lost of a distinct multiplet structure.

The operators corresponding to each of the three single-ion interactions are calculated as matrices using $LS$-coupling basis states, which are labelled according to a group theoretical classification. This is to facilitate the use of coefficients of fractional parentage (cfp) in the calculation of the matrix elements. These coefficients allow states of given configuration $nl^\nu$ to be determined by the states of a parent configuration $nl^{\nu-1}$ with one fewer electron. That is, we express the $\nu^{\mathrm{th}}$ electron state $\vert\Omega\rangle$ as a sum, $\sum_{\bar{\Omega},\omega} \langle\bar{\Omega};\omega\vert\Omega\rangle\vert\omega\rangle$, over all possible ways of adding a single electron $\vert\omega\rangle$ to a the $(\nu-1)^{\mathrm{th}}$ state, $\vert\bar{\Omega}\rangle$, where $\langle\bar{\Omega};\omega\vert\Omega\rangle$ is the cfp. We also use the results that a single-particle operator, $\langle\Omega\vert V\vert\Omega'\rangle$, acting on $\nu$ electrons is just $\nu$ times the same operator acting on one electron, $\langle\omega\vert v\vert\omega'\rangle$, to expand the states (and matrix element) in terms of the cfp's,

$$\langle\Omega\vert V\vert\Omega'\rangle = \nu \sum_{\bar{\Omega},... ...omega\vert v\vert\omega'\rangle \langle\bar{\Omega};\omega'\vert\Omega'\rangle$$

Once the cfp are known, only the single electron matrix element $\langle\omega\vert v\vert\omega'\rangle$ is needed to be determined for each particular operator. This is generally a straigtforward task. The cfp are also dependent on the classification of the states for which they are constructed, but the advantage of the use of a classification based on group theory is that we can use Racah's theorem [41]¹⁷, to factorise the cfp into terms involving only pairs of quantum numbers labelling the states, rather than all the quantum numbers together. Because of the properties of the groups from which they derive, most of these factors do not have explicit algebraic formulae¹⁸. However they make a chain calculation of the cfp from a single electron state practicable. These calculations were carried out in the atomic spectroscopy community in the 1950's and '60's and tabulated in the works of Racah [41], Wybourne [42], Allison and McNulty [43]¹⁹, and Nielson and Koster [44].

The LS-coupling basis states are tabulated by Nielsen and Koster [44], and the formulae for the matrix elements are from Racah [45,46,47] and Judd [48].

Finally expressions based on the above equation for the matrix elements of the spin-orbit and crystal field operators were derived by Elliot, Judd and Runciman [49]²⁰. For the Coulomb interaction, we have chosen to use the operators defined by Racah [41]²¹. In the case of the rank-1 spin and orbital angular momentum operators, an explicit algebraic expression exists for the sum over cfp in the above equation [48]²², resulting in the equations quoted by Chan and Lam [51].

The interested reader will find much more details in the books by Judd [48] and Wybourne [52], and for the more mathematically inclined, in the lecture notes of Racah [53]. In addition, this group theory approach was used to determine the neutron cross-sections in a series of papers in the 1970's and 1980's by Balcar and Lovesey, and summarised in their books [30,54].

Finally, the Slater and spin-orbit integrals which parameterise the strength of the Coulomb and SO interactions (these are the radial integrals, whereas what we have calculated using the above techniques is effectively the angular parts of these operators), are tabulated in various books. In particular we have used the experimentally determined values given by Abragam Bleaney [55], Griffith [56] and Chakravarty [57] for the transition metal ions, and by Carnall et al. [58,59] for the rare earths. Other relevant literature are the famous papers by Racah [45,46,47].

We now list formulae for the matrix elements of the principal operators used in ic1ion as they are programmed in:

  $$\mathcal H = \mathcal H_{\mathrm{coulomb}} + \mathcal H_{\mathrm {spinorbit}} + \mathcal H_{\mathrm{crystalfield}} + \mathcal H_{\mathrm{zeman}}$$ (114)

[Coulomb electron electron interaction operator:]

$$\mathcal{H}_{\mathrm{coulomb}} = \sum_{i>j=1}^\nu \frac{e^2}{4\pi... ... \sum_{ij} {\mathbf T}_i^{(k)} \cdot {\mathbf T}_j^{(k)} = \sum_k F^k \hat f_k$$

The denominator $\vert\mathbf{r}_i-\mathbf{r}_j\vert$ was expanded as a multipole series in terms of spherical harmonic functions $Y_{kq}$, which are in turn related to the tensor operators ${\mathbf T}_i^{(k)}$ of the $i^{\mathrm th}$ electron. The radial dependence of the electron wavefunction is then embodied in the Slater integrals, $F^k$, whilst the angular part is given by the Slater operators $\hat f_k$, and $k$ runs through even integers from 0 to $2l$. In the case of $f$-electrons we use Racah's operators, $\hat e_k$ which are linear combinations of $\hat f_k$ which transform according to specific irreducible representations of the group SO(7)$\times G_2$. This means that their matrix elements are more straightforward to calculate. The relation between these operators and $\hat f_k$, and their matrix elements are given by Racah [41]. For $d$- and $p$-electrons the $\hat f_k$ operators are used directly, from tabulated matrix elements given by [44].
[Spin-orbit interaction operator $\mathcal H_{\mathrm {spinorbit}}=\zeta \sum_i ({\mathbf s}_i . {\mathbf l}_i)$:] Eqn 16 from Elliot, Judd and Runciman [49].

[Crystal field interaction:]

The crystal field interaction is given by the electrical potential produced according to Coulombs law by neighbouring charges $q_j$ at positions $\mathbf R_j$ acting on n electrons.

  $$\mathcal H_{\mathrm{crystalfield}}=\frac{1}{4\pi \epsilon_0}\sum_{i=1}^\nu \sum_j \frac{-\vert e\vert q_j}{\vert\mathbf R_j - \mathbf r_i\vert}$$ (115)

The crystal field is rewritten in terms of tesseral harmonics $Z_{kq}(\Omega)$ (see appendix F)

  $$\mathcal H_{\mathrm{crystalfield}}=-\vert e\vert\sum_{i=1}^\nu \s... ...R_j^{k+1}} \sum_{q=-k}^k \frac{1}{2k+1} %%@ Z_{kq}(\Omega_i)Z_{kq}(\Omega_j)$$ (116)

we define crystal field parameters

  $$\gamma_{kq}=\sum_j \frac{q_j}{2k+1}\frac{1}{\epsilon_0R_j^{k+1}}Z_{kq}(\Omega_j)$$ (117)

and rewrite the crystal field as

  $$\mathcal H_{\mathrm{crystalfield}}=-\vert e\vert \sum_{i=1}^\nu \sum_{k=0}^{\infty} \sum_{q=-k}^k r_i^k \gamma_{kq}Z_{kq}(\Omega_i)$$ (118)

Defining rescaled tesseral harmonics $z_{kq}$ as

	$$z_{kq}=\sqrt{\frac{8\pi}{2k+1}}Z_{kq} \dots q\neq 0$$			(119)
	$$z_{k0}=\sqrt{\frac{4\pi}{2k+1}}Z_{k0} \dots q=0$$			(120)

and Wybourne parameters as

	$$L_{kq}=-\vert e\vert\langle r^k \rangle \gamma_{kq}\sqrt{\frac{2k+1}{8\pi}} \dots q\neq 0$$			(121)
	$$L_{k0}=-\vert e\vert\langle r^k \rangle \gamma_{kq}\sqrt{\frac{2k+1}{4\pi}} \dots q= 0$$			(122)

the crystal field may be rewritten as

  $$\mathcal H_{\mathrm{crystalfield}}= \sum_{k=0}^{\inf} \sum_{q=-k}^k L_{kq} \sum_{i=1}^\nu z_{kq}(\Omega_i)$$ (123)

In order to calculate the crystal field Hamiltonian matrix, we replace the tesseral harmonic functions $\sum_i z_{kq}(\Omega_i)$ by the tensor operators $\hat{T}_{kq}$ which transform in the same way as the tesseral harmonics under the action of the rotation group SO(3). That is, the action of the operator $\hat{T}_{kq}$ on a wavefunction $\vert J,m_J\rangle$ is the same as that of the functions $z_{kq}$ on the coordinates $\Omega \equiv (\theta,\phi)$. The $\hat{T}_{kq}$ are expressed in terms of similar operators $\hat{C}_{kq}$ (which transform in the same way as spherical harmonics) as:

  $$\hat{T}_{k0} = \hat{C}_{k0}, \qquad \hat{T}_{k,\pm\vert q\vert... ...C}_{k,-\vert q\vert} \pm (-1)^{\vert q\vert} \hat{C}_{k,\vert q\vert} \right]$$ (124)

and the matrix elements of these operators are given by

  $$\begin{array}{l} \langle \theta J m_J \vert \hat{C}_{kq} \vert... ...\begin{array}{ccc} J' & k & J \\ L & S & L' \end{array} \right\} \end{array}$$ (125)

with the single electron matrix element

  $$\langle l \vert \hat{c}_k \vert l \rangle = (-1)^l (2l+1) \left( \begin{array}{ccc} l & k & l \\ 0 & 0 & 0 \end{array} \right)$$ (126)

[The magnetic moment operators, $L_x L_y L_z S_x S_y S_z$ and $\mathcal H_{\mathrm{zeman}}$ :]

The orbital and spin operators have the same dependence on $m_J$, due to the Wigner-Eckart theorem:

	$$\langle \theta Jm_J\vert(L,S)_{x,y}\vert\theta'J'm_J'\rangle$$	$=$	$$\frac{(-1)^{J-m_J}}{\sqrt{\pm 2}} \left[ \left( \begin{array}{ccc... ...m \left( \begin{array}{ccc} J&-1&J' \\ -m_J&1&m'_J \end{array} \right) \right]$$	(127)
			$$\qquad \times \langle \theta J\vert\vert(L,S)\vert\vert\vert\theta'J'\rangle$$
	$$\langle \theta Jm_J\vert(L,S)_z\vert\theta'J'm_J'\rangle$$	$=$	$$(-1)^{J-m_J} \left( \begin{array}{ccc} J&1&J' \\ -m_J&0&m'_J \end{array} \right) \langle \theta J\vert\vert(L,S)\vert\vert\theta'J'\rangle$$	(128)

where the reduced matrix elements are:

			$$\langle \theta J\vert\vert(L,S)\vert\vert\theta'J'\rangle = \delta_{\theta,\theta'} (-1)^{S+L+(J,J')} \sqrt{(2J+1)(2J'+1)}$$
			$$\;\;\times \sqrt{(L,S)((L,S)+1)(2(L,S)+1)} \left\{ \begin{array}{ccc} J'&1&J \\ (L,S) & (S,L) & (L,S)' \end{array} \right\}$$	(129)

Using this the Zeman interaction can be written as

  $$\mathcal H_{\mathrm{zeman}}= -\mu_B (2\mathbf S^n+\mathbf L^n) {\mathbf H}$$ (130)