# Spin-orbital correlations from complex orbital order in MgV2O4

H. Lane, P. M. Sarte, K. Guratinder, Angel Arevalo Lopez, R. S. Perry, E. C.

Hunter, T. Weber, B. Roessli, A. Stunault, Y. Su, et al.

► To cite this version:

H. Lane, P. M. Sarte, K. Guratinder, Angel Arevalo Lopez, R. S. Perry, et al. Spin-orbital correlations from complex orbital order in MgV2O4. Physical Review Research, 2023, Physical Review Research, 5 (4), 10.1103/PhysRevResearch.5.043146. [https://doi.org/10.1103/PhysRevResearch.5.043146](https://doi.org/10.1103/PhysRevResearch.5.043146)

HAL Id: hal-04443893

https://hal.univ-lille.fr/hal-04443893

Submitted on 7 Feb 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Spin-orbital correlations from complex orbital order in MgV2O4

H. Lane, P. M. Sarte, K. Guratinder, Angel Arevalo Lopez, R. S. Perry, E. C.

Hunter, T. Weber, B. Roessli, A. Stunault, Y. Su, et al.

► To cite this version:

H. Lane, P. M. Sarte, K. Guratinder, Angel Arevalo Lopez, R. S. Perry, et al. Spin-orbital correlations from complex orbital order in MgV2O4. Physical Review Research, 2023, Physical Review Research, 5 (4), 10.1103/PhysRevResearch.5.043146. [https://doi.org/10.1103/PhysRevResearch.5.043146](https://doi.org/10.1103/PhysRevResearch.5.043146)

HAL Id: hal-04443893

https://lilloa.hal.science/hal-04443893v1

Submitted on 7 Feb 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
Spin-orbital correlations from complex orbital order in MgV$_{2}$O$_{4}$

H. Lane,¹,²,³,⁴ P. M. Sarte,⁵ K. Guratinder,¹ A. M. Arevalo-Lopez,⁶ R. S. Perry,³,⁷ E. C. Hunter,¹ T. Weber,⁸ B. Roessli,⁹

A. Stunault,⁸ Y. Su,¹⁰ R. A. Ewings,³ S. D. Wilson,¹¹ P. Böni,¹² J. P. Attfield,² and C. Stock¹

¹School of Physics and Astronomy, University of Edinburgh, Edinburgh EH9 3JZ, United Kingdom

²School of Chemistry and Centre for Science at Extreme Conditions, University of Edinburgh, Edinburgh EH9 3FJ, United Kingdom

³ISIS Pulsed Neutron and Muon Source, STFC Rutherford Appleton Laboratory, Harwell Campus, Didcot, Oxon OX11 0QX, United Kingdom

⁴School of Physics, Georgia Institute of Technology, Atlanta, Georgia 30332, USA

⁵Materials Department, University of California, Santa Barbara, California 93106, USA

⁶Univ. Lille, CNRS, Centrale Lille, Univ. Artois, UMR 8181 - UCCS - Unite de Catalyse et Chimie du Solide, F-59000 Lille, France

⁷London Centre for Nanotechnology and Department of Physics and Astronomy, University College London, London WC1E 6BT, United Kingdom

⁸Institute Laue-Langevin, 71 Avenue des Martyrs, CS 20156, 38042 Grenoble cedex 9, France

⁹Laboratory for Neutron Scattering and Imaging, Paul Scherrer Institut (PSI), 5232 Villigen PSI, Switzerland

¹⁰Jülich Centre for Neutron Science (JCNS) at Heinz Maier-Leibnitz Zentrum (MLZ), Forschungszentrum Jülich GmbH, Lichtenbergstraße 1, 85748 Garching, Germany

¹¹Materials Department and California Nanosystems Institute, University of California, Santa Barbara, California 93106, USA

¹²Physik-Department E21, Technische Universität München (TUM), James-Franck-Str. 1, 85748 Garching, Germany

<figure><img src="image_2.png" /></figure>

(Received 9 May 2023; accepted 22 September 2023; published 13 November 2023)

MgV$_{2}$O$_{4}$ is a spinel based on magnetic V$^{3+}$ ions, which host both spin ($S$ = 1) and orbital ($l_{eff}$ = 1) moments. Owing to the underlying pyrochlore coordination of the magnetic sites, the spins in MgV$_{2}$O$_{4}$ only antiferromagnetically order once the frustrating interactions imposed by the $Fd\bar{3}m$ lattice are broken through an orbitally-driven structural distortion at $T_s$ ≃ 60 K. Consequently, a Néel transition occurs at $T_N$ ≃ 40 K. Low-temperature spatial ordering of the electronic orbitals is fundamental to both the structural and magnetic properties; however, considerable discussion on whether it can be described by complex or real orbital ordering is ambiguous. We apply neutron spectroscopy to resolve the nature of the orbital ground state and characterize hysteretic spin-orbital correlations using x-ray and neutron diffraction. Neutron spectroscopy finds multiple excitation bands and we parametrize these in terms of a multilevel (or excitonic) theory based on the orbitally degenerate ground state. Meaningful for the orbital ground state, we report an “optical-like” mode at high energies that we attribute to a crystal-field-like excitation from the spin-orbital $j_{eff}$ = 2 ground-state manifold to an excited $j_{eff}$ = 1 energy level. We parametrize the magnetic excitations in terms of a Hamiltonian with spin-orbit coupling and local crystalline electric field distortions resulting from deviations from perfect octahedra surrounding the V$^{3+}$ ions. We suggest that this provides compelling evidence for complex orbital order in MgV$_{2}$O$_{4}$. We then apply the consequences of this model to understand hysteretic effects in the magnetic diffuse scattering where we propose that MgV$_{2}$O$_{4}$ displays a high-temperature orbital memory of the low-temperature spin order.

DOI: 10.1103/PhysRevResearch.5.043146

# I. INTRODUCTION

Atomic orbitals of magnetic ions provide a link between the crystallographic structure and local magnetic moments, and hence an avenue to couple structural and magnetic degrees of freedom [1–7]. A central parameter to controlling orbital order in materials is spin-orbit coupling, which exists in magnetic ions with an inherit single-ion orbital degeneracy.

Given that spin-orbit coupling in multielectron atoms scales with the atomic number squared ($\lambda \propto Z^2$) [8], first-row transition metals provide an opportunity to study magnetism when spin-orbit coupling is of a similar energy scale to the magnetic exchange and also energy scales of local structural distortions away from a perfect, for example, octahedral environment. These energy scales have consequences to the nature of the real-space magnetic orbitals, which spatially order at low temperatures, breaking any degeneracy. In particular is the situation of complex orbitals and real orbitals illustrated in Fig. 1. In this paper we investigate the interplay between spin-orbital physics in MgV$_{2}$O$_{4}$ and apply neutron spectroscopy to obtain information on the underlying orbital order pointing to ordering of complex orbitals.

<footnote>Published by the American Physical Society under the terms of the Creative Commons Attribution 4.0 International license. Further distribution of this work must maintain attribution to the author(s) and the published article's title, journal citation, and DOI.</footnote>
Spin-orbital correlations from complex orbital order in MgV$_{2}$O$_{4}$

H. Lane,<sup>1,2,3,4</sup> P. M. Sarte,<sup>5</sup> K. Guratinder,<sup>1</sup> A. M. Arevalo-Lopez,<sup>6</sup> R. S. Perry,<sup>3,7</sup> E. C. Hunter,<sup>1</sup> T. Weber,<sup>8</sup> B. Roessli,<sup>9</sup>

A. Stunault,<sup>8</sup> Y. Su,<sup>10</sup> R. A. Ewings,<sup>3</sup> S. D. Wilson,<sup>11</sup> P. Böni,<sup>12</sup> J. P. Attfield,<sup>2</sup> and C. Stock<sup>1</sup>

<sup>1</sup>School of Physics and Astronomy, University of Edinburgh, Edinburgh EH9 3JZ, United Kingdom

<sup>2</sup>School of Chemistry and Centre for Science at Extreme Conditions, University of Edinburgh, Edinburgh EH9 3FJ, United Kingdom

<sup>3</sup>ISIS Pulsed Neutron and Muon Source, STFC Rutherford Appleton Laboratory, Harwell Campus, Didcot, Oxon OX11 0QX, United Kingdom

<sup>4</sup>School of Physics, Georgia Institute of Technology, Atlanta, Georgia 30332, USA

<sup>5</sup>Materials Department, University of California, Santa Barbara, California 93106, USA

<sup>6</sup>Univ. Lille, CNRS, Centrale Lille, Univ. Artois, UMR 8181 - UCCS - Unite de Catalyse et Chimie du Solide, F-59000 Lille, France

<sup>7</sup>London Centre for Nanotechnology and Department of Physics and Astronomy, University College London, London WC1E 6BT, United Kingdom

<sup>8</sup>Institute Laue-Langevin, 71 Avenue des Martyrs, CS 20156, 38042 Grenoble cedex 9, France

<sup>9</sup>Laboratory for Neutron Scattering and Imaging, Paul Scherrer Institut (PSI), 5232 Villigen PSI, Switzerland

<sup>10</sup>Jülich Centre for Neutron Science (JCNS) at Heinz Maier-Leibnitz Zentrum (MLZ), Forschungszentrum Jülich GmbH, Lichtenbergstraße 1, 85748 Garching, Germany

<sup>11</sup>Materials Department and California Nanosystems Institute, University of California, Santa Barbara, California 93106, USA

<sup>12</sup>Physik-Department E21, Technische Universität München (TUM), James-Franck-Str. 1, 85748 Garching, Germany

<figure><img src="image_2.png" /></figure>

(Received 9 May 2023; accepted 22 September 2023; published 13 November 2023)

MgV$_{2}$O$_{4}$ is a spinel based on magnetic V$^{3+}$ ions, which host both spin ($S$ = 1) and orbital ($l_{eff}$ = 1) moments. Owing to the underlying pyrochlore coordination of the magnetic sites, the spins in MgV$_{2}$O$_{4}$ only antiferromagnetically order once the frustrating interactions imposed by the $Fd\bar{3}m$ lattice are broken through an orbitally-driven structural distortion at $T_s$ ≃ 60 K. Consequently, a Néel transition occurs at $T_N$ ≃ 40 K. Low-temperature spatial ordering of the electronic orbitals is fundamental to both the structural and magnetic properties; however, considerable discussion on whether it can be described by complex or real orbital ordering is ambiguous. We apply neutron spectroscopy to resolve the nature of the orbital ground state and characterize hysteretic spin-orbital correlations using x-ray and neutron diffraction. Neutron spectroscopy finds multiple excitation bands and we parametrize these in terms of a multilevel (or excitonic) theory based on the orbitally degenerate ground state. Meaningful for the orbital ground state, we report an “optical-like” mode at high energies that we attribute to a crystal-field-like excitation from the spin-orbital $j_{eff}$ = 2 ground-state manifold to an excited $j_{eff}$ = 1 energy level. We parametrize the magnetic excitations in terms of a Hamiltonian with spin-orbit coupling and local crystalline electric field distortions resulting from deviations from perfect octahedra surrounding the V$^{3+}$ ions. We suggest that this provides compelling evidence for complex orbital order in MgV$_{2}$O$_{4}$. We then apply the consequences of this model to understand hysteretic effects in the magnetic diffuse scattering where we propose that MgV$_{2}$O$_{4}$ displays a high-temperature orbital memory of the low-temperature spin order.

DOI: 10.1103/PhysRevResearch.5.043146

## I. INTRODUCTION

Atomic orbitals of magnetic ions provide a link between the crystallographic structure and local magnetic moments, and hence an avenue to couple structural and magnetic degrees of freedom [1–7]. A central parameter to controlling orbital order in materials is spin-orbit coupling, which exists in magnetic ions with an inherit single-ion orbital degeneracy.

Given that spin-orbit coupling in multielectron atoms scales with the atomic number squared ($\lambda \propto Z^2$) [8], first-row transition metals provide an opportunity to study magnetism when spin-orbit coupling is of a similar energy scale to the magnetic exchange and also energy scales of local structural distortions away from a perfect, for example, octahedral environment. These energy scales have consequences to the nature of the real-space magnetic orbitals, which spatially order at low temperatures, breaking any degeneracy. In particular is the situation of complex orbitals and real orbitals illustrated in Fig. 1. In this paper we investigate the interplay between spin-orbital physics in MgV$_{2}$O$_{4}$ and apply neutron spectroscopy to obtain information on the underlying orbital order pointing to ordering of complex orbitals.

<footnote>Published by the American Physical Society under the terms of the Creative Commons Attribution 4.0 International license. Further distribution of this work must maintain attribution to the author(s) and the published article's title, journal citation, and DOI.</footnote>
<figure><img src="image_3.png" /><figcaption>FIG. 1. Graphical representation of the real $|d_{αβ}|$ and imaginary $|\hat{l}_z = \pm 1\rangle$ orbital wavefunctions. The surface represents the absolute magnitude of the wavefunction whereas the color represents the phase.</figcaption></figure>

Consisting of orbitally degenerate $d^2$ trivalent $V^{3+}$ ($S = 1, l_{eff} = 1$) decorating a pyrochlore [9] sublattice of corner-sharing tetrahedra [Fig. 2(a)], the spinel vanadates $AV_2O_4$ with $A$ = divalent cation, provide a platform for the interplay of the underlying geometric constraints of a frustrated lattice [10] populated with strongly correlated electrons, each with charge, spin, and orbital degrees of freedom [1,11]. The structural, electronic, and magnetic properties of the $AV_2O_4$

series are based on two broadly different classes with the A site either being magnetic [as is the case for (Mn,Fe,Co)V$_2$O$_4$] or not [for example in (Mg,Zn)V$_2$O$_4$] [12]. These two classes of materials display some notable differences with the size of the ordered moment on the $V^{3+}$ reportedly stronger when the A site of the spinel structure is magnetic [13]. Also, magnetic A-site vanadate spinels display noncollinear magnetic order in contrast to their nonmagnetic counterparts [14]. For simplicity towards investigating the orbital ground state of the $V^{3+}$ ion in this class of compounds, we discuss in this paper the case where the A site is nonmagnetic. Given the orbital degeneracy on the $V^{3+}$ site, (Mg,Zn)V$_2$O$_4$ compounds exhibit a Jahn-Teller distortion with a contraction along the c axis accompanying a cubic to tetragonal structural phase transition at $T_S$. This distortion relieves the underlying geometric magnetic frustration, and hence is followed by the formation of long-range collinear antiferromagnetic order at $T_N < T_S$. The ordered moments in the Néel phase are reduced from expectations of a spin-only moment of $gS = 2 \mu_B$ with values of, for example, $\approx 0.5 \mu_B$ reported for MgV$_2$O$_4$ [15].

The specific case of MgV$_2$O$_4$ exhibits an orbitally driven Jahn-Teller structural transition at $T_S \approx 60$ K and a Néel transition at $T_N \approx 40$ K. The exact nature of the low-temperature orbital order is ambiguous because two types of low-temperature orbital order have been proposed for this compound. The first, termed real orbital ordering (ROO), corresponds to orbital ordering of the real $t_g$ and $e_g$ orbitals, which are eigenstates in the limit that the crystal field imposed on the

<figure><img src="image_4.png" /><figcaption>FIG. 2. (a) The low-temperature nuclear and magnetic structure of MgV$_2$O$_4$ obtained from our neutron diffraction results. (b) Our coaligned single crystals. (c) Comparison of powder x-ray and neutron results illustrating the structural and magnetic transitions. The shaded region indicates where both tetragonal and cubic phases coexist. (d) Single-crystal neutron diffraction refinements of the nuclear structure in the high-temperature $T = 140$ K cubic phase and low-temperature $T = 50$ K tetragonal phases. The calculated Bragg peak structure factor squared ($F_{calc}^2$) is plotted against the measured structure factor squared ($F_{obs}^2$) as measured on D9 (ILL).</figcaption></figure>
<figure><img src="image_3.png" /><figcaption>FIG. 1. Graphical representation of the real $|d_{αβ}|$ and imaginary $|\hat{l}_z = \pm 1\rangle$ orbital wavefunctions. The surface represents the absolute magnitude of the wavefunction whereas the color represents the phase.</figcaption></figure>

Consisting of orbitally degenerate $d^2$ trivalent $V^{3+}$ ($S = 1, l_{eff} = 1$) decorating a pyrochlore [9] sublattice of corner-sharing tetrahedra [Fig. 2(a)], the spinel vanadates $AV_2O_4$ with $A$ = divalent cation, provide a platform for the interplay of the underlying geometric constraints of a frustrated lattice [10] populated with strongly correlated electrons, each with charge, spin, and orbital degrees of freedom [1,11]. The structural, electronic, and magnetic properties of the $AV_2O_4$

series are based on two broadly different classes with the A site either being magnetic [as is the case for (Mn,Fe,Co)V$_2$O$_4$] or not [for example in (Mg,Zn)V$_2$O$_4$] [12]. These two classes of materials display some notable differences with the size of the ordered moment on the $V^{3+}$ reportedly stronger when the A site of the spinel structure is magnetic [13]. Also, magnetic A-site vanadate spinels display noncollinear magnetic order in contrast to their nonmagnetic counterparts [14]. For simplicity towards investigating the orbital ground state of the $V^{3+}$ ion in this class of compounds, we discuss in this paper the case where the A site is nonmagnetic. Given the orbital degeneracy on the $V^{3+}$ site, (Mg,Zn)V$_2$O$_4$ compounds exhibit a Jahn-Teller distortion with a contraction along the c axis accompanying a cubic to tetragonal structural phase transition at $T_S$. This distortion relieves the underlying geometric magnetic frustration, and hence is followed by the formation of long-range collinear antiferromagnetic order at $T_N < T_S$. The ordered moments in the Néel phase are reduced from expectations of a spin-only moment of $gS = 2 \mu_B$ with values of, for example, ≃ 0.5 $\mu_B$ reported for MgV$_2$O$_4$ [15].

The specific case of MgV$_2$O$_4$ exhibits an orbitally driven Jahn-Teller structural transition at $T_S \approx 60$ K and a Néel transition at $T_N \approx 40$ K. The exact nature of the low-temperature orbital order is ambiguous because two types of low-temperature orbital order have been proposed for this compound. The first, termed real orbital ordering (ROO), corresponds to orbital ordering of the real $t_g$ and $e_g$ orbitals, which are eigenstates in the limit that the crystal field imposed on the

<figure><img src="image_4.png" /><figcaption>FIG. 2. (a) The low-temperature nuclear and magnetic structure of MgV$_2$O$_4$ obtained from our neutron diffraction results. (b) Our coaligned single crystals. (c) Comparison of powder x-ray and neutron results illustrating the structural and magnetic transitions. The shaded region indicates where both tetragonal and cubic phases coexist. (d) Single-crystal neutron diffraction refinements of the nuclear structure in the high-temperature $T = 140$ K cubic phase and low-temperature $T = 50$ K tetragonal phases. The calculated Bragg peak structure factor squared ($F_{calc}^2$) is plotted against the measured structure factor squared ($F_{obs}^2$) as measured on D9 (ILL).</figcaption></figure>
$V^{3+}$ ion from the surrounding oxygen octahedra is large. This model has been advocated based on x-ray and electron beam diffraction data [16]. An alternate suggestion [17] has been made for ordering of the complex basis of the $\hat{L}_z$ observable, which are complex linear combinations of the real $d$ orbitals. This is denoted as complex orbital ordering (COO) and is the basis state used in the weak-intermediate crystal-field limit for transition metal ions. Based on powder neutron spectroscopy and diffraction, it has also been suggested that the orbital ordering could possibly be intermediate between the two ROO and COO extremes [15]. A graphical representation of the two different extremes is illustrated in Fig. 1.

A theoretical study outlined in Ref. [18] suggested neutron spectroscopy as a means of distinguishing between COO and ROO orders at low temperatures. In particular the study noted the importance of spin-orbit coupling with COO giving multiple magnetic branches in the neutron response and also a larger magnetic zone-center gap than would be expected in the ROO model. Motivated by this, we investigate the magnetic neutron response in single crystals of $MgV_2O_4$. We study the hysteretic critical dynamics using diffraction, measure the dynamic response, and develop a theory to describe the magnetic excitations from the spin-orbital ground state. Given the spatially diffuse and correlated nature of the spin-orbital states used for our theory, we term this an excitonic approach.

We apply this excitonic approach to model our neutron spectroscopy results using the single-ion states and then treating the single-ion Hamiltonian (including spin-orbit coupling) and exchange energetics equally [19–22]. This model combined with the data favors the ordering of complex orbitals (COO) in $MgV_2O_4$ and illustrates the importance of using a complex orbital basis for understanding the properties of the first-row transition metal ions.

This paper is divided into four sections. First, we outline the materials preparation methodology and experimental techniques used to probe the magnetic fluctuations and critical scattering. Second, we outline the spectroscopic experimental results of both the low-energy spin-wave excitations and a higher-energy optic-like mode. Third, we present a theory for the magnetic excitations including spin-orbit coupling and the orbitally degenerate ground state. Finally, we investigate the hysteretic effects resulting from the orbital degeneracy using energy integrated neutron diffuse scattering.

## II. EXPERIMENTAL INFORMATION

*Materials preparation.* Single crystals [Fig. 2(b)] of $MgV_2O_4$ were grown using the floating zone technique with details provided in the Appendix. Given the extreme sensitivity of the magnetic properties to stoichiometry and chemical order in $MgV_2O_4$ [23], we have characterized our single crystals using both thermodynamic and scattering probes with neutrons and x-rays. The diffraction results are summarized in Fig. 2(c), with high-temperature structural ($T_s \simeq 60$ K) and magnetic [$T_N \simeq 40$ K from the magnetic $\vec{Q} = (1,1,0)$ Bragg peak] transitions consistent with the published literature [15,24].

*Characterization.* Figure 2(c) illustrates powder diffraction measurements using synchrotron x-rays (Diamond, I11) and neutrons (PSI, RITA-II). Powder synchrotron data shows a

first-order transition at $T_s \simeq 60$ K from a high-temperature cubic to a low-temperature tetragonal phase. This is confirmed on single-crystal neutron diffraction data following the $\vec{Q} = (0,0,2)$ Bragg peak on warming and cooling. A large hysteresis in the neutron intensity is observed between warming and cooling over the same temperature range where synchrotron x-rays measure a coexistence of tetragonal and cubic phases. This is further discussed in the Appendix where this region was found to extend from ~55–60 K. We have further confirmed that the structural properties of single crystals are consistent with the published structural phases through single-crystal neutron diffraction data (D9, ILL). The refinement results are summarized in Fig. 2(d) in both the cubic ($T = 140$ K) and tetragonal ($T = 50$ K) phases with the structural parameters listed in the Appendix. Figure 2(d) illustrates the calculated Bragg peak intensities ($\propto F^2$) as a function of the measured intensities with a straight line indicative that the models describe the results well. The slope of the line is a calibration factor scaling experiment and calculated values. We have confirmed the magnetic structure to be consistent with that outlined in the literature for $ZnV_2O_4$ [25] using polarized neutrons (DNS, MLZ). The results of this are schematically shown in Fig. 2(a) with the magnetic $V^{3+}$ magnetic moments pointing along the tetragonal c axis.

*Critical magnetic fluctuations.* In this paper we discuss energy-integrated polarized diffuse scattering measurements sensitive to the magnetic critical scattering and neutron inelastic data probing the spin-orbital dynamics. To investigate the magnetic critical scattering as a function of temperature, we studied the polarized neutron cross section using an xyz geometry provided by the DNS diffractometer (MLZ, Munich) [26]. Further details of the polarized beam measurements and analysis are provided in the Appendix.

*Neutron spectroscopy.* To probe the magnetic dynamics sensitive to the spin-orbital ground state, neutron spectroscopy was performed on four different spectrometers. To study the high-energy excitations, which represent spin-orbital excitations we used the MAPS spectrometer. The sample was oriented so that $\vec{k}_i$ was oriented along the $c^*$ axis and the intensity along $L$ integrated as typically done with time of flight direct geometry spectrometers to measure one- or two-dimensional spin excitations. To obtain an overview of the low-energy and low-dimensional dispersive dynamics, the multi rep rate option on MERLIN [27] was exploited combined with rotating the crystal. Further data was taken on the EIGER triple-axis spectrometer [28] to characterize weakly three-dimensional dispersive excitations sensitive to interactions between chains. Finally, to probe the low-energy gapped excitations, sensitive to local single-ion anisotropy, we used the cold neutrons on the RITA-II spectrometer [29].

## III. EXPERIMENTAL RESULTS

### A. Low-energy dispersive dynamics

The magnetic excitations characterizing the underlying spin-orbital ground state below the Néel temperature $T_N$ are summarized in Fig. 3. The overall dispersion is displayed in Fig. 3(a) taken on the MERLIN spectrometer with $E_i = 49$ meV. This constant momentum slice illustrates a strong
$V^{3+}$ ion from the surrounding oxygen octahedra is large. This model has been advocated based on x-ray and electron beam diffraction data [16]. An alternate suggestion [17] has been made for ordering of the complex basis of the $\hat{L}_z$ observable, which are complex linear combinations of the real $d$ orbitals. This is denoted as complex orbital ordering (COO) and is the basis state used in the weak-intermediate crystal-field limit for transition metal ions. Based on powder neutron spectroscopy and diffraction, it has also been suggested that the orbital ordering could possibly be intermediate between the two ROO and COO extremes [15]. A graphical representation of the two different extremes is illustrated in Fig. 1.

A theoretical study outlined in Ref. [18] suggested neutron spectroscopy as a means of distinguishing between COO and ROO orders at low temperatures. In particular the study noted the importance of spin-orbit coupling with COO giving multiple magnetic branches in the neutron response and also a larger magnetic zone-center gap than would be expected in the ROO model. Motivated by this, we investigate the magnetic neutron response in single crystals of $MgV_2O_4$. We study the hysteretic critical dynamics using diffraction, measure the dynamic response, and develop a theory to describe the magnetic excitations from the spin-orbital ground state. Given the spatially diffuse and correlated nature of the spin-orbital states used for our theory, we term this an excitonic approach.

We apply this excitonic approach to model our neutron spectroscopy results using the single-ion states and then treating the single-ion Hamiltonian (including spin-orbit coupling) and exchange energetics equally [19–22]. This model combined with the data favors the ordering of complex orbitals (COO) in $MgV_2O_4$ and illustrates the importance of using a complex orbital basis for understanding the properties of the first-row transition metal ions.

This paper is divided into four sections. First, we outline the materials preparation methodology and experimental techniques used to probe the magnetic fluctuations and critical scattering. Second, we outline the spectroscopic experimental results of both the low-energy spin-wave excitations and a higher-energy optic-like mode. Third, we present a theory for the magnetic excitations including spin-orbit coupling and the orbitally degenerate ground state. Finally, we investigate the hysteretic effects resulting from the orbital degeneracy using energy integrated neutron diffuse scattering.

## II. EXPERIMENTAL INFORMATION

*Materials preparation.* Single crystals [Fig. 2(b)] of $MgV_2O_4$ were grown using the floating zone technique with details provided in the Appendix. Given the extreme sensitivity of the magnetic properties to stoichiometry and chemical order in $MgV_2O_4$ [23], we have characterized our single crystals using both thermodynamic and scattering probes with neutrons and x-rays. The diffraction results are summarized in Fig. 2(c), with high-temperature structural ($T_s \simeq 60$ K) and magnetic [$T_N \simeq 40$ K from the magnetic $\vec{Q} = (1,1,0)$ Bragg peak] transitions consistent with the published literature [15,24].

*Characterization.* Figure 2(c) illustrates powder diffraction measurements using synchrotron x-rays (Diamond, I11) and neutrons (PSI, RITA-II). Powder synchrotron data shows a

first-order transition at $T_s \simeq 60$ K from a high-temperature cubic to a low-temperature tetragonal phase. This is confirmed on single-crystal neutron diffraction data following the $\vec{Q} = (0,0,2)$ Bragg peak on warming and cooling. A large hysteresis in the neutron intensity is observed between warming and cooling over the same temperature range where synchrotron x-rays measure a coexistence of tetragonal and cubic phases. This is further discussed in the Appendix where this region was found to extend from ~55–60 K. We have further confirmed that the structural properties of single crystals are consistent with the published structural phases through single-crystal neutron diffraction data (D9, ILL). The refinement results are summarized in Fig. 2(d) in both the cubic ($T = 140$ K) and tetragonal ($T = 50$ K) phases with the structural parameters listed in the Appendix. Figure 2(d) illustrates the calculated Bragg peak intensities ($\propto F^2$) as a function of the measured intensities with a straight line indicative that the models describe the results well. The slope of the line is a calibration factor scaling experiment and calculated values. We have confirmed the magnetic structure to be consistent with that outlined in the literature for $ZnV_2O_4$ [25] using polarized neutrons (DNS, MLZ). The results of this are schematically shown in Fig. 2(a) with the magnetic $V^{3+}$ magnetic moments pointing along the tetragonal c axis.

*Critical magnetic fluctuations.* In this paper we discuss energy-integrated polarized diffuse scattering measurements sensitive to the magnetic critical scattering and neutron inelastic data probing the spin-orbital dynamics. To investigate the magnetic critical scattering as a function of temperature, we studied the polarized neutron cross section using an xyz geometry provided by the DNS diffractometer (MLZ, Munich) [26]. Further details of the polarized beam measurements and analysis are provided in the Appendix.

*Neutron spectroscopy.* To probe the magnetic dynamics sensitive to the spin-orbital ground state, neutron spectroscopy was performed on four different spectrometers. To study the high-energy excitations, which represent spin-orbital excitations we used the MAPS spectrometer. The sample was oriented so that $\vec{k}_i$ was oriented along the $c^*$ axis and the intensity along $L$ integrated as typically done with time of flight direct geometry spectrometers to measure one- or two-dimensional spin excitations. To obtain an overview of the low-energy and low-dimensional dispersive dynamics, the multi rep rate option on MERLIN [27] was exploited combined with rotating the crystal. Further data was taken on the EIGER triple-axis spectrometer [28] to characterize weakly three-dimensional dispersive excitations sensitive to interactions between chains. Finally, to probe the low-energy gapped excitations, sensitive to local single-ion anisotropy, we used the cold neutrons on the RITA-II spectrometer [29].

## III. EXPERIMENTAL RESULTS

### A. Low-energy dispersive dynamics

The magnetic excitations characterizing the underlying spin-orbital ground state below the Néel temperature $T_N$ are summarized in Fig. 3. The overall dispersion is displayed in Fig. 3(a) taken on the MERLIN spectrometer with $E_i = 49$ meV. This constant momentum slice illustrates a strong
<figure><img src="image_5.png" /><figcaption>FIG. 3. Panlel (a) illustrates the magnetic excitations, taken from MERLIN ($E_i$ = 49 meV) up to ~40 meV showing a band of excitations. This cut is taken along the chain direction. Panel (b) illustrates MERLIN ($E_i$ = 24 meV) displaying a comparatively weaker dispersion perpendicular to the chain direction. Panel (c) displays high resolution data taken on RITA-II displaying an excitation gap reflecting underlying anisotropy. Panels (d) and (f) show data taken from EIGER characterizing the dispersion perpendicular to the chains.</figcaption></figure>

band of magnetic excitations, which extends up to ~35 meV. Perpendicular to this direction, a much more weakly dispersive mode is illustrated in Fig. 3(b). This reflects the underlying strongly one-dimensional coordination [30] of the $V^{3+}$ ions. The dispersion is further studied through a series of constant momentum scans using the EIGER thermal triple-axis spectrometer in Figs. 3(d) and 3(f). Analyzing the low-energy fluctuations with the cold triple-axis spectrometer RITA-II in Fig. 3(c), illustrates an energy gap of ~6 meV in the fluctuation spectrum. These results demonstrate the presence of strongly one-dimensional magnetic fluctuations with an energetic gap due to crystalline anisotropy. We discuss the origin of these below when we present an excitonic model for the spin-orbital excitations where we also discuss the relative magnetic exchange energies.

## B. Higher-energy gapped mode

One of the distinctions between COO and ROO orbital-ordering scenarios is the presence of an optical-like excitation at higher energies [18] corresponding to a transition between

<figure><img src="image_6.png" /><figcaption>FIG. 4. Panel (a) displays a constant momentum slice at $E = [50,55]$ meV taken using the MAPS spectrometer. Panel (b) illustrates a background-corrected constant momentum cut showing a well-defined peak in energy indicative of a second higher-energy magnetic band. The red line shows a fitted Gaussian on an exponentially decaying background peaked at 51 meV with a FWHM of 3 meV.</figcaption></figure>

the $j_{eff} = 2$ and $j_{eff} = 1$ manifolds (as theoretically outlined below). A higher-energy band of magnetic excitations is investigated using the MAPS spectrometer. The results are summarized in Fig. 4. Figure 4(a) displays a constant energy slice showing weakly correlated excitations near $\tilde{Q} = (2,0,0)$ and equivalent positions. A background (taken from large momentum detectors) corrected constant momentum cut is shown in Fig. 4(b) displaying a well-defined peak at ~50 meV. We discuss the origin of this additional excitation and connect it with the low-energy response below by investigating the spin-orbital neutron response theoretically.

## IV. THEORY

In $MgV_2O_4$ the magnetic $V^{3+}$ ($3d^2$) ions form a pyrochlore lattice with the $V^{3+}$ ions on the spinel B sites surrounded by an octahedral coordination of oxygen $O^{2-}$ [Figs. 5(b)
<figure><img src="image_5.png" /><figcaption>FIG. 3. Panlel (a) illustrates the magnetic excitations, taken from MERLIN ($E_i$ = 49 meV) up to ~40 meV showing a band of excitations. This cut is taken along the chain direction. Panel (b) illustrates MERLIN ($E_i$ = 24 meV) displaying a comparatively weaker dispersion perpendicular to the chain direction. Panel (c) displays high resolution data taken on RITA-II displaying an excitation gap reflecting underlying anisotropy. Panels (d) and (f) show data taken from EIGER characterizing the dispersion perpendicular to the chains.</figcaption></figure>

band of magnetic excitations, which extends up to ~35 meV. Perpendicular to this direction, a much more weakly dispersive mode is illustrated in Fig. 3(b). This reflects the underlying strongly one-dimensional coordination [30] of the V³⁺ ions. The dispersion is further studied through a series of constant momentum scans using the EIGER thermal triple-axis spectrometer in Figs. 3(d) and 3(f). Analyzing the low-energy fluctuations with the cold triple-axis spectrometer RITA-II in Fig. 3(c), illustrates an energy gap of ~6 meV in the fluctuation spectrum. These results demonstrate the presence of strongly one-dimensional magnetic fluctuations with an energetic gap due to crystalline anisotropy. We discuss the origin of these below when we present an excitonic model for the spin-orbital excitations where we also discuss the relative magnetic exchange energies.

## B. Higher-energy gapped mode

One of the distinctions between COO and ROO orbital-ordering scenarios is the presence of an optical-like excitation at higher energies [18] corresponding to a transition between

<figure><img src="image_6.png" /><figcaption>FIG. 4. Panel (a) displays a constant momentum slice at $E = [50,55]$ meV taken using the MAPS spectrometer. Panel (b) illustrates a background-corrected constant momentum cut showing a well-defined peak in energy indicative of a second higher-energy magnetic band. The red line shows a fitted Gaussian on an exponentially decaying background peaked at 51 meV with a FWHM of 3 meV.</figcaption></figure>

the $j_{eff} = 2$ and $j_{eff} = 1$ manifolds (as theoretically outlined below). A higher-energy band of magnetic excitations is investigated using the MAPS spectrometer. The results are summarized in Fig. 4. Figure 4(a) displays a constant energy slice showing weakly correlated excitations near $\tilde{Q} = (2,0,0)$ and equivalent positions. A background (taken from large momentum detectors) corrected constant momentum cut is shown in Fig. 4(b) displaying a well-defined peak at ~50 meV. We discuss the origin of this additional excitation and connect it with the low-energy response below by investigating the spin-orbital neutron response theoretically.

## IV. THEORY

In $MgV_2O_4$ the magnetic V³⁺ ($3d^2$) ions form a pyrochlore lattice with the V³⁺ ions on the spinel B sites surrounded by an octahedral coordination of oxygen O²⁻ [Figs. 5(b)
<figure><img src="image_7.png" /><figcaption>FIG. 5. (a) Tanabe-Sugano diagram for a $d^2$ ion in an octahedral crystal field. The Racah parameters have been chosen such that $C/B = 4.43$ [31]. For the discussion of the crystal-field parameters, we have followed Ref. [32] (Table 7.3) and used $B = 0.11$ eV and $Dq = 0.22$ eV, indicated by the dashed line. The VO$_6$ octahedron as viewed from (b) an isometric viewpoint and (c) the c axis. The octahedron is tetragonally distorted, as evidenced by the shorted bond lengths along c. A further trigonal distortion is present along [111] leading to two different bond lengths in the a-b plane. (d) Eigenvalues of the single-ion Hamiltonian for an octahedrally coordinated V$^{3+}$ ion ($S = 1, l = 1$) subject to a crystallographic distortion as described in Sec. IV A. Red lines indicate the ground state and excited states for which dipole-allowed transitions exist.</figcaption></figure>

and 5(c)), which determines the orbital ground state. Neighboring V$^{3+}$ ions occupy edge-sharing octahedra with non-magnetic tetrahedrally-coordinated Mg filling the voids between VO$_6$ octahedra. The magnetic interion interactions are likely governed by a combination of direct d-d orbital overlap and oxygen mediated superexchange via the ~90° V-O-V bonds.

In the rare earth pyrochlores [33], where the B site is occupied by the heavier 4f ions, spin-orbit coupling ($\lambda \sim Z^2$) [8] is the dominant energy scale and motivates the projection onto a ground-state manifold $J \equiv L + S$. The crystalline electric

field further splits the ground state with the resulting energy distribution dependent on the species of ion and local environment. In the case that the ground state is a Kramers doublet, a projection onto an effective $S = 1/2$ pseudospin [34–36] has been utilized and provided the motivation for seeking out quantum spin liquid [37] behavior in these materials.

The physics of the 3d pyrochlores is unlike that of their rare earth cousins owing to a different hierarchy of single-ion energy scales. Particularly, in 3d ions, the spin-orbit coupling ($\lambda$) is smaller than the crystalline electric field (Dq) [32,38,39]. With typical energy scales of Dq ~ 0.1 eV and $\lambda$ ~ 10 meV,
<figure><img src="image_7.png" /><figcaption>FIG. 5. (a) Tanabe-Sugano diagram for a *d*² ion in an octahedral crystal field. The Racah parameters have been chosen such that *C/B* = 4.43 [31]. For the discussion of the crystal-field parameters, we have followed Ref. [32] (Table 7.3) and used *B* = 0.11 eV and *Dq* = 0.22 eV, indicated by the dashed line. The VO₆ octahedron as viewed from (b) an isometric viewpoint and (c) the *c* axis. The octahedron is tetragonally distorted, as evidenced by the shorted bond lengths along *c*. A further trigonal distortion is present along [111] leading to two different bond lengths in the *a-b* plane. (d) Eigenvalues of the single-ion Hamiltonian for an octahedrally coordinated V³⁺ ion (*S* = 1, *l* = 1) subject to a crystallographic distortion as described in Sec. IV A. Red lines indicate the ground state and excited states for which dipole-allowed transitions exist.</figcaption></figure>

and 5(c)), which determines the orbital ground state. Neighboring V³⁺ ions occupy edge-sharing octahedra with non-magnetic tetrahedrally-coordinated Mg filling the voids between VO₆ octahedra. The magnetic interion interactions are likely governed by a combination of direct *d*-d orbital overlap and oxygen mediated superexchange via the ~90° V-O-V bonds.

In the rare earth pyrochlores [33], where the B site is occupied by the heavier 4f ions, spin-orbit coupling (λ ~ Z²) [8] is the dominant energy scale and motivates the projection onto a ground-state manifold J ≡ L + S. The crystalline electric

field further splits the ground state with the resulting energy distribution dependent on the species of ion and local environment. In the case that the ground state is a Kramers doublet, a projection onto an effective *S* = 1/2 pseudospin [34–36] has been utilized and provided the motivation for seeking out quantum spin liquid [37] behavior in these materials.

The physics of the 3d pyrochlores is unlike that of their rare earth cousins owing to a different hierarchy of single-ion energy scales. Particularly, in 3d ions, the spin-orbit coupling (λ) is smaller than the crystalline electric field (Dq) [32,38,39]. With typical energy scales of Dq ~ 0.1 eV and λ ~ 10 meV,
spin-orbit coupling is a perturbation on the crystalline electric field Hamiltonian and the ground state is defined by the orbital angular momentum L. Within this ground-state manifold, the spin-orbit coupling, magnetic exchange, and Jahn-Teller energy scales are comparable. In the 3d transition metal compounds it is therefore necessary to consider both single-ion spin-orbital energy scales ($\mathcal{H}_{SI}$) and the corresponding magnetic exchange interactions ($\mathcal{H}_{exch}$) [38],

$$ \mathcal{H} = \mathcal{H}_{SI} + \mathcal{H}_{exch}. $$

We will begin by discussing the single-ion physics of the 3d² ions to understand the magnetic ground state before considering the magnetic interion exchange in vanadium spinels. Using a Green's function formalism, the dynamical structure factor measured with neutrons will then be calculated using the random phase approximation (RPA) applied to MgV₂O₄, treating the single-ion Hamiltonian, which produces a quantized multilevel ground state explicitly.

### A. $\mathcal{H}_{SI}$ - Single-ion physics

The single-ion physics discussed in this section is schematically shown in Fig. 5, which determines the quantized spin-orbital ground state of the magnetic V³⁺ ions. We consider first the free-ion V³⁺ ground state followed by the effects of the crystalline electric field, spin-orbit coupling, distortions from an octahedral environment, and finally the Zeeman-like molecular field originating from magnetic ordering in the low-temperature Néel phase.

#### 1. $^3F$ - Free ion ground state

For the case of a free V³⁺ ions with two electrons in the five degenerate d orbitals, the ground state is determined by Hund's rules, which defines L = 3 and S = 1. This fixes the orbital ground state of the free-ion V³⁺ state to be $^3F$ in spectroscopic notation.

#### 2. $\mathcal{H}_{CEF}$ - Crystalline electric field

As discussed above, the dominant energy scale for magnetic 3d is the crystalline electric field originating from the O²⁻ ions forming an octahedron around the V³⁺ ion. Application of an octahedral crystalline electric field in terms of Steven's operators $O_l^m$ gives the following:

$$ \mathcal{H}_{CEF} = B_4(O_4^0 + 5O_4^4) $$

where $B_4 < 0$ for d² ions [32]. The energy spectrum resulting from this crystalline electric-field Hamiltonian is schematically displayed in Fig. 5(a) (Tanabe-Sugano diagram) and parameterized in terms of the crystal-field strength ($Dq$) and the Racah parameter $B$, which physically corresponds to the energy cost associated with the Coulomb repulsion. The limit $Dq/B \to 0$ corresponds to the weak crystal-field limit where Hund's rules across all d orbitals applies. In the large $Dq/B$ limit, the crystal-field energy scale is dominant and in some transition metal ions (such as Co²⁺) can result in spin transitions.

Figure 5(a) shows that the octahedral crystal field splits the orbital $^3F$ ground state into two triplets ($^3T_{1,2}$) and a singlet ($^3A_2$), with the $^3T_1$ triplet being the orbital ground

state. The splitting between $^3T_1$ and the first excited triplet is $8Dq = 480B_4 \approx 1.8$ eV [32] (in the limit of small $Dq$). Each of the triplets under the octahedral crystal field forms an irreducible representation of the octahedral double group $O'$ [40]. By the orthogonality of different irreps of the same group, for an operator that is invariant under all octahedral symmetry operations, matrix elements between triplets are zero [40]. Assuming further symmetry lowering terms away from an octahedral field are small and the large energy separation in comparison to the temperatures of interest, one can justifiably work in the ground-state $^3T_1$ multiplet, neglecting the excited states.

Considering energy levels beyond the $^3F$ manifold in Fig. 5(a), the octahedral field mixes in the excited $^3P$ state [Fig. 5(a) when $Dq/B \to 0$]. The lowest energy-orbital state of the $^3P$ manifold has the same symmetry being a $^3T_1$ triplet [32] and therefore the overall ground state is a linear combination of the two triplet states from the $^3F$ and $^3P$ manifolds,

$$ |\psi(^3T_1)\rangle = \epsilon |^3F(^3T_1)\rangle + \tau|^3P(^3T_1)\rangle \quad (1) $$

where $\epsilon^2 + \tau^2 = 1$ [32]. One can represent the ground-state orbital triplet using the fictitious orbital angular momentum operator $\hat{L}$ via the transformation $\hat{L} = \alpha\hat{I}$, where the projection factor $\alpha < 0$ can be read off the block describing the ground-state manifold after projecting the matrices $\hat{L}_\alpha$ onto the space spanned by the eigenvectors of $\mathcal{H}_{CEF}$. For the $|^3F(^3T_1)\rangle$ block, $\alpha = -\frac{3}{2}$ [41], whilst for the $|^3P(^3T_1)\rangle$ state, $\alpha = +1$. The admixture of these two states leads to an effective projection factor that varies between these two extremal values, $\alpha = \frac{5}{2}\tau^2 - \frac{3}{2}$.

It is instructive to discuss the orbital triplet ground state in terms of the electronic orbital basis states, which have the advantage over the eigenstates of the observables $L^2$ and $L_z$ of being real and therefore give a greater connection to the microscopic picture of electron hopping through chemical bonding. In the 3d ions, the outer electrons partially fill d orbitals, which are split in an octahedral crystal field into the triply degenerate $t_{2g}$ level ($d_{xy}$, $d_{yz}$, $d_{xz}$) and doubly degenerate $e_g$ level ($d_{x^2-y^2}$, $d_{3z^2-r^2}$). The energy splitting between $t_{2g}$ and $e_g$ levels is set by $10Dq$, and therefore this parameter is key in distinguishing the weak (small $Dq$) from the strong (large $Dq$) crystal-field limits. For weak crystal fields it is more natural to use a complex orbital basis that are eigenstates of the observables $L^2$ and $L_z$. In the large crystal-field limit, the real $t_{2g}$ and $e_g$ states are the convenient basis.

Neglecting the spin-orbit coupling and further distortions, the ground state for two electrons in the 3d orbitals is triply degenerate with each member of ground-state manifold having two of the $t_{2g}$ levels occupied. In terms of these real basis states one can represent the ground-state triplet with the admixture caused by the crystal field [Eq. (1)] equivalently as

$$ |\psi(^3T_1)\rangle = \cos\theta |t_{2g}, t_{2g}\rangle - \sin\theta |t_{2g}, e_g\rangle. \quad (2) $$

The first of these two basis states has two occupied $t_{2g}$ levels as expected in an octahedral crystal field, and the second has an occupied $e_g$ level but has the same symmetry as $|t_{2g}, t_{2g}\rangle$. By diagonalizing the energy matrix for a d² ion with both a Coulomb term and an octahedral crystal field [32]
spin-orbit coupling is a perturbation on the crystalline electric field Hamiltonian and the ground state is defined by the orbital angular momentum L. Within this ground-state manifold, the spin-orbit coupling, magnetic exchange, and Jahn-Teller energy scales are comparable. In the 3d transition metal compounds it is therefore necessary to consider both single-ion spin-orbital energy scales ($\mathcal{H}_{SI}$) and the corresponding magnetic exchange interactions ($\mathcal{H}_{exch}$) [38],

$$ \mathcal{H} = \mathcal{H}_{SI} + \mathcal{H}_{exch}. $$

We will begin by discussing the single-ion physics of the 3d² ions to understand the magnetic ground state before considering the magnetic interion exchange in vanadium spinels. Using a Green's function formalism, the dynamical structure factor measured with neutrons will then be calculated using the random phase approximation (RPA) applied to MgV₂O₄, treating the single-ion Hamiltonian, which produces a quantized multilevel ground state explicitly.

### A. $\mathcal{H}_{SI}$ - Single-ion physics

The single-ion physics discussed in this section is schematically shown in Fig. 5, which determines the quantized spin-orbital ground state of the magnetic V³⁺ ions. We consider first the free-ion V³⁺ ground state followed by the effects of the crystalline electric field, spin-orbit coupling, distortions from an octahedral environment, and finally the Zeeman-like molecular field originating from magnetic ordering in the low-temperature Néel phase.

#### 1. $^3F$ - Free ion ground state

For the case of a free V³⁺ ions with two electrons in the five degenerate d orbitals, the ground state is determined by Hund's rules, which defines L = 3 and S = 1. This fixes the orbital ground state of the free-ion V³⁺ state to be $^3F$ in spectroscopic notation.

#### 2. $\mathcal{H}_{CEF}$ - Crystalline electric field

As discussed above, the dominant energy scale for magnetic 3d is the crystalline electric field originating from the O²⁻ ions forming an octahedron around the V³⁺ ion. Application of an octahedral crystalline electric field in terms of Steven's operators $O_l^m$ gives the following:

$$ \mathcal{H}_{CEF} = B_4(O_4^0 + 5O_4^4) $$

where $B_4 < 0$ for d² ions [32]. The energy spectrum resulting from this crystalline electric-field Hamiltonian is schematically displayed in Fig. 5(a) (Tanabe-Sugano diagram) and parameterized in terms of the crystal-field strength ($Dq$) and the Racah parameter $B$, which physically corresponds to the energy cost associated with the Coulomb repulsion. The limit $Dq/B \to 0$ corresponds to the weak crystal-field limit where Hund's rules across all d orbitals applies. In the large $Dq/B$ limit, the crystal-field energy scale is dominant and in some transition metal ions (such as Co²⁺) can result in spin transitions.

Figure 5(a) shows that the octahedral crystal field splits the orbital $^3F$ ground state into two triplets ($^3T_{1,2}$) and a singlet ($^3A_2$), with the $^3T_1$ triplet being the orbital ground

state. The splitting between $^3T_1$ and the first excited triplet is $8Dq = 480B_4 \approx 1.8$ eV [32] (in the limit of small $Dq$). Each of the triplets under the octahedral crystal field forms an irreducible representation of the octahedral double group $O'$ [40]. By the orthogonality of different irreps of the same group, for an operator that is invariant under all octahedral symmetry operations, matrix elements between triplets are zero [40]. Assuming further symmetry lowering terms away from an octahedral field are small and the large energy separation in comparison to the temperatures of interest, one can justifiably work in the ground-state $^3T_1$ multiplet, neglecting the excited states.

Considering energy levels beyond the $^3F$ manifold in Fig. 5(a), the octahedral field mixes in the excited $^3P$ state [Fig. 5(a) when $Dq/B \to 0$]. The lowest energy-orbital state of the $^3P$ manifold has the same symmetry being a $^3T_1$ triplet [32] and therefore the overall ground state is a linear combination of the two triplet states from the $^3F$ and $^3P$ manifolds,

$$ |\psi(^3T_1)\rangle = \epsilon |^3F(^3T_1)\rangle + \tau|^3P(^3T_1)\rangle \quad (1) $$

where $\epsilon^2 + \tau^2 = 1$ [32]. One can represent the ground-state orbital triplet using the fictitious orbital angular momentum operator $\hat{L}$ via the transformation $\hat{L} = \alpha\hat{I}$, where the projection factor $\alpha < 0$ can be read off the block describing the ground-state manifold after projecting the matrices $\hat{L}_\alpha$ onto the space spanned by the eigenvectors of $\mathcal{H}_{CEF}$. For the $|^3F(^3T_1)\rangle$ block, $\alpha = -\frac{3}{2}$ [41], whilst for the $|^3P(^3T_1)\rangle$ state, $\alpha = +1$. The admixture of these two states leads to an effective projection factor that varies between these two extremal values, $\alpha = \frac{5}{2}\tau^2 - \frac{3}{2}$.

It is instructive to discuss the orbital triplet ground state in terms of the electronic orbital basis states, which have the advantage over the eigenstates of the observables $L^2$ and $L_z$ of being real and therefore give a greater connection to the microscopic picture of electron hopping through chemical bonding. In the 3d ions, the outer electrons partially fill d orbitals, which are split in an octahedral crystal field into the triply degenerate $t_{2g}$ level ($d_{xy}$, $d_{yz}$, $d_{xz}$) and doubly degenerate $e_g$ level ($d_{x^2-y^2}$, $d_{3z^2-r^2}$). The energy splitting between $t_{2g}$ and $e_g$ levels is set by $10Dq$, and therefore this parameter is key in distinguishing the weak (small $Dq$) from the strong (large $Dq$) crystal-field limits. For weak crystal fields it is more natural to use a complex orbital basis that are eigenstates of the observables $L^2$ and $L_z$. In the large crystal-field limit, the real $t_{2g}$ and $e_g$ states are the convenient basis.

Neglecting the spin-orbit coupling and further distortions, the ground state for two electrons in the 3d orbitals is triply degenerate with each member of ground-state manifold having two of the $t_{2g}$ levels occupied. In terms of these real basis states one can represent the ground-state triplet with the admixture caused by the crystal field [Eq. (1)] equivalently as

$$ |\psi(^3T_1)\rangle = \cos\theta |t_{2g}, t_{2g}\rangle - \sin\theta |t_{2g}, e_g\rangle. \quad (2) $$

The first of these two basis states has two occupied $t_{2g}$ levels as expected in an octahedral crystal field, and the second has an occupied $e_g$ level but has the same symmetry as $|t_{2g}, t_{2g}\rangle$. By diagonalizing the energy matrix for a d² ion with both a Coulomb term and an octahedral crystal field [32]
[Fig. 5(a)] one finds $\tan 2\theta = 12B/(9B + 10Dq)$ and using the correspondence between the two descriptions [Eqs. (1) and (2)] one can quantify the fraction of the $^3P$ level in the ground-state orbital triplet, $\tau = \frac{1}{\sqrt{5}}(\cos\theta - 2\sin\theta)$. With $B = 0.11$ eV and $Dq = 0.22$ eV, as expected for a free $V^{3+}$ ion, we have $\tau \approx 0.27$. We therefore expect a projection factor of $\alpha \approx -1.32$ in comparison to a value of -1.5 in the absence of mixing. We discuss this parameter later in context of the reported ordered magnetic moment and the possibility for quantum fluctuations.

### 3. $H_{SO}$ - Spin-orbit coupling

The next largest term [illustrated in Fig. 5(d)] in the single-ion Hamiltonian for 3d ions is the spin-orbit coupling

$$H_{SO} = \lambda \hat{L} \cdot \hat{S} = \alpha \lambda \hat{I} \cdot \hat{S}$$

where $\lambda > 0$ for $d^2$ systems [42]. This term splits the ground-state ninefold spin-orbital manifold into three $\hat{j} = \hat{I} + \hat{S}$ levels according to the Landé interval rule with $j_{eff} = 0, 1, 2$ [Fig. 5(d)]. In $d^2$ ions, the ground state is the $j_{eff} = 2$, separated from the $j_{eff} = 1$ state by $|2\alpha\lambda|$ [32]. In the language of the strong crystal-field $t_{2g}$ levels, the effect of the inclusion of spin-orbit coupling is to mix the $d_{xz}$ and $d_{yz}$ orbitals in complex combinations such that three basis states in which the Hamiltonian is diagonal become $|l_z = 0\rangle = |d_{xy}\rangle$ and $|l_z = \pm 1\rangle = \frac{1}{\sqrt{2}}(|d_{xz}\rangle \pm i|d_{yz}\rangle)$ [43].

### 4. $H_{dis}$-Distortion

In MgV$_2$O$_4$ at the low temperatures of interest here, the symmetry of the local VO$_6$ octahedra is not $O_h$, since a Jahn-Teller distortion, originating from the orbital degeneracy of the $d^2$ vanadium ion, occurs on the transition from the high-temperature cubic phase to the low-temperature tetragonal phase [15]. This distortion results in the octahedra being subtly compressed along the fourfold $\hat{z}$ axis [illustrated in Fig. 5(b)]. The distortion can be modeled by the following term in the Hamiltonian:

$$H_{dis} = \Gamma (\hat{l}_z^2 - \frac{2}{3}) \quad (3)$$

where $\Gamma < 0$ for a compression. As displayed in Fig. 5(d) the effect of this term is to break the fivefold orbital ground state $j_{eff} = 2$ degeneracy into three levels with a ground-state doublet, an excited-state doublet, and an excited-state singlet.

Considering the effect of this distortion only on orbitals (without considering spin or spin-orbit coupling and in the limit of $|\Gamma| \gg |\lambda|$) and in terms of the strong crystal-field basis of real orbitals, the effect of this axial distortion is to break the ground-state $t_{2g}$ orbital triplet degeneracy, yielding a ground-state orbital doublet and excited singlet. The distortion lowers the energy of the $d_{xy}$ orbital relative to the $d_{xz}$ and $d_{yz}$ orbital. If we populate these two levels with two d electrons applying Pauli's exclusion principle, this results in a doubly degenerate ground state with a hole in either the $d_{xz}$ or $d_{yz}$ orbital.

We note that in addition to the primary tetragonal compression driven by orbital degeneracy, the VO$_6$ octahedra are trigonal distorted [Fig. 5(c)]—a compression along the threefold [111] axis (even within the high-temperature cubic

<table><caption>TABLE I. Orbital wavefunction coefficients in common orbital-ordering regimes for electron I, electron II. Coefficients correspond to Eq. (4). Column two states the dominant energy scale in each case.</caption><thead><tr><th>OO</th><th>Dom.</th><th>α</th><th>β</th><th>γ</th><th>θ</th><th>φ</th></tr></thead><tbody><tr><td>COO</td><td>λ</td><td>1,0</td><td>0, 1/√2</td><td>0, 1/√2</td><td>0,0</td><td>0, π/2</td></tr><tr><td>ROO</td><td>Γ[001] > 0</td><td>1,0</td><td>0, 1/√2</td><td>0, 1/√2</td><td>0,0</td><td>0,0</td></tr><tr><td>Trig. COO</td><td>Γ[111]</td><td>1/√3, 1/√3</td><td>1/√3, 1/√3</td><td>1/√3, 1/√3</td><td>0, 2π/3</td><td>0, -2π/3</td></tr></tbody></table>

phase). Within the manifold of the ground-state $t_{2g}$ orbitals, the multiplet structure under this distortion is the same as the tetragonal distortion [Eq. (3)] [43]. In this paper, we will take advantage of the projection onto the ground-state triplet, $\hat{L} = \alpha\hat{I}$ and treat the mixing of the $e_g$ levels as a small perturbation and so we can collect all contributions to the distortion Hamiltonian into a single distortion parameter $\Gamma$. Under a dominant trigonal distortion, the orbital spectrum has the form $|a_{1g}\rangle = \frac{1}{\sqrt{3}}(|d_{xy}\rangle + |d_{xz}\rangle + |d_{yz}\rangle)$ with the excited doublet $|e_g^\pi\rangle = \pm \frac{1}{\sqrt{3}}(|d_{xy}\rangle + e^{\pm 2\pi i/3} |d_{xz}\rangle + e^{\mp 2\pi i/3} |d_{yz}\rangle)$ [43]. From crystallographic considerations, the low-temperature tetragonal distortion discussed above is expected to be dominant; however, one might expect a small difference in orbital coefficients arising due to the subleading trigonal distortion. Such a situation maybe supported by the small ~8° canting of the spin towards the apical oxygen reported in some diffraction studies [15].

The orbital state for each electron can be written in general as

$$|\psi\rangle = \alpha |d_{xy}\rangle + \beta e^{i\theta} |d_{xz}\rangle + \gamma e^{i\phi} |d_{yz}\rangle \quad (4)$$

where $\alpha^2 + \beta^2 + \gamma^2 = 1$ ensures normalization. Based on the single-ion physics of the $V^{3+}$ ion ($3d^2$) alone, the orbital order in MgV$_2$O$_4$ is expected to be intermediate between the regimes of validity of real orbital order (ROO), which occurs when the tetragonal distortion dominates and complex orbital order (COO) where the spin-orbit coupling dominates, since the energy scales of the distortion and the spin-orbit coupling in the 3d ions are typically nearly comparable [19,41,44]. The weightings of the $d_{\alpha\beta}$ orbitals for each electron within each regime are summarized in Table I. We will discuss the nature of the orbital order further later in the paper based on the insight gained through our analysis.

### 5. $H_{MF}$ - Molecular field below $T_N$

The final term present in the single-ion Hamiltonian is the molecular field. This results from a mean-field decoupling of the exchange interaction between coupled ions, and is required such that the single-ion ground state about which one expands is correct. This term breaks time-reversal symmetry

$$H_{MF} = h_{MF} \hat{S}_z$$

as is consistent with the establishment of long-range magnetic order. The magnitude of $h_{MF}$ will be discussed further once the interion spin Hamiltonian $H_{exch}$ has been introduced.

In terms of the single-ion Hamiltonian, the ground-state order is expected to be ($S = 1, l = +1$) with positive $+l$
[Fig. 5(a)] one finds $\tan 2\theta = 12B/(9B + 10Dq)$ and using the correspondence between the two descriptions [Eqs. (1) and (2)] one can quantify the fraction of the $^3P$ level in the ground-state orbital triplet, $\tau = \frac{1}{\sqrt{5}}(\cos\theta - 2\sin\theta)$. With $B = 0.11$ eV and $Dq = 0.22$ eV, as expected for a free $V^{3+}$ ion, we have $\tau \approx 0.27$. We therefore expect a projection factor of $\alpha \approx -1.32$ in comparison to a value of -1.5 in the absence of mixing. We discuss this parameter later in context of the reported ordered magnetic moment and the possibility for quantum fluctuations.

### 3. $\mathcal{H}_{SO}$ - Spin-orbit coupling

The next largest term [illustrated in Fig. 5(d)] in the single-ion Hamiltonian for 3d ions is the spin-orbit coupling

$$ \mathcal{H}_{SO} = \lambda \hat{\mathbf{L}} \cdot \hat{\mathbf{S}} = \alpha \lambda \hat{\mathbf{I}} \cdot \hat{\mathbf{S}} $$

where $\lambda > 0$ for $d^2$ systems [42]. This term splits the ground-state ninefold spin-orbital manifold into three $\hat{\mathbf{j}} = \hat{\mathbf{I}} + \hat{\mathbf{S}}$ levels according to the Landé interval rule with $j_{\text{eff}} = 0, 1, 2$ [Fig. 5(d)]. In $d^2$ ions, the ground state is the $j_{\text{eff}} = 2$, separated from the $j_{\text{eff}} = 1$ state by $|2\alpha\lambda|$ [32]. In the language of the strong crystal-field $t_{2g}$ levels, the effect of the inclusion of spin-orbit coupling is to mix the $d_{xz}$ and $d_{yz}$ orbitals in complex combinations such that three basis states in which the Hamiltonian is diagonal become $|l_z = 0\rangle = |d_{xy}\rangle$ and $|l_z = \pm 1\rangle = \frac{1}{\sqrt{2}}(|d_{xz}\rangle \pm i|d_{yz}\rangle)$ [43].

### 4. $\mathcal{H}_{dis}$-Distortion

In MgV$_2$O$_4$ at the low temperatures of interest here, the symmetry of the local VO$_6$ octahedra is not $O_h$, since a Jahn-Teller distortion, originating from the orbital degeneracy of the $d^2$ vanadium ion, occurs on the transition from the high-temperature cubic phase to the low-temperature tetragonal phase [15]. This distortion results in the octahedra being subtly compressed along the fourfold $\hat{z}$ axis [illustrated in Fig. 5(b)]. The distortion can be modeled by the following term in the Hamiltonian:

$$ \mathcal{H}_{dis} = \Gamma(\hat{l}_z^2 - \frac{2}{3}) \quad (3) $$

where $\Gamma < 0$ for a compression. As displayed in Fig. 5(d) the effect of this term is to break the fivefold orbital ground state $j_{\text{eff}} = 2$ degeneracy into three levels with a ground-state doublet, an excited-state doublet, and an excited-state singlet.

Considering the effect of this distortion only on orbitals (without considering spin or spin-orbit coupling and in the limit of $|\Gamma| \gg |\lambda|$) and in terms of the strong crystal-field basis of real orbitals, the effect of this axial distortion is to break the ground-state $t_{2g}$ orbital triplet degeneracy, yielding a ground-state orbital doublet and excited singlet. The distortion lowers the energy of the $d_{xy}$ orbital relative to the $d_{xz}$ and $d_{yz}$ orbital. If we populate these two levels with two d electrons applying Pauli's exclusion principle, this results in a doubly degenerate ground state with a hole in either the $d_{xz}$ or $d_{yz}$ orbital.

We note that in addition to the primary tetragonal compression driven by orbital degeneracy, the VO$_6$ octahedra are trigonal distorted [Fig. 5(c)]—a compression along the threefold [111] axis (even within the high-temperature cubic

<table><caption>TABLE I. Orbital wavefunction coefficients in common orbital-ordering regimes for electron I, electron II. Coefficients correspond to Eq. (4). Column two states the dominant energy scale in each case.</caption><thead><tr><th>OO</th><th>Dom.</th><th>α</th><th>β</th><th>γ</th><th>θ</th><th>φ</th></tr></thead><tbody><tr><td>COO</td><td>λ</td><td>1,0</td><td>0, 1/√2</td><td>0, 1/√2</td><td>0,0</td><td>0, π/2</td></tr><tr><td>ROO</td><td>Γ<sub>[001]</sub> > 0</td><td>1,0</td><td>0, 1/√2</td><td>0, 1/√2</td><td>0,0</td><td>0,0</td></tr><tr><td>Trig. COO</td><td>Γ<sub>[111]</sub></td><td>1/√3, 1/√3</td><td>1/√3, 1/√3</td><td>1/√3, 1/√3</td><td>0, 2π/3</td><td>0, -2π/3</td></tr></tbody></table>

phase). Within the manifold of the ground-state $t_{2g}$ orbitals, the multiplet structure under this distortion is the same as the tetragonal distortion [Eq. (3)] [43]. In this paper, we will take advantage of the projection onto the ground-state triplet, $\hat{\mathbf{L}} = \alpha\hat{\mathbf{I}}$ and treat the mixing of the $e_g$ levels as a small perturbation and so we can collect all contributions to the distortion Hamiltonian into a single distortion parameter $\Gamma$. Under a dominant trigonal distortion, the orbital spectrum has the form $|a_{1g}\rangle = \frac{1}{\sqrt{3}}(|d_{xy}\rangle + |d_{xz}\rangle + |d_{yz}\rangle)$ with the excited doublet $|e_g^\pi\rangle = \pm \frac{1}{\sqrt{3}}(|d_{xy}\rangle + e^{\pm 2\pi i/3} |d_{xz}\rangle + e^{\mp 2\pi i/3} |d_{yz}\rangle)$ [43]. From crystallographic considerations, the low-temperature tetragonal distortion discussed above is expected to be dominant; however, one might expect a small difference in orbital coefficients arising due to the subleading trigonal distortion. Such a situation maybe supported by the small ~8° canting of the spin towards the apical oxygen reported in some diffraction studies [15].

The orbital state for each electron can be written in general as

$$ |\psi\rangle = \alpha |d_{xy}\rangle + \beta e^{i\theta} |d_{xz}\rangle + \gamma e^{i\phi} |d_{yz}\rangle \quad (4) $$

where $\alpha^2 + \beta^2 + \gamma^2 = 1$ ensures normalization. Based on the single-ion physics of the $V^{3+}$ ion ($3d^2$) alone, the orbital order in MgV$_2$O$_4$ is expected to be intermediate between the regimes of validity of real orbital order (ROO), which occurs when the tetragonal distortion dominates and complex orbital order (COO) where the spin-orbit coupling dominates, since the energy scales of the distortion and the spin-orbit coupling in the 3d ions are typically nearly comparable [19,41,44]. The weightings of the $d_{\alpha\beta}$ orbitals for each electron within each regime are summarized in Table I. We will discuss the nature of the orbital order further later in the paper based on the insight gained through our analysis.

### 5. $\mathcal{H}_{MF}$ - Molecular field below $T_N$

The final term present in the single-ion Hamiltonian is the molecular field. This results from a mean-field decoupling of the exchange interaction between coupled ions, and is required such that the single-ion ground state about which one expands is correct. This term breaks time-reversal symmetry

$$ \mathcal{H}_{MF} = h_{MF} \hat{S}_z $$

as is consistent with the establishment of long-range magnetic order. The magnitude of $h_{MF}$ will be discussed further once the interion spin Hamiltonian $\mathcal{H}_{exch}$ has been introduced.

In terms of the single-ion Hamiltonian, the ground-state order is expected to be ($S = 1, l = +1$) with positive $+l$
<figure><img src="image_8.png" /><figcaption>FIG. 6. Magnetic structure of MgV$_{2}$O$_{4}$ as viewed from an iso- metric viewpoint. The structure comprises chains in the a-b plane with the nearest-neighbor bond J$_{2}$. Interchain bond J$_{1}$ is frustrated and couples chains along [101] and [011]. A longer range bond J$_{3}$ couples perpendicular chains.</figcaption></figure>

(instead of − l ) selected due to the negative spin-orbit coupling constant in d 2 ions, which promotes alignment of the spin and orbital moments. This expected ground state gives rise to a magnetic moment µ = µ$_{B}$ ( ˆ L + 2 ˆ S ) = µ$_{B}$ ( − 1 . 32 + 2) = 0 . 68 µ$_{B}$ , which agrees well with the observed (reduced) mag- netic moment, µ = 0 . 47 µ$_{B}$ [15]. The discrepancy between the observed and calculated magnetic moment suggests the presence of quantum fluctuations, S − Δ S ≈ 0 . 9, as is to be expected given the small value of S and reduced dimension- ality originating from the frustrated lattice. This reduction in the ordered moment in turn contributes to a reduction in the molecular mean field over what would be expected for the full spin value.

B. Interion coupling

We now discuss the interion interactions present in MgV$_{2}$O$_{4}$. In the high-temperature cubic phase, the V 3+ ions lie on an ideal pyrochlore lattice, with each nearest-neighbor V-V bond equivalent by symmetry. In the case of antiferro- magnetic nearest-neighbor interactions, the spin arrangement on a pyrochlore lattice is geometrically frustrated [45], open- ing the possibility of spin liquid states [37,46] resulting from the large ground-state degeneracy. Pyrochlore systems with S = 1 have been studied both theoretically [47,48] and exper- imentally [49–51], particularly in search of an understanding of the crossover from quantum to classical spin liquid behav- ior [48,51]. The vanadium spinels have also attracted great interest [17,18,52–54] since the orbital degree of freedom al- lows for Jahn-Teller distortions, which break the ground-state degeneracy and permit magnetic order. This behavior is ob- served in MgV$_{2}$O$_{4}$ at T$_{S}$ ≃ 60 K, where the system becomes tetragonal before ordering antiferromagnetically at T$_{N}$ ≃ 40 K [15]. The cubic to tetragonal structural transition renders the nearest-neighbor bonds along [110] and [101] inequivalent (Fig. 6).

The orbital order plays an important role in the magnetic interion exchange since both the d - d and d - p orbital overlap are dependent on the V 3+ orbital ground state, and ultimately where the electron hole lies. This provides a mechanism for a

breaking of the equivalence of magnetic bonds both in terms of the amplitude of the direct and super-exchange interactions [43]. The question of the magnetic exchange in the vanadium spinels has been discussed at length by Di Matteo et al. [52] and Perkins et al. [18]. We will now briefly summarize the findings of the aforementioned references [18,52] insofar as they are relevant to the magnetic exchange in MgV$_{2}$O$_{4}$.

Di Matteo et al. considered the limit where the tetragonal distortion is dominant (ROO, real electronic orbitals discussed above) and where the spin-orbit coupling is dominant (COO, complex basis described above when the crystalline electric field Dq is weak) [52]. Both limits give rise to multiple mag- netic and orbital configurations dependent on the underlying physical parameters [52] and described by the general Heisen- berg spin Hamiltonian

$$
\mathcal{H}_{\text{exch}} = \frac{1}{2} \sum_{ij} \mathcal{J}_{ij} \hat{\mathbf{S}}_i \cdot \hat{\mathbf{S}}_j,
$$

where the factor of 1 2 accounts for double counting. In the case of ROO, only one possible orbital ground state is consistent with the experimentally observed magnetic structure [15]. This ROO ground state comprises strong antiferromagnetic (AFM) bonds in the xy plane and a weaker ferromagnetic cou- pling along [101] and [011]. The ROO model thus suggests that MgV$_{2}$O$_{4}$ is comprised of strongly-coupled AFM chains in the xy plane with weak FM coupling between chains. In each tetrahedron, two FM bonds are satisfied and two FM bonds are frustrated.

In contrast, when the spin-orbit coupling is dominant (COO), there are two ground states consistent with the observed magnetic structure. Both possible orbital ground states have strongly coupled AFM chains in the xy plane; however, one of these states has two weak AFM and two weak FM bonds per tetrahedron—thus no bond is frus- trated. In the other possible COO ground state, the interchain bonds are weak AFM bonds, two of which are frustrated. Based on the parameters obtained from spectroscopy data for vanadium perovskites [55,56] and typical strengths of the spin-orbit coupling for V 3+ reported in the literature [17,32], Di Matteo et al. concluded that the latter of the two COO states is likely the ground state [52].

Finally we note that the predicted magnetic ground states in the ROO and COO limits differ only in the sign of the interchain coupling, with both interchain couplings of a sim- ilar strength [52]. It is interesting to speculate as to whether these two limits are continuously connected via a phase with zero interchain coupling and hence might be tuned through with the application of strain, using the tetragonal distortion Γ as a tuning parameter. It is also interesting to note that a metal-insulator has been predicted in MgV$_{2}$O$_{4}$ based on ab initio calculations [57].

C. Neutron scattering intensity calculation

In order to calculate the neutron scattering response, we employ a spin-orbital exciton model in terms of Green’s func- tions as applied previously to describe multilevel systems [19–21,41,58,59]. A full derivation for a general single-Q multilevel magnetic system is presented in Ref. [59], here we
<figure><img src="image_8.png" /><figcaption>FIG. 6. Magnetic structure of MgV$_{2}$O$_{4}$ as viewed from an iso- metric viewpoint. The structure comprises chains in the a-b plane with the nearest-neighbor bond J$_{2}$. Interchain bond J$_{1}$ is frustrated and couples chains along [101] and [011]. A longer range bond J$_{3}$ couples perpendicular chains.</figcaption></figure>

(instead of − l ) selected due to the negative spin-orbit coupling constant in d 2 ions, which promotes alignment of the spin and orbital moments. This expected ground state gives rise to a magnetic moment µ = µ$_{B}$ ( ˆ L + 2 ˆ S ) = µ$_{B}$ ( − 1 . 32 + 2) = 0 . 68 µ$_{B}$ , which agrees well with the observed (reduced) mag- netic moment, µ = 0 . 47 µ$_{B}$ [15]. The discrepancy between the observed and calculated magnetic moment suggests the presence of quantum fluctuations, S − Δ S ≈ 0 . 9, as is to be expected given the small value of S and reduced dimension- ality originating from the frustrated lattice. This reduction in the ordered moment in turn contributes to a reduction in the molecular mean field over what would be expected for the full spin value.

B. Interion coupling

We now discuss the interion interactions present in MgV$_{2}$O$_{4}$. In the high-temperature cubic phase, the V 3+ ions lie on an ideal pyrochlore lattice, with each nearest-neighbor V-V bond equivalent by symmetry. In the case of antiferro- magnetic nearest-neighbor interactions, the spin arrangement on a pyrochlore lattice is geometrically frustrated [45], open- ing the possibility of spin liquid states [37,46] resulting from the large ground-state degeneracy. Pyrochlore systems with S = 1 have been studied both theoretically [47,48] and exper- imentally [49–51], particularly in search of an understanding of the crossover from quantum to classical spin liquid behav- ior [48,51]. The vanadium spinels have also attracted great interest [17,18,52–54] since the orbital degree of freedom al- lows for Jahn-Teller distortions, which break the ground-state degeneracy and permit magnetic order. This behavior is ob- served in MgV$_{2}$O$_{4}$ at T$_{S}$ ≃ 60 K, where the system becomes tetragonal before ordering antiferromagnetically at T$_{N}$ ≃ 40 K [15]. The cubic to tetragonal structural transition renders the nearest-neighbor bonds along [110] and [101] inequivalent (Fig. 6).

The orbital order plays an important role in the magnetic interion exchange since both the d - d and d - p orbital overlap are dependent on the V 3+ orbital ground state, and ultimately where the electron hole lies. This provides a mechanism for a

breaking of the equivalence of magnetic bonds both in terms of the amplitude of the direct and super-exchange interactions [43]. The question of the magnetic exchange in the vanadium spinels has been discussed at length by Di Matteo et al. [52] and Perkins et al. [18]. We will now briefly summarize the findings of the aforementioned references [18,52] insofar as they are relevant to the magnetic exchange in MgV$_{2}$O$_{4}$.

Di Matteo et al. considered the limit where the tetragonal distortion is dominant (ROO, real electronic orbitals discussed above) and where the spin-orbit coupling is dominant (COO, complex basis described above when the crystalline electric field Dq is weak) [52]. Both limits give rise to multiple mag- netic and orbital configurations dependent on the underlying physical parameters [52] and described by the general Heisen- berg spin Hamiltonian

$$
\mathcal{H}_{\text{exch}} = \frac{1}{2} \sum_{ij} \mathcal{J}_{ij} \hat{S}_i \cdot \hat{S}_j,
$$

where the factor of 1 2 accounts for double counting. In the case of ROO, only one possible orbital ground state is consistent with the experimentally observed magnetic structure [15]. This ROO ground state comprises strong antiferromagnetic (AFM) bonds in the xy plane and a weaker ferromagnetic cou- pling along [101] and [011]. The ROO model thus suggests that MgV$_{2}$O$_{4}$ is comprised of strongly-coupled AFM chains in the xy plane with weak FM coupling between chains. In each tetrahedron, two FM bonds are satisfied and two FM bonds are frustrated.

In contrast, when the spin-orbit coupling is dominant (COO), there are two ground states consistent with the observed magnetic structure. Both possible orbital ground states have strongly coupled AFM chains in the xy plane; however, one of these states has two weak AFM and two weak FM bonds per tetrahedron—thus no bond is frus- trated. In the other possible COO ground state, the interchain bonds are weak AFM bonds, two of which are frustrated. Based on the parameters obtained from spectroscopy data for vanadium perovskites [55,56] and typical strengths of the spin-orbit coupling for V 3+ reported in the literature [17,32], Di Matteo et al. concluded that the latter of the two COO states is likely the ground state [52].

Finally we note that the predicted magnetic ground states in the ROO and COO limits differ only in the sign of the interchain coupling, with both interchain couplings of a sim- ilar strength [52]. It is interesting to speculate as to whether these two limits are continuously connected via a phase with zero interchain coupling and hence might be tuned through with the application of strain, using the tetragonal distortion Γ as a tuning parameter. It is also interesting to note that a metal-insulator has been predicted in MgV$_{2}$O$_{4}$ based on ab initio calculations [57].

C. Neutron scattering intensity calculation

In order to calculate the neutron scattering response, we employ a spin-orbital exciton model in terms of Green’s func- tions as applied previously to describe multilevel systems [19–21,41,58,59]. A full derivation for a general single-Q multilevel magnetic system is presented in Ref. [59], here we
<table><caption>TABLE II. Summary of labeling convention for indices.</caption><thead><tr><th>Index</th><th>Description</th></tr></thead><tbody><tr><td>γ, γ'</td><td>V<sup>3+</sup> sites</td></tr><tr><td>p, q</td><td>spin-orbital crystal-field states</td></tr><tr><td>α, β, μ, ν</td><td>Cartesian coordinates</td></tr></tbody></table>

quote only the key results. The neutron scattering intensity is proportional to the structure factor

$$S(\mathbf{q}, \omega) = g_L^2 f^2(\mathbf{q}) \sum_{\alpha\beta} (\delta_{\alpha\beta} - \hat{q}_\alpha \hat{q}_\beta) S^{\alpha\beta}(\mathbf{q}, \omega),$$

where $g_L$ is the Landé g factor, the V$^{3+}$ form factor is $f(\mathbf{q})$, and $S^{\alpha\beta}(\mathbf{q}, \omega)$ is the partial dynamical structure factor (with Cartesian indices $\alpha, \beta$, Table II), defined as

$$S^{\alpha\beta}(\mathbf{q}, \omega) = \frac{1}{2\pi} \int dt e^{i\omega t} \langle \hat{S}^{\alpha}(\mathbf{q}, t) \hat{S}^{\beta}(-\mathbf{q}, 0) \rangle.$$

The factor preceding the partial dynamical structure factor is the polarization factor, which picks out the component of the structure factor perpendicular to the scattering wavevector, to which the scattered neutrons are sensitive. We note that neutron scattering is sensitive to magnetic correlations through interacting with the spin. It is through the presence of spin-orbit coupling in the magnetic Hamiltonian described above, neutrons are sensitive to orbital effects. The structure factor can be related to the Green's function by way of the fluctuation-dissipation theorem,

$$S^{\alpha\beta}(\mathbf{q}, \omega) = -\frac{1}{\pi} \frac{1}{1 - \exp(-\omega/k_B T)} \Re G^{\alpha\beta}(\mathbf{q}, \omega).$$

The Green's function for multilevel spin-systems can be written as a Dyson equation where the propagator describes the dynamics of the single ion at mean-field level with the interion interaction determining the self energy $\underline{\mathcal{J}}(\mathbf{q})$. The single-ion Green's function can be written as

$$g_{\tilde{\gamma}\tilde{\gamma}'}^{\alpha\beta}(\omega) = \sum_{qp} \frac{S_{\alphaqp}^{\tilde{\gamma}} S_{\beta pq}^{\tilde{\gamma}'}}{\omega - (\omega_p - \omega_q)}, \quad (5)$$

where $\omega_p$ is the single-ion eigenvalue of state $|p\rangle$ and $S_{\alphaqp}^{\tilde{\gamma}} = \langle p|\hat{S}_{\alpha}^{\tilde{\gamma}}|q\rangle$. $\phi_{qp} = (f_q - f_p)$ where $f_q$ is the Bose occupation factor of level $q$. The site indices are denoted as $\tilde{\gamma}, \tilde{\gamma}'$ (Table II). We perform a coordinate transformation onto the rotating frame so that the molecular mean field is nonoscillatory. In the rotating frame, the single-ion of each of the magnetic ions is independent of the moment direction and varies only between crystallographically inequivalent sites. In the rotating frame the Green's function is defined as

$$\tilde{G}_{\tilde{\gamma}\tilde{\gamma}'}^{\alpha\beta}(\mathbf{q}, \omega) = g_{\tilde{\gamma}\tilde{\gamma}'}^{\alpha\beta}(\omega) \delta_{\tilde{\gamma}\tilde{\gamma}'} + \sum_{\gamma'}^{\mu\nu} \tilde{J}_{\tilde{\gamma}\tilde{\gamma}'}^{\mu\nu}(\mathbf{q}) g_{\tilde{\gamma}\tilde{\gamma}'}^{\alpha\mu}(\omega) \tilde{G}_{\gamma'\gamma'}^{\nu\beta}(\mathbf{q}, \omega), \quad (6)$$

which can be solved as a matrix equation. Equation (6) contains the Fourier transformed exchange coupling in the rotating frame $\tilde{J}_{\tilde{\gamma}\tilde{\gamma}'}^{\mu\nu}(\mathbf{q})$, which takes the form of a matrix of dimension $3N \times 3N$, where $N$ is the number of sites in the

<table><caption>TABLE III. V$^{3+}$ ion positions in the unit cell.</caption><thead><tr><th>Index</th><th>Moment orientation</th><th>Position vector</th></tr></thead><tbody><tr><td>1</td><td>↑</td><td>(0.125,0.625,0.125)</td></tr><tr><td>2</td><td>↓</td><td>(0.375,0.875,0.125)</td></tr><tr><td>3</td><td>↑</td><td>(0.375,0.625,0.375)</td></tr><tr><td>4</td><td>↓</td><td>(0.125,0.875,0.375)</td></tr><tr><td>5</td><td>↑</td><td>(0.375,0.375,0.625)</td></tr><tr><td>6</td><td>↓</td><td>(0.625,0.625,0.625)</td></tr><tr><td>7</td><td>↓</td><td>(0.375,0.125,0.875)</td></tr><tr><td>8</td><td>↑</td><td>(0.125,0.375,0.875)</td></tr></tbody></table>

unit cell. At this point it is useful to examine the structure of MgV$_2$O$_4$. Below $T_N \simeq 40$ K, long-range antiferromagnetic magnetic order is established, with $\mathbf{Q} = (0, 0, 1)$ [15]. The full crystallographic and magnetic structure can be described by an eight-site unit cell as considered in Ref. [18].

In terms of this unit cell (Table III), no rotation between neighboring unit cells is required and, assuming Heisenberg coupling, in the rotating frame one has

$$\tilde{J}_{\tilde{\gamma}\tilde{\gamma}'}^{\mu\nu}(\mathbf{q}) = X'(\tilde{J}_{\tilde{\gamma}\tilde{\gamma}'}^{\mu\nu}(\mathbf{q}) \otimes \mathbb{I}_3)X \quad (7)$$

where $X = \text{diag}(R_1, ..., R_n)$ rotates the spins within the unit cell. The matrix $R_n$ is the $3 \times 3$ rotation matrix, which rotates the spin on site $n$ onto the $\hat{z}$ axis. Similarly $X' = \text{diag}(R_1^T, ..., R_n^T)$. Finally, the Green's function in the laboratory frame is found by rotating the Green's function back using the matrices, $X, X'$,

$$G_{\tilde{\gamma}\tilde{\gamma}'}^{\alpha\beta}(\mathbf{q}, \omega) = X \tilde{G}_{\tilde{\gamma}\tilde{\gamma}'}^{\alpha\beta}(\mathbf{q}, \omega) X'$$

We now apply this theory to the neutron scattering data presented in Sec. III.

### D. Application to MgV$_2$O$_4$

The excitonic dispersion relation was found according to the model presented in the previous section by finding the poles of the Green's function $G_{\tilde{\gamma}\tilde{\gamma}'}^{\alpha\beta}(\mathbf{q}, \omega)$. In Fig. 5(d), we illustrate dipole allowed transitions through the red lines when investigating the effect of $\mathcal{H}_{MF}$. Based on this, we would expect two branches with a lower strongly dispersive branch corresponding to transitions within the ground-state $j_{eff} = 2$ manifold and another transition from the ground state to a $j_{eff} = 1$ level. We might expect that given that this is a “crystal-field-like” transition, this second higher energy and dipole allowed transition will be weakly dispersive in momentum akin to crystal-field transitions commonly observed for rare-earth ions.

The experimentally observed dispersion was found both along the V-V chain direction and perpendicular to the chain by fitting Gaussian peaks to one-dimensional constant energy cuts through the neutron scattering data. These data points were then fitted using the dispersion derived according to the excitonic model in order to extract physical parameters. The data in comparison to theory is summarized through constant momentum and energy slices in Figs. 7–9. The excitonic theory qualitatively shows two sets of modes with
<table><caption>TABLE II. Summary of labeling convention for indices.</caption><thead><tr><th>Index</th><th>Description</th></tr></thead><tbody><tr><td>γ, γ'</td><td>V<sup>3+</sup> sites</td></tr><tr><td>p, q</td><td>spin-orbital crystal-field states</td></tr><tr><td>α, β, μ, ν</td><td>Cartesian coordinates</td></tr></tbody></table>

quote only the key results. The neutron scattering intensity is proportional to the structure factor

$$S(\mathbf{q}, \omega) = g_L^2 f^2(\mathbf{q}) \sum_{\alpha\beta} (\delta_{\alpha\beta} - \hat{q}_\alpha \hat{q}_\beta) S^{\alpha\beta}(\mathbf{q}, \omega),$$

where $g_L$ is the Landé g factor, the V$^{3+}$ form factor is $f(\mathbf{q})$, and $S^{\alpha\beta}(\mathbf{q}, \omega)$ is the partial dynamical structure factor (with Cartesian indices $\alpha, \beta$, Table II), defined as

$$S^{\alpha\beta}(\mathbf{q}, \omega) = \frac{1}{2\pi} \int dt e^{i\omega t} \langle \hat{S}^{\alpha}(\mathbf{q}, t) \hat{S}^{\beta}(-\mathbf{q}, 0) \rangle.$$

The factor preceding the partial dynamical structure factor is the polarization factor, which picks out the component of the structure factor perpendicular to the scattering wavevector, to which the scattered neutrons are sensitive. We note that neutron scattering is sensitive to magnetic correlations through interacting with the spin. It is through the presence of spin-orbit coupling in the magnetic Hamiltonian described above, neutrons are sensitive to orbital effects. The structure factor can be related to the Green's function by way of the fluctuation-dissipation theorem,

$$S^{\alpha\beta}(\mathbf{q}, \omega) = -\frac{1}{\pi} \frac{1}{1 - \exp(-\omega/k_B T)} \Re G^{\alpha\beta}(\mathbf{q}, \omega).$$

The Green's function for multilevel spin-systems can be written as a Dyson equation where the propagator describes the dynamics of the single ion at mean-field level with the interion interaction determining the self energy $\underline{\mathcal{J}}(\mathbf{q})$. The single-ion Green's function can be written as

$$g_{\tilde{\gamma}\tilde{\gamma}'}^{\alpha\beta}(\omega) = \sum_{qp} \frac{S_{\alphaqp}^{\tilde{\gamma}} S_{\beta pq}^{\tilde{\gamma}'}}{\omega - (\omega_p - \omega_q)}, \quad (5)$$

where $\omega_p$ is the single-ion eigenvalue of state $|p\rangle$ and $S_{\alphaqp}^{\tilde{\gamma}} = \langle p|\hat{S}_{\alpha}^{\tilde{\gamma}}|q\rangle$. $\phi_{qp} = (f_q - f_p)$ where $f_q$ is the Bose occupation factor of level $q$. The site indices are denoted as $\tilde{\gamma}, \tilde{\gamma}'$ (Table II). We perform a coordinate transformation onto the rotating frame so that the molecular mean field is nonoscillatory. In the rotating frame, the single-ion of each of the magnetic ions is independent of the moment direction and varies only between crystallographically inequivalent sites. In the rotating frame the Green's function is defined as

$$\tilde{G}_{\tilde{\gamma}\tilde{\gamma}'}^{\alpha\beta}(\mathbf{q}, \omega) = g_{\tilde{\gamma}\tilde{\gamma}'}^{\alpha\beta}(\omega) \delta_{\tilde{\gamma}\tilde{\gamma}'} + \sum_{\gamma'}^{\mu\nu} \tilde{J}_{\tilde{\gamma}\tilde{\gamma}'}^{\mu\nu}(\mathbf{q}) g_{\tilde{\gamma}\tilde{\gamma}'}^{\alpha\mu}(\omega) \tilde{G}_{\gamma'\gamma'}^{\nu\beta}(\mathbf{q}, \omega), \quad (6)$$

which can be solved as a matrix equation. Equation (6) contains the Fourier transformed exchange coupling in the rotating frame $\tilde{J}_{\tilde{\gamma}\tilde{\gamma}'}^{\mu\nu}(\mathbf{q})$, which takes the form of a matrix of dimension $3N \times 3N$, where $N$ is the number of sites in the

<table><caption>TABLE III. V$^{3+}$ ion positions in the unit cell.</caption><thead><tr><th>Index</th><th>Moment orientation</th><th>Position vector</th></tr></thead><tbody><tr><td>1</td><td>↑</td><td>(0.125,0.625,0.125)</td></tr><tr><td>2</td><td>↓</td><td>(0.375,0.875,0.125)</td></tr><tr><td>3</td><td>↑</td><td>(0.375,0.625,0.375)</td></tr><tr><td>4</td><td>↓</td><td>(0.125,0.875,0.375)</td></tr><tr><td>5</td><td>↑</td><td>(0.375,0.375,0.625)</td></tr><tr><td>6</td><td>↓</td><td>(0.625,0.625,0.625)</td></tr><tr><td>7</td><td>↓</td><td>(0.375,0.125,0.875)</td></tr><tr><td>8</td><td>↑</td><td>(0.125,0.375,0.875)</td></tr></tbody></table>

unit cell. At this point it is useful to examine the structure of MgV$_2$O$_4$. Below $T_N \simeq 40$ K, long-range antiferromagnetic magnetic order is established, with $\mathbf{Q} = (0, 0, 1)$ [15]. The full crystallographic and magnetic structure can be described by an eight-site unit cell as considered in Ref. [18].

In terms of this unit cell (Table III), no rotation between neighboring unit cells is required and, assuming Heisenberg coupling, in the rotating frame one has

$$\tilde{J}_{\tilde{\gamma}\tilde{\gamma}'}^{\mu\nu}(\mathbf{q}) = X'(\tilde{J}_{\tilde{\gamma}\tilde{\gamma}'}^{\mu\nu}(\mathbf{q}) \otimes \mathbb{I}_3)X \quad (7)$$

where $X = \text{diag}(R_1, \dots, R_n)$ rotates the spins within the unit cell. The matrix $R_n$ is the $3 \times 3$ rotation matrix, which rotates the spin on site $n$ onto the $\hat{z}$ axis. Similarly $X' = \text{diag}(R_1^T, \dots, R_n^T)$. Finally, the Green's function in the laboratory frame is found by rotating the Green's function back using the matrices, $X, X'$,

$$G_{\tilde{\gamma}\tilde{\gamma}'}^{\alpha\beta}(\mathbf{q}, \omega) = X \tilde{G}_{\tilde{\gamma}\tilde{\gamma}'}^{\alpha\beta}(\mathbf{q}, \omega) X'$$

We now apply this theory to the neutron scattering data presented in Sec. III.

### D. Application to MgV$_2$O$_4$

The excitonic dispersion relation was found according to the model presented in the previous section by finding the poles of the Green's function $G_{\tilde{\gamma}\tilde{\gamma}'}^{\alpha\beta}(\mathbf{q}, \omega)$. In Fig. 5(d), we illustrate dipole allowed transitions through the red lines when investigating the effect of $\mathcal{H}_{MF}$. Based on this, we would expect two branches with a lower strongly dispersive branch corresponding to transitions within the ground-state $j_{\text{eff}} = 2$ manifold and another transition from the ground state to a $j_{\text{eff}} = 1$ level. We might expect that given that this is a “crystal-field-like” transition, this second higher energy and dipole allowed transition will be weakly dispersive in momentum akin to crystal-field transitions commonly observed for rare-earth ions.

The experimentally observed dispersion was found both along the V-V chain direction and perpendicular to the chain by fitting Gaussian peaks to one-dimensional constant energy cuts through the neutron scattering data. These data points were then fitted using the dispersion derived according to the excitonic model in order to extract physical parameters. The data in comparison to theory is summarized through constant momentum and energy slices in Figs. 7–9. The excitonic theory qualitatively shows two sets of modes with
<figure><img src="image_9.png" /><figcaption>FIG. 7. Calculation of the component of the structure factor (a) parallel and (b) perpendicular to the chain direction. Overplotted is the extracted dispersion from neutron scattering. The green data points are from the EIGER data set and the red and blue points are from the MERLIN $E_i = 24$ meV and 49 meV data sets respectively. The white solid line indicates the approximate peak position of the weak in intensity high-energy mode and the dashed lines indicate the approximate width.</figcaption></figure>

low-energy fluctuations corresponding to transitions within the ground-state $j_{\text{eff}} = 2$ manifold and another very weak, and comparatively disperisonless, high-energy mode corresponding to transitions to the excited $j_{\text{eff}} = 1$ spin-orbit manifold. We have plotted the weak mode on the same linear intensity scale in Fig. 7 to illustrate the strong difference in intensity between these two branches of excitations. We note that the excitonic theory predicts a very weak and comparatively flat mode at ~10 meV. Due to the presence of overlapping domains and a comparable magnetic zone center excitation gap, we were not able to reliably resolve this mode from our data.

The interchain bond $J_1$ was determined to be negligible since the inclusion of this bond (in the presence of finite $J_3$) gives rise to a splitting of the lower dispersive mode at

approximately $Q \geq \pm 0.1$ (r.l.u), which is not seen in the data. We therefore determine the magnitude of $J_1$ to be negligible given the absence of any of these signatures, but note that frustrated nature of this bond makes its determination difficult experimentally. The neglect of $J_1$ is further justified given that previous calculations have suggested that this bond is weakly ferromagnetic or weakly antiferromagnetic in the limits of ROO and COO respectively [52]. In fact, the splitting caused by $J_1$ is a hybridization with the lifted zero-energy mode (which acquires some dispersion due to finite $J_3$) and the low-energy mode due to the intrachain coupling. As pointed out by [18], the presence of a zero-energy mode corresponds to the degeneracy associated with the rotation of the spins of the chains in the $xy$ plane. This degeneracy is lifted by the anisotropy provided by spin-orbit coupling and crystallo-

<figure><img src="image_10.png" /><figcaption>FIG. 8. (a) Constant energy slice through the intermultiplet mode at $E \sim 50$ meV. (b) One-dimensional cut through $Q = (2, 0)$, with multilevel spin wave calculation plotted in red. The Lorentzian half-width has been chosen to approximately match the energy resolution on MAPS, $\epsilon = 1.5$. The calculation has been integrated along the $c^*$ direction across the full width of the Brillouin zone. The presence of domain coexistence and the coupling of energy transfer and the wavevector along $c^*$ have been neglected. The differences between calculations and data at low energies originates from an oversubtraction of the data with the chosen background.</figcaption></figure>

<figure><img src="image_11.png" /><figcaption>FIG. 9. Constant energy slice at $10 \pm 1$ meV, showing the comparison between the measured data from MERLIN and the excitonic model with all three domains summed as described in the text. The intensities have been integrated along the $c^*$ direction ($L = -1 \pm 0.2$ r.l.u).</figcaption></figure>
<figure><img src="image_9.png" /><figcaption>FIG. 7. Calculation of the component of the structure factor (a) parallel and (b) perpendicular to the chain direction. Overplotted is the extracted dispersion from neutron scattering. The green data points are from the EIGER data set and the red and blue points are from the MERLIN $E_i = 24$ meV and 49 meV data sets respectively. The white solid line indicates the approximate peak position of the weak in intensity high-energy mode and the dashed lines indicate the approximate width.</figcaption></figure>

low-energy fluctuations corresponding to transitions within the ground-state $j_{\text{eff}} = 2$ manifold and another very weak, and comparatively disperisonless, high-energy mode corresponding to transitions to the excited $j_{\text{eff}} = 1$ spin-orbit manifold. We have plotted the weak mode on the same linear intensity scale in Fig. 7 to illustrate the strong difference in intensity between these two branches of excitations. We note that the excitonic theory predicts a very weak and comparatively flat mode at ~10 meV. Due to the presence of overlapping domains and a comparable magnetic zone center excitation gap, we were not able to reliably resolve this mode from our data.

The interchain bond $J_1$ was determined to be negligible since the inclusion of this bond (in the presence of finite $J_3$) gives rise to a splitting of the lower dispersive mode at

approximately $Q \geq \pm 0.1$ (r.l.u), which is not seen in the data. We therefore determine the magnitude of $J_1$ to be negligible given the absence of any of these signatures, but note that frustrated nature of this bond makes its determination difficult experimentally. The neglect of $J_1$ is further justified given that previous calculations have suggested that this bond is weakly ferromagnetic or weakly antiferromagnetic in the limits of ROO and COO respectively [52]. In fact, the splitting caused by $J_1$ is a hybridization with the lifted zero-energy mode (which acquires some dispersion due to finite $J_3$) and the low-energy mode due to the intrachain coupling. As pointed out by [18], the presence of a zero-energy mode corresponds to the degeneracy associated with the rotation of the spins of the chains in the $xy$ plane. This degeneracy is lifted by the anisotropy provided by spin-orbit coupling and crystallo-

<figure><img src="image_10.png" /><figcaption>FIG. 8. (a) Constant energy slice through the intermultiplet mode at $E \sim 50$ meV. (b) One-dimensional cut through $Q = (2, 0)$, with multilevel spin wave calculation plotted in red. The Lorentzian half-width has been chosen to approximately match the energy resolution on MAPS, $\epsilon = 1.5$. The calculation has been integrated along the $c^*$ direction across the full width of the Brillouin zone. The presence of domain coexistence and the coupling of energy transfer and the wavevector along $c^*$ have been neglected. The differences between calculations and data at low energies originates from an oversubtraction of the data with the chosen background.</figcaption></figure>

<figure><img src="image_11.png" /><figcaption>FIG. 9. Constant energy slice at $10 \pm 1$ meV, showing the comparison between the measured data from MERLIN and the excitonic model with all three domains summed as described in the text. The intensities have been integrated along the $c^*$ direction ($L = -1 \pm 0.2$ r.l.u).</figcaption></figure>
<table><caption>TABLE IV. Spin-orbit exciton parameters.</caption><thead><tr><th>Parameter</th><th>Bond distance (Å)</th><th>Value (meV)</th></tr></thead><tbody><tr><td>J<sub>1</sub></td><td>2.885</td><td>≈0</td></tr><tr><td>J<sub>2</sub></td><td>2.885</td><td>19.04 (±2.86)</td></tr><tr><td>J<sub>3</sub></td><td>4.997</td><td>-0.173 (±0.026)</td></tr><tr><td>αλ</td><td></td><td>-8.8 (±1.3)</td></tr><tr><td>Γ</td><td></td><td>-42.2 (±6.33)</td></tr><tr><td>Θ</td><td></td><td>0.8334 (±0.1250)</td></tr></tbody></table>

graphic distortions. A small, but finite, $J_3$ is expected based on the dispersive excitation seen in Fig. 3(b), which is consistent with coupling of parallel chains. We note that there exists a further neighbor bond $J_4$ with bond length marginally longer than $J_3$. Inclusion of this bond again splits the lower dispersive mode but does not appreciably change the bandwidth or gap, provided that it is small. In the absence of any observed splitting, and in the spirit of writing down a minimal model, we neglect this term.

For the lower-energy dispersive mode, the data were fitted by taking one-dimensional constant energy cuts through the neutron data to extract a dispersion curve to which the poles of the Green's function were fit as the parameters were varied (Fig. 7). To fit the higher-energy mode at ~50 meV (Fig. 8), a constant q cut was made at (2,0,0) and subtracted as a background. This method is justified at high-energy transfers hear the ~50 meV mode; however, it cuts through magnetic intensity at lower energies and therefore is an overestimate of the background. This is reflected in Fig. 7, which shows an overestimate of the intensity by the exciton model. The fitted parameters are listed in Table IV. Stated uncertainties represent a nominal 15% error to reflect the inherent uncertainties associated with the experimental measurements (such as spectrometer resolution) rather than the error bars associated with the numerical fitting procedure, which are unphysically small. The approximate magnitude of the intrachain coupling $J_2$ is consistent with that theoretically predicted in Ref. [18] using a Kugel-Khomskii model [60]. A reduction of the molecular mean-field strength from that expected for an $S = 1$ antiferromagnet

$$ \mathcal{H}_{MF} = h_{MF} \hat{S}_z = \Theta(-2J_2 + 4J_3), \quad (8) $$

where $\Theta < 1$, is present due to quantum fluctuations. A reduction factor of $\Theta \approx 0.9$ is expected based on the observed magnetic moment [15] and the orbital projection factor based on the free-ion values [32]. The calculated magnetic moment is sensitive to the projection factor $\alpha$, which is tied to the single-ion parameters $Dq$ and Racah parameter $B$. There is, however, a degree of uncertainty about the single-ion parameters, $B$, $C$, and $Dq$ in the crystal environment, hence $\Theta$ was refined in our analysis of the neutron spectroscopy results. The refined value of $\Theta$ implies that the true value of the orbital projection factor is reduced from that expected in a free ion. One possible origin of this is the presence of a nephelauxetic effect, which reduces the electron Coulomb repulsion and hence reduces the value of the Racah $B$ parameter. We note that reductions of ~10–20% have been reported [61],

consistent with the reduction that would be required to explain the value of $\Theta$ fitted here.

The calculated structure factor is plotted in Figs. 7(a)–7(c). Overplotted are the extracted peak positions from the neutron scattering data taken on EIGER and MERLIN. The calculation well describes both the dispersion and the intensity variation throughout the Brillouin zone.

In the experimental section above it was noted that single-crystal samples should be expected to exhibit domain coexistence as was observed in Ref. [16]. This occurs below the structural transition and results from the overall crystal structure retaining the average high-temperature cubic symmetry in the absence of strain. We now briefly note the consequences of domain coexistence on the neutron scattering spectrum. A rotation matrix can be defined for each domain, which transforms between the crystallographic basis vectors of each domain $d$,

$$ \underline{R}_d \hat{Q}_0 = \hat{Q}_d \quad (9) $$

with the $\hat{Q}_0$ defining the nominal reference basis plotted in Fig. 5. The crystallographic structure of $MgV_2O_4$ comprises chains along [110] and so one can define the rotation matrices

$$ \underline{R}_0 = \begin{pmatrix} 1 & 0 & 0 \\ 0 & 1 & 0 \\ 0 & 0 & 1 \end{pmatrix}, \quad (10a) $$

$$ \underline{R}_1 = \begin{pmatrix} 1 & 0 & 0 \\ 0 & 0 & 1 \\ 0 & 1 & 0 \end{pmatrix}, \quad (10b) $$

$$ \underline{R}_2 = \begin{pmatrix} 0 & 0 & 1 \\ 0 & 1 & 0 \\ 1 & 0 & 0 \end{pmatrix}, \quad (10c) $$

corresponding to three inequivalent domains. Figure 9 shows a constant energy slice with the scattering for all three domains superposed. The diagonal modes originate from the domains with basis vectors $\hat{Q}_1$ and $\hat{Q}_2$ with the horizontal and vertical lines originating from the $\hat{Q}_0$ domain. A small difference in the position of the diagonal modes likely originates from the fact that the unit cell is tetragonal with $c/a = 0.994$ [15]. In the calculation, the domain population has been assumed to be equal; however, the relative intensities of the streaks in the data indicate that this is not the case. This is in agreement with previous diffraction data [16].

### E. Orbital order

There has been a great deal of intrigue in the nature of the orbital order in $MgV_2O_4$. Much of this has centered around whether the system exhibits real or complex orbital order. We will briefly summarize both orbital-ordering schemes, remarking on how our results influence this discussion.

We begin by considering the isolated single ion in an octahedral environment, neglecting the interion coupling. An orbitally degenerate crystal-field ground state can lower its energy via an octahedral distortion [62] (the Jahn-Teller effect). An elongation of the octahedron lowers the $d_{xz}$ and $d_{yz}$ orbitals by $\frac{1}{2}E_{JT}$ and raises $d_{xy}$ by $E_{JT}$ [43]. The lower-energy doublet is then split by spin-orbital coupling into $|l_z = \pm 1\rangle = \frac{1}{\sqrt{2}}(d_{xz} \pm id_{yz})$ by $\pm \frac{1}{2}\alpha\lambda$. For a $d^2$ ion, assuming the spin-orbit
<table><caption>TABLE IV. Spin-orbit exciton parameters.</caption><thead><tr><th>Parameter</th><th>Bond distance (Å)</th><th>Value (meV)</th></tr></thead><tbody><tr><td>J<sub>1</sub></td><td>2.885</td><td>≈0</td></tr><tr><td>J<sub>2</sub></td><td>2.885</td><td>19.04 (±2.86)</td></tr><tr><td>J<sub>3</sub></td><td>4.997</td><td>-0.173 (±0.026)</td></tr><tr><td>αλ</td><td></td><td>-8.8 (±1.3)</td></tr><tr><td>Γ</td><td></td><td>-42.2 (±6.33)</td></tr><tr><td>Θ</td><td></td><td>0.8334 (±0.1250)</td></tr></tbody></table>

graphic distortions. A small, but finite, $J_3$ is expected based on the dispersive excitation seen in Fig. 3(b), which is consistent with coupling of parallel chains. We note that there exists a further neighbor bond $J_4$ with bond length marginally longer than $J_3$. Inclusion of this bond again splits the lower dispersive mode but does not appreciably change the bandwidth or gap, provided that it is small. In the absence of any observed splitting, and in the spirit of writing down a minimal model, we neglect this term.

For the lower-energy dispersive mode, the data were fitted by taking one-dimensional constant energy cuts through the neutron data to extract a dispersion curve to which the poles of the Green's function were fit as the parameters were varied (Fig. 7). To fit the higher-energy mode at ~50 meV (Fig. 8), a constant q cut was made at (2,0,0) and subtracted as a background. This method is justified at high-energy transfers hear the ~50 meV mode; however, it cuts through magnetic intensity at lower energies and therefore is an overestimate of the background. This is reflected in Fig. 7, which shows an overestimate of the intensity by the exciton model. The fitted parameters are listed in Table IV. Stated uncertainties represent a nominal 15% error to reflect the inherent uncertainties associated with the experimental measurements (such as spectrometer resolution) rather than the error bars associated with the numerical fitting procedure, which are unphysically small. The approximate magnitude of the intrachain coupling $J_2$ is consistent with that theoretically predicted in Ref. [18] using a Kugel-Khomskii model [60]. A reduction of the molecular mean-field strength from that expected for an $S = 1$ antiferromagnet

$$ \mathcal{H}_{MF} = h_{MF} \hat{S}_z = \Theta(-2J_2 + 4J_3), \quad (8) $$

where $\Theta < 1$, is present due to quantum fluctuations. A reduction factor of $\Theta \approx 0.9$ is expected based on the observed magnetic moment [15] and the orbital projection factor based on the free-ion values [32]. The calculated magnetic moment is sensitive to the projection factor $\alpha$, which is tied to the single-ion parameters $Dq$ and Racah parameter $B$. There is, however, a degree of uncertainty about the single-ion parameters, $B$, $C$, and $Dq$ in the crystal environment, hence $\Theta$ was refined in our analysis of the neutron spectroscopy results. The refined value of $\Theta$ implies that the true value of the orbital projection factor is reduced from that expected in a free ion. One possible origin of this is the presence of a nephelauxetic effect, which reduces the electron Coulomb repulsion and hence reduces the value of the Racah $B$ parameter. We note that reductions of ~10–20% have been reported [61],

consistent with the reduction that would be required to explain the value of $\Theta$ fitted here.

The calculated structure factor is plotted in Figs. 7(a)–7(c). Overplotted are the extracted peak positions from the neutron scattering data taken on EIGER and MERLIN. The calculation well describes both the dispersion and the intensity variation throughout the Brillouin zone.

In the experimental section above it was noted that single-crystal samples should be expected to exhibit domain coexistence as was observed in Ref. [16]. This occurs below the structural transition and results from the overall crystal structure retaining the average high-temperature cubic symmetry in the absence of strain. We now briefly note the consequences of domain coexistence on the neutron scattering spectrum. A rotation matrix can be defined for each domain, which transforms between the crystallographic basis vectors of each domain $d$,

$$ \underline{R}_d \hat{Q}_0 = \hat{Q}_d \quad (9) $$

with the $\hat{Q}_0$ defining the nominal reference basis plotted in Fig. 5. The crystallographic structure of $MgV_2O_4$ comprises chains along [110] and so one can define the rotation matrices

$$ \underline{R}_0 = \begin{pmatrix} 1 & 0 & 0 \\ 0 & 1 & 0 \\ 0 & 0 & 1 \end{pmatrix}, \quad (10a) $$

$$ \underline{R}_1 = \begin{pmatrix} 1 & 0 & 0 \\ 0 & 0 & 1 \\ 0 & 1 & 0 \end{pmatrix}, \quad (10b) $$

$$ \underline{R}_2 = \begin{pmatrix} 0 & 0 & 1 \\ 0 & 1 & 0 \\ 1 & 0 & 0 \end{pmatrix}, \quad (10c) $$

corresponding to three inequivalent domains. Figure 9 shows a constant energy slice with the scattering for all three domains superposed. The diagonal modes originate from the domains with basis vectors $\hat{Q}_1$ and $\hat{Q}_2$ with the horizontal and vertical lines originating from the $\hat{Q}_0$ domain. A small difference in the position of the diagonal modes likely originates from the fact that the unit cell is tetragonal with $c/a = 0.994$ [15]. In the calculation, the domain population has been assumed to be equal; however, the relative intensities of the streaks in the data indicate that this is not the case. This is in agreement with previous diffraction data [16].

### E. Orbital order

There has been a great deal of intrigue in the nature of the orbital order in $MgV_2O_4$. Much of this has centered around whether the system exhibits real or complex orbital order. We will briefly summarize both orbital-ordering schemes, remarking on how our results influence this discussion.

We begin by considering the isolated single ion in an octahedral environment, neglecting the interion coupling. An orbitally degenerate crystal-field ground state can lower its energy via an octahedral distortion [62] (the Jahn-Teller effect). An elongation of the octahedron lowers the $d_{xz}$ and $d_{yz}$ orbitals by $\frac{1}{2}E_{JT}$ and raises $d_{xy}$ by $E_{JT}$ [43]. The lower-energy doublet is then split by spin-orbital coupling into $|l_z = \pm 1\rangle = \frac{1}{\sqrt{2}}(d_{xz} \pm id_{yz})$ by $\pm \frac{1}{2}\alpha\lambda$. For a $d^2$ ion, assuming the spin-orbit
coupling does not lift $|l_z = -1\rangle$ above the $d_{xy}$ orbital, then both $|l_z = -1\rangle$ and $|l_z = 1\rangle$ are occupied and the resultant ground state has quenched orbital angular momentum.

Conversely, if the octahedron is compressed, the Jahn-Teller spectrum is inverted and the $d_{xy}$ level is lowered, with the $d_{xz}$ and $d_{yz}$ levels lying at $+\frac{1}{2}E_{JT}$. The ground state thus has one electron in the $d_{xy}$ orbital and one in the $|l_z = +1\rangle$ level. In this case the orbital angular momentum is unquenched owing to the orbital degeneracy, which is broken by spin-orbit coupling, lowering the ground-state energy. Based on single-ion physics alone, an elongation is favored if $|\alpha\lambda| > E_{JT} = \frac{2}{3}\Gamma$. From the fitted values in Table IV it is clear that the compressed structure is energetically favored by the comparatively large value of $\Gamma$, hence the parameters extracted from the neutron scattering measurements are consistent with the refined crystallographic structure.

We now turn to the effect of orbital order on the interion exchange in MgV$_2$O$_4$. Di Matteo et al. and Perkins et al. discussed this in Refs. [18,52] by way of a Kugel-Khomskii model [60] to model the superexchange. Here we only quote the key results. Two states were found to be consistent with the observed magnetic and crystallographic structure. The energy of the COO state is $E_{COO} = -\tilde{J}_1 - \frac{5}{2}\tilde{J}_2 - \lambda$, where $\tilde{J}_1 = J(1-\eta)/(1-3\eta)$, $\tilde{J}_2 = J(1+\eta)/(1+2\eta)$, $J = t^2/U_1$, and $\eta = J_H/U_1$. In terms of these variables, $J_2 = \tilde{J}_2$ and, assuming we have COO, $J_1 = \frac{1}{4}(\tilde{J}_2 - 2\tilde{J}_0)$, with $\tilde{J}_0 = \eta J/(1-3\eta)$. We thus find that our extracted exchange parameters suggest $\eta \approx 0.19$ and $J \approx 22.1$ meV, which are close to values suggested by photoemission spectroscopy ($\eta_{exp} = 0.11$, $J_{exp} = 20.4$ meV) [18,63]. In contrast, for ROO one has $J_1 = -\tilde{J}_0$, such that for $J_1 \approx 0$ one expects $\eta \approx 0$. Such a situation is clearly unphysical as it suggests either a vanishing Hund's coupling or extremely large Coulomb term. For reasonable values of $\eta \approx 0.1$ and with $J \approx 20$ meV, the ROO model predicts $|J_2| \approx 3$ meV, which would give rise to a splitting of the low-energy dispersive mode that would likely be resolvable in neutron scattering measurements.

As shown by Perkins et al. in Ref. [18], qualitative differences exist between the spectra expected for ROO and COO, namely an optical mode is present in the case of COO whereas ROO gives rise only to an acoustic mode. This optical mode finds its origin in the unquenched orbital angular momentum, which permits intermultiplet transitions between the $j_{eff} = 2$ and $j_{eff} = 1$ levels [Fig. 5(d)]. Our neutron scattering experiments confirm the existence of this optical mode at $E \sim 50$ meV, which we were able to model within a multi-level excitonic formalism. We find that the resultant single-ion ground state, consistent based on the fitted parameters (Table IV), is the $|l_z = 1, S_z = 1\rangle$ state presented in Ref. [52]. In this state, the effective orbital angular momentum aligns with the spin moment on all $V^{3+}$ sites. It should be noted that due to the negative projection factor, this corresponds to the antialignment of the spin and true orbital moments, which gives rise to the reduction in the magnetic moment measured in diffraction [15]. By treating the spin-orbit coupling, tetragonal distortion and molecular field on the same level using an excitonic model, we have demonstrated that MgV$_2$O$_4$ can be understood to behave according to the COO picture presented in Refs. [18,52] though the distortion, spin-orbit coupling and molecular mean-field strengths are

comparable. This analysis helps to explain the apparent discrepancy between the previous neutron scattering results [15] and the predicted spectrum for COO [18], which considered $|\Gamma| \ll |\lambda|$.

We note that Di Matteo et al. [52] propose several types of COO ranging from all $V^{3+}$ being in a complex orbital state and a mixture between real and complex orbital orders on different sites. In our theoretical analysis discussed above, we have considered the case that all $V^{3+}$ ions are in a complex orbital state. Qualitatively, a mixture of real and complex orbital states would decrease and likely damp the high-energy spin-orbit exciton corresponding to excitations from the ground-state $j_{eff} = 2$ manifold to a higher energy $j_{eff} = 1$ spin-orbit manifold. Given how weak this mode is this would likely not be observable. Another point against a mixture of complex and real orbitals is crystal lattice symmetry with the ground state being inconsistent with such a situation. We therefore consider here the case of full complex orbital order (COO) in MgV$_2$O$_4$.

## V. HYSTERETIC MAGNETIC CORRELATIONS

There are three different spin-orbital phases in MgV$_2$O$_4$ on cooling from high temperature. (1) At high temperatures exceeding $T_c$, the nuclear unit cell is cubic and the magnetism is paramagnetic. (2) In the temperature regime $T_N < T < T_S$, the structural unit cell distorts to tetragonal. (3) For low temperature $T < T_N$, antiferromagnetic order sets in. As displayed in Fig. 5(d), the single-ion spin-orbital ground state is different in each one of these phases transitioning from a degenerate $j_{eff} = 2$ state in the high-temperature cubic phase to an orbital doublet below $T_S$, and then finally this degeneracy being split by a molecular field in the low-temperature antiferromagnetic phase. This temperature dependence opens up the possibility of hysteretic effects on cooling and warming through the intermediate phase given the orbital degeneracy present, which is broken when cooling through the antiferromagnetic transition by a Zeeman-like molecular field. We now investigate these spin-orbital ground states using energy-integrated magnetic diffuse scattering.

To investigate the low-energy critical magnetic fluctuations as a function of temperature in MgV$_2$O$_4$, we studied the magnetic cross section using the DNS polarized diffractometer. The combined magnetic intensities from the diffuse scattering measurements are displayed in Fig. 10. We note these measurements display a hysteresis depending on if the sample is cooled below the magnetic ordering temperature, $T_N \approx 40$ K. Cooling to $T = 75$ K from $120$ K [Figs. 10(a) → 10(b)] one sees the development of magnetic diffuse scattering, consistent with a spin-frustrated pyrochlore system above the structural transition. On further cooling to $T = 50$ K [Figs. 10(b) → 10(c)], the magnetic scattering changes into chain-like rods below the structural transition. Upon heating from $T = 50$ K back to $75$ K [Figs. 10(c) → 10(b)] one recovers these short-ranged magnetic correlations. However, if the sample is cooled below the long-range magnetic transition $T_N$ and then warmed up to $T = 75$ K [Figs. 10(a) → 10(d)] a different behavior is observed. Whilst at $T = 50$ K [Figs. 10(d) → 10(e)] one still observes chain-like scattering, above the structural transition [Figs. 10(e) → 10(f)], $T_S$ the
coupling does not lift $|l_z = -1\rangle$ above the $d_{xy}$ orbital, then both $|l_z = -1\rangle$ and $|l_z = 1\rangle$ are occupied and the resultant ground state has quenched orbital angular momentum.

Conversely, if the octahedron is compressed, the Jahn-Teller spectrum is inverted and the $d_{xy}$ level is lowered, with the $d_{xz}$ and $d_{yz}$ levels lying at $+\frac{1}{2}E_{JT}$. The ground state thus has one electron in the $d_{xy}$ orbital and one in the $|l_z = +1\rangle$ level. In this case the orbital angular momentum is unquenched owing to the orbital degeneracy, which is broken by spin-orbit coupling, lowering the ground-state energy. Based on single-ion physics alone, an elongation is favored if $|\alpha\lambda| > E_{JT} = \frac{2}{3}\Gamma$. From the fitted values in Table IV it is clear that the compressed structure is energetically favored by the comparatively large value of $\Gamma$, hence the parameters extracted from the neutron scattering measurements are consistent with the refined crystallographic structure.

We now turn to the effect of orbital order on the interion exchange in MgV$_2$O$_4$. Di Matteo *et al.* and Perkins *et al.* discussed this in Refs. [18,52] by way of a Kugel-Khomskii model [60] to model the superexchange. Here we only quote the key results. Two states were found to be consistent with the observed magnetic and crystallographic structure. The energy of the COO state is $E_{COO} = -\tilde{J}_1 - \frac{5}{2}\tilde{J}_2 - \lambda$, where $\tilde{J}_1 = J(1-\eta)/(1-3\eta)$, $\tilde{J}_2 = J(1+\eta)/(1+2\eta)$, $J = t^2/U_1$, and $\eta = J_H/U_1$. In terms of these variables, $J_2 = \tilde{J}_2$ and, assuming we have COO, $J_1 = \frac{1}{4}(\tilde{J}_2 - 2\tilde{J}_0)$, with $\tilde{J}_0 = \eta J/(1-3\eta)$. We thus find that our extracted exchange parameters suggest $\eta \approx 0.19$ and $J \approx 22.1$ meV, which are close to values suggested by photoemission spectroscopy ($\eta_{exp} = 0.11$, $J_{exp} = 20.4$ meV) [18,63]. In contrast, for ROO one has $J_1 = -\tilde{J}_0$, such that for $J_1 \approx 0$ one expects $\eta \approx 0$. Such a situation is clearly unphysical as it suggests either a vanishing Hund's coupling or extremely large Coulomb term. For reasonable values of $\eta \approx 0.1$ and with $J \approx 20$ meV, the ROO model predicts $|J_2| \approx 3$ meV, which would give rise to a splitting of the low-energy dispersive mode that would likely be resolvable in neutron scattering measurements.

As shown by Perkins *et al.* in Ref. [18], qualitative differences exist between the spectra expected for ROO and COO, namely an optical mode is present in the case of COO whereas ROO gives rise only to an acoustic mode. This optical mode finds its origin in the unquenched orbital angular momentum, which permits intermultiplet transitions between the $j_{eff} = 2$ and $j_{eff} = 1$ levels [Fig. 5(d)]. Our neutron scattering experiments confirm the existence of this optical mode at $E \sim 50$ meV, which we were able to model within a multi-level excitonic formalism. We find that the resultant single-ion ground state, consistent based on the fitted parameters (Table IV), is the $|l_z = 1, S_z = 1\rangle$ state presented in Ref. [52]. In this state, the effective orbital angular momentum aligns with the spin moment on all $V^{3+}$ sites. It should be noted that due to the negative projection factor, this corresponds to the antialignment of the spin and true orbital moments, which gives rise to the reduction in the magnetic moment measured in diffraction [15]. By treating the spin-orbit coupling, tetragonal distortion and molecular field on the same level using an excitonic model, we have demonstrated that MgV$_2$O$_4$ can be understood to behave according to the COO picture presented in Refs. [18,52] though the distortion, spin-orbit coupling and molecular mean-field strengths are

comparable. This analysis helps to explain the apparent discrepancy between the previous neutron scattering results [15] and the predicted spectrum for COO [18], which considered $|\Gamma| \ll |\lambda|$.

We note that Di Matteo *et al.* [52] propose several types of COO ranging from all $V^{3+}$ being in a complex orbital state and a mixture between real and complex orbital orders on different sites. In our theoretical analysis discussed above, we have considered the case that all $V^{3+}$ ions are in a complex orbital state. Qualitatively, a mixture of real and complex orbital states would decrease and likely damp the high-energy spin-orbit exciton corresponding to excitations from the ground-state $j_{eff} = 2$ manifold to a higher energy $j_{eff} = 1$ spin-orbit manifold. Given how weak this mode is this would likely not be observable. Another point against a mixture of complex and real orbitals is crystal lattice symmetry with the ground state being inconsistent with such a situation. We therefore consider here the case of full complex orbital order (COO) in MgV$_2$O$_4$.

## V. HYSTERETIC MAGNETIC CORRELATIONS

There are three different spin-orbital phases in MgV$_2$O$_4$ on cooling from high temperature. (1) At high temperatures exceeding $T_c$, the nuclear unit cell is cubic and the magnetism is paramagnetic. (2) In the temperature regime $T_N < T < T_S$, the structural unit cell distorts to tetragonal. (3) For low temperature $T < T_N$, antiferromagnetic order sets in. As displayed in Fig. 5(d), the single-ion spin-orbital ground state is different in each one of these phases transitioning from a degenerate $j_{eff} = 2$ state in the high-temperature cubic phase to an orbital doublet below $T_S$, and then finally this degeneracy being split by a molecular field in the low-temperature antiferromagnetic phase. This temperature dependence opens up the possibility of hysteretic effects on cooling and warming through the intermediate phase given the orbital degeneracy present, which is broken when cooling through the antiferromagnetic transition by a Zeeman-like molecular field. We now investigate these spin-orbital ground states using energy-integrated magnetic diffuse scattering.

To investigate the low-energy critical magnetic fluctuations as a function of temperature in MgV$_2$O$_4$, we studied the magnetic cross section using the DNS polarized diffractometer. The combined magnetic intensities from the diffuse scattering measurements are displayed in Fig. 10. We note these measurements display a hysteresis depending on if the sample is cooled below the magnetic ordering temperature, $T_N \approx 40$ K. Cooling to $T = 75$ K from $120$ K [Figs. 10(a) → 10(b)] one sees the development of magnetic diffuse scattering, consistent with a spin-frustrated pyrochlore system above the structural transition. On further cooling to $T = 50$ K [Figs. 10(b) → 10(c)], the magnetic scattering changes into chain-like rods below the structural transition. Upon heating from $T = 50$ K back to $75$ K [Figs. 10(c) → 10(b)] one recovers these short-ranged magnetic correlations. However, if the sample is cooled below the long-range magnetic transition $T_N$ and then warmed up to $T = 75$ K [Figs. 10(a) → 10(d)] a different behavior is observed. Whilst at $T = 50$ K [Figs. 10(d) → 10(e)] one still observes chain-like scattering, above the structural transition [Figs. 10(e) → 10(f)], $T_S$ the
<figure><img src="image_12.png" /><figcaption>FIG. 10. The energy integrated magnetic cross section extracted from the scattering of polarized neutrons using the DNS diffractometer (FRM2). The arrows at the top and bottom of the figure indicate the cooling and heating history. Upon cooling to 3 K and reheating the diffuse scattering evolves hysteretically with no return to the short range ordered phase at T = 75 K. No such hysteresis is found if the sample is only cooled to 50 K and reheated (black arrow).</figcaption></figure>

short-ranged correlations are absent and no magnetic scatter-
ing is observable. This indicates that either the magnetic fluc-
tuations have moved outside our energy window determined
by kinematics on DNS with fixed Eᵢ = 4.64 meV, or the criti-
cal fluctuations are highly extended in momentum and energy.

We suggest that the explanation for this hysteresis tied
to cooling through the magnetic ordering transition at T_N
lies in the orbital order present in MgV₂O₄. At T = 75 K,
above the structural transition, in terms of a real orbital basis
the ground state has equal occupation of |xy⟩, |xz⟩ and |yz⟩
orbitals (Fig. 5). In this state the strength of the superexchange
is the same for both inter- and intrachain bonds, J₁ = J₂ and
MgV₂O₄ is a canonical frustrated pyrochlore. As one cools
below the structural transition, the local octahedra tetrago-
nally compress leading to a doubly degenerate state with l_z =
±1. The average occupancy of both the |xz⟩ and |yz⟩ orbitals is
½ and so the strength of the interchain coupling J₁ is reduced,
rendering the system quasi-one-dimensional, which gives rise
to rods in the magnetic diffuse scattering. At this point, long-
range magnetic order has yet to be established and so the
orbital ground state is doubly degenerate. Below the anti-
ferromagnetic ordering temperature, three-dimensional order
is established thanks to the cooperative effect of the third-
nearest-neighbor coupling J₃ and the anisotropy provided by
the distortion and spin-orbit coupling.

One key effect of longer-range coupling (beyond J₁ and J₂) is to reduce the degree to which the spin moment is suppressed due to fluctuations, both thermal and due to the logarithmic divergence of spin fluctuations in the one-dimensional quantum antiferromagnet (as per the Mermin-Wagner theorem [64]).

Upon establishment of long-range order, the orbital doublet
ground state is split by the molecular mean field, establish-
ing the COO ground state of Refs. [18,52]. As the sample
is then reheated, at T > T_N long-range magnetic order is
lost as fluctuations build, but we speculate that the orbital
moments remain frozen in their COO configuration as the
tetragonal symmetry prevents them from reorienting. Finally
as one reaches the structural transition from below, it should
be expected that the orbital moments reorient. However, given
the frozen nature of the orbital moments in the ground state
established by magnetic order, this is prevented by a potential
barrier established by the low-temperature molecular fields.
To overcome this requires a thermal energy equivalent to
the molecular field of order ~10 meV, or above ~100 K.
We speculate that through the Zeeman energy applied by the
molecular field induced through magnetic ordering, MgV₂O₄
displays an orbital memory of the low-temperature spin order.

VI. CONCLUSIONS

In conclusion, we have mapped out the spin-orbital fluctu-
ations in MgV₂O₄ using neutron spectroscopy where we have
observed two different branches. We have parameterized these
in terms of a spin-orbital excitonic theory where the single-ion
Hamiltonian is used to determine the quantized ground state,
which is then coupled on a lattice using RPA. The results are
strongly supportive of a COO ground state. We then used this
model to understand hysteretic magnetic fluctuations through
the three different spin-orbital phases.
<figure><img src="image_12.png" /><figcaption>FIG. 10. The energy integrated magnetic cross section extracted from the scattering of polarized neutrons using the DNS diffractometer (FRM2). The arrows at the top and bottom of the figure indicate the cooling and heating history. Upon cooling to 3 K and reheating the diffuse scattering evolves hysteretically with no return to the short range ordered phase at T = 75 K. No such hysteresis is found if the sample is only cooled to 50 K and reheated (black arrow).</figcaption></figure>

short-ranged correlations are absent and no magnetic scatter-
ing is observable. This indicates that either the magnetic fluc-
tuations have moved outside our energy window determined
by kinematics on DNS with fixed Eᵢ = 4.64 meV, or the criti-
cal fluctuations are highly extended in momentum and energy.

We suggest that the explanation for this hysteresis tied
to cooling through the magnetic ordering transition at T_N
lies in the orbital order present in MgV₂O₄. At T = 75 K,
above the structural transition, in terms of a real orbital basis
the ground state has equal occupation of |xy⟩, |xz⟩ and |yz⟩
orbitals (Fig. 5). In this state the strength of the superexchange
is the same for both inter- and intrachain bonds, J₁ = J₂ and
MgV₂O₄ is a canonical frustrated pyrochlore. As one cools
below the structural transition, the local octahedra tetrago-
nally compress leading to a doubly degenerate state with l_z =
±1. The average occupancy of both the |xz⟩ and |yz⟩ orbitals is
½ and so the strength of the interchain coupling J₁ is reduced,
rendering the system quasi-one-dimensional, which gives rise
to rods in the magnetic diffuse scattering. At this point, long-
range magnetic order has yet to be established and so the
orbital ground state is doubly degenerate. Below the anti-
ferromagnetic ordering temperature, three-dimensional order
is established thanks to the cooperative effect of the third-
nearest-neighbor coupling J₃ and the anisotropy provided by
the distortion and spin-orbit coupling.

One key effect of longer-range coupling (beyond J₁ and J₂) is to reduce the degree to which the spin moment is suppressed due to fluctuations, both thermal and due to the logarithmic divergence of spin fluctuations in the one-dimensional quantum antiferromagnet (as per the Mermin-Wagner theorem [64]).

Upon establishment of long-range order, the orbital doublet
ground state is split by the molecular mean field, establish-
ing the COO ground state of Refs. [18,52]. As the sample
is then reheated, at T > T_N long-range magnetic order is
lost as fluctuations build, but we speculate that the orbital
moments remain frozen in their COO configuration as the
tetragonal symmetry prevents them from reorienting. Finally
as one reaches the structural transition from below, it should
be expected that the orbital moments reorient. However, given
the frozen nature of the orbital moments in the ground state
established by magnetic order, this is prevented by a potential
barrier established by the low-temperature molecular fields.
To overcome this requires a thermal energy equivalent to
the molecular field of order ~10 meV, or above ~100 K.
We speculate that through the Zeeman energy applied by the
molecular field induced through magnetic ordering, MgV₂O₄
displays an orbital memory of the low-temperature spin order.

VI. CONCLUSIONS

In conclusion, we have mapped out the spin-orbital fluctu-
ations in MgV₂O₄ using neutron spectroscopy where we have
observed two different branches. We have parameterized these
in terms of a spin-orbital excitonic theory where the single-ion
Hamiltonian is used to determine the quantized ground state,
which is then coupled on a lattice using RPA. The results are
strongly supportive of a COO ground state. We then used this
model to understand hysteretic magnetic fluctuations through
the three different spin-orbital phases.
ACKNOWLEDGMENT

The authors thank E. Chan and M. Mourigal for useful dis-
cussions. H.L. was cofunded by the ISIS facility development
studentship programme and the EPSRC. C.S. acknowledges
support from the Carnegie Trust for the Universities of Scot-
land, EPSRC, and the STFC. P.M.S. acknowledges support
from the California NanoSystems Institute through the Elings
Fellowship program. S.D.W. acknowledges financial support
from the US Department of Energy (DOE), Office of Ba-
sic Energy Sciences, Division of Materials Sciences and
Engineering under Grant No. DE-SC0017752. The authors
gratefully acknowledge the financial support provided by the
Jülich Centre for Neutron Science (JCNS) to perform the
neutron scattering measurements at the Heinz Maier-Leibnitz
Zentrum (MLZ), Garching, Germany.

APPENDIX: SAMPLE CHARACTERIZATION

In this section we outline additional experimental information regarding crystal growth and characterization of the nuclear and magnetic properties of our MgV~2~O~4~ samples.

1. Additional experimental information

*Materials preparation.* Single crystals of MgV~2~O~4~ were grown using a four mirror image furnace (Crystal Systems Inc.). Dry magnesium oxide MgO (Sigma Aldrich, 99.995%) and vanadium (III) oxide V~2~O~3~ (Sigma Aldrich, 99.99%) were ground in an agate mortar and pestle with an excess of 1.5% V~2~O~3~. The excess V~2~O~3~ was added to combat the evaporation of vanadium oxide during the growth. For single-crystal growth, the ground powders were pressed in a rod and sintered in the image furnace on a lower power in a flowing argon atmosphere (10 bar pressure, 0.1 L/min). The crystal growth conditions were 5 mm/hr growth speed, 35 rpm rotation, 10 bar argon with flow rate ~0.1 L/min. The partial oxygen pressure was monitored on the input and output of the gas supply (Cambridge Sensotec) and the flow rate was adjusted to keep the partial pressure below ~10–15 ppm. A radiation-heated molybdenum getter was installed in the furnace to help reduce the oxygen partial pressure and keep it stable over the growth. Powders of MgV~2~O~4~ were synthesized using the same procedure (and samples) outlined in Ref. [65].

*Polarized neutron diffraction.* To investigate the momentum-broadened magnetic correlations characterizing the critical properties in the paramagnetic phases, polarized diffuse scattering measurements were performed on the cold neutron diffractometer DNS at the Heinz Maier-Leibnitz Zentrum (MLZ) with an incident neutron wavelength of 4.2 Å selected with a double-focusing pyrolytic graphite monochromator. While DNS has time-of-flight capabilities, data presented was collected in *non* time-of-flight mode, providing the integrated scattering intensity with energy transfers up to an *E*ᵢ = 4.64 meV. The flipping ratio was measured to be 20 and a weak 50 Gauss guide field was applied at the sample for the preservation of neutron polarization. A vanadium standard was used to account for variations of the detector response and the solid-angle coverage. Corrections for the polarization efficiency of the supermirror analyzers were performed by using the scattering

<figure><img src="image_13.png" /><figcaption>FIG. 11. The heat capacity and magnetic susceptibility measured on a piece of the single crystals used for our neutron studies and displayed in Fig. 1. The data clearly show the presence of two transitions, which we will show below are a cubic to tetragonal structural transition (≈ 60 K) and a lower antiferromagnetic transition (40 K) based on neutron and x-ray diffraction.</figcaption></figure>

from a NiCr alloy. Each measurement was performed over a
period of 8–12 hours after a sample thermalization period of
1 hour.

*Unpolarized neutron diffraction.* To characterize the struc-
tural properties of our single crystals in both the cubic and
tetragonal phases, unpolarized single-crystal neutron diffrac-
tion data were taken on the D9 diffractometer (ILL, France)
with a fixed incident wavelength of $\lambda = 0.842$ Å.

*Synchrotron x-ray diffraction.* Temperature-dependent x-ray diffraction on powders were performed to confirm structural properties and to search for any hysteretic effects in the lattice properties. Experiments were performed on the I11 high resolution beamline [66] at the Diamond Light Source (Didcot, UK) using the MAC detector system. A wavelength of $\lambda = 0.8258$ Å was used.

2. Thermodynamics

Heat capacity and magnetic susceptibility measurements (taken using Quantum Design MPMS and PPMS systems) on a small piece cut from the crystals [Fig. 2(b)] used for our neutron studies are displayed in Fig. 11. The data display two transitions at ≈ 60 K and 40 K. We discuss below based on neutron diffraction on pieces of these crystals that these two transitions represent a high-temperature structural transition from cubic to tetragonal unit cells (T~S~ ≈ 60 K) and a lower-temperature transition to an antiferromagnetically ordered phase (T~N~ ≈ 40 K).
ACKNOWLEDGMENT

The authors thank E. Chan and M. Mourigal for useful dis-
cussions. H.L. was cofunded by the ISIS facility development
studentship programme and the EPSRC. C.S. acknowledges
support from the Carnegie Trust for the Universities of Scot-
land, EPSRC, and the STFC. P.M.S. acknowledges support
from the California NanoSystems Institute through the Elings
Fellowship program. S.D.W. acknowledges financial support
from the US Department of Energy (DOE), Office of Ba-
sic Energy Sciences, Division of Materials Sciences and
Engineering under Grant No. DE-SC0017752. The authors
gratefully acknowledge the financial support provided by the
Jülich Centre for Neutron Science (JCNS) to perform the
neutron scattering measurements at the Heinz Maier-Leibnitz
Zentrum (MLZ), Garching, Germany.

APPENDIX: SAMPLE CHARACTERIZATION

In this section we outline additional experimental information regarding crystal growth and characterization of the nuclear and magnetic properties of our MgV~2~O~4~ samples.

1. Additional experimental information

*Materials preparation.* Single crystals of MgV~2~O~4~ were grown using a four mirror image furnace (Crystal Systems Inc.). Dry magnesium oxide MgO (Sigma Aldrich, 99.995%) and vanadium (III) oxide V~2~O~3~ (Sigma Aldrich, 99.99%) were ground in an agate mortar and pestle with an excess of 1.5% V~2~O~3~. The excess V~2~O~3~ was added to combat the evaporation of vanadium oxide during the growth. For single-crystal growth, the ground powders were pressed in a rod and sintered in the image furnace on a lower power in a flowing argon atmosphere (10 bar pressure, 0.1 L/min). The crystal growth conditions were 5 mm/hr growth speed, 35 rpm rotation, 10 bar argon with flow rate ~0.1 L/min. The partial oxygen pressure was monitored on the input and output of the gas supply (Cambridge Sensotec) and the flow rate was adjusted to keep the partial pressure below ~10–15 ppm. A radiation-heated molybdenum getter was installed in the furnace to help reduce the oxygen partial pressure and keep it stable over the growth. Powders of MgV~2~O~4~ were synthesized using the same procedure (and samples) outlined in Ref. [65].

*Polarized neutron diffraction.* To investigate the momentum-broadened magnetic correlations characterizing the critical properties in the paramagnetic phases, polarized diffuse scattering measurements were performed on the cold neutron diffractometer DNS at the Heinz Maier-Leibnitz Zentrum (MLZ) with an incident neutron wavelength of 4.2 Å selected with a double-focusing pyrolytic graphite monochromator. While DNS has time-of-flight capabilities, data presented was collected in *non* time-of-flight mode, providing the integrated scattering intensity with energy transfers up to an *E*ᵢ = 4.64 meV. The flipping ratio was measured to be 20 and a weak 50 Gauss guide field was applied at the sample for the preservation of neutron polarization. A vanadium standard was used to account for variations of the detector response and the solid-angle coverage. Corrections for the polarization efficiency of the supermirror analyzers were performed by using the scattering

<figure><img src="image_13.png" /></figure>

FIG. 11. The heat capacity and magnetic susceptibility measured on a piece of the single crystals used for our neutron studies and displayed in Fig. 1. The data clearly show the presence of two transitions, which we will show below are a cubic to tetragonal structural transition (≈ 60 K) and a lower antiferromagnetic transition (40 K) based on neutron and x-ray diffraction.

from a NiCr alloy. Each measurement was performed over a
period of 8–12 hours after a sample thermalization period of
1 hour.

*Unpolarized neutron diffraction.* To characterize the struc-
tural properties of our single crystals in both the cubic and
tetragonal phases, unpolarized single-crystal neutron diffrac-
tion data were taken on the D9 diffractometer (ILL, France)
with a fixed incident wavelength of $\lambda = 0.842$ Å.

*Synchrotron x-ray diffraction.* Temperature-dependent x-ray diffraction on powders were performed to confirm structural properties and to search for any hysteretic effects in the lattice properties. Experiments were performed on the I11 high resolution beamline [66] at the Diamond Light Source (Didcot, UK) using the MAC detector system. A wavelength of $\lambda = 0.8258$ Å was used.

**2. Thermodynamics**

Heat capacity and magnetic susceptibility measurements (taken using Quantum Design MPMS and PPMS systems) on a small piece cut from the crystals [Fig. 2(b)] used for our neutron studies are displayed in Fig. 11. The data display two transitions at ≈ 60 K and 40 K. We discuss below based on neutron diffraction on pieces of these crystals that these two transitions represent a high-temperature structural transition from cubic to tetragonal unit cells (T~S~ ≈ 60 K) and a lower-temperature transition to an antiferromagnetically ordered phase (T~N~ ≈ 40 K).
### 3. Temperature-dependent diffraction

We discuss the two phase transitions present in our single crystals using temperature-dependent powder and single-crystal diffraction. In terms of single-crystal neutron diffraction, Fig. 2(c) in the main paper illustrates the (0,0,2) intensity on cooling and warming displaying a structural change at ≈60 K with a significant difference in the intensity on heating and cooling. We speculate this originates from domains and strain in the sample resulting in changes in extinction on cooling/warming. We further discuss this point with regards to synchrotron x-ray measurements below. Figure 2(c) also plots the onset of magnetic order at $T_N$ ≈ 40 K probed with neutron diffraction by monitoring the (1,1,0) magnetic Bragg peak. The results are consistent with published data and illustrative of the sample quality [15, 23, 24, 67].

We further studied the structural transition using synchrotron x-ray diffraction and Fig. 2(c) displays the lattice parameters as a function of temperature. This further demonstrates that the high-temperature transition is structural. We have confirmed the unit cells and structural parameters in the literature in both the high- and low-temperature phases and representative diffraction patterns are shown in Fig. 12.

An interesting result in the context of our hysteretic studies outlined in the main text is shown in Fig. 13 where we plot representative powder data of the (4,0,0) Bragg peak (in the cubic setting). At high temperatures, a single peak is displayed representative of the cubic unit cell and a splitting is seen at low temperatures indicative of a structural distortion to a tetragonal unit cell. However, for a significant temperature region below the structural transition, a coexistence is observed. The relative fractions of the two structural phases, indicated by the relative peak intensities, is different dependent on the temperature history. This supports that the hysteresis in intensity displayed in Fig. 2 originates from relative domain populations.

### 4. Single-crystal nuclear structure

To characterize the structural properties in our single crystals we applied neutron diffraction on a small piece from our crystals used for spectroscopy. Owing to the large uncertainty related to possible mixing between the Mg and V sites in single crystalline $MgV_2O_4$, it was important to measure the nuclear structure above and below $T_S$ = 60 K. As demonstrated in Ref. [23], the site mixing has strong effects on the magnetic and structural transitions.

In agreement with previous studies [16], the refinements of the $T = 140$ K and 5 K temperature (Fig. 2) data sets confirm the high-temperature $Fd\bar{3}m$ and low-temperature $I4_1/amd$ crystal structures for $MgV_2O_4$, with $a = 8.42$ Å (cubic phase), and $a = 5.77$ Å and $c = 8.16$ Å in the tetragonal phase. The Mg and O sites were refined to be fully occupied within error, demonstrating the quality of the single crystal. The refined structural parameters (including domains, which are further discussed below) are listed in Table V.

There has been some discussion in the literature about the possibility of different structures based on the observation of a predicted extinct $\tilde{Q} = (0,0,2)$ Bragg peak [15]. We have observed this Bragg peak in our D9 unpolarized experiment. However, an azimuthal scan keeping the momentum transfer

<figure><img src="image_14.png" /><figcaption>FIG. 12. Representative refinements of $MgV_2O_4$ powders used to characterize the structural transition in the (a) cubic phase at $T = 100$ K and in the (b) tetragonal phase at $T = 6$ K.</figcaption></figure>

fixed found significant intensity variations indicating that this peak is likely due to multiple scattering processes. We therefore consider the low-temperature structure to be $I4_1/amd$ for the purposes of this paper.

The high-temperature cubic $Fd\bar{3}m$ space groups has a fourfold symmetry along each crystallographic axis. On cooling a macroscopic crystal through the Jahn-Teller transition to a $I4_1/amd$ tetragonal space group with a twofold and a fourfold symmetry axis. In the absence of an external symmetry breaking field such as strain this symmetry must, on average, be preserved, resulting in structural domains at low temperatures. There are therefore three possible low-temperature domains associated with a distortion along one of the axes of the high-temperature unit cell. These domains were refined not to be of equal population as indicated by the results in Table VI. This unequal domain population is also consistent with the findings of Ref. [16] and may be domain and environment dependent.
### 3. Temperature-dependent diffraction

We discuss the two phase transitions present in our single crystals using temperature-dependent powder and single-crystal diffraction. In terms of single-crystal neutron diffraction, Fig. 2(c) in the main paper illustrates the (0,0,2) intensity on cooling and warming displaying a structural change at ≈60 K with a significant difference in the intensity on heating and cooling. We speculate this originates from domains and strain in the sample resulting in changes in extinction on cooling/warming. We further discuss this point with regards to synchrotron x-ray measurements below. Figure 2(c) also plots the onset of magnetic order at $T_N$ ≈ 40 K probed with neutron diffraction by monitoring the (1,1,0) magnetic Bragg peak. The results are consistent with published data and illustrative of the sample quality [15, 23, 24, 67].

We further studied the structural transition using synchrotron x-ray diffraction and Fig. 2(c) displays the lattice parameters as a function of temperature. This further demonstrates that the high-temperature transition is structural. We have confirmed the unit cells and structural parameters in the literature in both the high- and low-temperature phases and representative diffraction patterns are shown in Fig. 12.

An interesting result in the context of our hysteretic studies outlined in the main text is shown in Fig. 13 where we plot representative powder data of the (4,0,0) Bragg peak (in the cubic setting). At high temperatures, a single peak is displayed representative of the cubic unit cell and a splitting is seen at low temperatures indicative of a structural distortion to a tetragonal unit cell. However, for a significant temperature region below the structural transition, a coexistence is observed. The relative fractions of the two structural phases, indicated by the relative peak intensities, is different dependent on the temperature history. This supports that the hysteresis in intensity displayed in Fig. 2 originates from relative domain populations.

### 4. Single-crystal nuclear structure

To characterize the structural properties in our single crystals we applied neutron diffraction on a small piece from our crystals used for spectroscopy. Owing to the large uncertainty related to possible mixing between the Mg and V sites in single crystalline $MgV_2O_4$, it was important to measure the nuclear structure above and below $T_S$ = 60 K. As demonstrated in Ref. [23], the site mixing has strong effects on the magnetic and structural transitions.

In agreement with previous studies [16], the refinements of the $T = 140$ K and 5 K temperature (Fig. 2) data sets confirm the high-temperature $Fd\bar{3}m$ and low-temperature $I4_1/amd$ crystal structures for $MgV_2O_4$, with $a = 8.42$ Å (cubic phase), and $a = 5.77$ Å and $c = 8.16$ Å in the tetragonal phase. The Mg and O sites were refined to be fully occupied within error, demonstrating the quality of the single crystal. The refined structural parameters (including domains, which are further discussed below) are listed in Table V.

There has been some discussion in the literature about the possibility of different structures based on the observation of a predicted extinct $\tilde{Q} = (0,0,2)$ Bragg peak [15]. We have observed this Bragg peak in our D9 unpolarized experiment. However, an azimuthal scan keeping the momentum transfer

<figure><img src="image_14.png" /><figcaption>FIG. 12. Representative refinements of $MgV_2O_4$ powders used to characterize the structural transition in the (a) cubic phase at $T = 100$ K and in the (b) tetragonal phase at $T = 6$ K.</figcaption></figure>

fixed found significant intensity variations indicating that this peak is likely due to multiple scattering processes. We therefore consider the low-temperature structure to be $I4_1/amd$ for the purposes of this paper.

The high-temperature cubic $Fd\bar{3}m$ space groups has a fourfold symmetry along each crystallographic axis. On cooling a macroscopic crystal through the Jahn-Teller transition to a $I4_1/amd$ tetragonal space group with a twofold and a fourfold symmetry axis. In the absence of an external symmetry breaking field such as strain this symmetry must, on average, be preserved, resulting in structural domains at low temperatures. There are therefore three possible low-temperature domains associated with a distortion along one of the axes of the high-temperature unit cell. These domains were refined not to be of equal population as indicated by the results in Table VI. This unequal domain population is also consistent with the findings of Ref. [16] and may be domain and environment dependent.
<figure><img src="image_15.png" /><figcaption>FIG. 13. Representative scans of the structural cubic (4,0,0) Bragg peak near the structural transition temperature. Different temperature cycles are displayed. The temperature cycle is (a) → (b) → (c) → (d).</figcaption></figure>

<table><caption>TABLE V. Atomic parameters for MgV$_2$O$_4$ at T = 140 K (cubic) and 50 K (tetragonal) from single-crystal neutron diffraction. Positions in the cubic phase are 8a (1/8, 1/8, 1/8), 16d (1/2, 1/2, 1/2), and 32e (x, x, x) for Mg, V, and O respectively, and in the 4a (0, 3/4, 1/8), 8d (0, 0, 1/2), and 16h (0, y, z) for the tetragonal phase. V was fixed to full occupancy. * Denotes parameters fixed from the T = 140 K refinement. Three different domains were observed and fractions included in the T = 50 K refinement yield 88(1)%, 8(1)%, and 4(1)% populations.</caption><thead><tr><th></th><th>T = 140 K</th><th>T = 50 K*</th></tr></thead><tbody><tr><td>Space Group</td><td>Fd<span style="text-decoration: overline">3</span>m</td><td>I4<sub>1</sub>/amd</td></tr><tr><td>x(O)</td><td>0.26000(6)</td><td>0.0</td></tr><tr><td>y(O)</td><td></td><td>0.47929(14)</td></tr><tr><td>z(O)</td><td></td><td>0.26002(8)</td></tr><tr><td>Mg B<sub>iso</sub> (Å<sup>2</sup>)</td><td>0.51(4)</td><td>0.16(3)</td></tr><tr><td>O B<sub>iso</sub> (Å<sup>2</sup>)</td><td>0.45(4)</td><td>0.24(2)</td></tr><tr><td>Mg occupancy (%)</td><td>97(7)</td><td>97*</td></tr><tr><td>O occupancy (%)</td><td>94(6)</td><td>94*</td></tr><tr><td>No. reflections</td><td>741</td><td>827</td></tr><tr><td>No. indep. Ref.</td><td>74</td><td>188</td></tr><tr><td>R<sub>f</sub> %</td><td>2.55</td><td>4.63</td></tr><tr><td>R<sub>w</sub> %</td><td>4.83</td><td>6.29</td></tr><tr><td>χ<sup>2</sup></td><td>2.15</td><td>3.14</td></tr></tbody></table>

## 5. Magnetic structure

To confirm the magnetic structure, we have used polarized neutron scattering measurements at DNS (MLZ) in the low-temperature ordered phase with the sample aligned in

<table><caption>TABLE VI. Atomic parameters for MgV$_2$O$_4$ at T = 100 K (cubic) and 7 K (tetragonal) from the powder sample (synchrotron). Positions in the cubic phase are 8a (1/8, 1/8, 1/8), 16d (1/2, 1/2, 1/2), and 32e (x, x, x) for Mg, V, and O respectively, and in the 4a (0, 3/4, 1/8), 8d (0, 0, 1/2), and 16h (0, y, z) for the tetragonal phase.</caption><thead><tr><th></th><th>T = 100 K</th><th>T = 7 K</th></tr></thead><tbody><tr><td>Space Group</td><td>Fd<span style="text-decoration: overline">3</span>m</td><td>I4<sub>1</sub>/amd</td></tr><tr><td>x(O)</td><td>0.25985(4)</td><td></td></tr><tr><td>y(O)</td><td></td><td>0.01833(2)</td></tr><tr><td>z(O)</td><td></td><td>0.26108(11)</td></tr><tr><td>Mg B<sub>iso</sub> (Å<sup>2</sup>)</td><td>0.386(5)</td><td>0.240(11)</td></tr><tr><td>O B<sub>iso</sub> (Å<sup>2</sup>)</td><td>0.578(6)</td><td>0.360(11)</td></tr><tr><td>V B<sub>iso</sub> (Å<sup>2</sup>)</td><td>0.366(1)</td><td>0.234(3)</td></tr><tr><td>R<sub>p</sub> %</td><td>5.17</td><td>8.96</td></tr><tr><td>R<sub>wp</sub> %</td><td>7.19</td><td>12.7</td></tr><tr><td>R<sub>B</sub> %</td><td>2.28</td><td>5.44</td></tr><tr><td>R<sub>f</sub> %</td><td>1.57</td><td>4.15</td></tr></tbody></table>
<figure><img src="image_15.png" /><figcaption>FIG. 13. Representative scans of the structural cubic (4,0,0) Bragg peak near the structural transition temperature. Different temperature cycles are displayed. The temperature cycle is (a) → (b) → (c) → (d).</figcaption></figure>

<table><caption>TABLE V. Atomic parameters for MgV$_2$O$_4$ at T = 140 K (cubic) and 50 K (tetragonal) from single-crystal neutron diffraction. Positions in the cubic phase are 8a (1/8, 1/8, 1/8), 16d (1/2, 1/2, 1/2), and 32e (x, x, x) for Mg, V, and O respectively, and in the 4a (0, 3/4, 1/8), 8d (0, 0, 1/2), and 16h (0, y, z) for the tetragonal phase. V was fixed to full occupancy. * Denotes parameters fixed from the T = 140 K refinement. Three different domains were observed and fractions included in the T = 50 K refinement yield 88(1)%, 8(1)%, and 4(1)% populations.</caption><thead><tr><th></th><th>T = 140 K</th><th>T = 50 K*</th></tr></thead><tbody><tr><td>Space Group</td><td>Fd<span style="text-decoration: overline">3</span>m</td><td>I4<sub>1</sub>/amd</td></tr><tr><td>x(O)</td><td>0.26000(6)</td><td>0.0</td></tr><tr><td>y(O)</td><td></td><td>0.47929(14)</td></tr><tr><td>z(O)</td><td></td><td>0.26002(8)</td></tr><tr><td>Mg B<sub>iso</sub> (Å<sup>2</sup>)</td><td>0.51(4)</td><td>0.16(3)</td></tr><tr><td>O B<sub>iso</sub> (Å<sup>2</sup>)</td><td>0.45(4)</td><td>0.24(2)</td></tr><tr><td>Mg occupancy (%)</td><td>97(7)</td><td>97*</td></tr><tr><td>O occupancy (%)</td><td>94(6)</td><td>94*</td></tr><tr><td>No. reflections</td><td>741</td><td>827</td></tr><tr><td>No. indep. Ref.</td><td>74</td><td>188</td></tr><tr><td>R<sub>f</sub> %</td><td>2.55</td><td>4.63</td></tr><tr><td>R<sub>w</sub> %</td><td>4.83</td><td>6.29</td></tr><tr><td>χ<sup>2</sup></td><td>2.15</td><td>3.14</td></tr></tbody></table>

### 5. Magnetic structure

To confirm the magnetic structure, we have used polarized neutron scattering measurements at DNS (MLZ) in the low-temperature ordered phase with the sample aligned in

<table><caption>TABLE VI. Atomic parameters for MgV$_2$O$_4$ at T = 100 K (cubic) and 7 K (tetragonal) from the powder sample (synchrotron). Positions in the cubic phase are 8a (1/8, 1/8, 1/8), 16d (1/2, 1/2, 1/2), and 32e (x, x, x) for Mg, V, and O respectively, and in the 4a (0, 3/4, 1/8), 8d (0, 0, 1/2), and 16h (0, y, z) for the tetragonal phase.</caption><thead><tr><th></th><th>T = 100 K</th><th>T = 7 K</th></tr></thead><tbody><tr><td>Space Group</td><td>Fd<span style="text-decoration: overline">3</span>m</td><td>I4<sub>1</sub>/amd</td></tr><tr><td>x(O)</td><td>0.25985(4)</td><td></td></tr><tr><td>y(O)</td><td></td><td>0.01833(2)</td></tr><tr><td>z(O)</td><td></td><td>0.26108(11)</td></tr><tr><td>Mg B<sub>iso</sub> (Å<sup>2</sup>)</td><td>0.386(5)</td><td>0.240(11)</td></tr><tr><td>O B<sub>iso</sub> (Å<sup>2</sup>)</td><td>0.578(6)</td><td>0.360(11)</td></tr><tr><td>V B<sub>iso</sub> (Å<sup>2</sup>)</td><td>0.366(1)</td><td>0.234(3)</td></tr><tr><td>R<sub>p</sub> %</td><td>5.17</td><td>8.96</td></tr><tr><td>R<sub>wp</sub> %</td><td>7.19</td><td>12.7</td></tr><tr><td>R<sub>B</sub> %</td><td>2.28</td><td>5.44</td></tr><tr><td>R<sub>f</sub> %</td><td>1.57</td><td>4.15</td></tr></tbody></table>
the (HHL) scattering plane. The measured magnetic Bragg peaks and intensity differences are consistent with the magnetic structure reported in ZnV₂O₄ [25] and is illustrated in Fig. 2(a). To guide the analysis of the magnetic structure, we have used the program SARAh [68] to calculate the irreducible representations and basis vectors for the space group. There are two common irreducible representations for the space group I4₁/amd with propagation vector (001). The polarized data at low temperatures shows a uniaxial structure within experimental error consistent with only one irreducible representation being involved in the magnetic transition. The involvement of basis functions from only one irreducible representation is consistent with the Tₙ ≈ 40 K magnetic transition being second order [Fig. 2(c)] as predicted from Landau theory. The magnetic moment was estimated by comparing magnetic $\vec{Q}$ = (110) and nuclear (111) Bragg peaks to be $\mu$ = 0.70 ± 0.15 $\mu_B$. This is in reasonable agreement with the *refined* magnetic moment value found by neutron diffraction in Ref. [15], $\mu_z$ = 0.47 $\mu_B$ (with an ~ 8° canting from the z axis). The magnetic structure in Figs. 2(a) and 6 is based upon interpenetrating chains of magnetic V³⁺, which carry a magnetic moment. The chain-like structure originates from the tetragonal structural distortion, which occurs at Tₛ that breaks the frustrating magnetic interactions present in the high-temperature cubic unit cell as discussed in the main text.

## 6. Polarization analysis

Polarized neutron spin-flip (SF) and non-spin-flip (NSF) scattering using an xyz geometry was performed on the DNS diffractometer. Magnetic neutron scattering obeys the selection rule that it is sensitive only to components of the magnetic vector perpendicular to the momentum transfer. In Ref. [69], this projected vector is denoted as $\vec{S}_⊥$ and polarized neutron scattering affords a means of separating out the various components of $\vec{S}_⊥$. The components of $\vec{S}_⊥$, which are perpendicular to the neutron polarization appear in the SF channel, while the component parallel to the neutron polarization will produce NSF scattering.

A problem with performing polarized neutron experiments on vanadium based materials is the fact that vanadium has

a strong incoherent cross section, which translates into an overall momentum-independent background. In polarized experiments, the incoherent cross section in the SF channel is twice that in the NSF channel, which complicates comparing the channels directly. To circumvent this, we considered the xyz channels in the SF and NSF channels separately. In our experiment, the x channel was defined as when the neutron polarization was parallel to the average momentum transfer $\vec{Q}$; y was set to be perpendicular to x, but within the horizontal scattering plane; and z was defined as vertical being perpendicular to x, y, and z is normal to the horizontal scattering plane.

Given the selection rules by neutron scattering $S_{⊥,x}$ is always absent, the problem is then to extract $S_{⊥,y}$ and $S_{⊥,z}$ reliably. Using these combined selection rules we note the following for the magnetic neutron intensities (I) in the SF channel including incoherent background ($B_{incoh}$):

$$ 
\begin{aligned}
I_{SF,x} &= S_{⊥,y} + S_{⊥,z} + B_{\text{incoh,SF}}, \\
I_{SF,y} &= S_{⊥,z} + B_{\text{incoh,SF}}, \\
I_{SF,z} &= S_{⊥,y} + B_{\text{incoh,SF}},
\end{aligned}
\quad (A1) $$

and in the NSF channel,

$$ 
\begin{aligned}
I_{\text{NSF},x} &= B_{\text{incoh,NSF}}, \\
I_{\text{NSF},y} &= S_{⊥,y} + B_{\text{incoh,NSF}}, \\
I_{\text{NSF},z} &= S_{⊥,z} + B_{\text{incoh,NSF}}.
\end{aligned}
\quad (A2) $$

Since the magnetic diffuse scattering does not overlap with the nuclear scattering we have not considered feedthrough from one channel to another due to the incomplete flipping ratio or misalignment of the spin with respect to the momentum transfer. We have removed the scattering around the nuclear Bragg peaks in this analysis. In an ideal experimental setup, $2 \times B_{incoh,NSF} = B_{incoh,SF}$; however, we have treated SF and NSF scattering separately to extract the individual cross section for $S_{⊥,y}$ and $S_{⊥,z}$. For example, from the SF channel, $S_{⊥,y} = I_{SF,x} - I_{SF,y}$ and $S_{⊥,z} = I_{SF,x} - I_{SF,z}$. From the NSF channel, $S_{⊥,y} = I_{NSF,y} - I_{NSF,x}$ and $S_{⊥,z} = I_{NSF,z} - I_{NSF,x}$. Therefore, the required components for $S_{⊥,y}$ and $S_{⊥,z}$ can be extracted independently from both the SF and NSF channels and combined. This method provides us a means of subtracting off the incoherent background from vanadium.

[1] S. H. Lee, H. Takagi, D. Louca, M. Matsuda, S. Ji, H. Ueda, T. Katsufuji, J. H. Chung, S. Park, S. W. Cheong, and C. Broholm, Frustrated magnetism and cooperative phase transitions in spinels, J. Phys. Soc. Jpn. 79, 011004 (2010).

[2] D. I. Khomskii, Review-orbital physics: Glorious past, bright future, ECS J. Solid State Sci. Technol. 11, 054004 (2022).

[3] D. I. Khomskii and S. V. Streltsov, Orbital effects in solids: Basics, recent progress, and opportunities, Chem. Rev. 121, 2992 (2021).

[4] S. V. Streltsov and D. I. Khomskii, Jahn-Teller effect and spin-orbit coupling: Friends or foes? Phys. Rev. X 10, 031043 (2020).

[5] S. V. Streltsov and D. I. Khomskii, Orbital physics in transition metal compounds: New trends, Phys. Usp. 60, 1121 (2017).

[6] S. V. Streltsov and D. I. Khomskii, Covalent bonds against magnetism in transition metal compounds, Proc. Natl. Acad. Sci. USA 113, 10491 (2016).

[7] T. Tokura and N. Nagaosa, Orbital physics in transition-metal oxides, Science 288, 462 (2000).

[8] L. D. Landau and E. M. Lifshitz, Quantum Mechanics (Non-Relativistic Theory) (Elsevier Science, Oxford, 1977).

[9] J. S. Gardner, M. J. P. Gingras, and J. E. Greedan, Magnetic pyrochlore oxides, Rev. Mod. Phys. 82, 53 (2010).

[10] R. Moessner, Magnets with strong geometric frustration, Can. J. Phys. 79, 1283 (2001).
the (HHL) scattering plane. The measured magnetic Bragg peaks and intensity differences are consistent with the magnetic structure reported in ZnV₂O₄ [25] and is illustrated in Fig. 2(a). To guide the analysis of the magnetic structure, we have used the program SARAh [68] to calculate the irreducible representations and basis vectors for the space group. There are two common irreducible representations for the space group I4₁/amd with propagation vector (001). The polarized data at low temperatures shows a uniaxial structure within experimental error consistent with only one irreducible representation being involved in the magnetic transition. The involvement of basis functions from only one irreducible representation is consistent with the Tₙ ≈ 40 K magnetic transition being second order [Fig. 2(c)] as predicted from Landau theory. The magnetic moment was estimated by comparing magnetic $\vec{Q} = (110)$ and nuclear (111) Bragg peaks to be $\mu = 0.70 \pm 0.15 \mu_B$. This is in reasonable agreement with the refined magnetic moment value found by neutron diffraction in Ref. [15], $\mu_z = 0.47 \mu_B$ (with an ~ 8° canting from the z axis). The magnetic structure in Figs. 2(a) and 6 is based upon interpenetrating chains of magnetic V³⁺, which carry a magnetic moment. The chain-like structure originates from the tetragonal structural distortion, which occurs at Tₛ that breaks the frustrating magnetic interactions present in the high-temperature cubic unit cell as discussed in the main text.

## 6. Polarization analysis

Polarized neutron spin-flip (SF) and non-spin-flip (NSF) scattering using an xyz geometry was performed on the DNS diffractometer. Magnetic neutron scattering obeys the selection rule that it is sensitive only to components of the magnetic vector perpendicular to the momentum transfer. In Ref. [69], this projected vector is denoted as $\vec{S}_⊥$ and polarized neutron scattering affords a means of separating out the various components of $\vec{S}_⊥$. The components of $\vec{S}_⊥$, which are perpendicular to the neutron polarization appear in the SF channel, while the component parallel to the neutron polarization will produce NSF scattering.

A problem with performing polarized neutron experiments on vanadium based materials is the fact that vanadium has

a strong incoherent cross section, which translates into an overall momentum-independent background. In polarized experiments, the incoherent cross section in the SF channel is twice that in the NSF channel, which complicates comparing the channels directly. To circumvent this, we considered the xyz channels in the SF and NSF channels separately. In our experiment, the x channel was defined as when the neutron polarization was parallel to the average momentum transfer $\vec{Q}$; y was set to be perpendicular to x, but within the horizontal scattering plane; and z was defined as vertical being perpendicular to x, y, and z is normal to the horizontal scattering plane.

Given the selection rules by neutron scattering $S_{⊥,x}$ is always absent, the problem is then to extract $S_{⊥,y}$ and $S_{⊥,z}$ reliably. Using these combined selection rules we note the following for the magnetic neutron intensities (I) in the SF channel including incoherent background ($B_{incoh}$):

$$ 
\begin{aligned}
I_{SF,x} &= S_{⊥,y} + S_{⊥,z} + B_{\text{incoh,SF}}, \\
I_{SF,y} &= S_{⊥,z} + B_{\text{incoh,SF}}, \\
I_{SF,z} &= S_{⊥,y} + B_{\text{incoh,SF}},
\end{aligned} 
\quad (A1) $$

and in the NSF channel,

$$ 
\begin{aligned}
I_{\text{NSF},x} &= B_{\text{incoh,NSF}}, \\
I_{\text{NSF},y} &= S_{⊥,y} + B_{\text{incoh,NSF}}, \\
I_{\text{NSF},z} &= S_{⊥,z} + B_{\text{incoh,NSF}}.
\end{aligned} 
\quad (A2) $$

Since the magnetic diffuse scattering does not overlap with the nuclear scattering we have not considered feedthrough from one channel to another due to the incomplete flipping ratio or misalignment of the spin with respect to the momentum transfer. We have removed the scattering around the nuclear Bragg peaks in this analysis. In an ideal experimental setup, $2 \times B_{\text{incoh,NSF}} = B_{\text{incoh,SF}}$; however, we have treated SF and NSF scattering separately to extract the individual cross section for $S_{⊥,y}$ and $S_{⊥,z}$. For example, from the SF channel, $S_{⊥,y} = I_{SF,x} - I_{SF,y}$ and $S_{⊥,z} = I_{SF,x} - I_{SF,z}$. From the NSF channel, $S_{⊥,y} = I_{NSF,y} - I_{NSF,x}$ and $S_{⊥,z} = I_{NSF,z} - I_{NSF,x}$. Therefore, the required components for $S_{⊥,y}$ and $S_{⊥,z}$ can be extracted independently from both the SF and NSF channels and combined. This method provides us a means of subtracting off the incoherent background from vanadium.

[1] S. H. Lee, H. Takagi, D. Louca, M. Matsuda, S. Ji, H. Ueda, T. Katsufuji, J. H. Chung, S. Park, S. W. Cheong, and C. Broholm, Frustrated magnetism and cooperative phase transitions in spinels, J. Phys. Soc. Jpn. 79, 011004 (2010).

[2] D. I. Khomskii, Review-orbital physics: Glorious past, bright future, ECS J. Solid State Sci. Technol. 11, 054004 (2022).

[3] D. I. Khomskii and S. V. Streltsov, Orbital effects in solids: Basics, recent progress, and opportunities, Chem. Rev. 121, 2992 (2021).

[4] S. V. Streltsov and D. I. Khomskii, Jahn-Teller effect and spin-orbit coupling: Friends or foes? Phys. Rev. X 10, 031043 (2020).

[5] S. V. Streltsov and D. I. Khomskii, Orbital physics in transition metal compounds: New trends, Phys. Usp. 60, 1121 (2017).

[6] S. V. Streltsov and D. I. Khomskii, Covalent bonds against magnetism in transition metal compounds, Proc. Natl. Acad. Sci. USA 113, 10491 (2016).

[7] T. Tokura and N. Nagaosa, Orbital physics in transition-metal oxides, Science 288, 462 (2000).

[8] L. D. Landau and E. M. Lifshitz, Quantum Mechanics (Non-Relativistic Theory) (Elsevier Science, Oxford, 1977).

[9] J. S. Gardner, M. J. P. Gingras, and J. E. Greedan, Magnetic pyrochlore oxides, Rev. Mod. Phys. 82, 53 (2010).

[10] R. Moessner, Magnets with strong geometric frustration, Can. J. Phys. 79, 1283 (2001).
[11] P. G. Radaelli, Orbital ordering in transition-metal spinels, *New J. Phys.* **7**, 53 (2005).

[12] V. Tsurkan, H.-A. K. von Nidda, J. Deisenhofer, P. Lunkenheimer, and A. Loidl, On the complexity of spinels: Magnetic, electronic, and polar ground states, *Phys. Rep.* **926**, 1 (2021).

[13] V. O. Garlea, R. Jin, D. Mandrus, B. Roessli, Q. Huang, M. Miller, A. J. Schultz, and S. E. Nagler, Magnetic and orbital ordering in the spinel MnV₂O₄, *Phys. Rev. Lett.* **100**, 066404 (2008).

[14] J. Krishna, N. Singh, S. Shallcross, J. K. Dewhurst, E. K. U. Gross, T. Maitra, and S. Sharma, Complete description of the magnetic ground state in spinel vanadates, *Phys. Rev. B* **100**, 081102(R) (2019).

[15] E. M. Wheeler, B. Lake, A. T. M. N. Islam, M. Reehuis, P. Steffens, T. Guidi, and A. H. Hill, Spin and orbital order in the vanadium spinel MgV₂O₄, *Phys. Rev. B* **82**, 140406(R) (2010).

[16] S. Niitaka, H. Ohsumi, K. Sugimoto, S. Lee, Y. Oshima, K. Kato, D. Hashizume, T. Arima, M. Takata, and H. Takagi, A-type antiferro-orbital ordering with I4₁/a symmetry and geometrical frustration in the spinel vanadate MgV₂O₄, *Phys. Rev. Lett.* **111**, 267201 (2013).

[17] O. Tchernyshyov, Structural, orbital, and magnetic order in vanadium spinels, *Phys. Rev. Lett.* **93**, 157206 (2004).

[18] N. B. Perkins and O. Sikora, Magnetic excitations in vanadium spinels, *Phys. Rev. B* **76**, 214434 (2007).

[19] P. M. Sarte, M. Songvilay, E. Pachoud, R. A. Ewings, C. D. Frost, D. Prabhakaran, K. H. Hong, A. J. Browne, Z. Yamani, J. P. Attfield, E. E. Rodriguez, S. D. Wilson, and C. Stock, Spin-orbit excitons in CoO, *Phys. Rev. B* **100**, 075143 (2019).

[20] P. M. Sarte, C. Stock, B. R. Ortiz, K. H. Hong, and S. D. Wilson, Van Vleck excitons in Ca₂RuO₄, *Phys. Rev. B* **102**, 245119 (2020).

[21] W. J. L. Buyers, T. M. Holden, and A. Perreault, Temperature dependence of magnetic excitations in singlet-ground-state systems. II. Excited-state spin waves near the Curie temperature in Pr₃Ti, *Phys. Rev. B* **11**, 266 (1975).

[22] P. M. Sarte, S. D. Wilson, J. P. Attfield, and C. Stock, Magnetic fluctuations and the spin-orbit interaction in Mott insulating CoO, *J. Phys.: Condens. Matter* **32**, 374011 (2020).

[23] A. T. M. N. Islam, E. M. Wheeler, M. Reehuis, K. Siemensmeyer, M. Tovar, B. Klemke, K. Kiefer, A. H. Hill, and B. Lake, Growth and magnetic properties of stoichiometric and site-disordered single crystalline MgV₂O₄, *Phys. Rev. B* **85**, 024203 (2012).

[24] H. Mamiya, M. Onoda, T. Furubayashi, J. Tang, and I. Nakatani, Structural and magnetic studies on vanadium spinel MgV₂O₄, *J. Appl. Phys.* **81**, 5289 (1997).

[25] M. Reehuis, A. Krimmel, N. Buttgen, A. Loidl, and A. Prokofiev, Crystallographic and magnetic structure of ZnV₂O₄, *Eur. Phys. J. B* **35**, 311 (2003).

[26] H. Maier-Leibnitz Zentrum *et al.*, DNS: Diffuse scattering neutron time-of-flight spectrometer, *J. Large-Scale Res. Facilities* **1**, A27 (2015).

[27] R. I. Bewley, R. S. Eccleston, K. A. McEwen, S. M. Hayden, M. T. Dove, S. M. Bennington, J. R. Treadgold, and R. L. S. Coleman, MERLIN, a new high count rate spectrometer at ISIS, *Phys. B: Condens. Matter* **385-386**, 1029 (2006).

[28] U. Stuhr, B. Roessli, S. Gvasaliya, H. M. Rønnow, U. Filges, D. Graf, A. Bollhalder, D. Hohl, R. Bürge, M. Schild *et al.*,

The thermal triple-axis-spectrometer EIGER at the continuous spallation source SINQ, Nucl. Instrum. Methods Phys. Res., Sect. A **853**, 16 (2017).

[29] K. Lefmann, D. F. McMorrow, H. M. Rønnow, K. Nielsen, K. N. Clausen, B. Lake, and G. Aeppli, Added flexibility in triple axis spectrometers: The two RITAs at Risø, *Phys. B: Condens. Matter* **283**, 343 (2000).

[30] S.-H. Lee, D. Louca, H. Ueda, S. Park, T. J. Sato, M. Isobe, Y. Ueda, S. Rosenkranz, P. Zschack, J. Íñiguez, Y. Qiu, and R. Osborn, Orbital and spin chains in ZnV₂O₄, *Phys. Rev. Lett.* **93**, 156407 (2004).

[31] D. S. McClure, *Electronic Spectra of Molecules and Ions in Crystals, Part II. Spectra of Ions in Crystals* (Academic Press, New York, 1959), pp. 399–525.

[32] A. Abragam and B. Bleaney, *Electron Paramagnetic Resonance of Transition Ions* (Dover Publications, New York, 1986).

[33] J. G. Rau and Michel J. P. Gingras, Frustrated quantum rare-earth pyrochlores, *Annu. Rev. Condens. Matter Phys.* **10**, 357 (2019).

[34] K. A. Ross, L. Savary, B. D. Gaulin, and L. Balents, Quantum excitations in quantum spin ice, *Phys. Rev. X* **1**, 021002 (2011).

[35] L. Savary, K. A. Ross, B. D. Gaulin, J. P. C. Ruff, and L. Balents, Order by quantum disorder in Er₂Ti₂O₇, *Phys. Rev. Lett.* **109**, 167201 (2012).

[36] J. D. Thompson, P. A. McClarty, D. Prabhakaran, I. Cabrera, T. Guidi, and R. Coldea, Quasiparticle breakdown and spin hamiltonian of the frustrated quantum pyrochlore Yb₂Ti₂O₇ in a magnetic field, *Phys. Rev. Lett.* **119**, 057203 (2017).

[37] L. Savary and L. Balents, Quantum spin liquids: A review, *Rep. Prog. Phys.* **80**, 016502 (2017).

[38] P. M. Sarte, R. A. Cowley, E. E. Rodriguez, E. Pachoud, D. Le, V. García-Sakai, J. W. Taylor, C. D. Frost, D. Prabhakaran, C. MacEwen, A. Kitada, A. J. Browne, M. Songvilay, Z. Yamani, W. J. L. Buyers, J. P. Attfield, and C. Stock, Disentangling orbital and spin exchange interactions for Co²⁺ on a rocksalt lattice, *Phys. Rev. B* **98**, 024415 (2018).

[39] P. M. Sarte, A. M. Arévalo-López, M. Songvilay, D. Le, T. Guidi, V. García-Sakai, S. Mukhopadhyay, S. C. Capelli, W. D. Ratcliff, K. H. Hong, G. M. McNally, E. Pachoud, J. P. Attfield, and C. Stock, Ordered magnetism in the intrinsically decorated $j_{eff} = \frac{1}{2}\alpha$-CoV₃O₈, *Phys. Rev. B* **98**, 224410 (2018).

[40] M. Tinkham, *Group Theory and Quantum Mechanics*, Dover Books on Chemistry Series (Dover Publications, Mineola, NY, 2003).

[41] H. Lane, E. Pachoud, J. A. Rodriguez-Rivera, M. Songvilay, G. Xu, P. M. Gehring, J. P. Attfield, R. A. Ewings, and C. Stock, Two-dimensional ferromagnetic spin-orbital excitations in honeycomb VI₃, *Phys. Rev. B* **104**, L020411 (2021).

[42] K. Yosida, *Theory of Magnetism* (Springer, New York, 1996).

[43] D. I. Khomskii, *Transition Metal Compounds* (Cambridge University Press, Cambridge, 2014).

[44] R. A. Cowley, W. J. L. Buyers, P. Martel, and R. W. H. Stevenson, Magnetic excitations and magnetic critical scattering in cobalt fluoride, *J. Phys. C* **6**, 2997 (1973).

[45] S. T. Bramwell and M. J. Harris, Frustration in Ising-type spin models on the pyrochlore lattice, *J. Phys.: Condens. Matter* **10**, L215 (1998).

[46] R. Moessner and J. T. Chalker, Properties of a classical spin liquid: The Heisenberg pyrochlore antiferromagnet, *Phys. Rev. Lett.* **80**, 2929 (1998).
[11] P. G. Radaelli, Orbital ordering in transition-metal spinels, *New J. Phys.* **7**, 53 (2005).

[12] V. Tsurkan, H.-A. K. von Nidda, J. Deisenhofer, P. Lunkenheimer, and A. Loidl, On the complexity of spinels: Magnetic, electronic, and polar ground states, *Phys. Rep.* **926**, 1 (2021).

[13] V. O. Garlea, R. Jin, D. Mandrus, B. Roessli, Q. Huang, M. Miller, A. J. Schultz, and S. E. Nagler, Magnetic and orbital ordering in the spinel MnV₂O₄, *Phys. Rev. Lett.* **100**, 066404 (2008).

[14] J. Krishna, N. Singh, S. Shallcross, J. K. Dewhurst, E. K. U. Gross, T. Maitra, and S. Sharma, Complete description of the magnetic ground state in spinel vanadates, *Phys. Rev. B* **100**, 081102(R) (2019).

[15] E. M. Wheeler, B. Lake, A. T. M. N. Islam, M. Reehuis, P. Steffens, T. Guidi, and A. H. Hill, Spin and orbital order in the vanadium spinel MgV₂O₄, *Phys. Rev. B* **82**, 140406(R) (2010).

[16] S. Niitaka, H. Ohsumi, K. Sugimoto, S. Lee, Y. Oshima, K. Kato, D. Hashizume, T. Arima, M. Takata, and H. Takagi, A-type antiferro-orbital ordering with *I*4₁/a symmetry and geometrical frustration in the spinel vanadate MgV₂O₄, *Phys. Rev. Lett.* **111**, 267201 (2013).

[17] O. Tchernyshyov, Structural, orbital, and magnetic order in vanadium spinels, *Phys. Rev. Lett.* **93**, 157206 (2004).

[18] N. B. Perkins and O. Sikora, Magnetic excitations in vanadium spinels, *Phys. Rev. B* **76**, 214434 (2007).

[19] P. M. Sarte, M. Songvilay, E. Pachoud, R. A. Ewings, C. D. Frost, D. Prabhakaran, K. H. Hong, A. J. Browne, Z. Yamani, J. P. Attfield, E. E. Rodriguez, S. D. Wilson, and C. Stock, Spin-orbit excitons in CoO, *Phys. Rev. B* **100**, 075143 (2019).

[20] P. M. Sarte, C. Stock, B. R. Ortiz, K. H. Hong, and S. D. Wilson, Van Vleck excitons in Ca₂RuO₄, *Phys. Rev. B* **102**, 245119 (2020).

[21] W. J. L. Buyers, T. M. Holden, and A. Perreault, Temperature dependence of magnetic excitations in singlet-ground-state systems. II. Excited-state spin waves near the Curie temperature in Pr₃Ti, *Phys. Rev. B* **11**, 266 (1975).

[22] P. M. Sarte, S. D. Wilson, J. P. Attfield, and C. Stock, Magnetic fluctuations and the spin-orbit interaction in Mott insulating CoO, *J. Phys.: Condens. Matter* **32**, 374011 (2020).

[23] A. T. M. N. Islam, E. M. Wheeler, M. Reehuis, K. Siemensmeyer, M. Tovar, B. Klemke, K. Kiefer, A. H. Hill, and B. Lake, Growth and magnetic properties of stoichiometric and site-disordered single crystalline MgV₂O₄, *Phys. Rev. B* **85**, 024203 (2012).

[24] H. Mamiya, M. Onoda, T. Furubayashi, J. Tang, and I. Nakatani, Structural and magnetic studies on vanadium spinel MgV₂O₄, *J. Appl. Phys.* **81**, 5289 (1997).

[25] M. Reehuis, A. Krimmel, N. Buttgen, A. Loidl, and A. Prokofiev, Crystallographic and magnetic structure of ZnV₂O₄, *Eur. Phys. J. B* **35**, 311 (2003).

[26] H. Maier-Leibnitz Zentrum *et al.*, DNS: Diffuse scattering neutron time-of-flight spectrometer, *J. Large-Scale Res. Facilities* **1**, A27 (2015).

[27] R. I. Bewley, R. S. Eccleston, K. A. McEwen, S. M. Hayden, M. T. Dove, S. M. Bennington, J. R. Treadgold, and R. L. S. Coleman, MERLIN, a new high count rate spectrometer at ISIS, *Phys. B: Condens. Matter* **385-386**, 1029 (2006).

[28] U. Stuhr, B. Roessli, S. Gvasaliya, H. M. Rønnow, U. Filges, D. Graf, A. Bollhalder, D. Hohl, R. Bürge, M. Schild *et al.*,

The thermal triple-axis-spectrometer EIGER at the continuous spallation source SINQ, Nucl. Instrum. Methods Phys. Res., Sect. A **853**, 16 (2017).

[29] K. Lefmann, D. F. McMorrow, H. M. Rønnow, K. Nielsen, K. N. Clausen, B. Lake, and G. Aeppli, Added flexibility in triple axis spectrometers: The two RITAs at Risø, *Phys. B: Condens. Matter* **283**, 343 (2000).

[30] S.-H. Lee, D. Louca, H. Ueda, S. Park, T. J. Sato, M. Isobe, Y. Ueda, S. Rosenkranz, P. Zschack, J. Íñiguez, Y. Qiu, and R. Osborn, Orbital and spin chains in ZnV₂O₄, *Phys. Rev. Lett.* **93**, 156407 (2004).

[31] D. S. McClure, *Electronic Spectra of Molecules and Ions in Crystals, Part II. Spectra of Ions in Crystals* (Academic Press, New York, 1959), pp. 399–525.

[32] A. Abragam and B. Bleaney, *Electron Paramagnetic Resonance of Transition Ions* (Dover Publications, New York, 1986).

[33] J. G. Rau and Michel J. P. Gingras, Frustrated quantum rare-earth pyrochlores, *Annu. Rev. Condens. Matter Phys.* **10**, 357 (2019).

[34] K. A. Ross, L. Savary, B. D. Gaulin, and L. Balents, Quantum excitations in quantum spin ice, *Phys. Rev. X* **1**, 021002 (2011).

[35] L. Savary, K. A. Ross, B. D. Gaulin, J. P. C. Ruff, and L. Balents, Order by quantum disorder in Er₂Ti₂O₇, *Phys. Rev. Lett.* **109**, 167201 (2012).

[36] J. D. Thompson, P. A. McClarty, D. Prabhakaran, I. Cabrera, T. Guidi, and R. Coldea, Quasiparticle breakdown and spin hamiltonian of the frustrated quantum pyrochlore Yb₂Ti₂O₇ in a magnetic field, *Phys. Rev. Lett.* **119**, 057203 (2017).

[37] L. Savary and L. Balents, Quantum spin liquids: A review, *Rep. Prog. Phys.* **80**, 016502 (2017).

[38] P. M. Sarte, R. A. Cowley, E. E. Rodriguez, E. Pachoud, D. Le, V. García-Sakai, J. W. Taylor, C. D. Frost, D. Prabhakaran, C. MacEwen, A. Kitada, A. J. Browne, M. Songvilay, Z. Yamani, W. J. L. Buyers, J. P. Attfield, and C. Stock, Disentangling orbital and spin exchange interactions for Co²⁺ on a rocksalt lattice, *Phys. Rev. B* **98**, 024415 (2018).

[39] P. M. Sarte, A. M. Arévalo-López, M. Songvilay, D. Le, T. Guidi, V. García-Sakai, S. Mukhopadhyay, S. C. Capelli, W. D. Ratcliff, K. H. Hong, G. M. McNally, E. Pachoud, J. P. Attfield, and C. Stock, Ordered magnetism in the intrinsically decorated $j_{eff} = \frac{1}{2}\alpha$-CoV₃O₈, *Phys. Rev. B* **98**, 224410 (2018).

[40] M. Tinkham, *Group Theory and Quantum Mechanics*, Dover Books on Chemistry Series (Dover Publications, Mineola, NY, 2003).

[41] H. Lane, E. Pachoud, J. A. Rodriguez-Rivera, M. Songvilay, G. Xu, P. M. Gehring, J. P. Attfield, R. A. Ewings, and C. Stock, Two-dimensional ferromagnetic spin-orbital excitations in honeycomb VI₃, *Phys. Rev. B* **104**, L020411 (2021).

[42] K. Yosida, *Theory of Magnetism* (Springer, New York, 1996).

[43] D. I. Khomskii, *Transition Metal Compounds* (Cambridge University Press, Cambridge, 2014).

[44] R. A. Cowley, W. J. L. Buyers, P. Martel, and R. W. H. Stevenson, Magnetic excitations and magnetic critical scattering in cobalt fluoride, *J. Phys. C* **6**, 2997 (1973).

[45] S. T. Bramwell and M. J. Harris, Frustration in Ising-type spin models on the pyrochlore lattice, *J. Phys.: Condens. Matter* **10**, L215 (1998).

[46] R. Moessner and J. T. Chalker, Properties of a classical spin liquid: The Heisenberg pyrochlore antiferromagnet, *Phys. Rev. Lett.* **80**, 2929 (1998).
[47] Y. Iqbal, T. Müller, P. Ghosh, M. J. P. Gingras, H. O. Jeschke, S. Rachel, J. Reuther, and R. Thomale, Quantum and classical phases of the pyrochlore Heisenberg model with competing interactions, Phys. Rev. X **9**, 011005 (2019).

[48] S. Zhang, H. J. Changlani, K. W. Plumb, O. Tchernyshyov, and R. Moessner, Dynamical structure factor of the three-dimensional quantum spin liquid candidate $NaCaNi_2F_7$, Phys. Rev. Lett. **122**, 167203 (2019).

[49] J. S. Gardner, B. D. Gaulin, S.-H. Lee, C. Broholm, N. P. Raju, and J. E. Greedan, Glassy statics and dynamics in the chemically ordered pyrochlore antiferromagnet $Y_2Mo_2O_7$, Phys. Rev. Lett. **83**, 211 (1999).

[50] J. W. Krizan and R. J. Cava, $NaCaNi_2F_7$: A frustrated high-temperature pyrochlore antiferromagnet with $S = 1$ $Ni^{2+}$, Phys. Rev. B **92**, 014406 (2015).

[51] K. W. Plumb, H. J. Changlani, A. Scheie, S. Zhang, J. W. Krizan, J. A. Rodriguez-Rivera, Y. Qiu, B. Winn, R. J. Cava, and C. L. Broholm, Continuum of quantum fluctuations in a three-dimensional $S = 1$ Heisenberg magnet, Nat. Phys. **15**, 54 (2019).

[52] S. Di Matteo, G. Jackeli, and N. B. Perkins, Orbital order in vanadium spinels, Phys. Rev. B **72**, 020408(R) (2005).

[53] Y. Motome and H. Tsunetsugu, Spin frustration and orbital order in vanadium spinels, J. Phys. Soc. Jpn. **74**, 208 (2005).

[54] G.-W. Chern, N. Perkins, and Z. Hao, Quantum 120° model on pyrochlore lattice: Orbital ordering in $MnV_2O_4$, Phys. Rev. B **81**, 125127 (2010).

[55] H. Tsunetsugu and Y. Motome, Magnetic transition and orbital degrees of freedom in vanadium spinels, Phys. Rev. B **68**, 060405(R) (2003).

[56] T. Mizokawa and A. Fujimori, Electronic structure and orbital ordering in perovskite-type 3d transition-metal oxides studied by Hartree-Fock band-structure calculations, Phys. Rev. B **54**, 5368 (1996).

[57] D. Baldomir, V. Pardo, S. Blanco-Canosa, F. Rivadulla, J. Rivas, P. Pineiro, and J. E. Arias, Pressure-induced metal-insulator transition in $MgV_2O_4$, Phys. B: Condens. Matter **403**, 1639 (2008).

[58] H. Lane, E. E. Rodriguez, H. C. Walker, Ch. Niedermayer, U. Stuhr, R. I. Bewley, D. J. Voneshen, M. A. Green, J. A. Rodriguez-Rivera, P. Fouquet, S.-W. Cheong, J. P. Attfield, R. A. Ewings, and C. Stock, Metastable antiphase boundary ordering in $CaFe_2O_4$, Phys. Rev. B **104**, 104404 (2021).

[59] H. Lane, M. Songvilay, R. A. Ewings, and C. Stock, Ex-citonic transverse and amplitude fluctuations in noncollinear and charge-ordered $RbFe^{2+}Fe^{3+}F_6$, Phys. Rev. B **106**, 054431 (2022).

[60] K. I. Kugel' and D. I. Khomskï, The Jahn-Teller effect and magnetism: transition metal compounds, Sov. Phys. Usp. **25**, 231 (1982).

[61] C. E. Housecroft and A. G. Sharpe, *Inorganic Chemistry* (Pearson, London, 2012).

[62] H. A. Jahn and E. Teller, Stability of polyatomic molecules in degenerate electronic states. I. Orbital degeneracy, Proc. R. Soc. London A **161**, 220 (1937).

[63] K. Takubo, J.-Y. Son, T. Mizokawa, H. Ueda, M. Isobe, Y. Matsushita, and Y. Ueda, Electronic structure of $AV_2O_4$ ($A = Li, Zn, and Cd$) studied by x-ray photoemission spectroscopy, Phys. Rev. B **74**, 155103 (2006).

[64] N. D. Mermin and H. Wagner, Absence of ferromagnetism or antiferromagnetism in one- or two-dimensional isotropic Heisenberg models, Phys. Rev. Lett. **17**, 1133 (1966).

[65] A. J. Browne and J. P. Attfield, Orbital molecules in vanadium oxide spinels, Phys. Rev. B **101**, 024112 (2020).

[66] S. P. Thompson, J. E. Parker, J. Potter, T. P. Hill, A. Birt, T. M. Cobb, F. Yuan, and C. C. Tang, Beamline I11 at Diamond: A new instrument for high resolution powder diffraction, Rev. Sci. Instr. **80**, 075107 (2009).

[67] T. Watanabe, T. Ishikawa, S. Hara, A. T. M. N. Islam, E. M. Wheeler, and B. Lake, Multiple lattice instabilities resolved by magnetic-field and disorder sensitivities in $MgV_2O_4$, Phys. Rev. B **90**, 100407(R) (2014).

[68] A. S. Wills, A new protocol for the determination of magnetic structures using simulated annealing and representational analysis (SARAh), Phys. B: Condens. Matter **276-278**, 680 (2000).

[69] G. Shirane, S. M. Shapiro, and J. M. Tranquada, *Neutron Scattering with a Triple-Axis Spectrometer: Basic Techniques* (Cambridge University Press, Cambridge, 2002).
[47] Y. Iqbal, T. Müller, P. Ghosh, M. J. P. Gingras, H. O. Jeschke, S. Rachel, J. Reuther, and R. Thomale, Quantum and classical phases of the pyrochlore Heisenberg model with competing interactions, Phys. Rev. X **9**, 011005 (2019).

[48] S. Zhang, H. J. Changlani, K. W. Plumb, O. Tchernyshyov, and R. Moessner, Dynamical structure factor of the three-dimensional quantum spin liquid candidate $NaCaNi_2F_7$, Phys. Rev. Lett. **122**, 167203 (2019).

[49] J. S. Gardner, B. D. Gaulin, S.-H. Lee, C. Broholm, N. P. Raju, and J. E. Greedan, Glassy statics and dynamics in the chemically ordered pyrochlore antiferromagnet $Y_2Mo_2O_7$, Phys. Rev. Lett. **83**, 211 (1999).

[50] J. W. Krizan and R. J. Cava, $NaCaNi_2F_7$: A frustrated high-temperature pyrochlore antiferromagnet with $S = 1$ $Ni^{2+}$, Phys. Rev. B **92**, 014406 (2015).

[51] K. W. Plumb, H. J. Changlani, A. Scheie, S. Zhang, J. W. Krizan, J. A. Rodriguez-Rivera, Y. Qiu, B. Winn, R. J. Cava, and C. L. Broholm, Continuum of quantum fluctuations in a three-dimensional $S = 1$ Heisenberg magnet, Nat. Phys. **15**, 54 (2019).

[52] S. Di Matteo, G. Jackeli, and N. B. Perkins, Orbital order in vanadium spinels, Phys. Rev. B **72**, 020408(R) (2005).

[53] Y. Motome and H. Tsunetsugu, Spin frustration and orbital order in vanadium spinels, J. Phys. Soc. Jpn. **74**, 208 (2005).

[54] G.-W. Chern, N. Perkins, and Z. Hao, Quantum 120° model on pyrochlore lattice: Orbital ordering in $MnV_2O_4$, Phys. Rev. B **81**, 125127 (2010).

[55] H. Tsunetsugu and Y. Motome, Magnetic transition and orbital degrees of freedom in vanadium spinels, Phys. Rev. B **68**, 060405(R) (2003).

[56] T. Mizokawa and A. Fujimori, Electronic structure and orbital ordering in perovskite-type 3d transition-metal oxides studied by Hartree-Fock band-structure calculations, Phys. Rev. B **54**, 5368 (1996).

[57] D. Baldomir, V. Pardo, S. Blanco-Canosa, F. Rivadulla, J. Rivas, P. Pineiro, and J. E. Arias, Pressure-induced metal-insulator transition in $MgV_2O_4$, Phys. B: Condens. Matter **403**, 1639 (2008).

[58] H. Lane, E. E. Rodriguez, H. C. Walker, Ch. Niedermayer, U. Stuhr, R. I. Bewley, D. J. Voneshen, M. A. Green,

J. A. Rodriguez-Rivera, P. Fouquet, S.-W. Cheong, J. P. Attfield, R. A. Ewings, and C. Stock, Metastable antiphase boundary ordering in $CaFe_2O_4$, Phys. Rev. B **104**, 104404 (2021).

[59] H. Lane, M. Songvilay, R. A. Ewings, and C. Stock, Excitonic transverse and amplitude fluctuations in noncollinear and charge-ordered $RbFe^{2+}Fe^{3+}F_6$, Phys. Rev. B **106**, 054431 (2022).

[60] K. I. Kugel' and D. I. Khomskï, The Jahn-Teller effect and magnetism: transition metal compounds, Sov. Phys. Usp. **25**, 231 (1982).

[61] C. E. Housecroft and A. G. Sharpe, *Inorganic Chemistry* (Pearson, London, 2012).

[62] H. A. Jahn and E. Teller, Stability of polyatomic molecules in degenerate electronic states. I. Orbital degeneracy, Proc. R. Soc. London A **161**, 220 (1937).

[63] K. Takubo, J.-Y. Son, T. Mizokawa, H. Ueda, M. Isobe, Y. Matsushita, and Y. Ueda, Electronic structure of $AV_2O_4$ ($A = Li, Zn, and Cd$) studied by x-ray photoemission spectroscopy, Phys. Rev. B **74**, 155103 (2006).

[64] N. D. Mermin and H. Wagner, Absence of ferromagnetism or antiferromagnetism in one- or two-dimensional isotropic Heisenberg models, Phys. Rev. Lett. **17**, 1133 (1966).

[65] A. J. Browne and J. P. Attfield, Orbital molecules in vanadium oxide spinels, Phys. Rev. B **101**, 024112 (2020).

[66] S. P. Thompson, J. E. Parker, J. Potter, T. P. Hill, A. Birt, T. M. Cobb, F. Yuan, and C. C. Tang, Beamline I11 at Diamond: A new instrument for high resolution powder diffraction, Rev. Sci. Instr. **80**, 075107 (2009).

[67] T. Watanabe, T. Ishikawa, S. Hara, A. T. M. N. Islam, E. M. Wheeler, and B. Lake, Multiple lattice instabilities resolved by magnetic-field and disorder sensitivities in $MgV_2O_4$, Phys. Rev. B **90**, 100407(R) (2014).

[68] A. S. Wills, A new protocol for the determination of magnetic structures using simulated annealing and representational analysis (SARAh), Phys. B: Condens. Matter **276-278**, 680 (2000).

[69] G. Shirane, S. M. Shapiro, and J. M. Tranquada, *Neutron Scattering with a Triple-Axis Spectrometer: Basic Techniques* (Cambridge University Press, Cambridge, 2002).