# Topological Green’s function zeros in an exactly solved model and beyond

Steffen Bollmann  Max-Planck Institute for Solid State Research, 70569
Stuttgart, Germany  Chandan Setty  Department of Physics and Astronomy, Rice
Center for Quantum Materials, Rice University, Houston, Texas 77005, USA
Urban F. P. Seifert  Kavli Institute for Theoretical Physics, University of
California, Santa Barbara, CA 93106, USA  Elio J. König  Max-Planck Institute
for Solid State Research, 70569 Stuttgart, Germany

###### Abstract

The interplay of topological electronic band structures and strong
interparticle interactions provides a promising path towards the constructive
design of robust, long-range entangled many-body systems. As a prototype for
such systems, we here study an exactly integrable, local model for a
fractionalized topological insulator. Using a controlled perturbation theory
about this limit, we demonstrate the existence of topological bands of zeros
in the exact fermionic Green’s function and show that in this model they do
affect the topological invariant of the system, but not the quantized
transport response. Close to (but prior to) the Higgs transition signaling the
breakdown of fractionalization, the topological bands of zeros acquire a
finite “lifetime”. We also discuss the appearance of edge states and edge
zeros at real space domain walls separating different phases of the system.
This model provides a fertile ground for controlled studies of the
phenomenology of Green’s function zeros and the underlying exactly solvable
lattice gauge theory illustrates the synergetic cross-pollination between
solid-state theory, high-energy physics and quantum information science.

The spectrum of low-lying excitations is central to quantum-many body physics
and reflected in poles of the electronic Green’s function (GF) [1].
Interpreted as a propagator, the GF encodes quantum mechanical transition
amplitudes. From this perspective, it is a natural and tantalizing question to
also study nodes of the amplitudes indicating destructive many-body
interference. This results in GF zeros which are key to Fermi surface
reconstruction and pseudogap physics — a major open problem in quantum
materials research [2, 3, 4, 5, 6, 7, 8, 9]. Just like their counterpart, the
poles, GF zeros defined by
$\mathbf{G}^{R}(\epsilon_{\emptyset}(\mathbf{k}),\mathbf{k})\psi(\mathbf{k})=0$
trace a bandstructure $\epsilon_{\emptyset}(\mathbf{k})$ in the Brillouin zone
(BZ). The associated “Bloch”-like eigenstates $\psi(\mathbf{k})$ may be
topological [10, 11, 12, 13, 14, 15] leading to, e.g., topological boundary
zeros [16]. Beyond numerics, analytically exact results for Hatsugai-Kohmoto-
type models with ultra-non-local interactions[17, 18, 19, 20, 21] demonstrate
non-trivial band representations of poles and zeros, alike[11, 13, 14, 22].

Figure 1: a) Schematics of the model. b) Phase diagram as a function of $h/J$
and $w_{2}/w_{0}$. Acronyms: I = insulator, TI = topological insulator, OTI =
orthogonal TI, ZTI = zeros TI, SM = semimetal, OSM = orthogonal SM. In the OTI
phase, the poles of the GF are topological, but the zeros of the GF are
trivial. In the ZTI phase, both are topological. The shaded region close to
the horizontal confinement-deconfinement transition at $h/J=3$ indicates that
GF zeros acquire a “lifetime”.

In this letter, we study the topological properties of Bloch-bands of zeros by
investigating a local many-body Hamiltonian which is exactly soluble, and in
which fermions couple to an emergent $\mathbb{Z}_{2}$ gauge theory [23, 24,
25, 26, 27, 4, 5, 28, 29, 30, 31]. In this theory, correlators of physical
fermions are given by convolutions of propagators of fractionalized
excitations. Equipped with this fundamental structure we can exactly prove in
the soluble limit (i.e. with a static $\mathbb{Z}_{2}$ gauge field) that the
physical GFs possess zeros. This provides us with a starting point for a
controlled perturbative expansion about the integrable limit, in which we
demonstrate that dispersing GF zeros form a topologically non-trivial
bandstructure, but that their existence is not directly related to quantized
topological transport coefficients. An infinite order resummation of leading
Feynman diagrams allows to access the exact long-wavelength GF in the vicinity
of the deconfinement-confinement Higgs transition. As the transition is
approached, the mass of GF zeros increases and hits the continuum of
excitations (i.e. poles) in the system. Beyond this point, the energy of GF
zeros becomes complex.

A GF based topological invariant [32, 33, 34, 10, 35, 36] may be employed to
attribute topology to strongly interacting systems, e.g., for Altland-
Zirnbauer class A in two dimensions

$N_{3}[\mathbf{G}]=\int\frac{\text{d}\epsilon\text{d}^{2}k}{24\pi^{2}}\epsilon^{\mu\nu\rho}\tr[\mathcal{A}_{\mu}\mathcal{A}_{\nu}\mathcal{A}_{\rho}],$
(1)

where
$\mathcal{A}_{\mu}=\mathbf{G}^{-1}(i\epsilon,\mathbf{k})\partial_{K_{\mu}}\mathbf{G}(i\epsilon,\mathbf{k})$.
In the non-interacting limit $N_{3}$ is identified with the quantized
anomalous Hall conductance $\sigma_{H}$, but the equivalence may break down in
the presence of strong interactions [37, 12, 38, 39]. Note that GF zeros enter
$N_{3}$ analogously to poles [10] and topological transitions in $N_{3}$ may
occur through gap closing of either poles or zeros. Related to this, the free
Dirac fermion at the non-interacting topological transition may experience
symmetric mass generation [40, 41, 42, 43] at strong coupling to account for
the topological transition in the GF invariant without requiring a gap closure
in the physical fermionic spectrum.

Model. In this letter, we consider a system in which Haldane’s model for the
quantum Hall effect without Landau levels is coupled to the toric code (TC) on
a honeycomb lattice, cf. Fig. 1 a), $H=H_{\rm QAH}+H_{\rm TC}+H_{J}$ where
[44, 25]

$\displaystyle H_{\rm QAH}=$ $\displaystyle-
w_{0}\sum_{\mathbf{x}}(-1)^{\mathbf{x}}c^{\dagger}_{\mathbf{x}}c_{\mathbf{x}}-w_{1}\sum_{\langle\mathbf{x},\mathbf{x}^{\prime}\rangle}c^{\dagger}_{\mathbf{x}}Z_{\langle\mathbf{x},\mathbf{x}^{\prime}\rangle}c_{\mathbf{x}^{\prime}}$
$\displaystyle-
iw_{2}\sum_{\llangle\mathbf{x},\mathbf{x}^{\prime}\rrangle}(-1)^{\mathbf{x}}[c^{\dagger}_{\mathbf{x}}\prod_{\mathbf{b}\in\gamma_{\mathbf{x},\mathbf{x}^{\prime}}}Z_{\mathbf{b}}c_{\mathbf{x}^{\prime}}-H.c.]$
$\displaystyle H_{\rm TC}=$ $\displaystyle-
K_{\varhexagon}\sum_{\varhexagon}B_{\varhexagon}-h\sum_{\mathbf{x}}Q_{\mathbf{x}},$
(2a) $\displaystyle H_{J}=$ $\displaystyle-J\sum_{\mathbf{b}}Z_{\mathbf{b}}.$
(2b)

Throughout the letter, $X_{\mathbf{b}},Z_{\mathbf{b}}$ denote Pauli matrices
residing on bonds and
$B_{\varhexagon}=\prod_{\mathbf{b}\in\partial\varhexagon}Z_{\mathbf{b}},Q_{\mathbf{x}}=(-1)^{\hat{n}_{\mathbf{x}}}\prod_{\mathbf{b}\in\partial\mathbf{x}}X_{\mathbf{b}}$
are plaquette and star operators respectively. We consider $h>0$ and, given
that $m-$particles are entirely irrelevant for the present discussion, send
$K_{\varhexagon}\rightarrow\infty$. Note that we added a staggered sublattice
potential $w_{0}$, which contrary to the Haldane next-nearest neighbor hopping
$w_{2}$ opens a trivial gap near the K, K’ points. Both $w_{0},w_{2}$ have
opposite contributions in opposite sublattices, as indicated by the notation
$(-1)^{\mathbf{x}}$. The path $\gamma_{\mathbf{x},\mathbf{x}^{\prime}}$ is the
shortest connection of next-nearest neighbors. Planck’s constant and the speed
of light are set to unity throughout this letter.

Exact solutions. This model is exactly soluble [25] in two limits. First, when
$h/J=0$, all qubits of the TC are frozen and the ground state is

$\ket{GS}_{h/J=0}=\ket{+}_{TC}\otimes\ket{FS}_{c}.$ (2ca) The state
$\ket{+}_{TC}$ indicates that $\langle
Z_{\mathbf{b}}\rangle=1\forall\mathbf{b}$ and $\ket{FS}_{c}$ is the filled
Fermi sea at half filling of the free fermion model defined by the Bloch
Hamiltonian in sublattice space
$h(\mathbf{k})=\left(\begin{array}[]{cc}w_{0}+w_{2}m(\mathbf{k})&-w_{1}s(\mathbf{k})\\\
-w_{1}s^{*}(\mathbf{k})&-w_{0}-w_{2}m(\mathbf{k})\end{array}\right),$ (2cb)
where
$s(\mathbf{k})=1+2e^{-i{\sqrt{3}k_{y}}/{2}}\cos\left({k_{x}}/{2}\right),m(\mathbf{k})=4\sin\left({k_{x}}/{2}\right)[\cos\left({k_{x}}/{2}\right)-\cos\left({\sqrt{3}k_{y}}/{2}\right)]$.
The Hamiltonian Eq. (2cb) is topological (topologically trivial) for
$w_{2}/w_{0}>1/\sqrt{27}$ ($w_{2}/w_{0}<1/\sqrt{27}$). We denote this phase as
TI (I) in Fig. 1 b), it clearly displays $\sigma_{H}=N_{3}=1$
($\sigma_{H}=N_{3}=0$) and an edge state at the TI-I interface, Fig. 2 d).

By virtue of the fact that all $Q_{\mathbf{x}}$ and $B_{\varhexagon}$ are
integrals of motion when $J=0$, Eq. (2) is also exactly soluble in the limit
$h/J\rightarrow\infty$, in this case

$\ket{GS}_{h/J\rightarrow\infty}=\prod_{\mathbf{x}}\frac{1+Q_{\mathbf{x}}}{2}\Big{(}\ket{+}_{TC}\otimes\ket{FS}_{c}\Big{)}.$
(2ce)

It is important to emphasize that Eq. (2ca) is the direct product of a product
state in the TC sector and a Slater determinant of a fully filled band.
Therefore it is short range entangled and strongly distinct from (2ce) which
involves topologically order stemming from the TC sector.

Figure 2: a) Real space representation of the GF, Eq. (2d), and the current-
current correlator. In the deconfining phase, the gapped $e$-particles (red,
wiggly propagator) imply an exponentially decaying fermionic GF even when the
propagator of orthogonal fermions is gapless. However, $e$-particle
propagators do not enter the current-current correlator (pink). b) In the
limit $J=0$ the spectral function is exactly momentum independent (color plot
for $h=w_{1}$ and $w_{0}=w_{2}=0$ along a BZ cut). Perturbative inclusion of
$J$ leads to moderate dispersion of Hubbard bands (not shown), while GF zeros
acquire a graphene-like dispersion of bandwidth $\mathcal{O}(Jw_{1}/h)$ (red).
c) The additional inclusion of the staggered potential leads to a trivial gap
$\mathcal{O}(w_{0})$ in the dispersion of zeros (top circle). In contrast,
Haldane-like next-nearest neighbor hopping opens a topological gap
$\mathcal{O}(J^{2}w_{2}/h^{2})$ in the Dirac nodes of the zeros (bottom
circle). At $w_{0}\sim J^{2}w_{2}/h^{2}$ the gap vanishes in one of the
K/K’-valleys signaling a topological phase transition of GF zeros (central
circle). d) Real space diagram of edge states/zeros with color coding of
phases as in Fig. 1 b).

Propagators of fractionalized theory: $h\gg J$. The TC coupled to fermions,
Eq. (2), is related to $\mathbb{Z}_{2}$ slave spin theories [45, 46] of
Hubbard-like models by gauge fixing [47]. In this context $h$ is the local
Hubbard repulsion, while $w_{1,2}$ and $J$ correspond to the hopping
amplitudes of fractionalized fermionic and slave spin excitations,
respectively. In particular, the gapped $e$-particles of the TC are given by
flipped star operators (i.e. sites at which $Q_{\mathbf{x}}=-1$ instead of
$Q_{\mathbf{x}}=1$). In the integrable limit $J=0$ the Euclidean time GF
$D(\tau;\mathbf{x},\mathbf{x}^{\prime})=e^{-2h|\tau|}\delta_{\mathbf{x},\mathbf{x}^{\prime}}$
of $e$-particles is ultralocal, this can be interpreted as a vanishing kinetic
energy of slave-spin-flips. The fermionic excitations on top of the filled
Fermi sea $\ket{FS}_{c}$ with propagator
$G(i\epsilon,\mathbf{k})=[i\epsilon-h(\mathbf{k})]^{-1}$ will be referred to
as “orthogonal fermions”.

As the creation of the physical fermion in the fractionalized theory involves
the simultaneous excitation in the fermionic and TC (slave spin) sectors
alike, the total two-point function takes the form of a convolution, Fig. 2
a):

$\mathbf{G}(\tau;\mathbf{x},\mathbf{x}^{\prime})\equiv-\langle\mathcal{T}[c_{\mathbf{x}}(\tau)c_{\mathbf{x}^{\prime}}^{\dagger}(0)]\rangle=D(\tau;\mathbf{x},\mathbf{x}^{\prime})G(\tau;\mathbf{x},\mathbf{x}^{\prime}).$
(2d)

We now go away from the solvable limit $J=0$ and use perturbation theory in
$J/h$ to obtain $D$ and $\mathbf{G}$. The perturbative evaluation [48] of the
fermionic propagator up to second order in $J/h$ is conveniently written in
the momentum- and frequency domain as
$\mathbf{G}(i\epsilon,\mathbf{k})=-\left[i\epsilon+\mathcal{C}_{0}h_{\emptyset}(i\epsilon,\mathbf{k})h/w_{1}\right]/{4h^{2}},$
where

$\displaystyle h_{\emptyset}(i\epsilon,\mathbf{k})$
$\displaystyle=\left(\begin{array}[]{cc}w_{0}+w_{2}\frac{J^{2}}{h^{2}}\mathcal{C}_{2}m(\mathbf{k})&-w_{1}\frac{J}{h}\mathcal{C}_{1}s(\mathbf{k})\\\
-w_{1}\frac{J}{h}\mathcal{C}_{1}s^{*}(\mathbf{k})&-w_{0}-w_{2}\frac{J^{2}}{h^{2}}\mathcal{C}_{2}m(\mathbf{k})\end{array}\right),$
(2g)

and the matrix form results from the sublattice-resolved structure. In the
limit $w_{0}\ll w_{2}\ll w_{1}\ll h,\epsilon\ll h$, the functions
$\mathcal{C}_{0,1,2}$ acquire a finite, positive value of order unity.

After analytical continuation to real energies, the following observations are
noteworthy: First, for $J=0$, the ultralocal, gapped behavior of
$D(\tau;\mathbf{x},\mathbf{x}^{\prime})$ implies a $\mathbf{k}$ independent
spectral weight, Fig. 2 b), with all spectral weight pushed to “Hubbard” bands
at energies larger than $2h$. Second, inside the gap there are flat bands of
zeros at $E=\pm w_{0}$ (not shown)[29]. One may readily check that in the
limit $h/J\rightarrow\infty$, $N_{3}=0$. While this behavior persists
throughout, for $w_{2}/w_{0}<1/\sqrt{27}$ ($w_{2}/w_{0}>1/\sqrt{27}$) i.e.
when $h(\mathbf{k})$ is topologically trivial (topological), we use the name
orthogonal (topological) insulator, abbreviated as OI (OTI) in Fig. 1 b).
Third, perturbation theory in $J$ leads to slow $\mathbf{k}$ dependence of
“Hubbard” bands (not shown) and importantly lifts the degeneracy of the flat
bands of zeros, see red contours in Fig. 2 b). Note that nearest neighbor
hopping, $s(\mathbf{k})$, enters $\mathbf{G}(i\epsilon,\mathbf{k})$ to first
order, while the Haldane-like next-nearest neighbor hopping, $m(\mathbf{k})$,
appears in the GF only at the second order in $J/h$. Fourth, by comparing
Haldane mass and staggered potential of $h_{\emptyset}$, the band structure of
zeros is topological for $h/J\lesssim\sqrt{w_{2}/w_{0}}$, Fig. 2 c). We call
this phase “zeros topological insulator” (ZTI) and explicitly demonstrate
below that $N_{3}=1$ in the ZTI, just as in the TI phase. This equality
follows from $h_{\emptyset}(\mathbf{k})$ entering the numerator of
$\mathbf{G}(i\epsilon,\mathbf{k})$ with the opposite sign than $h(\mathbf{k})$
in the denominator of $G(i\epsilon,\mathbf{k})$. This also implies that the
Bloch wave functions $\psi(\mathbf{k})$ associated to zeros are not
wannierizable, i.e. do not correspond to elementary band representations.
Sixth, the topological phase transition between ZTI and OTI involves a Dirac-
like gap-closing of zeros, central panel of Fig. 2 c), which is depicted by a
dot-dashed orange line in Fig. 1 b). Finally, in the limit of slowly spatially
varying $w_{2}(\mathbf{x})$, semiclassical arguments imply boundary zeros in
real space at real space domain walls between OTI and ZTI, Fig. 2 d), bottom
panel.

Exact continuum GF for $h\gg\sqrt{h^{2}-3hJ}\gg J^{3}/h^{2}$. The inclusion of
$J$-induced $e$-particle hopping can be resummed to yield dispersive
$e$-particle bands, of which the lower can be approximated by
$\omega(\mathbf{q})=\sqrt{r^{2}+v^{2}q^{2}}$ at lowest energies $r\ll h$, Fig.
3 a). Here $v=\sqrt{hJ}$, $r^{2}=4h[h-3J]$. The closing of the gap indicates
condensation of $e$-particles, i.e. the Higgs transition. Below the
transition, the fermionic GF, Eq. (2d), acquires a finite quasiparticle weight
[4, 5].

In the present case of 2+1 dimension, the immediate vicinity of the transition
is governed by non-Gaussian fluctuations ($e$-particle self-interaction), but
extending standard arguments by Ginzburg and Levanyuk [49, 50, 48], one may
neglect those as long as $r\gg J^{3}/h^{2}$ [48, 5]. Moreover, in the present
model, the interaction between fermion and $e$-particles is absent to leading
order [5]. Such terms might be generated due to fluctuations near the
confinement transition, but are renormalization group irrelevent in view of
the fermionic gap and therefore disregarded throughout.

For excitations close to the K point we use a continuum Hamiltonian for
fermions, too,

$h_{\rm eff}(\mathbf{p})=c\mathbf{p}\cdot\bm{\sigma}+m\sigma_{z}.$ (2h)

Here, $m=w_{0}+\sqrt{27}w_{2},c=\sqrt{3}w_{1}/2$, and $\mathbf{p}$ is the
momentum relative to the K-point, $\bm{\sigma}=(\sigma_{x},\sigma_{y})$, where
$\sigma_{x,y,z}$ are Pauli matrices in sublattice space. An analogous equation
holds close to the K’ point, but there the mass is
$\bar{m}=w_{0}-\sqrt{27}w_{2}$.

Using the effective continuum theory with the above assumption, the convoluted
GF, Eq. (2d), can be calculated exactly in the simplified case $c=v$ of
identical speed for $e$-particles and Dirac fermions

$\displaystyle\mathbf{G}(i\epsilon,\mathbf{p})$
$\displaystyle=-{\frac{\sqrt{3}}{J}}\frac{\arctan\left(\frac{K}{r+m}\right)}{4\pi
K}$
$\displaystyle\times\left[F_{r,m}(K)(i\epsilon+c\mathbf{p}\cdot\bm{\sigma})+2m\sigma_{z}\right]$
(2ia) where $K=\sqrt{\epsilon^{2}+c^{2}p^{2}}$ and $\displaystyle F_{r,m}(K)$
$\displaystyle=1-\frac{r-m}{r+m}\left[\frac{r+m}{K\arctan\left(\frac{K}{r+m}\right)}-\frac{(r+m)^{2}}{K^{2}}\right].$
(2ib)

This result for a non-Fermi liquid GF is central to this work and has multiple
implications. First, we remind the reader of the minimal energy associated to
the creation of a physical fermion, which fractionalizes into $e$-particle
(mass $r$) and orthogonal fermion (mass $m$). This is technically reflected in
the branch cuts of the inverse trigonometric function entering the
analytically continued $\mathbf{G}^{R}(E,\mathbf{p})$ and physically implies a
continuum in the spectral weight for energies $E>r+m$, but a spectral gap
below this energy, Fig. 3 b). Second, as can be readily guessed from the $r=m$
limit of the GF numerator of Eq. (2i), there are zeros in the numerator of
zeros with dispersion
$\epsilon_{\emptyset}(\mathbf{p})=\pm\sqrt{m_{\emptyset}^{2}+c^{2}p^{2}}$.
Third, the case $r=m$ is very special, as the mass of zeros
$m_{\emptyset}=r+m=2r$ coincides with the minimal energy of the continuum of
poles (the branch cut). For $r>m$ the minimal energy $m_{\emptyset}$ of GF
zeros decreases and vanishes when $r/m\rightarrow\infty$, Fig. 3 b). However,
for negative $r-m$, the cost to create an orthogonal fermion exceeds that of
an $e$-particle, and the solutions to $\det[\mathbf{G}^{R}(E,0)]=0$ become
complex[48]. By invoking semiclassics, this also implies the absence of
fermionic bound states or boundary zeros at smooth domain walls in $r$
describing ZTI-TI interfaces in real space. Fourth, to obtain the equivalent
of Eq. (2i) for the K’ point on needs to replace $m\rightarrow\bar{m}$ in
front of $\sigma_{z}$ in the numerator and $m\rightarrow|\bar{m}|$ anywhere
else. At $\bar{m}=0$, this implies a direct OI-ZTI transition with a Dirac
node of zeros and semiclassically a boundary zero at an OI-ZTI interface, Fig.
3 d). Fifth, it is worthwhile to highlight a very different behavior in the
spectral weight $-\text{Im}[\mathbf{G}^{R}(E,0)]$ for $r>m$ and $r<m$,
respectively. In particular, when $r>m$, i.e. the $e$-particle gap exceeds the
fermionic gap, the spectral weight displays a rapid suppression as the energy
is lowered from outside the spectral gap towards $|E|=r+m$. In contrast, for
$r<m$ the spectral gap at at the K-point features “coherence” peak structures,
which can be interpreted as a precursor of the previously mentioned finite
quasiparticle weight on the confining side of the phase transition (i.e. when
$r^{2}$ becomes negative). These features are in principle visible in the
color coding of Fig. 3 b), but exemplified more clearly for two cross cuts in
Fig. 3 c). Finally, one may use Eq. (2ia) to evaluate the topological GF
invariant, Eq. (1). Summation over valleys yields $N_{3}=1$ ($N_{3}=0$) in the
topological (non-topological) case $w_{2}/w_{0}>1/\sqrt{27}$
($w_{2}/w_{0}<1/\sqrt{27}$) corresponding to a ZTI (OI). Crucially, the
contribution to $N_{3}$ entirely stems from the bands of zeros, while
“Hubbard” bands of poles drop out of $N_{3}$ [48, 38, 16]. As neither poles
nor zeros display a gap closing as $h/J$ is increased, it follows that
$N_{3}=1$ throughout the ZTI phase.

Figure 3: a) Dispersion of low-lying $e$-particles (red), approximated by a
hyperbolic dispersion relation, and of low-lying fermions (blue) with
effective Hamiltonian Eq. (2h). b) Color plot: Spectral weight obtained from
Eq. (2i). The equation $\det[\mathbf{G}^{R}(E,\mathbf{p}=0)]=0$ defines the
lowest energy of a GF zero (red). Importantly, this energy acquires an
imaginary part (a “lifetime”) once the solution hits the continuum (i.e. once
$r$ drops below $m$). c) Top: (Bottom:) Spectral weight at $r-m=[r+m]/2$
($r-m=-[r+m]/2$). d) Planar tunneling setup in which the ZTI system is
sandwiched between two metallic leads.

Transport signatures. The present integrable model readily allows to access
quantized transport coefficients (here the Hall conductance). Importantly,
unlike the Hatsugai-Kohmoto approach, the Hamiltonian Eq. (2) is local in real
space and thus straightforwardly coupled to external U(1) gauge fields.
Lattice models of the type of Eq. (2) encode “orthogonal” metals, semimetals
or (topological) insulators[23, 24, 25, 51]. A defining characteristic of the
latter is that two-particle GF (e.g. polarization bubble or current-current
correlator) are exactly the same as in the corresponding non-fractionalized
phase, despite the fact that the single particle propagator (i.e. the two-
point function) resembles a Mott insulator. Another way to state this is that
diagrams associated to electromagnetic response are evaluated using
$G(i\epsilon,\mathbf{p})$ rather than $\mathbf{G}(i\epsilon,\mathbf{p})$. To
understand why, for simplicity consider the hopping term along a single link
$c^{\dagger}_{\mathbf{x}}Z_{\mathbf{x},\mathbf{x}^{\prime}}c_{\mathbf{x}^{\prime}}+H.c.$.
Peierls substitution and derivative with respect to an external U(1) gauge
field leads to a current operator $j\sim
ic^{\dagger}_{\mathbf{x}}Z_{\mathbf{x},\mathbf{x}^{\prime}}c_{\mathbf{x}^{\prime}}+H.c.$.
This current operator commutes with $Q_{\mathbf{x}_{0}}$ for any
$\mathbf{x}_{0}$, and, by consequence, gapped $e$-particle propagators are not
entering the current-current correlators, Fig. 2 a). As a corollary, it
follows that the Hall conductance is [25]

$\sigma_{H}=\Theta(\sqrt{27}w_{2}-w_{0}),$ (2j)

irrespective of the value of $h/J$. Throughout the letter, including Fig. 1
b), we express $\sigma_{H}$ in units of the von-Klitzing constant.

What is the relationship between $N_{3}[\mathbf{G}]$, Eq. (1), and
$\sigma_{H}$? By applying the Widom-Streda formula[52, 53] it has been
shown[38, 39, 37, 12] on general grounds that
$\sigma_{H}[\mathbf{G}]=N_{3}[\mathbf{G}]+2\pi\partial_{B}I_{L}[\mathbf{G}]|_{B=0}/e$,
where $I_{L}[\mathbf{G}]$ is the Luttinger integral[48]. Note that in the
context of this proof, $\sigma_{H}[\mathbf{G}]$ is evaluated with the full GF
$\mathbf{G}$, even though, as argued above, for the present system the
physical Hall response $\sigma_{H}$ is evaluated using the free fermion GF
$G$. We explicitly evaluate the Widom-Streda formula at $J=0$ to deduce
$\sigma_{H}[\mathbf{G}]=\Theta(\sqrt{27}w_{2}-w_{0})$[48]. Given that we
independently calculated $N_{3}[\mathbf{G}]=0$ in this regime, one may
conclude
$2\pi\partial_{B}I_{L}[\mathbf{G}]|_{B=0}/e=\Theta(\sqrt{27}w_{2}-w_{0})$ and
that $\sigma_{H}[\mathbf{G}]=\sigma_{H}$.

Discussion and Outlook. In summary, we have studied an exactly soluble, local
model of a fractionalized topological insulator [25, 54]. Analytically
controlled perturbation about the integrable limit reveals topological bands
of GF zeros, which are sharp deep in the deconfining phase, but acquire a
“lifetime” near the Higgs condensation.

As an outlook devoted to the experimental detection of (topological) GF zeros,
we discuss momentum conserving tunneling which has previously served as a
crucial tool in determining many-body Fermi and Luttinger liquid physics in 2D
(“planar tunneling”[55, 56]) and 1D [57], respectively. In a device with
momentum conservation such as Fig. 3 d), the decay rate $\Gamma(\epsilon^{\rm
lead}_{\mathbf{p}},\mathbf{p})$ of electrons in the upper lead determines the
current across the junction. We use the notation $\epsilon^{\rm
lead}_{\mathbf{p}},G^{R}_{\rm lead}$ for dispersion and GF of lead electrons,
respectively, consider identical top and bottom leads and, for simplicity,
scalar Green’s functions. As we argue in the following [48], the planar
tunneling rate and, by consequence, the tunneling conductance $dI/dV$ display
characteristic dips in the presence of GF zeros.

To lowest order $\Gamma(\epsilon^{\rm
lead}_{\mathbf{p}},\mathbf{p})=-t^{2}\text{Im}\mathbf{G}^{R}(\epsilon_{\mathbf{p}}^{\rm
lead},\mathbf{p})$ which vanishes for energies below the spectral gap of the
central layer of the device. In this case, the leading elastic contribution is
$\Gamma(\epsilon_{\mathbf{p}},\mathbf{p})=-t^{4}[\text{Re}\mathbf{G}^{R}(\epsilon_{\mathbf{p}}^{\rm
lead},\mathbf{p})]^{2}\text{Im}G^{R}_{\rm lead}(\epsilon^{\rm
lead}_{\mathbf{p}},\mathbf{p})$, which however also vanishes when the lead
electron energy and momentum $(\epsilon^{\rm lead}_{\mathbf{p}},\mathbf{p})$
coincide with energy and momentum of zeros in the central layer. For 2D planar
tunneling, these kinematic constraints might admittedly be challenging to
fine-tune experimentally. However, 1D momentum conserving tunneling in a
sandwich (chiral edge state – ZTI edge zero – chiral edge state) is practical,
since it can be expected that TI edge state and ZTI edge zero dispersions
cross each other at a given energy $E_{*}$. While momentum-conserving co-
tunneling in the absence of edge zeros implies a finite differential
conductance at $\mathcal{O}(t^{4})$, the same experiment shows a suppressed
signal at the characteristic voltage $E_{*}/e$ in the presence of topological
edge zeros.

Note added: During the completion of this manuscript we became aware of a
related work [58]; where there is overlap our results agree.

Acknowledgements: It is a pleasure to thank A. Blason, J. Cano, A. Chubukov,
L. Classen, M. Fabrizio, F. Pollmann, M. Vojta, S. Sur for useful discussions
on the topic. The authors are particularly indebted to G. Sangiovanni, N.
Wagner for sharing their unpublished results prior to posting them on the
arXiv. CS and EJK acknowledge hospitality by the Kavli Institute for
Theoretical Physics, where this work was initiated. This research was
supported in part by grants NSF PHY-1748958 and PHY-2309135 to the Kavli
Institute for Theoretical Physics (KITP). UFPS was supported by the Deutsche
Forschungsgemeinschaft (DFG, German Research Foundation) through a Walter
Benjamin fellowship, Project ID 449890867, and the DOE office of BES, through
award number DE-SC0020305.

## References

  * Abrikosov  _et al._ [2012] A. Abrikosov, L. Gorkov, and I. Dzyaloshinski,  _Methods of Quantum Field Theory in Statistical Physics_ , Dover Books on Physics (Dover Publications, 2012).
  * Dzyaloshinskii [2003] I. Dzyaloshinskii, Some consequences of the Luttinger theorem: The Luttinger surfaces in non-Fermi liquids and Mott insulators, Phys. Rev. B 68, 085113 (2003).
  * Sakai  _et al._ [2016] S. Sakai, M. Civelli, and M. Imada, Hidden fermionic excitation boosting high-temperature superconductivity in cuprates, Phys. Rev. Lett. 116, 057003 (2016).
  * Gazit  _et al._ [2020] S. Gazit, F. F. Assaad, and S. Sachdev, Fermi surface reconstruction without symmetry breaking, Physical Review X 10, 041057 (2020).
  * König  _et al._ [2020] E. J. König, P. Coleman, and A. M. Tsvelik, Soluble limit and criticality of fermions in z2 gauge theories, Phys. Rev. B 102, 155143 (2020).
  * Fabrizio [2022] M. Fabrizio, Emergent quasiparticles at luttinger surfaces, Nature communications 13, 1561 (2022).
  * Fabrizio [2023] M. Fabrizio, Spin-liquid insulators can be landau’s fermi liquids, Phys. Rev. Lett. 130, 156702 (2023).
  * Vasiliou  _et al._ [2023] K. Vasiliou, Y. He, and N. Bultinck, Electrons interacting with goldstone modes and the rotating frame, arxiv:2307.05699 (2023).
  * Nikolaenko  _et al._ [2023] A. Nikolaenko, J. von Milczewski, D. G. Joshi, and S. Sachdev, Spin density wave, fermi liquid, and fractionalized phases in a theory of antiferromagnetic metals using paramagnons and bosonic spinons, Phys. Rev. B 108, 045123 (2023).
  * Gurarie [2011] V. Gurarie, Single-particle green’s functions and interacting topological insulators, Phys. Rev. B 83, 085426 (2011).
  * Setty  _et al._ [2023a] C. Setty, S. Sur, L. Chen, F. Xie, H. Hu, S. Paschen, J. Cano, and Q. Si, Symmetry constraints and spectral crossing in a mott insulator with green’s function zeros, arxiv:2301.13870 (2023a).
  * Setty  _et al._ [2023b] C. Setty, F. Xie, S. Sur, L. Chen, M. G. Vergniory, and Q. Si, Electronic properties, correlated topology and green’s function zeros, arxiv:2309.14340 (2023b).
  * Setty  _et al._ [2023c] C. Setty, F. Xie, S. Sur, L. Chen, S. Paschen, M. G. Vergniory, J. Cano, and Q. Si, Topological diagnosis of strongly correlated electron systems, arXiv preprint arXiv:2311.12031 (2023c).
  * Mai  _et al._ [2023] P. Mai, J. Zhao, B. E. Feldman, and P. W. Phillips, 1/4 is the new 1/2 when topology is intertwined with mottness, Nature communications 14, 5999 (2023).
  * Manning-Coe and Bradlyn [2023] D. Manning-Coe and B. Bradlyn, Ground state stability, symmetry, and degeneracy in mott insulators with long range interactions, arxiv:2306.00221 (2023).
  * Wagner  _et al._ [2023a] N. Wagner, L. Crippa, A. Amaricci, P. Hansmann, M. Klett, E. König, T. Schäfer, D. Di Sante, J. Cano, A. Millis, _et al._ , Mott insulators with boundary zeros, Nat Commun 14, 7531 (2023a).
  * Hatsugai and Kohmoto [1992] Y. Hatsugai and M. Kohmoto, Exactly solvable model of correlated lattice electrons in any dimensions, Journal of the Physical Society of Japan 61, 2056 (1992).
  * Phillips  _et al._ [2018] P. W. Phillips, C. Setty, and S. Zhang, Absence of a charge diffusion pole at finite energies in an exactly solvable interacting flat-band model in d dimensions, Physical Review B 97, 195102 (2018).
  * Setty [2021a] C. Setty, Superconductivity from luttinger surfaces: Emergent sachdev-ye-kitaev physics with infinite-body interactions, Physical Review B 103, 014501 (2021a).
  * Setty [2021b] C. Setty, Dilute magnetic moments in an exactly solvable interacting host, arXiv preprint arXiv:2105.15205 (2021b).
  * Phillips  _et al._ [2020] P. W. Phillips, L. Yeo, and E. W. Huang, Exact theory for superconductivity in a doped mott insulator, Nature Physics 16, 1175 (2020).
  * Jabłonowski  _et al._ [2023] K. Jabłonowski, J. Skolimowski, W. Brzezicki, K. Byczuk, and M. M. Wysokiński, Topological mott insulator in the odd-integer filled anderson lattice model with hatsugai-kohmoto interactions, Phys. Rev. B 108, 195145 (2023).
  * Nandkishore  _et al._ [2012] R. Nandkishore, M. A. Metlitski, and T. Senthil, Orthogonal metals: The simplest non-fermi liquids, Phys. Rev. B 86, 045128 (2012).
  * Zhong  _et al._ [2012] Y. Zhong, K. Liu, Y.-Q. Wang, and H.-G. Luo, Correlated metallic state in honeycomb lattice: Orthogonal dirac semimetal, Phys. Rev. B 86, 165134 (2012).
  * Zhong  _et al._ [2013] Y. Zhong, Y.-F. Wang, and H.-G. Luo, ${Z}_{2}$ fractionalized chern/topological insulators in an exactly soluble correlated model, Phys. Rev. B 88, 045109 (2013).
  * Maciejko and Rüegg [2013] J. Maciejko and A. Rüegg, Topological order in a correlated chern insulator, Phys. Rev. B 88, 241101 (2013).
  * Gazit  _et al._ [2017] S. Gazit, M. Randeria, and A. Vishwanath, Emergent dirac fermions and broken symmetries in confined and deconfined phases of Z2 gauge theories, Nature Physics 13, 484 (2017).
  * Borla  _et al._ [2022] U. Borla, B. Jeevanesan, F. Pollmann, and S. Moroz, Quantum phases of two-dimensional Z2 gauge theory coupled to single-component fermion matter, Phys. Rev. B 105, 075132 (2022).
  * Roy and König [2023] K. Roy and E. J. König, ZN lattice gauge theories with matter fields, arxiv:2308.13083 (2023).
  * Golterman and Shamir [2023] M. Golterman and Y. Shamir, Propagator zeros and lattice chiral gauge theories, arxiv:2311.12790 (2023).
  * Emonts  _et al._ [2023] P. Emonts, A. Kelman, U. Borla, S. Moroz, S. Gazit, and E. Zohar, Finding the ground state of a lattice gauge theory with fermionic tensor networks: A $2+1\mathrm{D}$ Z2 demonstration, Phys. Rev. D 107, 014505 (2023).
  * Ishikawa and Matsuyama [1986] K. Ishikawa and T. Matsuyama, Magnetic field induced multi-component QED 3 and quantum hall effect, Zeitschrift für Physik C Particles and Fields 33, 41 (1986).
  * Volovik and Yakovenko [1989] G. Volovik and V. M. Yakovenko, Fractional charge, spin and statistics of solitons in superfluid 3He film, Journal of Physics: Condensed Matter 1, 5263 (1989).
  * Volovik [2003] G. E. Volovik,  _The universe in a helium droplet_ , Vol. 117 (OUP Oxford, 2003).
  * Essin and Gurarie [2011] A. M. Essin and V. Gurarie, Bulk-boundary correspondence of topological insulators from their respective green’s functions, Phys. Rev. B 84, 125132 (2011).
  * Raghu  _et al._ [2008] S. Raghu, X.-L. Qi, C. Honerkamp, and S.-C. Zhang, Topological mott insulators, Phys. Rev. Lett. 100, 156401 (2008).
  * Zhao  _et al._ [2023] J. Zhao, P. Mai, B. Bradlyn, and P. Phillips, Failure of topological invariants in strongly correlated matter, Phys. Rev. Lett. 131, 106601 (2023).
  * Blason and Fabrizio [2023] A. Blason and M. Fabrizio, Unified role of green’s function poles and zeros in correlated topological insulators, Phys. Rev. B 108, 125115 (2023).
  * Peralta Gavensky  _et al._ [2023] L. Peralta Gavensky, S. Sachdev, and N. Goldman, Connecting the many-body chern number to luttinger’s theorem through Streda’s Formula, Phys. Rev. Lett. 131, 236601 (2023).
  * Manmana  _et al._ [2012] S. R. Manmana, A. M. Essin, R. M. Noack, and V. Gurarie, Topological invariants and interacting one-dimensional fermionic systems, Phys. Rev. B 86, 205119 (2012).
  * Slagle  _et al._ [2015] K. Slagle, Y.-Z. You, and C. Xu, Exotic quantum phase transitions of strongly interacting topological insulators, Phys. Rev. B 91, 115121 (2015).
  * You  _et al._ [2018] Y.-Z. You, Y.-C. He, C. Xu, and A. Vishwanath, Symmetric fermion mass generation as deconfined quantum criticality, Phys. Rev. X 8, 011026 (2018).
  * Wang and You [2022] J. Wang and Y.-Z. You, Symmetric mass generation, Symmetry 14, 1475 (2022).
  * Haldane [1988] F. D. M. Haldane, Model for a quantum hall effect without landau levels: Condensed-matter realization of the ”parity anomaly”, Phys. Rev. Lett. 61, 2015 (1988).
  * de’Medici  _et al._ [2005] L. de’Medici, A. Georges, and S. Biermann, Orbital-selective mott transition in multiband systems: Slave-spin representation and dynamical mean-field theory, Phys. Rev. B 72, 205124 (2005).
  * Rüegg  _et al._ [2010] A. Rüegg, S. D. Huber, and M. Sigrist, Z2-slave-spin theory for strongly correlated Fermions, Phys. Rev. B 81, 155118 (2010).
  * Tupitsyn  _et al._ [2010] I. S. Tupitsyn, A. Kitaev, N. V. Prokof’ev, and P. C. E. Stamp, Topological multicritical point in the phase diagram of the toric code model and three-dimensional lattice gauge higgs model, Phys. Rev. B 82, 085114 (2010).
  * [48] See supplementary information for detailed calculations.
  * Levanyuk [1959] A. Levanyuk, Contribution to the theory of light scattering near the second-order phase-transition points, Sov. Phys. JETP 9, 571 (1959).
  * Ginzburg [1960] V. Ginzburg, Some remarks on second order phase transitions and microscopic theory of ferroelectrics, Fiz. Tverd. Tela 2, 2031 (1960).
  * Hohenadler and Assaad [2018] M. Hohenadler and F. F. Assaad, Fractionalized metal in a falicov-kimball model, Phys. Rev. Lett. 121, 086601 (2018).
  * Widom [1982] A. Widom, Thermodynamic derivation of the hall effect current, Physics Letters A 90, 474 (1982).
  * Streda [1982] P. Streda, Theory of quantised hall conductivity in two dimensions, Journal of Physics C: Solid State Physics 15, L717 (1982).
  * Stern [2016] A. Stern, Fractional topological insulators: A pedagogical review, Annual Review of Condensed Matter Physics 7, 349 (2016).
  * Murphy  _et al._ [1995] S. Q. Murphy, J. P. Eisenstein, L. N. Pfeiffer, and K. W. West, Lifetime of two-dimensional electrons measured by tunneling spectroscopy, Phys. Rev. B 52, 14825 (1995).
  * Britnell  _et al._ [2012] L. Britnell, R. V. Gorbachev, R. Jalil, B. D. Belle, F. Schedin, M. I. Katsnelson, L. Eaves, S. V. Morozov, A. S. Mayorov, N. M. Peres, _et al._ , Electron tunneling through ultrathin boron nitride crystalline barriers, Nano letters 12, 1707 (2012).
  * Barak  _et al._ [2010] G. Barak, H. Steinberg, L. N. Pfeiffer, K. W. West, L. Glazman, F. Von Oppen, and A. Yacoby, Interacting electrons in one dimension beyond the luttinger-liquid limit, Nature Physics 6, 489 (2010).
  * Wagner  _et al._ [2023b] N. Wagner, D. Guerci, A. Millis, and G. Sangiovanni, Edge zeros and boundary spinons in topological mott insulators, arXiv:2312.13226 (2023b).

Supplementary materials on

”Topological Green’s function zeros in an exactly solved model and beyond”
Steffen Bollmann1, Chandan Setty 2, Urban F. P. Seifert 3, Elio J. König1 1
Max-Planck Institute for Solid State Research, 70569 Stuttgart, Germany
2Department of Physics and Astronomy, Rice Center for Quantum Materials, Rice
University, Houston, Texas 77005, USA 3Kavli Institute for Theoretical
Physics, University of California, Santa Barbara, CA 93106, USA

These supplementary materials contain a summary of all definitions employed in
this manuscript, Sec. S1, details on the derivations of the Green’s function
in two different limits of the fractionalized phase, Sec. S2, a discussion of
the Widom-Streda formula, Sec. S3, of the Green’s function based topological
invariant Sec. S4, and about momentum conserving tunneling probes Sec. S5.

All citations refer to the bibliography of the main text.

##  S1 Definitions and Conventions

###  S1.1 The local lattice limit

We consider the thermodynamic limit and treat the lattice momentum k as a
continuous variable. This means, in particular, the sum

$\lim_{N\rightarrow\infty}\frac{1}{N}\sum_{\textbf{k}}\rightarrow\frac{1}{\text{Vol
B.Z.}}\int_{\text{B.Z.}}\text{d}\textbf{k}$ (S1)

becomes an integral. We chose the primitive lattice vectors

$\textbf{n}_{1}=\frac{1}{2}\begin{pmatrix}1\\\
\sqrt{3}\end{pmatrix}\quad\text{and}\quad\textbf{n}_{2}=\frac{1}{2}\begin{pmatrix}-1\\\
\sqrt{3}\end{pmatrix}$ (S2)

to be normalized. Thus the volume of the first Brillouin zone evaluates to
$\text{Vol B.Z.}=\frac{2}{\sqrt{3}}(2\pi)^{2}$. The k integration becomes
$\frac{\sqrt{3}}{2}\int(\text{d}k)$ where we use the shorthand notation
$(\text{d}k)=\frac{\text{d}\textbf{k}}{(2\pi)^{2}}$.

In the main text, we use the notation x for vertices of the honeycomb lattice.
Here we use $\textbf{x}=(\textbf{r},s)$ where $\mathbf{r}$ is a Bravais
lattice vector and $s\in\\{A,B\\}$ describe the sublattice degrees of freedom.
The Fourier transformation of the annihilation operator is defined as

$c_{\textbf{r},A/B}=\sqrt{\frac{\sqrt{3}}{2}}\int(\text{d}k)\text{e}^{i\textbf{k}\cdot\textbf{r}}c_{A/B}(\textbf{k}).$
(S3)

from which directly follows that

$G_{s,s^{\prime}}(\tau,\textbf{r})=\frac{\sqrt{3}}{2}\int(\text{d}k)\text{e}^{i\textbf{k}\cdot\textbf{r}}G_{s,s^{\prime}}(\tau,\textbf{k}),$
(S4)

where
$G_{s,s^{\prime}}(\tau,\textbf{r})=-\langle\mathcal{T}\left[c^{\phantom{\dagger}}_{\textbf{r},s}(\tau)c^{\dagger}_{0,s^{\prime}}(\tau)\right]\rangle$
is the imaginary time ordered GF and $s,s^{\prime}\in\\{A,B\\}$.

###  S1.2 The Model

For further references, a more detailed representation of the model
Hamiltonian (2), $H=H_{\rm TC}+H_{\rm f}+\delta H$ is worthwhile:

$\displaystyle H_{\rm TC}$
$\displaystyle=-K\sum_{\varhexagon}B_{\varhexagon}-h\sum_{\mathbf{r}}[Q_{\mathbf{r}}^{A}+Q_{\mathbf{r}}^{B}],$
(S5a) $\displaystyle H_{\rm f}$
$\displaystyle=-w\sum_{\mathbf{r}}\left[\sum_{\mathbf{r}^{\prime}\in\\{\mathbf{r},\mathbf{r}-\hat{n}_{1},\mathbf{r}-\hat{n}_{2}\\}}c^{\dagger}_{\mathbf{r},A}Z_{\mathbf{b}_{(\mathbf{r},A),(\mathbf{r}^{\prime},B)}}c_{\mathbf{r}^{\prime},B}+H.c.\right]+w_{0}\sum_{\mathbf{r}}\sum_{s}(-1)^{s}c^{\dagger}_{\mathbf{r},s}c_{\mathbf{r},s}$
$\displaystyle-
iw_{2}\sum_{\mathbf{r}}\left[\sum_{\begin{subarray}{c}\mathbf{r}^{\prime}\in\\\
\\{\mathbf{r}+\hat{n}_{1}-\hat{n}_{2},\mathbf{r}-\hat{n}_{1}.,\mathbf{r}+\hat{n}_{2}\\}\end{subarray}}c^{\dagger}_{\mathbf{r},A}\prod_{\begin{subarray}{c}\mathbf{b}\in\\\
\gamma_{(\mathbf{r},A),(\mathbf{r}^{\prime},A)}\end{subarray}}Z_{\mathbf{b}}c_{\mathbf{r}^{\prime},A}-\sum_{\begin{subarray}{c}\mathbf{r}^{\prime}\in\\\
\\{\mathbf{r}+\hat{n}_{1}-\hat{n}_{2},\mathbf{r}-\hat{n}_{1}.,\mathbf{r}+\hat{n}_{2}\\}\end{subarray}}c^{\dagger}_{\mathbf{r},B}\prod_{\begin{subarray}{c}\mathbf{b}\in\\\
\gamma_{(\mathbf{r},B),(\mathbf{r}^{\prime},B)}\end{subarray}}Z_{\mathbf{b}}c_{\mathbf{r}^{\prime},B}-H.c.\right],$
(S5b) $\displaystyle\delta H$
$\displaystyle=-J\sum_{\mathbf{b}}Z_{\mathbf{b}},$ (S5c)

where
$B_{\varhexagon}=\prod_{\mathbf{b}\in\partial\varhexagon}Z_{\mathbf{b}},Q_{\mathbf{r}}^{s}=(-1)^{\hat{n}_{\mathbf{r}}^{s}}\prod_{\mathbf{b}\in\partial(\mathbf{r},s)}X_{\mathbf{b}}$

and $X_{\textbf{b}},Z_{\textbf{b}}$ are Pauli Matrices living on the bond b.
The path $\gamma_{(\mathbf{r},s),(\mathbf{r}^{\prime},s^{\prime})}$ connects
two points on the shortest path.

Note that each unit cell contains three qubits and three integrals of motion
(at $J=0$, that is), i.e.
$B_{\varhexagon},Q_{\mathbf{r}}^{A},Q_{\mathbf{r}}^{B}$.

###  S1.3 Ground state in the exactly solvable limit $J=0$

As shown in [5], the ground state of a TC coupled to fermions can be
constructed by projecting the tensor product of a fully polarized gauge field
state with the fermionic sector.

To explicitly construct the Fermi sea, we need to switch to Fourier space.
Starting from the gauge $Z_{\mathbf{b}}=1$ everywhere, we find

$\displaystyle H_{\rm f}=\int_{\rm
BZ}(\text{d}k)\left(c^{\dagger}_{A}(\mathbf{k}),c_{B}^{\dagger}(\mathbf{k})\right)h(\mathbf{k})\left(\begin{array}[]{c}c_{A}(\mathbf{k})\\\
c_{B}(\mathbf{k})\end{array}\right)$ (S6c) Here, we used the Fourier transform
as explained in detail in App. S1.1. The single-particle Hamiltonian has the
form
$h(\mathbf{k})=\left(\begin{array}[]{cc}w_{2}m(\mathbf{k})+w_{0}&-w_{1}s(\mathbf{k})\\\
-w_{1}s^{*}(\mathbf{k})&-w_{2}m(\mathbf{k})-w_{0}\end{array}\right)$ (S6d)
with $\displaystyle s(\mathbf{k})$
$\displaystyle=[e^{-i\mathbf{k}\cdot\mathbf{n}_{1}}+e^{-i\mathbf{k}\cdot\mathbf{n}_{2}}+1]$
(S6g) $\displaystyle m(\mathbf{k})$
$\displaystyle=\left[\sin(\mathbf{k}\cdot\mathbf{n}_{1-2})-\sin(\mathbf{k}\cdot\mathbf{n}_{1})+\sin(\mathbf{k}\cdot\mathbf{n}_{2})\right]$
(S6h)

We will also need the energy

$\epsilon(\mathbf{k})=\sqrt{w_{1}^{2}|s(\mathbf{k})|^{2}+(w_{2}m(\mathbf{k})+w_{0})^{2}}.$
(S7)

##  S2 Fermionic Green’s Function in the presence of topological order

According to the discussion in the main text and a comprehensive discussion in
Ref. [5], the full GF for fermionic excitations in the OTI and ZTI is given as
the component-wise product of the fermionic GF
$G_{s,s^{\prime}}(\tau,\mathbf{r},\mathbf{r}^{\prime})$ and the $e$-particle
GF $D_{s,s^{\prime}}(\tau,\mathbf{r},\mathbf{r}^{\prime})$, that is

$\displaystyle\textbf{G}_{s,s^{\prime}}(\tau,\mathbf{r},\mathbf{r}^{\prime})=G_{s,s^{\prime}}(\tau,\mathbf{r},\mathbf{r}^{\prime})D_{s,s^{\prime}}(\tau,\mathbf{r},\mathbf{r}^{\prime})$
(S8)

which translates to the convolution

$\textbf{G}_{s,s^{\prime}}(i\epsilon,\textbf{k})=\frac{\sqrt{3}}{2}\int(\text{d}Q)D_{s,s^{\prime}}(i\omega,q)G_{s,s^{\prime}}(i(\epsilon-\omega),\textbf{k}-\textbf{q})$
(S9)

in momentum/Matsubara frequency space, and $(dQ)=d\omega d^{2}q/(2\pi)^{3}$.

###  S2.1 Propagator of orthogonal fermions

As explained in the main text, the GF of the orthogonal fermions in momentum
space is given as

$G(i\epsilon,\mathbf{k})=[i\epsilon-h(\mathbf{k})]^{-1}=-\frac{i\epsilon+h(\mathbf{k})}{\epsilon^{2}+\epsilon(\mathbf{k})^{2}}.$
(S10)

###  S2.2 Propagator of $e$-particles

In the limit of vanishing $J$, the $e$-particle can not move over the lattice,
and the GF becomes ultra-local and takes the simple shape

$D^{(0)}_{s,s^{\prime}}(\tau,\mathbf{r},\mathbf{r}^{\prime})=\text{e}^{-2h|\tau|}\delta_{s,s^{\prime}}\delta_{\mathbf{r},\mathbf{r}^{\prime}}$
(S11)

or

$D^{(0)}_{s,s^{\prime}}(i\omega,\mathbf{q})=\frac{4h}{\omega^{2}+4h^{2}}.$
(S12)

If non-Gaussian fluctuations of the gauge field are subleading (an estimate is
given in Sec. S2.4, below), the $e$-particle propagator can be approximated by
the recursion formula

$\displaystyle D_{s_{f},s_{i}}\left(\mathbf{r}_{f},i\omega\right)$
$\displaystyle=D^{(0)}(i\omega)\delta_{s_{f},s_{i}}\delta_{\mathbf{r}_{f},\mathbf{r}_{i}}$
$\displaystyle+J\sum_{\langle(\mathbf{r},s),(\mathbf{r}_{f},s_{F})\rangle}D_{s,s_{i}}(\mathbf{r},i\omega)D^{(0)}(i\omega),$
(S13)

where the brackets $\langle...\rangle$ stand for nearest neighbors. It
incorporates all possible paths the particle can take from its initial
position $(r_{i},s_{i})$ to the final position $(r_{f},s_{f})$ that do not
intersect themselves. A more detailed discussion can be found in Ref. [5, 29].

In momentum space, the recursion equation (S13) takes the shape

$D(i\omega,\mathbf{q})=D^{(0)}\left[\bm{1}+J\begin{pmatrix}0&s(\mathbf{q})\\\
s^{*}(\mathbf{q})&0\end{pmatrix}D(i\omega,\mathbf{q})\right]$ (S14)

and we obtain a formal solution

$\begin{split}\,&D(i\omega,\mathbf{q})=4h\left[\omega^{2}+4h^{2}-4hJ\left(\begin{array}[]{cc}0&s(\mathbf{q})\\\
s^{*}(\mathbf{q})&0\end{array}\right)\right]^{-1}\\\
&=\frac{4h}{\omega^{2}+4h^{2}-(4hJ)^{2}|s(\mathbf{q})|^{2}}\begin{pmatrix}\omega^{2}+4h^{2}&4hJs(\mathbf{q})\\\
4hJs^{*}(\mathbf{q})&(\omega^{2}+4h^{2})\end{pmatrix}\end{split}$ (S15)

for the $e$-particle propagator. The matrices act in sublattice space.

In the following, we present the calculation of $\mathbf{G}$, first, in a
perturbative expansion for $h/J\gg 1$ and, second, in a continuum limit when
$J\sim h$ where we, however, stay far enough from the Higgs transition where
non-Gaussian gauge field fluctuations become important.

###  S2.3 Pertuabative GF for $h/J\gg 1$

We now expand Eq. (S15) up to the second order in $J$ and insert it order by
order in the convolution (S9).

For $w_{0}=0$, the dispersion relation $\epsilon(\mathbf{k})$ would be $C_{6}$
symmetric, simplifying further calculations. Any finite $w_{0}$ reduces the
symmetry to $C_{3}$. However, if we stay within the hierarchy of limits $h\gg
w_{1}\gg w_{2}\gg w_{0}$ (in the top left of the phase diagram 1b) and
$h\gg\epsilon$ we can expand the dispersion relation

$\epsilon(\mathbf{k})=\sqrt{\underbrace{w_{1}^{2}|s(\mathbf{k})|^{2}+w_{2}^{2}m(\mathbf{k})^{2}+w_{0}^{2}}_{\xi^{2}(\mathbf{k})}+2w_{0}w_{2}m(\mathbf{k})}.$
(S16)

Note that the term $2w_{0}w_{2}m(k)$ is, for any $\mathbf{k}$, small compared
to $\xi^{2}(\mathbf{k})$. Thus, we can expand $\epsilon(\mathbf{k})$ in powers
of $w_{0}w_{2}m(\mathbf{k})$ and obtain up to first order

$\epsilon(\mathbf{k})\simeq\xi(\mathbf{k})+\frac{w_{0}w_{2}m(\mathbf{k})}{\xi(\mathbf{k})}.$
(S17)

We split $\epsilon(\mathbf{k})$ into a dominating $C_{6}$ symmetric and a
subleading $C_{6}$ anti-symmetric part.

It will be useful to introduce the function

$g(x,y)=\frac{1}{x^{2}+y^{2}}.$ (S18)

####  S2.3.1 GF to order $J^{0}$

The integral

$\begin{split}\mathbf{G}(i\epsilon)|_{J^{0}}&=-\frac{\sqrt{3}}{2}\int(\text{d}q)(\text{d}\omega)D^{(0)}(i\omega)\\\
&\times\frac{i(\epsilon-\omega)+(w_{2}m(\mathbf{k}-\mathbf{q})+w_{0})\sigma_{z}}{(\epsilon-\omega)^{2}+\epsilon^{2}(\mathbf{k}-\mathbf{q})}\end{split}$
(S19)

gives the zeroth order contribution to the GF. Performing the $\omega$
integration and applying the shift
$\mathbf{q}\rightarrow-\mathbf{q}+\mathbf{k}$ yields

$\begin{split}\mathbf{G}(i\epsilon)|_{J^{0}}&=-\frac{\sqrt{3}}{2}\int(\text{d}q)\Bigg{(}i\epsilon
g(\epsilon,\epsilon(\mathbf{q})+2h)\\\
&+\left.\left(1+\frac{2h}{\epsilon(\mathbf{q})}\right)g(\epsilon,\epsilon(\mathbf{q})+2h)(w_{2}m(\mathbf{q})+w_{0})\sigma_{z}\right).\end{split}$
(S20)

Now, we use that (in the assumed limit) $h$ is the dominant energy scale, and
approximate the function
$g(\epsilon,\epsilon(\mathbf{k}))\rightarrow\frac{1}{4h^{2}}$. This
approximation will be used throughout the rest of this section. Thus, the
first part simply contributes $\frac{-i\epsilon}{4h^{2}}$. In the second part,
we only keep the term in leading order in $1/h$.

Also, we expand the function
$g(\epsilon,\epsilon(\mathbf{q})+2h)/\epsilon(\mathbf{q})$ in the subleading
$C_{6}$ anti-symmetric part, according to the procedure outlined in the last
section. We obtain

$\begin{split}\frac{g(\epsilon,\epsilon(\mathbf{q}+2h))}{\epsilon(\mathbf{q})}&\approx\frac{g(\epsilon,\xi(\mathbf{q})+2h)}{\xi(\mathbf{q})}\\\
&-\frac{w_{0}w_{2}g(\epsilon,\xi(\mathbf{q}))}{\xi^{3}(\mathbf{q})}m(\mathbf{q})\\\
&-2\frac{w_{0}w_{2}g^{2}(\epsilon,\xi(\mathbf{q})+2h)}{\xi^{2}(\mathbf{q})}\left(1+\frac{2h}{\xi(\mathbf{q})}\right)m(\mathbf{q}).\end{split}$
(S21)

All anti-symmetric terms (i.e., those are proportional to $m(\mathbf{q})$) are
either of higher order in $1/h$ or are suppressed by additional factors
$w_{0}/w_{1}$ or $w_{2}/w_{1}$ compared to the symmetric part. Therefore, we
also ignore these terms. Applying the approximations yields

$\begin{split}\mathbf{G}|_{J^{0}}(i\epsilon)&=-\frac{i\epsilon}{4h^{2}}-\frac{\sqrt{3}}{2}\int(\text{d}q)\frac{2hg(\epsilon,\xi(\mathbf{q})+2h)}{\xi(\mathbf{q})}w_{0}\sigma_{z}\\\
&\approx-\frac{i\epsilon}{4h^{2}}-\frac{1}{hw_{1}}\underbrace{\sqrt{3}\int(\text{d}q)\frac{1}{|s(\mathbf{q})|}}_{\mathcal{C}_{0}}\sigma_{z}\\\
&=-\frac{1}{4h^{2}}(i\epsilon+\frac{h}{w_{1}}\mathcal{C}_{0}\,w_{0}\sigma_{z}),\end{split}$
(S22)

where $\mathcal{C}_{0}\approx 1.8$. Here, we first used that, due to symmetry,

$\int(\text{d}q)g(\epsilon,\xi(\mathbf{q})+2h)m(\mathbf{q})=0$ (S23)

and, second, that in the limit $w_{1}\gg w_{2}\gg w_{0}$ the scaling of
$\xi(\mathbf{q})$ is essentially given by $w_{1}$. Therefore, we can
approximate $\xi(\mathbf{q})\approx w_{1}|s(\mathbf{q})|$.

####  S2.3.2 GF to order $J^{1}$

The integral

$\begin{split}\mathbf{G}(i\epsilon,\mathbf{k})|_{J^{1}}&=-J\frac{\sqrt{3}}{2}\int(\text{d}q)(\text{d}\omega)(D^{(0)}(i\omega))^{2}\\\
&\times\begin{pmatrix}0&\frac{-w_{1}s(\mathbf{q})s(\mathbf{k}-\mathbf{q})}{(\epsilon-\omega)^{2}+\epsilon^{2}(\mathbf{k}-\mathbf{q})}\\\
\frac{s^{*}(\mathbf{q})s^{*}(\mathbf{k}-\mathbf{q})}{{(\epsilon-\omega)^{2}+\epsilon^{2}(\mathbf{k}-\mathbf{q})}}&0\end{pmatrix}\end{split}$
(S24)

gives the first order contribution to the GF. Using the identity

$(D^{(0)}(i\omega))^{2}=-\frac{h}{2}\partial_{h}\frac{1}{h}D^{(0)}(i\omega)$
(S25)

the frequency integration is analogous to the one for the zeroth order.
Performing the integration and only keeping terms that are leading in $1/h$ we
obtain

$\begin{split}\mathbf{G}(i\epsilon,\mathbf{k})|_{J^{1}}&=Jh\frac{\sqrt{3}}{2}\partial_{h}\int(\text{d}q)\frac{g(\epsilon,\xi(\mathbf{q})+2h)}{\xi(\mathbf{q})}\\\
&\times\begin{pmatrix}0&-w_{1}s(\mathbf{q})s(\mathbf{k}-\mathbf{q})\\\
-w_{1}s^{*}(\mathbf{q})s^{*}(\mathbf{k}-\mathbf{q})&0\end{pmatrix}\\\
&=Jh\frac{\sqrt{3}}{2}\begin{pmatrix}0&-w_{1}s(\mathbf{k})\\\
-w_{1}s^{*}(\mathbf{k})\end{pmatrix}\\\
&\times\partial_{h}\int(\text{d}q)\frac{g(\epsilon,\xi(\mathbf{q})+2h)}{\xi(\mathbf{q})}(1+2\cos(\mathbf{q}\cdot\mathbf{n}_{1/2}))\\\
&\approx-\frac{1}{4h^{2}}J\mathcal{C}_{0}\mathcal{C}_{1}\begin{pmatrix}0&-s(\mathbf{k})\\\
-s^{*}(\mathbf{k}),\end{pmatrix}\end{split}$ (S26)

where

$\mathcal{C}_{0}\mathcal{C}_{1}\equiv\frac{\sqrt{3}}{4}\int(\text{d}q)\frac{1+2\cos(\mathbf{q}\cdot\mathbf{n}_{1/2})}{|s(\mathbf{q})|}\approx
0.26.$ (S27)

To achieve the second equality in Eq. (S26), we expand the product
$s(\mathbf{q})s(\mathbf{k}-\mathbf{q})$ and use the $C_{6}$ symmetry of the
integrand. $\mathbf{n}_{1/2}$ denotes that both $\mathbf{n}_{1}$ and
$\mathbf{n}_{2}$ can be used since they yield the same numerical value (again
due to the $C_{6}$ symmetry).

####  S2.3.3 GF to order $J^{2}$

The integral

$\begin{split}\mathbf{G}(i\epsilon,\mathbf{k})|_{J^{2}}&=-\frac{\sqrt{3}}{2}J^{2}\int(\text{d}q)(\text{d}\omega)(D^{(0)}(i\omega))^{3}\\\
&\times|s(q)|^{2}\frac{w_{2}m(\mathbf{k}-\mathbf{q})+w_{0}}{(\epsilon-\omega)^{2}+\epsilon(\mathbf{k}-\mathbf{q})}\end{split}$
(S28)

gives the second-order contribution to the GF. Using a similar trick as for
the first-order, that is

$(D^{0}(i\omega))^{3}=\underbrace{\frac{1}{8}\left(h\partial_{h}^{2}-\partial_{h}\right)}_{\mathcal{D}_{h}}\frac{1}{h}D^{(0)}(i\omega)$
(S29)

enables a straightforward frequency integration. We obtain

$\begin{split}\mathbf{G}(i\epsilon,\mathbf{k})|_{J^{2}}&=-\frac{\sqrt{3}}{2}\mathcal{D}_{h}\frac{1}{h}\int(\text{d}q)|s(\mathbf{k}-\mathbf{q})|^{2}m(\mathbf{q})\\\
&\times\frac{g(\epsilon,\xi(\mathbf{q})+2h)}{\xi(\mathbf{q})}\sigma_{z}\\\
&\approx-\frac{\sqrt{3}}{{2}}\mathcal{D}_{h}\frac{m(\mathbf{k})}{2h^{2}}\int(\text{d}q)\\\
&\times\frac{2\sin^{2}(\mathbf{q}\cdot\mathbf{n}_{1/2})\sin^{2}(\frac{\mathbf{q}\cdot\mathbf{n}_{2/1}}{2})}{w_{1}|s(\mathbf{q})|}\sigma_{z}\\\
&=-\frac{1}{4h^{2}}\frac{J^{2}}{hw_{1}}\mathcal{C}_{0}\mathcal{C}_{2}m(\mathbf{k})\sigma_{z},\end{split}$
(S30)

where

$\mathcal{C}_{0}\mathcal{C}_{2}=4\sqrt{3}\int(\text{d}q)\frac{\sin^{2}(\mathbf{q}\cdot\mathbf{n}_{1/2})}{|s(\mathbf{q})|}\sin^{2}\left(\frac{\mathbf{q}\cdot\mathbf{n}_{2/1}}{2}\right)\approx
2.4.$ (S31)

Here, we used trigonometric identities again to expand
$|s(\mathbf{k}-\mathbf{q})|^{2}m(\mathbf{q})$ and the $C_{6}$ symmetry to
simplify the expression under the integral sign.

####  S2.3.4 Full Result

Adding all calculated orders yields the final expression for the GF:

$\begin{split}\mathbf{G}(i\epsilon,\mathbf{k})&=\mathbf{G}(i\epsilon,\mathbf{k})|_{J^{0}}+\mathbf{G}(i\epsilon,\mathbf{k})|_{J^{1}}+\mathbf{G}(i\epsilon,\mathbf{k})|_{J^{2}}\\\
&=\frac{-1}{4h}\left(i\epsilon+\mathcal{C}_{0}\frac{h}{w_{1}}h_{\emptyset}\right),\end{split}$
(S32a) where
$h_{\emptyset}=\begin{pmatrix}w_{0}+\frac{J^{2}}{h^{2}}w_{2}\mathcal{C}_{2}m(\mathbf{k})&-\frac{J}{h}w_{1}\mathcal{C}_{1}s(\mathbf{k})\\\
-\frac{J}{h}w_{1}\mathcal{C}_{1}s^{*}(\mathbf{k})&-w_{0}-\frac{J^{2}}{h^{2}}\mathcal{C}_{2}w_{2}m(\mathbf{k})\end{pmatrix}.$
(S32b)

This is the origin of equation (2g).

###  S2.4 Estimate of Ginzburg-Levanyuk number

For the next section, in which we employ a continuum approximation of both
e-particle and fermionic propagators it is useful to discuss the limits of
applicability of the present theories without non-Gaussian fluctuaions.

The confinement transition is described by an effective $\phi^{4}$ theory

$\mathcal{L}=\phi[-\partial_{\tau}^{2}-\nabla^{2}+r^{2}]\phi+\lambda\phi^{4}$
(S33)

Note that $[r]=E,[\phi]=E^{-1+D/2}$, $[\lambda]=E^{4-D}$ One may estimate the
regime of negligible non-Gaussian fluctuations as

$r>\lambda^{-\frac{1}{D-4}}.$ (S34)

Our bare estimates of coupling constants are[5] are $r=\sqrt{h^{2}-3hJ}$,
$\lambda=J^{3}/h^{2}$, leading to

$\sqrt{h^{2}-3hJ}>J^{3}/h^{2}.$ (S35)

###  S2.5 $e$-particle and orthogonal fermion GFs in continuum limit

In the case where $J\sim h$ we can find an exact expression of
$\textbf{G}_{s,s^{\prime}}(i\epsilon,\mathbf{k})$ for the low energy/momentum
limit. First, we diagonalize the $e$-particle propagator (S15) and obtain

$\begin{split}D(i\omega,\mathbf{q})&=\frac{4h}{\omega^{2}+4h^{2}-4hJ|s(q)|}\begin{pmatrix}1&0\\\
0&0\end{pmatrix}\\\
&+\frac{4h}{\omega^{2}+4h^{2}+4hJ|s(q)|}\begin{pmatrix}0&0\\\
0&1\end{pmatrix}\\\
&\approx\frac{4h}{\omega^{2}+4h^{2}-4hJ|s(q)|}\begin{pmatrix}1&0\\\
0&0\end{pmatrix},\end{split}$ (S36)

where we only keep the lower e-particle band in the second line. The
transformation matrix reads

$T=\frac{1}{\sqrt{2}}\begin{pmatrix}\sqrt{\frac{s(q)}{s^{*}(q)}}&-\sqrt{\frac{s(q)}{s^{*}(q)}}\\\
1&1\end{pmatrix}\approx\frac{1}{\sqrt{2}}\begin{pmatrix}1&-1\\\
1&1\end{pmatrix}.$ (S37)

As we will be interested at the lowest $e-$particle excitaitons (near the
$\Gamma$ point), we expanded the transformation matrix near $\mathbf{q}=0$.
Transforming the propagator for the lower band back into the sublattice basis
using $T$ yields

$D(i\omega,\mathbf{q})\simeq\frac{4h}{\omega^{2}+4h^{2}-4hJ|s(\mathbf{q})|}\begin{pmatrix}1&1\\\
1&1\end{pmatrix}.$ (S38)

Since we multiply the $e$-particle and fermion GF component-wise, the matrix
structure of $D(i\omega,\mathbf{q}),$ allows us to treat it as a scalar in the
subsequent calculation.

In the long wavelength limit, we can approximate $|s(\mathbf{q})|$ as a
parabola (see Fig. 3 of the main text) and obtain for the $e$-particle
propagator

$D(i\omega,\mathbf{q})=\frac{4h}{\omega^{2}+r^{2}+v^{2}q^{2}},$ (S39)

where $v=\sqrt{hJ}$, $r^{2}=4h(h-3J)$ and $q=|\mathbf{q}|$.

The GF for the fermions in the long wavelength limit can be calculated
straightforwardly using the continuum limit of (S6d) and reads

$G(i\epsilon,\mathbf{p})=(i\epsilon-
h_{\text{eff}}(\mathbf{p}))^{-1}=-\frac{i\epsilon+h_{\text{eff}}(\mathbf{p})}{\epsilon^{2}+m^{2}+c^{2}p^{2}},$
(S40)

where $m=w_{0}+3\sqrt{3}w_{2}$, $c=\frac{\sqrt{3}}{2}w_{1}$,$p=|\mathbf{p}|$
and $\mathbf{p}$ is the momentum relative to the K point.

###  S2.6 Full GF in continuum’s limit

The full GF in the continuum limit is best calculated in real space. It is
straightforward to show that the $e$-particle propagator in real space reads
(we here use the notation $\mathbf{x}$ for the continuum coordinate to avoid
confusion with the e-particle gap $r$)

$D(\tau,\mathbf{x})=\frac{4h}{v^{2}}Y_{r}(X_{v}),$ (S41)

where

$Y_{r}(X_{v})=\frac{\text{e}^{-rX_{v}}}{4\pi
X_{v}}=\int(\text{d}\omega)(\text{d}q)\frac{v^{2}\text{e}^{i\mathbf{q}\cdot\mathbf{x}}\text{e}^{-i\omega\tau}}{\omega^{2}+r^{2}+v^{2}q^{2}}$
(S42)

and $X_{v}=\sqrt{\tau^{2}+|\textbf{r}|^{2}/v^{2}}$.

The fermion GF can be similarly calculated:

$\begin{split}G(\tau,\mathbf{x})&=-\int(dp)(d\epsilon)e^{i\mathbf{p}\cdot\mathbf{x}}e^{-i\epsilon\tau}\frac{i\epsilon+h(\mathbf{p})}{\epsilon^{2}+c^{2}p^{2}+m^{2}}\\\
&=-\int(dp)(d\epsilon)e^{i(\mathbf{p}\cdot\mathbf{x}+\epsilon\tau)}\frac{-i\epsilon+\mathbf{p}\cdot\bm{\sigma}+m\sigma_{z}}{\epsilon^{2}+m^{2}+c^{2}p^{2}}\\\
&=\frac{[\partial_{\tau}+ic\nabla\cdot\sigma-m\sigma_{z}]}{c^{2}}\int(dp)(d\epsilon)\frac{c^{2}e^{i(\mathbf{p}\cdot\mathbf{x}+\epsilon\tau)}}{\epsilon^{2}+m^{2}+c^{2}p^{2}}\\\
&=\frac{[\partial_{\tau}+ic\nabla\cdot\sigma-m\sigma_{z}]}{c^{2}}Y_{m}(X_{c}).\end{split}$
(S43)

.

If restrict too the special case $v=c=1$ (i.e.
$X_{v}=X_{c}=X=\sqrt{\tau^{2}+|\mathbf{x}|^{2}}$) the full GF in real space
reads

$\begin{split}\textbf{G}(\tau,\mathbf{x})&=2\sqrt{3}hY_{r}(X)\left[\partial_{\tau}+i\nabla\bm{\sigma}-m\sigma_{z}\right]Y_{m}(X)\\\
&=\sqrt{3}h\left[\mathcal{D}_{\tau,\mathbf{x}}+(r-m)\frac{\tau+i\mathbf{x}\cdot\bm{\sigma}}{X}-2m\sigma_{z}\right]\Upsilon_{r+m}(X)\end{split}$
(S44)

where $\Upsilon_{r+m}(X)=\frac{\text{e}^{-(r+m)X}}{(4\pi X)^{2}}$ and
$\mathcal{D}_{\tau,\mathbf{x}}=\partial_{\tau}+i\nabla\cdot\bm{\sigma}$.

Eq. (S44) can be Fourier transformed back into the momentum/frequency space.
First, we calculate

$\begin{split}\Upsilon_{r+m}(K)&=\int\text{d}\mathbf{x}\text{d}\tau\Upsilon_{r+m}(X)\text{e}^{-i\mathbf{x}\cdot\mathbf{p}}\text{e}^{i\epsilon\tau}\\\
&=\frac{1}{4\pi}\int_{0}^{\infty}\text{d}X\text{e}^{-(r+m)X}\frac{\sin(KX)}{KX}\\\
&=\frac{\arctan(K/(r+m))}{4\pi K}\end{split},$ (S45)

where $K=\sqrt{\epsilon^{2}+p^{2}}$. We can use this result to calculate
further

$\begin{split}\,&\int\text{d}\mathbf{x}\text{d}\tau\frac{(\tau+i\mathbf{x}\cdot\bm{\sigma})\text{e}^{-i\mathbf{x}\cdot\mathbf{p}}\text{e}^{i\epsilon\tau}}{X}\Upsilon_{r+m}(X)\\\
&=-\underbrace{(i\partial_{\epsilon}+\nabla_{\mathbf{p}}\cdot\bm{\sigma})}_{\mathcal{D}_{\epsilon,\mathbf{p}}}\int\text{d}\mathbf{x}\text{d}\tau\frac{\text{e}^{-i\mathbf{x}\cdot\mathbf{p}}\text{e}^{i\epsilon\tau}}{X}\Upsilon_{r+m}(X)\\\
&=-\frac{1}{4\pi}\int_{0}^{\infty}\text{d}X\mathcal{D}_{\epsilon,\mathbf{p}}\left(\frac{\sin(KX)}{KX^{2}}\right)e^{-(r+m)X}\\\
&=-\frac{\mathcal{D}_{\epsilon,\mathbf{p}}K}{4\pi}\int_{0}^{\infty}\text{d}X\frac{1}{X}\partial_{K}\left(\frac{\sin(KX)}{KX}\right)e^{-(r+m)X}\\\
&=-\frac{i\epsilon+\mathbf{p}\cdot\bm{\sigma}}{4\pi
K^{2}}\int_{0}^{\infty}\text{d}X\partial_{X}\left(\frac{\sin(KX)}{KX}\right)e^{-(r+m)X}\\\
&=-\frac{i\epsilon+\mathbf{p}\cdot\bm{\sigma}}{4\pi
K^{2}}\left(\frac{r+m}{K}\arctan(K/(r+m))-1\right)\\\
&=\frac{i\epsilon+\mathbf{p}\cdot\sigma}{r+m}f(K/(r+m))\Upsilon_{r+m}(K),\end{split}$
(S46)

where $f=(\frac{1}{x\arctan(x)}-\frac{1}{x^{2}})$. In the 4th equality we used
that
$\frac{1}{X}\partial_{K}\left(\frac{\sin(KX)}{KX}\right)=\frac{1}{K}\partial_{X}\left(\frac{\sin(KX)}{KX}\right)$
and
$\mathcal{D}_{\epsilon,\mathbf{k}}K=\frac{i\epsilon-\mathbf{k}\cdot\bm{\sigma}}{K}$.

Combining the equations (S44), (S45) and (S46) yields the exact GF in
frequency space:

$\displaystyle\mathbf{G}(i\epsilon,\mathbf{p})$
$\displaystyle=-{\frac{\sqrt{3}}{J}}\frac{\arctan\left(\frac{K}{r+m}\right)}{4\pi
K}$
$\displaystyle\times\left[F_{r,m}(K)(i\epsilon+\mathbf{p}\cdot\bm{\sigma})+2m\sigma_{z}\right]$
(S47a) where $\displaystyle F_{r,m}(K)$
$\displaystyle=1-\frac{r-m}{r+m}\left[\frac{r+m}{K\arctan\left(\frac{K}{r+m}\right)}-\frac{(r+m)^{2}}{K^{2}}\right].$
(S47b)

Note that we replaced $h$ by $1/J$ since we set $v=\sqrt{hJ}=1$. This is the
origin of Eq. (2i)

###  S2.7 Spectral properties of GF in continuum limit

In this section, we analytically continue the Green’s function (2i) to obtain
the Green’s function zeros from the retarded GF. We start analytical
continuing $K$, which is a well-behaved function of $i\epsilon$

$\lim_{i\epsilon\rightarrow
E+i0^{+}}K=i\underbrace{\sqrt{E^{2}-p^{2}}}_{\mathcal{K}},$ (S48)

where we assumed that $E>p$. Thus, the retarded Green’s function reads

$\begin{split}\mathbf{G}^{R}(E,\mathbf{p})&=-\frac{\sqrt{3}}{J}\frac{\text{artanh}(\frac{\mathcal{K}}{r+m})}{4\pi\mathcal{K}}\\\
&\times\left[F^{R}_{r,m}(E+\mathbf{p}\cdot\bm{\sigma})+2m\sigma_{z}\right],\end{split}$
(S49)

where

$F^{R}_{r,m}=1-\frac{r-m}{r+m}\left[\frac{(r+m)^{2}}{\mathcal{K}^{2}}-\frac{r+m}{\mathcal{K}\text{artanh}(\frac{\mathcal{K}}{r+m})}\right].$
(S50)

####  S2.7.1 Determinant of the GF

The spectrum of the zero modes can be found by taking the determinant of the
retarded GF and solving for 0, that is

$\begin{split}\det\mathbf{G}^{R}(\mathcal{K})&=\frac{3}{J^{2}}\frac{\text{artanh}^{2}[\mathcal{K}/(r+m)]}{(4\pi\mathcal{K})^{2}}\\\
&\times\left[(F^{R}_{r,m})^{2}(\mathcal{K})\mathcal{K}^{2}-4m^{2}\right]=0.\end{split}$
(S51)

This equation can be recasted in the condition

$\mathcal{K}=\frac{2m}{F_{r,m}(\mathcal{K})}.$ (S52)

We identify the solutions $\mathcal{K}_{0}$ to this equation with the restmass
of zeros $\mathcal{K}_{0}\equiv m_{\emptyset}$, since the spectrum of zeros
near K, K’ is given by

$E=\sqrt{p^{2}+m_{\emptyset}^{2}}.$ (S53)

Figure S1: Real and complex part of the zero’s energy mass gap as a function
of $r-m$.

####  S2.7.2 Approximate solution near $r=m$

We iteratively solve the equation $\det\mathbf{G}^{R}(\mathcal{K})=0$ near
$r=m$ and $\mathcal{K}=1$. Using units in which $r+m=1$ and $r-m\equiv\Delta$,
the defining equation becomes to leading order in $\Delta$

$\mathcal{K}=1+\frac{2\Delta}{\ln\left(\frac{1-\mathcal{K}}{2}\right)}=0.$
(S54)

We see that the term with a large log in the denominator is subleading, so we
solve the equation iteratively and we obtain

$\mathcal{K}=1+2\frac{\Delta}{\log\left(-\frac{\Delta}{\log\left(-\frac{\Delta}{\log\left(-\frac{\Delta}{\log\left(\dots\right)}\right)}\right)}\right)}.$
(S55)

The real and imaginary part of the GF zero mass as a function of $r-m$ at the
K point is plotted in Fig. S1. We can clearly see that the imaginary part of
the mass suddenly becomes finite once the mass of the orthogonal fermions
exceeds the e-particle mass, that is, $m>r$.

##  S3 Widom-Streda Formula

In this section we demonstrate an independent calculation of $\sigma_{H}$ (in
units of the von-Klitzing constant) using the Widom Streda formula

$\sigma_{H}=\left.\frac{2\pi}{e}\frac{\partial n}{\partial B}\right|_{B=0}.$
(S56)

where the density is calculated with the full propagator $\mathbf{G}$. To be
concrete, in the presence of a magnetic field it is convenient to express the
density using the full real space GF

$n[\mathbf{G}]=\frac{1}{\Omega}\int(\text{d}\epsilon)e^{i\epsilon
0^{+}}\sum_{\mathbf{r}}\tr[\mathbf{G}(i\epsilon;\mathbf{r},\mathbf{r})],$
(S57)

where $\Omega$ is the area of the system (total number of lattice sites). It
is emphasized in the main text that the physical density $n$ is evaluated
using $G$ in the present case, and the identification of $n[\mathbf{G}]$ with
the density is a priori not guaranteed.

Using this expression, it was shown in Refs. [38, 39] that

$\sigma_{H}[\mathbf{G}]=N_{3}[\mathbf{G}]+\frac{2\pi}{e}\left.\frac{\partial
I_{L}[\mathbf{G}]}{\partial B}\right|_{B=0},$ (S58)

where the Luttinger integral is

$I_{L}[\mathbf{G}]=\frac{1}{\Omega}\int(\text{d}\epsilon)\text{Tr}[\mathbf{G}(i\epsilon)\frac{\partial\Sigma(i\epsilon)}{i\partial\epsilon}].$
(S59)

Here, Tr includes the trace over internal degrees of freedom as well as
space/momentum space.

In the present model, a direct evaluation of $I_{L}$ is not natural, given
that $\mathbf{G}$ is defined as a convolution and
$\Sigma=G^{-1}-\mathbf{G}^{-1}$ is never explicitly evaluated. Instead, we
here directly calculate $\sigma_{H}[\mathbf{G}]$ using the Widom Streda
formula. Given that this result is independent from the previous section we
thereby implicitly demonstrate that Eq. (S58) holds for the case of our system
if $\sigma_{H}[\mathbf{G}]$ and $\sigma_{H}$ are identified.

We use for arbitrary eigenbasis
$[\hat{H}\psi_{n}](\mathbf{r})=E_{n}\psi_{n}(\mathbf{r})$ to evaluate Eq.
(S57) at $J=0$

$\displaystyle\mathbf{G}(i\epsilon;\mathbf{r},\mathbf{r}^{\prime})$
$\displaystyle=\sum_{n}\int(\text{d}\omega)D(i\omega)\sum_{n}\frac{\Psi_{n}(\mathbf{r})\Psi_{n}(\mathbf{r}^{\prime})}{i\epsilon-
E_{n}}$
$\displaystyle=\sum_{n}\frac{\Psi_{n}(\mathbf{r})\Psi_{n}(\mathbf{r}^{\prime})}{i\epsilon-\text{sign}(E_{n})[|E_{n}|+2h]},$
(S60)

so that

$\displaystyle n[\mathbf{G}]$
$\displaystyle=\frac{1}{\Omega}\int(\text{d}\epsilon)e^{i\epsilon
0^{+}}\sum_{n}\frac{1}{i\epsilon-\text{sign}(E_{n})[|E_{n}|+2h]}$ (S61)
$\displaystyle=\frac{1}{\Omega}\sum_{n}\Theta(-E_{n}).$ (S62)

Note that $n[\mathbf{G}]$ at $h/J=\infty$ and $n$ are thus equal.

Finally, for the evaluation of $\sigma_{H}$ we calculate the difference of
$\sigma_{H}$ at large and small $w_{2}/w_{0}$ and exploit that for $w_{2}=0$,
$\sigma_{H}[\mathbf{G}]=0$ by symmetry. We may readily use the continuum
approximation Eq. (2h) using canonical momenta $\Pi_{\mu}=p_{\mu}+eA_{\mu}$
and $m=w_{0}+\sqrt{27}w_{2}$, $\bar{m}=w_{0}-\sqrt{27}w_{2}$ We further use
for $e>0,B>0$ that

$\displaystyle b=\frac{1}{\sqrt{2eB}}i[\Pi_{x}-i\Pi_{y}],$ (S63)

where $b,b^{\dagger}$ are regular bosonic ladder operators. Thus, using the
notation $\ell_{B}=1/\sqrt{eB}$, we find

$\displaystyle h_{K}$
$\displaystyle=\left(\begin{array}[]{cc}m&-i\sqrt{2}\frac{c}{\ell_{B}}b\\\
-i\sqrt{2}\frac{c}{\ell_{B}}b^{\dagger}&-m\end{array}\right),$ (S66)
$\displaystyle h_{K^{\prime}}$
$\displaystyle=\left(\begin{array}[]{cc}\bar{m}&i\sqrt{2}\frac{c}{\ell_{B}}b^{\dagger}\\\
-i\sqrt{2}\frac{c}{\ell_{B}}b&-\bar{m}\end{array}\right).$ (S69)

The spectrum is ($\Omega_{c}^{2}=\sqrt{2}{c}/{\ell_{B}}$)

$\displaystyle E_{n}^{(K)}$ $\displaystyle=\begin{cases}-m,&n=0,\\\
\pm\sqrt{m^{2}+\Omega_{c}^{2}n^{2}},&n\neq 0,\end{cases}$ (S70) $\displaystyle
E_{n}^{(K^{\prime})}$ $\displaystyle=\begin{cases}\bar{m},&n=0,\\\
\pm\sqrt{\bar{m}^{2}+\Omega_{c}^{2}n^{2}},&n\neq 0.\end{cases}$ (S71)

Note that each Landau level has a degeneracy of $|eB|\Omega/2\pi$. At the
topological transition $\bar{m}$ changes sign and the zeroth Landau level of
the K’ point plunges under the Fermi level leading to an increase in density

$n_{\bar{m}<0}[\mathbf{G}]-n_{\bar{m}>0}[\mathbf{G}]=\frac{eB}{2\pi},$ (S72)

i.e. $\sigma_{H}[\mathbf{G}]=\Theta(-\bar{m})$.

##  S4 Green’s function based topological formula

In this section we consider the contribution of a given valley to the
quantized topological invariant, Eq. (1) of the main text,

$\displaystyle N_{3}$
$\displaystyle=\frac{\pi}{3}\int(\text{d}K)q[\mathbf{G}(\vec{K})],$ (S73a)
$\displaystyle q[\mathbf{G}(\vec{K})]$
$\displaystyle=\epsilon^{\mu\nu\rho}\tr\\{[\mathbf{G}^{-1}(\vec{K})\partial_{K_{\mu}}\mathbf{G}(\vec{K})]$
$\displaystyle\times[\mathbf{G}^{-1}(\vec{K})\partial_{K_{\nu}}\mathbf{G}(\vec{K})]$
$\displaystyle\times[\mathbf{G}^{-1}(\vec{K})\partial_{K_{\rho}}\mathbf{G}(\vec{K})]\\}.$
(S73b)

Here $(\text{d}K)=\text{d}^{2}kd\epsilon/(2\pi)^{3}$ and
$\vec{K}=(\epsilon,k_{x},k_{y})$, and $\mu,\nu,\rho\in\\{0,1,2\\}$, and the
totally antisymmetric Levi-Civita tensor $\epsilon^{\mu\nu\rho}$ is uniquely
defined by setting $\epsilon^{012}=1$.

Given GF (2i), we rewrite the GF in the form

$\displaystyle\mathbf{G}(\vec{K})$
$\displaystyle=\mathcal{G}(K)\tilde{G}(i\epsilon,\mathbf{k})$ (S74a)
$\displaystyle\tilde{G}(i\epsilon,\mathbf{k})$
$\displaystyle=i\sigma_{z}\underbrace{[\vec{K}\cdot(\sigma_{z},\sigma_{y},-\sigma_{x})^{T}+im(K)]}_{\equiv\tilde{g}(\vec{K})},$
(S74b) where
$\mathcal{G}(\vec{K})=-\frac{\sqrt{3}}{J}\Upsilon_{r+m}(K)F_{r+m}(K)$ (S74c)
and $m(K)=\frac{2m}{F_{r+m}(K)}$ (S74d)

Note that for constant $m\neq m(K)$, $\tilde{g}(\vec{K})$ has the same
topology as the non-fractionalized counterpart of our theory.

We next use that for an arbitrary product of GFs

$\begin{split}q[\mathbf{G}]&=q[\mathcal{G}\tilde{G}]=q[\mathcal{G}]+q[\tilde{G}]\\\
&+2\epsilon_{\mu\nu\rho}\partial_{\mu}\tr[(\mathcal{G}^{-1}\partial_{\nu}\mathcal{G})(\tilde{G}^{-1}\partial_{\rho}\tilde{G})]\\\
&=q[\tilde{G}]=q[\tilde{g}].\end{split}$ (S75)

Since $\mathcal{G}(K)$ is a function instead of a matrix, $q[\mathcal{G}]$
vanishes. The second equality follows because the surface term vanishes under
the integral sign. This is shown in the following. We start by noting

$\tr[(\mathcal{G}^{-1}\partial_{\nu}\mathcal{G})(\tilde{G}^{-1}\partial_{\rho}\tilde{G})]=(\mathcal{G}^{-1}\partial_{\nu}\mathcal{G})\tr[\tilde{g}^{-1}\partial_{\rho}\tilde{g}].$
(S76)

The factor
$(\mathcal{G}^{-1}\partial_{\nu}\mathcal{G})=(\mathcal{G}^{-1}\partial_{K}\mathcal{G})\hat{K}_{\nu}$
behaves asymptotically at infinity as $1/K$ and is proportional to
$\hat{K}=\vec{K}/K$. We next calculate (here
$\tilde{\bm{\sigma}}=(\sigma_{z},\sigma_{y},-\sigma_{x})$, note it has the
correct SU(2) algbra)

$\begin{split}\tilde{g}^{-1}\partial_{\rho}\tilde{g}&=\frac{\vec{K}\cdot\tilde{\bm{\sigma}}+im(K)}{K^{2}+m^{2}(K)}\partial_{\rho}[\vec{K}\cdot\tilde{\bm{\sigma}}-im(K)]\\\
&=\frac{\vec{K}\cdot\tilde{\bm{\sigma}}+im(K)}{K^{2}+m^{2}(K)}[\tilde{\sigma}_{\rho}-im^{\prime}(K)\hat{K}_{\rho}]\\\
&=-i\frac{M_{\rho\mu}}{K^{2}+m^{2}(K)}\tilde{\sigma}_{\mu}+h(K)\hat{K}_{\rho}\mathbf{1},\end{split}$
(S77)

where

$M_{\rho\mu}=\epsilon^{\rho\nu\mu}K_{\nu}+m(K)\delta_{\rho\mu}-\frac{K_{\rho}K_{\mu}}{K}m^{\prime}(K),$
(S78)

and

$h(K)=\frac{K+m(K)m^{\prime}(K)}{K^{2}+m^{2}(K)}.$ (S79)

In this equation, a ′ indicates a $K$-derivative. Only those terms
proportional to $\mathbf{1}$ survive under the trace (S76). Therefore, both
terms $(\mathcal{G}^{-1}\nabla\mathcal{G})\propto\vec{K}$, and
$\tr[\tilde{g}^{-1}\nabla\tilde{g}]\propto\vec{K}$ are proportional to
$\vec{K}$. Therefore, their cross product vanishes which ultimately proves the
validity of equation (S75).

Now, we can calculate $N_{3}$ according to equation (S73). Note that the part
of $\tilde{g}^{-1}\partial_{\mu}\tilde{g}$ which is proportional to
$\mathbf{1}$ immediately drops out of $q[\tilde{g}]$. Thus,

$\displaystyle N_{3}$
$\displaystyle=i\frac{\pi}{3}\int(\text{d}K)\epsilon_{\mu\nu\lambda}\frac{M_{\mu\rho}M_{\nu\sigma}M_{\lambda\tau}}{(\text{K}^{2}+m^{2}(K))^{3}}\tr[\sigma_{\rho}\sigma_{\sigma}\sigma_{\tau}]$
$\displaystyle=-\frac{2\pi}{3}\int(\text{d}K)\epsilon_{\mu\nu\lambda}\epsilon_{\rho\sigma\tau}\frac{M_{\mu\rho}M_{\nu\sigma}M_{\lambda\tau}}{(K^{2}+m^{2}(K))^{3}}$
$\displaystyle=-{4\pi}\int(\text{d}K)\frac{1}{(K^{2}+m^{2}(K))^{3}}$
$\displaystyle\times[M_{xx}M_{yy}M_{zz}+M_{xy}M_{yz}M_{zx}+M_{xz}M_{yx}M_{zy}$
$\displaystyle-M_{yx}M_{xy}M_{zz}-M_{yy}M_{xz}M_{zx}-M_{yz}M_{xx}M_{zy}]$
$\displaystyle=4\pi\int(\text{d}K)\frac{m(K)-Km^{\prime}(K)}{(m(K)^{2}+K^{2})^{2}}$
$\displaystyle=\frac{2}{\pi}\int_{0}^{\infty}\text{d}KK^{2}\frac{m(K)-Km^{\prime}(K)}{(m(K)^{2}+K^{2})^{2}}$
$\displaystyle=\text{sign}(m(K))\frac{2}{\pi}\int_{0}^{\infty}\text{d}xx^{2}\frac{1}{(1+x^{2})^{2}}$
$\displaystyle=\text{sign}(m(K))\frac{1}{2}$ (S80)

In the second to last equation, we introduced $x=K/|m(K)|$ (in our model,
$m(K)/m$ is positive definite). This proofs that $N_{3}$ is the same for
fractionalized and non-fractionalized continuum models.

##  S5 Momentum conserving tunneling

In this section we present details about the momentum conserving tunneling
proposed as a probe of zeros in the discussion section of the main text.

Consider a tunneling Hamiltonian

$H_{\rm
tun}=-t\sum_{l\in{t/b}}\int(dk)c^{\dagger}(\mathbf{k})p_{l}(\mathbf{k})+H.c.,$
(S81)

where $p_{l}(\mathbf{k})$ describe lead electrons with lead index $l=t/b$ for
the top or bottom lead. For simplicity we consider scalar Green’s functions in
this appendix – this explicitly includes the case of 1D planar tunneling
between quantum anomalous Hall edge states and zeros.

We want to determine the transition rate of charge between the two metallic
layers. To this purpose, we study the self-energy of electrons and electron-
multiplets in the top lead order by order in $t$, see Fig. S2, in which non-
interacting lead electrons from top/bottom lead are represented by purple/pink
and the exact Green’s function of electrons in the strongly central layer of
Fig. 3 d) are represented by thick boldface arrows. Most of the discussion in
this section is general, i.e. it is independent of the soluble model discussed
in the main text, yet the diagrammatic notation is only suggestive: generally
the interaction vertices entering Figs. S2 c), d) involve multiple pairs of
internal particle-hole lines (not drawn, but implied). As an illustration for
a non-trival effective interaction for the model Eq. (2), see the inset of
Fig. S2 c).

Figure S2: Contributions to the self-energy of the electrons in the upper lead
of Fig. 3 d). All explanations in the main text of the supplement.

By means of the optical theorem, cuts through the self-energies (illustrated
by dashed gray lines in Fig. S2) can be used to determine the transition rate
$\Gamma$ discussed in the main text.

First, consider the leading order self-energy, i.e. $\mathcal{O}(t^{2})$.
Clearly, the imaginary part of this diagram vanishes when there is a spectral
gap in the exact fermion Green’s function of the interacting system (the
transition rate for charge transfer vanishes).

Next, consider the leading diagrams to fourth order in $t$. The upper diagram
in Fig. S2 b) is not one-particle reducible – it is already taken into account
by considering the self-energy discussed in panel a). The imaginary part of
the lower diagram of Fig. S2 b) describes the rate of tunneling from top
metallic lead to bottom metallic lead through virtual excitations of the
central correlated system. Since the imaginary part of the Green’s function of
the central layer vanishes, the transition rate is given by
$\Gamma(\epsilon_{\mathbf{p}},\mathbf{p})=-t^{4}[\text{Re}\mathbf{G}^{R}(\epsilon_{\mathbf{p}}^{\rm
lead},\mathbf{p})]^{2}\text{Im}G^{R}_{\rm lead}(\epsilon^{\rm
lead}_{\mathbf{p}},\mathbf{p})$, as quoted in the main text.

Next, we consider inelastic tunneling and co-tunneling, which are of
particular interest for the fractionalized theory associated with Eq. (2),
because certain limits of four-point functions of fermions behave the same as
in Fermi liquids. The co-tunneling process Fig. S2 c) describes a charge 2
transfer between the two leads, but is subleading in orders of $t$. The same
is true for the exemplary inelastic process creating a particle-hole pair in
the lower lead.

Inelastic processes which are of the same order as the elastic process
discussed in Fig. S2 c) involve creating an excitation in the correlated
layer. Fig. S2 e) illustrates such a process leaving an electron-hole pair in
the correlated layer (this generally can include spin-excitations). Whether
such particle-hole pairs are gapped or gapless depends on details of the
model.

To illustrate this, the inelastic processes allowed to $\mathcal{O}(t^{4})$
within the model, Eq. (2), of this letter are presented in Fig. S2 f) g). The
excitation presented in panel f) is gapped: it requires the energy of the
incoming particle (the voltage) to be higher than twice the energy gap of
e-particles. Such a process is thus not contributing at the lowest energies
(of the order of the bands of dispersing Green’s function zeros). In contrast,
Fig. S2 g) illustrates a process where a particle-hole pair of “orthogonal
fermions” is left behind, which allows for low energy inelastic charge
tunneling.

In the bulk, for OTI and ZTI phases one may expect that the excitation energy
of this particle-hole exceeds the (topological) gap of zeros. This implies
that the bulk probe of zeros perisists, i.e. the dip in the tunneling
conductance is still clearly observable.

At the boundary ZTI-OTI we expect boundary zeros but no edge states (poles) at
the lowest energy. Thus the inelastic processes Fig. S2 f), g) are gapped at
lowest energies while boundary zeros show up in momentum-conserving tunneling.
In such a situation(Green’s function boundary zeros are present, and boundary
poles are absent), momentum-conserving tunneling probes of topological edge
zeros should work best.

Finally, at the boundary ZTI-OI we expect both zeros and edge states. We
expect that inelastic processes creating particle-hole pairs in the edge
display characteristic temperature dependence which can be discriminated from
the elastic contribution discussed in the main text and leave their study to a
separate work.