# Spin and Orbital Metallic Magnetism in Rhombohedral Trilayer Graphene

Chunli Huang  Theoretical Division, T-4, Los Alamos National Laboratory, Los
Alamos, New Mexico 87545, USA  Tobias M. R. Wolf  Wei Qin  Nemin Wei  Igor V.
Blinov  Allan H. MacDonald  Department of Physics, University of Texas at
Austin, Austin TX 78712

(August 27, 2024)

###### Abstract

We provide a complete theoretical interpretation of the metallic broken
spin/valley symmetry states recently discovered in ABC trilayer graphene (ABC)
perturbed by a large transverse displacement field. Our conclusions combine
insights from ABC trilayer graphene electronic structure models and mean field
theory, and are guided by precise magneto-oscillation Fermi-surface-area
measurements. We conclude that the physics of ABC trilayer graphene is shaped
by the principle of momentum-space condensation, which favors Fermi surface
reconstructions enabled by broken spin/valley flavor symmetries when the
single-particle bands imply thin annular Fermi seas. We find one large outer
Fermi surface enclosed majority-flavor states and one or more small inner
hole-like Fermi surfaces enclosed minority-flavor states that are primarily
responsible for nematic order. The smaller surfaces can rotate along a ring of
van-Hove singularities or reconstruct into multiple Fermi surfaces with little
cost in energy. We propose that the latter property is responsible for the
quantum oscillation frequency fractionalization seen experimentally in some
regions of the carrier-density/displacement-field phase diagram.

Figure 1: Zero temperature phase diagram of ABC trilayer graphene. Fermi
surface topologies in the various phases are indicated schematically. The
phase boundaries mark magnetic transitions (MT) and Lifshitz transitions, both
of which change Fermi surface topology. We distinguish two types of Lifshitz
transitions, ALTs (annular Lifshitz transitions) that occur when the Fermi
level crosses $\bm{k}=0$ energy band maxima to form an electron pocket inside
an annular hole Fermi sea, and SPLTs (saddle point Lifshitz transitions) that
occur when the Fermi level crosses $\bm{k}\neq 0$ energy band saddle points to
break the annulus into three separate pockets. At the first MT the
paramagnetic annular Fermi sea supports a pair of Shubnikov-de Haas quantum
oscillations with normalized frequencies $\bm{f}_{\nu}^{*}=(0.06,0.31)$. The
interlayer potential $U$ is proportional to the displacement field $D$ in the
experiment. Superconductivity occurs in region of the phase diagram shaded
yellow close to the high hole density first MT. We used the dielectric
constant $\epsilon_{r}=4$.

Introduction:– Multilayer graphene systems continue to surprise researchers
with unexpected states of matter. Zhou. et.al [1, 2, 3] have recently
uncovered rich phase diagrams in both AB bilayer and rhombohedral ABC trilayer
graphene containing half-metal, quarter-metal, partially-isopsin polarized
(PIP) metal, spin-triplet superconductor, and spin-singlet superconductor
states in a three dimensional parameter space spanned by the electric
displacement field $D$, the carrier density $n_{e}$, and the magnetic field
$B$, see also Ref. [4, 5]. The origin of superconductivity has been studied in
Refs. [6, 7, 8, 9, 10]. In this Letter, we explain the pattern of spontaneous
symmetry breaking and Fermi surface reconstructions observed in the normal
metal states of ABC trilayer graphene, and comment on what they say about
electron-interaction physics in two-dimensional multilayer graphene electron
gases.

Figure 2: Mean-field valence bands and Fermi surfaces of competing states at
$n_{e}=-1.4\cdot 10^{12}\,\text{cm}^{-2}$ and $U=60$ meV (yellow star in Fig.
1): (a) the paramagnetic metal, (b–c) the valley-XY nematic metal and (d) the
spin nematic metal. The yellow, green, and blue regions are respectively
occupied by zero, one, and two holes. In the ordered region all magnetic
states gain roughly the same amount of total energy relative to the
paramagnetic state, ($\Delta E_{\text{tot}}^{\text{p.h.}}\sim 200\mu$eV for
the illustrated case) by lowering exchange energy $\Delta
E_{\text{ex}}^{\text{p.h.}}\sim 1500\mu$eV per hole at a cost in band energy.
Their energy differences between the different ordered states ($\sim 1\to
20\mu$eV per hole) is mainly determined by the characteristics of the inner
nematic Fermi sea. The magnetic anisotropic energy in the $\mathrm{SU}(4)$
spin–valley space first minimizes the area of the inner hole Fermi sea,
measured by the quantum oscillation frequencies $\bm{f}$, and then breaks
rotational symmetry to make the inner surface more compact. The arrows in (c)
suggest the sliding degree-of-freedom of the crescent-shaped Fermi sea. This
mean-field calculation was performed with an accuracy of $\lesssim 1\mu$eV per
hole, using $\epsilon_{r}=4.5$ and $D=50$nm.

At low carrier densities $|n_{e}|\lesssim 2\times 10^{12}\text{cm}^{-2}$ and
large electric displacement fields $D\gtrsim 0.2\text{V/nm}$, holes
($n_{e}<0$) in ABC graphene occupy low velocity Bloch states [11, 12] with
momenta near one of two inequivalent honeycomb lattice Brillouin-zone corners.
These regions of momentum space are referred to below as valleys, and endow
electrons with a valley pseudospin in addition to spin. States in the two
valleys are related by time-reversal symmetry. Because the spin-
($s=\uparrow,\downarrow$) and valley ($\tau=K,K^{\prime}$) projected densities
of states are identical for each spin/valley, the paramagnetic state occupies
the four flavors equally. The broken spin/valley symmetry states seen
experimentally occupy the four flavors unequally, just as magnetic metals
occupy majority and minority spins unequally, and are expected because the
Fermi energy at low carrier densities is small compared to the Coulomb energy
per particle. Because of the high quality of the ABC graphene devices studied
by Zhou et al., it has been possible to measure magnetic oscillations at weak
magnetic fields, and in this way to accurately measure the areas enclosed by
most of the Fermi surfaces present in most phases. Magneto-oscillations do not
on their own distinguish different symmetry-broken states and do not provide
the information about Fermi surface shapes. These we infer by combining
insights from ABC trilayer graphene electronic structure models and mean field
theory calculations [13] that are corroborated by features that appear in
resistance vs. $n_{e}$ and $D$ maps, in addition to the thermodynamic
compressibility data [1, 2, 3]. We detail our band structure model and mean-
field calculations in the Supplemental Material [13].

We seek to explain the evolution of the ground state as the hole density
decreases and interactions increase in importance. As shown in Fig. 1, at
large interlayer potential $U$ the first transition is an annular Liftshitz
transition (ALT) at which an electron pocket emerges in the interior of the
hole Fermi sea. This annular Fermi sea undergoes a first magnetic transition
(1st MT) to a valley ordered state accompanied by a strong Fermi surface
reconstruction; i.e., the number of Fermi surfaces (per spin) is reduced from
four to two and the annulus near $k=0$ is removed. This state is susceptible
to nematic broken symmetries that give rise to small crescent-shaped Fermi
surfaces, which explain the small quantum oscillation frequencies reported in
Ref. [1, 3]. We identify such a valley ordered state as a valley-XY nematic
metal. 111Note that other works sometimes refer to valley-XY ordered states as
‘intervalley coherent states‘. When the hole density is further lowered, the
crescent-shaped (nematic) Fermi liquid shrinks and leads to a negative
electronic compressibility as reported in Fig. 1d of Ref. [1]. The crescent-
shaped Fermi liquid disappears at a second magnetic transition (2nd MT) to a
spin-polarized half-metal. As the hole densities is lowered even further, a
second ALT changes the topology of the spin-polarized half-metal to an annular
shape. Both the paramagnetic phase ALT (blue line) and the spin-polarized
state ALT are prominent experimentally because they lead to easily identified
peaks in resistance along the lines from $(n_{e}\,,\,D)=(-1.2,0)$ to
$(-1.8,0.28)$ and from $(-0.6,0)$ to $(-1.1,0.5)$ in Fig. 1e of Ref. [1],
respectively. At small $U$, there is only one magnetic transition because the
paramagnetic metal avoids the transition into a valley-XY nematic metal. In
what follows, we first address the instability of a paramagnetic metal with an
annular Fermi sea [15], then the properties of magnetic metals.

Annular Fermi sea:– Depending on $n_{e}$ and $U$, the single-particle bands
can have one of three distinct Fermi surface topologies (see Fig. 1): a
simply-connected Fermi sea, a Fermi sea with three disjoint pockets, and an
annular Fermi sea. With decreasing hole density, the annular topology appears
at the ALT when the Fermi level intersects a shallow local valence band
minimum at $\bm{k}=0$, and persists to the saddle-point Lifshitz transition
(SPLT) at which the annulus splits into three pockets and the density-of-
states diverges logarithmically. We focus below on the experimentally
demonstrated instability of the $C_{3}$-distorted annular Fermi sea, which
occurs well before the SPLT is reached and is therefore not simply because of
the large density-of-states at the Fermi level.

The annular Fermi sea has an exterior hole-like Fermi surface and an interior
electron-like Fermi surface, as shown in Fig. 2(a) 222We overlay the Fermi
surfaces centered at $K,K^{\prime}$ point to the $\Gamma$ point for better
comparison to the symmetry broken states.. The enclosed areas
$A_{FS,\text{in}}$ and $A_{FS,\text{out}}$ lead to a pair of normalized
quantum oscillation frequencies
$\bm{f}_{\nu}=(f_{\nu,\text{in}}\,,\,f_{\nu,\text{out}})$:

$f_{\nu,\text{in}}=\frac{A_{FS,\text{in}}}{\sum_{\alpha}A_{h,\alpha}}\;,\;f_{\nu,\text{out}}=\frac{A_{FS,\text{out}}}{\sum_{\alpha}A_{h,\alpha}},$
(1)

where the sum in the denominator is over the flavor index
$\alpha=\\{K\uparrow,K\downarrow,K^{\prime}\uparrow,K^{\prime}\downarrow\\}$.
The paramagnetic metal partitions holes equally into the four flavors so that

$A_{h,\alpha}\equiv A_{h}^{0}=\frac{(2\pi)^{2}n_{e}}{4}\;\;\forall\,\alpha.$
(2)

For an annular Fermi sea, $A_{h}^{0}=A_{FS,\text{out}}-A_{FS,\text{in}}$.

Valley-XY nematic metal:– As the inner Fermi surface of the paramagnet becomes
larger, occupied hole momenta are more widely dispersed and this leads to a
reduction in the long-range Coulomb exchange energy [15]. We identify this
reduction as the driver of the first magnetic transition, which occurs near
$A_{FS,\text{in}}/A_{h}^{0}\approx 0.3$ and therefore with normalized quantum
oscillation frequencies

$\bm{f}_{\nu}^{*}\sim(0.06,0.31).$ (3)

The fact that these frequencies are nearly independent of $U$ emphasizes that
Fermi sea shape is more important in driving the instability than hole
density.

The 1st MT leads to strong flavor polarization and Fermi surface
reconstructions where the majority-flavor states have a large hole-like Fermi
surface with no annulus. Since the outer Fermi surface $A_{FS,\text{out}}$ has
a large Fermi velocity, it is energetically favorable to retain a small
minority-flavor hole pocket centered at $k\neq 0$ to limit its expansion. For
a given $n_{e}$, the exchange energy favors the valley-XY state over the spin-
polarized and valley-Ising states illustrated in Fig. 2 because the outer
Fermi surface of its valley-mixed quasiparticles has $C_{6}$ symmetry and can
therefore enclose a larger area for given Fermi radii. The relationship
between the spin–valley $\mathrm{SU}(4)$ magnetic anisotropy energy and the
area of the inner Fermi surface is illustrated in Fig. 2(b,c) for valley-XY
polarized states and in Fig. 2(d) for a spin-polarized state. The 1st MT leads
to the discontinuous changes in quantum oscillation frequencies from the
paramagnetic state values $\bm{f}_{\nu}^{*}$ to the ones shown in Fig. 2.
Although the candidate ordered states differ little in energy, they do have
distinct quantum oscillation frequencies. The values reported in Ref. [3]
support our identification of the valley-XY state as the ground state in this
region of the phase diagram. Note that the quantum oscillation data put a
severe constraint on the many possible magnetic states in the large
$\mathrm{SU}(4)$ Hilbert space, see Ref. [13]. In particular, an intervalley
coherent state (valley-XY) without nematicity cannot explain the sudden jump
in the oscillation frequencies at the phase-boundary.

We now explain the crescent shape of the reconstructed inner Fermi sea of the
valley-XY metal shown in Fig. 2(c). Note first that the valence bands have a
ring of maxima at $\bm{k}_{\text{VBM}}\sim k_{0}(\cos(\theta),\sin(\theta))$
where $k_{0}a\sim 0.05$ and $2\pi>\theta\geq 0$. When a point on the ring
rises above the Fermi level it nucleates a hole-like Fermi sea. Such Fermi sea
nucleation can be understood as a form of momentum-space condensation
introduced by Heisenberg [17]. Because the energy dispersion along the radial
direction is steep, due mainly to strong non-local exchange-splitting, whereas
the dispersion along the ring of valence band maximum is anomalously flat, the
nucleated Fermi sea has the form of a crescent. Furthermore, the crescent
Fermi sea wants to be as compact as possible to increase the exchange energy
from the long-range Coulomb interaction. Thus, although the reconstructed
inner Fermi seas in Fig. 2(b,c) have the same area, their total energy differ
because of the relative compactness of their Fermi seas. The single-crescent
state is the most compact and hence has the lowest energy. The crescent Fermi
sea leads to nematicity, which we characterize by the $2\times 2$ traceless
symmetric tensor

$\displaystyle\hat{N}_{ij}=\hat{D}_{ij}-\frac{1}{2}\text{tr}(\hat{D})\delta_{ij},$
(4)

where $i,j\in\\{x,y\\}$ and $\hat{D}$ is the Drude tensor given by

$\displaystyle\hat{D}_{ij}=\sum_{a=1,2...N_{F}}\int\left(\frac{v_{a,i}v_{a,j}}{|\bm{v}_{a}|}\right)\frac{dk_{\parallel,a}}{4\pi^{2}\hbar}.$
(5)

Here $\bm{v}_{a}(\bm{k})=\nabla\epsilon_{a}(\bm{k})$ is the Fermi velocity and
$a=1,2,\dots,N_{F}$ enumerates the Fermi surfaces. The eigenvectors of
$\hat{N}$ are the directors of the nematic Fermi sea, and the eigenvalues
$\pm\sqrt{\text{det}(\hat{N})}$ defines a nematic order parameter

$\mathcal{N}=\sqrt{|\text{det}(\hat{D})|}/\text{tr}(\hat{D}).$ (6)

The values for the single crescent and double crescent states are
$\mathcal{N}=4\%$ and $\mathcal{N}=6\%$, respectively. This nematicity can
lead to significant transport anisotropy below the transition temperature,
that we predict will be most prominent close to the phase boundary of the 1st
MT.

The crescent Fermi sea is very fragile in the sense that it costs very little
energy to deform and slide it around the ring of valence band maxima at
$\bm{k}_{\text{VBM}}$ as indicated in Fig. 2(c). We can control the center of
mass of the crescent Fermi sea on $\bm{k}_{\text{VBM}}$ by minimising
$\tilde{H}=H^{\text{MF}}-\mu_{\text{nem}}\,\rho(\bm{k}_{\text{VBM}})$ where
$\mu_{\text{nem}}$ is a Lagrange multiplier, $H^{\text{MF}}$ is the mean-field
Hamiltonian, and $\rho(\bm{k}_{\text{VBM}})$ is the density matrix evaluated
at $\bm{k}_{\text{VBM}}$ [13]. The pinning potential which favors particular
crescent orientations is $\lesssim 1$ $\mu$eV per hole, and cannot be resolved
by our finite $k$-grid mean-field calculations 333We used a triangular
$k$-grid with $140\times 140$ point-resolution in regions up to radii
$k_{\text{max}}a=0.2$. When the crescent Fermi sea becomes small, it costs
very little energy to split it into two or more smaller crescents. We
attribute the fractionalization of the small quantum oscillation frequencies
reported in Ref. [3] as the system advances deeper into the valley-XY state,
to this fragility.

Figure 3: Transition from a nematic valley-XY phase to the valley-XY half-
metal and to the spin-polarized half-metal as a function of decreasing hole
density $n_{e}$ at fixed interlayer potential $U=62$ meV. (a) The evolution of
the Fermi sea for valley-XY and spin ordered states, where the ground states
are shaded light blue. The valley-XY ordered state undergoes a Lifshitz
transition (LT) upon which it loses its inner pockets before going through the
2nd magnetic transition (MT). (b) Quantum oscillation frequencies $f_{\nu}$
(big markers) and chemical potential $\mu$ as the electron density increases.
This calculation was performed with model parameters $\epsilon_{r}=4.0$,
$d_{\text{gate}}=50$ nm, and a $k$-mesh with $N_{k}=10^{4}$ points.

2nd magnetic transition:– As the hole density decreases in the valley-XY metal
phase, the paramagnetic region enclosed by the crescent Fermi sea shrinks as
shown in Fig. 3(a) and increases the valley-XY order parameter. This leads a
decrease of Fermi energy ($\mu$) as shown in Fig. 3(b) and explains the
observed (slightly) negative inverse compressibility
$\kappa=\partial\mu/\partial n_{e}$ in the PIP phase of Ref. [1]. When the
nematic Fermi sea disappears completely, $\mathcal{N}=0$ and
$\bm{f}_{\nu}=(0,0.5)$. The valley-XY half-metal is then characterized by a
single compact Fermi surface with $C_{6}$ rotational symmetry (for each spin).
Although the $C_{6}$ Fermi surface is closer to a circle, compared to the
spin-polarized half-metal which only has $C_{3}$, the content of the
wavefunction in the valley-XY metal changes dramatically over the enclosed
area. The valley-XY half-metal eventually undergoes a 2nd magnetic transition
to a spin-polarized half-metal, as shown in Fig. 1 and indicated by the arrow
in Fig. 3(a). Although this is a first-order phase transition, there is no
change in $\bm{f}$ because both states are half-metals, and no change in
$\mathcal{N}$ because both states lack nematicity. Their phase boundary,
however, can easily be shifted by an in-plane spin-polarizing Zeeman field [1,
2, 3]. The transition to the spin-polarized half-metal can be understood as a
reset transition [1, 19, 20, 21], at which the minority spin Hilbert space
becomes inactive while there is still equipartition of density in the majority
spin Hilbert space:

$A_{h,K\uparrow}=A_{h,K^{\prime}\uparrow}=\frac{(2\pi)^{2}n_{e}}{2}\;\;,\;\;A_{h,K\downarrow}=A_{h,K^{\prime}\downarrow}=0.$
(7)

As the hole density is lowered, a Lifshitz transition creates an annular
topology in the (majority) $\uparrow$-spin Fermi sea and the same symmetry-
breaking pattern repeats in the $\uparrow$-spin Hilbert space as the density
of holes is further lowered.

Spin-fluctuation and spin-orbit interactions:– At weaker displacement fields,
the Fermi velocity of the outer Fermi surface is smaller at the first magnetic
phase transition than it is at large $U,|n_{e}|$, especially where it is close
to the three SPLT momenta. It is hence favorable to make a direct transition
from the paramagnetic state to the spin-polarized half-metal 444The
groundstate energy of a spin-polarized half metal and a valley-polarized half-
metal appear to be almost degenerate due to the $SU(4)$ invariant Coulomb
interaction. Proximity to a transition between a paramagnetic metal and a
ferromagnetic metal is known to suppress spin singlet superconductivity [23]
and valley-singlet superconductivity [24], explaining the termination of
singlet superconductivity at small $D$ seen in Ref. [3], which occurs despite
an increasing density-of-states. The direct transition from the paramagnetic
state to the spin-polarized half metal leads to a sudden jump of quantum
oscillation frequencies directly from $\bm{f}_{\nu}^{*}$ to
$\bm{f}_{\nu}=(0,0.5)$. In our mean-field calculations, this direct transition
continues to zero displacement field ($U=0$), as shown in Fig. 1. In
experiment, however, the magnetic phase fades away at small displacement field
and blends into the paramagnetic phase. We attribute this behavior to spin
fluctuations that become important when the electric displacement field is
small. The displacement field can suppress spin fluctuations because it
increases the strength of spin-orbit coupling [25, 26, 27], providing the
magnetic anisotropy needed for spin-order.

Summary and Discussions– The magnetic order recently discovered [1, 2, 3] in
metallic ABC trilayer graphene is unusual because it involves valley as well
as spin degrees of freedom, and because the quasiparticle bands imply that the
unordered paramagnetic states have annular Fermi seas. By combining mean-field
theory with a model that captures pertinent electronic structure details, we
are able to account for the phase transitions observed experimentally and shed
light on the nature of the different phases. In particular we identify the
experimental partially-isospin-polarized (PIP) state as a valley-XY nematic
metal. We find that, as in conventional itinerant electron magnets, the
magnetic condensation energy greatly exceeds a much smaller magnetic
anisotropy energy. The $\mathrm{SU}(4)$ spin–valley magnetic anisotropy energy
controls a competition between valley and spin-polarized states, which are
very close in energy throughout the phase diagram. We predict that the valley-
XY metal, which is stabilized at large interlayer potential, is nematic
because small minority-hole pockets, driven by the energetic expense of
expanding the outer Fermi surface, break rotational symmetries in order to
become more compact. The nematicity of these states should be observable in
transport measurements. When the minority-holes are absent at low hole
densities and at low displacement fields, spin order is preferred over valley
order. In closing, let us briefly mention the possibility for the small
nematic Fermi sea to crystallize into a Wigner crystal state that breaks both
translation and $C_{3}$ symmetry. In the spirit of Ref. [28], at sufficiently
low density of holes, such anisotropic Wigner crystal states can have lower
energy than the nematic fluid and conventional hexagonal Wigner crystal states
due to the ring of Van Hove singularities in k-space.

Note added: Recently, we found that Ref. [29] discusses momentum condensation
(“flocking” effect) in Bernal bilayer graphene with three disjoint Fermi
pockets whose microscopic origin is similar to the formation of nematic Fermi
surface we discuss here in ABC trilayer graphene.

###### Acknowledgements.

Acknowledgments:– We thank Anna Seiler, Andrea Young and Haoxin Zhou for
informative discussions. The work done at LANL was carried out under the
auspices of the US DOE NNSA under Contract No. 89233218CNA000001 through the
LDRD Program. This work was supported by the U.S. Department of Energy, Office
of Science, Basic Energy Sciences, under Award DE-SC0022106. T.M.R.W. is
grateful for the financial support from the Swiss National Science Foundation
(Postdoc.Mobility Grant No. 203152).

## References

  * Zhou  _et al._ [2021a] H. Zhou, T. Xie, A. Ghazaryan, T. Holder, J. R. Ehrets, E. M. Spanton, T. Taniguchi, K. Watanabe, E. Berg, M. Serbyn, and A. F. Young, Half and quarter metals in rhombohedral trilayer graphene, Nature 598, 429 (2021a), arXiv:2104.00653 .
  * Zhou  _et al._ [2021b] H. Zhou, Y. Saito, L. Cohen, W. Huynh, C. L. Patterson, F. Yang, T. Taniguchi, K. Watanabe, and A. F. Young, Isospin magnetism and spin-triplet superconductivity in Bernal bilayer graphene, (2021b), arXiv:2110.11317 [cond-mat] .
  * Zhou  _et al._ [2021c] H. Zhou, T. Xie, T. Taniguchi, K. Watanabe, and A. F. Young, Superconductivity in rhombohedral trilayer graphene, Nature 598, 434 (2021c), arXiv:2106.07640 .
  * Seiler  _et al._ [2021] A. M. Seiler, F. R. Geisenhof, F. Winterer, K. Watanabe, T. Taniguchi, T. Xu, F. Zhang, and R. T. Weitz, Quantum cascade of new correlated phases in trigonally warped bilayer graphene, arXiv preprint arXiv:2111.06413 (2021).
  * de la Barrera  _et al._ [2021] S. C. de la Barrera, S. Aronson, Z. Zheng, K. Watanabe, T. Taniguchi, Q. Ma, P. Jarillo-Herrero, and R. Ashoori, Cascade of isospin phase transitions in Bernal bilayer graphene at zero magnetic field, arXiv e-prints , arXiv:2110.13907 (2021), arXiv:2110.13907 [cond-mat.mes-hall] .
  * Chatterjee  _et al._ [2021] S. Chatterjee, T. Wang, E. Berg, and M. P. Zaletel, Inter-valley coherent order and isospin fluctuation mediated superconductivity in rhombohedral trilayer graphene, (2021), arXiv:2109.00002 [cond-mat] .
  * Ghazaryan  _et al._ [2021] A. Ghazaryan, T. Holder, M. Serbyn, and E. Berg, Unconventional superconductivity in systems with annular fermi surfaces: Application to rhombohedral trilayer graphene, Physical review letters 127, 247001 (2021).
  * You and Vishwanath [2021] Y.-Z. You and A. Vishwanath, Kohn-Luttinger Superconductivity and Inter-Valley Coherence in Rhombohedral Trilayer Graphene, arXiv e-prints , arXiv:2109.04669 (2021), arXiv:2109.04669 [cond-mat.str-el] .
  * Cea  _et al._ [2022] T. Cea, P. A. Pantaleón, V. T. Phong, and F. Guinea, Superconductivity from repulsive interactions in rhombohedral trilayer graphene: A Kohn-Luttinger-like mechanism, Phys. Rev. B 105, 075432 (2022), arXiv:2109.04345 [cond-mat.mes-hall] .
  * Qin  _et al._ [2022] W. Qin, C. Huang, T. M. R. Wolf, N. Wei, I. V. Blinov, and A. H. MacDonald, Functional Renormalization Group Study of Superconductivity in Rhombohedral Trilayer Graphene, arXiv e-prints , arXiv:2203.09083 (2022), arXiv:2203.09083 [cond-mat.supr-con] .
  * Zhang  _et al._ [2010] F. Zhang, B. Sahu, H. Min, and A. H. MacDonald, Band structure of A B C -stacked graphene trilayers, Phys. Rev. B 82, 035409 (2010).
  * McCann and Koshino [2013] E. McCann and M. Koshino, The electronic properties of bilayer graphene, Rep. Prog. Phys. 76, 056503 (2013), arXiv:1205.6953 .
  * [13] See Supplemental Material for details on the band structure model, the Coulomb interaction potential, the self-consistent Hartree-Fock approximation, and a discussion of intervalley coherent states without nematicity.
  * Note [1] Note that other works sometimes refer to valley-XY ordered states as ‘intervalley coherent states‘.
  * Jung  _et al._ [2015] J. Jung, M. Polini, and A. H. MacDonald, Persistent current states in bilayer graphene, Physical Review B 91, 155423 (2015).
  * Note [2] We overlay the Fermi surfaces centered at $K,K^{\prime}$ point to the $\Gamma$ point for better comparison to the symmetry broken states.
  * Heisenberg [1984] W. Heisenberg, Two lectures, in  _Scientific Review Papers, Talks, and Books Wissenschaftliche Übersichtsartikel, Vorträge und Bücher_, edited by W. Blum, H.-P. Dürr, and H. Rechenberg (Springer Berlin Heidelberg, Berlin, Heidelberg, 1984) pp. 443–455.
  * Note [3] We used a triangular $k$-grid with $140\times 140$ point-resolution in regions up to radii $k_{\text{max}}a=0.2$.
  * Zondiner  _et al._ [2020] U. Zondiner, A. Rozen, D. Rodan-Legrain, Y. Cao, R. Queiroz, T. Taniguchi, K. Watanabe, Y. Oreg, F. von Oppen, A. Stern, _et al._ , Cascade of phase transitions and dirac revivals in magic-angle graphene, Nature 582, 203 (2020).
  * Saito  _et al._ [2021] Y. Saito, F. Yang, J. Ge, X. Liu, T. Taniguchi, K. Watanabe, J. Li, E. Berg, and A. F. Young, Isospin pomeranchuk effect in twisted bilayer graphene, Nature 592, 220 (2021).
  * Wong  _et al._ [2020] D. Wong, K. P. Nuckolls, M. Oh, B. Lian, Y. Xie, S. Jeon, K. Watanabe, T. Taniguchi, B. A. Bernevig, and A. Yazdani, Cascade of electronic transitions in magic-angle twisted bilayer graphene, Nature 582, 198 (2020).
  * Note [4] The groundstate energy of a spin-polarized half metal and a valley-polarized half-metal appear to be almost degenerate due to the $SU(4)$ invariant Coulomb interaction.
  * Berk and Schrieffer [1966] N. Berk and J. Schrieffer, Effect of ferromagnetic spin correlations on superconductivity, Physical Review Letters 17, 433 (1966).
  * Huang  _et al._ [2021] C. Huang, N. Wei, W. Qin, and A. MacDonald, Pseudospin paramagnons and the superconducting dome in magic angle twisted bilayer graphene, arXiv preprint arXiv:2110.13351 (2021).
  * Konschuh  _et al._ [2012] S. Konschuh, M. Gmitra, D. Kochan, and J. Fabian, Theory of spin-orbit coupling in bilayer graphene, Physical Review B 85, 115423 (2012).
  * Huertas-Hernando  _et al._ [2006] D. Huertas-Hernando, F. Guinea, and A. Brataas, Spin-orbit coupling in curved graphene, fullerenes, nanotubes, and nanotube caps, Physical Review B 74, 155426 (2006).
  * Min  _et al._ [2006] H. Min, J. Hill, N. A. Sinitsyn, B. Sahu, L. Kleinman, and A. H. MacDonald, Intrinsic and rashba spin-orbit interactions in graphene sheets, Physical Review B 74, 165310 (2006).
  * Berg  _et al._ [2012] E. Berg, M. S. Rudner, and S. A. Kivelson, Electronic liquid crystalline phases in a spin-orbit coupled two-dimensional electron gas, Phys. Rev. B 85, 035116 (2012).
  * Dong  _et al._ [2021] Z. Dong, M. Davydova, O. Ogunnaike, and L. Levitov, Isospin ferromagnetism and momentum polarization in bilayer graphene, arXiv e-prints , arXiv:2110.15254 (2021), arXiv:2110.15254 [cond-mat.mes-hall] .
  * Koshino and McCann [2009] M. Koshino and E. McCann, Trigonal warping and Berry’s phase $n\pi$ in ABC-stacked multilayer graphene, Phys. Rev. B 80, 165409 (2009).

Supplementary Material: Pseudospin Paramagnons and the Superconducting Dome in
Magic Angle Twisted Bilayer Graphene

Supplemental Material:
Spin and Orbital Metallic Magnetism in Rhombohedral Trilayer Graphene

##  I Band structure model

Rhombohedral trilayer graphene (RTG) has three layers $\ell=1,2,3$ separated
by interlayer distance $d=3.4$ Å and two lattices $\sigma=A,B$ each, such that
the unit cell of the triangular lattice with lattice constant of $a=2.46$ Å
contains $6$ atoms. The layers of RTG are in a ABC stacking configuration.
Following Refs. [11, 30], we consider a single electronic orbital ($\pi$) per
atom, with spin $s=\pm 1/2$, and use the continuum model that can be derived
from a tight-binding Hamiltonian by expanding the low-energy dispersion around
each valley $\tau K$ with $\tau=\pm 1$ at the corners of the Brillouin zone.
We label the four flavor combinations of spin $s$ and valley $\tau$ with the
flavor index $\alpha$. In the basis $\psi_{\alpha}(\bm{k})=(\psi_{\alpha
A1}(\bm{k}),\psi_{\alpha B1}(\bm{k}),\psi_{\alpha A2}(\bm{k}),\psi_{\alpha
B2}(\bm{k}),\psi_{\alpha A3}(\bm{k}),\psi_{\alpha B3}(\bm{k}))$, the resulting
continuum Hamiltonian for given spin–valley $\alpha$ is [11]

$\displaystyle
h(\bm{k})=\begin{bmatrix}t(\bm{k})+U_{1}&t_{12}(\bm{k})&t_{13}\\\
t_{12}^{\dagger}(\bm{k})&t(\bm{k})+U_{2}&t_{12}(\bm{k})\\\
t_{13}^{\dagger}&t_{12}^{\dagger}(\bm{k})&t(\bm{k})+U_{3}\end{bmatrix}_{6\times
6},$ (S1)

where $\bm{k}=(k_{x},k_{y})$ is the Bloch momentum measured with respect to
the center of valley $\tau$. This Hamiltonian contains matrices with
intralayer hopping amplitudes ($t$), nearest-layer hopping amplitudes
($t_{12}$), and next-nearest-layer hopping amplitudes ($t_{13}$), as well as
possible potential differences between different layers and/or sublattices due
to external gates and/or broken symmetries (terms $U_{i}$). Explicitly, the
various hopping amplitudes are given by

$\displaystyle t(\bm{k})=\begin{bmatrix}0&v_{0}\pi^{\dagger}\\\
v_{0}\pi&0\end{bmatrix},\quad
t_{12}(\bm{k})=\begin{bmatrix}-v_{4}\pi^{\dagger}&v_{3}\pi\\\
\gamma_{1}&-v_{4}\pi^{\dagger}\end{bmatrix},\quad
t_{13}(\bm{k})=\begin{bmatrix}0&\gamma_{2}/2\\\ 0&0\end{bmatrix},$ (S2)

where $\pi=\tau k_{x}+ik_{y}$ is a linear momentum, and $\gamma_{i}$
($i=0,\dots,4$) are hopping amplitudes of the tight-binding model with
corresponding velocity parameters $v_{i}=(\sqrt{3}/2)a\gamma_{i}/\hbar$. We
use the model parameters listed in Table 1, which are chosen to match quantum
oscillation frequency signatures in Ref. [1].

Table 1: Tight-binding parameters (in eV) for rhombohedral trilayer graphene, see also Refs. [11, 30, 1]. $\gamma_{0}$ | $\gamma_{1}$ | $\gamma_{2}$ | $\gamma_{3}$ | $\gamma_{4}$ | $U$ | $\Delta$ | $\delta$
---|---|---|---|---|---|---|---
$3.160$ | $0.380$ | $-0.015$ | $-0.290$ | $0.141$ | $0.030$ | $-0.0023$ | $-0.0105$

As shown in Fig. S1(a), the electronic band structure of $h^{\alpha}(\bm{k})$
in Eq. S1 is strongly trigonally-warped (by interlayer hopping amplitudes) and
has van Hove singularities near charge neutrality. The dispersion is
semimetalic in absence of a interlayer potential $U$ but has a insulating gap
when it is finite. In our work, we focus on the hole-doped side of the
dispersion, in the electronic density range of $n_{e}\in[-1.8\cdot
10^{12},0]\;\text{cm}^{-2}$ accessible to recent experiments [1], and consider
interlayer potentials in the range $U\in[0,80]$ meV, corresponding to
displacement fields $D\in[0,0.5]\;\text{V}/\text{nm}$ after including Hartree
screening between the system and gates [12]. The latter is done implicitly by
matching the quantum oscillation frequency data [1] at finite interlayer
potential to band structure calculations. As shown in Fig. Fig. S1(b), the
topology and shape of the Fermi sea (FS) undergoes Lifshitz transitions when
the electron density $n_{e}$ and the potentials $U$ changes over the
experimentally accessible range: from a single simple FS to an annular FS, and
then to three or six FS pockets. Each transition is accompanied by a van-Hove
singularity that induces a discontinuity in the density of states: a jump in
the case of going from single FS to annular FS, and a saddle-point log-
divergence when tuning towards isolated pockets.

Figure S1: Non-interacting band structure of rhombohedral trilayer graphene.
(a) Electronic bands for valley $K$ (purple) and valley $K^{\prime}$ (green)
without and with interlayer potential $U=20$ meV, and corresponding density of
states. (b) Density of states (DOS) at the Fermi surface as function of
electronic density $n_{e}$ and interlayer potential $U$. Insets indicate the
Fermi sea (FS) shape and topology per valley–spin found in different regions
of the parameter space. The Lifshitz transition from single FS to annular FS
is accompanied by a jump in the DOS, while transitions from annular to
disjoint pockets is accompanied by a log-divergent saddle-point van-Hove
singularity.

##  II Coulomb interactions

### Hamiltonian

In this work, we investigate the ground state of $H=H_{0}+H_{\rm ee}$, where
$H_{0}$ is the electronic Hamiltonian of rhombohedral trilayer graphene as in
Eq. S1 and $H_{\rm ee}$ is the Coulomb interaction between electrons. In
second quantization, they are

$\displaystyle H=H_{0}+H_{\rm ee}$
$\displaystyle=\sum_{\bm{k},\alpha,ij}\psi_{\alpha
i}^{\dagger}(\bm{k})\,h_{ij}(\tau
k_{x},k_{y})\,\psi^{\phantom{\dagger}}_{\alpha
j}(\bm{k})+\frac{1}{2A}\sum_{\bm{k},\bm{k}^{\prime}\\!\\!,\bm{q}}\\!V_{\bm{q}}\,\psi_{\alpha
i}^{\dagger}(\bm{k}\\!+\\!\bm{q})\psi_{\beta
j}^{\dagger}(\bm{k}^{\prime}\\!\\!-\\!\bm{q})\psi^{\phantom{\dagger}}_{\beta
j}(\bm{k}^{\prime})\psi^{\phantom{\dagger}}_{\alpha i}(\bm{k}),$ (S3)

where $\bm{k}=(k_{x},k_{y})$ is the same Bloch momentum as before, $i,j$ run
through the six atoms in the unit cell, $\alpha,\beta$ run through the four
flavor combinations of spin and valley, $A$ is the system area, and
$V_{\bm{q}}$ is the electron–electron interaction potential in Fourier space.
Note that the interaction here only depends on the in-plane distance between
electrons, but is independent of layer, spin and valley degrees of freedom.
The layer dependence is small and was neglected, which does not qualitatively
affect our result. The spin–valley symmetry in our interaction is
$\mathrm{SU}(4)$ since we neglect intervalley Coulomb scattering, which would
reduce it to $\mathrm{SU}_{\text{spin}}(2)\times\mathrm{SU}_{\text{spin}}(2)$
(the same as that of the non-interacting Hamiltonian). As explained below, we
use the dual-gate-screened interaction potential

$\displaystyle V_{\bm{q}}=\frac{2\pi
k_{e}}{\epsilon_{r}}\frac{\tanh(|\bm{q}|D)}{|\bm{q}|}$ (S4)

with metallic gates at distance $D\sim 30,\dots,50$ nm, relative permeability
$\epsilon_{r}\sim 4$–$10$, and Coulomb constant $k_{e}=1.44$ eV nm. The
parameters $D$ and $\epsilon_{r}$ are treated as tuning parameters in our
work, and are chosen informed by the onset of symmetry breaking observed in
the experiment. Screening due to electron–electron interactions is not taken
into account. We argue that the internal screening plays a smaller role than
the gate screening, because the latter already removes the $\bm{q}\to 0$
divergence and exponentially cuts of the range in $k$-space of the interaction
potential.

### Gate screening

Gate screening is taken into account by considering the mirror charges induced
by the metallic surface of nearby gates, which in the dual-gated case leads to
an effective potential

$\displaystyle
V(\bm{r},\bm{r}^{\prime})=\frac{k_{e}}{\epsilon_{r}}\sum_{n=-\infty}^{\infty}(-1)^{n}\frac{1}{\sqrt{((\bm{r}_{\perp}-\bm{r}^{\prime}_{\perp})^{2}+\xi_{n}^{2}}},$
(S5)

where $\bm{r}_{\perp}^{(\prime)}$ denotes the in-plane position and
$\xi_{n}=r_{z}-(-1)^{n}(r_{z}^{\prime}-nD)$ is the vertical distance from a
test charge to other (mirror) charges. The Fourier transform of this potential
can be evaluated to

$\displaystyle V(q)$ $\displaystyle=\frac{2\pi
k_{e}}{\epsilon_{r}q}\sum_{n=-\infty}^{\infty}(-1)^{n}e^{-\xi_{n}q}=\frac{2\pi
k_{e}}{\epsilon_{r}q}\begin{cases}\frac{\cosh(Dq)-\cosh(dq)}{\sinh(Dq)}&\text{intralayer}\\\
\frac{\cosh((D-d)q)-1}{\sinh(Dq)}&\text{interlayer}\end{cases}$ (S6)

where $q=|\bm{q}_{\perp}|$ is the in-plane momentum, and simplifies to Eq. S4
if $D\gg d$.

##  III Self-consistent Hartree-Fock approximation

We seek to minimize the the total energy $E_{\text{tot}}$ of $H=H_{0}+H_{\rm
ee}$ within the subspace of Slater-determinant states. In this variational
ansatz, we have
$E_{\text{tot}}\approx\tr(\rho[h+h^{\text{SCMF}}(\rho)])/(2\tr(\rho))$ with a
density matrix $\rho$ that satisfies the self-consistent equations

$\displaystyle
h^{\text{SCMF}}(\rho,\bm{k})\ket{n\bm{k}}=\varepsilon_{n\bm{k}}\ket{n\bm{k}},\quad\rho\\!=\\!\frac{1}{N_{k}}\sum_{n\bm{k}}f_{n\bm{k}}\outerproduct{n\bm{k}}{n\bm{k}},$
(S7)

where $f_{n\bm{k}}$ are Fermi-Dirac weights and
$h^{\text{SCMF}}(\rho)=h+\Sigma^{H}(\rho)+\Sigma^{F}(\rho)$ is an effective
single-particle Hamiltonian. The latter contains a Hartree potential and a
Fock exchange term, i.e.,

$\displaystyle\Sigma^{H}_{ab}(\rho)=\delta_{ab}\frac{1}{A}\sum_{b}V^{ab}(\bm{q}=0)\,n_{b},$
$\displaystyle\Sigma^{F}_{ab}(\rho,\bm{k})=-\frac{1}{A}\sum_{\bm{k}^{\prime}}V_{ab}(\bm{k}-\bm{k}^{\prime})\,\rho_{ab}^{*}(\bm{k}^{\prime}),$
(S8)

that contribute distinctly to $E_{\text{tot}}=E_{\text{band}}+E_{H}+E_{F}$.
Hartree describes how the total electronic density in degree of freedom $b$
affects (repels if $V$ is positive) electrons in degree of freedom $a$, and is
just a classical electrostatic effect. The Fock term is a quantum mechanical
contribution due to the antisymmetry of Fermions under exchange. The chemical
potential $\mu$ is constrained through the density
$n_{e}A=N_{k}^{-1}\sum_{n\bm{k}}f_{n\bm{k}}$, where $N_{k}$ is the number of
discrete k-points.

Solving Eq. S7 requires iterative methods that relax an initial guess
$\rho_{0}$ into a fixpoint. We identify lowest-energy solutions by exploring
different initial states $\rho_{0}$ (e.g., random, valley-polarized, spin-
polarized ones). Substituting $h^{\text{SCMF}}$ with
$h^{\text{SCMF}}(\rho)-\lambda\,X$ in Eq. S7 imposes an additional soft
constraint: $\lambda$ is a (small) Lagrange multiplier that favors the order
parameter $\langle X\rangle$. In the main text, we used
$X=\rho(\bm{k}_{\text{VBM}})$ to demonstrate how the crescent Fermi surface
can be seeded at any point along the valence band maximum ring.

##  IV Intervalley coherent state without nematicity

In this section, we show that an intervalley coherent state without nematicity
cannot explain the sudden jump of Shubnikov-de Haas oscillation frequencies at
the first magnetic transitions. In order to understand this, we consider an
intervalley coherent state in Fig. S2 which has two hole-like Fermi surfaces
and one electron-like Fermi surface (the smallest one). The areas enclosed by
these Fermi surfaces lead to the following three normalized frequencies
$\bm{f}=(0.39,0.25,0.14)$. Note that the difference between the total area
enclosing holes ($0.39+0.25=0.64$) and the area enclosing electrons ($0.14$)
is half. However, such a state is not consistent with the experimental finding
which has only two prominent frequencies $\bm{f}_{\rm exp}\sim(0.42,0.08)$
suggesting that there are only two hole-like Fermi surfaces [3]. The electron-
like Fermi surface can be removed when electrons cluster at a point on the
ring of van Hove singularities (valence band maxima) in Fig. S2(a) and this
leads to the valley-XY nematic metal described in the main text.

Figure S2: A local minimum mean-field solution of the Hartree-Fock equation
which shows finite intervalley coherence $\langle\tau_{x}\rangle\neq 0$ but no
nematicity. (a) Intervalley-mixed quasiparticle bands and (b) corresponding
Fermi sea where the color scale indicates the magnitude of the order parameter
$\tau_{x}$. The characters of the three Fermi surfaces with increasing radii
are respectively electron-like, hole-like and hole-like. This state is
calculated with dielectric permeability $\epsilon_{r}=4.3$, gate distance
$D=50$nm at electron density $n_{e}=-1.75\cdot 10^{12}\text{cm}^{-2}$ and
interlayer potential $U=70$meV. While this is a local minimum of the Hartree-
Fock equation, the solution with finite nematicity we described in the main
text has a lower energy.


* [50] Grassi A., Mariño M., The complex side of the TS/ST correspondence, J. Phys. A: Math. Theor. 52 (2019), 055402, 22 pages, arXiv:1708.08642.
  * [51] Gu J., Sulejmanpasic T., High order perturbation theory for difference equations and Borel summability of quantum mirror curves, J. High Energy Phys. 2017 (2017), no. 12, 014, 36 pages, arXiv:1709.00854.
  * [52] Gutzwiller M.C., The quantum mechanical Toda lattice, Ann. Physics 124 (1980), 347–381.
  * [53] Gutzwiller M.C., The quantum mechanical Toda lattice. II, Ann. Physics 133 (1981), 304–331.
  * [54] Hatsuda Y., Mariño M., Exact quantization conditions for the relativistic Toda lattice, J. High Energy Phys. 2016 (2016), no. 5, 133, 35 pages, arXiv:1511.02860.
  * [55] Hatsuda Y., Mariño M., Moriyama S., Okuyama K., Non-perturbative effects and the refined topological string, J. High Energy Phys. 2014 (2014), no. 9, 168, 42 pages, arXiv:1306.1734.
  * [56] Huang M.-X., Klemm A., Poretschkin M., Refined stable pair invariants for E-, M- and $[p,q]$-strings, J. High Energy Phys. 2013 (2013), no. 11, 112, 117 pages, arXiv:1308.0619.
  * [57] Huang M.-X., Sun K., Wang X., Blowup equations for refined topological strings, J. High Energy Phys. 2018 (2018), no. 10, 196, 87 pages, arXiv:1711.09884.
  * [58] Humphreys J.E., Introduction to Lie algebras and representation theory, Graduate Texts in Mathematics, Vol. 9, Springer-Verlag, New York – Berlin, 1972.
  * [59] Iqbal A., Kozçaz C., Vafa C., The refined topological vertex, J. High Energy Phys. 2009 (2009), no. 10, 069, 58 pages, arXiv:hep-th/0701156.
  * [60] Ito K., Shu H., ODE/IM correspondence and the Argyres–Douglas theory, J. High Energy Phys. 2017 (2017), no. 8, 071, 22 pages.
  * [61] Jentschura U.D., Surzhykov A., Zinn-Justin J., Multi-instantons and exact results. III. Unification of even and odd anharmonic oscillators, Ann. Physics 325 (2010), 1135–1172, arXiv:1001.3910.
  * [62] Kozlowski K., Teschner J., TBA for the Toda chain, in New Trends in Quantum Integrable Systems, World Sci. Publ., Hackensack, NJ, 2011, 195–219, arXiv:1006.2906.
  * [63] Källén J., Mariño M., Instanton effects and quantum spectral curves, Ann. Henri Poincaré 17 (2016), 1037–1074, arXiv:1308.6485.
  * [64] Kashani-Poor A.-K., Quantization condition from exact WKB for difference equations, J. High Energy Phys. 2016 (2016), no. 6, 180, 34 pages, arXiv:1604.01690.
  * [65] Katz S., Klemm A., Vafa C., Geometric engineering of quantum field theories, Nuclear Phys. B 497 (1997), 173–195, arXiv:hep-th/9609239.
  * [66] Kharchev S., Lebedev D., Integral representation for the eigenfunctions of a quantum periodic Toda chain, Lett. Math. Phys. 50 (1999), 53–77, arXiv:hep-th/9910265.
  * [67] Klemm A., Lerche W., Mayr P., Vafa C., Warner N., Self-dual strings and $N=2$ supersymmetric field theory, Nuclear Phys. B 477 (1996), 746–764, arXiv:hep-th/9604034.
  * [68] Klemm A., Lerche W., Theisen S., Nonperturbative effective actions of $(N=2)$-supersymmetric gauge theories, Internat. J. Modern Phys. A 11 (1996), 1929–1973, arXiv:hep-th/9505150.
  * [69] Klemm A., Lerche W., Yankielowicz S., Theisen S., Simple singularities and $N=2$ supersymmetric Yang–Mills theory, Phys. Lett. B 344 (1995), 169–175, arXiv:hep-th/9411048.
  * [70] Konishi K., Paffuti G., Quantum mechanics: a new introduction, Oxford University Press, 2009.
  * [71] Landau L.D., Lifshitz E.M., Quantum mechanics: non-relativistic theory, Elsevier, 2013.
  * [72] Laptev A., Schimmer L., Takhtajan L.A., Weyl type asymptotics and bounds for the eigenvalues of functional-difference operators for mirror curves, Geom. Funct. Anal. 26 (2016), 288–305, arXiv:1510.00045.
  * [73] Losev A.S., Marshakov A.V., Nekrasov N.A., Small instantons, little strings and free fermions, in From Fields to Strings: Circumnavigating Theoretical Physics, Vol. 1, World Sci. Publ., Singapore, 2005, 581–621, arXiv:hep-th/0302191.
  * [74] Macdonald I.G., Symmetric functions and Hall polynomials, Oxford Mathematical Monographs, The Clarendon Press, Oxford University Press, New York, 1979.
  * [75] Mariño M., Lectures on non-perturbative effects in large $N$ gauge theories, matrix models and strings, Fortschr. Phys. 62 (2014), 455–540, arXiv:1206.6272.
  * [76] Mariño M., Instantons and large $N$. An introduction to non-perturbative methods in quantum field theory, Cambridge University Press, Cambridge, 2015.
  * [77] Mariño M., Spectral theory and mirror symmetry, in String-Math 2016, Proc. Sympos. Pure Math., Vol. 98, Amer. Math. Soc., Providence, RI, 2018, 259–294, arXiv:1506.07757.
  * [78] Mariño M., Zakany S., Exact eigenfunctions and the open topological string, J. Phys. A: Math. Theor. 50 (2017), 325401, 50 pages, arXiv:1606.05297.
  * [79] Mariño M., Zakany S., Wavefunctions, integrability, and open strings, arXiv:1706.07402.
  * [80] Martinec E.J., Warner N.P., Integrable systems and supersymmetric gauge theory, Nuclear Phys. B 459 (1996), 97–112, arXiv:hep-th/9509161.
  * [81] Matone M., Instantons and recursion relations in $N=2$ SUSY gauge theory, Phys. Lett. B 357 (1995), 342–348, arXiv:hep-th/9506102.
  * [82] Matsuyama A., Periodic Toda lattice in quantum mechanics, Ann. Physics 220 (1992), 300–334.
  * [83] Meneghelli C., Yang G., Mayer-cluster expansion of instanton partition functions and thermodynamic Bethe ansatz, J. High Energy Phys. 2014 (2014), no. 5, 112, 42 pages, arXiv:1312.4537.
  * [84] Mironov A., Morosov A., Nekrasov functions and exact Bohr–Sommerfeld integrals, J. High Energy Phys. 2010 (2010), no. 4, 040, 15 pages, arXiv:0910.5670.
  * [85] Mironov A., Morozov A., Nekrasov functions from exact Bohr–Sommerfeld periods: the case of ${\rm SU}(N)$, J. Phys. A: Math. Theor. 43 (2010), 195401, 11 pages, arXiv:0911.2396.
  * [86] Nekrasov N., Five-dimensional gauge theories and relativistic integrable systems, Nuclear Phys. B 531 (1998), 323–344, arXiv:hep-th/9609219.
  * [87] Nekrasov N.A., Seiberg–Witten prepotential from instanton counting, Adv. Theor. Math. Phys. 7 (2003), 831–864, arXiv:hep-th/0206161.
  * [88] Nekrasov N.A., Okounkov A., Seiberg–Witten theory and random partitions, in The Unity of Mathematics, Progr. Math., Vol. 244, Birkhäuser Boston, Boston, MA, 2006, 525–596, arXiv:hep-th/0306238.
  * [89] Nekrasov N.A., Shatashvili S.L., Quantization of integrable systems and four dimensional gauge theories, in XVIth International Congress on Mathematical Physics, World Sci. Publ., Hackensack, NJ, 2010, 265–289, arXiv:0908.4052.
  * [90] Pasquier V., Gaudin M., The periodic Toda chain and a matrix generalization of the Bessel function recursion relations, J. Phys. A: Math. Gen. 25 (1992), 5243–5252.
  * [91] Poghossian R., Deforming SW curve, J. High Energy Phys. 2011 (2011), no. 4, 033, 12 pages, arXiv:1006.4822.
  * [92] Sciarappa A., Bethe/gauge correspondence in odd dimension: modular double, non-perturbative corrections and open topological strings, J. High Energy Phys. 2016 (2016), no. 10, 014, 57 pages, arXiv:1606.01000.
  * [93] Sciarappa A., Exact relativistic Toda chain eigenfunctions from separation of variables and gauge theory, J. High Energy Phys. 2017 (2017), no. 10, 116, 87 pages, arXiv:1706.05142.
  * [94] Seara T.M., Sauzin D., Borel summation and the theory of resurgence, Butl. Soc. Catalana Mat. 18 (2003), 131–153.
  * [95] Seiberg N., Witten E., Electric-magnetic duality, monopole condensation, and confinement in ${\mathcal{N}}=2$ supersymmetric Yang–Mills theory, Nuclear Phys. B 426 (1994), 19–52, Erratum, Nuclear Phys. B 430 (1994), 485–486, arXiv:hep-th/9407087.
  * [96] Silverstone H.J., JWKB connection-formula problem revisited via Borel summation, Phys. Rev. Lett. 55 (1985), 2523–2526.
  * [97] Sklyanin E.K., The quantum Toda chain, in Nonlinear Equations in Classical and Quantum Field Theory (Meudon/Paris, 1983/1984), Lecture Notes in Phys., Vol. 226, Springer, Berlin, 1985, 196–233.
  * [98] Sulejmanpasic T., Ünsal M., Aspects of perturbation theory in quantum mechanics: the BenderWu Mathematica® package, Comput. Phys. Commun. 228 (2018), 273–289, arXiv:1608.08256.
  * [99] Sun K., Wang X., Huang M.-X., Exact quantization conditions, toric Calabi–Yau and non-perturbative topological string, J. High Energy Phys. 2017 (2017), no. 1, 061, 102 pages, arXiv:1606.07330.
  * [100] Tachikawa Y., A review on instanton counting and W-algebras, in New Dualities of Supersymmetric Gauge Theories, Springer, Cham, 2016, 79–120, arXiv:1412.7121.
  * [101] Taki M., Refined topological vertex and instanton counting, J. High Energy Phys. 2008 (2008), no. 3, 048, 22 pages, arXiv:0710.1776.
  * [102] Voros A., Spectre de l’équation de Schrödinger et méthode BKW, Publications Mathématiques d’Orsay 81, Vol. 9, Université de Paris-Sud, Orsay, 1982.
  * [103] Voros A., The return of the quartic oscillator: the complex WKB method, Ann. Inst. H. Poincaré Sect. A 39 (1983), 211–338.
  * [104] Voros A., Exact anharmonic quantization condition (in one dimension), in Quasiclassical Methods (Minneapolis, MN, 1995), IMA Vol. Math. Appl., Vol. 95, Springer, New York, 1997, 189–224.
  * [105] Voros A., Exact resolution method for general $1$D polynomial Schrödinger equation, J. Phys. A: Math. Gen. 32 (1999), 5993–6007, arXiv:math-ph/9903045.
  * [106] Wang X., Zhang G., Huang M.-X., New exact quantization condition for toric Calabi–Yau geometries, Phys. Rev. Lett. 115 (2015), 121601, 5 pages, arXiv:1505.05360.
  * [107] Wilczek F., Quantum time crystals, Phys. Rev. Lett. 109 (2012), 160401, 5 pages, arXiv:1202.2539.
  * [108] Yaris R., Bendler J., Lovett R.A., Bender C.M., Fedders P.A., Resonance calculations for arbitrary potentials, Phys. Rev. A 18 (1978), 1816–1825.
  * [109] Zinn-Justin J., Multi-instanton contributions in quantum mechanics. II, Nuclear Phys. B 218 (1983), 333–348.
  * [110] Zinn-Justin J., Jentschura U.D., Multi-instantons and exact results. I. Conjectures, WKB expansions, and instanton interactions, Ann. Physics 313 (2004), 197–267, arXiv:quant-ph/0501136.
  * [111] Zinn-Justin J., Jentschura U.D., Multi-instantons and exact results. II. Specific cases, higher-order effects, and numerical calculations, Ann. Physics 313 (2004), 269–325, arXiv:quant-ph/0501137.