# Absence of thermal Hall half-quantization by field-suppression of magnetic
order in a Kitaev-Heisenberg model

Aman Kumar  Department of Theoretical Physics, Tata Institute of Fundamental
Research, Homi Bhabha Road, Navy Nagar, Mumbai 400005, India  Vikram Tripathi
Department of Theoretical Physics, Tata Institute of Fundamental Research,
Homi Bhabha Road, Navy Nagar, Mumbai 400005, India

(August 27, 2024)

###### Abstract

We investigate thermal Hall conductivity $\kappa_{xy}$ of a $J$-$K$ Kitaev-
Heisenberg model and a Zeeman field in the $(111)$ direction in the light of
the recent debate surrounding the possible re-emergence of Ising topological
order (ITO) and half-quantized $\kappa_{xy}/T$ upon field-suppression of long-
range magnetic order in Kitaev materials. We use the purification-based finite
temperature Tensor Network approach making no prior assumptions about the
nature of the excitations: Majorana, visons or spin waves. For purely Kitaev
interactions and fields $h/K\gtrsim 0.02$ sufficient to degrade ITO, the peak
$\kappa_{xy}/T$ monotonously decreases from half-quantization associated with
lower fields - a behavior physically consistent with vison fluctuation
corrections and qualitatively different from spin wave treatments. In the
$J$-$K$ model (ferro-$K$ and antiferro-$J$) in the vicinity of field-
suppressed magnetic order, we found $\kappa_{xy}/T$ to be significant, with
peak magnitudes exceeding half-quantization followed by a decrease with
increasing $h.$ We thus conclude that half-quantized thermal Hall effect, if
observed in the vicinity of field suppressed magnetic order in our model, is a
fine-tuning effect and is not associated with a Majorana Hall state with ITO.

Thermal Hall conductivity has become a very important transport diagnostic for
emergent phenomena in insulating magnetic systems such as the underdoped
cuprates Grissonnanche  _et al._ , frustrated magnets such as pyrochlores
Hirschberger  _et al._ (2015a), kagome systems Hirschberger  _et al._ (2015b),
and more recently a Kitaev material $\alpha$-RuCl3 Kasahara  _et al._ (2018);
Yokoi  _et al._ (2021); Bruin  _et al._ (2021); Czajka  _et al._ (2021).
Specifically in the Kitaev material, one set of studies Yokoi  _et al._
(2021); Kasahara  _et al._ (2018); Bruin  _et al._ (2021) has reported the
half-quantized thermal Hall effect at finite magnetic fields comparable to the
magnetic ordering scale, suggesting a re-emergent Majorana Hall insulator
state originally predicted Kitaev (2006) for the honeycomb Kitaev model at
very small magnetic fields and low temperatures which is a deconfined phase
with robust Ising topological order (ITO). This view has been challenged
Czajka  _et al._ (2021) in recent experiments with an in-plane magnetic field
where half-quantization has been ruled out. Nevertheless $\kappa_{xy}/T$ is
found to be large, sometimes even exceeding the half-quantized value at
temperatures and fields comparable to the magnetic ordering scale, implying
half-quantized peaks/plateaus of $\kappa_{xy}/T,$ if observed near field
suppressed magnetic order, are a consequence of fine tuning and do not
represent re-emergence of a deconfined ITO phase. Alternate explanations based
on a spin wave picture are then offered. Theoretical studies have relied on
various interacting quasiparticle pictures, with spin wave approximations
Chern  _et al._ (2021); Zhang  _et al._ (2021); Koyama and Nasu (2021);
McClarty  _et al._ (2018); Joshi (2018) most commonly used. Since spin wave
approximations in frustrated spin-$1/2$ magnets are best justified at high
fields or in the presence of long-range magnetic order, the possibility of
half-quantization at low temperatures near field-suppressed magnetic order
(not seen in spin wave treatments) cannot be completely ruled out. Other
theoretical approaches begin from the fractionalized quasiparticles, and
suppression of half-quantized $\kappa_{xy}/T$ by external magnetic field or
finite temperatures is understood from the perspective of vison (gauge field)
fluctuations that couple to the Kitaev spinons Guo  _et al._ (2020), or from
thermal excitation of bulk Majorana fermions Nasu  _et al._ (2017). Although
the fractionalization analyses are based around expansion from the ITO phase,
even here it is not clear whether ITO can re-emerge through field-suppression
of long-range magnetic order. These challenges motivated us to take a
different route of numerically estimating $\kappa_{xy}$ without making _any_
quasiparticle approximation - whether fractionalized or spin waves - which
forms the basis of existing theoretical approaches. We obtain the thermal Hall
current and conductivity numerically using purification-based Tensor Network
techniquesHauschild  _et al._ (2018). To our knowledge, this is the first use
of the Tensor Network technique for thermal Hall calculation of an interacting
many-body system.

We consider for concreteness a ferromagnetic (FM) Kitaev model with a
competing antiferromagnetic (AFM) Heisenberg interaction ($J$-$K$ model),
subjected to a Zeeman field along the $\mathbf{c}=(1,1,1)$ direction,

$\displaystyle H$ $\displaystyle=-K\\!\\!\\!\\!\\!\sum_{\langle
ij\rangle;\gamma-\text{links}}\\!\\!\\!\\!\\!S_{i}^{\gamma}S_{j}^{\gamma}+J\sum_{\langle
ij\rangle}\mathbf{S}_{i}\cdot\mathbf{S}_{j}+h\sum_{i,\gamma}S_{i}^{\gamma}.$
(1)

where $\gamma=x,y,z$ denotes the bond dependent spin axis on honeycomb
lattice. Thermal Hall conductivity of this model has recently been studied
using spin wave approximations Koyama and Nasu (2021). In the pure Kitaev
model ($J/K=0$), there is a fairly large window of magnetic fields (e.g.
$h_{\parallel\mathbf{c}}\gtrsim 0.02K$) Zhu  _et al._ (2018); Gohlke  _et al._
(2018) where ITO is significantly degraded (see also S.M. for a similar
signature in the topological entanglement entropy $\gamma(h)$) but interacting
spin wave analysis McClarty  _et al._ (2018) shows that magnons are good
quasiparticles for $h/K\gtrsim 1$ although not at significantly lower fields
due to interaction broadening. Fields in the range $0.02\lesssim h/K\lesssim
1$ thus constitute an interpolating region between the ITO and spin wave
regimes not easily lending to quasiparticle-based analyses, and we focus on
this window. Likewise for the spin density wave (SDW) ordered states of the
$J$-$K$ model, we are interested in the range of magnetic fields starting from
field-suppressed SDW to deeper in the spin-polarized phase, where the nature
of excitations may be changing.

Our main findings are as follows. In the Kitaev limit $J/K=0,$ near the lower
end of our field range $h\approx 0.02K,$ $\kappa_{xy}/T$ peaks at small finite
temperatures, with peak heights smaller but comparable to the half-quantized
value $\pi/12$ associated with the ITO phase (see Fig. 2). With increasing $h$
as well as temperature $T,$ $\kappa_{xy}/T$ declines sharply as expected from
vison fluctuations Guo  _et al._ (2020). Thermal excitation of the bulk
majoranas also suppresses half-quantization within the ITO phase Nasu  _et
al._ (2017); however we are outside this regime. Our result qualitatively
differs with the spin wave calculations in Ref. McClarty  _et al._ (2018)
where a non-monotonous temperature dependence with a peak-dip feature was
reported. Another key difference from earlier spin wave or fractionalization
treatments is that in our case $\kappa_{xy}/T$ decreases to zero when either
$h$ or $T$ becomes large, showing no sign of saturation at high temperatures.
In the magnetically ordered phase, we studied $\kappa_{xy}/T$ at intermediate
fields where ITO has been predicted to re-emerge following suppression of
magnetic order. In this field regime, we found $\kappa_{xy}/T$ to large, with
peak values that can even exceed half-quantization. Further increasing the
magnetic field decreases the extrema, which monotonously approach zero at high
fields like in the Kitaev limit. We posit that half-quantized thermal Hall
effect appears in the ITO phase and possibly its immediate vicinity. Therefore
in the field-suppressed SDW regime, do not accord any physical significance to
any isomagnetic of $\kappa_{xy}/T$ showing extremum near half-quantization
since the isomagnetics commonly peak at much larger values for nearby fields.

_Strategy_ : Thermal currents can be readily expressed in terms of local
microscopic degrees of freedom (spin operators in our case) and do not require
any knowledge of quasiparticles, which was the approach followed in Kitaev’s
original work Kitaev (2006). Kitaev’s prescription for calculation of edge
thermal currents is however difficult to implement numerically in a general
nonintegrable situation because of the large temperature differences used to
mimic the bulk and vacuum regions. An alternate way of using the Kubo linear
response approach is also numerically challenging because it involves small
differences of two large contributions represented by the energy-current
correlators and the bound energy magnetization contributions Cooper  _et al._
(1997). Our strategy is similar to Kitaev’s except that the system is immersed
in a uniform temperature bath. Although there is no net transport Hall
current, it is because of cancellation of equal and opposite contributions of
the two edges. Then, for gapped systems that are sufficiently long, the
thermal Hall current associated with one edge is simply obtained by summing
over the Hall currents from the relevant edge to the center, where the Hall
current has decayed to very small values.

Figure 1: Cylindrical geometry for the honeycomb Kitaev model used in our
calculations, with periodic boundary conditions along the $y$-direction and
open zigzag edges with $L_{y}$ unit cells. Calculations have been performed
with $N_{y}=2L_{y}=6$ of sites along the $y$ direction and convergence checked
for $N_{y}=8.$ There are up to $N_{x}=L_{x}=20$ sites along the $x$-direction,
the cylinder’s axis. The dashed dual lattice lines separate adjacent layers
(ellipsoids) associated with Hamiltonian $H_{i}$, $H_{i+1}.$

For calculations, we put our model on a cylinder along the $x$-direction with
periodic boundary conditions along the $y$-direction, with a pair of open
zigzag edges, see Fig. 1. The Hamiltonian is re-expressed as a sum of layer
contributions, $H=\sum_{i}H_{i},$ connecting the two open edges. Links that
cut the dual curve (dashed line in Fig. 1) separating adjacent layers, make
equal contributions to the adjacent $H_{i}.$ The energy current across the
dual curve is obtained from the continuity equation,

$\displaystyle\frac{dH_{i}}{dt}\,=\,-i[H_{i},H]\,=\,\sum_{j}\hat{J}_{ij},\quad\hat{J}_{ij}=-i[H_{i},H_{j}].$
(2)

Here $\hat{J}_{ij}$ is the energy current from layer $j$ to $i,$ (see S.M. for
expression of interlayer thermal currents) and we consider it as a sum of
contributions $\hat{J}_{ij}(x)$ starting from one edge ($x=0$) to the other
($x=L_{x}$). In bulk gapped systems, for a sufficiently long cylinder,
$J_{ij}(x)=\langle\hat{J}_{ij}(x)\rangle_{T}$ will decay and ultimately vanish
deep in the bulk. In this situation, $J_{ij}$ is essentially the sum of two
(Hall) currents associated with the respective open ends, and like in Ref.
Kitaev (2006), we are interested in the Hall current associated with one of
the edges, e.g., $J_{ij}^{\text{H}}=\sum_{x\leq L_{x}/2}J_{ij}(x).$ The
averaging includes both quantum and thermal. In the steady state, the labels
$ij$ can also be dropped. Thermal Hall conductivity $\kappa_{xy}$ is obtained
by taking the temperature derivative of $J^{\text{H}}.$

We use the standard purification-based finite temperature Tensor Network
method Hauschild  _et al._ (2018) for performing the quantum and thermal
averaging, where imaginary time evolution of state is carried by applying
evolution operator as a matrix product operator(MPO) with W-II approximation,
described in reference Zaletel  _et al._ (2015) – see S.M. for more details.
Benchmarking of the finite temperature purification method was done (see S.M.)
against Exact Diagonalization (ED) calculation of one of the Kitaev plaquette
fluxes at finite temperature and a small field of $h/K=0.025,$ and we found
very good agreement between the two even at temperatures as low as $T/K=0.01.$

_Results_ : We consider first $\kappa_{xy}/T$ in the Kitaev limit $J/K=0,$ see
Fig. 2. Introduction of a magnetic field hybridizes the dispersing majoranas
with pairs of visons Zhang  _et al._ (2022) whose excitation energy at zero
field is $\Delta_{\text{pair}}\approx 0.065K.$ The hybridization results in
gapping of the bulk majoranas and appearance of chiral Majorana edge modes
Kitaev (2006), a finite dispersion for the static visons, and a suppression of
topological order Hauschild  _et al._ (2018); Zhu  _et al._ (2018) beyond
$h\sim 0.02K$ with a concomitant gapping of the edge states Zhu  _et al._
(2018). Fractionalization may however survive to higher fields, as evidenced
from quasiparticle stability analysis Kumar and Tripathi (2020), and also
experimentally from specific heat measurements Tanaka  _et al._ (2022). At
high fields $h/K\sim 1,$ fractionalization is lost and local spin flips
(magnons) constitute the elementary excitations. We have taken sufficiently
large cylinder lengths such that the edge currents have declined to very small
values in the center of the cylinder (see S.M.) thus ensuring that the left
and right edge currents do not overlap.

Figure 2 shows the temperature dependence of $\kappa_{xy}/T$ in the Kitaev
limit for different values of $h.$ Table 1 additionally shows the temperatures
$T^{*}$ corresponding to the peak values of $\kappa_{xy}/T,$ together with the
calculated bulk gap $\Delta_{\text{bulk}}$ 111The bulk gap is calculated using
DMRG on a torus for a 32-site system.. The $T^{*}$ all lie below
$\Delta_{\text{bulk}},$ signifying the low temperature peaks are associated
with gapped edge modes. For all fields except the lowest $h/K=0.03$ $T^{*}$ is
not very sensitive to the finite edge - we confirmed this by increasing the
edge size from $N_{y}=6$ to $N_{y}=8$ spins, and found no significant change
in the position of $T^{*},$ which is not surprising since edge modes get
gapped upon transitioning out of the ITO phase Zhu  _et al._ (2018). For the
lowest field shown in Fig. 2, $h/K=0.03,$ the peak value of $\kappa_{xy}/T$
reaches nearly $83\%$ of the half-quantized value associated with the ITO
phase, although topological order as measured by topological entanglement
entropy $\gamma$ is significantly degraded, $\gamma(h)/\gamma(0)\approx 0.14,$
(see S.M.). We attribute the reduction here to vison fluctuations as analyzed
in Ref. Guo  _et al._ (2020) since the peak is comparable to the half-
quantized value at lower fields, and the alternate, interacting spin wave
approximation becomes progressively inaccurate at low fields below $h/K\sim 1$
McClarty  _et al._ (2018). Further increase of $h$ monotonously decreases the
peak values, while simultaneously pushing $T^{*}$ to higher temperatures. At
high temperatures, the $\kappa_{xy}/T$ all tend to vanish, and show no signs
of saturation predicted from spin wave McClarty  _et al._ (2018) based
analyses.

Figure 2: Isomagnetic plots of $\kappa_{xy}/T$ for the ferromagnetic Kitaev
model in the field range $0.03\leq h/K\leq 0.3$ interpolating the Ising
topological order (ITO) and spin wave limits. The chosen system size
$N_{x}=10$ was enough for Hall current to decay in bulk. At low fields,
signatures of proximate ITO are visible in the form of $\kappa_{xy}/T$ peaks
that approach half-quantization at the lower end. Thermal Hall effect for this
range of fields appears more consistent with a picture of interacting
fractionalized quasiparticles Guo  _et al._ (2020) than spin waves McClarty
_et al._ (2018). Solid gray line describes the half-quantized thermal Hall
effect.

We now turn on the competing Heisenberg AFM interaction $J.$ At $h=0,$ the
Kitaev spin liquid phase is known to survive up to $J/K\approx 0.1$ Chaloupka
_et al._ (2010); Gohlke  _et al._ (2017), beyond which a stripy SDW phase
appears. The $(111)$ Zeeman field competes with both topological and magnetic
order Jiang  _et al._ (2011), and at large values yields a trivial polarized
phase. Recent field-theoretical phenomenology Zou and He (2020) as well as
thermal Hall measurements of $\alpha$-RuCl3 support the view that ITO re-
emerges upon field-suppression of SDW order. However other experimental Czajka
_et al._ (2021) and theoretical studies based on spin wave treatments McClarty
_et al._ (2018); Koyama and Nasu (2021) disagree, and even the robustness of
the reported half-quantization in this regime has been questioned. Before
examining $\kappa_{xy}/T,$ which depends on the excitations in the model, we
first checked if the ground state state shows re-emerged ITO by calculating
$\gamma(h)$ for different values of $J/K$ (see S.M.). For the magnetically
ordered phase, $J/K>0.1,$ we found that there is indeed a small partial
revival of $\gamma(h)$ from very small values to near $\gamma(h)$ associated
with the the pure Kitaev model - this happens for fields $h\sim J/3.$ However
the revival never takes $\gamma$ anywhere near the value $\ln(2)$ associated
with ITO (see S.M.). In the discussion below we focus our attention in the
vicinity of the fields where this partial revival of ITO takes place, since
$\kappa_{xy}/T$ also peaks here.

Figure 3 shows isomagnetic curves of $\kappa_{xy}/T$ for different values of
$J/K$ across the Kitaev spin liquid to stripy AFM transition. Consider the
curve corresponding to $(J/K,h/K)=(0.1,0.05).$ Apart from the sign reversal,
$\kappa_{xy}/T$ shows a deep minimum that exceeds the half-quantized value,
indicating we are well outside the ITO regime at such fields. This is also
supported by the strong degradation of $\gamma(h)$ (see S.M.) for the above
parameters. We noticed that sign reversal of $\kappa_{xy}/T$ happens at
$J/K=0.02$, observed from the computation of Hall current in ground state,
which is a zero temperature calculation. Higher values of $J/K=0.2,\,0.3$ take
us well within the stripy SDW phase at $h=0.$ We dial up the field to values
where ITO is partially revived, and calculate $\kappa_{xy}/T$ in its vicinity.
We find deep minima of comparable strength to the one seen for
$(J/K,h/K)=(0.1,0.05).$ Upon increasing the field, the minimum becomes
shallower, very reminiscent of the intermediate field behavior of the Kitaev
model where ITO has degraded. The temperature corresponding to the extremum,
$T^{*},$ lies well within the bulk gap (see Table 1) for these parameter
values, and it is not surprising that the edge currents are significant at
these temperatures. In regimes with strong magnetic order and low fields
($h/J\ll 1$), we found that the Hall currents do not decay in the bulk for the
length-scales we could study; consequently restricting us to the field-
suppressed SDW regime and beyond. However it is encouraging for us that
experimentally, the thermal Hall response at low fields in the SDW phase is
not observed to be significant Kasahara  _et al._ (2018), and this is also
observed in spin wave treatments Koyama and Nasu (2021) that are relevant for
magnetically ordered states.

Figure 3:  Isomagnetic plots of $\kappa_{xy}/T$ for the $J$-$K$ model, with Heisenberg interactions chosen both inside the Kitaev spin liquid phase ($J/K=0.1$) as well as in the SDW phase ($J/K=0.2,\,0.3$). For $J/K=0.1,$ the fields chosen exceed the region of ITO order, and for the SDW, the fields are in the vicinity of the partially revived ITO and beyond. Generically, $\kappa_{xy}/T$ is characterized by deep troughs at low fields that get shallower as the field increases. The minimum values in the vicinity of partially revived ITO typically go beyond half-quantization. For the field-suppressed SDW cases, we found no phase characterized by half-quantization and accompanying fully revived ITO. $(J/K,h/K)$ | $T^{*}/K$ | $\Delta_{\text{bulk}}/K$
---|---|---
(0,0.03) | 0.0149 | 0.0473
(0, 0.035) | 0.0182 | 0.0494
(0, 0.06) | 0.0278 | 0.0720
(0, 0.10) | 0.0482 | 0.1160
(0, 0.15) | 0.0725 | 0.1860
(0, 0.5) | 0.2101 | 0.77
(0.1,0.05) | 0.015 | 0.0529
(0.1,0.1) | 0.036 | 0.1039
(0.2,0.1) | 0.022 | 0.0738
(0.2,0.15) | 0.044 | 0.14
(0.3,0.15) | 0.0242 | 0.080
Table 1: Table showing the temperatures corresponding to the extrema $T^{*}$
of the isomagnetic plots of $\kappa_{xy}/T$ in Figs. 2 and 3, together with
the bulk gap $\Delta_{\text{bulk}}.$

To summarize, we studied the temperature dependence of thermal Hall
conductivity of the $J$-$K$ Kitaev-Heisenberg model for Zeeman fields in the
$(111)\parallel\mathbf{c}$ direction using a purification based finite
temperature Tensor Network technique that goes down to $K/T\approx 100.$ We
make no prior assumption about the nature of the quasiparticles, whether
fractionalized ones or spin waves. In the Kitaev limit ($J/K$=0,) we focused
on fields $h>0.02K$ where ITO has substantially degraded, extending the
understanding of thermal Hall conductivity to the topologically trivial phase
in the neighborhood of ITO. We found that $\kappa_{xy}/T$ features a low-
temperature peak whose value approaches half-quantization at the lower end of
our field range, and continuously declines with increasing $h,$ approaching
zero at larger fields. We associate the low-temperature peak with gapped edge
modes expected in the absence of Ising topological order. After peaking,
$\kappa_{xy}/T$ monotonously decreases with increasing temperature as
predicted for vison fluctuation effects Guo  _et al._ (2020), but unlike the
prediction from earlier works, shows no sign of saturation at high
temperatures. To understand the effect of Heisenberg perturbations we
considered the cases $J/K=0.1$ within the spin liquid phase (but close to the
SDW phase boundary), and $J/K=0.2,\,0.3$ that lie well within the stripy SDW
phase. Characteristically deep minima were seen in the isomagnetic curves of
$\kappa_{xy}/T$ versus $T,$ significantly exceeding half-quantization of the
ITO phase, reminiscent of the experimental observations in
$\alpha$-RuCl${}_{3}.$ Although many of our observations are consistent with
vison fluctuation effects, the large extrema in $\kappa_{xy}/T$ in the
vicinity of field-suppressed magnetic order are not. At higher fields, any
partially revived topological order seen after destruction of SDW order gets
suppressed again, resulting in the polarized phase with small thermal Hall
effects.

We did not find any intermediate field regime following field-suppression of
SDW order where robust half-quantized $\kappa_{xy}/T$ associated with ITO
exists. Although in the vicinity of the field-induced SDW suppression we do
note a weak partial revival of topological order (as measured by the
topological entanglement entropy $\gamma$), the revived $\gamma$ remain well
below $\gamma=\ln 2$ associated with ITO. Since the magnitude of
$\kappa_{xy}/T$ in the $J$-$K$ model may exceed the half-quantized value, it
is possible to fine tune the field to bring it near half-quantization but that
does not imply ITO. We believe that similar conclusions should also hold near
field-suppressed SDW order even in the presence of other competing
interactions like the anisotropic terms $\Gamma,\Gamma^{\prime}$ \- this is
ultimately to be decided in future studies.

###### Acknowledgements.

AK and VT acknowledge support of the Department of Atomic Energy, Government
of India, under Project Identification No. RTI 4002, and the Department of
Theoretical Physics, TIFR, for computational resources. AK thanks helpful
discussions with Sounak Biswas about finite temperature tensor network
methods.

##  Appendix A Purification method and benchmarking

For calculation of finite temperature quantities, we have used a purification
based tensor network method Hauschild  _et al._ (2018). The purification
strategy involves expressing the mixed state as a pure state in an enlarged
Hilbert space. We obtain the finite temperature (T=1/$\beta$) thermofield
state (TFD) $|\psi_{\beta}\rangle$ by action of imaginary time evolution
operator $\exp(-\frac{\beta}{2}H)$ on the infinite temperature state
$|\psi_{0}\rangle$, which is constructed by Bell entangled pair of physical
and ancilla spin $1/2$ degree of freedom at each sites. Here, $H$ acts on the
physical part of Hilbert space and we use the compact matrix product operator
(MPO) representation for the finite temperature evolution operator with W-II
approximation Zaletel  _et al._ (2015), provided in TenPy library Hauschild
and Pollmann (2018). Finite temperature expectation of a physical observable
$\hat{O}$ is obtained from
$\langle\hat{O}\rangle_{\beta}=\langle\psi_{\beta}|\hat{O}|\psi_{\beta}\rangle.$
The purification method works best at higher temperatures. Figure 4 shows the
benchmarking of the method against Exact Diagonalization (ED) calculation of
one of the Kitaev plaquette fluxes $\langle W_{p}\rangle$ at finite
temperatures and a small field of $h/K=0.025$ for N= 18 site system. We find
very good agreement between the two even at temperatures as low as $T/K=0.01.$
At higher temperatures, the discrepancy between purification and ED results is
actually on account of the latter being inadequate (insufficient number of
excited states taken into consideration).

Figure 4: Benchmarking done via calculation of flux value for an $18$-spin
system with periodic boundary conditions using ED and MPS based purification
methods. The bond dimension $\chi=400$ in the MPS calculation is adequate for
good agreement with ED at the low temperature end. At higher temperatures, the
small discrepancy is actually on account of ED - the number of excited states
used is insufficient for correct description of the high temperature behavior.

The accuracy of matrix product state (MPS) calculations requires enhancing the
bond dimension $\chi$ and reducing the Trotter step size $\Delta\beta$ used in
W-II approximation Zaletel  _et al._ (2015). Figure 5 shows the
$\kappa_{xy}/T$ (in the Kitaev limit as well as for the $J$-$K$ model) as a
function of the Trotter step size. Calculations in the Kitaev limit typically
require larger $\chi$ and smaller $\Delta\beta$ compared to the Heisenberg-
perturbed cases that we studied. It is evident from Fig. 5 that fine-graining
the Trotter steps from $0.04$ to $0.02$ does not affect the expectation
values, while in the Kitaev limit we still see some sensitivity, especially at
lower field values. We observed that for the bond dimensions that we work with
($\chi=400$ and $600$), Trotter errors are the main limitation on accuracy.
For this reason, we limit our calculations in the Kitaev limit to fields
greater than $h/K=0.03.$

Figure 5:  Plot (a) shows the $\kappa_{xy}/T$ vs. T/K with varying
($\Delta\beta$) and $\chi$ in pure Kitaev limit for different magnetic field
parameter range. $\kappa_{xy}/T$ is almost converged at high $\chi$ and lower
($\Delta\beta$). Plot (b) shows $\kappa_{xy}/T$ vs. T/K for J-K Kitaev-
Heisenberg model with $J/K=0.1$ for two magnetic field parameters
$h/K=0.05,0.075$ and they found to be converged, implying accuracy of
computation.

Figure 6: Plot (a) shows the decay of thermal Hall current along the the
length of cylinder which is measured at the peak of $\kappa_{xy}/T$ in pure
Kitaev limit with varying magnetic field strength. We show the current decay
up to half of the chosen cylinder length $(L_{x}=10)$ and in other half
current will be equal and opposite. It is evident from the plot that thermal
current vanishes deep in the bulk. Plot (b) show the current decay in field
suppressed SDW phase. Here, we choose the cylinder length $(L_{x}=20)$ to make
sure current decays well with in bulk gap.

##  Appendix B Thermal Hall current decay in the bulk

In gapped phases, the low-temperature behavior of thermal Hall effect is
governed by low-lying edge modes. Figure 6 shows the decay of the thermal Hall
current along the $x$-direction (cylinder axis) for in the Kitaev limit for
various values of the Zeeman field along the $(111)$ direction. The current
decays rapidly as we go deeper in the bulk allowing us to treat the left and
right edge currents independently. Figure 6 shows the current decay in the
$J$-$K$ Kitaev-Heisenberg model in the presence of a Zeeman field along the
$(111)$ direction. The ratios $J/K$ span both sides of the spin liquid
($J/K=0.1$) and stripy SDW phase ($J/K=0.2,\,0.3,$) boundary, and in the SDW
phase, the field is chosen to just suppress SDW order and partially revive
topological order (as measured by the topological entanglement entropy
$\gamma$, see Fig. 8). We observe that the current changes its sign relative
to the Kitaev limit and also decays rapidly into the bulk.

Figure 7: Ground state plaquette flux revival upon suppression of magnetic
order. We considered the system size of $N_{x}=10,N_{y}=8,$ where ground state
is computed using finite DMRG method. Figure 8: Field dependence of
topological entanglement entropy $\gamma.$ It is small in the magnetically
ordered phase as well as at high fields but there is a small, partial revival
at intermediate fields where magnetic field is comparable to magnetic ordering
scale.

##  Appendix C Flux revival and Topological entanglement entropy ($\gamma$) in
field suppressed SDW phase

Motivated by seen experimental large thermal Hall in field revived SDW phase,
we compute flux expectation value in ground state for J-K model with varying
Zeeman field strength. We see that flux gets revived at Zeeman field strength
$\sim J/3$ \- see figure 7. We also compute topological entanglement entropy
$\gamma$ using Kitaev and Preskill construction Kitaev and Preskill (2006). In
this intermediate field regime, we see that the $\gamma(h)$ gets revived,
albeit weakly, see Figure 8. We carried out this computation with system size
of $N_{x}=10,N_{y}=8$ using finite DMRG method.

## References

  * (1) G. Grissonnanche, S. Thériault, A. Gourgout, M. Boulanger, E. Lefrançois, A. Ataei, F. Laliberté, M. Dion, J. Zhou, S. Pyon, _et al._ , arXiv preprint arXiv:2003.00111 .
  * Hirschberger  _et al._ (2015a) M. Hirschberger, J. W. Krizan, R. Cava, and N. Ong, arXiv preprint arXiv:1502.02006 (2015a).
  * Hirschberger  _et al._ (2015b) M. Hirschberger, R. Chisnell, Y. S. Lee, and N. P. Ong, Phys. Rev. Lett. 115, 106603 (2015b).
  * Kasahara  _et al._ (2018) Y. Kasahara, K. Sugii, T. Ohnishi, M. Shimozawa, M. Yamashita, N. Kurita, H. Tanaka, J. Nasu, Y. Motome, T. Shibauchi, and Y. Matsuda, Phys. Rev. Lett. 120, 217205 (2018).
  * Yokoi  _et al._ (2021) T. Yokoi, S. Ma, Y. Kasahara, S. Kasahara, T. Shibauchi, N. Kurita, H. Tanaka, J. Nasu, Y. Motome, C. Hickey, _et al._ , Science 373, 568 (2021).
  * Bruin  _et al._ (2021) J. Bruin, R. Claus, Y. Matsumoto, N. Kurita, H. Tanaka, and H. Takagi, arXiv preprint arXiv:2104.12184 (2021).
  * Czajka  _et al._ (2021) P. Czajka, T. Gao, M. Hirschberger, P. Lampen-Kelley, A. Banerjee, J. Yan, D. G. Mandrus, S. E. Nagler, and N. Ong, Nature Physics 17, 915 (2021).
  * Kitaev (2006) A. Kitaev, Annals of Physics 321, 2 (2006).
  * Chern  _et al._ (2021) L. E. Chern, E. Z. Zhang, and Y. B. Kim, Phys. Rev. Lett. 126, 147201 (2021).
  * Zhang  _et al._ (2021) E. Z. Zhang, L. E. Chern, and Y. B. Kim, Phys. Rev. B 103, 174402 (2021).
  * Koyama and Nasu (2021) S. Koyama and J. Nasu, Phys. Rev. B 104, 075121 (2021).
  * McClarty  _et al._ (2018) P. A. McClarty, X.-Y. Dong, M. Gohlke, J. G. Rau, F. Pollmann, R. Moessner, and K. Penc, Phys. Rev. B 98, 060404 (2018).
  * Joshi (2018) D. G. Joshi, Phys. Rev. B 98, 060405 (2018).
  * Guo  _et al._ (2020) H. Guo, R. Samajdar, M. S. Scheurer, and S. Sachdev, Physical Review B 101, 195126 (2020).
  * Nasu  _et al._ (2017) J. Nasu, J. Yoshitake, and Y. Motome, Phys. Rev. Lett. 119, 127204 (2017).
  * Hauschild  _et al._ (2018) J. Hauschild, E. Leviatan, J. H. Bardarson, E. Altman, M. P. Zaletel, and F. Pollmann, Phys. Rev. B 98, 235163 (2018).
  * Zhu  _et al._ (2018) Z. Zhu, I. Kimchi, D. N. Sheng, and L. Fu, Phys. Rev. B 97, 241110 (2018).
  * Gohlke  _et al._ (2018) M. Gohlke, R. Moessner, and F. Pollmann, Phys. Rev. B 98, 014418 (2018).
  * Cooper  _et al._ (1997) N. R. Cooper, B. I. Halperin, and I. M. Ruzin, Phys. Rev. B 55, 2344 (1997).
  * Zaletel  _et al._ (2015) M. P. Zaletel, R. S. K. Mong, C. Karrasch, J. E. Moore, and F. Pollmann, Phys. Rev. B 91, 165112 (2015).
  * Zhang  _et al._ (2022) S.-S. Zhang, G. B. Halász, and C. D. Batista, Nature communications 13, 1 (2022).
  * Kumar and Tripathi (2020) A. Kumar and V. Tripathi, Phys. Rev. B 102, 100401 (2020).
  * Tanaka  _et al._ (2022) O. Tanaka, Y. Mizukami, R. Harasawa, K. Hashimoto, K. Hwang, N. Kurita, H. Tanaka, S. Fujimoto, Y. Matsuda, E.-G. Moon, _et al._ , Nature Physics 18, 429 (2022).
  * Note (1) The bulk gap is calculated using DMRG on a torus for a 32-site system.
  * Chaloupka  _et al._ (2010) J. c. v. Chaloupka, G. Jackeli, and G. Khaliullin, Phys. Rev. Lett. 105, 027204 (2010).
  * Gohlke  _et al._ (2017) M. Gohlke, R. Verresen, R. Moessner, and F. Pollmann, Phys. Rev. Lett. 119, 157203 (2017).
  * Jiang  _et al._ (2011) H.-C. Jiang, Z.-C. Gu, X.-L. Qi, and S. Trebst, Phys. Rev. B 83, 245104 (2011).
  * Zou and He (2020) L. Zou and Y.-C. He, Phys. Rev. Research 2, 013072 (2020).
  * Hauschild and Pollmann (2018) J. Hauschild and F. Pollmann, SciPost Phys. Lect. Notes , 5 (2018).
  * Kitaev and Preskill (2006) A. Kitaev and J. Preskill, Phys. Rev. Lett. 96, 110404 (2006).