# Intervalley coherence and intrinsic spin-orbit coupling in rhombohedral
trilayer graphene

Trevor Arp  These authors contributed equally: Trevor Arp and Owen Sheekey
Owen Sheekey  These authors contributed equally: Trevor Arp and Owen Sheekey
Haoxin Zhou  Current address: Department of Electrical Engineering and
Computer Sciences, University of California, Berkeley, California 94720, USA.
C.L. Tschirhart  Current address: Laboratory of Atomic and Solid-State
Physics, Cornell University, Ithaca, NY, USA Caitlin L. Patterson  H. M. Yoo
Ludwig Holleis  Evgeny Redekop  Grigory Babikyan  Tian Xie  Department of
Physics, University of California at Santa Barbara, Santa Barbara CA 93106,
USA  Jiewen Xiao  Yaar Vituri  Tobias Holder  Department of Condensed Matter
Physics, Weizmann Institute of Science, Rehovot, Israel  Takashi Taniguchi
International Center for Materials Nanoarchitectonics, National Institute for
Materials Science, 1-1 Namiki, Tsukuba 305-0044, Japan  Kenji Watanabe
Research Center for Functional Materials, National Institute for Materials
Science, 1-1 Namiki, Tsukuba 305-0044, Japan  Martin E. Huber  Departments of
Physics and Electrical Engineering, University of Colorado Denver, Denver, CO
80217, USA.  Erez Berg  Department of Condensed Matter Physics, Weizmann
Institute of Science, Rehovot, Israel  Andrea F. Young
<EMAIL_ADDRESS>Department of Physics, University of California at
Santa Barbara, Santa Barbara CA 93106, USA

Rhombohedral graphene multilayers provide a clean and highly reproducible
platform to explore the emergence of superconductivity[1, 2, 3, 4] and
magnetism[5, 6, 2, 7, 8, 9, 10, 11] in a strongly interacting electron system.
Here, we use electronic compressibility and local magnetometry to explore the
phase diagram of this material class in unprecedented detail. We focus on
rhombohedral trilayer in the quarter metal regime, where the electronic ground
state is characterized by the occupation of a single spin and valley isospin
flavor. Our measurements reveal a subtle competition between valley imbalanced
(VI) orbital ferromagnets and intervalley coherent (IVC) states in which
electron wave functions in the two momentum space valleys develop a
macroscopically coherent relative phase. Contrasting the in-plane spin
susceptibility of the IVC and VI phases reveals the influence of graphene’s
intrinsic spin-orbit coupling[12, 13, 14, 15], which drives the emergence of a
distinct correlated phase with hybrid VI and IVC character. Spin-orbit also
suppresses the in-plane magnetic susceptibility of the VI phase, which allows
us to extract the spin-orbit coupling strength of $\lambda\approx 50\mu eV$
for our hexagonal boron nitride-encapsulated graphene system. We discuss the
implications of finite spin-orbit coupling on the spin-triplet superconductors
observed in both rhombohedral and twisted graphene multilayers.

A defining feature of correlated low-energy electron physics in both moiré and
crystalline graphene is the large number of closely competing ground states.
This near-degeneracy is traceable to the approximate SU(4) symmetry within the
combined spin and valley ‘isospin’ space, which allows for a large number of
broken symmetry phases with nearly degenerate energies. In experiment, this
degeneracy may be lifted either spontaneously or by weak symmetry breaking
terms in the Hamiltonian, which compete to determine the ground state. Some
symmetry breaking effects arise at the single particle level, such as atomic
scale spin-orbit coupling[12, 13, 14]–though this effect has been presumed to
be negligible in the theoretical literature. Other symmetries are broken by
the inter-particle interactions themselves. For example, differences between
inter- and intravalley scattering reduces the symmetry of the Hamiltonian to
$SU(2)_{\text{spin}}\times U(1)_{K}\times U(1)_{K^{\prime}}$ representing spin
rotation and the independent conservation of charge in the K and K’ valleys,
respectively[16]. Finally, interactions may drive the spontaneous breaking of
the symmetries that do remain: for example through the formation of
intervalley coherent (IVC) ground states characterized by a macroscopically
coherent phase between the wave function in the two inequivalent valleys. At
the theoretical level, the strength of the numerous microscopic parameters
that govern breaking of the symmetries cannot be accurately determined from
first principles. As a result, experimental determination of the ground state
provides the primary method to constrain the microscopic Hamiltonian. This has
proven particularly challenging in twisted graphene multilayers, where
experimental irreproducibility has hampered efforts to reliably determine the
phase diagram.

Rhombohedral graphene multilayers provide a structurally simple and
experimentally reproducible platform, enabling a concrete connection between
precision measurements and many-body theory. In rhombohedral graphene alone,
experiments have revealed signatures of symmetry breaking states that include
nematics[17, 9, 4, 3], spin and orbital magnets[5, 6, 2, 7, 8, 18, 11, 19],
and superconductors[1, 2, 3, 4]. However, much of the observed
phenomenology—particularly the role of spin in the superconducting
phases—remains unexplained. A particularly striking puzzle is the stability of
superconductivity in hexagonal boron nitride encapsulated Bernal bilayer
graphene, which shows no superconductivity at zero magnetic field but in which
a spin-polarized superconducting state emerges above a threshold magnetic
field applied in the plane of the sample[2]. Spin-valley locked
superconductivity is also induced (at zero magnetic field) by supporting the
bilayer on a WSe2[3, 4] substrate, which is known to induce a large Ising-type
spin-orbit coupling[20]—further pointing to the role of spin in the
superconducting phase diagram. Of course, in the absence of a clear
understanding of what stabilizes the superconducting states, consensus has
also proved elusive regarding the underlying mechanism, with both all-
electronic[21, 22, 23, 24, 25, 26, 27, 28, 29] and phonon-mediated[30]
mechanisms remaining viable in light of current experimental data.

Figure 1:  Thermodynamics of rhombohedral trilayer graphene in the hole-doped
quarter metal regime. (A) Schematic of the measurement geometry. The bottom
gate is modulated at kHz frequencies with amplitude $\delta\mathrm{V}$. A
transistor amplifier connected to the top gate measures the modulated charge
$\delta Q$, from which the inverse compressibility, $\kappa$, is derived[31].
A scanning SQUID magnetometer measures the local magnetic field $\delta B$.
(B) $\kappa$, as a function of carrier density $n_{e}$ and applied
displacement field $D$ at $T\approx 20$ mK and $B_{\parallel}=0.5$ T. (C)
Schematic phase diagram showing competing phases, classified by isospin
polarization and Fermi surface topology. Uppercase labels denote the isospin
polarization, which can be valley imbalanced (VI), intervalley coherent (IVC),
spin imbalanced (SI), or isospin symmetric (Sym); lower-case labels denote
simple (s), annular (a), or disjoint (d) Fermi surface topology. The region
marks SC2 denotes the superconductor reported previously[1] near the edge of
the SI phase. As discussed in the text, the isospin order of the gray region
cannot be determined conclusively. (D) Spatial image of $\delta B$ taken $\sim
300$ nm above the rhombohedral trilayer. The image is acquired at
$B_{\perp}=15$ mT, $n_{e}=-0.67\times 10^{12}$ cm-2 and $D=0.55$ V/nm,
corresponding to a first order transition into an orbitally magnetized phase.
The black outline indicates the design geometry of the dual gated region of
the device. Scale bar is 1 $\mu$m. (E) $\delta B$ at $B_{\parallel}=0.5$T and
$T=300$ mK, acquired at the position marked in panel C. Lines indicate the
boundaries between phases from A, with solid lines indicating negative
compressibility peaks and dashed lines compressibility steps associated with
first order and Lifshitz transitions, respectively. The arrow at bottom right
shows the direction of the $\delta\mathrm{V}$ modulation. (F) $\delta B$
(blue) and $\kappa$ (red) along the dashed line shown in panel A. The regions
with VI and IVC order are separated by a first order phase transition (visible
in $\kappa$), as well as a change in magnetization indicated by a peak in
$\delta B$. Changes in Fermi surface topology occur at Lifshitz transitions
—steps in $\kappa$–and are not accompanied by features in $\delta B$.

Here we combine global charge sensing and local magnetometry measurements in
rhombohedral trilayer graphene to explore the nature of the isospin
ferromagnetic phases that emerge as the material is doped through the Van Hove
singularities bracketing the neutrality point. Fig. 1A shows a schematic of
our measurement geometry. Our rhombohedral trilayer graphene sample is
encapsulated between hexagonal boron nitride gate dielectrics and graphite
flakes which serve as electrostatic gates to control the total charge carrier
density $n_{e}$ and applied displacement field $D$ (details of sample
fabrication, and transport measurements from the same device, were reported
previously[6, 1]). We use a cryogenic transistor amplifier[32] connected to
the top gate to monitor the charge $\delta Q$ induced by the modulated bottom
gate voltage $\delta\mathrm{V}$. $\delta Q$ is proportional to the spatially
averaged inverse compressibility $\kappa=\partial\mu/\partial n_{e}$. In
addition, we use a scanning superconducting quantum interference device to
image the local magnetic field, $\delta B$, that arises in response to the
same modulated voltage. For purely out of plane magnetic moments, $\delta B$
is proportional to the gate-modulated change in magnetization.

##  I Intervalley coherent quarter metal

Figure 1B shows $\kappa$ measured for hole doping, with $B_{\parallel}=0.5$ T
magnetic field applied to ensure in-plane spin polarization. As shown
previously, in the low-$|n_{e}|$, low $D$ extreme of this range the electron
system is in a symmetric (‘Sym’) phase that preserves isospin symmetry. At low
$D$ and high-$|n_{e}|$, the system is a valley unpolarized, spin-polarized
half-metal[1] in which the two Fermi seas have annular topology. We denote
this phase SI in the schematic in Fig. 1C, as it is an example of a spin-
imbalanced phase. Spin-triplet superconductivity was reported at the
low-$|n_{e}|$ extreme of this ‘half-metal’ SI phase[1]. In the rest of the
phase space spanned by Figs. 1B,C measurements of quantum oscillations have
shown[1] evidence of valley polarization. In the regions rendered in red and
blue in Fig. 1B, quantum oscillations show a single Fermi sea with either
simple or annular topology—a quarter metal. Quantum oscillations in the
intermediate phase between the SI and quarter metal phases were interpreted as
indicating partial valley polarization[6, 1].

In the present work, we improve upon previous compressibility measurements
through the use of a thermal decoupling capacitor between sample and cryogenic
amplifier and an additional cryogenic amplification stage at 4K (see Methods).
The data in Fig. 1B reveals both ‘step’ features, associated with Lifshitz
transitions, and negative compressibility spikes, associated with first order
phase transitions, that were not previously reported. Most strikingly, a dip
in $\kappa$ bisects both the simple and annular quarter metal regime. As
indicated in Fig. 1C, we ascribe this feature to a transition between distinct
valley imbalanced (VI) and inter-valley coherent (IVC) quarter metal phases.

The contrast between these phases can be visualized directly using local
magnetometry. Fig. 1D shows a spatial map of $\delta B$ acquired $300$ nm
above the sample layer at a point in the $n_{e},D,B_{\perp}$ and
$B_{\parallel}$-tuned parameter space where the system is at a transition
between phases with contrasting valley population imbalance. The $\delta B$
image of Fig. 1D thus represents the change in fringe magnetic field
associated with a change in the sample’s orbital magnetic moment. We find the
spatial structure of $\delta B$ to be largely independent of $n_{e}$ and $D$,
indicating that the phase transitions observed in our global $\kappa$
measurements occur in the entire sample simultaneously. An exception to this
rule is the behavior of transitions associated with spin imbalance changes as
a function of $B_{\parallel}$ (Extended Data Fig. 5); specifically, we observe
the signatures of in-plane spin moments at these transitions when
$|B_{\parallel}|>|B_{\perp}|$.

Fig. 1E shows $\delta B$ measured at a single spatial position as a function
of $n_{e}$ and $D$ for $B_{\parallel}=0.5$ T. The position is chosen to
minimize sensitivity to the in-plane spin magnetic moment. The spatial
uniformity of the sample then allows us to interpret $\delta B$ as equivalent
to the modulated out-of-plane magnetic moment, $\delta m_{z}$, up to a
multiplicative constant. We assume the spins are fully polarized in the plane,
so that $\delta m_{z}$ arises only from the intrinsically out-of-plane orbital
magnetization.[33]

Comparing $\delta B$ (Figs. 1E) and $\kappa$ (Figs. 1B) reveals the nature of
the two quarter metal phases. Only transitions where the valley imbalance—and
thus orbital moment—changes show contrast in $\delta B$. Within the quarter
metal phase, the region to the left of the first order transition shows the
signatures of a finite orbital moment, and we assign it to a valley imbalanced
phase characterized by the complete- or near-complete polarization of the
electron system into a single valley. Additional transitions are visible in
$\kappa$ within the VI phase, which we assign to Lifshitz transitions in the
fermi sea topology that are not accompanied by any change in isospin symmetry.
Fermi sea topology may be simple ($s$), annular ($a$), or disjoint ($d$)–where
$d$ means there exist fermi pockets of different sizes but the same carrier
sign. These transitions manifest as steps in the electronic compressibility
without signature in $\delta B$, as shown in Fig. 1F.

Outside of this VI phase, none of the phases in the range of Figs. 1B,D have
finite orbital moment. This is expected for the Sym and SI phases, but is
novel in the case of the simple- and annular quarter metals, as well as the
region, indicated in gray in Fig. 1C separating the quarter metal and SI
phases. In the case of quarter metals only an intervalley coherent (IVC)
state, characterized by a coherent superposition of states in the two valleys,
is consistent with a vanishing orbital moment and full flavor polarization. We
do not make a definitive assignment for the phase marked in gray. This phase
is spin polarized and valley balanced (Extended Data Fig. 6), but it
disappears from the phase diagram by $B_{\perp}\approx 100$mT preventing
detailed study of quantum oscillations. While Hartree Fock calculations
(Extended Data Fig. 7) show that an intervalley coherent state may be
energetically favorable in this regime, the close energetic competition
between the IVC state and other spin-polarized states that break orbital
symmetries but preserve conservation of charge in the two valleys–for example,
nematics–cannot be ruled out.

Figure 2: Intervalley coherent quarter metal. (A) $\kappa$ as a function of
$n_{e}$ and $B_{\perp}$ along a trajectory that traverses the VI, IVC and Sym
phases at $B_{\parallel}=1.3$T and $T=20$ mK. (B) Higher resolution
measurement over the boxed region of panel A. (C) $R_{\text{xx}}(1/B_{\perp})$
(top) and its Fourier transform (bottom), plotted as a function of the
normalized frequency $f_{\nu}$. Curves measured in the VI and IVC phases both
show a single peak at $f_{\nu}=1$ indicating a single flavor polarized Fermi
surface. (D) Schematic depiction of the valley polarization of the IVC and VI
quarter metal phases (left) and their representation on a valley Block sphere
(right). For the IVC phase, (net) in-plane valley polarization, e.g.
$\langle\tau_{x}\rangle=1$, corresponds to a momentum dependent valley
polarization $\langle\tau_{z}\rangle$.

To substantiate the assignment of IVC order in Fig. 1E, we analyze $\kappa$,
measured across the phase boundary separating VI and IVC phases with simple
Fermi surface topology, as a function of magnetic field . The first order
nature of the transition implies a Clausius-Clapeyron-type relation between
the critical density $n_{e}^{*}$ and the magnetic moment difference
$\Delta\vec{m}$ between the two phases,

$\frac{\Delta\vec{m}}{\Delta\mu}=\frac{dn_{e}^{\ast}}{d\vec{B}},$ (1)

where $\Delta\mu$ is magnitude of the chemical potential jump between the two
phases and $\vec{m}\parallel\vec{B}$. Fig. 2A shows the IVC to VI transition
as a function of $B_{\perp}$. To isolate the effects of orbital magnetism, the
measurements are again performed in a large, fixed in-plane magnetic field
$B_{\parallel}=1.3T\gg B_{\perp}$. The boundary between the VI and IVC phase
evolves rapidly with $B_{\perp}$ in favor of the VI phase, consistent with
$m_{\perp}^{\text{VI}}>m_{\perp}^{\text{IVC}}$. The sharp cusp in the critical
density at $B_{\perp}=0$ implies a divergent orbital magnetic susceptibility
in the VI phase.

In contrast, the IVC-s/Sym phase boundary remains at fixed density as
$B_{\perp}$ is tuned, consistent with
$m_{\perp}^{\text{IVC}}=m_{\perp}^{\text{Sym}}=0$. Notably, quantum
oscillations show no break or phase shift across the IVC to VI transition
(Fig. 2B), implying consistent Fermi surface topology between the two phases.
This is confirmed in Fig. 2C, which shows $R_{\text{xx}}$ traces within the
two phases. Data is shown both as a function of $1/B_{\perp}$ and of
$f_{\nu}$, the frequency normalized to the total Luttinger volume. Both phases
show a single peak at $f_{\nu}=1$, corresponding to a single Fermi surface, a
behavior which extends across the entire simple quarter metal regime (see
Extended Data Fig. 8). Finally, Hall measurements show hysteretic anomalous
Hall effect only within the VI phase but not in the IVC phase (Extended Data
Fig. 9). Taken together, these data provide strong evidence for intervalley
coherent phases in which the ground state is an equal superposition of
wavefunctions in the two valleys (shown in Fig. 2D), leading to a cancellation
of the net orbital moment. We note that in real space, the IVC phases
described here are expected to show lattice-scale charge density wave order,
making them detectable by the same scanning tunneling microscopy techniques
that have identified intervalley coherent order in monolayer graphene at high
magnetic fields[34, 35] and moiré graphene systems[36, 37].

##  II Spin-orbit Coupling

Figure 3: Effects of spin-orbit coupling. (A) $\kappa$ across the VI to IVC
transition as a function of $B_{\parallel}$ at $T\approx 300$mK. The white
line shows the fit used to extract $\lambda\approx$ 46 $\mu$eV described in
the main text. (B) $\Delta\mu$ as determined by integrating the dip in
$\kappa$ (red) and $\Delta m_{\parallel}=\Delta\mu(\partial n_{e}/\partial
B_{\parallel})$ (blue). (C) In-plane spin magnetic moment of the VI phase in
the presence of spin-orbit coupling as a function of $B_{\parallel}$. Inset
shows the evolution of the VI phase in the Bloch sphere representation; the
spin cants as a function of $B_{\parallel}$ while the valley remains
polarized. (D) $\kappa$ across the VI to IVC transition at B=0, showing the
intermediate phase. (E) $\kappa$ in the hole-doped quarter metal regime at
$B_{\parallel}=B_{\perp}=0$. White overlay denotes range of SC2.[1] (F)
$\kappa$ across the VI to IVC transition as a function of $B_{\parallel}$
showing the suppression of the intermediate phase at low in-plane field. (G)
$\kappa$ as a function of $B_{\perp}$. The intermediate phase rapidly grows
with $B_{\perp}$, displacing the IVC phase completely by 50 mT. (H) Phase
diagram of a phenomenological model, described in the supplementary
information, that accounts for the competing effects of magnetic field,
exchange energy, Hund’s coupling, and spin-orbit effects. Insets show Fermi
surfaces, color coded by $\langle\tau_{z}\rangle$ as in Fig. 2D.  Figure 4:
Intervalley coherence for electron doping (A) $\kappa$ as a function of
carrier density $n_{e}$ and applied displacement field $D$ for electron doping
at $T\approx 20$ mK and $B_{\parallel}=0$. (B) Hartree-Fock phase diagram,
showing the ground state isospin polarization as a function of $n_{e}$ and the
interlayer potential $\Delta$, which is approximately proportional to $D$[38].
(C) Evolution of the Sym, VI, and IVC phases along the dashed line in panel A
as a function of $B_{\perp}$ and (D) $B_{\parallel}$.

The phase diagrams presented above are all acquired at large in-plane magnetic
field. However, many of the phase boundaries show non-trivial behavior as
$B_{\parallel}$ is reduced. Fig. 3A shows the VI to IVC boundary studied in
Figs. 1 and 2 as a function of $B_{\parallel}$. At high $B_{\parallel}$, the
transition is nearly $B_{\parallel}$ independent, but at low $B_{\parallel}$
the IVC phase is increasingly disfavored relative to the VI phase. As shown in
Fig. 3B, $\Delta\mu$ is independent of $B_{\parallel}$ over much of the range.
It follows from Eq. (1) that the evolution of $n_{e}^{\ast}$ implies a finite,
and $B_{\parallel}$ dependent, difference in $m_{\parallel}$ between the two
phases. Given the near complete layer polarization of the electron wave
function at the high $D$ of the experiment, $m_{\parallel}$ arises exclusively
from the electron spin. The finite $\Delta m_{\parallel}$ is thus unexpected,
as both the VI-s and IVC-s phases are in principle fully spin polarized, and
are expected, in the absence of spin-orbit coupling, to have identical and
divergent susceptibility to spin polarization in any direction.

Intrinsic spin-orbit coupling[12]—typically neglected in theoretical
treatments of graphene many body physics—provides a natural explanation for
the anomalous behavior of $\Delta m_{\parallel}$. Projected into the low-
energy bands near the K-points of the Brillouin zone, spin-orbit coupling in
graphene takes the form[12]
$H_{\text{KM}}=\frac{\lambda}{2}\sigma_{z}\tau_{z}s_{z}$ (where $\sigma_{z}$,
$\tau_{z}$, and $s_{z}$ are Pauli matrices in the sublattice, valley, and spin
space respectively). Theoretical estimates of $\lambda$ range over two orders
of magnitude, from 1 to 100 $\mu$eV [12, 13, 15, 14, 39]. Recent experiments
in mono- and bi-layer graphene, however, suggest that $\lambda$ falls
somewhere in the range 40 to 80 $\mu$eV [40, 41, 42].

Spin-orbit coupling has contrasting effects on the IVC and VI phases. At high
displacement field in the layer (and therefore sublattice) polarized limit,
Kane-Mele spin-orbit coupling reduces to a staggered effective Zeeman energy
which is opposite in sign in the two valleys (the ‘Ising’ spin-orbit coupling
familiar from other honeycomb systems with broken inversion symmetry). For
finite $\lambda$, a spin-polarized, valley-imbalanced phase may lower its
energy by aligning its spin in the spin-orbit favored, out-of-plane direction.
As $B_{\parallel}$ increases, the spins in the VI phase cant into the plane.
Taking $E_{z,\parallel}=g\mu_{B}B_{\parallel}$, the in-plane component of the
canted spin is given by
$\text{m}_{\parallel}^{\mathrm{VI}}=E_{z,\parallel}/[E_{z,\parallel}^{2}+\lambda^{2}]^{1/2}$,
plotted in Fig. 3C. In contrast, the energy of the valley balanced IVC phase
is independent of spin-orbit coupling to first order in $\lambda$. The IVC
phase may thus polarize in-plane for arbitrarily small $B_{\parallel}$, with a
divergent susceptibility and subsequent energy gain that is linear in the
applied field. A fit to the measured $n_{e}^{\ast}$, shown in overlay on Fig.
3A, gives $\lambda\approx 46\mu$eV, in agreement with the experimental
literature where this coupling was measured in different contexts and using
different experimental techniques[40, 41, 42].

Near $B_{\parallel}=0$, spin-orbit coupling produces even richer effects.
Measurements across the VI to IVC transition reveal an intermediate phase,
separated by compressibility dips from the IVC and VI phases (Figs. 3D-E).
Quantum oscillations in this phase indicate that it is also a quarter metal
with a single Fermi surface (Extended Data Fig. 8). As shown in Figs. 3F-G,
the intermediate phase has in- and out-of plane susceptibilities that are
intermediate between those of the IVC and VI phases. This is consistent with a
intervalley coherent phase with partial valley imbalance, endowing it with a
finite orbital moment but a smaller suppression of in-plane susceptibility
relative to the VI phase.

In fact, such a phase is also accounted for by the effect of spin-orbit
coupling on the spin-polarized IVC quarter metal phase. In the absence of
spin-orbit coupling, the spin magnetic moment can point in any direction.
Spin-orbit coupling induces either an easy-plane or an easy-axis spin
anisotropy within any spin polarized phase, including the IVC phases. If the
average spin moment is in-plane, the system will lower its energy by slightly
canting the spins out of the plane in an opposite direction in the two
valleys; we dub this phase IVCx. The energy gain associated with the canting
scales as $-\lambda^{2}/J_{H}$, where $J_{H}$ is the ferromagnetic inter-
valley Hund’s coupling. If the spin is pointing out of the plane, on the other
hand, the system lowers its energy by introducing a small imbalance between
the population of the two valleys. We dub this phase IVCz. The associated
energy gain is of the order of $-\chi_{v}\lambda^{2}$, where $\chi_{v}$ is the
valley susceptibility of the IVC phase in the absence of spin-orbit coupling.
As we approach the IVC to VI phase transition from the IVC side, and assuming
that the transition is not too strongly first order, $\chi_{v}$ is expected to
increase, and hence the out-of-plane polarized IVCz phase is favored over the
in-plane polarized IVCx phase. Fig. 3H shows the phase diagram of a simple
model that captures the relative energetics of the IVCx, IVCz, and VI phases
as a function of $B_{\parallel}$, $B_{\perp}$, and $\Delta_{\rm{V}}$, the
exchange splitting between VI and IVC phases. In the vicinity of $n_{e}^{*}$,
we take $\Delta_{\rm{V}}\propto n_{e}-n_{e}^{*}$, making $\Delta_{\rm{V}}$ a
reasonable proxy for the carrier density tuned in experiment. This assumption
is reproduced by Hartree-Fock calculations which show a large $\chi_{V}$ in
the vicinity of the VI to IVC phase transition (see Extended Data Fig. 10).
Despite its simplicity, the model reproduces all the key features of the
experimental data, and quantitatively reproduces the in- and out-of-plane
fields required to appreciably tune the balance between phases.

##  III Electron doping

Theoretically, the competition between VI and IVC states is governed by the
interplay of exchange and kinetic energy[26]; VI states effectively minimize
the exchange energy, but the momentum-dependent valley polarization of the IVC
phases reduces the kinetic energy in the presence of trigonal warping (see
Fig. 2D). For hole doping, where the Van Hove singularity occurs at a finite-
momentum saddle point, there is a large kinetic energy contribution from the
trigonally warped pockets near the band maximum, accounting for the dominance
of IVC order at lower hole density as observed in our experiment. This
contrasts with the situation for electron doping, where the bands have a
nearly quartic dispersion and the Van Hove singularity is very close to the
band edge for a large range of $D$. Within the electron band, the effects of
trigonal warping are most pronounced at high density, contrasting with the
hole band.

Fig. 4A shows $\kappa$ measured for electron doping ($n_{e}>0$). As at hole
doping, we observe a first-order phase boundary separating two phases with a
single Fermi surface, in addition to the previously reported[6] first-order
transitions separating the quarter metal region from an isospin symmetric
(Sym) phase at low $n_{e},D$ and a spin-polarized half metal state (SI) at
high $n_{e}$. Fig. 4B shows a phase diagram calculated within the Hartree-Fock
approximation, plotted as a function of $n_{e}$ and the interlayer potential
difference $\Delta$, which is proportional to $D$ (see Methods and
Supplementary information for further details). The simulations predict a
valley imbalanced quarter metal at low $n_{e}$ and a VI phase at intermediate
$n_{e}$, and we assign these labels to the two quarter metal phases in Fig.
4A.

To verify this hypothesis, we track the phase transitions along the dotted
line in Fig. 4A as a function of magnetic field. The VI phase has divergent
susceptibility to $B_{\perp}$, both as compared to the IVC phase and the Sym
phase (Fig. 4C). The behavior of the VI phase is consistent with a vanishing
in-plane susceptibility at $B_{\parallel}=0$,and the gradual canting of the
spin moment with increasing $B_{\parallel}$. This manifests at the VI to IVC
transition in Fig. 4D, where the curvature produces an estimate of
$\lambda\approx 32\mu$eV. The VI-Sym boundary is also consistent with this
mechanism, following a hyperbolic trajectory from which we extract
$\lambda\approx 65\mu$eV (a summary of the $\lambda$ fitting is presented in
the methods, and in Extended Data Fig. 11). These values of $\lambda$ are
consistent with each other and with literature values within the large
systematic errors (detailed in the Methods and supplementary information)
inherent to this fitting procedure. We note that the VI to IVC transition in
the electron-doped quarter metal was missed in the lower resolution
measurements of Ref. [6]; as a result, anomalous Hall measurements in the VI
phase and the divergent in-plane spin susceptibility of the IVC phase were
mistakenly attributed to a single, valley polarized phase with vanishing spin
anisotropy.

##  IV Discussion

Theoretically, intervalley coherent fluctuations may generate attractive
electronic interactions even in the absence of acoustic phonons. Concrete
scenarios for purely electronic superconductivity from this mechanism have
been proposed for both rhombohedral graphene multilayers[26, 28] and twisted
graphene multilayers[43]. In twisted bilayer and trilayer graphene IVC order
was recently detected[36, 37] in regions of the parameter space where
superconductivity also occurs[44, 45, 46]. While this scenario is appealing,
it has been in tension with the fact that superconductivity is so widespread
in graphene systems: superconductivity is observed in a broad range of twist
angles and densities in graphene multilayers[47], including far from where IVC
order has been detected[48, 49].

Our observations open up the possibility that intervalley coherence may be
just as widespread as superconductivity. Our measurements are able to detect
IVC order only in the quarter metal regime where the space of possible phases
is already significantly constrained. Outside this regime, we cannot
conclusively distinguish between intervalley coherence and other theoretically
predicted phases, for example, nematic orders with no net orbital
magnetization [50, 51]. of particular interest, both the phase represented in
gray in Fig. 1F, adjacent to the spin-triplet SC2, and the spin-unpolarized
phase adjacent to the Pauli-limited SC1, are compatible with intervalley
coherence given existing experimental constraints.

We also highlight the observation of a large valley susceptibility within the
IVC phase, despite the first order nature of the VI to IVC transition.
Superconductivity is observed on the disordered side of phase transitions to
potential IVC states[4]. Naively, the first order nature of these transitions
argues against the importance of fluctuations. However, our results suggest
that in rhombohedral graphene these transitions are sufficiently weakly first
order so as to allow for a large susceptibility.

We finally comment on the implications of finite spin-orbit coupling, which
has not been extensively discussed in the context of many-body physics in
graphene. A number of experimental reports in twisted graphene multilayers
have observed anomalously large effects of in-plane field on orbital
ferromagnetic order[52, 53]. This was presumed not to originate from intrinsic
spin-orbit coupling due to its negligible magnitude. However, our estimate for
$\lambda\approx 50\mu$eV is non-negligible when compared to the energy
differences between competing magnetic states—expected to be in the few 100
$\mu$eV range[26, 50, 54]—and is considerably larger than the superconducting
pairing energy inferred from the transition temperature in rhombohedral
multilayers[1, 2, 3, 4].

Spin-orbit coupling is likely to play a particularly important role in the
spin-triplet superconducting states observed in a variety of graphene
multilayers[1, 2, 55], selecting the direction of the spin of the Cooper pairs
and even forbidding certain phases. For example, in rhombohedral bi- and
trilayers triplet superconductivity arises from a spin-imbalanced phase,
identified from quantum oscillations as an annular half metal in the trilayer
and partially spin polarized phase in the bilayer. In the presence of spin-
orbit coupling, the phase diagram of these spin-imbalanced states is expected
to be similar to that of the IVC quarter metal shown in Fig. 3H, i.e. with
competing easy-axis and easy-plane spin ferromagnetic phases. Interestingly,
while easy-plane phases maintain valley balance, easy axis phases develop
finite valley polarization, which is strongly pair breaking for a
superconductor whose order parameter carries zero momentum. It is notable that
in bilayer graphene, superconductivity emerges at values of $B_{\parallel}$
where the spins become fully in-plane polarized, at the same field scale where
we find a IVCz to IVCx transition in the rhombohedral trilayer. It is possible
that superconductivity in graphene bilayers appears once the pair breaking
effect of spin-orbit coupling in phases with an out-of-plane spin moment is
quenched by the in-plane field. This hypothesis can be tested by measurements
of the in-plane and out of plane magnetic moments as a function of in-plane
field.

Note added: While finalizing this work, we became aware of a parallel
theoretical analysis predicting a transition between a valley-XY quarter metal
and a valley-Ising quarter metal – equivalent to the intervalley coherent and
valley imbalanced quarter metals discussed in this work, respectively. Ref.
[56] predicts a slope of the phase transition in $B_{\parallel}$ in rough
agreement with what we find experimentally at low in-plane fields, $\partial
n^{\ast}/\partial B_{\parallel}\approx-0.5\times 10^{11}$ cm-2 T-1.

Acknowledgements. The authors would like to acknowledge discussions with L.
Levitov, Z. Dong, A. Macdonald, C. Huang, M. Zaletel, S. Chatterjee, O. Vafek.
This project was primarily funded by the Department of Energy under DE-
SC0020043. A.F.Y. acknowledges the support of the Gordon and Betty Moore
Foundation under award GBMF9471 and the Packard Foundation under award
2016-65145 for general group activities. This research was supported in part
by the National Science Foundation under Grants No. NSF PHY-1748958 and
PHY-2309135. A portion of this work was performed in the UCSB Nanofabrication
Facility, an open access laboratory. This work made use of facilities funded
by Enabling Quantum Leap: Convergent Accelerated Discovery Foundries for
Quantum Materials Science, Engineering and Information (Q-AMASE-i) award
number DMR-1906325 from the National Science Foundation. K.W. and T.T.
acknowledge support from the Elemental Strategy Initiative conducted by the
MEXT, Japan (Grant Number JPMXP0112101001) and JSPSKAKENHI (Grant Numbers
19H05790, 20H00354 and 21H05233). E.B. was supported by the European Research
Council (ERC) under grant HQMAT (Grant Agreement No. 817799), the Israel-US
Binational Science Foundation (BSF), and by CRC 183 of the Deutsche
Forschungsgemeinschaft (project C02). CLT acknowledges support from the Hertz
Foundation and from the National Science Foundation Graduate Research
Fellowship Program under grant 1650114.

Author Contributions. A.F.Y. conceived of and directed the project. H.Z. and
T.X. designed and fabricated the device. C.L.T., T.A., O.S., G.B. and E.R.
fabricated nanoSQUID tips for the measurements. K.W. and T.T. grew the hBN
crystals. M.E.H., C.L.T., T.A., O.S. and E.R. developed the nanoSQUID
microscope. T.A., O.S., C.L.T. performed nanoSQUID measurements and T.A.,
O.S., H.Z., H.M.Y., C.L.P., and L.H. performed capacitance measurements. J.X.,
Y.V., T.H., and E. B. contributed to the theoretical interpretation and
performed the numerical simulations. T.A., O.S., E.B., and A.F.Y. wrote the
paper with inputs from all other authors.

##  V Materials and Methods

#### Compressibility and electronic transport measurements:

Ref. [6] and ref. [1] studied the same device, and details of the sample
fabrication, electronic transport, and compressibility measurements are
provided in those references. In contrast to those references, however, data
in Figs. 1B,E, 3A,E,F,G, 4A,C,D, 11, 6 were taken in an improved
compressibility set-up that introduced a decoupling capacitor between the
transistor amplifier and the sample as well as a second stage cryogenic
amplifier located at the 4K stage. This setup resulted in lower electronic
temperature and provided improved signal to noise, allowing for smaller
amplitude modulation and consequently better resolution of fine features (see
Ref. [4] for details).

#### nanoSQUID on tip measurements:

The nanoSQUID on tip (nSOT) microscopy was performed with an indium SQUID
fabricated at the tip of a pulled quartz pipette[57]. The field period of 71
mT corresponds to an effective diameter of 193 nm. The nSOT signal was
acquired using a series SQUID array amplifier (SSAA) in feedback mode, where
the feedback voltage $V_{\text{SSAA}}$ is proportional to the current through
the nSOT[57]. The SQUID was mechanically controlled with a piezoelectrically
pumped quartz tuning fork in a phase locked loop and positioned over the
sample by piezoelectric nanopositioners. Unless otherwise noted the magnetic
data was taken in constant height mode 150 nm above the surface of the sample,
corresponding to a height $h\approx 300$ nm above the ABC graphene layer.

At 300mK, the nSOT had good sensitivity up to $H_{c,\perp}\geq 500$ mT. During
operation the magnetic field and nSOT voltage bias are held constant at a
particular working point determined by the applied bias and magnetic field,
with the bias selected to maximize the sensitivity at a given field. For the
measurements shown here, the sensitivity varied between 25-35
nT/$\sqrt{\mathrm{Hz}}$. A small out-of-plane field is typically necessary to
get reasonable sensitivity: for Fig. 1C and 5A,
$(B_{\perp},B_{\parallel})=(15\mathrm{mT},0)$, and for Fig. 1D,E
$(B_{\perp},B_{\parallel})=(19\mathrm{mT},500\mathrm{mT})$.

As described in the main text, magnetometry measurements are performed by
applying a finite frequency (typically 0.5-4 kHz, chosen to minimize noise)
excitation to the bottom gate (with peak-peak amplitude $15$ mV for Fig. 1B
and $35$ mV otherwise) and measuring the SQUID response at the same frequency.
Images thus show contrast at the onset (and disappearance) of magnetic phases
in the $n_{e}$-$D$ plane. Fig. 1C shows a real-space image of the device
acquired at $n_{e}=-0.67\times 10^{12}$ cm-1 and $D=0.55$ V/nm.

####  Determination of $\lambda$:

In the main text Fig. 3A and Extended Data Fig. 11 we fit two different
functional forms for the curvature of spin-orbit affected transitions as a
function of $B_{\parallel}$. Three free parameters were included in the fit,
labeled below as $n_{\text{offset}}$, $\alpha$, and $\lambda$.

$n_{e}^{\ast}(B_{\parallel},n_{\text{offset}},\alpha,\lambda)=\left\\{\begin{array}[]{ll}n_{\text{offset}}+\alpha\left(B_{\parallel}+\lambda/g\mu_{B}-\sqrt{B_{\parallel}^{2}+(\lambda/g\mu_{B})^{2}}\right)&\text{Fig.~{}\ref{fig:3}A,
Fig.~{}\ref{fig:ext_data_electron}A,B,D}\\\
n_{\text{offset}}+\alpha\left(\lambda/g\mu_{B}-\sqrt{B_{\parallel}^{2}+(\lambda/g\mu_{B})^{2}}\right)&\text{Fig.~{}\ref{fig:ext_data_electron}C}\end{array}\right.$
(2)

Here, $n_{\text{offset}}$ represents a net density offset to account for the
finite density location of transitions, $\lambda$ is the SO coupling strength,
and $\alpha=\pm\frac{|m_{spin}|}{\Delta\mu}$. A complete derivation of the
functional forms shown here can be found in the supplementary information. The
location of the transition was determined by fitting compressibility data to a
Gaussian peak. In Fig. 11F, we give a complete list of extracted $\lambda$
values for all four studied transitions. Near zero field our data is affected
by the presence of additional phases, particularly the IVCz, which prevents
effective determination of $n_{e}^{\ast}$. At high $B_{\parallel}$, meanwhile,
the transitions saturate to a fixed $n^{*}_{e}$, and fits become insensitive
to $\lambda$. The dominant source of error thus arises from the choice of
range for the fit. To capture this uncertainty, we fit both the whole range as
well as all sub-ranges between the high- and low $B_{\parallel}$ limits that
span at least fifteen times the single-point error in $n^{*}_{e}$ as
determined from the standard deviation of the Gaussian peak fit. In Extended
Data Fig. 11, we report both the best fit to the whole curve as well as the
maximum and minimum fits from the sub-ranges, $\lambda_{\text{max}}$ and
$\lambda_{\text{min}}$. The best fit for the four transitions give values
between 30-65 $\mu$eV.

#### Determination of $\Delta m$:

Knowing $\Delta\mu$ at a phase boundary allows, in principle, for precise
determination of $\Delta m$. To measure $\Delta\mu$, we integrate our
thermodynamic compressibility $\kappa=d\mu/dn_{e}$ over the density range of
the negative $\kappa$ peak corresponding to a phase transition. In practice,
integrating our experimental data is challenging due to the weakness of the
first order transitions and differences between the compressibility of the two
adjoining phases. We thus expect an overall systematic error; for example,
systematic underestimates of $\Delta\mu$ may be responsible for the smaller
than expected values of $\Delta m_{\parallel}$ reported in Fig. 3B. However,
we do not expect this type of systematic error to be magnetic field dependent,
allowing comparison of the field dependence of $\Delta\mu$ and $\partial
n^{\ast}/\partial B$. This comparison is crucial to determining whether the
curvature observed in a phase transition arises from a change in $\Delta
m$–curvature may also be generated by a field-dependent $\Delta\mu$.

In the data shown in Fig. 3B, we integrate compressibility data to produce the
plotted $\Delta\mu$. This data shows no systematic trend above our noise
floor, so we assume a constant $\Delta\mu=$ over the depicted range,
corresponding to the average measured value. $\partial n^{\ast}/\partial B$ is
determined by numerically differentiating the $n_{e}^{*}$ values extracted
from peak fitting as described in the previous section. $\Delta m_{\parallel}$
is then calculated by multiplying this number by the average $\Delta\mu$.

## References

  * Zhou  _et al._ [2021a] H. Zhou, T. Xie, T. Taniguchi, K. Watanabe, and A. F. Young, Nature 598, 434 (2021a).
  * Zhou  _et al._ [2022] H. Zhou, L. Holleis, Y. Saito, L. Cohen, W. Huynh, C. L. Patterson, F. Yang, T. Taniguchi, K. Watanabe, and A. F. Young, Science 375, 774 (2022).
  * Zhang  _et al._ [2023] Y. Zhang, R. Polski, A. Thomson, E. Lantagne-Hurtubise, C. Lewandowski, H. Zhou, K. Watanabe, T. Taniguchi, J. Alicea, and S. Nadj-Perge, Nature 613, 268 (2023).
  * Holleis  _et al._ [2023] L. Holleis, C. L. Patterson, Y. Zhang, H. M. Yoo, H. Zhou, T. Taniguchi, K. Watanabe, S. Nadj-Perge, and A. F. Young, Ising Superconductivity and Nematicity in Bernal Bilayer Graphene with Strong Spin Orbit Coupling (2023), arXiv:2303.00742 [cond-mat].
  * Shi  _et al._ [2020] Y. Shi, S. Xu, Y. Yang, S. Slizovskiy, S. V. Morozov, S.-K. Son, S. Ozdemir, C. Mullan, J. Barrier, J. Yin, A. I. Berdyugin, B. A. Piot, T. Taniguchi, K. Watanabe, V. I. Fal’ko, K. S. Novoselov, A. K. Geim, and A. Mishchenko, Nature 584, 210 (2020).
  * Zhou  _et al._ [2021b] 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, Nature 598, 429 (2021b).
  * de la Barrera  _et al._ [2022] S. C. de la Barrera, S. Aronson, Z. Zheng, K. Watanabe, T. Taniguchi, Q. Ma, P. Jarillo-Herrero, and R. Ashoori, Nature Physics 18, 771 (2022).
  * Seiler  _et al._ [2022] A. M. Seiler, F. R. Geisenhof, F. Winterer, K. Watanabe, T. Taniguchi, T. Xu, F. Zhang, and R. T. Weitz, Nature 608, 298 (2022).
  * Lin  _et al._ [2023] J.-X. Lin, Y. Wang, N. J. Zhang, K. Watanabe, T. Taniguchi, L. Fu, and J. I. A. Li, Spontaneous momentum polarization and diodicity in Bernal bilayer graphene (2023), arXiv:2302.04261 [cond-mat].
  * Seiler  _et al._ [2023] A. M. Seiler, M. Statz, I. Weimer, N. Jacobsen, K. Watanabe, T. Taniguchi, Z. Dong, L. S. Levitov, and R. T. Weitz, Interaction-driven (quasi-) insulating ground states of gapped electron-doped bilayer graphene (2023), arXiv:2308.00827 [cond-mat].
  * Han  _et al._ [2023a] T. Han, Z. Lu, G. Scuri, J. Sung, J. Wang, T. Han, K. Watanabe, T. Taniguchi, L. Fu, H. Park, and L. Ju, Orbital Multiferroicity in Pentalayer Rhombohedral Graphene (2023a), arXiv:2308.08837 [cond-mat].
  * Kane and Mele [2005] C. L. Kane and E. J. Mele, Physical Review Letters 95, 226801 (2005).
  * Min  _et al._ [2006] H. Min, J. E. Hill, N. A. Sinitsyn, B. R. Sahu, L. Kleinman, and A. H. MacDonald, Physical Review B 74, 165310 (2006).
  * Yao  _et al._ [2007] Y. Yao, F. Ye, X.-L. Qi, S.-C. Zhang, and Z. Fang, Physical Review B 75, 041401 (2007).
  * Konschuh  _et al._ [2012] S. Konschuh, M. Gmitra, D. Kochan, and J. Fabian, Physical Review B 85, 115423 (2012).
  * Cvetkovic  _et al._ [2012] V. Cvetkovic, R. E. Throckmorton, and O. Vafek, Physical Review B 86, 075467 (2012).
  * Mayorov  _et al._ [2011] A. S. Mayorov, D. C. Elias, M. Mucha-Kruczynski, R. V. Gorbachev, T. Tudorovskiy, A. Zhukov, S. V. Morozov, M. I. Katsnelson, A. K. Geim, and K. S. Novoselov, Science 333, 860 (2011).
  * Winterer  _et al._ [2023] F. Winterer, F. R. Geisenhof, N. Fernandez, A. M. Seiler, F. Zhang, and R. T. Weitz, Ferroelectric and anomalous quantum Hall states in bare rhombohedral trilayer graphene (2023), arXiv:2305.04950 [cond-mat].
  * Han  _et al._ [2023b] T. Han, Z. Lu, G. Scuri, J. Sung, J. Wang, T. Han, K. Watanabe, T. Taniguchi, H. Park, and L. Ju, Correlated Insulator and Chern Insulators in Pentalayer Rhombohedral Stacked Graphene (2023b), arXiv:2305.03151 [cond-mat].
  * Island  _et al._ [2019] J. O. Island, X. Cui, C. Lewandowski, J. Y. Khoo, E. M. Spanton, H. Zhou, D. Rhodes, J. C. Hone, T. Taniguchi, K. Watanabe, L. S. Levitov, M. P. Zaletel, and A. F. Young, Nature 571, 85 (2019).
  * You and Vishwanath [2022] Y.-Z. You and A. Vishwanath, Physical Review B 105, 134524 (2022).
  * Ghazaryan  _et al._ [2021] A. Ghazaryan, T. Holder, M. Serbyn, and E. Berg, Physical Review Letters 127, 247001 (2021).
  * Ghazaryan  _et al._ [2023] A. Ghazaryan, T. Holder, E. Berg, and M. Serbyn, Physical Review B 107, 104502 (2023).
  * Cea  _et al._ [2022] T. Cea, P. A. Pantaleón, V. T. Phong, and F. Guinea, Physical Review B 105, 075432 (2022).
  * Wagner  _et al._ [2023] G. Wagner, Y. H. Kwan, N. Bultinck, S. H. Simon, and S. A. Parameswaran, Superconductivity from repulsive interactions in Bernal-stacked bilayer graphene (2023), arXiv:2302.00682 [cond-mat].
  * Chatterjee  _et al._ [2022] S. Chatterjee, T. Wang, E. Berg, and M. P. Zaletel, Nature Communications 13, 6013 (2022).
  * Dong and Levitov [2021] Z. Dong and L. Levitov, arXiv:2109.01133 [cond-mat] (2021), arXiv: 2109.01133 [cond-mat].
  * Dong  _et al._ [2023] Z. Dong, P. A. Lee, and L. S. Levitov, Signatures of Cooper pair dynamics and quantum-critical superconductivity in tunable carrier bands (2023), arXiv:2304.09812 [cond-mat].
  * Dong  _et al._ [2022] Z. Dong, A. V. Chubukov, and L. Levitov, Spin-triplet superconductivity at the onset of isospin order in biased bilayer graphene (2022), arXiv:2205.13353 [cond-mat].
  * Chou  _et al._ [2021] Y.-Z. Chou, F. Wu, J. D. Sau, and S. D. Sarma, Acoustic-phonon-mediated superconductivity in rhombohedral trilayer graphene (2021), arXiv:2106.13231 [cond-mat].
  * Eisenstein  _et al._ [1992] J. P. Eisenstein, L. N. Pfeiffer, and K. W. West, Phys. Rev. Lett. 68, 674 (1992).
  * Ashoori  _et al._ [1992] R. C. Ashoori, H. L. Stormer, J. S. Weiner, L. N. Pfeiffer, S. J. Pearton, K. W. Baldwin, and K. W. West, Phys. Rev. Lett. 68, 3088 (1992).
  * Das and Huang [2023a] M. Das and C. Huang, Unconventional Metallic Magnetism: Non-analyticity and Sign-changing Behavior of Orbital Magnetization in ABC Trilayer Graphene (2023a), arXiv:2308.01996 [cond-mat].
  * Coissard  _et al._ [2022] A. Coissard, D. Wander, H. Vignaud, A. G. Grushin, C. Repellin, K. Watanabe, T. Taniguchi, F. Gay, C. B. Winkelmann, H. Courtois, H. Sellier, and B. Sacépé, Nature 605, 51 (2022).
  * Liu  _et al._ [2022] X. Liu, G. Farahi, C.-L. Chiu, Z. Papic, K. Watanabe, T. Taniguchi, M. P. Zaletel, and A. Yazdani, Science 375, 321 (2022).
  * Nuckolls  _et al._ [2023] K. P. Nuckolls, R. L. Lee, M. Oh, D. Wong, T. Soejima, J. P. Hong, D. Călugăru, J. Herzog-Arbeitman, B. A. Bernevig, K. Watanabe, T. Taniguchi, N. Regnault, M. P. Zaletel, and A. Yazdani, Nature 620, 525 (2023).
  * Kim  _et al._ [2023] H. Kim, Y. Choi, E. Lantagne-Hurtubise, C. Lewandowski, A. Thomson, L. Kong, H. Zhou, E. Baum, Y. Zhang, L. Holleis, K. Watanabe, T. Taniguchi, A. F. Young, J. Alicea, and S. Nadj-Perge, Imaging inter-valley coherent order in magic-angle twisted trilayer graphene (2023), arXiv:2304.10586 [cond-mat].
  * Zhang  _et al._ [2009] Y. Zhang, T.-T. Tang, C. Girit, Z. Hao, M. C. Martin, A. Zettl, M. F. Crommie, Y. R. Shen, and F. Wang, Nature 459, 820 (2009).
  * Ochoa  _et al._ [2012] H. Ochoa, A. H. Castro Neto, V. I. Fal’ko, and F. Guinea, Physical Review B 86, 245411 (2012).
  * Sichau  _et al._ [2019] J. Sichau, M. Prada, T. Anlauf, T. Lyon, B. Bosnjak, L. Tiemann, and R. Blick, Physical Review Letters 122, 046403 (2019).
  * Banszerus  _et al._ [2020] L. Banszerus, B. Frohn, T. Fabian, S. Somanchi, A. Epping, M. Müller, D. Neumaier, K. Watanabe, T. Taniguchi, F. Libisch, B. Beschoten, F. Hassler, and C. Stampfer, Physical Review Letters 124, 177701 (2020).
  * Kurzmann  _et al._ [2021] A. Kurzmann, Y. Kleeorin, C. Tong, R. Garreis, A. Knothe, M. Eich, C. Mittag, C. Gold, F. K. de Vries, K. Watanabe, T. Taniguchi, V. Fal’ko, Y. Meir, T. Ihn, and K. Ensslin, Nature Communications 12, 6004 (2021).
  * Lee  _et al._ [2019] J. Y. Lee, E. Khalaf, S. Liu, X. Liu, Z. Hao, P. Kim, and A. Vishwanath, Nature Communications 10, 5333 (2019), number: 1 Publisher: Nature Publishing Group.
  * Cao  _et al._ [2018] Y. Cao, V. Fatemi, S. Fang, K. Watanabe, T. Taniguchi, E. Kaxiras, and P. Jarillo-Herrero, Nature 556, 43 (2018).
  * Park  _et al._ [2021] J. M. Park, Y. Cao, K. Watanabe, T. Taniguchi, and P. Jarillo-Herrero, Nature 590, 249 (2021).
  * Hao  _et al._ [2021] Z. Hao, A. M. Zimmerman, P. Ledwith, E. Khalaf, D. H. Najafabadi, K. Watanabe, T. Taniguchi, A. Vishwanath, and P. Kim, Science 371, 1133 (2021).
  * Arora  _et al._ [2020] H. S. Arora, R. Polski, Y. Zhang, A. Thomson, Y. Choi, H. Kim, Z. Lin, I. Z. Wilson, X. Xu, J.-H. Chu, K. Watanabe, T. Taniguchi, J. Alicea, and S. Nadj-Perge, Nature 583, 379 (2020).
  * Zhang  _et al._ [2022] Y. Zhang, R. Polski, C. Lewandowski, A. Thomson, Y. Peng, Y. Choi, H. Kim, K. Watanabe, T. Taniguchi, J. Alicea, F. von Oppen, G. Refael, and S. Nadj-Perge, Science 377, 1538 (2022).
  * Park  _et al._ [2022] J. M. Park, Y. Cao, L.-Q. Xia, S. Sun, K. Watanabe, T. Taniguchi, and P. Jarillo-Herrero, Nature Materials 21, 877 (2022).
  * Koh  _et al._ [2023] J. M. Koh, J. Alicea, and E. Lantagne-Hurtubise, Correlated Phases in Spin-Orbit-Coupled Rhombohedral Trilayer Graphene (2023), arXiv:2306.12486 [cond-mat].
  * Zhumagulov  _et al._ [2023] Y. Zhumagulov, D. Kochan, and J. Fabian, Emergent correlated phases in rhombohedral trilayer graphene induced by proximity spin-orbit and exchange coupling (2023), arXiv:2305.14277 [condmat].
  * Sharpe  _et al._ [2021] A. L. Sharpe, E. J. Fox, A. W. Barnard, J. Finney, K. Watanabe, T. Taniguchi, M. A. Kastner, and D. Goldhaber-Gordon, Nano Letters 21, 4299 (2021).
  * Kuiri  _et al._ [2022] M. Kuiri, C. Coleman, Z. Gao, A. Vishnuradhan, K. Watanabe, T. Taniguchi, J. Zhu, A. H. MacDonald, and J. Folk, Nature Communications 13, 6468 (2022).
  * Huang  _et al._ [2023] C. Huang, T. M. R. Wolf, W. Qin, N. Wei, I. V. Blinov, and A. H. MacDonald, Physical Review B 107, L121405 (2023).
  * Cao  _et al._ [2021] Y. Cao, J. M. Park, K. Watanabe, T. Taniguchi, and P. Jarillo-Herrero, Nature 595, 526 (2021).
  * Das and Huang [2023b] M. Das and C. Huang, Quarter-Metal Phases in Multilayer Graphene: Ising-XY and Annular Lifshitz Transitions (2023b), arXiv:2310.10759 [cond-mat].
  * Anahory  _et al._ [2020] Y. Anahory, H. R. Naren, E. O. Lachman, S. B. Sinai, A. Uri, L. Embon, E. Yaakobi, Y. Myasoedov, M. E. Huber, R. Klajn, and E. Zeldov, Nanoscale 12, 3174 (2020).
  * Zhang  _et al._ [2010] F. Zhang, B. Sahu, H. Min, and A. H. MacDonald, Physical Review B 82, 035409 (2010).
  * Kudin  _et al._ [2002] K. N. Kudin, G. E. Scuseria, and E. Cancès, The Journal of Chemical Physics 116, 8255 (2002).

Extended Data Figures

Figure 5: Half metal ferromagnetism. (A) The nSOT phase diagram of the hole
side versus density and displacement field. Most nSOT data in this work
focuses on magnetism in the vicinity of the VI phase, in this figure we focus
on magnetism near the SI transition, marked by black circle. (B) Optical image
of device. Scale bar is 3 $\mu$m. (C) Spatial images of $\delta B$ for
different directions of in-plane and out-of-plane field with $n_{e}$ and D set
on the half metal transition, marked with a back circle in A. Marked location
indicates where the the data in Fig. 1D was taken. Scale bars are 1 $\mu$m.
(D) Simulated magnetism from a sheet of dipoles in the shape of the sample
with magnetic moment oriented in the direction of the fields in D and the
signal normalized to correspond to the total density $n_{e}$ with each carrier
having $1\mu_{B}$ of magnetic moment. The qualitative and approximate
quantitative agreement of the data in D and the simple dipole model in E
indicates that magnetism arises due to the spins of the carriers aligning with
the external field. Scale bars are 1 $\mu$m.  Figure 6: Electronic Phases in
the vicinity of SC2 (A) Inverse compressibilty as a function of $B_{\perp}$
and $n_{e}$ in the vicinity of the quarter metal at $B_{\perp}$ = 0. (B)
Schematic phase diagram indicating the electronic phases of A as a function of
$n_{e}$. Solid lines between different colored areas represent first order
phase transitions in $\kappa$ between electronic phases: Sym, VI, IVC and SI.
Dashed lines indicate Lifshitz transitions between disjointed (d), simple (s)
and annular (a) Fermi surfaces as marked by lower case labels. There is a
superconducting pocket (SC2) at the boundary between the IVC and SI phases as
marked.[1] (C) Inverse compressibility, $\kappa$, as a function of carrier
density $n_{e}$ and versus $B_{\parallel}$ along the dashed line in A.

Figure 7:  Hartree Fock phase diagram for hole doping. (A) Phase diagram
without SOC, showing the ground state isospin polarization as a function of
$n_{e}$ and the interlayer potential $\Delta_{1}$. (B) Phase diagram of Fermi
surface topology from the same calculation. Regions that break the C3 symmetry
(i.e., nematics) are indicated by white dots. “Disjoint” refers to any case
where there are multiple Fermi surfaces with different areas, excluding
annuli. (C) Flavor occupation of the Hartree-Fock ground state at the points
labeled in panels A, B. The color indicates occupation expectation values. The
isospin order, Fermi surface topology, and nematicity in these eight points is
summarized in the table below. Note that in point no. 5 there are two distinct
Fermi surfaces with different areas, and hence the FS topology is “disjoint”
in our scheme.

Point | Isospin order | Topology | Nematic
---|---|---|---
1 | VI, SI | Annular | ✗
2 | SI, IVC | Annular | ✗
3 | SI, IVC | Annular | ✗
4 | VI, SI | Disjoint | ✗
5 | SI, IVC | Disjoint | ✓
6 | SI | Annular | ✗
7 | SI | Simple | ✗
8 | SI, IVC | Disjoint | ✗

Figure 8: Longitudinal Resistance of Quarter Metals. (A) $R_{xx}$ for hole
doping at $B_{\parallel}=0$T at T=20mK. The IVC and VI phases are not distinct
in longitudinal resistivity. (B) $R_{xx}$ compared to $\kappa$ as a function
of $n_{e}$ and $D$ along the trajectory indicated in A. The first order phase
transitions between VI and IVC do not generate contrast in $R_{xx}$. (C)
$R_{xx}$ along the trajector in panel A, as a function of $B_{\perp}$ with
$B_{\parallel}=0$. As explained in the main text, this range largely spans the
VI and $IVC_{z}$ phases. (D) FFT of the data outlined by the white box in B,
with the FFT taken as a function of $1/B_{\perp}$ at each value of $n_{e}$.
$f_{\nu}$ is the FFT frequency normalized to $n_{e}$.  Figure 9: Anomalous
Hall effect. (A) Resistance measured in an R${}_{\text{XY}}$-like
configuration at different points in the quarter metal phase for $T\approx 20$
mK and B${}_{\parallel}=1$ T. Anomalous Hall effect is observed only in the VI
phase. (B) Resistance in the quarter metal phase at $T\approx 20$ mK and
B${}_{\parallel}=0$ T. Both VI and $IVC_{z}$ phases show anomalous Hall
signal.  Figure 10:  Effects of spin-orbit coupling in the IVC phase from
perturbation theory. (A) Spin-valley susceptibility for the IVCx and IVCz
phases in the absence of SOC, as a function of carrier density at
$\Delta_{1}=17\,\rm{meV}$. (B) Stability of the IVCx and IVCz phases within
the IVC region ‘2’ in panel A. Between IVCx and IVCz, the state with larger
$\chi_{\tau_{z}s_{z}}$ (panel B) has lower energy in the presence of weak but
non-zero SOC.  Figure 11: Estimation of $\lambda$ from different transitions.
(A) SI to VI transition at hole doping, at $D=0.46\,\text{V}\,\text{nm}^{-1}$.
White points represent extracted transition location; the dark line indicates
the fit described in the Methods. (B) fit for the VI to IVC transition at hole
doping, (C) fit for the Sym to VI transition at electron doping, for
$D=-0.23\,\text{V}\,\text{nm}^{-1}$, and (D) VI to IVC transition for electron
doping for $D=-0.23\,\text{V}\,\text{nm}^{-1}$. The extracted
$(\lambda_{min},\lambda_{fit},\lambda_{max})$, in $\mu eV$, for the
transitions are A, $(25,46,60)$; B, $(39,56,98)$; C, $(54,65,174)$; D,
$(29,32,52)$.

Supplementary information

###  .1 Determination of $\lambda$

To elucidate the effects of spin-orbit coupling on the in-plane spin
polarization of a VI phase we consider a toy model given by

$H=\frac{1}{2}\lambda\tau_{z}s_{z}-\frac{1}{2}E_{z,\parallel}s_{x}$ (S1)

where $\tau$ and $s$ are the valley and spin Pauli matrices, respectively,
$E_{z,\parallel}=g\mu_{B}B_{\parallel}$, and we take $g=2$ for the spin
gyromagnetic ratio. The resulting in-plane magnetization, $m_{\parallel}$, of
a valley polarized phase where $\langle\tau_{z}\rangle=1$ is given by

$m_{\parallel}=|m_{spin}|\frac{E_{z}}{\sqrt{E_{z}^{2}+\lambda^{2}}}$ (S2)

In our experiment, we measure the magnetic field evolution of the density,
$n_{e}^{*}$, and displacement field, $D^{*}$, at which we observe first order
phase transitions. The condition for a first order phase transition between
two phases with contrasting in-plane magnetic moments is
$E_{1}(n_{e}^{*},D^{*})+m_{1,\parallel}B_{\parallel}=E_{2}(n_{e}^{*},D^{*})+m_{2,\parallel}B_{\parallel}$.
Experiments are performed by sweeping both $n$ and $D$ along a trajectory
parameterized by $t-t_{0}=n+D/a$, where $t_{0}$ and $a$ are constants. It
follows that, at a given $B_{\parallel}$, the magnetic moment difference

$\Delta
m_{\parallel}=\frac{\Delta\mu_{1,2}+a\,\partial_{D}(E_{1}-E_{2})}{\Big{[}\frac{\partial
B_{\parallel}}{\partial n^{\ast}}+a\,\frac{\partial B_{\parallel}}{\partial
D^{\ast}}\Big{]}}$ (S3)

For the three phase boundaries studied in Fig. 11A, C, and D, $a=0$. For the
transition studied in Fig. 11B, we make the simplifying assumption that the
two phases have similar polarizability $\partial E/\partial D$, and that
$\partial B_{\parallel}/\partial D^{*}$ is negligible. As described in the
main text, $\Delta\mu$ and $n_{e}^{*}(B_{\parallel})$ may both be measured
independently. We find that $\Delta\mu$ is independent of $B_{\parallel}$ to
within our experimental resolution for all four transitions, and so we may
take it as a constant in our analysis. It follows that

$\displaystyle\frac{dn_{e}^{*}}{dB_{\parallel}}$ $\displaystyle=\frac{\Delta
m_{\parallel}}{\Delta\mu}$ (S4) $\displaystyle\Delta n_{e}^{*}(B_{\parallel})$
$\displaystyle=\frac{1}{\Delta\mu}\int_{0}^{B_{\parallel}}\Delta
m_{\parallel}dB_{\parallel}$ (S5)

To fit the phase boundaries between VI and SI or IVC phases (Figs. 11A,B, and
D), we assume $m_{\parallel}=m_{spin}$ for IVC and SI phases. Then

$\displaystyle\Delta n_{e}^{*}$
$\displaystyle=\pm\frac{|m_{spin}|/(g\mu_{B})}{\Delta\mu}\int_{0}^{B_{\parallel}}\left(1-\frac{E_{z}}{\sqrt{E_{Z}^{2}+\lambda^{2}}}\right)dE_{z}$
(S6) $\displaystyle\Delta n_{e}^{*}$
$\displaystyle=\pm\frac{|m_{spin}|/(g\mu_{B})}{\Delta\mu}\left(E_{Z}+\lambda-\sqrt{E_{Z}^{2}+\lambda^{2}}\right)$
(S7)

To fit the phase boundaries between VI and Sym phase (Fig. 11B), we assume
$m_{\parallel}=0$ for the Sym phase. Then

$\displaystyle\Delta n_{e}^{*}$
$\displaystyle=\pm\frac{|m_{spin}|/(g\mu_{B})}{\Delta\mu}\int_{0}^{B_{\parallel}}\left(-\frac{E_{z}}{\sqrt{E_{Z}^{2}+\lambda^{2}}}\right)dE_{z}$
(S8) $\displaystyle\Delta n_{e}^{*}$
$\displaystyle=\pm\frac{|m_{spin}|/(g\mu_{B})}{\Delta\mu}\left(\lambda-\sqrt{E_{Z}^{2}+\lambda^{2}}\right)$
(S9)

###  .2 Band structure and interactions

For the band structure of rhombohedral trilayer graphene, we use the six band
model of Ref. [58]:

$H_{0}(\bm{k},\xi)=\left(\begin{matrix}\Delta_{1}+\Delta_{2}+\delta&\frac{1}{2}\gamma_{2}&v_{0}\pi^{*}&v_{4}\pi^{*}&v_{3}\pi&0\\\
\frac{1}{2}\gamma_{2}&\Delta_{2}-\Delta_{1}+\delta&0&v_{3}\pi^{*}&v_{4}\pi&v_{0}\pi\\\
v_{0}\pi&0&\Delta_{1}+\Delta_{2}&\gamma_{1}&v_{4}\pi^{*}&0\\\
v_{4}\pi&v_{3}\pi&\gamma_{1}&-2\Delta_{2}&v_{0}\pi^{*}&v_{4}\pi^{*}\\\
v_{3}\pi^{*}&v_{4}\pi^{*}&v_{4}\pi&v_{0}\pi&-2\Delta_{2}&\gamma_{1}\\\
0&v_{0}\pi^{*}&0&v_{4}\pi&\gamma_{1}&\Delta_{2}-\Delta_{1}\end{matrix}\right),$
(S10)

where $\pi=\xi k_{x}+ik_{y}$ ($\xi=\pm$ corresponds to valleys $K$ and
$K^{\prime}$) and the Hamiltonian is written in the basis
$(A_{1},B_{3},B_{1},A_{2},B_{2},A_{3})$, where $A_{i}$ and $B_{i}$ label the
two sublattices at layer $i$. The velocities $v_{i}$ ($i=0,3,4$) are related
to the microscopic hopping parameters $\gamma_{i}$ by
$v_{i}=\sqrt{3}a_{0}\gamma_{i}/2$, where $a_{0}=2.46$ Åis the lattice constant
of monolayer graphene. $\Delta_{1}$ is a potential difference between outer
layers, which is approximately proportional to the applied displacement field,
while $\Delta_{2}$ is the potential difference between the middle layer and
the average potential of the outer layers. Finally, $\delta$ is an on-site
potential on $A_{1}$ and $B_{3}$. We have used the following parameters:
$\gamma_{0}=3.1\,{\rm eV},$$\gamma_{1}=0.38\,{\rm eV},$
$\gamma_{2}=-0.015\,{\rm eV}$, $\gamma_{3}=-0.29\,{\rm eV}$,
$\gamma_{4}=-0.141\,{\rm eV}$, $\delta=-0.0105\,{\rm eV}$,
$\Delta_{2}=-0.0023\,{\rm eV}$. See Ref. [6] for the procedure used to fix
these parameters.

The total Hamiltonian used in our theoretical analysis is

$\mathcal{H}=\sum_{\bm{k},\xi=K,K^{\prime},s=\uparrow,\downarrow}\psi_{\bm{k},\xi,s}^{\dagger}H_{0}(\bm{k},\xi)\psi_{\bm{k},\xi,s}+\mathcal{H}_{{\rm
C}}+\mathcal{H}_{{\rm Hunds}}+\mathcal{H}_{{\rm SOC}},$ (S11)

where $\psi_{\bm{k},\xi,s}$ is a six component spinor of annihilation
operators $\psi_{\bm{k},\xi,s,\sigma_{i}}$ that annihilate electrons with
momentum $\bm{k}$, valley index $\xi$, spin $s$, and sublattice index
$\sigma_{i}=A_{i},B_{i}$ on layer $i$. $H_{{\rm C}}$ is the long-ranged
Coulomb interaction,

$\mathcal{H}_{{\rm
C}}=\frac{1}{A}\sum_{\bm{q}}V_{\bm{q}}:\hat{n}_{\bm{q}}\hat{n}_{-\bm{q}}:$

where
$\hat{n}_{\bm{q}}=\sum_{\bm{k},\xi,s}\psi^{\dagger}_{\bm{k+q},\xi,s}\psi_{\bm{k},\xi,s}$.
For the screened Coulomb interaction, we use $V_{\bm{q}}=2\pi
e^{2}\tanh(qd)/(\epsilon q)$. Here $d=40\rm{nm}$ is the distance to the gates,
and for the dielectric constant we employ $\epsilon=24$, which incorporates
both the dielectric constant of the BN substrate and the effects of screening
in the graphene trilayer [22, 26, 54, 50]. The spatially short-ranged Hund’s
coupling is given by

$\mathcal{H}_{{\rm
Hunds}}=\frac{J_{H}}{A}\sum_{\bm{k},\bm{q}}:\bm{\hat{s}}_{\bm{q},+}\cdot\bm{\hat{s}}_{-\bm{q},-}:,$
(S12)

where the spin density operator in valley $\xi=\pm 1$ is given by
$\bm{\hat{s}}_{\bm{q},\xi}=\sum_{\bm{k},s_{1},s_{2}}\psi^{\dagger}_{\bm{k+q},\xi,s_{1}}\bm{s}_{s_{1},s_{2}}\psi_{\bm{k},\xi,s_{2}}$
($\bm{s}_{s_{1},s_{2}}$ are Pauli matrices spin space). As mentioned in the
main text, the intrinsic spin-orbit coupling $\mathcal{H}_{{\rm
SOC}}=\sum_{\bm{k},\xi,s}\psi_{\bm{k},\xi,s}^{\dagger}H_{\rm
KM}\psi_{\bm{k},\xi,s}$ is of the Kane-Mele type, $H_{{\rm
KM}}=\lambda\sigma_{z}\tau_{z}s_{z}/2$, where the Pauli matrices denote
sublattice ($\sigma_{z}$), spin ($s_{z}$) and valley ($\tau_{z}$).

###  .3 Hartree-Fock analysis

We seek mean-field solutions for the interacting Hamiltonian Eq. (S11) which
lift the degeneracy between the four flavors spanned by the spin and valley
degrees of freedom within one given band, either on the electron or hole side.
To this end, we project the interactions onto a single band nearest to the
Fermi level, followed by a self-consistent determination of a momentum-
dependent density matrix elements in the $4\times 4$ flavor subspace [26].
Upon projection, the form factors of the band appear in the interaction
Hamiltonian.

The Hartree-Fock calculations are performed using a rectangular momentum grid
with size $61^{2}$, with spacing $0.003\cdot 2\pi/a$, at zero temperature. The
self-consistency loop is implemented using an adaptive-step algorithm
following Ref. [59], which uses for each iteration step the linear
interpolation between the previous solution and the updated solution and
determines the point with minimal energy.

For spin and valley polarized states, different solutions can be enforced by
selectively disabling certain matrix elements of the density matrix. In
contrast, intervalley coherent states always compete with polarized solutions,
and are thus more fragile numerically. The typical energy difference between a
flavor symmetric state and flavor polarized ones is of the order of $1\,{\rm
meV}$ per electron, while the IVC state can additionally gain around
$0.05\,{\rm meV}$ over the polarized solutions. Lastly, within the IVC phase,
the energy differences between the in-plane and out-of-plane spin
polarizations (IVCx vs. IVCz) are of order $1\,{\rm\mu eV}$ for
$\lambda=40\,{\rm\mu eV}$, which is close to the energy resolution of our
calculation (determined by how well converged the Hartree-Fock solution is).
For this reason, we additionally calculate the spin-valley polarization
susceptibilities for both IVCx and IVCz, defined as
$\chi_{\tau^{z}s^{z}}\langle\tau_{z}s_{z}\rangle/\lambda$. For small
$\lambda$, the energy gain from SOC in the two phases is then $\Delta
E_{a=X,Z}=-\frac{1}{2}\chi_{\tau^{z}s^{z}}(\lambda/2)^{2}$, and the phase
whose susceptibility is larger is favored. We found that this method is more
accurate than calculating the energy directly. An additional benefit of this
method is that it is not necessary to specify a value for $\lambda$, as long
as $\lambda$ is sufficiently small. We find that both the free energy and the
spin-valley susceptibility lead to the same overall trends in the resulting
phase diagrams, with only small quantitative differences.

For the determination of the phase diagram, it is necessary to first fix the
value of $J_{H}$. By testing a few values of the Hund’s coupling starting from
$J_{H}=0$, up to $J_{H}=-0.02V_{\bm{q}=0}$ while $\lambda$ is kept zero, we
find that already a small $J_{H}=-0.011V_{0}$ significantly favors the
formation of IVC states. Using this value for the Hund’s coupling, one can
then determine the spin-valley susceptibility inside the IVC region of the
phase diagram.

The numerical results of this calculation for the electron and the hole side
are shown in the main text Fig. 4B and extended data Fig. 7 respectively. The
regions of stability of the IVCx and IVCz within the IVC quarter metal phase
are shown in Fig. 10. Also shown is the spin-valley susceptibilities of the
two phases along a cut through the phase diagram, showing the large
enhancement of the susceptibility in the IVCZ state when the phase boundary to
the VI is approached. In the following section we show, using a
phenomenological ansatz for the free energy, how the various trends in these
phase diagrams can be understood.

###  .4 Model for IVCx to IVCz transition in quarter metal

To understand the sequence of transitions that occurs between the IVC and the
VI states within the quarter metal phase, it is useful to consider a
simplified model that includes the interplay between exchange, Hund’s rule,
and spin-orbit interactions. Within the Hartree-Fock approximation, each
occupied $\bm{k}$ orbital is polarized in a certain “direction” in spin and
valley space. In the simplified model, we assume that this isospin
polarization is nearly $\bm{k}-$independent. We may then write the system’s
energy as a function of a single state $|\Psi\rangle$ in spin and valley
space,

$|\Psi\rangle=\alpha|K\uparrow\rangle+\beta|K\downarrow\rangle+\gamma|K^{\prime}\uparrow\rangle+\delta|K^{\prime}\downarrow\rangle,$
(S13)

parameterized by the four-component spinor
$\psi^{T}=(\alpha,\beta,\gamma,\delta)$.

In general, the energy density $E$ in the quarter metal can be expressed as a
function of $\psi$. We use the following phenomenological model for $E[\psi]$:

$E[\psi]=\Delta_{\rm{V}}\langle\tau^{z}\rangle^{2}+\Delta^{\prime}_{\rm{V}}\langle\tau^{z}\rangle^{4}-J_{H}\langle\bm{s}_{+}\rangle\cdot\langle\bm{s}_{-}\rangle-\frac{\lambda}{2}\langle
s^{z}\tau^{z}\rangle-\frac{g\mu_{B}\bm{B}}{2}\cdot\langle\bm{s}\rangle-\frac{g_{{\rm
v}}B_{\perp}}{2}\langle\tau^{z}\rangle.$ (S14)

Here $\langle O\rangle\equiv\psi^{\dagger}O\psi$, where $O$ is a $4\times 4$
matrix in spin- and valley space. $\bm{\tau}$ and $\bm{s}$ are Pauli matrices
in valley and spin space, respectively, and
$\bm{s}_{\pm}=\frac{1\pm\tau^{z}}{2}\bm{s}$ are the spin operators projected
to the two valleys.

The parameters $\Delta_{\rm{V}}$, $\Delta^{\prime}_{\rm{V}}$ represent the
anisotropy of the energy in valley space, proportional to $n_{e}$. These
parameters come from the kinetic energy and the long-ranged part of the
Coulomb interactions. They are symmetric under
$\rm{U_{+}(2)}\times\rm{U_{-}(2)}$ transformations corresponding to the
approximate conservation of charge and spin in each valley separately. $J_{H}$
is the Hund’s coupling between the spins in the two valleys (assumed to be
ferromagnetic, $J_{H}\geq 0$). $\lambda$ is the strength of the Ising spin-
orbit coupling. $J_{H}$ and $\lambda$ reduce the symmetry down to
$\rm{U_{c}}(1)\times\rm{U_{v}}(1)\times\rm{U}_{s^{z}}(1)$ (corresponding to
the conservation of the charge, valley, and $z$ component of the total spin).
$\bm{B}$ is the magnetic field measured in energy units. The field couples
both to the spin via the Zeeman effect, and to the perpendicular valley
orbital moment (proportional to $\langle\tau^{z}\rangle$) with a valley
$g-$factor, $g_{\rm{v}}$.

We first set $\lambda=0$, $\bm{B}=0$, and consider the phase diagram as a
function of $\Delta_{\rm{V}}$ for fixed $\Delta^{\prime}_{\rm{V}}$ and
$J_{H}$. For $\Delta^{\prime}_{\rm{V}}<0$, we find a first-order phase
transition from the fully valley and spin polarized (VI) phase to an IVC phase
at a critical value of $\Delta_{\rm{V}}$ given by
$\Delta_{\rm{V},c}=-\Delta^{\prime}_{\rm{V}}-J_{H}/4$. The IVC is stable for
$\Delta_{\rm{V}}>\Delta_{\rm{V},c}$. For $\Delta^{\prime}_{\rm{V}}>0$ the
transition is of second order. In order to describe the experiment, where the
IVC to VI transition is observed to be first order, we focus on
$\Delta^{\prime}_{\rm{V}}\leq 0$ here.

Next, we consider the effect of $\lambda$, assumed for simplicity to be the
smallest energy scale in the problem. In the VI phase, $\lambda$ locks the
spin to be perpendicular to the plane, gaining an energy $-|\lambda|$. In the
IVC regime, the spin can be oriented either in the plane (which we dub the
IVCx phase) or out of the plane (the IVCz phase). The IVCx phase gains an
energy $\Delta E_{x}=-(\lambda/2)^{2}/(2J_{H})$ by canting the spins slightly
out of the plane. The IVCz phase gains an energy by slightly imbalancing the
valleys, $\langle\tau^{z}\rangle\neq 0$. For sufficiently small $\lambda$ such
that this valley imbalance is small, the energy gain is given by $\Delta
E_{z}=-\frac{1}{2}\chi_{v}(\lambda/2)^{2}$, where

$\chi_{v}=\frac{1}{2(|\Delta^{\prime}_{\rm{V}}|+\Delta_{\rm{V}}-\Delta_{1,c})}$
(S15)

is the valley susceptibility in the IVC regime (defined as
$\chi_{v}=\langle\tau^{z}\rangle/h_{v}$, where $h_{v}$ is a small valley
Zeeman field that couples to $\tau^{z}$, and $\langle\tau^{z}\rangle$ is the
resulting valley polarization).

Comparing the energies of the IVCx and IVCz phases, we find that the IVCz is
stable in a region
$0<\Delta_{\rm{V}}-\Delta_{\rm{V},c}<-\left|\Delta^{\prime}_{\rm{V}}\right|+\frac{J_{H}}{2}$
near the IVC to VI transition, assuming that
$J_{H}>2|\Delta^{\prime}_{\rm{V}}|$ (if $J_{H}<2|\Delta^{\prime}_{\rm{V}}|$,
IVCx is the ground state for all $\Delta_{\rm{V}}>\Delta_{1,c}$). Physically,
the intermediate IVCz is stable if the IVC to VI transition is sufficiently
weakly first order (i.e., $|\Delta^{\prime}_{\rm{V}}|$ is sufficiently small),
such that $\chi_{v}$ becomes sufficiently large upon approaching the
transition from the IVC side.

Fig. S1 shows the phase diagram of the model in the plane
$(\Delta_{\rm{V}},B_{\parallel})$ for $\lambda=0.2$ meV, $J_{H}=0.5$ meV, and
different values of $\Delta^{\prime}_{\rm{V}}$, obtained by minimizing the
energy (S14) over $\psi$. The color indicates $\langle\tau^{z}\rangle$. For
$\Delta^{\prime}_{\rm{V}}=0$, all the transitions are continuous, whereas for
$\Delta^{\prime}_{\rm{V}}<0$, the transitions from the VI to either IVCx or
IVCz is first order. The transition from IVCx to IVCz is always continuous.
Note how the intermediate IVCz phase expands as $|\Delta^{\prime}_{\rm{V}}|$
decreases; this is due to the enhxwancement of valley susceptibility
$\chi_{v}$ in the IVC side of the IVCz transition, which favors spin
polarization along $z$. Figure 3H in the main text shows phase diagrams as a
function of both in-plane and out-of-plane field and $\Delta_{\rm{V}}$
calculated from this model, for $\lambda=0.2$ meV,
$\Delta^{\prime}_{\rm{V}}=-0.025$ meV, $g_{v}=5g$, $J_{H}=0.5$ meV.

Figure S1: Phase diagram of the theoretical model as a function of
$\Delta_{\rm{V}}$ and $B_{\parallel}$ for different values of
$\Delta^{\prime}_{\rm{V}}$. The color shows $\langle\tau_{z}\rangle$. In the
VI phase, $|\langle\tau^{z}\rangle|=1$. In the IVCz phase,
$\langle\tau^{+}\rangle\neq 0$, $0<|\langle\tau^{z}\rangle|<1$ and $|\langle
s^{z}\rangle|\neq 0$. In the IVCx phase, $\langle\tau^{+}\rangle$ and $\langle
s^{+}\rangle$ are non-zero, whereas $\langle\tau^{z}\rangle=\langle
s^{z}\rangle=0$.

###  .5 Competition between IVC and valley imbalanced states

Both in the experimental data and in the Hartree-Fock calculations, the order
of phases as a function of increasing carrier density is dissimilar between
hole doping and electron doping. Namely, on the hole side the symmetric phase
at the lowest densities is followed first by an IVC, which in turn gives way
to a valley imbalanced phase. On the electron side, the symmetric phase is
instead first replaced by the valley imbalanced phase, and only then an IVC
sets in. As long as the Hund’s coupling is not too large, these differences
can be understood based on the different dispersions in the electron and hole
sides, as illustrated in Fig. S2. The IVC phase is favored by the difference
in kinetic energy between the two valleys at a given momentum (after folding
the two valleys on top of each other), which allows for a kinetic energy gain
by creating a momentum-dependent valley polarization. The VI phase is favored
by the exchange energy, which is maximized for electrons in the same valley,
due to the different form factors of the two valleys.

The dispersion near the bottom of the conduction band is significantly flatter
than the dispersion near the top of the valence band (Fig. S2). In particular,
the crossing between the dispersions of the two valleys, where the IVC gap
opens, happens closer in energy to the band edge on the hole side compared to
the electron side. These differences can explain why on the hole side, the IVC
phase quarter metal is favored at low hole density and the VI phase at higher
hole density (where the Berry curvature enclosed by the Fermi surface is
larger), whereas on the electron side the VI phase is favored at lower density
(where the bands are very flat), and the IVC phase is favored at higher
density (where the bands are more dispersive).

Figure S2:  Comparison of the non-interacting band structures of hole and
electron bands. The single-particle dispersions of valleys $K$ (blue) and
$K^{\prime}$ (red) for both (A) hole and (B) electron side, for $k_{y}=0$.
Indicated by a dashed gray line is the crossing region, which leads to the
formation of an IVC phase if the Fermi hlevel resides in the vicinity.