††thanks: Current address: Department of Physics, Northeastern University,
Boston, MA 02115

# Pulling order back from the brink of disorder: Observation of a nodal line
spin liquid and fluctuation stabilized order in K2IrCl6

Qiaochu Wang  Department of Physics, Brown University, Providence, Rhode
Island 02912, United States  Alberto de la Torre  Department of Physics, Brown
University, Providence, Rhode Island 02912, United States  Jose A. Rodriguez-
Rivera  NIST Center for Neutron Research, National Institute of Standards and
Technology, Gaithersburg, MD 20899, USA  Department of Materials Science and
Engineering, University of Maryland, College Park, MD 20742, USA  Andrey A.
Podlesnyak  Wei Tian  Adam A. Aczel  Masaaki Matsuda  Neutron Scattering
Division, Oak Ridge National Laboratory, Oak Ridge, Tennessee 37831, USA
Philip J. Ryan  Jong-Woo Kim  X-ray Science Division, Argonne National
Laboratory, Argonne, Illinois 60439, USA  Jeffrey G. Rau  Department of
Physics, University of Windsor, Windsor, Ontario, N9B 3P4, Canada  Kemp W.
Plumb<EMAIL_ADDRESS>Department of Physics, Brown University,
Providence, Rhode Island 02912, United States

###### Abstract

Competing interactions in frustrated magnets can give rise to highly
degenerate ground states from which correlated liquid-like states of matter
often emerge. The scaling of this degeneracy influences the ultimate ground
state, with extensive degeneracies potentially yielding quantum spin liquids,
while sub-extensive or smaller degeneracies yield static orders. A
longstanding problem is to understand how ordered states precipitate from this
degenerate manifold and what echoes of the degeneracy survive ordering. Here,
we use neutron scattering to experimentally demonstrate a new “nodal line”
spin liquid, where spins collectively fluctuate within a sub-extensive
manifold spanning one-dimensional lines in reciprocal space. Realized in the
spin-orbit coupled, face-centered cubic iridate K2IrCl6, we show that the sub-
extensive degeneracy is robust, but remains susceptible to fluctuations or
longer range interactions which cooperate to select a magnetic order at low
temperatures. Proximity to the nodal line spin liquid influences the ordered
state, enhancing the effects of quantum fluctuations and stabilizing it
through the opening of a large spin-wave gap. Our results demonstrate quantum
fluctuations can act counter-intuitively in frustrated materials: instead of
destabilizing ordering, at the brink of the nodal spin liquid they can act to
stabilize it and dictate its low-energy physics.

The promise of discovering novel states of matter governed by quantum and
thermal fluctuations drives a sustained interest in highly frustrated magnetic
materials. One of the oldest and most important examples is the Heisenberg
antiferromagnet (AFM) on the face-centered cubic (FCC) lattice. Anderson [2]
showed that classically its ground state degeneracy scales sub-extensively in
the system size, and includes co-planar spin spirals characterized by all
magnetic wavevectors spanning one-dimensional lines in reciprocal space. As
temperature is lowered, fluctuations of the spins become more confined to
states on these lines, forming a correlated state we term a “nodal line spin
liquid”. This spin-liquid is distinguished from a conventional paramagnet by
having growing correlations along two spatial directions as temperature is
lowered, but only short-range correlations along the third. Among correlated
paramagnets, including spiral spin-liquids, where the degenerate states span
surfaces in momentum space [6, 17], or classical spin liquids [58, 48, 41, 25,
44], where degenerate states span a volume, the nodal line spin liquid has the
least degeneracy, and thus sits closest to conventional magnetic order.

Due to its sub-extensive degeneracy, the nodal line spin liquid is susceptible
to order via quantum or thermal fluctuations – “order-by-disorder" [59, 23,
45, 24]. For the FCC Heisenberg AFM, thermal [19, 65] and quantum [55]
fluctuations each favor a sub-extensive set of collinear Néel states [2, 13],
both ultimately selecting so-called “type I” magnetic order characterized by a
$(1,0,0)$ propagation vector. However, the selection of this type I order
among the collinear states is extremely weak – with an energy difference of
only $\sim\\!0.25\%$ of the Heisenberg exchange, $J$, and many other
metastable collinear states lying within $\sim\\!1\%$ of $J$ in energy [55].
While this classical degeneracy and fluctuation-induced selection are robust
to nearest-neighbor anisotropic exchange interactions, it can be lifted
energetically by further neighbor couplings. For example, next nearest-
neighbor (NNN) Heisenberg exchange selects “type I” order when ferromagnetic,
but selects the “type III” state, a collinear Néel order with propagation
vector $(1,0.5,0)$ when antiferromagnetic [62].

Figure 1:  Dynamic correlations of the nodal line spin liquid. (a), (b),
Constant energy neutron scattering intensity at $T=6$ K integrated over 0.5 to
1 meV in (a) and 1 to 2 meV in (b). Lower right of (a), and (b), shows the
result from classical spin dynamics of our best fit FCC Heisenberg-Kitaev
model with exchange parameters $J=0.74$ meV and $K=0.15$ meV. (c), Momentum
and energy resolved magnetic excitation spectra of the nodal line spin liquid
at $T=6$ K. Momentum directions are integrated over $\pm 0.15$ r.l.u. . (d),
Classical spin dynamics at a temperature twice the model’s ordering
temperature computed as described in [56].

In each case, ordering out of the nodal-line-spin-liquid proceeds through a
multi-step process. First, order-by-disorder singles out the sub-manifold of
discrete collinear states from the continuous manifold of the nodal line spin
liquid. Despite their near degeneracy, each of these local minima exhibits
conventional magnetic order with gapped spin-wave excitations. Second, a
specific collinear state among the local minima is selected by either
fluctuations or energetic effects. While the ordered state is ultimately
determined by the second stage, the low-energy physics – including the
excitations and thermodynamic properties – are dictated by the fluctuation-
induced selection arising in the first stage. This leads to a counter-
intuitive result: while the final order-by-disorder selection of the ground
state is very weak, each of the collinear states is stable, with a large spin-
wave gap generated by strong quantum fluctuations along the nodal line
directions of the parent spin liquid. This large gap then suppresses
fluctuations of the ordered moment. Thus, despite the difficulty of
fluctuations in selecting a unique ground state, they act to stabilize
whichever particular ordering is singled out by sub-leading perturbations, and
ultimately dictates the physics of the material.

In this article, we show that such physics is realized in the FCC iridate
K2IrCl6. First, we establish that the magnetic correlations form one-
dimensional rods in reciprocal space characteristic of a nodal line spin
liquid as the magnetic transition is approached from above. We then present a
comprehensive series of neutron scattering measurements in the ordered state
that, when combined with theoretical modeling of the dynamic spin
correlations, conclusively demonstrate that K2IrCl6 is a nearly ideal
realization of the Heisenberg-Kitaev model on the FCC lattice. This model is
highly frustrated, with a specific ordering only weakly selected by NNN
interactions or quantum fluctuations. Despite the appearance of static
magnetic order, as observed by diffraction measurements, a hard excitation
gap, and long-lived spin waves, linear spin wave theory (LSWT) theory
completely fails as a description of magnetic excitations in K2IrCl6. Instead,
we find that the physics of the excitation spectrum is dominated by
fluctuation effects that generate a large, $\sim 0.7$ meV, excitation gap –
about 30 % of the overall $\sim 2.5$ meV magnon bandwidth. We thus show how
proximity to the nodal line spin liquid enhances quantum fluctuations and,
contrary to conventional lore, stabilizes magnetic order by quenching out the
zero-point corrections to the ordered moment.

While much theoretical and experimental effort has been made on FCC magnets
[4, 20, 21, 34, 57, 38, 1, 5], the antifluorite-type K2IrCl6 has sustained
more than 70 years of interest [43, 4, 18, 11, 29, 20, 21, 34, 22, 57, 27, 39,
36, 61, 42, 30, 47, 7]. Despite such effort, the magnetic Hamiltonian of
K2IrCl6 has yet to be resolved due to an absence of any inelastic neutron
scattering data. In addition to the presence of geometric frustration, the
magnetic degrees of freedom in K2IrCl6 are spin-orbital entangled $j=1/2$
states [28], supporting bond anisotropic Kitaev interactions, which are
believed to be significant relative to the Heisenberg exchange [11, 30]. Early
magnetic neutron diffraction experiments revealed a type III antiferromagnetic
order below $\sim\\!3$ K [27]. Based on this observation and a consistent
explanation of the bulk magnetic properties, K2IrCl6 has long been considered
to represent a near ideal FCC lattice with weak antiferromagnetic NNN
interactions [30, 47]. However, recent electron spin resonance measurements
indicate NNN terms may only play a perturbative role and order-by-disorder
physics may be essential [7] in explaining its properties. To date, there are
no reported measurements of momentum resolved dynamic correlations
characterizing the excitation spectrum of this highly frustrated anisotropic
magnet.

##  I Nodal line spin liquid

We first present inelastic neutron scattering results that demonstrate two-
dimensional long-range correlations above the Néel transition at $T=6$ K, well
above our observed ordering temperature $T_{\rm N}=3.2$ K, through
characteristic rods of intensity in the dynamical structure factor averaged
over a broad energy window [Fig. 1]. Such highly anisotropic, two dimensional,
correlations arising from an underlying isotropic, three dimensional, magnetic
Hamiltonian are the direct experimental signature of a nodal-line spin liquid.

The rods of intensity extend along equivalent $[1,q,0]$ reciprocal space
directions characteristic of the set of co-planar spin spirals of the nodal
line spin liquid. Momentum-energy slices shown in Fig. 1(c) reveal a continuum
of scattering extending from these rods, with intensity concentrated near
$\bm{Q}=(1,0,0)$ momentum transfers and extending to 2 meV. There is no
evidence for a build-up of critical intensity near the $\bm{Q}=(1,0.5,0)$
wavevector of the type III magnetic ordering that appears below $T_{\rm N}$.
Instead, dynamic magnetic correlations are localized along lines in reciprocal
space, identifying K2IrCl6 as a _nodal-line spin-liquid_ state in this
temperature regime.

We compare these correlations with expectations from the classical spin
dynamics of the FCC Heisenberg-Kitaev (HK) model. As shown in Fig. 1(a), the
low energy momentum dependent correlations we observe are well captured by
classical spin dynamics of the FCC HK model within its paramagnetic phase that
we will further establish below as the appropriate Hamiltonian for K2IrCl6.
Although the classical calculation can reproduce much of the observed momentum
dependent scattering, there is a notable discrepancy for energy transfers
above 1 meV where the data reveals an intensity maxima near $\bm{Q}=(1,0,0)$
not seen in the theoretical calculation [Fig. 1(b,d)]. This discrepancy
reflects strong quantum fluctuations in K2IrCl6 that are not captured in the
classical spin dynamics and redistributes magnetic fluctuations to higher
energies.

Since the degeneracy of the nodal line spin liquid scales subextensively,
either fluctuations, small additional magnetic interactions, or both are
expected to select a long range ordered state as temperature is lowered.
Consistent with these expectations, K2IrCl6 undergoes a magnetic ordering
transition at $T_{\rm N}=3.2$ K. However, we will show below that, despite a
well formed type III Néel order, the excitations out of the ordered state
cannot be captured on even a qualitative level by semi-classical LSWT starting
from the ordered state and, as in the paramagnetic phase, large corrections
from quantum fluctuations are required to capture the physics of this
material.

Figure 2: Spin waves in the ordered phase. Measured (left) and calculated
(right) energy–momentum slices along (a), $[1,k,0]$ and (b), $[0,k,0]$
integrated over $h\pm 0.15$ r.l.u. and $l\pm 0.15$ r.l.u.. (c)-(f), Constant
energy slices across the magnon spectra integrated over $\omega\pm 0.15$ meV.
The measured spectrum is shown in the upper left and the calculated spectrum
is shown in the lower right. (g), Constant $\bm{Q}$ cuts at high symmetry
wavevectors indicated by colored arrows at the top of (a) and integrated over
$\pm~{}0.1$ r.l.u.. Data in panels (a)–(g) was collected at $T=300$ mK and
errorbars represent one standard deviation. Spectrum calculated at
$O(1/S^{2})$ using self-consistent non-linear spin-wave theory for an FCC
Heisenberg-Kitaev model with type III magnetic order (domain averaged) using
best-fit exchange parameters $J=0.74$ meV and $K=0.15$ meV.

##  II Fluctuation induced spin wave gap

Fig. 2 shows an overview of the magnetic neutron spectra at $T=300$ mK, well
below the onset of type III $(0.5,0,1)$ magnetic order. There are sharp spin
waves present throughout the Brillouin zone spanning a bandwidth of $2.5$ meV
with a $0.7\pm 0.05$ meV gap at the $(1,0,0)$ dispersion minimum. Although the
strong neutron absorption from Ir places limitations on the signal-to-noise
and achievable energy resolution of our experiments, the spectrum we observe
is dominated by resolution limited single magnon excitations, and any multi-
magnon scattering was not measurable. The observed spectrum is surprising for
two reasons. First, the dispersion minimum does not occur at the wavevector of
static magnetic order $(1,0.5,0)$, but instead near $(1,0,0)$, and second for
the generic symmetry allowed spin model relevant to K2IrCl6 LSWT predicts a
gapless pseudo-Goldstone mode in the type III phase [46, 56].

The presence of a large excitation gap could potentially be explained within
LSWT by the lowering of the crystal symmetry from cubic to (say) tetragonal,
as might arise from a weak structural distortion at $T_{\rm N}$ [1, 60].
Although all previous diffraction measurements have reported K2IrCl6 as cubic
for $T<T_{\rm N}$ [30, 47], there is spectroscopic evidence that K2IrCl6 may
exhibit local structural distortions [33]. Furthermore, recent high-resolution
single crystal x-ray diffraction measurements reveal a 0.02% compression of
the c-axis accompanied by a 0.014% expansion of $a$ and $b$ axes coinciding
with the magnetic transition [60]. To conclusively identify this large spin-
wave gap as fluctuation induced, we first determine whether these distortions
could potentially account for it’s observed magnitude.

Since $j=1/2$ pseudospins have both a spin and orbital contribution to their
magnetic moment, they can couple more easily with lattice deformations of the
appropriate symmetry [35]. Due to this coupling, the type III ordering will
necessarily generate structural deformations through a pseudo-Jahn-Teller
mechanism, which is consistent with the observed tetragonal distortion [60]
and the reported anisotropic magnetization density in K2IrCl6 [37]. Such
magneto-elastic couplings will also give rise to spatially anisotropic
exchange interactions, $D$, that can induce a magnon gap $\propto\sqrt{JD}$
[35], where $J$ is the (dominant) Heisenberg exchange. We estimate the
magnetic anisotropy from the measured tetragonal strain $\varepsilon\approx
2\times 10^{-4}$ [60] and spin-lattice coupling constant
$\tilde{g}\simeq\left(t^{2}/U\right)3/2\lambda\left(J_{\rm H}/U\right)$ where
$t$, $U$, $J_{\rm H}$, and $\lambda$ are the hopping amplitude, Coulomb
repulsion, Hund’s coupling, and spin-orbit coupling, respectively [35]. Based
on similar estimates in Sr2IrO4 [35], but with a nearest-neighbor exchange
constant $J\sim t^{2}/U$ in K2IrCl6 that is a factor of $\sim 100$ smaller, we
estimate $\tilde{g}\approx 0.25$ meV and thus
$D\sim\tilde{g}\varepsilon\approx 0.025\ \mu$eV. This yields a distortion-
induced gap $\sim 0.05$ $\mu$eV – _four-orders of magnitude too small_ to
account for the observed 0.7 meV value [56]. As an additional check, we have
carried out a LSWT calculation including all symmetry-allowed terms with only
tetragonal symmetry [56]. This includes different exchanges within ($J$, $K$)
and perpendicular to ($J_{\perp},K_{\perp}$) the basal plane (perpendicular to
the tetragonal axis). Fitting a tetragonal model that can account for the
measured gap yields significant differences in the exchanges within the plane
and in the perpendicular direction; the Heisenberg exchange yielding a modest
$J_{\perp}/J=1.29$, but a much larger $K^{\prime}_{\perp}/K=1.93$ for the
Kitaev exchange [56]. Given the smallness of the observed tetragonal
distortion and our estimate of the magneto-elastic coupling, we find this
level of spatial anisotropy in the exchanges implausible and rule it out as an
explanation for the large spin-wave gap in our data.

Given that magnetic excitations in the ordered state of K2IrCl6 cannot be
explained qualitatively by LSWT when cubic symmetry is maintained, we must
then consider the effects of quantum fluctuations beyond LSWT. Including
interactions between the spin-waves can lift the pseudo-Goldstone modes,
inducing a gap through quantum fluctuations [54, 46]. To account for these
fluctuation induced gaps, we have computed the dynamical structure factor for
the FCC HK model (as measured by INS) about the type III ordering using
$O(1/S^{2})$ self-consistent non-linear spin wave theory (NLSWT) [56]. As
shown in Fig. 2, we find excellent quantitative agreement with our
experimental data – including the large gap – using only nearest-neighbor
Heisenberg, $J=0.74$ meV, and Kitaev, $K=0.15$ meV ($K/J=0.2$), exchanges.
These exchange parameters are comparable to DFT predictions [30] and estimates
from ESR measurements [7]. Neither the symmetric anisotropy $\Gamma$ nor NNN
interactions are required to describe the inelastic neutron data for K2IrCl6,
suggesting they are small relative to $J$ and $K$.

However, despite providing a quantitative description of the excitations in
K2IrCl6, without further neighbor interactions our best-fit model does not
clearly predict a type III ground state. Even including fluctuations to
$O(1/S^{2})$ the energy difference between type I and III states is expected
to be incredibly small [54]. This close competition between type I and type
III states indicates the mechanism selecting type III order in K2IrCl6 is
delicate and the ordering transition in K2IrCl6 deserves closer scrutiny.

##  III Competing type I and III magnetic order

Heat capacity measurements on our samples display a single sharp peak at $3.2$
K, in agreement with the critical temperature extracted from the $(1,0.5,0)$
magnetic order parameter [Fig. 3(a,b)]. However, a broader reciprocal space
survey on the same samples revealed additional weak magnetic Bragg reflections
at $\bm{Q}=(1,0,0)$, corresponding to type I order, with an onset temperature
of $T_{\rm N}^{\rm I}=3.4$ K [Fig. 3(b)]. The magnetic origin of $(1,0,0)$
reflections was confirmed through energy and polarization-dependent resonant
x-ray diffraction that indicates an absence of any charge signal [Fig.
3(d-f)]. The fine momentum space resolution x-ray measurements also revealed
that $(0.5,0,5)$ magnetic Bragg peaks, corresponding to type III Néel order,
are resolution limited, while the $(0,0,3)$ magnetic peak, corresponding to
type I order, has a finite correlation length of 877 Å. Although there is no
heat capacity anomaly visible at $T_{\rm N}^{\rm I}$, we find substantial
magnetic entropy extending up to $T\sim 15$ K that is consistent with the
onset of short-range type I order at $3.5$ K.

Figure 3:  Coexistence of type I and type III collinear magnetic orders in
K2IrCl6. (a), Measured zero-field specific heat $C_{p}$ and estimated lattice
contribution $C_{\rm ph}$. (b), Temperature dependence of $(0.5,1,0)$ (type
III) and $(0,1,0)$ (type I) intensities obtained from neutron diffraction.
(c), Polarized neutron diffraction of a type III peak at 1.5 K indicating the
absence of spin component along the $[0,0,1]$ direction. (d), Energy dependent
elastic X-ray scattering demonstrating that both peaks arise from a magnetic
resonance. Resonant x-ray rocking curves on (e), type III and (f), type I
peaks at 2 K measured on the same crystal as for (a). Errorbars in all panels
represent one standard deviation.

A simultaneous appearance of $(1,0.5,0)$ and $(1,0,0)$ magnetic peaks is
consistent with either a coexistence of collinear magnetic phases across
multiple domains, or a multi-$Q$ order. Polarized neutron diffraction
presented in Fig. 3(c) rules out the possibility of non-coplanar multi-$Q$
states. The spin-flip neutron intensity of the $\bm{Q}=(1,-0.5,0)$ Bragg
reflection is sensitive to both in and out-of-plane components of the
magnetization for guide fields $\bm{P}$ parallel to $\bm{Q}$, and in-plane
components of the magnetization, along $[1,2,0]$, for $\bm{P}\perp\bm{Q}$.
While the non-spin flip neutron intensity is only sensitive to out-of-plane
components of the magnetization, along $[0,0,1]$, when $\bm{P}\perp\bm{Q}$. We
observed identical $(1,-0.5,0)$ peak intensities for $\bm{P}\parallel\bm{Q}$
and $\bm{P}\perp\bm{Q}$ spin-flip channels and no magnetic signal in the non-
spin flip component, ruling out any out-of-plane component, as expected for a
collinear type III Néel state [Fig. 3(c)]. We also find that the polarized
diffraction from $\bm{Q}=(1,0,0)$ magnetic peaks is most consistent with a
collinear arrangement of magnetic moments that are parallel to the $[0,0,1]$
direction, expected for a type I Néel state [56]. Thus, the diffraction data
support a coexistence of type I and type III magnetic orders across different
spatial regions of the sample. This is consistent with expectations from
order-by-disorder, which strongly favor collinear single-$Q$ states. From the
integrated neutron diffraction Bragg intensities, we extract a relative phase
fraction $N_{\rm I}/N_{\rm III}\approx 11\%$. This phase fraction did not
change between measurements on two different samples from different crystal
growths indicating that the small fraction of type I phase is likely intrinsic
and does not arise from quenched disorder in our samples.

The appearance of both type I and III order near $T_{\rm N}$ in K2IrCl6
implicates weak and competing selection effects from fluctuations or further-
neighbor exchange for these nearly degenerate states. To better understand the
role of these selection mechanisms near $T_{\rm N}$, we consider energetic
selection via a NNN exchange $J_{2}$ in the classical model. When $J_{2}/J=0$,
we find that thermal order by disorder selects a type I order for the FCC HK
model [19, 55, 56]. A small $J_{2}/J\sim 0.02$ selects type III order instead,
onsetting at a nearly identical Néel temperature to the $J_{2}=0$ case [56].
Together, our data and modeling demonstrate that K2IrCl6 uniquely displays a
competition between energetic and entropic mechanisms for order selection.
Thermal fluctuations initially favor type I order, but as temperature is
lowered and the influence of entropy diminishes, quantum fluctuations or weak
further neighbor interactions select type III order. Such a competition
between weak selection mechanisms for the magnetic order in K2IrCl6 is
consistent with nearly negligible NNN interactions and underscores the
dominant role played by quantum fluctuations in this material.

##  IV Discussion

The conventional lore of frustrated quantum magnetism is that degenerate and
competing ground states suppress the tendency towards magnetic long range
order. In the quantum limit of $S=1/2$, quantum fluctuations among degenerate
ground states persist as $T\rightarrow 0$ and act to reduce the static ordered
moment size. Experimentally, this manifests as a transfer of neutron spectral
weight from the magnetic Bragg peak (order parameter) to the excitations [32].
K2IrCl6 presents an important counterpoint to these expectations where despite
the appearance of a well formed static Néel state and long lived magnon
excitations, higher order (in $1/S$) corrections from quantum fluctuations are
essential to describe the spectra on even a qualitative level. The ordered
moment reduction can be quantified through the energy distribution of neutron
intensity. The momentum and energy integrated intensity gives the total moment
$S(S+1)$, while the elastic intensity for a fully ordered moment is $S^{2}$.
Longitudinal fluctuations will reduce the static moment so that elastic
intensity is $(S-\Delta S)^{2}$ and the ratio of elastic to total integrated
intensity is $R=(S-\Delta S)^{2}/S(S+1)$. For K2IrCl6, we find $R\approx
0.31$, or $\Delta S\rightarrow 0$, which is consistent with the small moment
reduction predicted by NLSWT and demonstrates that by generating a large gap
for longitudinal magnetic fluctuations, quantum fluctuations in K2IrCl6 act to
stabilize the ordered moment.

Above $T_{\rm N}$, K2IrCl6 is a unique realization of a new classical spin
liquid: the nodal line spin liquid. This correlated paramagnet is similar to
the spiral-spin liquid recently observed in MnSc2S4 [17]. Both have a
subextensive degeneracy and are susceptible to thermal order-by-disorder, but
fluctuations in the nodal line liquid are confined to one dimensional lines
instead of spiral surfaces. Despite this reduced degeneracy, the role of
quantum fluctuations is significantly enhanced in K2IrCl6. Magnetic ordering
in MnSc2S4 is controlled by further neighbor, dipolar, and anisotropic
exchange interactions, the excitation spectrum is fully described by LSWT [31,
16], and there is no experimental support for an order by disorder transition
[31, 15, 16]. This difference likely arises from the small $j=1/2$ moments
(Mn2+ is $S=5/2$) and negligible dipolar interactions in K2IrCl6.

It is the scale of quantum fluctuations in K2IrCl6 that distinguish it as an
important frustrated magnet. These can be quantified by the relative size of
the fluctuation induced magnon gap compared to the bandwidth that exceeds 30 %
in K2IrCl6. For comparison, quantum order-by-disorder spin wave gaps have been
found in the garnet Fe2Ca3(GeO4)3 [8], the pyrochlore Er2Ti2O7 [9, 51, 52,
50], and the honeycomb lattice CoTiO3 [63, 14]. In all of these, the spin wave
gap to bandwidth ratio is less than 10% and nearly all details of the magnon
spectra besides the fluctuation induced gaps are reproduced by linear spin
wave theory, while this fails qualitatively for K2IrCl6.

The scale of quantum fluctuations in K2IrCl6 is further reflected by the
dynamics of the nodal line spin liquid phase that are only partially
reproduced by a classical spin dynamics calculation. This is distinct from
other correlated paramagnets. For example even the excitation spectrum of the
highly frustrated Heisenberg model on the pyrochlore lattice can be
qualitatively captured by classical spin dynamics calculations renormalized by
a classical to quantum correspondence factor [64]. Departures from
renormalized classical spin dynamics calculations in the nodal line spin
liquid realized by K2IrCl6 calls for further work to understand the quantum to
classical crossover in classical spin liquids.

## Acknowledgments

Work at Brown University was supported by the U.S. Department of Energy,
Office of Basic Energy Sciences, under Grant No. DE-SC0021223. A portion of
this research used resources at the Spallation Neutron Source and High Flux
Isotope Reactor, the DOE Office of Science User Facilities operated by the Oak
Ridge National Laboratory. Access to MACS was provided by the Center for High
Resolution Neutron Scattering, a partnership between the National Institute of
Standards and Technology and the National Science Foundation under Agreement
No. DMR-2010792. This research used resources of the Advanced Photon Source, a
U.S. Department of Energy (DOE) Office of Science user facility operated for
the DOE Office of Science by Argonne National Laboratory under Contract No.
DE-AC02-06CH11357. Work at the University of Windsor was funded by the Natural
Sciences and Engineering Research Council of Canada (NSERC) (Funding Reference
No. RGPIN-2020-04970). We acknowledge the use of computational resources
provided by Digital Research Alliance of Canada.

##  Appendix A Methods

###  A.1 Sample Synthesis and Characterization

Commercially available powder of K2IrCl6 was obtained from Fisher Scientific.
K2IrCl6 crystals were prepared by slow evaporation from a saturated solution
of K2IrCl6 in dilute hydrochloric acid (pH $\sim 3$). The evaporation process
was controlled to maintain a temperature of 32 degrees Celsius. Larger
crystals were obtained by successive seeding of solution growths.

Specific heat and magnetization were measured in a Quantum Design Physical
Properties Measurement System (PPMS). Specific heat was measured upon warming
in zero magnetic field after a zero-field cooling. Additional sample
characterization is provided in the Supplemental information [56].

###  A.2 Neutron Scattering

Inelastic neutron scattering measurements were performed on a sample comprised
of 11 co-aligned crystals with total mass 0.3 g using the Cold Neutron Chopper
Spectrometer (CNCS) at the Spallation Neutron Source (SNS) at Oak Ridge
National Laboratory (ORNL) and the MACS spectrometer [49] at the NIST Center
for Neutron Research. The neutron momentum transfer is indexed using the
Miller indices of the cubic unit cell,
$(h,k,\ell)\\!=\\!(2\pi/a,2\pi/a,2\pi/a)$, where $a\\!=\\!9.66$ Å and all
inelastic neutron scattering data has been corrected for the energy dependent
neutron absorption from Ir [56].

On CNCS, we used an incident neutron energy of $E_{i}=3.32$ meV in the high-
flux configuration with a chopper frequency of 180 Hz to give an elastic line
energy resolution of 0.17 meV (FWHM). A He cryostat and He3 insert was used to
cool the sample to a base temperature of $T=0.3$ K. Data were collected for a
360∘ rotation of the sample with 1 degree step size. Measured neutron count
rates were placed into absolute units of the neutron scattering cross-section
using incoherent elastic scattering from the sample. The scale factor for
conversion to absolute units was additionally cross-checked against the
integrated intensity for (200) and (220) nuclear Bragg reflections. All data
reduction and analysis was carried out using the Mantid software suite [3].

Measurements on MACS were conducted with the sample oriented in the
$(h,h,\ell)$ scattering plane using a double-focusing configuration and fixed
final energy of 3.7 meV with a BeO filter after the sample and no incident
beam filter. The data was corrected for contamination from high-order
harmonics in the incident beam neutron monitor.

Unpolarized neutron diffraction measurements were carried out on the HB-1A
diffractometer of the High Flux Isotope Reactor(HFIR) at ORNL with
collimations of 40’-40’-40’-80’ and fixed incident energy $E_{i}=14.5$ meV.
Polarized neutron diffraction measurements were carried out on the HB-1 triple
axis spectrometer of HFIR at ORNL with collimations of 48’-80’-60’-open and
fixed incident energy $E_{i}=13.5$ meV. Both experiments on HB-1A and HB-1
were conducted on the same single crystal with a mass of 42.5 mg aligned in
the $(h,k,0)$ scattering plane.

###  A.3 X-ray Scattering

Resonant elastic x-ray scattering measurements at the Ir L3 edge,
$E_{i}=11.215$ keV, were performed on Beamline 6-ID-B at the Advanced Photon
Source. A Pilatus 100K detector with 487$\times$195 pixels of 175$\times$175
micron size, was utilized for the measurement. Polarization analysis was
performed using the $(008)$ reflection of a pyrolytic graphite analyzer and
the Cyberstar scintillator detector. A Joule-Thompson displex cryostat was
used to reach a base temperature of 2 K.

###  A.4 Heisenberg-Kitaev model on the FCC lattice

We model the $j=1/2$ doublets of the Ir4+ atoms as effective $S=1/2$ spins,
$\bm{S}_{i}$, on an FCC lattice with an isotropic $g$-factor. We consider a
model with nearest neighbor Heisenberg exchange and Kitaev exchanges [29, 10,
1]

$\sum_{\braket{ij}_{\gamma}}\left[J\bm{S}_{i}\cdot\bm{S}_{j}+KS^{\gamma}_{i}S^{\gamma}_{j}\right]$
(A1)

where we have divided the bonds of the FCC lattice into three types: $x$, $y$
and $z$, depending on whether they lie in the $yz$, $zx$ or $xy$ planes. When
NNN interactions are included we consider only a Heisenberg exchange
$J_{2}\sum_{\braket{\braket{ij}}}\bm{S}_{i}\cdot\bm{S}_{j}$. Our best fit
parameters (as described in the main text) correspond to $J=0.74$ meV and
$K=0.15$ meV.

###  A.5 Self-consistent non-linear spin-wave theory

We consider a semi-classical expansion about the type III ordered state,
starting from the $S\rightarrow\infty$ limit. We use the Holstein-Primakoff
representation [26] of the spin

$\displaystyle\bm{S}_{i}\equiv$
$\displaystyle\sqrt{S}\left[\left(1-\frac{n_{\bm{r}\alpha}}{2S}\right)^{1/2}{a}^{\phantom{\dagger}}_{\bm{r}\alpha}\bm{\hat{e}}_{\alpha,-}+{a}^{\dagger}_{\bm{r}\alpha}\left(1-\frac{n_{\bm{r}\alpha}}{2S}\right)^{1/2}\bm{\hat{e}}_{\alpha,+}\right]$
$\displaystyle+\left(S-n_{\bm{r}\alpha}\right)\bm{\hat{e}}_{\alpha,0},$ (A2)

where
$n_{\bm{r}\alpha}\equiv{a}^{\dagger}_{\bm{r}\alpha}{a}^{\phantom{\dagger}}_{\bm{r}\alpha}$
and $i\equiv\bm{r},\alpha$ labels the unit cell and sublattice of a spin. For
the collinear type-III orders a four-sublattice unit cell is sufficient for
each of the domains. The vectors
$\bm{\hat{e}}_{\alpha,\pm}\equiv(\bm{\hat{x}}_{\alpha}\pm
i\bm{\hat{y}}_{\alpha})/\sqrt{2}$ and
$\bm{\hat{e}}_{\alpha,0}\equiv\bm{\hat{z}}_{\alpha}$ define a local frame with
$\bm{\hat{z}}_{\alpha}$ being along the ordering directions of the type III
state. Expanding in powers of $1/S$ then yields a semi-classical expansion.
Including terms up to $O(1)$ we find

$H=NS(S+1)\epsilon_{\rm cl}+H_{2}+H_{4},$

where the $O(S^{2})$ part, $\epsilon_{\rm cl}$, is the classical ground state
energy density and we define the $O(S)$ and $O(1)$ parts in symmetrized form
as

$\displaystyle{H}_{2}=\frac{1}{2}\sum_{\alpha\beta}\sum_{\bm{k}}$
$\displaystyle\left[{A}_{\bm{k}}^{\alpha\beta}{a}^{\dagger}_{\bm{k}\alpha}{a}^{\phantom{\dagger}}_{\bm{k}\beta}+{A}_{-\bm{k}}^{\beta\alpha}{a}^{\phantom{\dagger}}_{-\bm{k}\alpha}{a}^{\dagger}_{-\bm{k}\beta}+\right.$
(A3a)
$\displaystyle\left.\left({B}^{\alpha\beta}_{\bm{k}}{a}^{\dagger}_{\bm{k}\alpha}{a}^{\dagger}_{-\bm{k}\beta}+\bar{B}^{\alpha\beta}_{\bm{k}}{a}^{\phantom{\dagger}}_{-\bm{k}\beta}{a}^{\phantom{\dagger}}_{\bm{k}\alpha}\right)\right],$
$\displaystyle{H}_{4}=\frac{1}{N_{c}}\sum_{\alpha\beta\mu\nu}\sum_{\bm{k}\bm{k}^{\prime}\bm{q}}$
$\displaystyle\left[\frac{1}{(2!)^{2}}{V}_{\bm{k}\bm{k}[\bm{q}]}^{\alpha\beta\mu\nu}{a}^{\dagger}_{\bm{k}+\bm{q},\alpha}{a}^{\dagger}_{\bm{k}^{\prime}-\bm{q},\beta}{a}^{\phantom{\dagger}}_{\bm{k}^{\prime}\mu}{a}^{\phantom{\dagger}}_{\bm{k}\nu}+\right.$
(A3b)
$\displaystyle\left.\frac{1}{3!}\left({D}_{\bm{k}\bm{k}^{\prime}\bm{q}}^{\alpha\beta\mu\nu}{a}^{\dagger}_{\bm{k}\alpha}{a}^{\dagger}_{\bm{k}^{\prime}\beta}{a}^{\dagger}_{\bm{q}\mu}{a}^{\phantom{\dagger}}_{\bm{k}+\bm{k}^{\prime}+\bm{q},\nu}+\textrm{h.c.}\right)\right].$

In terms of the exchange matrices expressed in these local frames,
$\mathcal{J}^{\mu\mu^{\prime}}_{\bm{\delta},\alpha\alpha^{\prime}}\equiv{\bm{\hat{e}}}^{\intercal}_{\alpha,\mu}{\bm{J}}_{\bm{\delta},\alpha\alpha^{\prime}}{\bm{\hat{e}}}^{\phantom{\dagger}}_{\alpha^{\prime},\mu^{\prime}}$,
the coefficients are

$\displaystyle{A}_{\bm{k}}^{\alpha\beta}$
$\displaystyle=S\left(\mathcal{J}^{+-}_{\bm{k},\alpha\beta}-\delta_{\alpha\beta}\sum_{\mu}\mathcal{J}^{00}_{\bm{0},\alpha\mu}\right),$
(A4a) $\displaystyle{B}_{\bm{k}}^{\alpha\beta}$
$\displaystyle=S\mathcal{J}^{++}_{\bm{k},\alpha\beta},$ (A4b)
$\displaystyle{V}^{\alpha\beta\mu\nu}_{\bm{k}\bm{k}^{\prime}[\bm{q}]}$
$\displaystyle=\left(\delta_{\alpha\mu}\delta_{\beta\nu}\mathcal{J}^{00}_{\bm{k}-\bm{k}^{\prime}+\bm{q},\alpha\beta}\delta_{\alpha\nu}\delta_{\beta\mu}\mathcal{J}^{00}_{\bm{q},\alpha\beta}\right)$
(A4c)
$\displaystyle-\left(\delta_{\mu\nu}\delta_{\mu\beta}\mathcal{J}^{+-}_{\bm{k}+\bm{q},\alpha\nu}+\delta_{\alpha\beta}\delta_{\alpha\mu}\mathcal{J}^{+-}_{\bm{k},\alpha\nu}\right),$
$\displaystyle{D}^{\alpha\beta\mu\nu}_{\bm{k}\bm{k}^{\prime}\bm{q}}$
$\displaystyle=-\frac{3}{4}\left(\delta_{\alpha\mu}\delta_{\alpha\nu}\mathcal{J}^{++}_{\bm{k}^{\prime},\beta\alpha}+\delta_{\mu\beta}\delta_{\nu\beta}\mathcal{J}^{++}_{\bm{k},\alpha\beta}\right),$
(A4d)

where the interaction vertices have been left unsymmetrized for brevity.

At leading order in perturbation theory, these magnon interaction terms
renormalize the LSWT spectrum, giving corrections to $\delta{\bm{A}}_{\bm{k}}$
and $\delta{\bm{B}}_{\bm{k}}$ to ${\bm{A}}_{\bm{k}}$ and ${\bm{B}}_{\bm{k}}$ .
We treat these corrections self-consistently by writing

$\displaystyle\delta
A^{\alpha\beta}_{\bm{k}}=\frac{1}{N_{c}}\sum_{\bm{q}\mu\nu}$
$\displaystyle\Big{[}V^{\alpha\mu\nu\beta}_{\bm{k}\bm{q}[\bm{0}]}\braket{{a}^{\dagger}_{\bm{q}\mu}{a}^{\phantom{\dagger}}_{\bm{q}\nu}}_{\textrm{MF}}+$
$\displaystyle\frac{1}{2}\left(D^{\alpha\mu\nu\beta}_{\bm{k},-\bm{q},\bm{q}}\braket{{a}^{\dagger}_{-\bm{q}\mu}{a}^{\dagger}_{\bm{q}\nu}}_{\textrm{MF}}+\bar{D}^{\beta\mu\nu\alpha}_{\bm{k},\bm{q},-\bm{q}}\braket{{a}^{\phantom{\dagger}}_{\bm{q}\mu}{a}^{\phantom{\dagger}}_{-\bm{q}\nu}}_{\textrm{MF}}\right)\Big{]},$
$\displaystyle\delta
B^{\alpha\beta}_{\bm{k}}=\frac{1}{N_{c}}\sum_{\bm{q}\mu\nu}$
$\displaystyle\left[D^{\mu\alpha\beta\nu}_{\bm{q},\bm{k},-\bm{k}}\braket{{a}^{\dagger}_{\bm{q}\mu}{a}^{\phantom{\dagger}}_{\bm{q}\nu}}_{\textrm{MF}}+\frac{1}{2}V^{\alpha\beta\nu\mu}_{\bm{q},-\bm{q},[\bm{k}-\bm{q}]}\braket{{a}^{\phantom{\dagger}}_{\bm{q}\mu}{a}^{\phantom{\dagger}}_{-\bm{q}\nu}}_{\textrm{MF}}\right].$

where $N_{c}$ is the number of unit cells and the averages
$\braket{\cdots}_{\textrm{MF}}$ are with respect to the quadratic Hamiltonian,
$H_{2}+\delta H_{2}$, where we have replaced
${\bm{A}}_{\bm{k}}\rightarrow{\bm{A}}_{\bm{k}}+\delta{\bm{A}}_{\bm{k}}$ and
${\bm{B}}_{\bm{k}}\rightarrow{\bm{B}}_{\bm{k}}+\delta{\bm{B}}_{\bm{k}}$. These
equations are solved self-consistently (via iteration) to obtain an
(effective) corrected quadratic Hamiltonian.

###  A.6 Calculation of inelastic neutron scattering intensity

The inelastic neutron scattering intensity is computed theoretically in terms
of the spin-spin-correlation function

$\mathcal{S}_{\mu\nu}(\bm{k},\omega)\equiv\frac{1}{2\pi N}\int dte^{-i\omega
t}\braket{M^{\mu}_{-\bm{k}}M^{\nu}_{\bm{k}}(t)}$ (A6)

where $\bm{M}_{\bm{k}}\equiv
g\mu_{B}\sum_{i}e^{-i\bm{k}\cdot\bm{r}_{i}}\bm{S}_{i}$ is the magnetization
operator at wave-vector $\bm{k}$ and $N$ is the total number of spins. The
observed intensity is given by

$\displaystyle I(\bm{k},\omega)$ $\displaystyle\propto
F(k)^{2}\sum_{\mu\nu}\left(\delta_{\mu\nu}-\hat{k}_{\mu}\hat{k}_{\nu}\right)\mathcal{S}_{\mu\nu}(\bm{k},\omega)$
(A7)

where $F(k)$ is the Ir4+ magnetic form factor.

Within linear spin-wave theory the dynamical structure factor can be expressed
in terms of the transverse-transverse part of the spin-spin correlation
function. This remains true in our non-linear spin-wave theory (up to an
overall intensity renormalization) and $I(\bm{k},\omega)$ can be computed
using the self-consistently determined $H_{2}+\delta H_{2}$ as in LSWT. For
more details, including the averaging needed to emulate the binning of the
experimental data, we refer the reader to the Supplemental Material [56].

###  A.7 Classical spin dynamics simulations

In the paramagnetic phase we consider the FCC HK model with $J=0.74$ meV and
$K=0.15$ meV in the classical limit where the spins are unit length vectors
$|\bm{S}_{i}|^{2}=1$. The spins obey conventional Landau-Lifshitz dynamics

$\frac{d\bm{S}_{i}}{dt}=-\bm{S}_{i}\times\frac{\partial
H}{\partial\bm{S}_{i}},$ (A8)

where $-\partial H/\partial\bm{S}_{i}\equiv\bm{B}_{i}$ is the (local) exchange
field due to the neighboring spins. The initial conditions $\bm{S}_{i}(0)$ are
drawn from a thermal distribution at temperature $T$ using Monte Carlo
sampling. Once a sample of trajectories $\bm{S}_{i}(t)$ are obtained, the part
of the dynamical structure factor relevant for inelastic neutron scattering is
given by

$\mathcal{S}_{\rm
cl}(\bm{k},\omega)=\sum_{\mu\nu}\left(\delta_{\mu\nu}-\hat{k}_{\mu}\hat{k}_{\nu}\right)\braket{\bar{S}^{\mu}_{\bm{k}}(\omega)S^{\nu}_{\bm{k}}(\omega)},$
(A9)

where $\bm{S}_{\bm{k}}(\omega)$ is the Fourier transform of $\bm{S}_{i}(t)$ in
both space and time. To account for K2IrCl6 being $S=1/2$ we rescale
frequencies by a factor of $S$, ensuring that in the low temperature limit the
classical and quantum spin-wave frequencies agree. To partially account for
quantum effects we also multiply by an energy dependent correction factor [64,
53]

$F_{\rm qu}(\omega)\equiv\beta\omega\left(1+n_{B}(\omega)\right),$ (A10)

where $n_{B}(\omega)=1/(e^{\beta\omega}-1)$ is the Bose distribution.

Thermal samples are generated using Monte Carlo with single-site heat-bath
[40] and over-relaxation [12] updates, annealing down from high temperature
$O(10J)$ to the temperature of interest and then thermalizing for the same
number of sweeps. For Fig. 1 in the paramagnetic phase only a small number of
sweeps, typically $O(10^{3})$, are necessary to reach equilibrium and the
system size used was $N=L^{3}$ with $L=32$. For each initial state we solve
the coupled non-linear ordinary differential equations using an adaptive
fourth-order Runge-Kutta method. At the temperatures of interest $O(10^{2})$
samples are sufficient to reach convergence in both energy resolved and energy
averaged quantities. The adaptive time stepping was performed with (absolute
and relative) error tolerances of $10^{-8}$.

Simulations used to establish the dependence of the ordering type and ordering
temperature on $J_{2}/J$ used a larger number of sweeps, $O(10^{5})$, on
smaller system sizes with a conventional cell $N=4L^{3}$ where $L=8$, $10$,
$12$, as well as parallel tempering to aid in reaching equilibrium. The phase
transition temperature was inferred from the location of a sharp maximum in
the heat capacity and the onset of a type III order parameter.

## References

  * [1] A. A. Aczel, A. M. Cook, T. J. Williams, S. Calder, A. D. Christianson, G.-X. Cao, D. Mandrus, Y. Kim, and A. Paramekanti (2016-06) Highly anisotropic exchange interactions of jeff = $\frac{1}{2}$ iridium moments on the fcc lattice in La2BIrO6 ( B = Mg , Zn ).  Phys. Rev. B 93 (21),  pp. 214426.  External Links: Document Cited by: §A.4, §II, Pulling order back from the brink of disorder: Observation of a nodal line spin liquid and fluctuation stabilized order in K2IrCl6.
  * [2] P. W. Anderson (1950-08) Generalizations of the weiss molecular field theory of antiferromagnetism.  Phys. Rev. 79 (4),  pp. 705–710.  External Links: Document Cited by: Pulling order back from the brink of disorder: Observation of a nodal line spin liquid and fluctuation stabilized order in K2IrCl6, Pulling order back from the brink of disorder: Observation of a nodal line spin liquid and fluctuation stabilized order in K2IrCl6.
  * [3] O. Arnold, J.C. Bilheux, J.M. Borreguero, A. Buts, S.I. Campbell, L. Chapon, M. Doucet, N. Draper, R. Ferraz Leal, M.A. Gigg, V.E. Lynch, A. Markvardsen, D.J. Mikkelson, R.L. Mikkelson, R. Miller, K. Palmen, P. Parker, G. Passos, T.G. Perring, P.F. Peterson, S. Ren, M.A. Reuter, A.T. Savici, J.W. Taylor, R.J. Taylor, R. Tolchenov, W. Zhou, and J. Zikovsky (2014) Mantid—Data analysis and visualization package for neutron scattering and $\mu$SR experiments.  Nucl. Instrum. Methods Phys. Res. Sect. Accel. Spectrometers Detect. Assoc. Equip. 764,  pp. 156–166.  External Links: Document Cited by: §A.2.
  * [4] C. A. Bailey and P. L. Smith (1959-05) Specific Heats of Ammonium, Potassium, and Sodium Chloroiridates.  Phys. Rev. 114 (4),  pp. 1010–1010.  External Links: Document Cited by: Pulling order back from the brink of disorder: Observation of a nodal line spin liquid and fluctuation stabilized order in K2IrCl6.
  * [5] P. Balla, Y. Iqbal, and K. Penc (2020-11) Degenerate manifolds, helimagnets, and multi-Q chiral phases in the classical Heisenberg antiferromagnet on the face-centered-cubic lattice.  Phys. Rev. Res. 2 (4),  pp. 043278.  External Links: Document Cited by: Pulling order back from the brink of disorder: Observation of a nodal line spin liquid and fluctuation stabilized order in K2IrCl6.
  * [6] D. Bergman, J. Alicea, E. Gull, S. Trebst, and L. Balents (2007-07) Order-by-disorder and spiral spin-liquid in frustrated diamond-lattice antiferromagnets.  Nat. Phys. 3 (7),  pp. 487–491.  External Links: Document Cited by: Pulling order back from the brink of disorder: Observation of a nodal line spin liquid and fluctuation stabilized order in K2IrCl6.
  * [7] L. Bhaskaran, A. N. Ponomaryov, J. Wosnitza, N. Khan, A. A. Tsirlin, M. E. Zhitomirsky, and S. A. Zvyagin (2021-11) Antiferromagnetic resonance in the cubic iridium hexahalides (NH4)2IrCl6 and K2IrCl6.  Phys. Rev. B 104 (18),  pp. 184404.  External Links: Document Cited by: §II, Pulling order back from the brink of disorder: Observation of a nodal line spin liquid and fluctuation stabilized order in K2IrCl6.
  * [8] Th. Brueckel, B. Dorner, A. G. Gukasov, V. P. Plakhty, W. Prandl, E. F. Shender, and O. P. Smirnow (1988-12) Dynamical interaction of antiferromagnetic subsystems: a neutron scattering study of the spinwave spectrum of the garnet Fe2Ca3(GeO4)3.  Z. Für Phys. B Condens. Matter 72 (4),  pp. 477–485.  External Links: Document Cited by: §IV.
  * [9] J. D. M. Champion, M. J. Harris, P. C. W. Holdsworth, A. S. Wills, G. Balakrishnan, S. T. Bramwell, E. Č Čiž žmár, T. Fennell, J. S. Gardner, J. Lago, D. F. McMorrow, M. Orendáč č, A. Orendáč čová, D. McK. Paul, R. I. Smith, M. T. F. Telling, and A. Wildes (2003-07) Er2Ti2O7: Evidence of quantum order by disorder in a frustrated antiferromagnet.  Phys. Rev. B 68 (2),  pp. 020401.  External Links: Document Cited by: §IV.
  * [10] A. M. Cook, S. Matern, C. Hickey, A. A. Aczel, and A. Paramekanti (2015-07) Spin-orbit coupled jeff=$\frac{1}{2}$ iridium moments on the geometrically frustrated fcc lattice.  Phys. Rev. B 92 (2),  pp. 020417.  External Links: Document Cited by: §A.4.
  * [11] Cooke (1959-02) Exchange interactions in antiferromagnetic salts of iridium II. Magnetic susceptibility measurements.  Proc. R. Soc. Lond. Ser. Math. Phys. Sci. 250 (1260),  pp. 97–109.  External Links: Document Cited by: Pulling order back from the brink of disorder: Observation of a nodal line spin liquid and fluctuation stabilized order in K2IrCl6.
  * [12] M. Creutz (1987-07) Overrelaxation and monte carlo simulation.  Phys. Rev. D 36,  pp. 515–519.  External Links: Document, Link Cited by: §A.7.
  * [13] A. Danielian (1961-06) Ground state of an ising face-centered cubic lattice.  Phys. Rev. Lett. 6 (12),  pp. 670–671.  External Links: Document Cited by: Pulling order back from the brink of disorder: Observation of a nodal line spin liquid and fluctuation stabilized order in K2IrCl6.
  * [14] M. Elliot, P. A. McClarty, D. Prabhakaran, R. D. Johnson, H. C. Walker, P. Manuel, and R. Coldea (2021-12) Order-by-Disorder from Bond-Dependent Exchange and Intensity Signature of Nodal Quasiparticles in a Honeycomb Cobaltate.  Nat. Commun. 12 (1),  pp. 3936.  External Links: Document Cited by: §IV.
  * [15] S. Gao, K. Guratinder, U. Stuhr, J. S. White, M. Mansson, B. Roessli, T. Fennell, V. Tsurkan, A. Loidl, M. Ciomaga Hatnean, G. Balakrishnan, S. Raymond, L. Chapon, V. O. Garlea, A. T. Savici, A. Cervellino, A. Bombardi, D. Chernyshov, Ch. Rüegg, J. T. Haraldsen, and O. Zaharko (2018-04) Manifolds of magnetic ordered states and excitations in the almost Heisenberg pyrochlore antiferromagnet MgCr2O4.  Phys. Rev. B 97 (13),  pp. 134430.  External Links: Document Cited by: §IV.
  * [16] S. Gao, H. D. Rosales, F. A. Gómez Albarracín, V. Tsurkan, G. Kaur, T. Fennell, P. Steffens, M. Boehm, P. Čermák, A. Schneidewind, E. Ressouche, D. C. Cabra, C. Rüegg, and O. Zaharko (2020) Fractional antiferromagnetic skyrmion lattice induced by anisotropic couplings.  Nature 586 (7827),  pp. 37–41.  External Links: Document Cited by: §IV.
  * [17] S. Gao, O. Zaharko, V. Tsurkan, Y. Su, J. S. White, G. S. Tucker, B. Roessli, F. Bourdarot, R. Sibille, D. Chernyshov, T. Fennell, A. Loidl, and C. Rüegg (2017-02) Spiral spin-liquid and the emergence of a vortex-like state in MnSc2S4.  Nat. Phys. 13 (2),  pp. 157–161.  External Links: Document Cited by: §IV, Pulling order back from the brink of disorder: Observation of a nodal line spin liquid and fluctuation stabilized order in K2IrCl6.
  * [18] J. H. E. Griffiths, J. Owen, J. G. Park, M. F. Partridge, and B. Bleaney (1959) Exchange interactions in antiferromagnetic salts of iridium. I. Paramagnetic resonance experiments.  Proc. R. Soc. Lond. Ser. Math. Phys. Sci. 250 (1260),  pp. 84–96.  External Links: Document Cited by: Pulling order back from the brink of disorder: Observation of a nodal line spin liquid and fluctuation stabilized order in K2IrCl6.
  * [19] M. V. Gvozdikova and M. E. Zhitomirsky (2005-03) A Monte Carlo study of the first-order transition in a Heisenberg FCC antiferromagnet.  J. Exp. Theor. Phys. Lett. 81 (5),  pp. 236–240.  External Links: Document Cited by: §III, Pulling order back from the brink of disorder: Observation of a nodal line spin liquid and fluctuation stabilized order in K2IrCl6.
  * [20] D. T. Haar and M. E. Lines (1962) A molecular-field theory of anisotropic ferromagnetica.  Phil. Trans. R. Soc. A 254,  pp. 521–555.  External Links: Document Cited by: Pulling order back from the brink of disorder: Observation of a nodal line spin liquid and fluctuation stabilized order in K2IrCl6.
  * [21] D. T. Haar and M. E. Lines (1962) A spin-wave theory of anisotropic antiferromagnetica.  Phil. Trans. R. Soc. A 255 (1049),  pp. 1–30.  External Links: 73109, Document Cited by: Pulling order back from the brink of disorder: Observation of a nodal line spin liquid and fluctuation stabilized order in K2IrCl6.
  * [22] Harris (1965-12) The magnetic properties of three exchange coupled iridium ions.  Proc. R. Soc. Lond. Ser. Math. Phys. Sci. 289 (1416),  pp. 122–139.  External Links: Document Cited by: Pulling order back from the brink of disorder: Observation of a nodal line spin liquid and fluctuation stabilized order in K2IrCl6.
  * [23] C. L. Henley (1987-04) Ordering by disorder: Ground-state selection in fcc vector antiferromagnets.  J. Appl. Phys. 61 (8),  pp. 3962–3964.  External Links: Document Cited by: Pulling order back from the brink of disorder: Observation of a nodal line spin liquid and fluctuation stabilized order in K2IrCl6.
  * [24] C. L. Henley (1989-04) Ordering due to disorder in a frustrated vector antiferromagnet.  Phys. Rev. Lett. 62 (17),  pp. 2056–2059.  External Links: Document Cited by: Pulling order back from the brink of disorder: Observation of a nodal line spin liquid and fluctuation stabilized order in K2IrCl6.
  * [25] C. L. Henley (2010) The “coulomb phase” in frustrated systems.  Annu. Rev. Condens. Matter Phys. 1 (1),  pp. 179–210.  External Links: Document Cited by: Pulling order back from the brink of disorder: Observation of a nodal line spin liquid and fluctuation stabilized order in K2IrCl6.
  * [26] T. Holstein and H. Primakoff (1940-12) Field dependence of the intrinsic domain magnetization of a ferromagnet.  Phys. Rev. 58,  pp. 1098–1113.  External Links: Document, Link Cited by: §A.5.
  * [27] M. T. Hutchings and C. G. Windsor (1967-08) The magnetic structure of K2IrCl6.  Proc. Phys. Soc. 91 (4),  pp. 928–932.  External Links: Document Cited by: Pulling order back from the brink of disorder: Observation of a nodal line spin liquid and fluctuation stabilized order in K2IrCl6.
  * [28] G. Jackeli and G. Khaliullin (2009-01) Mott insulators in the strong spin-orbit coupling limit: From Heisenberg to a quantum compass and Kitaev models.  Phys. Rev. Lett. 102 (1),  pp. 017205.  External Links: Document Cited by: Pulling order back from the brink of disorder: Observation of a nodal line spin liquid and fluctuation stabilized order in K2IrCl6.
  * [29] Judd (1959-02) Exchange interactions in antiferromagnetic salts of iridium III. Theory of the interaction.  Proc. R. Soc. Lond. Ser. Math. Phys. Sci. 250 (1260),  pp. 110–120.  External Links: Document Cited by: §A.4, Pulling order back from the brink of disorder: Observation of a nodal line spin liquid and fluctuation stabilized order in K2IrCl6.
  * [30] N. Khan, D. Prishchenko, Y. Skourski, V. G. Mazurenko, and A. A. Tsirlin (2019-04) Cubic symmetry and magnetic frustration on the fcc spin lattice in K2IrCl6.  Phys. Rev. B 99 (14),  pp. 144425.  External Links: Document Cited by: §II, §II, Pulling order back from the brink of disorder: Observation of a nodal line spin liquid and fluctuation stabilized order in K2IrCl6.
  * [31] A. Krimmel, M. Mücksch, V. Tsurkan, M. M. Koza, H. Mutka, C. Ritter, D. V. Sheptyakov, S. Horn, and A. Loidl (2006-01) Magnetic ordering and spin excitations in the frustrated magnet MnSc2S4.  Phys. Rev. B 73 (1),  pp. 014413.  External Links: Document Cited by: §IV.
  * [32] S.-H. Lee, C. Broholm, M. F. Collins, L. Heller, A. P. Ramirez, Ch. Kloc, E. Bucher, R. W. Erwin, and N. Lacevic (1997-10) Less than 50% sublattice polarization in an insulating S=3/2 kagomé antiferromagnet at T$\approx$0.  Phys. Rev. B 56 (13),  pp. 8091–8097.  External Links: Document Cited by: §IV.
  * [33] S. Lee, B. H. Kim, M. Seong, and K. Choi (2022-05) Noncubic local distortions and spin-orbit excitons in K2IrCl6.  Phys. Rev. B 105 (18),  pp. 184433.  External Links: Document Cited by: §II.
  * [34] M. E. Lines (1963-01) Antiferromagnetism in the face-centred cubic lattice by a spin-wave method.  Proc. R. Soc. Lond. A 271 (1344),  pp. 105–119.  External Links: Document Cited by: Pulling order back from the brink of disorder: Observation of a nodal line spin liquid and fluctuation stabilized order in K2IrCl6.
  * [35] H. Liu and G. Khaliullin (2019-02) Pseudo-Jahn-Teller Effect and Magnetoelastic Coupling in Spin-Orbit Mott Insulators.  Phys. Rev. Lett. 122 (5),  pp. 057203.  External Links: Document Cited by: §II.
  * [36] J. W. Lynn, G. Shirane, and M. Blume (1976-07) Covalency effects in the magnetic form factor of Ir in K2IrCl6.  Phys. Rev. Lett. 37 (3),  pp. 154–157.  External Links: Document Cited by: Pulling order back from the brink of disorder: Observation of a nodal line spin liquid and fluctuation stabilized order in K2IrCl6.
  * [37] J.W. Lynn, H.H. Patterson, G. Shirane, and R.G. Wheeler (1978-09) Soft rotary mode and structural phase transitions in K2ReCl6.  Solid State Commun. 27 (9),  pp. 859–862.  External Links: Document Cited by: §II.
  * [38] M. Matsuura, Y. Endoh, H. Hiraka, K. Yamada, A. S. Mishchenko, N. Nagaosa, and I. V. Solovyev (2003-09) Classical and quantum spin dynamics in the fcc antiferromagnet NiS2 with frustration.  Phys. Rev. B 68 (9),  pp. 094409.  External Links: Document Cited by: Pulling order back from the brink of disorder: Observation of a nodal line spin liquid and fluctuation stabilized order in K2IrCl6.
  * [39] V.J. Minkiewicz, G. Shirane, B.C. Frazer, R.G. Wheeler, and P.B. Dorain (1968) Neutron diffraction study of magnetic ordering in K2IrCl6, K2ReBr6 and K2ReCl6.  J. Phys. Chem. Solids 29 (6),  pp. 881–884.  External Links: Document Cited by: Pulling order back from the brink of disorder: Observation of a nodal line spin liquid and fluctuation stabilized order in K2IrCl6.
  * [40] Y. Miyatake, M. Yamamoto, J. Kim, M. Toyonaga, and O. Nagai (1986) On the implementation of the’heat bath’algorithms for Monte Carlo simulations of classical Heisenberg spin systems.  J. Phys. C Solid State Phys. 19 (14),  pp. 2539.  External Links: Document Cited by: §A.7.
  * [41] R. Moessner and J. T. Chalker (1998-03) Properties of a classical spin liquid: The heisenberg pyrochlore antiferromagnet.  Phys. Rev. Lett. 80 (13),  pp. 2929–2932.  External Links: Document Cited by: Pulling order back from the brink of disorder: Observation of a nodal line spin liquid and fluctuation stabilized order in K2IrCl6.
  * [42] D. Moses, M. Sutton, R. L. Armstrong, and P. P. M. Meincke (1979) Specific Heat Measurements at Low Temperatures in K2MCl6 Antifluorite Crystals.  J. Low Temp. Phys. 36,  pp. 587–597.  External Links: Document Cited by: Pulling order back from the brink of disorder: Observation of a nodal line spin liquid and fluctuation stabilized order in K2IrCl6.
  * [43] J. Owen (1959) Exchange interactions in antiferromagnetic salts of iridium.  J. Phys. Radium 20 (2-3),  pp. 138–140.  External Links: Document Cited by: Pulling order back from the brink of disorder: Observation of a nodal line spin liquid and fluctuation stabilized order in K2IrCl6.
  * [44] 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 (2019) Continuum of quantum fluctuations in a three dimensional S = 1 Heisenberg magnet.  Nat. Phys. 15,  pp. 54–59.  External Links: Document Cited by: Pulling order back from the brink of disorder: Observation of a nodal line spin liquid and fluctuation stabilized order in K2IrCl6.
  * [45] E. Rastelli and A. Tassi (1987-05) Order produced by quantum disorder in the Heisenberg rhombohedral antiferromagnet.  J. Phys. C Solid State Phys. 20 (15),  pp. L303.  External Links: Document Cited by: Pulling order back from the brink of disorder: Observation of a nodal line spin liquid and fluctuation stabilized order in K2IrCl6.
  * [46] J. G. Rau, P. A. McClarty, and R. Moessner (2018-12) Pseudo-Goldstone Gaps and Order-by-Quantum Disorder in Frustrated Magnets.  Phys. Rev. Lett. 121 (23),  pp. 237201.  External Links: Document Cited by: §II, §II.
  * [47] D. Reig-i-Plessis, T. A. Johnson, K. Lu, Q. Chen, J. P. C. Ruff, M. H. Upton, T. J. Williams, S. Calder, H. D. Zhou, J. P. Clancy, A. A. Aczel, and G. J. MacDougall (2020-12) Structural, electronic, and magnetic properties of nearly ideal Jeff=$\frac{1}{2}$ iridium halides.  Phys. Rev. Mater. 4 (12),  pp. 124407.  External Links: Document Cited by: §II, Pulling order back from the brink of disorder: Observation of a nodal line spin liquid and fluctuation stabilized order in K2IrCl6.
  * [48] J. N. Reimers, A. J. Berlinsky, and A.-C. Shi (1991-01) Mean-field approach to magnetic ordering in highly frustrated pyrochlores.  Phys. Rev. B 43 (1),  pp. 865–878.  External Links: Document Cited by: Pulling order back from the brink of disorder: Observation of a nodal line spin liquid and fluctuation stabilized order in K2IrCl6.
  * [49] J. A. Rodriguez, D. M. Adler, P. C. Brand, C. Broholm, J. C. Cook, C. Brocker, R. Hammond, Z. Huang, P. Hundertmark, J. W. Lynn, N. C. Maliszewskyj, J. Moyer, J. Orndorff, D. Pierce, T. D. Pike, G. Scharfstein, S. A. Smee, and R. Vilaseca (2008-01) MACS—a new high intensity cold neutron spectrometer at NIST.  Meas. Sci. Technol. 19 (3),  pp. 034023.  External Links: Document Cited by: §A.2.
  * [50] K. A. Ross, Y. Qiu, J. R. D. Copley, H. A. Dabkowska, and B. D. Gaulin (2014-02) Order by disorder spin wave gap in the XY pyrochlore magnet Er2Ti2O7.  Phys. Rev. Lett. 112 (5),  pp. 057201.  External Links: Document Cited by: §IV.
  * [51] J. P. C. Ruff, J. P. Clancy, A. Bourque, M. A. White, M. Ramazanoglu, J. S. Gardner, Y. Qiu, J. R. D. Copley, M. B. Johnson, H. A. Dabkowska, and B. D. Gaulin (2008-10) Spin waves and quantum criticality in the frustrated XY pyrochlore antiferromagnet Er2Ti2O7.  Phys. Rev. Lett. 101 (14),  pp. 147205.  External Links: Document Cited by: §IV.
  * [52] L. Savary, K. A. Ross, B. D. Gaulin, J. P. C. Ruff, and L. Balents (2012-10) Order by quantum disorder in Er2Ti2O7.  Phys. Rev. Lett. 109 (16),  pp. 167201.  External Links: Document Cited by: §IV.
  * [53] A. Scheie, O. Benton, M. Taillefumier, L. D. C. Jaubert, G. Sala, N. Jalarvo, S. M. Koohpayeh, and N. Shannon (2022-11) Dynamical scaling as a signature of multiple phase competition in Yb2Ti2O7.  Phys. Rev. Lett. 129 (21),  pp. 217202.  External Links: Document Cited by: §A.7.
  * [54] R. Schick, O. Götze, T. Ziman, R. Zinke, J. Richter, and M. E. Zhitomirsky (2022-09) Ground-state selection by magnon interactions in a fcc antiferromagnet.  Phys. Rev. B 106 (9),  pp. 094431.  External Links: Document Cited by: §II, §II.
  * [55] R. Schick, T. Ziman, and M. E. Zhitomirsky (2020-12) Quantum versus thermal fluctuations in the fcc antiferromagnet: Alternative routes to order by disorder.  Phys. Rev. B 102 (22),  pp. 220405.  External Links: Document Cited by: §III, Pulling order back from the brink of disorder: Observation of a nodal line spin liquid and fluctuation stabilized order in K2IrCl6.
  * [56] See supplemental material at [url will be inserted by publisher] for further details of sample characterization, additional supporting neutron scattering data, details of classical spin dynamics and non-linear spin wave theory, and considerations of the effects of a tetragonal distortion..  Cited by: §A.1, §A.2, §A.6, Figure 1, §II, §II, §II, §III, §III.
  * [57] R.A. Tahir-Kheli, H.B. Callen, and H. Jarrett (1966) Magnetic ordering in cubic crystals with first and second neighbor exchange.  J. Phys. Chem. Solids 27 (1),  pp. 23–32.  External Links: Document Cited by: Pulling order back from the brink of disorder: Observation of a nodal line spin liquid and fluctuation stabilized order in K2IrCl6.
  * [58] J. Villain (1979-03) Insulating spin glasses.  Z. Für Phys. B Condens. Matter 33 (1),  pp. 31–42.  External Links: Document Cited by: Pulling order back from the brink of disorder: Observation of a nodal line spin liquid and fluctuation stabilized order in K2IrCl6.
  * [59] Villain, J., Bidaux, R., Carton, J.-P., and Conte, R. (1980) Order as an effect of disorder.  J. Phys. France 41 (11),  pp. 1263–1272.  External Links: Document Cited by: Pulling order back from the brink of disorder: Observation of a nodal line spin liquid and fluctuation stabilized order in K2IrCl6.
  * [60] Q. Wang, A. De La Torre, W. Tian, A. A. Aczel, J. C. Neuefeind, M. M. Everett, Y. Feng, M. H. Upton, J. Kim, P. J. Ryan, J. Kim, and K. W. Plumb (2024) In preparation.  Cited by: §II, §II.
  * [61] H. W. Willemsen, R. L. Armstrong, and P. P. M. Meincke (1977) Thermal Expansion Measurements Near the Antiferromagnetic Phase Transitions in K2ReCI6 and K2IrCI6.  J. Low Temp. Phys. 26,  pp. 11.  External Links: Document Cited by: Pulling order back from the brink of disorder: Observation of a nodal line spin liquid and fluctuation stabilized order in K2IrCl6.
  * [62] T. Yildirim, A. B. Harris, and E. F. Shender (1998-08) Frustration and quantum fluctuations in Heisenberg fcc antiferromagnets.  Phys. Rev. B 58 (6),  pp. 3144–3159.  External Links: Document Cited by: Pulling order back from the brink of disorder: Observation of a nodal line spin liquid and fluctuation stabilized order in K2IrCl6.
  * [63] B. Yuan, I. Khait, G. Shu, F. C. Chou, M. B. Stone, J. P. Clancy, A. Paramekanti, and Y. Kim (2020-03) Dirac magnons in a honeycomb lattice quantum XY magnet CoTiO3.  Phys. Rev. X 10 (1),  pp. 011062.  External Links: Document Cited by: §IV.
  * [64] S. Zhang, H. J. Changlani, K. W. Plumb, O. Tchernyshyov, and R. Moessner (2019-04) Phys. Rev. Lett. 122 (16),  pp. 167203.  External Links: Document Cited by: §A.7, §IV.
  * [65] M. E. Zhitomirsky, M. V. Gvozdikova, P. C. W. Holdsworth, and R. Moessner (2012-08) Quantum order by disorder and accidental soft mode in Er2Ti2O7.  Phys. Rev. Lett. 109 (7),  pp. 077204.  External Links: Document Cited by: Pulling order back from the brink of disorder: Observation of a nodal line spin liquid and fluctuation stabilized order in K2IrCl6.

Supplemental Material for “Pulling order back from the brink of disorder:
Observation of a nodal line spin liquid and fluctuation stabilized order in
K2IrCl6 ”

###### Contents

  1. I Nodal line spin liquid
  2. II Fluctuation induced spin wave gap
  3. III Competing type I and III magnetic order
  4. IV Discussion
  5. A Methods
    1. A.1 Sample Synthesis and Characterization
    2. A.2 Neutron Scattering
    3. A.3 X-ray Scattering
    4. A.4 Heisenberg-Kitaev model on the FCC lattice
    5. A.5 Self-consistent non-linear spin-wave theory
    6. A.6 Calculation of inelastic neutron scattering intensity
    7. A.7 Classical spin dynamics simulations
  6. A Experiment
    1. A.1 Sample characterization
    2. A.2 Details of Absorption Corrections for Inelastic Neutron Scattering Data
    3. A.3 Additional Inelastic Neutron Scattering Data
    4. A.4 Additional Polarized Neutron Scattering Data
    5. A.5 Survey of order by disorder gaps in various materials
  7. B Theory
    1. B.1 Model
    2. B.2 Classical Ground States
      1. B.2.1 Luttinger-Tisza
      2. B.2.2 Ground State Manifold
      3. B.2.3 Effect of Second-neighbor Exchange
    3. B.3 Linear Spin Wave Theory
      1. B.3.1 Spin Waves at $O(1/S)$
      2. B.3.2 Order by Quantum Disorder
      3. B.3.3 Dynamical Structure Factor
    4. B.4 Non-Linear Spin Wave Theory
      1. B.4.1 Self-consistent Spin Waves
      2. B.4.2 Spin Waves at $O(1/S^{2})$
      3. B.4.3 Corrections to the Dynamical Structure Factor
    5. B.5 Effect of Tetragonal Distortion
      1. B.5.1 Effect on Magnon Gap in Type-III State
      2. B.5.2 Best Fit to Tetragonal Model
    6. B.6 Estimating the Magneto-elastic Coupling
    7. B.7 Dynamics in Paramagnetic Phase
    8. B.8 Effect of $J_{2}$ on Néel temperature
    9. B.9 Constraints on Size of $\Gamma$ Exchange
    10. B.10 Details of Comparison Between Theoretical and Experimental Dynamical Structure Factor
      1. B.10.1 Energy Resolution
      2. B.10.2 Wave-vector Averaging
      3. B.10.3 Iridium Form Factor
      4. B.10.4 Interpolation of Theoretical Results

##  Appendix A Experiment

###  A.1 Sample characterization

Powder X-ray diffraction measurements were carried out on the 11-BM
diffractometer at the Advanced Photon Source using an X-ray wavelength of
0.45893 Å. Measurements were conducted at at $T=90$ K and 300 K. Reitveld
refinement of powder data was carried out using the FullProf software [42].

Fig. S1 shows the diffraction data and Rietveld refinement at 90 K (a) and 300
K (b). K2IrCl6 maintains a cubic crystal structure at all measurement
temperatures with no observable symmetry lowering down to 90 K. Table S1
summarizes the resulting crystallographic parameters and atomic positions. All
refined parameters are consistent with published reports [41, 25].

Table S1: Crystallographic parameters and details of the structure refinement for K2IrCl6.  $T$ | 90 K | 300 K
---|---|---
Space group | $Fm\overline{3}m$ | $Fm\overline{3}m$
$a=b=c$ (Å) | 9.6973(2) | 9.7777(3)
$\alpha=\beta=\gamma$ (∘) | 90∘ | 90∘
Atomic parameters |  |
K | $x/a=0.25,y/b=0.25,z/c=0.25$ | $x/a=0.25,y/b=0.25,z/c=0.25$
| $B_{\rm iso}=1.119(39)$ | $B_{\rm iso}=2.779(68)$
| ${\rm Occ}=1.003$ | ${\rm Occ}=0.994$
Ir | $x/a=0,y/b=0,z/c=0$ | $x/a=0,y/b=0,z/c=0$
| $B_{\rm iso}=0.260(8)$ | $B_{\rm iso}=0.886(10)$
| ${\rm Occ}=1.002$ | ${\rm Occ}=1.001$
Cl | $x/a=0.2394(2),y/b=0,z/c=0$ | $x/a=0.2373(2),y/b=0,z/c=0$
| $B_{\rm iso}=0.840(28)$ | $B_{\rm iso}=2.203(47)$
| ${\rm Occ}=0.991$ | ${\rm Occ}=0.998$
Refinement |  |
$R_{\rm p}$ | 9.32 | 8.59
$R_{\rm exp}$ | 4.90 | 5.57
$\chi^{2}$ | 6.72 | 3.52

To complement the diffraction measurements and verify a local cubic structure
in K2IrCl6, $T=300$ K x-ray pair distribution function (PDF) measurements were
carried out at the Advanced Photon Source, 11-ID-B with an X-ray wavelength of
0.2116 Å. We used GSAS-II [50] and PDFgui [18] for data processing and
refinement. The measured PDF and refinement are shown in Fig. S2 with
positions for distances of Ir-Cl, Cl-Cl and Ir-Ir indicated out. The data
reveal a single nearest neighbor Ir-Cl distance as expected for a local cubic
environment. No deviation from the cubic structure is observed at room
temperature.

Heat capacity measurements were conducted using a Quantum Design PPMS in zero
field during a warming process. The phonon contribution was estimated using
the empirical model by N. Khan et al [25], which incorporates one Debye-type
and three Einstein-type terms. The Debye term accounts for 24% of the total
modes, characterized by a Debye temperature $\Theta_{\rm D}=91.1$ K. The three
Einstein terms have Einstein temperatures $\Theta_{\rm E}=459,~{}242,~{}136$
K, respectively.

Magnetization measurements were carried out using a Quantum Design PPMS
vibrating sample magnetometer (VSM). Fig. S3(a) shows the temperature
dependent magnetic susceptibility measured under a $[1,0,0]$ oriented $H=0.1$
T field. Data was collected on warming after zero-field cooling. A Curie-Weiss
temperature of $\theta_{\rm CW}=-32.98$ K and paramagnetic moment of 1.59
$\mu_{B}$ were extracted from a Curie-Weiss fit to the data over the
temperature range $100<T<280$ K. These parameters are in agreement with
published characterization [25]. Fig. S3(b) shows the field dependent
magnetization at 2 K, 6 K, 60 K with the field along $[1,0,0]$ and $[1,1,1]$
crystal directions.

Figure S1:  High-resolution synchrotron powder X-ray diffraction on K2IrCl6.
Results of a Rietveld refinement(black line) are displayed overtop of the
high-resolution powder x-ray diffraction (red circle) at (a) 90 K and (b) 300
K . The expected reflections are indexed with green ticks and the blue line is
the fit residual. Insets show the zoom-in of the high-$Q$ diffraction pattern.
No structural distortion or impurity phases are discernible. Figure S2:  X-ray
pair distribution function (PDF) for K2IrCl6 measured at 300 K. The measured
PDF (red circles) on top of the refinement (blue line) with the difference
(green line). The positions for distances of Ir-Cl, Cl-Cl and Ir-Ir are
pointed out. No deviation from the cubic structure is observed at room
temperature. Figure S3:  Magnetic properties of K2IrCl6. (a), Temperature
dependence of magnetic susceptibility (black curve) and inverse magnetic
susceptibility (blue curve) for K2IrCl6 measured under a magnetic field of
$H=0.1$ T. The red line is the Curie-Weiss fit to the inverse susceptibility
data at $T>100$ K. (b), Field dependence of the magnetic moment at 2 K (red),
6 K (blue) and 60 K (black) along $[1,0,0]$ (solid) and $[1,1,1]$ (dashed)
directions of the crystal lattice.

###  A.2 Details of Absorption Corrections for Inelastic Neutron Scattering
Data

Iridium has a large neutron absorption cross-section of 425 barn. The
absorption lengths of K2IrCl6 are 4.2 and 8.2 cm-1 at neutron energies 14.5
and 3.32 meV, respectively. Thus an absorption correction is necessary in
processing the neutron scattering data. The absorption correction for CNCS
data was conducted using the built-in Mantid algorithm [6]. The Mantid
algorithm relies on a numerical integration method, which effectively
calculates the attenuation factors arising from sample absorption based on its
given material properties. In this study, an approximate cylindrical crystal
array with a radius of 0.12 cm and a height of 3 cm was utilized for the
measurements. The cylindrical crystal array was discretized into small cubes,
each measuring $1\times 1\times 1$ mm in size. The cubes whose centers lie
within the sample make up the set of integration elements. Path lengths
through the sample were then computed for each center-point of the selected
integration elements. Finally, a numerical integration technique was employed
over the volume elements using these calculated path lengths to determine the
necessary corrections [6].

###  A.3 Additional Inelastic Neutron Scattering Data

Additional inelastic neutron scattering measurements carried out on the MACS
spectrometer with the K2IrCl6 crystal aligned in the $[h,h,l]$ scattering
plane are shown in Fig. S2. The observed spin wave dispersion and gap are
consistent with the CNCS data and non-linear spin wave modeling of the nearest
neighbor Heisenberg-Kitaev model described below and in the main text.

Figure S4: Additional inelastic neutron scattering data. Magnetic excitation
spectra measured in the $[h,h,l]$ scattering plane on MACS at 1.8 K.

###  A.4 Additional Polarized Neutron Scattering Data

Polarized neutron scattering measurement on the type-I magnetic peak [1,0,0]
are shown in Fig. S5. The spin-flip intensity of [1,0,0] reflection is
sensitive to both in and out-of-plane components, i.e. [0,0,1] direction of
the magnetization for guide fields $\bm{P}\parallel\bm{Q}$, and only in-plane
components for $\bm{P}\perp\bm{Q}$. Although the [1,0,0] peak is extremely
weak, we still observe a stronger intensity in the $\bm{P}\parallel\bm{Q}$
channel compared with the $\bm{P}\perp\bm{Q}$ channel. Based on a finite
flipping ratio of 15.2 which gives rise to the intensity in the
$\bm{P}\perp\bm{Q}$ channel, a total background of 0.126 counts/s is estimated
in the $\bm{P}\parallel\bm{Q}$ channel, and thus an intensity of $\sim$ 0.008
counts/s is left for the out-of-plane component of magnetization from type-I
magnetic order.

Figure S5: Polarized neutron diffraction measurement on the type-I peak at 1.5
K. The spin-flip channel was measured for neutron spins parallel to the
scattering vector ($\bm{P}\parallel\bm{Q}$) and perpendicular to the
scattering plane ($\bm{P}\perp\bm{Q}$). The difference between the two
channels corresponds to the out-of-plane component.

###  A.5 Survey of order by disorder gaps in various materials

In order to place our observations and the scale of quantum fluctuations in
K2IrCl6 in the broader context of order by quantum disorder candidate
materials, below we summarize magnon bandwidths and excitation gaps that have
been attributed to quantum order by disorder.

Table S2: Scale of quantum order by disorder gap relative to bandwidth in candidate materials.  Material | Spin |  $W$ (Bandwidth) |  $\Delta$ (Gap size) | $\Delta/W$ | Ref.
---|---|---|---|---|---
Ca3Fe2Ge3O12 | $S=5/2$ |  $\sim$ 1.95 meV | 136 $\mu$eV | 7.0% |  [10]
CoTiO3 | $S_{\rm eff}=1/2$ |  $\sim$ 12 meV | 1 meV | 8.3% |  [52, 17]
ErTi2O7 | $S_{\rm eff}=1/2$ |  $\sim$ 0.45 meV | 43 $\mu$eV | 9.6% |  [43, 12, 55, 44]
K2IrCl6 | $j=1/2$ | 2.5 meV | 0.7 meV | 30% | This study

##  Appendix B Theory

###  B.1 Model

We consider a minimal model for the (effective) $j=1/2$ doublets [23] of the
Ir4+ ions on an FCC lattice relevant for K2IrCl6. The degrees of freedom of
these doublets are described by an effective $S=1/2$ spin, $\bm{S}_{i}$, with
an isotropic $g$-factor. From electron spin resonance measurements [7] we
expect the $g$-factor to be $g\approx 1.8$, not far from $g=2$ expected for an
ideal $j_{\rm eff}=1/2$ doublet [37].

The cubic symmetry strongly constrains the symmetry-allowed exchange
interactions. At nearest-neighbor level, this includes nearest neighbor
Heisenberg exchange as well as bond-dependent Kitaev and $\Gamma$ exchanges
[24, 13, 2]

$\sum_{\braket{ij}_{\alpha\beta(\gamma)}}\left[J\bm{S}_{i}\cdot\bm{S}_{j}+KS^{\gamma}_{i}S^{\gamma}_{j}+(-1)^{\sigma^{\alpha\beta}_{ij}}\Gamma\left(S^{\alpha}_{i}S^{\beta}_{j}+S^{\alpha}_{i}S^{\beta}_{j}\right)\right]$
(S1)

We have divided the bonds of the lattice into three types: $x$, $y$ and $z$,
depending on whether they lie in the $yz$, $zx$ or $xy$ planes. The sign of
the $\Gamma$ term is determined by the bond direction
$\bm{d}_{ij}\equiv\bm{r}_{j}-\bm{r}_{i}$ as
$\sigma^{\alpha\beta}_{ij}\equiv{\rm sgn}(d^{\alpha}_{ij}d^{\beta}_{ij})$.
Second neighbor interactions along the cubic axes are restricted to only a
Heisenberg or Kitaev exchange [15]

$J_{2}\sum_{\braket{\braket{ij}}}\bm{S}_{i}\cdot\bm{S}_{j}+K_{2}\sum_{\braket{\braket{ij}}_{\gamma}}S^{\gamma}_{i}S^{\gamma}_{j}$
(S2)

where $\gamma$ is the bond direction. For simplicity, we will restrict our
discussion to only isotropic second neighbor exchange setting $K_{2}=0$.

###  B.2 Classical Ground States

####  B.2.1 Luttinger-Tisza

Figure S6: Phase diagram for $J$-$K$-$\Gamma$ model [Eq. (S1)] with $J>0$ from
the Luttinger-Tisza approach. For each phase we have indicated the minimal
wave-vector, with $X\equiv[1,0,0]$, $L\equiv[0.5,0.5,0.5]$,
$K\equiv[\tfrac{3}{2},\tfrac{3}{2},0]$, $\Gamma\equiv[0,0,0]$ and
$U\equiv[1,0.25,0.25]$. The notation $X|U$ denotes a minimum wave-vector along
the line connecting $X$ and $U$ (for example). The region with $|\Gamma|<K/2$,
$K>0$ does not have a unique minimum, but instead a degenerate line of minima
along $X|W$ where $W\equiv[0.5,1,0]$. Note that we have distinguished phases
only by wave-vector, not by moment direction.

To understand the qualitative features of the classical phase diagram, we
begin by revisiting the Luttinger-Tisza method [31, 29] that has been used in
previous works [26, 13, 25, 15]. We find a rich phase diagram including
multiple commensurate and incommensurate magnetic phases, as well as a broad
region with a degeneracy along one-dimensional wave-vector manifolds – “spiral
lines” along $[q,1,0]$, $[0,q,1]$, $[1,0,q]$ and equivalents. This is
illustrated in Fig. S6. The commensurate phase with wave-vector $X$ is a type
I phase, while the type III phase would correspond to a wave-vector of $W$. We
note that the incommensurate phases do not always satisfy the spin-length
constraint and thus do not always represent true classical ground states. [38]
In these regions direct simulation via parallel tempering Monte Carlo [22, 49,
34] using heat-bath updates [35] and iterative minimization [48] yield a
complex set of multi-$Q$ incommensurate spirals whose wave-vectors
qualitatively track the Luttinger-Tisza result.

####  B.2.2 Ground State Manifold

Figure S7: Illustration of a generic “stacked” ground state (along
$\bm{\hat{x}}$) of the nearest-neighbor model [Eq. (S1)] when $J>0$ and
$|\Gamma|<K/2$,. The moment direction would be parallel to stacking direction,
$\bm{\hat{x}}$, with filled and open circles denoting the sign.

In the region with $K>0$ and $|\Gamma|<K/2$ containing the “spiral lines” the
Luttinger-Tisza energy is $E_{0}=-(2J+K)S^{2}$, representing a lower-bound on
the true classical ground state energy of the system. One can find classical
ground states that saturate this bound and thus represent true ground states
of the model. We first construct a sub-extensive number of collinear states
with this energy. These can be understood as an arbitrary stacking of Néel
planes along the cubic directions, with two choices for the Néel state in each
plane. With three choice of stacking direction ($x$, $y$ or $z$) and two
choices per plane, we have $3\cdot 2^{2L}$ states where $L$ is the number of
cubic unit cells along the stacking direction. More explicitly we can define
three families of states

$\displaystyle\bm{S}^{x}_{i}$
$\displaystyle=S(-1)^{n_{1}+n_{2}}\sigma_{n_{2}+n_{3}}^{x}\bm{\hat{x}},$
$\displaystyle\bm{S}^{y}_{i}$
$\displaystyle=S(-1)^{n_{2}+n_{3}}\sigma_{n_{3}+n_{1}}^{y}\bm{\hat{y}},$
$\displaystyle\bm{S}^{z}_{i}$
$\displaystyle=S(-1)^{n_{3}+n_{1}}\sigma_{n_{1}+n_{2}}^{z}\bm{\hat{z}}$

where we have expressed the position
$\bm{r}_{i}=n_{1}\bm{a}_{1}+n_{2}\bm{a}_{2}+n_{3}\bm{a}_{3}$ in terms of the
primitive lattice vectors $\bm{a}_{1}$, $\bm{a}_{2}$ and $\bm{a}_{3}$. The
$\sigma^{\mu}_{n}$ are Ising variables equal to $\pm 1$ that encode the choice
Néel plane. For each of these states the inter-plane couplings cancel and only
the plane perpendicular to the stacking direction contributes to the classical
energy. For $|\Gamma|<K/2$ this is minimized by a Néel state (satisfying the
Heisenberg part) with moment parallel to the stacking direction (satisfying
the Kitaev part), with $\Gamma$ dropping out. This gives the required
$E=-(2J+K)S^{2}$. These “stacked” states are identical those that make up the
sub-extensive discrete degeneracy of the FCC Ising anti-ferromagnet [32].

This discrete degeneracy does not exhaust the ground state manifold. As these
stacked states are collinear and can be oriented along perpendicular axes,
linear combinations of these states are also ground states. More explicitly,
if we define

$\bm{S}_{i}\equiv
S\left[\alpha_{x}(-1)^{n_{1}+n_{2}}\sigma_{n_{2}+n_{3}}^{x}\bm{\hat{x}}+\alpha_{y}(-1)^{n_{2}+n_{3}}\sigma_{n_{3}+n_{1}}^{y}\bm{\hat{y}}+\alpha_{z}(-1)^{n_{3}+n_{1}}\sigma_{n_{1}+n_{2}}^{z}\bm{\hat{z}}\right]$

where $\sum_{\mu}\alpha_{\mu}^{2}=1$, then we can see that these are
normalized $|\bm{S}_{i}|^{2}=1$. The cross terms in the energy vanish for
distinct stackings, yielding the same $-(2J+K)S^{2}$ energy. We have thus
established a very large _continuous_ manifold of ground states: for each
discrete ground state (say) stacked along the $\bm{\hat{z}}$ direction, we can
“rotate” it into two arbitrary ground states stacked in the $\bm{\hat{x}}$ and
$\bm{\hat{y}}$ directions, producing a new (non-collinear) ground state.

This continuous manifold is a subset of the continuous manifold found in the
Heisenberg limit when $J>0$ and $K=\Gamma=0$ [5, 4]. In the Heisenberg limit,
this manifold is considerably larger. For example, the stacked Néel planes can
be mixed with states stacked along the same direction, but with orthogonal
moment directions.

####  B.2.3 Effect of Second-neighbor Exchange

\begin{overpic}[width=138.76157pt]{type-III-x-
crop.pdf}\put(42.0,-5.0){$X$}\end{overpic}\begin{overpic}[width=99.73074pt]{type-
III-y-
crop.pdf}\put(42.0,-7.0){$Y$}\end{overpic}\begin{overpic}[height=138.76157pt]{type-
III-z-crop.pdf}\put(42.0,-5.0){$Z$}\end{overpic}

Figure S8: Illustration of the sign structure of type-III anti-ferromagnetic
states favored by anti-ferromagnetic second-neighbor exchange, $J_{2}>0$.
Filled circles indicate moments along the direction of the stacked Néel
planes, open circles indicate moments against the stacking direction.

\begin{overpic}[width=86.72267pt]{type-I-100-crop.pdf}\put(42.0,-10.0){$[100]$}\end{overpic}\begin{overpic}[width=86.72267pt]{type-I-010-crop.pdf}\put(42.0,-10.0){$[010]$}\end{overpic}\begin{overpic}[width=86.72267pt]{type-I-001-crop.pdf}\put(42.0,-10.0){$[001]$}\end{overpic}

Figure S9: Illustration of the sign structure of the type-I anti-ferromagnetic
states favored by ferromagnetic second-neighbor exchange, $J_{2}<0$. Filled
circles indicate moments along the direction of the stacked Néel planes, open
circles indicate moments against the stacking direction.

Adding a finite second neighbor coupling immediately lifts most, but not all,
of the large accidental degeneracy present in the nearest-neighbor model with
$K>0$ and $|\Gamma|<K/2$.

Antiferromagnetic second neighbor exchange, $J_{2}>0$, favors a collinear
type-III anti-ferromagnet containing wave-vectors $(0.5,1,0)$, $(0,0.5,1)$,
$(1,0,0.5)$ or equivalents. This can be viewed as a period four stacking of
the form $(\pm,\pm,\mp,\mp)$ or equivalent of Néel planes, with the moment
direction parallel to the stacking direction. There are four choices for the
periodicity and three stacking directions and thus 12 distinct collinear
states. The sign structure of this state is illustrated in Fig. S8.

For $J_{2}<0$ it favors a collinear type-I anti-ferromagnet with wave-vectors
$(1,0,0)$, $(0,1,0)$ or $(0,0,1)$. These can be viewed as stacks of
_ferromagnetic_ planes with the moments oriented perpendicular to the stacking
direction. There are 12 distinct collinear type I states, corresponding to
three stacking directions, two choices of moment direction and two choices for
the alternation pattern. In terms of stacking of _Néel_ planes, these can be
viewed as being composed of the four stackings of the form $(\pm,\pm)$ or
$(\pm,\mp)$ along each of the three cubic directions with the moment parallel
to the stacking direction. The sign structure of this state is illustrated in
Fig. S9.

In both cases there is a remnant of the larger degeneracy of the nearest-
neighbor model: states stacked along different directions can be continuously
mixed. For $J_{2}>0$ one can still mix a type-III state with any other type-
III state with perpendicular moment, and for $J_{2}<0$ you can smoothly mix a
type-I state with another type-I state so along as their moments are
orthogonal. This degeneracy cannot be lifted by any exchange interactions
bilinear in the spins that respect the cubic symmetry of the crystal. For each
case, three orthogonal type I or type III states can be grouped into a order
parameter $\bm{m}$ that transforms as a vector under the cubic symmetries.
Since the only bilinear in $\bm{m}$ that respects cubic symmetry is
$\propto|\bm{m}|^{2}$ we see that an accidental $O(3)$ will remain regardless
of the addition of further (symmetry preserving) exchange interactions [44].

###  B.3 Linear Spin Wave Theory

We now consider a semi-classical expansion about the one of ground states
described above. The spin operators can expressed in terms of Holstein-
Primakoff bosons as

$\bm{S}_{\bm{r}\alpha}\equiv\sqrt{S}\left[\left(1-\frac{n_{\bm{r}\alpha}}{2S}\right)^{1/2}{a}^{\phantom{\dagger}}_{\bm{r}\alpha}\bm{\hat{e}}_{\alpha,-}+{a}^{\dagger}_{\bm{r}\alpha}\left(1-\frac{n_{\bm{r}\alpha}}{2S}\right)^{1/2}\bm{\hat{e}}_{\alpha,+}\right]+\left(S-n_{\bm{r}\alpha}\right)\bm{\hat{e}}_{\alpha,0},$
(S3)

where
$n_{\bm{r}\alpha}\equiv{a}^{\dagger}_{\bm{r}\alpha}{a}^{\phantom{\dagger}}_{\bm{r}\alpha}$
and $\bm{r},\alpha$ denotes the unit cell and sublattice index of the spin.
For the collinear type-I or type-III orders a four-sublattice unit cell is
sufficient, while for other stacked states or for non-collinear mixtures of
type-I and type-III larger unit cells are necessary. The vectors
$\bm{\hat{e}}_{\alpha,\pm}$, $\bm{\hat{e}}_{\alpha,0}$ define a local frame of
reference; in a more conventional Cartesian basis one defines
$\bm{\hat{e}}_{\alpha,\pm}\equiv(\bm{\hat{x}}_{\alpha}\pm
i\bm{\hat{y}}_{\alpha})/\sqrt{2}$ and
$\bm{\hat{e}}_{\alpha,0}\equiv\bm{\hat{z}}_{\alpha}$. It is useful to write
the exchange matrix in the frame aligned with these axes

$\mathcal{J}^{\mu\mu^{\prime}}_{\bm{\delta},\alpha\alpha^{\prime}}\equiv{\bm{\hat{e}}}^{\intercal}_{\alpha,\mu}{\bm{J}}_{\bm{\delta},\alpha\alpha^{\prime}}{\bm{\hat{e}}}^{\phantom{\dagger}}_{\alpha^{\prime},\mu^{\prime}}.$
(S4)

Expanding in powers of $1/S$ then yields a semi-classical expansion about the
ordered state defined by $\bm{\hat{e}}_{\alpha,0}$, typically chosen to be
along the classical ordering direction.

At order $O(S)$ in the Holstein-Primakoff operators we have

$\bm{S}_{\bm{r}\alpha}\approx\sqrt{S}\left[{a}^{\phantom{\dagger}}_{\bm{r}\alpha}\bm{\hat{e}}_{\alpha,-}+{a}^{\dagger}_{\bm{r}\alpha}\bm{\hat{e}}_{\alpha,+}\right]+\left(S-n_{\bm{r}\alpha}\right)\bm{\hat{e}}_{\alpha,0},$
(S5)

Inserting this into our spin Hamiltonian and keeping only terms to $O(S)$
yields

$H=NS(S+1)\epsilon_{\rm
cl}+\frac{1}{2}\sum_{\bm{k}}\left({[{\bm{a}}^{\dagger}_{\bm{k}}]}^{\intercal}\
{\bm{a}}^{\intercal}_{-\bm{k}}\right)\left(\begin{array}[]{cc}{\bm{A}}_{\bm{k}}&{\bm{B}}_{\bm{k}}\\\
\bar{{\bm{B}}}_{-\bm{k}}&\bar{{\bm{A}}}_{-\bm{k}}\end{array}\right)\left(\begin{array}[]{c}{\bm{a}}^{\phantom{\dagger}}_{\bm{k}}\\\
{\bm{a}}^{\dagger}_{-\bm{k}}\end{array}\right)+O(S^{1/2}),$ (S6)

where $N_{s}$ is the number of sublattices, $N$ is the total number of sites
and we have defined the Fourier transforms of the bosons as
$a_{\bm{k}\alpha}\equiv
N_{c}^{-1/2}\sum_{\bm{r}}e^{-i\bm{k}\cdot\bm{r}}a_{\bm{r}\alpha}$ where
$N=N_{c}N_{s}$. The classical energy per site is defined as

$\epsilon_{\rm
cl}\equiv\frac{1}{2N_{s}}\sum_{\alpha\alpha^{\prime}}\sum_{\bm{\delta}}\mathcal{J}^{00}_{\bm{\delta},\alpha\alpha^{\prime}},$
(S7)

and the matrices ${\bm{A}}_{\bm{k}}$ and ${\bm{B}}_{\bm{k}}$ are given by

$\displaystyle A^{\alpha\alpha^{\prime}}_{\bm{k}}$
$\displaystyle=S\left(\mathcal{J}^{+-}_{\bm{k},\alpha\alpha^{\prime}}-\delta_{\alpha\alpha^{\prime}}\sum_{\mu}\mathcal{J}^{00}_{\bm{0},\alpha\mu}\right),$
(S8a) $\displaystyle B^{\alpha\alpha^{\prime}}_{\bm{k}}$
$\displaystyle=S\mathcal{J}^{++}_{\bm{k},\alpha\alpha^{\prime}},$ (S8b)

where we have defined the Fourier transforms of the local exchange matrices as

$\mathcal{J}^{\mu\mu^{\prime}}_{\bm{k},\alpha\alpha^{\prime}}\equiv\sum_{\bm{\delta}}\mathcal{J}^{\mu\mu^{\prime}}_{\bm{\delta},\alpha\alpha^{\prime}}e^{i\bm{k}\cdot\bm{\delta}}.$
(S9)

The linear spin-wave Hamiltonian [Eq. (S6)] can be diagonalized by a
Bogoliubov transformation. To do this one diagonalizes the modified matrix [8]

$\displaystyle\left(\begin{array}[]{cc}{\bm{A}}_{\bm{k}}&{\bm{B}}_{\bm{k}}\\\
-\bar{{\bm{B}}}_{-\bm{k}}&-\bar{{\bm{A}}}_{-\bm{k}}\end{array}\right)\equiv{\bm{\sigma}}_{3}{\bm{M}}_{\bm{k}}$
(S12)

This yields pairs of eigenvectors $\bm{V}_{\bm{k}\alpha}$ and
$\bm{W}_{-\bm{k},\alpha}={\bm{\sigma}}_{1}\bar{\bm{V}}_{-\bm{k}\alpha}$ with
eigenvalues $+\epsilon_{\bm{k}\alpha}$ and $-\epsilon_{-\bm{k},\alpha}$. These
vectors can be normalized such that [8]

$\displaystyle{\bm{V}}^{\dagger}_{\bm{k}\alpha}{\bm{\sigma}}_{3}{\bm{V}}^{\phantom{\dagger}}_{\bm{k}\alpha^{\prime}}$
$\displaystyle=+\delta_{\alpha\alpha^{\prime}},$
$\displaystyle{\bm{W}}^{\dagger}_{-\bm{k}\alpha}{\bm{\sigma}}_{3}{\bm{W}}^{\phantom{\dagger}}_{-\bm{k}\alpha^{\prime}}$
$\displaystyle=-\delta_{\alpha\alpha^{\prime}},$
$\displaystyle{\bm{W}}^{\dagger}_{-\bm{k}\alpha}{\bm{\sigma}}_{3}{\bm{V}}^{\phantom{\dagger}}_{\bm{k}\alpha^{\prime}}$
$\displaystyle=0.$ (S13)

One can then write the Hamiltonian in terms of diagonalized bosons,
$\gamma_{\bm{k}\alpha}$, as [8]

$\displaystyle H\equiv NS(S+1)\epsilon_{\rm cl}+N\epsilon_{\rm
qu}+\sum_{\bm{k}\alpha}\epsilon_{\bm{k}\alpha}{\gamma}^{\dagger}_{\bm{k}\alpha}{\gamma}^{\phantom{\dagger}}_{\bm{k}\alpha}+O(S^{1/2}),$
(S14)

where we have identified the energy per site from quantum zero-point motion,
$\epsilon_{\rm qu}$, as

$\displaystyle\epsilon_{\rm qu}$
$\displaystyle\equiv\frac{1}{2N}\sum_{\bm{k}\alpha}\epsilon_{\bm{k}\alpha}$
(S15)

####  B.3.1 Spin Waves at $O(1/S)$

Figure S10:  Linear spin-wave spectrum of the $S=1/2$ nearest-neighbor model
[Eq. (S1)] in the type III antiferromagnet along a high-symmetry path with
$J>0$ and $\Gamma=0$. Spectra for several values of $K/J$ and for three
symmetry related domains are shown. We see that some of the line nodes present
when $K/J=0$ are lifted by finite $K/J$, but those that are perpendicular to
the stacking direction remain, along with several pseudo-Goldstone modes.

The type III ordering observed in K2IrCl6 can be described using a four-site
magnetic unit cell using the same four atoms of the conventional cubic cell.
We write

$\displaystyle\bm{r}_{1}$ $\displaystyle=\bm{0}$ $\displaystyle\bm{r}_{2}$
$\displaystyle=\bm{a}_{1}$ $\displaystyle\bm{r}_{3}$
$\displaystyle=\bm{a}_{2}$ $\displaystyle\bm{r}_{4}$
$\displaystyle=\bm{a}_{3}$

where $\bm{a}_{1}=a(\bm{\hat{y}}+\bm{\hat{z}})/2$,
$\bm{a}_{2}=a(\bm{\hat{x}}+\bm{\hat{z}})/2$ and
$\bm{a}_{3}=(\bm{\hat{x}}+\bm{\hat{y}})/2$ are the usual primitive lattice
vectors. The type III domains can be described by ordering directions on each
sublattice $\bm{\hat{z}}_{1}$, $\bm{\hat{z}}_{2}$, $\bm{\hat{z}}_{3}$,
$\bm{\hat{z}}_{4}$ and new lattice translations $\bm{A}_{1}$, $\bm{A}_{2}$,
$\bm{A}_{3}$, as given below

$\begin{array}[]{c|cccc|ccc}\textrm{Domain}&\bm{\hat{z}}_{1}&\bm{\hat{z}}_{2}&\bm{\hat{z}}_{3}&\bm{\hat{z}}_{4}&\bm{A}_{1}&\bm{A}_{2}&\bm{A}_{3}\\\
\hline\cr
X&\pm\bm{\hat{x}}&\mp\bm{\hat{x}}&\pm\bm{\hat{x}}&\mp\bm{\hat{x}}&a\bm{\hat{x}}+\bm{a}_{1}&a\bm{\hat{y}}&a\bm{\hat{z}}\\\
&\pm\bm{\hat{x}}&\mp\bm{\hat{x}}&\mp\bm{\hat{x}}&\pm\bm{\hat{x}}&a\bm{\hat{x}}+\bm{a}_{1}&a\bm{\hat{y}}&a\bm{\hat{z}}\\\
Y&\pm\bm{\hat{y}}&\mp\bm{\hat{y}}&\mp\bm{\hat{y}}&\pm\bm{\hat{y}}&a\bm{\hat{x}}&a\bm{\hat{y}}+\bm{a}_{2}&a\bm{\hat{z}}\\\
&\pm\bm{\hat{y}}&\pm\bm{\hat{y}}&\mp\bm{\hat{y}}&\mp\bm{\hat{y}}&a\bm{\hat{x}}&a\bm{\hat{y}}+\bm{a}_{2}&a\bm{\hat{z}}\\\
Z&\pm\bm{\hat{z}}&\pm\bm{\hat{z}}&\mp\bm{\hat{z}}&\mp\bm{\hat{z}}&a\bm{\hat{x}}&a\bm{\hat{y}}&a\bm{\hat{z}}+\bm{a}_{3}\\\
&\pm\bm{\hat{z}}&\mp\bm{\hat{z}}&\pm\bm{\hat{z}}&\mp\bm{\hat{z}}&a\bm{\hat{x}}&a\bm{\hat{y}}&a\bm{\hat{z}}+\bm{a}_{3}\\\
\end{array}$

Note that we have provided four sign structures for the
$\bm{\hat{z}}_{\alpha}$ of each of the $X$, $Y$ and $Z$ domains, giving in
total twelve domains (as expected). The four sign structures correspond to
different type III states related by translational or time-reversal symmetries
(leaving the spectrum and scattering intensity unchanged), while the $X$, $Y$
and $Z$ families are related by a three-fold rotation.

For each choice of domain and sign structure we can construct the Fourier
transforms of the local exchange matrices,
$\mathcal{J}^{\mu\nu}_{\bm{k},\alpha\alpha^{\prime}}$, as outlined in Sec. B.3
and from those the ${\bm{A}}_{\bm{k}}$ and ${\bm{B}}_{\bm{k}}$ matrices. Each
of these are four by four matrices, leading to a final eight by eight
eigenproblem encoded in ${\bm{M}}_{\bm{k}}$ to determine the spin-wave
energies and intensities. Practically, we tabulate a list of nearest-neighbor
bonds originating within our unit cell and construct these matrices
numerically for each wave-vector of interest.

Results for the Heisenberg-Kitaev limit ($\Gamma/J=0$) at $S=1/2$ are shown in
Fig. S10, with three inequivalent domains indicated. We see that for $K/J=0$
one has nodal lines along the $[1,0,0]$, $[0,1,0]$ and $[0,0,1]$ directions
(or equivalents). As $K/J$ is rendered finite we see that the nodal lines
perpendicular to the stacking are preserved, while those along the stacking
direction are lifted. Further, in addition to these nodal lines we see several
pseudo-Goldstone modes at wave-vectors characteristic of the type III order
(such as $(1,0.5,0)$ and equivalents). The presence of finite $\Gamma$ with
$|\Gamma|<K/2$ does affect the linear spin-wave spectrum, but only away from
high-symmetry lines (for details see Sec. B.9), and the nodal lines and
pseudo-Goldstone modes are preserved. We thus see the linear spin-wave
spectrum is qualitatively different than what is observed experimentally due
to the presence of a large number of zero modes.

As the accidental degeneracy that gives rise to these zero modes can be lifted
by a finite second neighbor exchange $J_{2}$, we should expect that the nodal
lines will be gapped for $J_{2}/J>0$ (stabilizing the type III order). The
effect of finite $J_{2}/J$ is shown in Fig. S11 for several values of
$J_{2}/J$. We see that indeed the remaining nodal lines are lifted, but
several of the pseudo-Goldstone modes associated with the type III ordering
wave-vectors remain. This is a consequence of the accidental degeneracy of the
of the three type III domains that is not lifted by finite $J_{2}/J$. As
discussed in Sec. B.2.3, this degeneracy is robust and is not lifted by _any_
bilinear exchanges consistent with the symmetry of the model. We thus expect
these pseudo-Goldstone modes to be persistent at the linear spin-wave level
even if additional anisotropic exchanges (e.g. a second neighbor Kitaev
coupling) or further neighbor interactions are included. As with the
Heisenberg-Kitaev model itself, the presence of these low-lying bands and
gapless pseudo-Goldstone modes disagrees qualitatively with what is observed
experimentally.

Figure S11:  Linear spin-wave spectrum of the $S=1/2$ nearest-neighbor model
in the type III antiferromagnet along a high-symmetry path with $J>0$,
$\Gamma=0$ and finite $J_{2}/J>0$. Spectra for several values of $J_{2}/J$
with $K/J=0.1$ and for three symmetry related domains are shown. We see that
some of the line nodes present when $J_{2}/J=0$ are lifted by finite
$J_{2}/J$, but several pseudo-Goldstone modes remain.

####  B.3.2 Order by Quantum Disorder

In the Heisenberg limit, quantum fluctuations at $O(1/S)$ are expected to
select the collinear type III ordered from the stacked states described in
Sec. B.2.2. This was explored in Ref. [47] for type I and III states, as well
as a family of incommensurate spirals. We have further confirmed that at
$O(1/S)$ that among all collinear stacked states up to period-12 in the
stacking direction that type III remains the minimum. We have also confirmed
that non-collinear states obtained by continuously mixing different type I
states, different type III states as well as type I and type III states do not
yield a different selection. These results hold true at finite $K/J$ and for
small $\Gamma$ as well.

Superficially this would suggest order-by-quantum-disorder selection at
$O(1/S)$ could be responsible for the appearance of type III order in K2IrCl6.
However, as shown in Ref. 46, $O(1/S^{2})$ and higher corrections likely
change this result in the Heisenberg limit. How this result extends to finite
$K/J$ or finite $\Gamma/J$ is unclear. Calculations at $O(1/S^{2})$ following
the strategy of Ref. [46] are inconclusive when $K/J\neq 0$, with the type I
state developing negative energy modes at $S=1/2$. We leave exploration of the
$O(1/S^{2})$ order-by-quantum-disorder selection in the degenerate phase of
the nearest-neighbor model [Eq. (S1)] to future work.

####  B.3.3 Dynamical Structure Factor

The inelastic neutron scattering intensity per spin is determined by the spin-
spin correlation function

$\mathcal{S}_{\mu\nu}(\bm{k},\omega)\equiv\frac{1}{2\pi N}\int dte^{-i\omega
t}\braket{M^{\mu}_{-\bm{k}}M^{\nu}_{\bm{k}}(t)}$

where $\bm{M}_{\bm{k}}\equiv
g\mu_{B}\sum_{i}e^{-i\bm{k}\cdot\bm{r}_{i}}\bm{S}_{i}$ is the magnetization
operator at wave-vector $\bm{k}$ and $N$ is the total number of spins.
Explicitly the intensity per spin is given by

$\displaystyle I(\bm{k},\omega)$
$\displaystyle=\frac{k_{f}}{k_{i}}\left(\frac{\gamma
r_{0}}{2\mu_{B}}\right)^{2}F(k)^{2}\sum_{\mu\nu}\left(\delta_{\mu\nu}-\hat{k}_{\mu}\hat{k}_{\nu}\right)\mathcal{S}_{\mu\nu}(\bm{k},\omega)$

where $r_{0}=\mu_{0}e^{2}/(4\pi m_{e})\approx 2.818\cdot 10^{-15}$ m is the
classical electron radius, $\mu_{B}$ is the Bohr magneton, $\gamma=1.913$,
$F(k)$ is the magnetic form factor and $k_{i}$, $k_{f}$ are the incoming and
outgoing neutron wave-vectors. It is useful to write $(\gamma
r_{0})^{2}\approx 0.291$ barns where $1$ barn $=10^{-28}$ m.

Within linear spin-wave theory the dynamical structure factor can be expressed
in terms of the transverse-transverse part of the spin-spin correlation
function. At zero temperature the dynamical structure factor with linear spin-
wave theory takes the form [19]

$\mathcal{S}_{\mu\nu}(\bm{k},\omega)=\sum_{\bm{k},n}W^{\mu\nu}_{\bm{k},n}\delta(\omega-\omega_{\bm{k},n})$

where we defined the weights for each spin-wave mode as
$W^{\mu\nu}_{\bm{k},n}\equiv\left(S\mu_{B}g^{2}/N_{s}\right)\Phi^{n,\mu}_{\bm{k}}\bar{\Phi}^{n,\mu^{\prime}}_{\bm{k}}$
where $N_{s}$ is the number magnetic sublattices. The quantities
$\Phi^{n,\mu}_{\bm{k}}$ are defined as

$\Phi^{n,\mu}_{\bm{k}}\equiv\sum_{\alpha}e^{-i\bm{k}\cdot\bm{r}_{\alpha}}\left(X^{\alpha}_{\bm{k},n}\hat{e}^{\mu}_{\alpha,-}+Y^{\alpha}_{\bm{k},n}\hat{e}^{\mu}_{\alpha,+}\right)$

The vectors ${\bm{X}}_{\bm{k}}$ and ${\bm{Y}}_{\bm{k}}$ are blocks of the
eigenvectors of ${\bm{M}}_{\bm{k}}$ with
$\bm{V}_{\bm{k},n}=(\bm{X}_{\bm{k},n},\bm{Y}_{\bm{k},n})$.

Note that the gap observed experimentally is much larger than the sample
temperature, $\Delta\sim 0.65$ meV compared to $T\sim 0.25$ K and so we are
justified in taking the zero temperature limit. Corrections due finite
temperature effects will be bounded by the Bose factor $n_{B}(\Delta)$, which
is entirely negligible at $0.25$ K.

###  B.4 Non-Linear Spin Wave Theory

Going to next order in $1/S$, one obtains a more complex Hamiltonian for the
magnon excitations. This can be written as a sum of the usual two-magnon
terms, as well as generally three- and four-magnon interactions. Following the
insights from Sec. B.3 we will consider only the case where $\Gamma=0$. With
this restriction the three-magnon interaction terms vanish and the spin-wave
interactions simplify considerably. One obtains

$H=NS(S+1)\epsilon_{\rm cl}+H_{2}+H_{4}$ (S16)

where we define the individual pieces in symmetrized form as

$\displaystyle{H}_{2}$
$\displaystyle=\frac{1}{2}\sum_{\alpha\beta}\sum_{\bm{k}}\left[{A}_{\bm{k}}^{\alpha\beta}{a}^{\dagger}_{\bm{k}\alpha}{a}^{\phantom{\dagger}}_{\bm{k}\beta}+{A}_{-\bm{k}}^{\beta\alpha}{a}^{\phantom{\dagger}}_{-\bm{k}\alpha}{a}^{\dagger}_{-\bm{k}\beta}+\left({B}^{\alpha\beta}_{\bm{k}}{a}^{\dagger}_{\bm{k}\alpha}{a}^{\dagger}_{-\bm{k}\beta}+\bar{B}^{\alpha\beta}_{\bm{k}}{a}^{\phantom{\dagger}}_{-\bm{k}\beta}{a}^{\phantom{\dagger}}_{\bm{k}\alpha}\right)\right],$
$\displaystyle{H}_{4}$
$\displaystyle=\frac{1}{N_{c}}\sum_{\alpha\beta\mu\nu}\sum_{\bm{k}\bm{k}^{\prime}\bm{q}}\left[\frac{1}{(2!)^{2}}{V}_{\bm{k}\bm{k}[\bm{q}]}^{\alpha\beta\mu\nu}{a}^{\dagger}_{\bm{k}+\bm{q},\alpha}{a}^{\dagger}_{\bm{k}^{\prime}-\bm{q},\beta}{a}^{\phantom{\dagger}}_{\bm{k}^{\prime}\mu}{a}^{\phantom{\dagger}}_{\bm{k}\nu}+\frac{1}{3!}\left({D}_{\bm{k}\bm{k}^{\prime}\bm{q}}^{\alpha\beta\mu\nu}{a}^{\dagger}_{\bm{k}\alpha}{a}^{\dagger}_{\bm{k}^{\prime}\beta}{a}^{\dagger}_{\bm{q}\mu}{a}^{\phantom{\dagger}}_{\bm{k}+\bm{k}^{\prime}+\bm{q},\nu}+\textrm{h.c.}\right)\right].$

In terms of the local exchange matrices [Eq. (S4)] one can write

$\displaystyle{A}_{\bm{k}}^{\alpha\beta}$
$\displaystyle=S\left(\mathcal{J}^{+-}_{\bm{k},\alpha\beta}-\delta_{\alpha\beta}\sum_{\mu}\mathcal{J}^{00}_{\bm{0},\alpha\mu}\right),$
(S18a) $\displaystyle{B}_{\bm{k}}^{\alpha\beta}$
$\displaystyle=S\mathcal{J}^{++}_{\bm{k},\alpha\beta},$ (S18b)
$\displaystyle{V}^{\alpha\beta\mu\nu}_{\bm{k}\bm{k}^{\prime}[\bm{q}]}$
$\displaystyle=\left(\delta_{\alpha\mu}\delta_{\beta\nu}\mathcal{J}^{00}_{\bm{k}-\bm{k}^{\prime}+\bm{q},\alpha\beta}+\delta_{\alpha\nu}\delta_{\beta\mu}\mathcal{J}^{00}_{\bm{q},\alpha\beta}\right)-\left(\delta_{\mu\nu}\delta_{\mu\beta}\mathcal{J}^{+-}_{\bm{k}+\bm{q},\alpha\nu}+\delta_{\alpha\beta}\delta_{\alpha\mu}\mathcal{J}^{+-}_{\bm{k},\alpha\nu}\right),$
(S18c) $\displaystyle{D}^{\alpha\beta\mu\nu}_{\bm{k}\bm{k}^{\prime}\bm{q}}$
$\displaystyle=-\frac{3}{4}\left(\delta_{\alpha\mu}\delta_{\alpha\nu}\mathcal{J}^{++}_{\bm{k}^{\prime},\beta\alpha}+\delta_{\mu\beta}\delta_{\nu\beta}\mathcal{J}^{++}_{\bm{k},\alpha\beta}\right),$
(S18d)

where the four-magnon vertices have been left unsymmetrized for brevity. At
leading order in perturbation theory, these magnon interaction terms
renormalize the linear spectrum, giving corrections to ${\bm{A}}_{\bm{k}}$ and
${\bm{B}}_{\bm{k}}$

$\displaystyle\Delta A^{\alpha\beta}_{\bm{k}}$
$\displaystyle=\frac{1}{N_{c}}\sum_{\bm{q}}\sum_{\mu\nu}\left[V^{\alpha\mu\nu\beta}_{\bm{k}\bm{q}[\bm{0}]}\braket{{a}^{\dagger}_{\bm{q}\mu}{a}^{\phantom{\dagger}}_{\bm{q}\nu}}_{0}+\frac{1}{2}\left(D^{\alpha\mu\nu\beta}_{\bm{k},-\bm{q},\bm{q}}\braket{{a}^{\dagger}_{-\bm{q}\mu}{a}^{\dagger}_{\bm{q}\nu}}_{0}+\bar{D}^{\beta\mu\nu\alpha}_{\bm{k},\bm{q},-\bm{q}}\braket{{a}^{\phantom{\dagger}}_{\bm{q}\mu}{a}^{\phantom{\dagger}}_{-\bm{q}\nu}}_{0}\right)\right],$
$\displaystyle\Delta B^{\alpha\beta}_{\bm{k}}$
$\displaystyle=\frac{1}{N_{c}}\sum_{\bm{q}}\sum_{\mu\nu}\left[D^{\mu\alpha\beta\nu}_{\bm{q},\bm{k},-\bm{k}}\braket{{a}^{\dagger}_{\bm{q}\mu}{a}^{\phantom{\dagger}}_{\bm{q}\nu}}_{0}+\frac{1}{2}V^{\alpha\beta\nu\mu}_{\bm{q},-\bm{q},[\bm{k}-\bm{q}]}\braket{{a}^{\phantom{\dagger}}_{\bm{q}\mu}{a}^{\phantom{\dagger}}_{-\bm{q}\nu}}_{0}\right].$

where $\braket{\dots}_{0}$ is an average with respect to the linear spin-wave
Hamiltonian, $H_{2}$.

####  B.4.1 Self-consistent Spin Waves

Due to the large number of zero modes present in the linear spin-wave
spectrum, divergences in the interactions can make the leading perturbative
result unreliable [47]. One way to resolve these divergences is by including
the interaction terms self-consistently [46]. With this in mind we consider a
adding and subtracting a quadratic piece to the Hamiltonian, writing

$H_{2}+H_{4}=\left(H_{2}+\delta H_{2}\right)+\left(H_{4}-\delta H_{2}\right).$
(S20)

where we define $\delta H_{2}$ as

$\delta
H_{2}\equiv\sum_{\alpha\beta}\sum_{\bm{k}}\left[\delta{A}_{\bm{k}}^{\alpha\beta}{a}^{\dagger}_{\bm{k}\alpha}{a}^{\phantom{\dagger}}_{\bm{k}\beta}+\frac{1}{2}\left(\delta{B}^{\alpha\beta}_{\bm{k}}{a}^{\dagger}_{\bm{k}\alpha}{a}^{\dagger}_{-\bm{k}\beta}+\delta\bar{B}^{\alpha\beta}_{\bm{k}}{a}^{\phantom{\dagger}}_{-\bm{k}\beta}{a}^{\phantom{\dagger}}_{\bm{k}\alpha}\right)\right].$
(S21)

We now consider $H_{2}+\delta H_{2}$ as the linear spin-wave problem,
evaluating the energies and eigenvectors with respect to this Hamiltonian. We
treat $\delta H_{2}$ perturbatively; it appears once at leading order and thus
simply produces an additive contribution to ${\bm{A}}_{\bm{k}}$ and
${\bm{B}}_{\bm{k}}$. We define our self-consistent mean-field theory by
choosing the perturbations $\delta{\bm{A}}_{\bm{k}}$ and
$\delta{\bm{B}}_{\bm{k}}$ so that they _cancel_ the corrections produced by
the interactions. This can be done via an iterative process

$\displaystyle\delta A^{\alpha\beta}_{\bm{k}}$
$\displaystyle\leftarrow\frac{1}{N_{c}}\sum_{\bm{q}}\sum_{\mu\nu}\left[V^{\alpha\mu\nu\beta}_{\bm{k}\bm{q}[\bm{0}]}\braket{{a}^{\dagger}_{\bm{q}\mu}{a}^{\phantom{\dagger}}_{\bm{q}\nu}}_{\textrm{MF}}+\frac{1}{2}\left(D^{\alpha\mu\nu\beta}_{\bm{k},-\bm{q},\bm{q}}\braket{{a}^{\dagger}_{-\bm{q}\mu}{a}^{\dagger}_{\bm{q}\nu}}_{\textrm{MF}}+\bar{D}^{\beta\mu\nu\alpha}_{\bm{k},\bm{q},-\bm{q}}\braket{{a}^{\phantom{\dagger}}_{\bm{q}\mu}{a}^{\phantom{\dagger}}_{-\bm{q}\nu}}_{\textrm{MF}}\right)\right],$
$\displaystyle\delta B^{\alpha\beta}_{\bm{k}}$
$\displaystyle\leftarrow\frac{1}{N_{c}}\sum_{\bm{q}}\sum_{\mu\nu}\left[D^{\mu\alpha\beta\nu}_{\bm{q},\bm{k},-\bm{k}}\braket{{a}^{\dagger}_{\bm{q}\mu}{a}^{\phantom{\dagger}}_{\bm{q}\nu}}_{\textrm{MF}}+\frac{1}{2}V^{\alpha\beta\nu\mu}_{\bm{q},-\bm{q},[\bm{k}-\bm{q}]}\braket{{a}^{\phantom{\dagger}}_{\bm{q}\mu}{a}^{\phantom{\dagger}}_{-\bm{q}\nu}}_{\textrm{MF}}\right].$

where $\braket{\dots}_{\textrm{MF}}$ is evaluated with respect to the shifted
Hamiltonian $H_{2}+\delta H_{2}$ and depends on the current values of
$\delta{\bm{A}}_{\bm{k}}$ and $\delta{\bm{B}}_{\bm{k}}$. Once converged, one
has an effective quadratic model of the spin-wave spectrum, with interaction
effects encoded in $\delta{\bm{A}}_{\bm{k}}$ and $\delta{\bm{B}}_{\bm{k}}$.

####  B.4.2 Spin Waves at $O(1/S^{2})$

We follow the self-consistent scheme discussed in Sec. B.4.1. As for the
linear spin-wave case, we tabulate a list of nearest-neighbor bonds without
our unit cell and construct the matrices and interaction vertices numerically
for each wave-vector needed in the self-consistent sums.

Practically, we iterate this loop for each wave-vector until the corrections
have stopped changing to an absolute tolerance of size $10^{-8}$ in all matrix
elements. To resolve some of the divergences due to zero modes present when
interactions are absent, the $\delta{\bm{A}}_{\bm{k}}$ and
$\delta{\bm{B}}_{\bm{k}}$ matrices were initialized to non-zero values at the
start of the iterative loop with $\delta{\bm{A}}_{\bm{k}}=\mu{\bm{1}}$ and
$\delta{\bm{B}}_{\bm{k}}={\bm{0}}$ with $\mu>0$. We have confirmed that the
converged result is independent of the precise choice of $\mu$.

The resulting $\delta{\bm{A}}_{\bm{k}}$ and $\delta{\bm{B}}_{\bm{k}}$ yield a
quadratic spin-wave Hamiltonian with ${\bm{M}}^{\rm
eff}_{\bm{k}}={\bm{M}}_{\bm{k}}+\delta{\bm{M}}_{\bm{k}}$ which can be
diagonalized as in Sec. B.3 to yield spin-wave energies and wave-functions.

Figure S12:  Self-consistent non-linear spin-wave spectrum of the $S=1/2$
nearest-neighbor model [Eq. (S1)] in the type III antiferromagnet along a
high-symmetry path with $J>0$ and $\Gamma=0$. Spectra for several values of
$K/J$ and for three symmetry related domains are shown. We see that the
remaining line nodes present when $K/J=0$ in linear spin-wave theory are
lifted by interaction corrections, as are all remaining pseudo-Goldstone
modes.

Figure S13:  Linear and self-consistent non-linear spin-wave spectrum of the
$S=1/2$ nearest-neighbor model [Eq. (S1)] in the type III antiferromagnet
along a high-symmetry path for the best fit parameters $J=0.74$ meV and
$K=0.15$ meV. Three symmetry related domains are shown. We see that the
presence of spin-wave interactions introduces a large gap to all of the nodal
lines and pseudo-Goldstone modes present in the linear spin-wave spectrum.

####  B.4.3 Corrections to the Dynamical Structure Factor

The higher order terms in the Holstein-Primakoff expansion also induce
corrections to the spin-spin correlation functions that appear in the
dynamical structure factor. When magnon decay is forbidden, this generally
involves modifications to the intensities from the transverse spin-spin
correlators, as well as contributions to the longitudinal spin-spin
correlators from the two-magnon continuum [36].

The renormalization of the transverse part depends on the normal and anomalous
on-site averages [54, 36]

$\displaystyle n_{\alpha}$
$\displaystyle\equiv\braket{{a}^{\dagger}_{\bm{r}\alpha}a_{\bm{r}\alpha}}_{0}$
$\displaystyle\delta_{\alpha}$
$\displaystyle\equiv\braket{{a}_{\bm{r}\alpha}a_{\bm{r}\alpha}}_{0}$

For the nearest-neighbor model with $\Gamma=0$ we have that $n_{\alpha}\equiv
n$ which is independent of sublattice and $\delta_{\alpha}=0$. In this
scenario the dynamical structure factor in the self-consistent theory can then
be computed from $H_{2}+\delta H_{2}$ as in standard linear-spin-wave theory
(see Sec. B.3.3), save with the an overall intensity correction

$I(\bm{k},\omega)\rightarrow\left(1-\frac{n}{2S}\right)I(\bm{k},\omega)$

We have included this correction in our comparisons to our experimental data.
We note that for our best fit parameters the size of $n$ is not particular
large, $n\sim 0.1$, and this correction only amounts to a $10\%$ reduction in
scattering intensity.

The contribution from the two-magnon continuum was found to be insignificant
in our self-consistent spin-wave theory, with the spectral weight amounting to
a small fraction of the contribution of the transverse parts. We have
therefore not included these contributions in our comparisons to the
experimental data. Further, we note that these contributions would only modify
the intensities for $\omega\geq 2\Delta\sim 1.3\ \textrm{meV}$ where
$\Delta\sim 0.65\ \textrm{meV}$ is the one-magnon gap. This intensity also
only increases slowly from that minimum, following roughly the two-magnon
density of states.

###  B.5 Effect of Tetragonal Distortion

If the symmetry of the crystal structure is lowered from cubic to tetragonal,
additional exchange interactions are allowed in the spin Hamiltonian. For
illustration we consider the allowed exchange interactions in $I4/mmm$ (#139)
which is one of the larger tetragonal subgroups of $Fm\bar{3}m$ (#225). A
standard analysis of the space group symmetry (taking $\bm{\hat{z}}$ to be the
tetragonal axis) yields

$\displaystyle H=\sum_{\braket{ij}_{z}}$
$\displaystyle\left[J_{\perp}\bm{S}_{i}\cdot\bm{S}_{j}+\left(K_{\perp}+A_{\perp}\right)S^{z}_{i}S^{z}_{j}+(-1)^{\sigma^{xy}_{ij}}\Gamma_{\perp}\left(S^{x}_{i}S^{y}_{j}+S^{y}_{i}S^{x}_{j}\right)\right]+$
$\displaystyle\sum_{\gamma=x,y}\sum_{\braket{ij}_{\gamma}}$
$\displaystyle\left[J\bm{S}_{i}\cdot\bm{S}_{j}+KS^{\gamma}_{i}S^{\gamma}_{j}+(-1)^{\sigma^{\alpha\beta}_{ij}}\Gamma\left(S^{\alpha}_{i}S^{\beta}_{j}+S^{\alpha}_{i}S^{\beta}_{j}\right)+A_{\perp}S^{z}_{i}S^{z}_{j}\right]$
(S23)

This model includes different $J$, $K$ and $\Gamma$ interactions in the
perpendicular to the tetragonal axis, as well as global XXZ anisotropy encoded
in $A_{\perp}$. Given a small lattice distortion we should loosely expect each
of these terms to be proportional to some component of the strain $\epsilon$,
as well as a spin-lattice coupling constant. The presence of these
anisotropies breaks the continuous degeneracy between the type III states and
other stacked classical ground states.

####  B.5.1 Effect on Magnon Gap in Type-III State

At the linear spin-wave level the magnon energies at the high-symmetry wave-
vectors $(1,0,0)$, $(0,1,0)$ and $(0,0,1)$ for the $Z$-domain of the type-III
order gives two distinct energies for each wave-vector

$\displaystyle\omega_{[001],1}$
$\displaystyle=2\sqrt{(K_{\perp}+A_{\perp})(2(J_{\perp}-J)+K_{\perp}-K+A_{\perp})},$
$\displaystyle\omega_{[001],2}$
$\displaystyle=2\sqrt{(K_{\perp}+A_{\perp})(2(J_{\perp}+J)+K_{\perp}+K+A_{\perp})},$
$\displaystyle\omega_{[100],1}=\omega_{[010],1}$
$\displaystyle=2\sqrt{(K_{\perp}-K+A_{\perp})(2J_{\perp}+K_{\perp}+A_{\perp}))},$
$\displaystyle\omega_{[100],2}=\omega_{[010],2}$
$\displaystyle=2\sqrt{(K+K_{\perp}+A_{\perp})(2J_{\perp}+K_{\perp}+A_{\perp})}.$

Note that the symmetric anisotropies have dropped out. For these wave-vectors
we also see that the effect of $K_{\perp}$ and $A_{\perp}$ cannot be
distinguished, with only $K_{\perp}^{\prime}\equiv K_{\perp}+A_{\perp}$
appearing. We see that, as expected, all zero modes have are gapped out when
tetragonal distortions such as $J_{\perp}\neq J$, $K_{\perp}\neq K$ or
$A_{\perp}\neq 0$ are included.

####  B.5.2 Best Fit to Tetragonal Model

Given the observed spectrum we can make reasonable inferences on which of
these modes would correspond to which experimental features. First, the
$\omega_{(001),1}$, $\omega_{(100),1}$ and $\omega_{(010),1}$ modes become
gapless when $J_{\perp}=J$, $K^{\prime}_{\perp}=K$. Following how these are
expected acquire gaps in Fig. S13 for the $Z$ domain, we identify
$\omega_{(001),1}\approx 0.65$ meV and
$\omega_{(100),1}=\omega_{(010),1}\approx 0.9$ meV. The $\omega_{(001),2}$
mode would then correspond to the higher-lying mode at $(0,0,0)$ which was
observed in ESR near $\omega_{(001),2}\approx 1.2\ \textrm{meV}$. The final
mode would then correspond to the high intensity mode near
$\omega_{(010),2}=\omega_{(100),2}\approx 1.6\ \textrm{meV}$ at $(1,1,0)$,
which also appears at $(1,0,0)$ and equivalents for this domain. These
experimental constraints do not uniquely determine the four spatially
anisotropic exchange parameters. In fact, this linear spin-wave model cannot
reproduce all four energies due to satisfying the relation
$\omega_{(001),1}^{2}+\omega_{(001),2}^{2}=\omega_{(100),1}^{2}+\omega_{(100),2}^{2}$
which is not satisfied for the experimentally determined energies.

To resolve this we take two steps: first, we only consider the two low lying
modes ($0.65$ meV and $0.9$ meV) and the highest lying mode ($1.6$ meV). To
find a unique solution we additionally require that the Curie-Weiss constant
$\Theta_{\rm CW}=-29.3$ K [41] is reproduced, with
$(2J+J_{\perp})+(2K+K_{\perp})/3\approx 2.52$ meV. This yields the solution

$\displaystyle J$ $\displaystyle=0.724\ \textrm{meV},$ $\displaystyle K$
$\displaystyle=0.106\ \textrm{meV},$ $\displaystyle J_{\perp}$
$\displaystyle=0.935\ \textrm{meV},$ $\displaystyle K_{\perp}^{\prime}$
$\displaystyle=0.203\ \textrm{meV}.$

These parameters reproduce the three stated spin-wave energies as well as the
Curie-Weiss temperature, but (necessarily) give an incorrect
$\omega_{(001),2}\approx 1.71$ meV instead of $1.2$ meV for the excluded ESR-
visible mode. This represents a moderate increase in the Heisenberg coupling
$J$ of about $29\%$ for the perpendicular plane $J_{\perp}/J=1.29$, but a near
doubling of the effective Kitaev coupling $K_{\perp}^{\prime}/K=1.93$. Note
that this could be partitioned arbitrarily into a smaller increase $K$, but a
corresponding non-zero value of $A_{\perp}$. The spectrum of this best fit
using linear spin-wave theory is shown in Fig. S14.

Figure S14:  Linear spin-wave spectrum of the nearest-neighbor model with
tetragonal distortion [Eq. (B.5)] in the type III antiferromagnet along a
high-symmetry path for parameters $J=0.724$ meV, $K=0.106$ meV,
$J_{\perp}=0.935$ meV and $K_{\perp}=0.203$ meV with
$\Gamma=\Gamma_{\perp}=A_{\perp}=0$. Three symmetry related domains are shown
(with the tetragonal distortion following the domain direction). We see that
the anisotropy introduces a large gap to all of the nodal lines and pseudo-
Goldstone modes present in the undistorted case.

###  B.6 Estimating the Magneto-elastic Coupling

As structural distortions due to magneto-elastic couplings can also
(effectively) lower the symmetry and gap out the pseudo-Goldstone modes, it is
worthwhile to estimate (at least roughly) the expected energy scale for such
effects in K2IrCl6. We will consider two approaches: first, we will estimate
the change in exchange constants due to the observed tetragonal distortion
below $T_{\rm N}$. Next we will directly estimate the distortion and induced
exchanges based on more microscopic arguments. For both we find the effects of
magneto-elastic coupling are far too small to account for the size of the gap
observed in K2IrCl6.

We follow the treatment of Liu and Khaliullin [30] that was developed to
explain features of another iridate compound, Sr2IrO4, with similar on-site
and exchange physics [39]. Consider a magneto-elastic coupling of the form

$\tilde{g}\sum_{ij}\sum_{\mu\nu}\epsilon_{\mu\nu}\sum_{\mu^{\prime}\nu^{\prime}}C^{\mu\nu,\mu^{\prime}\nu^{\prime}}_{ij}S^{\mu^{\prime}}_{i}S^{\nu^{\prime}}_{j}$

where $\epsilon_{\mu\nu}$ is the (symmetric) strain, $\tilde{g}$ is the
overall strength of the magneto-elastic couplings and
$C^{\mu\nu,\mu^{\prime}\nu^{\prime}}_{ij}$ represent dimensionless factors
that form the correct spin quadrupoles out of the bond operators
$S^{\mu^{\prime}}_{i}S^{\nu^{\prime}}_{j}$ to couple to the
$\epsilon_{\mu\nu}$ strain component. We expect $\tilde{g}$ to be decomposable
into $\tilde{g}=\kappa g$ where $g$ is an atomic energy scale and $\kappa$
relates to the energy scale of exchange interactions. For the Ir4+ we use an
estimate [30] of $g\sim 5$ eV, while we assume $\kappa$ scales with the
exchange $J$. Since the exchange energy scale in K2IrCl6 is a factor of
$10^{2}$ smaller than in Sr2IrO4 [23, 30], we will take $\kappa\sim 5\cdot
10^{-5}$, thus yielding $\tilde{g}\sim 0.25$ meV.

The observed splitting of the structural Bragg peaks (see Fig. S3) below
$T_{\rm N}$ implies a strain of roughly $\epsilon\sim 10^{-4}$. For this value
of strain the modification of the magnetic exchanges would then be of order
$\tilde{g}\epsilon\sim 0.025$ $\mu$eV. This is negligible relative to all
magnetic energy scales in K2IrCl6. Roughly, one would expect for type I
pseudo-Goldstone modes [40] that the induced gap would scale as
$\Delta\sim\sqrt{JD}$ with $D\sim\tilde{g}\epsilon$. This yields a gap of
$\sim 5$ $\mu$eV; as expected, far too small to account for the observed gap
in K2IrCl6.

Consider now the expected distortion due to $\tilde{g}$. This is given by
$\epsilon\sim\tilde{g}/\mathcal{K}$ where $\mathcal{K}$ defines the elastic
energy $\sim\frac{1}{2}\mathcal{K}\epsilon^{2}$. If we estimate that the
elastic constants in K2IrCl6 are comparable to that in Sr2IrO4, then we would
have an induced exchange $D\sim\tilde{g}^{2}/\mathcal{K}$ and thus gaps of
size $\Delta\sim\sqrt{JD}=\tilde{g}\sqrt{J/\mathcal{K}}$. With $J$ and
$\tilde{g}$ two orders of magnitude smaller in K2IrCl6, this yields a factor
of $10^{3}$ in the gap relative to the estimate of Liu and Khaliullin [30],
i.e. $\Delta\sim O(\mu{\rm eV})$ as above. One might expect a smaller elastic
constant in K2IrCl6, e.g. due to some evidence for proximity to a structural
distortion [51], but several orders of magnitude difference in $\mathcal{K}$
would be required to render the estimate for $\Delta$ significant. Note that
for $\mathcal{K}\sim 100$ GPa yields an estimate of $\epsilon\sim 10^{-6}$ or
$10^{-5}$; not too far from the observed $10^{-4}$ distortion.

###  B.7 Dynamics in Paramagnetic Phase

Accessing dynamical properties in the paramagnetic regime when $T>T_{N}$ is
difficult when quantum effects are included. We therefore adopt a classical
approach for the spin dynamics in the regime, including renormalizations in
the energy scale and intensities to partly mimic some quantum features.

We consider the Hamiltonian Eq. (S1) augmented with second-neighbor exchange
[Eq. (S2)] in the classical limit where the spins are unit length vectors
$|\bm{S}_{i}|^{2}=1$. We assume they follow conventional Landau-Lifshitz
dynamics

$\frac{d\bm{S}_{i}}{dt}=-\bm{S}_{i}\times\frac{\partial
H}{\partial\bm{S}_{i}},$ (S25)

where $-\partial H/\partial\bm{S}_{i}\equiv\bm{B}_{i}$ is the (local) exchange
field due to the neighboring spins. The initial conditions $\bm{S}_{i}(0)$ for
this equation are drawn from a thermal distribution at temperature $T$ using
Monte Carlo sampling. Once a sample of trajectories $\bm{S}_{i}(t)$ are
obtained, the part of the dynamical structure factor relevant for inelastic
neutron scattering is given by

$\mathcal{S}_{\rm
cl}(\bm{k},\omega)=\sum_{\mu\nu}\left(\delta_{\mu\nu}-\hat{k}_{\mu}\hat{k}_{\nu}\right)\braket{\bar{S}^{\mu}_{\bm{k}}(\omega)S^{\nu}_{\bm{k}}(\omega)},$

where $\bm{S}_{\bm{k}}(\omega)$ is the Fourier transform of $\bm{S}_{i}(t)$ in
both space and time. To account for our spins being $S=1/2$ we simply rescale
frequencies by a factor of $S$, ensuring (for example) that in the low
temperature limit the classical and quantum spin-wave frequencies agree.

To partially account quantum effects we also multiply by an energy dependent
correction factor [53, 45]

$F_{\rm qu}(\omega)\equiv\beta\omega\left(1+n_{B}(\omega)\right),$

where $n_{B}(\omega)=1/(e^{\beta\omega}-1)$ is the Bose factor. This factor
enforces the quantum version of the fluctuation-dissipation theorem relating
positive and negative frequencies [11], as $F_{\rm
qu}(-\omega)=e^{-\beta\omega}F_{\rm qu}(\omega)$.

Practically, we generate our thermal samples using Monte Carlo with single-
site heat-bath [35] and over-relaxation [14] updates, annealing down from high
temperature $O(10J)$ to the temperature of interest and then thermalizing for
the same number of sweeps. At the high temperatures of interest here only a
small number of sweeps, typically $O(10^{3})$, are necessary to reach
equilibrium.

For each sample we solve the coupled non-linear ordinary differential
equations of Eq. (S25) using an adaptive fourth-order Runge-Kutta method [16]
as implemented in the odeint library [3]. Fourier transforms in both space and
time, as implemented in the FFTW library [20, 9], are then used to obtain
$\bm{S}_{\bm{k}}(\omega)$ and then compute the dynamical structure factor. At
the temperatures of interest $O(10^{2})$ samples are sufficient to reach
convergence in both energy resolved and energy averaged quantities. The
adaptive time stepping was performed with (absolute and relative) error
tolerances of $10^{-8}$. The final trajectories were evaluated by
interpolation on a grid to obtain a frequency spacing of $\Delta\omega=0.04$
meV and a maximum frequency of $\omega_{\rm max}=3$ meV, to allow comparison
to the full measured spectrum.

For the simulations presented in the main text we used our best fit
parameters, $J=0.74$ meV, $K=0.15$ meV and $\Gamma=0$, A small finite $J_{2}$
can be included to classically stabilize a type III order (see Sec. B.8)
though it does not significantly affect any of our results. To make meaningful
comparisons to experiment, we performed simulations at $T=2T_{N}$ where
$T_{N}$ is the (classical) Néel temperature to match the 6 K used
experimentally. Each spectrum and slice is averaged in the same way as the
experimental data, as described in more detail in Sec. B.10.

###  B.8 Effect of $J_{2}$ on Néel temperature

Monte Carlo simulations for the static properties of our best model fit
(augmented with finite $J_{2}/J$) were carried out using parallel tempering
Monte Carlo [22, 49, 34]. Each temperature was annealed from the highest
temperature in our grid then thermalized before samples were taken. Single-
site heat-bath [35] and over-relaxation [14] updates were used throughout and
proved sufficient to equilibriate our system when combined with the parallel
tempering. Typical simulations were for systems of $L^{3}$ cubic cells
($N=4L^{3}$ spins) with $L=8,10,12$ and used $4\cdot 10^{4}$ sweeps for
annealing and thermalization and a further $10^{5}$ sweeps for sample
production. Each sweep consists of $N$ heat bath and $N$ over-relaxation
updates followed by a sweep attempting swaps of (random) neighboring
temperatures for each replica.

The classical transition temperature $T_{N}/J\sim 0.67$ was located via the
position of the maximum in heat capacity and in the onset of the static
structure factor at $(1,0.5,0)$ and equivalents. We confirmed that the
$(1,0,0)$ type I order parameter remained small through the transition. Note
that with $J_{2}/J=0$ the transition remains near $T/J\sim 0.66$ but the order
switches to type I due to the effects of order-by-thermal-disorder as has been
established in the Heisenberg case [21]. While, small finite $J_{2}/J>0$ has
little effect on $T_{N}/J$ the entire phase below $T_{N}$ switches to type III
almost immediately.

###  B.9 Constraints on Size of $\Gamma$ Exchange

At the level of linear spin-waves the effect of finite $\Gamma$ is limited.
For the spectra, slices and cuts shown in Fig. 2 of the main text the
contributions to ${\bm{A}}_{\bm{k}}$ and ${\bm{B}}_{\bm{k}}$ that are
$\propto\Gamma$ all vanish. We thus do not expect significant changes to our
results if small $\Gamma\neq 0$ is included.

To determine whether $\Gamma$ is not necessary to explain our data, we have
examined $I(\bm{k},\omega)$ at wave-vector $\bm{k}=[0.5,0.5,0.5]$ where
$\Gamma$ would expected to contribute maximally at the level of linear spin-
wave theory. We find no significant differences between our non-linear spin-
waves with $\Gamma=0$ and the experimental data. Finally, we note that ab-
initio calculations [25] estimate $\Gamma/J\approx 0.08$, small relative to
$J$ and a factor of three smaller than our best fit $K/J=0.2$.

Going beyond linear spin-waves a finite $\Gamma$ can induce spontaneous magnon
decay [54] when the one- and two-magnon excitations overlap in energy. As
discussed in Sec. B.4.3, the two-magnon intensity would be expected to begin
at $2\Delta\sim 1.3\ \textrm{meV}$ where $\Delta\sim 0.65\ \textrm{meV}$ is
the one-magnon gap. Any spontaneous decay would thus be limited to that
regime. The two magnon density of states grows slowly starting from its
minimum energy, further suppressing this decay channel. Given the resolution
of the inelastic data presented here it is unlikely this relatively small
level of decay could be resolved.

###  B.10 Details of Comparison Between Theoretical and Experimental Dynamical
Structure Factor

To perform a detailed comparison between the theoretical calculation and
experimental data we must model not only $I(\bm{k},\omega)$, but also the
finite energy resolution, averaging of the data in $\bm{k}$ and the magnetic
form factor of Ir4+.

####  B.10.1 Energy Resolution

Due to the finite energy resolution of the experimental measurement, we
convolve our theoretical result with a Gaussian lineshape with width
$\sigma_{\omega}=0.07$ meV, chosen to match the peak widths observed in the
individual $\bm{k}$ cuts [Main text, Fig. 2(g)]. The resulting intensity is