Van Vleck excitons in Ca₂RuO₄

corresponding to the individual contributions from the octahedral crystalline electric field $\hat{\mathcal{H}}_{CEF}$, spin-orbit $\hat{\mathcal{H}}_{SO}$, the structural distortion $\hat{\mathcal{H}}_{dis}$ away from ideal octahedral coordination, and a mean molecular field $\hat{\mathcal{H}}_{MF}$ stemming from long range magnetic order. Additional contributions such as hyperfine nuclear splitting exhibit much weaker energy scales ($\sim\mu$eV), and are thus neglected in the current discussion. Chatterji and Schneider [2009] We will now discuss each contributing term to this single-ion Hamiltonian, and thus the terms that ultimately parametrize the uncoupled single site susceptibilities $g(\omega)$.

where $\lambda$ is the spin-orbit coupling constant and expected to increase with atomic number as $\sim Z^{2}$. Landau and Lifshitz [1977] The inclusion of spin-orbit coupling in the magnetic Hamiltonian yields a non-zero $[S_{z},\mathcal{H}]$, and as a result, the expectation value of $S_{z}$ is not conserved. This lack of conservation is in stark contrast with an exclusively Heisenberg magnetic Hamiltonian, and enables the possibility of an amplitude, or longitudinal $zz$, mode to exist. Pekker and Varma [2015] In the case of Ca₂RuO₄, $|\lambda|\sim 0.1$ eV is an order of magnitude smaller than $10Dq$ Das et al. [2018]; Mizokawa et al. [2001]; Fatuzzo et al. [2015], and thus spin-orbit coupling can be considered as a perturbation to the crystal field states defined by $\hat{\mathcal{H}}_{CEF}$ in the $\{|t\rangle,|e\rangle\}$ basis.

consisting of new orbital angular momentum operators that act on the new $|l=1,m_{l}\rangle$ basis, accompanied by a scalar projection factor $\alpha^{\prime}$. In contrast to magnets located in the weak-intermediate field regime, the determination of the scalar $\alpha^{\prime}$ for Ca₂RuO₄ is not particularly straightforward. In the case of the weak crystal field limit, the projection is between bases with good quantum numbers $L$ ($m_{L}$) or $l$ (m_l), and $S$ ($m_{S}$). Possessing fixed orbital and spin angular momentum values, the projection is amenable to methods based on the matrix representation of the angular momentum operators, significantly simplifying the process for determining $\alpha^{\prime}$ as was done for CoO. Sarte et al. [2019] For magnetic ions located in the strong crystal field regime such as is the case for Ca₂RuO₄, the value of $L$ ($m_{L}$) was addressed by Griffith Griffith [1960] using the $T$-$P$ equivalence relation

Valid in the $\Delta\gg\lambda$ regime where minimal mixing between $t$ and $e$ states occur, the projection of the $L=2$ $t$ onto the $L=1$ $p$ states greatly simplifies the process of calculating $\alpha^{\prime}$ directly using representation theory. In the case of the orbital triplet ground state of the $d^{4}$ configuration, Griffith Griffith [1960] and Moffitt Moffitt et al. [1959] determined a value of $-\frac{1}{2}$ for $\alpha^{\prime}$ in the pure $L$-$S$ Russell-Saunders coupling scheme.

confirms a singlet ground state, separated from triplet and pentet excited manifolds by $\Delta$(singlet$\rightarrow$triplet) = $37.5$ meV $\equiv\alpha^{\prime}\lambda$, and $\Delta$(singlet$\rightarrow$pentet) = $112.5$ meV $\equiv 3\alpha^{\prime}\lambda$ (Fig. 1(b)), in agreement with Eq. 15. Confirmation of the assignment of $j\rm{{}_{eff}}=0$, 1, and 2 to the singlet, triplet, and pentet manifolds, respectively, is accomplished by the projection of the components of the effective total angular momentum operator $\hat{\mathbf{j}}=\hat{\mathbf{l}}+\hat{\mathbf{S}}$ onto the subspaces that are spanned by the three individual manifolds defined by $\hat{\mathcal{H}}_{SO}$. Such a projection corresponds to a rotation of the individual angular momentum operators from the original $|l,m_{l},S,m_{S}\rangle$ basis to a basis $|\phi{\rm{{}_{SO}}}\rangle$, consisting of the individual eigenvectors of $\hat{\mathcal{H}}_{SO}$. In the matter of a generic angular momentum operator $\hat{\mathcal{O}}$, this particular rotation is achieved by

corresponding to the matrix multiplication of $\hat{\mathcal{O}}$ by a transformation matrix $\mathcal{C}$ (and its inverse) which consists of the individual eigenvectors $\phi\rm{{}_{SO}}$ that are arranged in order of increasing energy. In the case of $z$-component $\hat{j_{z}}$, the rotation given by Eq. 16 yields

corresponding to a matrix whose top 1 $\times$ 1, middle 3 $\times$ 3, and bottom 5 $\times$ 5 block matrices are identical to the $\hat{J_{z}}$ operator in $|j{\rm{{}_{eff}}}=0,m_{j,{\rm{eff}}}\rangle$, $|j{\rm{{}_{eff}}}=1,m_{j,{\rm{eff}}}\rangle$, and $|j{\rm{{}_{eff}}}=2,m_{j,{\rm{eff}}}\rangle$ bases, respectively. By performing the same projection for the $\hat{j_{x}}$ and $\hat{j_{y}}$ operators, it can be shown that these block matrices satisfy the canonical commutation relations of angular momentum $\hat{\mathbf{j}}=i\hat{\mathbf{j}}\times\hat{\mathbf{j}}$, thus confirming that these block matrices do indeed correspond to valid angular momentum operators.

The Stevens parameter $B^{m}_{l}$ prefactor is a numerical coefficient given by Hutchings [1964]; Bauer and Rotter [2010]

where $\Theta_{l}$ corresponds to projection constants accompanying the conversion from Cartesian coordinates to angular momentum operators via the Wigner-Eckart theorem, and is commonly denoted as $\alpha$ and $\beta$ for $l=2$ and $l=4$, respectively. $\langle r^{l}\rangle$ denotes the expectation values of the radial wavefunction and correspond to 3.319$a^{2}_{0}$ and 20.22$a^{4}_{0}$ for $l=2$ and $l=4$, respectively. Abragam and Bleaney [1986]; Hotta [2006] $p^{m}_{l}$ are the scalar coefficients of the corresponding Tesseral functions $Z^{m}_{l}$. In the point charge approximation, where the charge density $\rho(\mathbf{r})\propto\delta(\mathbf{r})$, the Tesseral functions are incoroporated into the $\gamma^{m}_{l}$ term given by

where $q_{i}$ denotes the $i\rm{{}^{th}}$ charge.

and parametrized by $\Gamma$, whose sign and magnitude are determined by $B^{0}_{2}$.

As illustrated in Fig. 1(b), $\hat{\mathcal{H}}_{MF}$ corresponds to a Zeeman-like term that removes time-reversal symmetry, splitting originally degenerate $j\rm{{}_{eff}}$ levels. In the case where $|J|\rightarrow|\lambda|$ (or equivalently $|H{\rm{{}_{MF}}}|/|\alpha^{\prime}\lambda|\rightarrow 1$), the corresponding increase in splitting induced by $\hat{\mathcal{H}}_{MF}$ results in the significant entanglement between individual $j\rm{{}_{eff}}$ levels.Sarte et al. [2018a, b] The presence of such a strong admixture renders the modeling of the these systems difficult using psuedo-bosonic (Holstein-Primakoff) approaches that are based on conventional linear spin-wave theory, ultimately making it necessary to employ alternative methods such as the multi-level spin-orbit exciton model discussed here.

and

for $J_{i,i+k}(\mathbf{Q})$ in Eq. 12, where labels $\mathbf{d}$ and $\mathbf{s}$ denote different $i\neq i+k$ and same $i=i+k$ site indices, respectively.

and

where $\hbar$ has been set to 1, and thus $\omega$ being relabeled as $E=\hbar\omega$. We note here that we have assumed the different single-ion levels are coupled with the same $J(\mathbf{Q})$. We make this assumption for simplicity and test this below against data in the next section. Solving these four coupled equations yields:

where $G^{-+}(\mathbf{Q},E)$ has the same form as $G^{+-}(\mathbf{Q},E)$ with indices $+$ $\longleftrightarrow$ $-$. Here, the energy $E$ was redefined as $E+i\delta$, where $\delta$ is a positive infinitesimal to ensure analyticity of $g^{\alpha\beta}$ (Eq. 9). Its value was set to 50% of the experimental elastic resolution width (HWHM) on the ARCS time-of-flight neutron spectrometer (SNS, ORNL) set to the experimental parameters employed by Jain et al. Jain et al. [2017]