# Spin-orbit enabled unconventional Stoner magnetism

Yue Yu1, Tatsuya Shishidou1, Shuntaro Sumita2,3,4, Michael Weinert1, and
Daniel F. Agterberg1 1Department of Physics, University of
Wisconsin–Milwaukee, Milwaukee, Wisconsin 53201, USA  2Department of Basic
Science, The University of Tokyo, Meguro, Tokyo 153-8902, Japan  3Komaba
Institute for Science, The University of Tokyo, Meguro, Tokyo 153-8902, Japan
4Condensed Matter Theory Laboratory, RIKEN CPR, Wako, Saitama 351-0198, Japan

###### Abstract

The Stoner instability remains a cornerstone for understanding metallic
ferromagnets. This instability captures the interplay of Coulomb repulsion,
Pauli exclusion, and two-fold fermionic spin degeneracy. In materials with
spin-orbit coupling, fermionic spin is generalized to a two-fold degenerate
pseudospin. Here we identify a new fermionic pseudospin with a symmetry that
forbids it to couple to a Zeeman field. This spinless pseudospin exists in
five non-symmorphic space groups and appears at the Brillouin zone (BZ)
boundary. With Coulomb repulsion, this spinless pseudopsin gives rise to
Stoner instabilities into magnetic states that are not usual ferromagnets.
These spinless-pseudopsin ferromagnets break time-reversal symmetry, have a
vanishing magnetization, and are generally non-collinear. The non-collinear
magnetization enables a scalar spin chirality driven Berry curvature. In
addition, much like altermagnets, these pseudospin ferromagnets exhibit energy
band spin-splittings that vanish by symmetry along lines in the BZ, allowing
for drumhead surfaces states. We discuss candidate materials.

Spin orbit coupling (SOC) is an essential interaction in quantum materials. It
underlies realizations of the Kitaev model Jackeli and Khaliullin (2009),
topological insulators Schnyder et al. (2008), electric field control of
ferromagnets and antiferromagnets Baltz et al. (2018), and high critical
fields in Ising superconductors Wickramaratne and Mazin (2023). SOC typically
does not alter the symmetry properties of the underlying fermionic Bloch
spins. In particular, when time-reversal ($T$) and inversion ($I$) symmetries
are present, the two-fold Kramers’ degeneracy, or pseudopsin, that appears
when SOC is present typically behaves as usual spin-1/2 under rotations. An
immediate consequence is that a Zeeman magnetic field couples to this
pseudospin and can be used to control the quantum states discussed above. Here
we find a remarkable fermionic Bloch pseudospin that does not couple to a
Zeeman field, no matter which direction this field is applied.

As we will show, spinless pseudopsin is a consequence of two perpendicular
glide mirror symmetries. These two mirrors preserve Bloch momentum on lines in
momentum space. The magnetic response of Fermi surfaces near these momentum
lines are consequently dominated by spinless pseudospin. Hence, when we
examine Stoner instabilities Stoner (1938) of such Fermi surfaces, we find
that spinless-pseudopsin ferromagnets are energetically much more favorable
than usual ferromagnets.

Unlike usual ferromagnets, these pseudospin ferromagnets have no net
magnetization, are non-collinear in real space, and exhibit altermagnetic-like
spin-textures in the BZ Šmejkal et al. (2022). We find they generically occur
for the orthorhombic space groups Pbcm, Pbcn, Pbca, Pnma, and the cubic space
group Pa$\overline{3}$. All the pseudospin ferromagnetic states we find share
two common features. The first is that the non-collinear magnetism gives rise
to a non-zero scalar-spin chirality
$\chi_{ijk}=\vec{S}_{i}\cdot(\vec{S}_{j}\times\vec{S}_{k})$ Wen et al. (1989)
which implies a spin gauge flux of $\pm\pi/2$ for an electron traveling around
any closed triangle of three neighboring sites. The second is that symmetry
requires energy band spin-splittings that vanish along Weyl lines in the BZ.
These Weyl fermion lines are protected by mirror symmetry and allow drumhead
surface states to appear at mirror-invariant surfaces Chiu and Schnyder
(2014); Chan et al. (2016); Belopolski et al. (2019). Symmetries of pseudospin
Pauli Matrices: Our key finding is that in the five space groups listed above,
SOC implies a spinless pseudospin that does not couple to a Zeeman field.
Prior to discussing our general results, we highlight the symmetry differences
between usual pseudospin and the spinless pseudospin found here.

Class |  $C_{2v}$ Symmetry | $\tilde{\sigma}_{i}$
---|---|---
spinless pseudospin | $A_{1}$ | $\tilde{\sigma}_{x}=k_{z}k_{y}\sigma_{x}$
| $A_{1}$ | $\tilde{\sigma}_{y}=k_{z}k_{x}\sigma_{y}$
| $A_{1}$ | $\tilde{\sigma_{z}}=k_{x}k_{y}\sigma_{z}$
usual pseudospin | $B_{2}$ | $\tilde{\sigma}_{x}=\sigma_{x}$
| $B_{1}$ | $\tilde{\sigma}_{y}=\sigma_{y}$
| $A_{2}$ | $\tilde{\sigma}_{z}=\sigma_{z}$
Table 1: Symmetry distinction between spinless pseudospin Pauli matrices and
usual pseudospin Pauli matices. Since all spinless pseudospin components
belong to a different symmetry than usual spin, it cannot couple to a Zeeman
field.

Central to spinless pseudopsin are lines in the BZ that are invariant under
the point group $C_{2v}$, which consists of the identity $E$, two mirror
symmetries $M_{x}$ and $M_{y}$, and their product $C_{2z}$. At each ${\bf k}$
on this line, labeling the Kramers’ degenerate states as $|1\rangle$ and
$|2\rangle$, we can define pseudospin Pauli matrices
$\tilde{\sigma}_{x}=|1\rangle\langle 2|+|2\rangle\langle 1|$,
$\tilde{\sigma}_{y}=i|1\rangle\langle 2|-i|2\rangle\langle 1|$, and
$\tilde{\sigma}_{z}=|1\rangle\langle 1|-|2\rangle\langle 2|$. We find that on
these $C_{2v}$ invariant lines, in addition to usual pseudopsin, spinless
pseudospin exists. These two pseudospins differ in their symmetry
transformation properties under $C_{2v}$. As summarized in Table 1, for usual
pseudospin, the pseudospin Pauli matrices belong to the $A_{2}$, $B_{1}$, and
$B_{2}$ representations of the group $C_{2v}$. These are the same
representations as usual spin-1/2 Pauli matrices (here denoted by
$\sigma_{i}$). For spinless pseudospin, all pseudospin Pauli matrices
$\tilde{\sigma}_{i}$ belong to the $A_{1}$ representation of $C_{2v}$. It is
this symmetry property that ensures spinless pseudospin cannot couple to a
Zeeman field. To construct spinless pseudospin Pauli matrices in terms of
usual spin-1/2 Pauli matrices necessarily requires including some momentum
dependent functions. In Table 1, we carry out such a representative expansion
valid near the $\Gamma$ point.

SOC and origin of spinless pseudospin: A common feature of the five space
groups mentioned above is the large multiplicity of Wyckoff positions: there
are (at least) 4 equivalent sites in the primitive unit cell. This large sub-
lattice degree of freedom is a direct consequence of the non-symmorphic
symmetries and enables spinless pseudospin, which is not possible in single-
or double-site systems. In the Supplemental Materials (SM) we give general
symmetry arguments that demonstrate that spinless pseudopsin must emerge along
the following momentum lines: in Pbcm (SG 57) along $(k_{x},\pi,\pi)$; in Pbcn
(SG 60) along$(\pi,k_{y},\pi)$; in Pnma (SG 62) along$(\pi,\pi,k_{z})$; in
Pbca (SG 61) along $(k_{x},\pi,\pi),(\pi,k_{y},\pi)$, and $(\pi,\pi,k_{z})$;
and in Pa$\overline{3}$ (SG 205) along $(k_{x},\pi,\pi),(\pi,k_{y},\pi)$, and
$(\pi,\pi,k_{z})$. In the following we use a minimal microscopic model to
demonstrate how multifold sub-lattice degrees of freedom, non-symmorphic
symmetry, and SOC give rise to spinless pseudospin.

Here we focus our discussion on the space group Pbca with electronic
$s$-orbitals on the four 4a Wyckoff positions. Here these four Wyckoff
positions given by the labels
$\\{a_{1},a_{2},b_{1},b_{2}\\}=\\{(0,0,0),(0,1/2,1/2),(1/2,0,1/2),(1/2,1/2,0)\\}$.
We use the Pauli matrices $\tau_{i}$ and $\rho_{i}$ to differentiate between
the $(a,b)$ and $(1,2)$ sublattices respectively, and the Pauli matrices
$\sigma_{i}$ for the fermion spin. Without SOC, the tight-binding Hamiltonian
contains the three nearest neighbor inter-sublattice hopping terms:

$\begin{split}&H_{0}=\varepsilon_{\bf
k}+t_{ab}\cos\frac{k_{x}}{2}\cos\frac{k_{z}}{2}\tau_{x}+t_{12}\cos\frac{k_{y}}{2}\cos\frac{k_{z}}{2}\rho_{x}\\\
&+t_{a1b2}\cos\frac{k_{x}}{2}\cos\frac{k_{y}}{2}\tau_{x}\rho_{x}.\end{split}$
(1)

An important property of Eq. 1 is that along any linear edge of the BZ
boundary, for example $k_{x}=k_{y}=\pi$, there is an eight-fold degeneracy.
This general property is a consequence of the glide mirror symmetries
$\widetilde{M_{x}}=\\{M_{x}|1/2,1/2,0\\}$ and
$\widetilde{M_{y}}=\\{M_{y}|0,1/2,1/2\\}$ and the resultant eight-fold
degeneracy is known as a ‘double nodal line’ Takahashi et al. (2017); Li et
al. (2021). As we show below, these double nodal lines generally give rise to
spinless pseudospin.

We now include the following symmetry allowed nearest neighbor SOC terms

$\begin{split}&H_{SOC}=\lambda_{1}\sin\frac{k_{x}}{2}\sin\frac{k_{y}}{2}\tau_{y}\rho_{x}\sigma_{x}\\\
&+\lambda_{2}\sin\frac{k_{y}}{2}\sin\frac{k_{z}}{2}\rho_{y}\sigma_{y}+\lambda_{3}\sin\frac{k_{x}}{2}\sin\frac{k_{z}}{2}\tau_{y}\rho_{z}\sigma_{z}\end{split}$
(2)

On the $k_{x}=k_{y}=\pi$ line, Eq.2 is non-vanishing and consists of two
commuting terms ($\tau_{y}\rho_{x}\sigma_{x}$, $\rho_{y}\sigma_{y}$) and their
product. They break the eight-fold degeneracy into four two-fold Kramers’
degenerate pairs. Interestingly, these two commuting terms are proportional to
the mirror operators $\widetilde{M}_{x}$ and $\widetilde{M}_{y}$. Each of
these Kramers’ pairs is hence labelled by eigenvalues of $\widetilde{M}_{x}$
and $\widetilde{M}_{y}$. Both of the Kramers’ degenerate states have the same
mirror eigenvalues, hence Pauli matrices made from these two states must be
invariant under both mirror operations, implying spinless pseudospin. Indeed
any operator that changes either mirror eigenvalue must be inter-band. This
can be readily verified for the spin operators $\sigma_{i}$, which anticommute
with at least one SOC term.

Figure 1: Intra-band $g$-factors in UCoGe calculated along the $k$ line
$(\pi+k,\pi+k,0.2\pi)$. These $g$-factors quantify the energy splitting
between Kramers’ doublets for fields applied along the $x$, $y$, and $z$
directions. The value $g=1$ corresponds to usual spin-1/2. Inset: Fermi
surface of UCoGe. The $g$ calculations are done along the arrow and for the
bands producing the inner (pale red) Fermi surface.

Vanishing $g$-factors: application to UCoGe: The response of pseudospin to a
Zeeman field is generally described by the Hamiltonian
$H_{Z}=\mu_{B}\sum_{i,j}g_{ij}H_{i}\tilde{\sigma}_{j}$. For spinless
pseudospin, symmetry requires that $g_{i,j}=0$. This result has implications
for pseudospin-singlet superconductors occurring on a Fermi surface
surrounding a double nodal line. In particular, since spinless pseudopsin
dominates the spin physics on such a Fermi surface, there is nearly no Pauli
suppression of superconductivity for the field applied in any direction. This
generalizes the well known lack of Pauli suppression for in-plane fields
observed in Ising superconductors Suh et al. (2023); Wickramaratne and Mazin
(2023).

One relevant material is UCoGe. UCoGe has space group Pnma and exhibits
critical fields that exceed the Pauli limit for all field directions Aoki et
al. (2019). DFT shows that UCoGe exhibits quasi-2D cylindrical Fermi surfaces
near the R-S line consistent with quantum oscillation measurements Bastien et
al. (2016). Along the R-S line we predict spinless pseudospin. Figure 1 shows
the DFT calculated $g$-factors near the R-S line and reveals that they all
vanish on this line. Furthermore, they are small near this line, consistent
with the large observed Pauli fields. UCoGe is often argued to be a
pseudospin-triplet superconductor based on the large observed critical fields.
Our results allow the possibility that this a pseudospin-singlet
superconductor and further necessitates a more careful analysis of spin-
fluctuations in such superconductors.

Stoner spinless-pseudospin ferromagnets: An important finding of our work is
that spinless pseudospin leads to Stoner instabilities into novel magnetic
states. To understand this, it is useful to revisit the Stoner instability of
a single band into a ferromagnet. Central to this is the local Coulomb
repulsion on each lattice site $i$, $Un_{i,\uparrow}n_{i,\downarrow}$ (where
$n_{i,s}$ gives the number of electrons with spin $s$). Treating this
interaction within mean field theory leads to the usual Stoner criterion for a
ferromagnetic instability, $UN(0)>1$ where $N(0)$ is the density of states at
the chemical potential. The key observation for more general pseudospin is
that the local Coulomb repulsion is statistical since it originates from the
Pauli exclusion of the two degenerate spin states that comprise the band.
Consequently, the same argument implies a Stoner instability into a pseudospin
ferromagnet. When the pseudopsin symmetry differs from that of usual spin-1/2,
the resulting Stoner driven magnetic state must have a different symmetry than
a usual ferromagnet.

To verify this argument and to understand the resultant spinless-pseudospin
ferromagnet states, we have treated an on-site Coulomb repulsive interaction
within the random phase approximation (RPA) for the three space groups Pbca,
Pnma, and Pa$\overline{3}$. For space groups Pbca and Pa$\overline{3}$ this
was carried out for single Kramers’ degenerate band with a spherical Fermi
surface around the R-point $(\pi,\pi,\pi)$. For space group Pnma, this was
carried out for a cylindrical Fermi surface around the R-S line. We included
all possible on-site charge and magnetic orders which, together with the four
sub-lattice degrees of freedom, allows for sixteen instability channels. The
details of this RPA analysis are given in the SM. Among these instability
channels, we generically find only three channels for which the critical value
of $U$ satisfies $U_{c}N(0)\approx O(1)$, similar to the usual Stoner
criterion. All other instability channels, including usual ferromagnetic
instabilities, yield unreasonably large $U_{c}$, with $U_{c}N(0)\gg 1$.

Figure 2: Non-collinear structure of a spinless-pseudospin ferromagnetic for
space group Pbca Momma and Izumi (2011). This structure was calculated for a
spherical Fermi surface around the R-point. A similar non-collinear structure
was found for space group Pnma for a cylindrical Fermi surface surrounding the
R-S line.

The resultant leading magnetic instability for space group Pbca is shown in
Fig. 2. Remarkably, very unlike a usual ferromagnetic Stoner instability, this
instability is to a non-collinear magnetic state. The moments are given by
$\vec{S}_{a1}=(m_{x},m_{y},m_{z})$, $\vec{S}_{a2}=(-m_{x},m_{y},-m_{z})$,
$\vec{S}_{b1}=(-m_{x},-m_{y},m_{z})$, and
$\vec{S}_{b2}=(m_{x},-m_{y},-m_{z})$. Similar non-collinear magnetic
structures are found for the leading three instabilities for all the space
groups considered here. The resultant magnetic space groups are all magnetic
type 1 Federov groups for which all space group symmetries remain and only $T$
symmetry is broken.

Given that these different space groups all generate similar magnetic states
from the same mechanism, it is reasonable to ask if they share any common
physical properties. In the following we discuss two such properties: $\pi/2$
spin gauge flux generated from the non-collinear magnetic order and symmetry
required Weyl line nodes.

Scalar spin chirality: Non-collinear magnetic structures can give rise to
scalar spin chirality
$\chi_{ijk}=\vec{S}_{i}\cdot(\vec{S}_{j}\times\vec{S}_{k})$ which is closely
linked to orbital magnetism and Berry curvature Wen et al. (1989); Shindou and
Nagaosa (2001); Taguchi et al. (2001); Machida et al. (2010). For all the
magnetic groups we find here, $\chi_{ijk}=\pm 4m_{x}m_{y}m_{z}$ for any
elementary triangle of nearest neighbor spins. This non-collinear magnetic
structure implies spin-gauge flux for a spin-1/2 electron that travels each
such triangle. We find this spin gauge flux is $\pm\pi/2$ for all the
spinless-pseudopsin ferromagnetic states. Such spin-gauge fluxes are known to
generate Berry curvature and an anomalous Hall effect (AHE). However, due to
the presence of mirror symmetries in spinless-fermion ferromagnets, the total
AHE will vanish. However, the application of a shear strain, which breaks a
mirror symmetry, will allow both an orbital magnetization and an AHE to
appear. This has been explored in distorted fcc lattices for a magnetic state
that is identical to that we have found for space group Pa$\overline{3}$
Shindou and Nagaosa (2001). Our findings imply that the results of Ref.
Shindou and Nagaosa, 2001 can be readily generalized to all the orthorhombic
space groups we consider here. We leave a detailed discussion to a future
publication.

Weyl lines: Spinless-pseudopsin ferromagnets generically have vanishing
energy-band spin-splitting along lines in momentum space. These Weyl lines
imply drumhead surface states on surfaces oriented along the mirror planes
Chiu and Schnyder (2014); Chan et al. (2016); Belopolski et al. (2019). Near
the $\Gamma$ point, symmetry requires that the Hamiltonian for spinless-
pseudopsin ferromagnets is (a related Hamiltonian is discussed in Ref.
Fernandes et al., 2023)

$H=\epsilon_{0}({\bm{k}})+\alpha_{x}k_{z}k_{y}\sigma_{x}+\alpha_{y}k_{z}k_{x}\sigma_{y}+\alpha_{z}k_{x}k_{y}\sigma_{z}$
(3)

where $\alpha_{i}$ are constants. This yields Weyl lines along the $x$, $y$,
and $z$ axes that are each topologically protected by two of the three
orthogonal mirror symmetries. On any mirror-symmetric loop around a Weyl line,
the Berry phase is $\pi$ (here the relevant mirror symmetry is orthogonal to
the Weyl line). Such a loop can be extended to the BZ edges, forming two
counter-propagating one-dimensional (1d) integrals that each yield a Zak phase
Chan et al. (2016). For these 1d integrals, the same mirror symmetry serves as
an effective inversion symmetry that quantizes the Zak phase to be zero or
$\pi$. Provided there are no additional Weyl lines that intersect these 1d
contours, the sum of the two Zak phases must be $\pi$. This implies one of the
two 1d contours is topological and yields drumhead surface states on crystal
surfaces that are mirror invariant. Key to the existence of these states is
the position of possible additional Weyl lines away from the $\Gamma$ point.
Such additional Weyl lines can imply non-topological, or zero, Zak phases
along the 1d contours across the BZ. Generically, using Elcoro et al. (2021),
we find that such Weyl lines exist on the BZ boundary for all the spinless-
pseudopsin ferromagnets discussed here. The positions of all Weyl lines are
shown in Fig. 3. This implies that drumhead surface states will appear on all
$x$, $y$, and $z$ surfaces for spinless-pseudospin ferromagnets in space
groups Pbca, Pbcn, Pnma, and Pa$\overline{3}$, and on the $y$ and $z$ surfaces
for space group Pbcm. In the SM we reveal these surface states for space group
Pa$\overline{3}$. In addition to the Weyl lines discussed above, the $A_{1}$
spinless-pseudospin ferromagnet in space group Pa$\overline{3}$ also hosts
three-fold degenerate fermions at the R point Cano et al. (2019).

Figure 3: BZs for (a) primitive cubic and (b) primitive orthorhombic space
groups. (c–f) Weyl fermion lines (in red) for the spinless-pseudospin
ferromagnetic space groups found here.

Materials realizations: The Magndata database Gallego et al. (2016) contains
37 entries that are candidates for the Stoner mechanism we discuss above. Two
material classes are of particular interest. One is a group of 11 perovskite
materials with space group Pnma. These have the form RTO3, where R is a rare
earth and T is a transition metal. DFT shows that many of these materials are
metallic in the normal state (see for example Hasan et al. (2022)). The second
is NiS2-xSex with space group Pa$\overline{3}$. This material is metallic
under doping Xu et al. (2014) and under pressure Friedemann et al. (2016). DFT
calculations reveal that NiS2 shows Fermi surfaces surrounding the R-point and
additional Fermi surfaces surrounding the double nodal lines Reiss et al.
(2022), making this materials class a strong candidate for the mechanism we
propose.

Discussion and conclusions: Our results emphasize Stoner instabilities of
Fermi surfaces near double nodal lines that support spinless pseudopsin. We
expect spinless-pseudospin ferromagnets occur more generally for space groups
Pbca, Pnma, and Pa$\overline{3}$. For these space groups, the pseudospin is
not unusual just on the double nodal lines but also on all the BZ zone faces
Suh et al. (2023). In particular, on the BZ faces, the pseudospin only couples
to a magnetic field perpendicular to the face Cavanagh et al. (2022); Suh et
al. (2023). This is energetically important for the formation of usual
ferromagnetism: for any ferromagnetic state with a moment oriented normal to
the a BZ face, the other two BZ zone faces must exhibit vanishing spin
splittings (see also Ref. Wilde et al., 2021). Since these Weyl fermion plane
nodes do not allow an energetic gain from the usual spin-splitting of
ferromagnetism, they disfavor a Stoner instability into a usual ferromagnetic
state. As Figure 3 shows, the spinless-pseudopsin ferromagnetic states
discussed here do not have such nodal planes. This implies larger spin-
splittings on the BZ faces for the spinless-pseudospin ferromagnetic states
compared to usual ferromagnets. Hence, for materials with Fermi surfaces near
all these BZ faces, the Stoner instability will likely be into spinless-
pseudospin ferromagnetic states.

In conclusion, we have identified a Bloch pseudospin that does not couple to a
Zeeman field and exists on double nodal lines in space groups Pbcm, Pbcn,
Pbca, Pnma, and Pa$\overline{3}$. We have shown that this spinless pseudospin
enables field-robust superconductivity and Stoner instabilities in non-
collinear magnetic states that have no net magnetization. Finally, we have
shown that these spinless-pseudopsin ferromagnetic states generically have
scalar spin chirality with $\pm\pi/2$ spin gauge flux and Weyl fermion lines
that give rise to drumhead surface states.

We thank Rafael Fernandes, Tomas Jungwirth, and Lian Li for useful
discussions. D.F.A., T.S., M.W., and Y.Y were supported by the National
Science Foundation Grant No. DMREF 2323857. S.S. was supported by JSPS KAKENHI
Grants No. JP23K03333 and No. JP23K13056, and JST CREST Grant No. JPMJCR19T2.

##  Appendix A Symmetry arguments for spinless pseudospin

Here, we will first review the symmetry arguments for the double nodal line in
the absence of SOC. We then include SOC to reveal how the eight-fold
degeneracy along these double nodal lines are split into four Kramers’
doublets.

The eight-fold band degeneracy requires four-fold orbital degeneracy. Consider
the symmetry elements that preserve a momentum point on the nodal line
$(\pi,\pi,k_{z})$. They are generated by
$\\{E,\widetilde{M_{x}},\widetilde{M_{y}},IT\\}$, where $\widetilde{M_{x}}$ is
the usual mirror symmetry $x\rightarrow-x$ followed by a translation
${\bf{t_{x}}}=(t_{x}^{x},t_{x}^{y},t_{x}^{z})$. The origin is chosen to be the
inversion center. We have
$\widetilde{M_{x}}^{2}=\\{0,2t^{y}_{x},2t^{z}_{x}\\}$,
$\widetilde{M_{x}}I=\\{2t^{x}_{x},2t^{y}_{x},2t^{z}_{x}\\}I\widetilde{M_{x}}$,
and $I^{2}=T^{2}=1$. Here $\\{...\\}$ denotes translation. The spin is ignored
in the absence of SOC. We show that the four states (below labelled as
$1\rightarrow 4$):

$|\psi\rangle,\ \widetilde{M_{x}}IT|\psi\rangle,\
\widetilde{M_{x}}\widetilde{M_{y}}|\psi\rangle,\
\widetilde{M_{y}}IT|\psi\rangle$ (4)

are degenerate orthogonal eigenstates. They follows from two two-fold Kramers’
degeneracies. The first Kramers’ degeneracy arises from $\widetilde{M_{x}}IT$.
This satisfies
$(\widetilde{M_{x}}IT)^{2}=\\{2t^{x}_{x},0,0\\}=\exp(2ik_{x}t_{x}^{x})=-1$,
(here we used the fact that the $T$ commutes with other symmetry operations).
This symmetry makes states (1,2) as well as (3,4) orthogonal to each other.
The second Kramers’ degeneracy arises from $\widetilde{M_{y}}IT$, which
satisfies $(\widetilde{M_{y}}IT)^{2}=\exp(2ik_{y}t_{y}^{y})=-1$. This makes
states (1,4) as well as (2,3) orthogonal to each other. For these two Kramers’
degeneracies to be independent, states (1,3): $|\psi\rangle$ and
$\widetilde{M_{x}}\widetilde{M_{y}}|\psi\rangle$ need to be orthogonal. This
requires $\\{\widetilde{M_{x}},\widetilde{M_{y}}\\}=0$, such that
$|\psi\rangle$ can be chosen as an eigenstate of $\widetilde{M_{x}}$, while
$\widetilde{M_{x}}\widetilde{M_{y}}|\psi\rangle$ is an orthogonal eigenstate
with the opposite $\widetilde{M_{x}}$ eigenvalue. In summary, without SOC, the
eight-fold band degeneracy requires:
$\exp(2ik_{x}t_{x}^{x})=\exp(2ik_{y}t_{y}^{y})=-1$ and
$\\{\widetilde{M_{x}},\widetilde{M_{y}}\\}=0$. The former condition gives
$t_{x}^{x}=t_{y}^{y}=1/2$ and the latter condition restricts
$(t_{x}^{y},t_{y}^{x})=(0,1/2)$ or $(1/2,0)$.

We now include the SOC. The time reversal operation now satisfies $T^{2}=-1$.
The bands have the usual Kramers’ degeneracy from $IT$. The mirror operators
$\widetilde{M_{x,y}}$ now include the spin part $i\sigma_{x,y}$. This leads to
$[\widetilde{M_{x}},\widetilde{M_{y}}]=0$. There is now no additional
degeneracy arising from other symmetry operators involving $T$, since
$\widetilde{M_{x}}$ and $\widetilde{M_{y}}$ commute with the Hamiltonian and
with each other, so define four non-degenerate orthogonal eigenstates. Hence
the eight-fold degeneracy splits into four Kramers’ pairs.

We now prove that each Kramers’ pair shares the same $\widetilde{M_{x}}$
eigenvalues. The proof for $\widetilde{M_{y}}$ is similar. Since
$\widetilde{M_{x}}^{2}=(i\sigma_{x})^{2}\\{0,2t^{y}_{x},2t^{z}_{x}\\}=-\exp(2ik_{y}t_{x}^{y}+2ik_{z}t_{x}^{z})$,
the eigenvalue of $\widetilde{M_{x}}$ is $\pm
i\exp(ik_{y}t_{x}^{y}+ik_{z}t_{x}^{z})$. Let $|\psi\rangle$ be an eigenstate
of $\widetilde{M_{x}}$, then consider its Kramers’ partner $TI|\psi\rangle$:

$\begin{split}&\widetilde{M_{x}}|\psi\rangle=ie^{ik_{y}t_{x}^{y}+ik_{z}t_{x}^{z}}|\psi\rangle\\\
\Rightarrow&\widetilde{M_{x}}TI|\psi\rangle=\\{2t_{x}^{x},2t_{x}^{y},2t_{x}^{z}\\}TI\widetilde{M_{x}}|\psi\rangle=e^{2ik_{x}t_{x}^{x}+2ik_{y}t_{x}^{y}+2ik_{z}t_{x}^{z}}TI(ie^{ik_{y}t_{x}^{y}+ik_{z}t_{x}^{z}}|\psi\rangle)\\\
=&(-e^{2ik_{x}t_{x}^{x}})ie^{ik_{y}t_{x}^{y}+ik_{z}t_{x}^{z}}TI|\psi\rangle=ie^{ik_{y}t_{x}^{y}+ik_{z}t_{x}^{z}}TI|\psi\rangle\end{split}$
(5)

Hence $TI|\psi\rangle$ is another eigenstate with the same $\widetilde{M_{x}}$
eigenvalue as $|\psi\rangle$. Since the two Kramers’ partners share the same
mirror eigenvalues, Pauli matrices made from these Kramers’ partners must be
invariant under the two mirror operations $\widetilde{M_{x}}$ and
$\widetilde{M_{y}}$ Cavanagh et al. (2022); Suh et al. (2023). This implies
that all the Kramers’ doublets form spinless pseudospin.

##  Appendix B Tight binding models used for space groups Pbca, Pnma, and
Pa$\overline{3}$

For all these space groups we considered $s$-orbitals on Wyckoff position 4a.
For space groups Pbca and Pa$\overline{3}$ these are given by
$(a_{1},a_{2},b_{1},b_{2})=\\{(0,0,0),(0,1/2,1/2),(1/2,0,1/2),(1/2,1/2,0)\\}$.
We use the Pauli matrices $\tau_{i}$ and $\rho_{i}$ to differentiate between
the $(a,b)$ and $(1,2)$ sublattices respectively, and the Pauli matrices
$\sigma_{i}$ for the spin of these fermions. The tight-binding Hamiltonian
that includes all nearest neighbor interactions is:

$\begin{split}H&=2t_{x}\cos k_{x}+2t_{y}\cos k_{y}+2t_{z}\cos k_{z}-\mu\\\
&+t_{ab}\cos\frac{k_{x}}{2}\cos\frac{k_{z}}{2}\tau_{x}+t_{12}\cos\frac{k_{y}}{2}\cos\frac{k_{z}}{2}\rho_{x}+t_{a1b2}\cos\frac{k_{x}}{2}\cos\frac{k_{y}}{2}\tau_{x}\rho_{x}+H_{SOC1}+H_{SOC2}\\\
H_{SOC1}&=\lambda_{1}\sin\frac{k_{x}}{2}\sin\frac{k_{y}}{2}\tau_{y}\rho_{x}\sigma_{x}+\lambda_{2}\sin\frac{k_{y}}{2}\sin\frac{k_{z}}{2}\rho_{y}\sigma_{y}+\lambda_{3}\sin\frac{k_{x}}{2}\sin\frac{k_{z}}{2}\tau_{y}\rho_{z}\sigma_{z}\\\
H_{SOC2}&=a_{1}\cos\frac{k_{y}}{2}\cos\frac{k_{z}}{2}\rho_{y}\sigma_{z}+a_{2}\cos\frac{k_{y}}{2}\cos\frac{k_{z}}{2}\tau_{z}\rho_{y}\sigma_{x}+a_{3}\cos\frac{k_{x}}{2}\cos\frac{k_{z}}{2}\tau_{y}\sigma_{y}\\\
&+a_{4}\cos\frac{k_{x}}{2}\cos\frac{k_{z}}{2}\tau_{y}\rho_{z}\sigma_{x}+a_{5}\cos\frac{k_{x}}{2}\cos\frac{k_{y}}{2}\tau_{x}\rho_{y}\sigma_{z}+a_{6}\cos\frac{k_{x}}{2}\cos\frac{k_{y}}{2}\tau_{y}\rho_{x}\sigma_{y}\end{split}$
(6)

For space group Pbca, the coefficients are generically different. For space
group Pa$\overline{3}$, the following constraints are imposed by symmetry:

$t_{x}=t_{y}=t_{z}=:t,\quad
t_{ab}=t_{12}=t_{a1b2}=:t^{\prime},\quad\lambda_{1}=\lambda_{2}=\lambda_{3}=:\lambda,\quad
a_{1}=a_{4}=a_{6}=:a,\quad a_{2}=a_{3}=a_{5}=:a^{\prime}.$ (7)

For space group Pnma, we denote the 4a Wyckoff positions as
$(a_{1},a_{2},b_{1},b_{2})=\\{(0,0,0),(1/2,0,1/2),(0,1/2,0),(1/2,1/2,1/2)\\}$.
The tight-binding Hamiltonian is:

$\begin{split}H&=2t_{x}\cos k_{x}+2t_{y}\cos k_{y}+2t_{z}\cos k_{z}-\mu\\\
&+t_{ab}\cos\frac{k_{y}}{2}\tau_{x}+t_{12}\cos\frac{k_{x}}{2}\cos\frac{k_{z}}{2}\rho_{x}+t_{a1b2}\cos\frac{k_{x}}{2}\cos\frac{k_{y}}{2}\cos\frac{k_{z}}{2}\tau_{x}\rho_{x}+H_{SOC1}+H_{SOC2}\\\
H_{SOC1}&=\lambda_{1}\sin\frac{k_{x}}{2}\sin\frac{k_{y}}{2}\cos\frac{k_{z}}{2}\tau_{x}\rho_{y}\sigma_{x}+\lambda_{2}\sin\frac{k_{y}}{2}\sin{k_{z}}\tau_{y}\sigma_{y}+\lambda_{3}\sin\frac{k_{x}}{2}\sin\frac{k_{z}}{2}\tau_{z}\rho_{y}\sigma_{z}\\\
H_{SOC2}&=a_{1}\cos\frac{k_{x}}{2}\cos\frac{k_{z}}{2}\tau_{0}\rho_{y}\sigma_{y}+a_{2}\cos\frac{k_{x}}{2}\cos\frac{k_{z}}{2}\tau_{z}\rho_{y}\sigma_{x}+a_{3}\cos\frac{k_{y}}{2}\tau_{y}\rho_{0}\sigma_{z}\\\
&+a_{4}\cos\frac{k_{y}}{2}\tau_{y}\rho_{z}\sigma_{x}+a_{5}\cos\frac{k_{x}}{2}\cos\frac{k_{y}}{2}\cos\frac{k_{z}}{2}\tau_{x}\rho_{y}\sigma_{y}+a_{6}\cos\frac{k_{x}}{2}\cos\frac{k_{y}}{2}\cos\frac{k_{z}}{2}\tau_{y}\rho_{x}\sigma_{z}\\\
&+a_{7}\cos\frac{k_{y}}{2}\sin\frac{k_{x}}{2}\sin\frac{k_{z}}{2}\tau_{y}\rho_{x}\sigma_{x}+a_{8}\cos\frac{k_{x}}{2}\sin\frac{k_{y}}{2}\sin\frac{k_{z}}{2}\tau_{x}\rho_{y}\sigma_{z}+a_{9}\cos\frac{k_{x}}{2}\sin\frac{k_{y}}{2}\sin\frac{k_{z}}{2}\tau_{y}\rho_{x}\sigma_{y}\end{split}$
(8)

##  Appendix C Random phase approximation

In this section, we provide the details for the RPA calculations. We consider
on-site repulsive Hubbard interactions on each site and all on-site charge and
magnetic vertices. On each site $i=1...4$, we include the four vertices
$V_{i\uparrow\uparrow}$, $V_{i\downarrow\downarrow}$,
$V_{i\uparrow\downarrow}$ and $V_{i\downarrow\uparrow}$. We treat this within
an RPA approach. Importantly, we formulate this within the original
sublattice/spin basis, a basis in which our electronic Green’s functions are
not diagonal. We then project the results onto the band basis for a Kramers’
doublet that crosses the chemical potential.

The one-loop, or RPA, corrections for determining these vertices can be
written in matrix form:

$\begin{split}\bf V=MV+V_{0}\end{split}$ (9)

Here, ${\bf V}$ is a vector containing the 16 vertices. ${\bf V_{0}}$ is the
vector containing the bare vertices and $\bf M$ is a $16\times 16$ matrix that
quantifies all the one-loop corrections. The dominant instability corresponds
to the eigenvector of $\bf M$ with the largest eigenvalue.

The RPA calculation involves evaluating the two types of one-loop diagrams
shown in Fig. 4. The left diagram is relevant for output vertices with
opposite spins, while the right diagram is for output vertices with the same
spins. As an example, we evaluate two elements of $\bf M$ for these two
Feynman diagrams. For a particular point on the Fermi surface, we have two
pseudospin eigenvectors $u^{I}$ and $u^{II}$. The two elements are

$\begin{split}M_{(j\uparrow,\uparrow),(i\uparrow,\downarrow)}&={UN(0)}\left\langle\sum_{a,b=I}^{II}u^{a}_{j\uparrow}u^{a*}_{i\uparrow}u^{b}_{i\downarrow}u^{b*}_{j\uparrow}\right\rangle\\\
M_{(j\uparrow,\uparrow),(i\uparrow,\uparrow)}&=-{UN(0)}\left\langle\sum_{a,b=I}^{II}u^{a}_{j\uparrow}u^{a*}_{i\downarrow}u^{b}_{i\downarrow}u^{b*}_{j\uparrow}\right\rangle\\\
\end{split}$ (10)

Here, $u^{a}_{j\uparrow}$ denotes the $j\uparrow$ component in the eigenvector
$u^{a}$. $U$ is the Hubbard interaction strength, and $N(0)$ is the density of
states at the Fermi level. $\langle...\rangle$ denotes the averaging over the
Fermi surface.

Figure 4: One-loop Feynman diagrams used here. (Left) For output vertex
$\uparrow\downarrow$ or $\downarrow\uparrow$. (Right) For output vertex
$\uparrow\uparrow$ or $\downarrow\downarrow$. i and j denote sublattice sites.

For a system with one sublattice and no SOC, the largest eigenvector of $\bf
M$ is the usual Stoner ferromagnetic state with eigenvalue equal to $UN(0)$.
The critical interaction strength, or the zero-temperature Stoners’ criteria
is then $U_{c}^{-1}=N(0)$. For a system with four sublattices and no SOC, the
largest eigenvector of $\bf M$ is still a ferromagnetic state with eigenvalue
equal to $UN(0)/4$. There is a factor of $1/4$ as $N(0)$ is defined as the
total DOS at the Fermi surface. The critical interaction strength is then
$U_{c0}^{-1}\equiv N(0)/4$.

In the calculation for space group Pbca, we take $t_{x}=1$, $t_{y}=1.05$,
$t_{z}=1.1$, $\mu=-6$, $t_{ab}=0.1$, $t_{12}=0.095$, $t_{a1b2}=0.105$
$\lambda_{i}=0.1$ and $a_{i}=0.1$. We consider the elliptical Fermi surface
from the band with the lowest energy. The dominant instability is the
altermagnetic-like state
$0.4839\rho_{z}\sigma_{z}+0.4663\tau_{z}\rho_{z}\sigma_{x}+0.0498\tau_{z}\sigma_{y}$
with critical interaction strength $U_{c}^{-1}=0.9964U_{c0}^{-1}$. The next
two leading instabilities are also altermagnetic-like states with similar
critical strengths: $U_{c}^{-1}=0.9961U_{c0}^{-1}$ and
$U_{c}^{-1}=0.9951U_{c0}^{-1}$. The fourth leading instability is a canted
ferromagnetic state
$0.1823\sigma_{z}+0.0449\tau_{z}\sigma_{x}+0.7728\tau_{z}\rho_{z}\sigma_{y}$.
This state requires a much larger critical strength
$U_{c}^{-1}=0.0024U_{c0}^{-1}$.

In the calculation for space group Pnma, we take $t_{x}=1=t_{y}=1$,
$t_{z}=-0.1$, $\mu=-3.9$, and $t_{ab}=t_{12}=t_{a1b2}=\lambda_{i}=a_{i}=0.01$.
We consider the elliptical Fermi surface from the band with the lowest energy.
The dominant instability is into the altermagnetic-like state
$0.1404\tau_{z}\sigma_{z}+0.185\tau_{z}\rho_{z}\sigma_{x}+0.6747\rho_{z}\sigma_{y}$
with critical interaction strength $U_{c}^{-1}=0.95U_{c0}^{-1}$. The next two
leading instabilities are also altermagnetic-like states with similar critical
strengths: $U_{c}^{-1}=0.82U_{c0}^{-1}$ and $U_{c}^{-1}=0.78U_{c0}^{-1}$. The
fourth leading instability is into a ferromagnetic state
$0.0139\rho_{z}\sigma_{z}+0.7353\sigma_{x}+0.2508\tau_{z}\sigma_{y}$. This
requires a much larger critical strength: $U_{c}^{-1}=0.18U_{c0}^{-1}$.

##  Appendix D Drumhead surface states

In this section, we demonstrate the drumhead surface states for space group
Pa$\overline{3}$. For this purpose, we use the tight-binding model in Eq. (6)
under the conditions (7). Furthermore, additional SOC terms allowed in the
space group symmetry are considered here:

$\displaystyle H_{\text{SOC4}}$ $\displaystyle=\lambda(\delta_{1}\cos
k_{x}+\delta_{2}\cos k_{y}+\delta_{3}\cos
k_{z})\sin\frac{k_{x}}{2}\sin\frac{k_{y}}{2}\tau_{y}\rho_{x}\sigma_{x}$
$\displaystyle+\lambda(\delta_{3}\cos k_{x}+\delta_{1}\cos
k_{y}+\delta_{2}\cos
k_{z})\sin\frac{k_{y}}{2}\sin\frac{k_{z}}{2}\rho_{y}\sigma_{y}$
$\displaystyle+\lambda(\delta_{2}\cos k_{x}+\delta_{3}\cos
k_{y}+\delta_{1}\cos
k_{z})\sin\frac{k_{x}}{2}\sin\frac{k_{z}}{2}\tau_{y}\rho_{z}\sigma_{z},$ (11)

which lift accidental degeneracies of the energy bands and enable us to
clearly see the existence of the surface states.

To discuss the altermagnetic-like state induced by the Stoner instability for
spinless pseudopsin, we treat the magnetic order as a molecular field:
$h(\rho_{z}\sigma_{z}+\tau_{z}\rho_{z}\sigma_{x}+\tau_{z}\sigma_{y})$. Figure
5 shows energy spectra on the slice with $k_{x}=\pi/2$. The blue solid lines
are obtained by considering periodic boundary conditions (PBCs) for all
directions, whereas the red dashed lines represent the spectra under an open
boundary condition (OBC) along the $z$ axis and PBCs along the other $x$ and
$y$ directions. There are several surface bands in the bulk gap, which
originate from the topology of the Weyl lines in Fig. 3(c) of the main text.
For example, the drumhead surface states emerging in the $-\pi<k_{y}<0$ region
of the upper bands (the inset of Fig. 5) are characterized by the Zak phase,

$\gamma_{n}(k_{x},k_{y})=\frac{1}{i}\int_{-\pi}^{\pi}\mathrm{d}k_{z}\langle
u_{n}(\bm{k})|\partial_{k_{z}}|u_{n}(\bm{k})\rangle\pmod{2\pi},$ (12)

which is quantized into a $Z_{2}$ variable by the mirror symmetry
$\widetilde{M}_{z}$. Here $n$ represents a band index. We calculate the Zak
phase $\gamma_{n}(k_{x},k_{y})$ by choosing $n$ as the top band, and confirm
that it takes a nontrivial value $\pi$ for $-\pi<k_{y}<0$, while it is zero
for $0<k_{y}<\pi$. These results correspond to the presence or absence of
drumhead surface states.

Figure 5: Weyl-line node driven topological surface states. Parameters are set
to
$(t,t^{\prime},\lambda,a,a^{\prime},\delta_{1},\delta_{2},\delta_{3},h,\mu)=(-0.2,0.2,1,1,-1,0.2,0.5,0.5,1,2.5)$.
The blue solid and red dashed lines represent the spectrum under the PBC and
OBC along the $z$ axis, respectively. The inset indicates the enlarged energy
spectrum in the orange dash–dotted box.

## References

  * Jackeli and Khaliullin (2009) G. Jackeli and G. Khaliullin, Phys. Rev. Lett. 102, 017205 (2009), URL https://link.aps.org/doi/10.1103/PhysRevLett.102.017205.
  * Schnyder et al. (2008) A. P. Schnyder, S. Ryu, A. Furusaki, and A. W. W. Ludwig, Phys. Rev. B 78, 195125 (2008), URL https://link.aps.org/doi/10.1103/PhysRevB.78.195125.
  * Baltz et al. (2018) V. Baltz, A. Manchon, M. Tsoi, T. Moriyama, T. Ono, and Y. Tserkovnyak, Rev. Mod. Phys. 90, 015005 (2018), URL https://link.aps.org/doi/10.1103/RevModPhys.90.015005.
  * Wickramaratne and Mazin (2023) D. Wickramaratne and I. I. Mazin, Applied Physics Letters 122 (2023), URL https://doi.org/10.1063%2F5.0153345.
  * Stoner (1938) E. C. Stoner, Proceedings of the Royal Society of London. Series A. Mathematical and Physical Sciences 165, 372 (1938).
  * Šmejkal et al. (2022) L. Šmejkal, J. Sinova, and T. Jungwirth, Phys. Rev. X 12, 040501 (2022), URL https://link.aps.org/doi/10.1103/PhysRevX.12.040501.
  * Wen et al. (1989) X. G. Wen, F. Wilczek, and A. Zee, Phys. Rev. B 39, 11413 (1989), URL https://link.aps.org/doi/10.1103/PhysRevB.39.11413.
  * Chiu and Schnyder (2014) C.-K. Chiu and A. P. Schnyder, Phys. Rev. B 90, 205136 (2014), URL https://link.aps.org/doi/10.1103/PhysRevB.90.205136.
  * Chan et al. (2016) Y.-H. Chan, C.-K. Chiu, M. Y. Chou, and A. P. Schnyder, Phys. Rev. B 93, 205132 (2016), URL https://link.aps.org/doi/10.1103/PhysRevB.93.205132.
  * Belopolski et al. (2019) I. Belopolski, K. Manna, D. S. Sanchez, G. Chang, B. Ernst, J. Yin, S. S. Zhang, T. Cochran, N. Shumiya, H. Zheng, et al., Science 365, 1278 (2019), eprint https://www.science.org/doi/pdf/10.1126/science.aav2327, URL https://www.science.org/doi/abs/10.1126/science.aav2327.
  * Takahashi et al. (2017) R. Takahashi, M. Hirayama, and S. Murakami, Phys. Rev. B 96, 155206 (2017), URL https://link.aps.org/doi/10.1103/PhysRevB.96.155206.
  * Li et al. (2021) X.-P. Li, B. Fu, D.-S. Ma, C. Cui, Z.-M. Yu, and Y. Yao, Phys. Rev. B 103, L161109 (2021), URL https://link.aps.org/doi/10.1103/PhysRevB.103.L161109.
  * Suh et al. (2023) H. G. Suh, Y. Yu, T. Shishidou, M. Weinert, P. M. R. Brydon, and D. F. Agterberg, _Superconductivity of anomalous pseudospin_ (2023), eprint 2301.11458.
  * Aoki et al. (2019) D. Aoki, K. Ishida, and J. Flouquet, Journal of the Physical Society of Japan 88, 022001 (2019), ISSN 1347-4073, URL http://dx.doi.org/10.7566/JPSJ.88.022001.
  * Bastien et al. (2016) G. Bastien, A. Gourgout, D. Aoki, A. Pourret, I. Sheikin, G. Seyfarth, J. Flouquet, and G. Knebel, Phys. Rev. Lett. 117, 206401 (2016), URL https://link.aps.org/doi/10.1103/PhysRevLett.117.206401.
  * Momma and Izumi (2011) K. Momma and F. Izumi, Journal of Applied Crystallography 44, 1272 (2011), URL https://doi.org/10.1107/S0021889811038970.
  * Shindou and Nagaosa (2001) R. Shindou and N. Nagaosa, Phys. Rev. Lett. 87, 116801 (2001), URL https://link.aps.org/doi/10.1103/PhysRevLett.87.116801.
  * Taguchi et al. (2001) Y. Taguchi, Y. Oohara, H. Yoshizawa, N. Nagaosa, and Y. Tokura, Science 291, 2573 (2001), eprint https://www.science.org/doi/pdf/10.1126/science.1058161, URL https://www.science.org/doi/abs/10.1126/science.1058161.
  * Machida et al. (2010) Y. Machida, S. Nakatsuji, S. Onoda, T. Tayama, and T. Sakakibara, Nature 463, 210 (2010).
  * Fernandes et al. (2023) R. M. Fernandes, V. S. de Carvalho, T. Birol, and R. G. Pereira, _Topological transition from nodal to nodeless zeeman splitting in altermagnets_ (2023), eprint 2307.12380.
  * Elcoro et al. (2021) L. Elcoro, B. J. Wieder, Z. Song, Y. Xu, B. Bradlyn, and B. A. Bernevig, Nature Communications 12 (2021), URL https://doi.org/10.1038%2Fs41467-021-26241-8.
  * Cano et al. (2019) J. Cano, B. Bradlyn, and M. G. Vergniory, APL Materials 7, 101125 (2019), URL https://doi.org/10.1063%2F1.5124314.
  * Gallego et al. (2016) S. V. Gallego, J. M. Perez-Mato, L. Elcoro, E. S. Tasci, R. M. Hanson, M. I. Aroyo, and G. Madariaga, Journal of Applied Crystallography 49, 1941 (2016), URL https://doi.org/10.1107/S1600576716015491.
  * Hasan et al. (2022) W. Hasan, A. M. Hossain, M. Rasheduzzaman, M. A. Rahman, M. M. Hossain, K. R. Mohammad, R. Chowdhury, K. M. Hossain, M. M. Hossen, and M. Z. Hasan, RSC Adv. 12, 27492 (2022), URL http://dx.doi.org/10.1039/D2RA04273H.
  * Xu et al. (2014) H. C. Xu, Y. Zhang, M. Xu, R. Peng, X. P. Shen, V. N. Strocov, M. Shi, M. Kobayashi, T. Schmitt, B. P. Xie, et al., Phys. Rev. Lett. 112, 087603 (2014), URL https://link.aps.org/doi/10.1103/PhysRevLett.112.087603.
  * Friedemann et al. (2016) S. Friedemann, H. Chang, M. Gamża, P. Reiss, X. Chen, P. Alireza, W. Coniglio, D. Graf, S. Tozer, and F. Grosche, Scientific Reports 6, 25335 (2016).
  * Reiss et al. (2022) P. Reiss, S. Friedemann, and F. M. Grosche, Phys. Rev. B 106, 205131 (2022), URL https://link.aps.org/doi/10.1103/PhysRevB.106.205131.
  * Cavanagh et al. (2022) D. C. Cavanagh, T. Shishidou, M. Weinert, P. M. R. Brydon, and D. F. Agterberg, Physical Review B 105 (2022), URL http://dx.doi.org/10.1103/PhysRevB.105.L020505.
  * Wilde et al. (2021) M. A. Wilde, M. Dodenhöft, A. Niedermayr, A. Bauer, M. M. Hirschmann, K. Alpin, A. P. Schnyder, and C. Pfleiderer, Nature 594, 374 (2021), URL https://doi.org/10.1038%2Fs41586-021-03543-x.


# Local and Nonlocal Transport Spectroscopy in Planar Josephson Junctions

A. Banerjee  Center for Quantum Devices, Niels Bohr Institute, University of
Copenhagen, Universitetsparken 5, 2100 Copenhagen, Denmark  O. Lesser
Department of Condensed Matter Physics, Weizmann Institute of Science,
Rehovot, Israel 76100  M. A. Rahman  Center for Quantum Devices, Niels Bohr
Institute, University of Copenhagen, Universitetsparken 5, 2100 Copenhagen,
Denmark  C. Thomas  Department of Physics and Astronomy, and Birck
Nanotechnology Center, Purdue University, West Lafayette, Indiana 47907 USA
T. Wang  Department of Physics and Astronomy, and Birck Nanotechnology Center,
Purdue University, West Lafayette, Indiana 47907 USA  M. J. Manfra  Department
of Physics and Astronomy, and Birck Nanotechnology Center, Purdue University,
West Lafayette, Indiana 47907 USA  School of Materials Engineering, and School
of Electrical and Computer Engineering, Purdue University, West Lafayette,
Indiana 47907 USA  E. Berg  Department of Condensed Matter Physics, Weizmann
Institute of Science, Rehovot, Israel 76100  Y. Oreg  Department of Condensed
Matter Physics, Weizmann Institute of Science, Rehovot, Israel 76100  Ady
Stern  Department of Condensed Matter Physics, Weizmann Institute of Science,
Rehovot, Israel 76100  C. M. Marcus  Center for Quantum Devices, Niels Bohr
Institute, University of Copenhagen, Universitetsparken 5, 2100 Copenhagen,
Denmark

###### Abstract

We report simultaneously acquired local and nonlocal transport spectroscopy in
a phase-biased planar Josephson junction based on an epitaxial InAs/Al hybrid
two-dimensional heterostructure. Quantum point contacts at the junction ends
allow measurement of the 2$\times$2 matrix of local and nonlocal tunneling
conductances as a function of magnetic field along the junction, phase
difference across the junction, and carrier density. A closing and reopening
of a gap was observed in both the local and nonlocal tunneling spectra as a
function of magnetic field. For particular tunings of junction density, gap
reopenings were accompanied by zero-bias conductance peaks (ZBCPs) in local
conductances. End-to-end correlation of gap reopening was strong, while
correlation of local ZBCPs was weak. A simple, disorder-free model of the
device shows comparable conductance matrix behavior associated with a
topological phase transition. Phase dependence helps distinguish possible
origins of the ZBCPs.

Topological materials obey a bulk-boundary correspondence, establishing a
connection between the topological index of the bulk and the number of
boundary modes Qi and Zhang (2011); Fu and Kane (2008). In one-dimensional
topological superconductors (1D-TSCs) Oreg  _et al._ (2010); Lutchyn  _et al._
(2010), these considerations imply that the bulk modes undergo a
characteristic closing and reopening of the superconducting gap whenever the
system is driven through a topological phase transition. In this situation,
the gap reopening is connected to the appearance of zero-energy Majorana modes
at the two ends Sau  _et al._ (2010); Stanescu  _et al._ (2011); Tewari  _et
al._ (2012); Rainis  _et al._ (2013). Several experimental works have reported
zero-bias conductance peaks (ZBCPs) at the ends of wires or 1D structures,
identified as signatures of 1D-TSCs Mourik  _et al._ (2012); Nadj-Perge  _et
al._ (2014); Das  _et al._ (2012); Deng  _et al._ (2016); Nichele  _et al._
(2017); Grivnin  _et al._ (2019); Fornieri  _et al._ (2019). However, in most
of these cases, though not all Vaitiekėnas  _et al._ (2018); Banerjee  _et
al._ (2022), an associated gap closing and reopening was not observed in
tunneling conductance.

An emerging method that allows simultaneous observation of end modes and bulk
gap behavior is nonlocal spectroscopy, where measurement of the tunneling
current between the ends of a device provides information about the bulk
Rosdahl  _et al._ (2018). This technique requires a three-terminal (3T)
configuration, and has been theoretically explored for nanowires Rosdahl  _et
al._ (2018); Pan  _et al._ (2021); Danon  _et al._ (2020); Hess  _et al._
(2021); Wang  _et al._ (2021). Nonlocal transport experiments, also in
nanowires, were used to probe symmetries of the conductance matrix arising
from current conservation and measure the local charge of Andreev bound states
Ménard  _et al._ (2020); Danon  _et al._ (2020). Experiments in short nanowire
segments identified end-to-end correlation between extended Andreev bound
states Anselmetti  _et al._ (2019). In long nanowire segments, local
conductance showed ZBCPs while the gap in the nonlocal spectrum remained
closed Puglia  _et al._ (2021), suggesting non-topological ZBCPs arising from
strong disorder Bagrets and Altland (2012); Prada  _et al._ (2012); Liu  _et
al._ (2017); Pan  _et al._ (2020); Pan and Das Sarma (2020). These experiments
demonstrated the importance of combining local and nonlocal transport to
differentiate trivial and potentially topological ZBCPs.

Figure 1:  Device and measurement setup. False-color micrograph of a
representative device showing three-terminal configuration. Meandering
perforations etched onto the superconducting leads allow partial depletion of
the semiconductor using gate voltage $V_{\rm sc}$. Al loop allows phase
biasing of the junction with a small out-of-plane magnetic field, $B_{\perp}$.
An in-plane magnetic field, $B_{\parallel}$ is applied parallel to the S-N
interfaces. Voltage biases $V_{\rm T}$ and $V_{\rm B}$ are applied to the top
and bottom ohmic contacts through the current amplifiers (CA). Gates $V_{\rm
top(bot)}$ and $V_{\rm t(b)qpc}$ form QPCs at the junction ends. $V_{1}$
controls carrier density in the junction. All connections to the device are
via $\sim 1-2$k$\Omega$ fridge wires and filters, see Supplemental Material
for details.

Planar Josephson junctions (PJJs) of superconductor-semiconductor hybrids have
recently emerged as a promising alternative platform for topological
superconductivity, providing several knobs that can control a possible
topological superconducting phase, including, notably, the novel control
parameter of the phase difference across the junction Hell  _et al._ (2017);
Pientka  _et al._ (2017); Setiawan  _et al._ (2019); Fornieri  _et al._
(2019); Dartiailh  _et al._ (2021); Ren  _et al._ (2019); Ke  _et al._ (2019);
Haxell  _et al._ (2022).

In this Letter, we extend previous work Banerjee  _et al._ (2022) by
investigating nonlocal conductance spectroscopy, measured simultaneously with
local spectroscopy, on a 3T PJJ device with quantum point contacts (QPCs) at
both ends, measured as a function of in-plane magnetic field along the
junction, phase difference across the junction, and carrier density within the
junction. Our main observation is a closing and reopening of the
superconducting gap in nonlocal conductance correlated with the appearance of
ZBCPs in local conductance. We find that the gap closing and reopening in both
local and nonlocal conductance is robust against variations of the junction
carrier density, but the observation of ZBCPs at one or both ends require
careful tuning of voltages on the junction and QPC gates.

To help interpret these results, we investigate numerically a simple model of
a PJJ Banerjee  _et al._ (2022). Within the model, a gap reopening in nonlocal
conductances (without ZBCP) appearing together with a ZBCP in local
conductance is characteristic of a topological phase transition.

Figure 1 shows a micrograph of one of the devices, along with a schematic
electrical circuit. The PJJ can be probed by a pair of integrated QPCs at the
ends of the junction, and phase-biased by applying a small ($\sim 0.1$ mT
scale) out-of-plane magnetic field through a superconducting loop.

The device was fabricated on a molecular-beam-epitaxy grown heterostructure
stack with a shallow InAs quantum well separated from a top Al layer by an
In0.75Ga0.25As barrier. A combination of wet etching of the Al layer and deep
wet etching of the semiconductor stack was used to define the superconducting
loop, the Josephson junction and the mesa with a U-shaped trench. A patch of
the mesa (with Al removed) within the loop was contacted by a layer of Ti/Au
to form an internal submicron ohmic contact to enable bottom-end tunneling
spectroscopy. A layer of HfO2, grown by atomic layer deposition (ALD) and
patterned in a rectangular shape, was used to isolate the Ti/Au layer from the
superconducting loop and the conducting mesa. A second layer of HfO2 was
deposited globally followed by the deposition of Ti/Au gates for electrostatic
control of the junction and the QPCs.

The carrier density in the normal barrier of the JJ (width $w_{n}=100$ nm,
length $l=1.6\,\mu$m) was controlled by gate voltage $V_{1}$. Gate voltage
$V_{\rm sc}$ controlled the carrier density in the semiconductor underneath
the superconducting leads. Split gates controlled by voltages $V_{\rm tqpc}$
and $V_{\rm bqpc}$ define QPCs at the top and bottom of the junction.
Additional gate voltages $V_{\rm top}$ and $V_{\rm bot}$ controlled densities
in the normal regions outside the QPCs, and were typically fixed at $\sim+100$
mV.

Figure 2:  Magnetic field dependence of the conductance matrix. (a) Local
differential conductance $G_{\rm TT}$ and (e) Nonlocal differential
conductance $G_{\rm BT}$ measured as a function of $V_{\rm T}$ and
$B_{\parallel}$. (c) Nonlocal differential conductance $G_{\rm TB}$ and (g)
local differential conductance $G_{\rm BB}$ measured as a function of $V_{\rm
B}$ and $B_{\parallel}$. The phase bias is set to $\Phi=0$. Line-cuts at
$B_{\parallel}=0.3$ T where (b) $G_{\rm TT}$ shows a ZBCP. (h) $G_{\rm BB}$
shows a ZBCP. (d) $G_{\rm TB}$ and (h) $G_{\rm BT}$ are strongly antisymmetric
at high DC biases, and zero in a finite range around zero DC bias. Gate
voltage settings used for this measurement were $V_{\rm tqpc}=-0.37$ V,
$V_{\rm bqpc}=-0.34$ V, $V_{1}=0.085$ V and $V_{\rm sc}=-3.6$ V. Figure 3:
Phase dependence of local and nonlocal conductance spectra. Local differential
conductance $G_{\rm TT}$ (left column) and nonlocal differential conductance
$G_{\rm TB}$ (right column) measured as a function of out-of-plane magnetic
field $B_{\perp}$ and source-drain bias, $V_{\rm T}$ and $V_{\rm B}$,
respectively at three values of the in-plane magnetic field. (a) and (b)
$B_{\parallel}=0$, showing a periodic modulation of the superconducting gap in
both the local and nonlocal spectrum. The nonlocal spectrum shows a larger
amplitude of the superconducting gap than the local conductance spectrum. (c)
and (d) $B_{\parallel}=0.15$ T, shows sub-gap states that are lowered in
energy. These states are phase-asymmetric and appear in both $G_{\rm TT}$ and
$G_{\rm TB}$. (e) and (f) $B_{\parallel}=0.3$ T shows a phase-independent ZBCP
in $G_{\rm BB}$, but not in $G_{\rm TB}$. $G_{\rm TB}$ displays a
superconducting gap that is modulated periodically with $B_{\perp}$. Gate
voltage settings used for this measurement were $V_{\rm tqpc}=-0.385$ V,
$V_{\rm bqpc}=-0.38$ V, $V_{1}=0.092$ V and $V_{\rm sc}=-3.6$ V.

The 3T measurement configuration is shown in Fig. 1(a). The top ohmic contact
is a region of InAs separated from the junction by the top QPC. The bottom
ohmic contact is formed by a Ti/Au electrode, separated by the bottom QPC. The
Al loop connecting the two sides of the junction provides the third contact,
held at ground. Low-frequency AC plus DC voltage biases $V_{\rm T(B)}$ are
applied through current amplifiers (denoted CA). The measured currents $I_{\rm
T(B)}$ then yield the 2$\times$2 conductance matrix, $G_{\rm ij}=dI_{\rm
i}/dV_{\rm j}$, with ${\rm i},{\rm j}=\rm{T,B}$ via standard AC lock-in
measurements (see additional details in Supplemental Material).

For measurements shown in Fig. 2, the conductance matrix was measured as a
function of in-plane magnetic field, $B_{\parallel}$, with $V_{\rm sc}=-3.6$
V, giving a hard superconducting gap in the leads, $V_{1}=+85$ mV, giving
ZBCPs in both top and bottom local conductances at $B_{\parallel}\sim 0.3$ T,
and QPCs set to $V_{\rm tqpc}=-0.37$ V, $V_{\rm bqpc}=-0.34$ V, to yield
sizable nonlocal conductances, $\sim 0.01\times 2e^{2}/h$ near the gap edge,
$|V_{\rm T,B}|\sim 150~{}\mu$V. To compensate any coupling of $B_{\parallel}$
through the superconducting loop controlling phase $\Phi$ across the junction,
a sweep of $B_{\perp}$ was made at each value of $B_{\parallel}$ and then
sliced along cuts of constant $\Phi$ numerically by following $\Phi$-dependent
lobe features (see Methods). This allowed us to obtain the $B_{\parallel}$
dependence of the conductance matrix at fixed flux, as shown, for instance, in
Fig. 2 for $\Phi=0$.

Local conductance spectra showed a finite superconducting gap around
$B_{\parallel}=0$ [Fig. 2(a, g)]. With increasing $B_{\parallel}$, a band of
resolvable discrete states moved towards zero bias, closing the gap at
$B_{\parallel}\sim$ 0.2 T followed by its reopening. Beyond the reopening (0.2
T $<B_{\parallel}<0.4$ T), but not before, ZBCPs were observed in both $G_{\rm
TT}$ and $G_{\rm BB}$ [Fig. 2(b, h)]. In this data set, the ZBCP at the top
end splits as $B_{\parallel}$ approaches 0.4 T, whereas the bottom end ZBCP
appears to remain at zero, but diminishes in amplitude. Additionally, the
ZBCPs observed at each end do not exhibit strong correlation with respect to
variations of $V_{1}$ (see Fig. S9). Both local conductances show a final gap
closure at $B_{\parallel}\sim~{}0.45$ T.

The corresponding behavior of the nonlocal conductance spectra is shown in
Figs. 2(c, e). A predominantly antisymmetric signal is observed throughout the
measured magnetic field range, with amplitude remaining roughly uniform. The
gap in the nonlocal spectrum undergoes a closure at $B_{\parallel}\sim 0.2$ T,
at the same magnetic field as the local conductance spectrum, and is visible
in both $G_{\rm TB}$ and $G_{\rm BT}$. Both nonlocal conductances remain
strongly antisymmetric around zero bias. The nonlocal gap then reopens
obtaining a maximum at $B_{\parallel}\sim$ 0.3 T, with line-cuts shown in
Figs. 2(d, f). Notably, no ZBCP is observed. Both nonlocal conductances
disappear in a finite window around zero bias before turning on sharply at
finite VT/B. The final closure of the nonlocal gap at $B_{\parallel}\sim$ 0.45
T is more pronounced, in terms of signal strength, than the closure at
$B_{\parallel}\sim$ 0.2 T.

Local and nonlocal conductance spectra are modulated by a small perpendicular
magnetic field, $B_{\perp}$, which threads flux through the $\sim 12\,\mu$m2
superconducting loop ($B_{\perp}\sim~{}$0.17 mT corresponds to $\Phi_{0}=h/2e$
through the loop), showing the same period in $B_{\perp}$ and in phase. Around
$B_{\parallel}=0$ [Figs. 3(a, b)], the local gap in $G_{\rm TT}$ appears
smaller than the nonlocal gap in $G_{\rm BT}$. At $B_{\parallel}\sim 0.15$ T
[Figs. 3(c),(d)], flux-dependent states are lowered in energy and fill the
sub-gap spectrum in both local and nonlocal conductances. Within each flux
lobe, states are asymmetric with respect to $\Phi$. At $B_{\parallel}=$0.3 T,
a phase-independent ZBCP is measured in $G_{\rm TT}$ [Fig. 3(e)] but absent in
$G_{\rm BT}$ [Fig. 3(f)]. At this field, $G_{\rm BT}$ remains zero until a
source-drain bias of $V_{\rm T}\sim 40~{}\mu$eV at $\Phi=0$ and closes at
$\Phi=\Phi_{0}/2$. $G_{\rm TB}$ and $G_{\rm BB}$, are qualitatively similar to
$G_{\rm BT}$ and $G_{\rm TT}$ respectively (see Fig. S8).

Figure 4:  Numerical simulation of local and nonlocal conductance. (a) Local
conductance $G_{\rm BB}$ evaluated at the bottom end of the junction as a
function of $B_{\parallel}$, and $\Phi=$ 0 including thermal broadening
equivalent to a temperature of 50 mK. At $B_{\parallel}\sim 0.2$ T the
junction undergoes a topological phase transition reflected by the closing and
reopening of the superconducting gap, followed by the appearance of a ZBCP.
(b) Nonlocal conductance $G_{\rm BT}$ shows a corresponding gap-reopening
transition at $B_{\parallel}\sim 0.2$ T without the formation of a ZBCP.

We also investigated nonlocal transport at gate settings where a ZBCP was
observed in the bottom local conductance, but not the top local conductance
(see Fig. S3 and S4). In this case, a gap-reopening signature was observed in
nonlocal conductance and the nonlocal gap remained finite in the reopened
state. In other devices where reasonably strong nonlocal conductance was
observed ($\sim 0.01\times 2e^{2}/h$), the nonlocal spectrum exhibited a gap-
reopening feature. In some cases, the sub-gap nonlocal conductance in the
reopened state remained finite, indicating a soft nonlocal gap (see Figs. S10
and S11).

To help interpret characteristic features of the observed conductance matrix,
we investigate a simple model of a PJJ using the Kwant software package Groth
_et al._ (2014), as described previously Banerjee  _et al._ (2022). Here, we
extend the model by tunnel coupling the system to metallic leads at the
junction ends. Top-bottom symmetry of the model ensures that $G_{\rm
TT}=G_{\rm BB}$ and $G_{\rm TB}=G_{\rm BT}$. Figure 4(a) shows the local
conductance spectrum undergoing a topological gap-reopening transition at
$B_{\parallel}\sim 0.2$ T, followed by a ZBCP arising from a Majorana zero
mode. The corresponding nonlocal conductance spectrum, shown in Fig. 4(b),
also shows a reopening of the gap, but no sub-gap structure once the gap
reopens.

These model results support the interpretation that nonlocal conductance is
mediated by a combination of co-tunneling and crossed-Andreev reflection of
quasiparticles carried by extended Andreev bound states. These states have a
finite probability density throughout the length of the junction, including
the two ends, and therefore also appear in the local conductances. Within this
picture, Majorana zero-modes appear as zero-bias peaks only in local
conductance, not in nonlocal conductance because of their localized
probability density. This is in contrast to extended Andreev bound states,
which are expected to appear in both local and nonlocal conductances Ménard
_et al._ (2020); Hess  _et al._ (2021).

Conductance matrix signatures obtained from our model are roughly consistent
with experiment. However, as opposed to the experiment, in the numerical
simulations we find a large symmetric component of the nonlocal conductance,
comparable in strength to the antisymmetric component. Their relative strength
depends on details used to model finite temperature, disorder, and tunneling
barriers Rosdahl  _et al._ (2018); Wang  _et al._ (2021), and may explain this
discrepancy. Finally, local ZBCPs are fully correlated in the model, but lack
such correlation in experiment.

Within a non-topological interpretation of our data, an inhomogeneous chemical
potential profile produces non-topological zero-energy Andreev bound states at
the two device ends. In the model these states are not stable at zero-energy
with respect to variation of $\Phi$, which may provide a distinguishing
signature. Another non-topological scenario is the case of strong disorder,
where a proliferation of low-energy sub-gap states prevents a topological
phase transition. This scenario can produce ZBCPs in local conductance, but
tends not to show a gap-reopening in the nonlocal conductance Pan  _et al._
(2021); Puglia  _et al._ (2021).

Within a topological interpretation, the presence of a finite nonlocal gap
without strong end-to-end ZBCP correlation may arise from disorder of weak-to-
intermediate strength Woods  _et al._ (2021). Such disorder creates disjointed
topological segments in the bulk of the junction, eliminating end-to-end ZBCP
correlation while preserving the nonlocal gap together and gap-reopening
signature. Further study of conductance matrix behavior with a greater degree
of the potential along the junction would help clarify the situation.

We thank Geoff Gardner and Sergei Gronin for contributions to materials
growth, and Asbjørn Drachmann for assistance with fabrication. We thank Lucas
Casparis, Tom Dvir, Karsten Flensberg, Max Geier, Esteban Martinez, and
Andreas Pöschl, for valuable discussions. We acknowledge a research grant
(Project 43951) from VILLUM FONDEN, support from the ERC under the Horizon
2020 Research and Innovation programme (LEGOTOP No. 788715 and HQMAT No.
817799), the DFG (CRC/Transregio 183, EI 519/7-1), the BSF and NSF (2018643),
the ISF Quantum Science and Technology (2074/19), and a research grant from
Irving and Cherna Moskowitz.

## References

  * Qi and Zhang (2011) X.-L. Qi and S.-C. Zhang, Rev. Mod. Phys. 83, 1057 (2011).
  * Fu and Kane (2008) L. Fu and C. L. Kane, Phys. Rev. Lett. 100, 096407 (2008).
  * Oreg  _et al._ (2010) Y. Oreg, G. Refael, and F. von Oppen, Phys. Rev. Lett. 105, 177002 (2010).
  * Lutchyn  _et al._ (2010) R. M. Lutchyn, J. D. Sau, and S. Das Sarma, Phys. Rev. Lett. 105, 077001 (2010).
  * Sau  _et al._ (2010) J. D. Sau, S. Tewari, R. M. Lutchyn, T. D. Stanescu, and S. Das Sarma, Phys. Rev. B 82, 214509 (2010).
  * Stanescu  _et al._ (2011) T. D. Stanescu, R. M. Lutchyn, and S. Das Sarma, Phys. Rev. B 84, 144522 (2011).
  * Tewari  _et al._ (2012) S. Tewari, J. D. Sau, V. W. Scarola, C. Zhang, and S. Das Sarma, Phys. Rev. B 85, 155302 (2012).
  * Rainis  _et al._ (2013) D. Rainis, L. Trifunovic, J. Klinovaja, and D. Loss, Phys. Rev. B 87, 024515 (2013).
  * Mourik  _et al._ (2012) V. Mourik, K. Zuo, S. M. Frolov, S. R. Plissard, E. P. A. M. Bakkers, and L. P. Kouwenhoven, Science 336, 1003 (2012).
  * Nadj-Perge  _et al._ (2014) S. Nadj-Perge, I. K. Drozdov, J. Li, H. Chen, S. Jeon, J. Seo, A. H. MacDonald, B. A. Bernevig, and A. Yazdani, Science 346, 602 (2014).
  * Das  _et al._ (2012) A. Das, Y. Ronen, Y. Most, Y. Oreg, M. Heiblum, and H. Shtrikman, Nature Physics 8, 887 (2012).
  * Deng  _et al._ (2016) M. Deng, S. Vaitiekėnas, E. B. Hansen, J. Danon, M. Leijnse, K. Flensberg, J. Nygård, P. Krogstrup, and C. M. Marcus, Science 354, 1557 (2016).
  * Nichele  _et al._ (2017) F. Nichele, A. C. C. Drachmann, A. M. Whiticar, E. C. T. O’Farrell, H. J. Suominen, A. Fornieri, T. Wang, G. C. Gardner, C. Thomas, A. T. Hatke, P. Krogstrup, M. J. Manfra, K. Flensberg, and C. M. Marcus, Phys. Rev. Lett. 119, 136803 (2017).
  * Grivnin  _et al._ (2019) A. Grivnin, E. Bor, M. Heiblum, Y. Oreg, and H. Shtrikman, Nature Comm. 10, 1 (2019).
  * Fornieri  _et al._ (2019) A. Fornieri, A. M. Whiticar, F. Setiawan, E. Portolés, A. C. Drachmann, A. Keselman, S. Gronin, C. Thomas, T. Wang, R. Kallaher, _et al._ , Nature 569, 89 (2019).
  * Vaitiekėnas  _et al._ (2018) S. Vaitiekėnas, M.-T. Deng, J. Nygård, P. Krogstrup, and C. M. Marcus, Phys. Rev. Lett. 121, 037703 (2018).
  * Banerjee  _et al._ (2022) A. Banerjee, O. Lesser, M. A. Rahman, H.-R. Wang, M.-R. Lee, A. Kringhøj, A. M. Whiticar, A. C. C. Drachmann, C. Thomas, T. Wang, G. C. Gardner, M. J. Manfra, E. Berg, Y. Oreg, A. Stern, and C. M. Marcus, ArXiv 2201.03453 (2022).
  * Rosdahl  _et al._ (2018) T. O. Rosdahl, A. Vuik, M. Kjaergaard, and A. R. Akhmerov, Phys. Rev. B 97, 045421 (2018).
  * Pan  _et al._ (2021) H. Pan, J. D. Sau, and S. Das Sarma, Phys. Rev. B 103, 014513 (2021).
  * Danon  _et al._ (2020) J. Danon, A. B. Hellenes, E. B. Hansen, L. Casparis, A. P. Higginbotham, and K. Flensberg, Phys. Rev. Lett. 124, 036801 (2020).
  * Hess  _et al._ (2021) R. Hess, H. F. Legg, D. Loss, and J. Klinovaja, Phys. Rev. B 104, 075405 (2021).
  * Wang  _et al._ (2021) G. Wang, T. Dvir, N. van Loo, G. P. Mazur, S. Gazibegovic, G. Badawy, E. P. Bakkers, L. P. Kouwenhoven, and G. de Lange, arXiv preprint arXiv:2110.05373 (2021).
  * Ménard  _et al._ (2020) G. C. Ménard, G. L. R. Anselmetti, E. A. Martinez, D. Puglia, F. K. Malinowski, J. S. Lee, S. Choi, M. Pendharkar, C. J. Palmstrøm, K. Flensberg, C. M. Marcus, L. Casparis, and A. P. Higginbotham, Phys. Rev. Lett. 124, 036802 (2020).
  * Anselmetti  _et al._ (2019) G. L. R. Anselmetti, E. A. Martinez, G. C. Ménard, D. Puglia, F. K. Malinowski, J. S. Lee, S. Choi, M. Pendharkar, C. J. Palmstrøm, C. M. Marcus, L. Casparis, and A. P. Higginbotham, Phys. Rev. B 100, 205412 (2019).
  * Puglia  _et al._ (2021) D. Puglia, E. A. Martinez, G. C. Ménard, A. Pöschl, S. Gronin, G. C. Gardner, R. Kallaher, M. J. Manfra, C. M. Marcus, A. P. Higginbotham, and L. Casparis, Phys. Rev. B 103, 235201 (2021).
  * Bagrets and Altland (2012) D. Bagrets and A. Altland, Phys. Rev. Lett. 109, 227005 (2012).
  * Prada  _et al._ (2012) E. Prada, P. San-Jose, and R. Aguado, Phys. Rev. B 86, 180503 (2012).
  * Liu  _et al._ (2017) C.-X. Liu, J. D. Sau, T. D. Stanescu, and S. Das Sarma, Phys. Rev. B 96, 075161 (2017).
  * Pan  _et al._ (2020) H. Pan, W. S. Cole, J. D. Sau, and S. Das Sarma, Phys. Rev. B 101, 024506 (2020).
  * Pan and Das Sarma (2020) H. Pan and S. Das Sarma, Phys. Rev. Research 2, 013377 (2020).
  * Hell  _et al._ (2017) M. Hell, M. Leijnse, and K. Flensberg, Phys. Rev. Lett. 118, 107701 (2017).
  * Pientka  _et al._ (2017) F. Pientka, A. Keselman, E. Berg, A. Yacoby, A. Stern, and B. I. Halperin, Phys. Rev. X 7, 021032 (2017).
  * Setiawan  _et al._ (2019) F. Setiawan, A. Stern, and E. Berg, Phys. Rev. B 99, 220506 (2019).
  * Dartiailh  _et al._ (2021) M. C. Dartiailh, W. Mayer, J. Yuan, K. S. Wickramasinghe, A. Matos-Abiague, I. Žutić, and J. Shabani, Phys. Rev. Lett. 126, 036802 (2021).
  * Ren  _et al._ (2019) H. Ren, F. Pientka, S. Hart, A. T. Pierce, M. Kosowsky, L. Lunczer, R. Schlereth, B. Scharf, E. M. Hankiewicz, L. W. Molenkamp, _et al._ , Nature 569, 93 (2019).
  * Ke  _et al._ (2019) C. T. Ke, C. M. Moehle, F. K. de Vries, C. Thomas, S. Metti, C. R. Guinn, R. Kallaher, M. Lodari, G. Scappucci, T. Wang, _et al._ , Nature Comm. 10, 1 (2019).
  * Haxell  _et al._ (2022) D. Z. Haxell, E. Cheah, F. Křížek, R. Schott, M. F. Ritter, M. Hinderling, W. Belzig, C. Bruder, W. Wegscheider, H. Riel, and F. Nichele, arXiv (2022), 2204.05619 .
  * Groth  _et al._ (2014) C. W. Groth, M. Wimmer, A. R. Akhmerov, and X. Waintal, New Journal of Physics 16, 063065 (2014).
  * Woods  _et al._ (2021) B. D. Woods, S. D. Sarma, and T. D. Stanescu, arXiv preprint arXiv:2103.06880 (2021).
  * Martinez  _et al._ (2021) E. A. Martinez, A. Pöschl, E. B. Hansen, M. A. Y. van de Poll, S. Vaitiekėnas, A. P. Higginbotham, and L. Casparis, arXiv:2104.02671 (2021).
  * Peierls (1933) R. Peierls, Zeitschrift für Physik 80, 763 (1933).
  * Tinkham (1996) M. Tinkham,  _Introduction to Superconductivity_ , 2nd ed., International Series in Pure and Applied Physics (McGraw-Hill, New York, 1996).

Supplemental Material

Wafer structure: The heterostructure stack used in this work consists of an
InAs two-dimensional quantum well in epitaxial contact with aluminum, grown by
molecular beam epitaxy. The wafer was grown on an insulating InP substrate and
comprises a 100-nm-thick In0.52Al0.48As matched buffer, a 1$\mu$m thick step-
graded buffer realized with alloy steps from In0.52Al0.48As to In0.89Al0.11As
(20 steps, 50 nm/step), a 58 nm In0.82Al0.18As layer, a 4 nm In0.75Ga0.25As
bottom barrier, a 7 nm InAs quantum well, a 10 nm In0.75Ga0.25As top barrier,
two monolayers of GaAs and a 7 nm film of epitaxially grown Al. The top Al
layer was grown without breaking the vacuum, in the same growth chamber. Hall
effect measurements performed in Hall bar devices of the same material, with
Al etched away, indicated a peak electron mobility peak $\mu=43,000$ cm2/V-s
at a carrier density of $n=8\times 10^{11}$ cm-2, corresponding to a peak
electron mean free path of $l_{e}\sim$ 600 nm. This suggests that our devices
are quasiballistic along the length $l\sim 3l_{e}$ and ballistic in the width
direction $w_{n}\sim l_{e}/6$. Transport characterization of an etched
aluminum Hall bar revealed an upper critical field $B_{\parallel}\sim$ 2.5 T,
indicating that the Al layer remains superconducting well above the collapse
of the induced superconducting gap at $B_{\parallel}\sim$ 0.5 T.

Device fabrication: Devices were fabricated using standard electron beam
lithography. Devices on the same chip were electrically isolated from one
another using a two-step mesa etch process, first by removing Al with Al
etchant Transene D, and then a standard III-V chemical wet etch using a
solution H2O : C6H8O7 : H3PO4 : H2O2 (220:55:3:3) to etch the mesa until a
depth of $\sim$300 nm. This step also defined the U-shaped trench and the
patch of mesa that eventually formed the submicron ohmic contact. In the next
lithography step, the Al layer was selectively removed leaving behind the
Josephson junction with a flux loop. Particularly, this step also removed Al
from the internal ohmic mesa patch. A patterned layer of dielectric,
comprising 15 nm thick HfO2 grown at 90∘C using atomic layer deposition (ALD),
was then deposited to galvanically isolate the ohmic Ti/Au layer from the rest
of the device. This was followed by the deposition of Ti/Au layers (5 nm/300
nm) for the inner ohmic contact. Next, a global layer of 15 nm thick HfO2 was
deposited over the entire sample to serve as the gate dielectric. Gates were
defined using electron beam lithography followed by e-beam evaporation of
Ti/Au layers of thickness (5 nm/20 nm) for finer structures and (5 nm/350 nm)
for the bonding pads.

Three-terminal electrical transport measurements: Electrical transport
measurements were performed in an Oxford Triton dilution refrigerator at a
base temperature of $\sim$20 mK using standard AC lock-in techniques. The
superconducting loop was grounded by connecting it to the fridge ground
through a $\sim$0.9 k$\Omega$ line. The top and the bottom ohmic contacts were
connected to low-impedance current-to-voltage converters, through
$\sim$1.5–2.0 k$\Omega$ line resistances. The grounding pin of each current-
to-voltage converter was biased through a combination of AC+DC voltages
$V_{\rm t(b)}+V_{\rm T(B)}$, at the top and bottom ends respectively. The AC
biases, $V_{\rm t}$ and $V_{\rm b}$, were generated by two lock-in amplifiers,
with the same excitation amplitude (3 $\mu$V), but different frequencies,
$f_{\rm t}=31.5$ Hz and $f_{\rm b}=77.5$ Hz respectively. The DC biases
$V_{\rm T(B)}$ were generated using two low-noise DC voltage sources. The
voltage output of each current-to-voltage converter was measured using two
lock-in amplifiers, operating at frequencies $f_{\rm t}$ and $f_{\rm b}$,
requiring a total of 4 lock-in amplifiers for each element of the 2$\times$2
conductance matrix.

We find experimentally that the conductance $G_{ij}$ depends strongly on
$V_{j}$ but only weakly on $V_{i}$ (unless $i=j$), as shown in Figs. S1 and
S2. This is expected theoretically Danon  _et al._ (2020). Therefore, for
measurement of $G_{\rm TT}$ and $G_{\rm BT}$, we sweep the DC bias $V_{\rm T}$
and set $V_{\rm B}=0$, and for measurement of $G_{\rm BB}$ and $G_{\rm TB}$,
we sweep the DC bias $V_{\rm B}$ and set $V_{\rm T}=0$.

Voltage-divider effects: In the three-terminal geometry, measurement circuit
effects arising from finite line impedances can lead to spurious voltage-
divider effects as discussed in Ref. Martinez  _et al._ (2021). DC and AC
voltage drops across line-resistances connected to the top ($R_{\rm T}$),
bottom ($R_{\rm B}$), and ground ($R_{\rm G}$) terminals of the devices
produce contributions to the local and nonlocal conductance signals, and can
be corrected using the method outlined in Ref. Martinez  _et al._ (2021).

Briefly, the $2\times 2$ conductance matrix $G_{ij}$ ($i$, $j$ = T, B)
measured as a function of voltages $V_{\rm T,B}$ applied outside the cryostat,
is converted to a corrected conductance matrix $G^{\prime}_{ij}$ obtained as a
function of $U_{\rm T,B}$ , the voltages estimated to be applied at the top
and bottom terminals of the device inside the cryostat, using the relation
${\bf G}^{\prime}={\bf G}[\bf{I}_{2}-{\bf Z}{\bf G}]^{-1}$. Here $\bf
Z=\begin{bmatrix}R_{\rm T}&R_{\rm T}+R_{\rm G}\\\ R_{\rm B}+R_{\rm G}&R_{\rm
B}\\\ \end{bmatrix}$ is the impedance matrix.

As an example, we depict in Fig. S12 the corrected conductance matrix for data
shown in Fig. 2. The line resistances were measured as $R_{\rm T}\simeq
2~{}{\rm k}\Omega$, $R_{\rm B}\simeq 2~{}{\rm k}\Omega$, and $R_{\rm G}\simeq
0.9~{}{\rm k}\Omega$. As shown in Fig. S12(b), the voltage-divider corrections
do not produce qualitative changes to the uncorrected data.

Phase-biasing and magnetic field alignment: Magnetic field to the sample is
applied using a three-axis ($B_{x}$, $B_{y}$, $B_{z}$)=(1 T, 1 T, 6 T) vector
rotate magnet. The sample is oriented with respect to the vector magnet such
that $B_{x}$ is nominally parallel to $B_{\perp}$ [Fig. 1(c)] and $B_{z}$ is
nominally parallel to $B_{\parallel}$ [Fig. 1(a)]. However, because of sample
misalignment, $B_{z}$ has a small but finite contribution to $B_{\perp}$,
which controls the flux through the phase-biasing loop. For a finite value of
$B_{z}$, is is therefore necessary to identify the value of $B_{x}$ that leads
to lines of constant flux. At zero $B_{z}$, the value of $B_{x}$ at which the
superconducting gap is maximised corresponds to zero and multiples of
$\Phi_{0}$, while distinct phase-slips appear at odd multiples of
$\Phi_{0}/2$. This allows calibration of the flux through the device at zero
$B_{z}$. At finite $B_{z}$, the superconducting gap acquires a phase-
asymmetric dispersion, and the maxima of the gap cannot be used to track lines
of constant flux. Instead, we use the phase-slips, which become more prominent
at high values of $B_{z}$, to identify lines of constant flux through the
device. The values of $B_{x}$ required for this compensation are typically
linear as a function of $B_{z}$.

Numerical simulations: To model our device, we use an extension of the
Hamiltonian proposed in Refs. Pientka  _et al._ (2017); Hell  _et al._ (2017)
utilizing a two-layer structure to account for finite thickness and include
orbital effects, as described in more detail in Ref. Banerjee  _et al._
(2022). Here we briefly describe the essentials of the model. The rectangular
device contains three parts: normal region in the middle with width $w_{n}$,
and superconducting regions of width $w$ on each side. In the Nambu basis
$(\psi_{\uparrow},\psi_{\downarrow},\psi_{\downarrow}^{\dagger},-\psi_{\uparrow}^{\dagger})$,
the Bogoliubov–de Gennes Hamiltonian is given by

$\displaystyle H$
$\displaystyle=\left[-\frac{\partial^{2}_{x}+\partial^{2}_{y}}{2m^{*}}-t_{\perp}\nu_{x}-\mu(y)+i\alpha(z)\left(\partial_{x}\sigma_{y}-\partial_{y}\sigma_{x}\right)\right]\tau_{z}$
(S1) $\displaystyle+\frac{g(y)\mu_{\rm
B}B_{\parallel}}{2}\sigma_{x}+\Delta(y,z)\tau^{+}+\Delta^{*}(y,z)\tau^{-}.$

Here $\sigma,\tau,\nu$ are Pauli matrices acting in spin, electron-hole, and
layer basis, respectively. $m^{*}$ is the effective mass of electrons in the
semiconductor, $\alpha(z)$ is the layer-dependent Rashba spin-orbit coupling
strength, $t_{\perp}$ is the inter-layer hopping amplitude, $B_{\parallel}$ is
the magnetic field applied along the junction, and $\mu_{\rm B}$ is the Bohr
magneton. The g-factor $g(y)$ and the chemical potential $\mu(y)$ are $g_{\rm
N},\mu_{\rm N}$ in the normal region and $g_{\rm S},\mu_{\rm S}$ in the
superconducting regions. $\Delta(y)$ is the superconducting pairing potential
which is non-zero only in the superconducting region and includes a phase
difference $\phi$:

$\Delta(y)=\left\\{\begin{array}[]{ll}\Delta
e^{i\phi/2}&\,\frac{w_{n}}{2}<y<w+\frac{w_{n}}{2}\\\
0&\,|y|<\frac{w_{n}}{2}\\\ \Delta
e^{-i\phi/2}&\,-w-\frac{w_{n}}{2}<y<-\frac{w_{n}}{2}.\end{array}\right.$ (S2)

Orbital effects due to the in-plane magnetic field are accounted for by adding
a complex phase factor to the inter-layer hopping amplitude $t_{\perp}$,
corresponding to the vector potential $\vec{A}=B_{\parallel}y\hat{z}$ Peierls
(1933). The parallel magnetic field also induces linear phase growth along the
junction’s cross section Tinkham (1996), which is modeled as an additional
modulation $\propto B_{\parallel}yz$ to the superconducting phase.

For the purpose of numerical simulations, we discretize the Hamiltonian to a
tight-binding model on a square lattice of spacing $a=10$ nm. Simulations are
performed with the following parameters: $m^{*}=0.026\,m_{e}$,
$\Delta=140~{}\mu$eV, $t_{\perp}=10$ meV, $l=10~{}\mu$m, $w_{n}=100$ nm,
$w=200$ nm, $w_{z}=10$ nm, $\mu_{\rm SC}=3.6$ meV, $\mu_{\rm N}=3.3$ meV. The
g-factors are taken to be $g_{\rm N}=8$ and $g_{\rm S}=4$, and the spin-orbit
coupling constants are $\alpha(0)=15$ meV nm, $\alpha(1)=-\alpha(0)/4$.
Transport simulations were performed using the Kwant software package Groth
_et al._ (2014). We attach two normal leads to the top and bottom of the
device, and add tunneling barriers of $25~{}{\rm meV}$ between the leads and
the device. We then numerically solve the scattering problem and extract the
conductance matrix. To account for thermal broadening, we perform a
convolution of the conductance with the derivative of the Fermi-Dirac
function. The temperature we used in these calculations is 50 mK.

Figure S1:  Conductance matrix as a function of bottom QPC voltages in Device
1 (a)-(b) Conductance matrix measured as function of bottom QPC gate voltage
$V_{\rm bqpc}$ and top end source-drain bias $V_{\rm T}$, with $V_{\rm B}=$ 0
and $V_{\rm tqpc}=-0.37$ V. (a) $G_{\rm TT}$ shows a hard superconducting gap
at the top end as a function of $V_{\rm T}$, but does not respond to $V_{\rm
bqpc}$. (b) $G_{\rm BT}$ exhibits anti-symmetry and is pinched off by $V_{\rm
bqpc}$. (c)-(d) Conductance matrix measured as function of bottom QPC gate
voltage, $V_{\rm bqpc}$ and bottom end source-drain bias, $V_{\rm B}$, with
$V_{\rm T}=0$ and $V_{\rm tqpc}=-0.37$ V. (c) $G_{\rm TB}$, measured at the
top end, is pinched off by $V_{\rm bqpc}$ at the bottom end, reflecting
nonlocality. (d) $G_{\rm BB}$ shows a hard superconducting gap and is pinched
off by $V_{\rm bqpc}$. Figure S2:  Conductance matrix as a function of top QPC
voltages in Device 1. (a)-(b) Conductance matrix measured as function of top
QPC gate voltage, $V_{\rm tqpc}$ and top end source-drain bias, $V_{\rm T}$,
with $V_{\rm B}=0$ and $V_{\rm bqpc}=-0.38$ V. (a) $G_{\rm TT}$ shows a hard-
superconducting gap and is pinched off by $V_{\rm tqpc}$. (b) $G_{\rm BT}$,
measured at the bottom end, is pinched off by $V_{\rm tqpc}$ at the top end,
reflecting nonlocality. (c)-(d) Conductance matrix measured as function of top
QPC gate voltage $V_{\rm tqpc}$ and bottom end source-drain bias $V_{\rm B}$,
with $V_{\rm T}=$ 0 and $V_{\rm bqpc}=-0.38$ V. (c) $G_{\rm TB}$, exhibits
antisymmetry with respect to zero-bias and is pinched off by $V_{\rm tqpc}$.
(d) $G_{\rm BB}$ shows a hard superconducting gap at the bottom end as a
function of $V_{\rm B}$, but does not respond to $V_{\rm tqpc}$. Figure S3:
Conductance matrix as a function of in-plane magnetic field in a regime
supporting ZBCP only at bottom end in Device 1. (a) Local differential
conductance $G_{\rm TT}$ and (e) Nonlocal differential conductance $G_{\rm
BT}$ measured as a function of $V_{\rm T}$ and $B_{\parallel}$. (c) Nonlocal
differential conductance $G_{\rm TB}$ and (g) local differential conductance
$G_{\rm BB}$ measured as a function of $V_{\rm B}$ and $B_{\parallel}$. The
phase bias is set to $\Phi=0$. Line-cuts at $B_{\parallel}=0.3$ T where (b)
$G_{\rm TT}$ is gapped. (h) $G_{\rm BB}$ shows a ZBCP. (d) $G_{\rm TB}$ and
(f) $G_{\rm BT}$ are strongly antisymmetric at high DC biases, and zero in a
finite range of DC bias around zero. Figure S4:  Phase dependence of the
conductance matrix in Device 1. Conductance matrix measured as a function of
out-of-plane magnetic field at different values of in-plane magnetic field.
(a)–(d) $B_{\parallel}$=0. All four conductance matrix elements exhibit
periodic modulation of the superconducting gap. (e)–(h) $B_{\parallel}$=0.15
T. Subgap states are lowered in energy. These states create a phase-asymmetric
modulation of the gap within each flux lobe. (i)–(l) $B_{\parallel}$=0.2 T.
The superconducting gap is closed at all values of $\Phi$, the closure being
most evident in $G_{\rm BT}$. (m)–(p) $B_{\parallel}$=0.33 T. The
superconducting gap reopens with a phase independent ZBCP visible in $G_{\rm
BB}$. Other conductance matrix elements show a phase-modulated superconducting
gap. Figure S5: Top gate voltage dependence of the conductance matrix in
Device 1. Conductance matrix measured as a function of top gate voltage
$V_{1}$ at different values of in-plane magnetic field. We set $\Phi=0$.
(a)–(d) $B_{\parallel}=0$. All four conductance matrix elements show a maximum
superconducting gap of size $\sim 150~{}\mu$V. $G_{\rm TT}$ shows a
suppression of the gap near $V_{1}\sim$ 5 mV, which is not observed in the
other conductance matrix elements. We interpret this as a local feature at the
top end. (e)–(h) $B_{\parallel}=0.15$ T. Subgap states are lowered in energy,
but the superconducting gap remains finite as a function of $V_{1}$. (i)–(l)
$B_{\parallel}=0.32$ T. Finite superconducting gap measured in the reopened
state, particularly visible in the nonlocal conductance matrix elements.
$G_{\rm BB}$ shows a ZBCP centered around $V_{1}=17$ mV. (m)–(p)
$B_{\parallel}=0.5$ T. Final closure of the gap, visible in all four
conductance matrix elements. Note that the nonlocal conductance signal is
strongly antisymmetric and concentrated around zero bias. Figure S6:
Conductance matrix as a function of bottom QPC pinch-off in Device 1.
Conductance matrix measured as a function of bottom QPC gate voltage $V_{\rm
bqpc}$ at $V_{1}$= 17 mV , $B_{\parallel}=$ 0.33 T and $\Phi=0$. (a) $G_{\rm
TT}$ shows a superconducting gap without subgap structure. The gap is
modulated weakly by $V_{\rm bqpc}$. (b) $G_{\rm TB}$ is strongly antisymmetric
and shows a superconducting gap. The signal is diminished and ultimately
pinched off by $V_{\rm bqpc}$. (c) $G_{\rm BT}$ is strongly antisymmetric and
shows a superconducting gap. The signal is diminished and ultimately pinched
off by $V_{\rm bqpc}$. (d) $G_{\rm BB}$ shows a sharp ZBCP that becomes
visible at $V_{\rm bqpc}\sim-0.31$ V and persists until pinch-off as a
function of $V_{\rm bqpc}$. Figure S7:  Conductance matrix as a function of
in-plane magnetic field in regimes supporting ZBCPs at both ends in Device 1
(a)–(d) $V_{1}~{}$= 85 mV and (e)–(h) $V_{1}~{}$= 90 mV. Both local
conductances show a gap-reopening transition at $B_{\parallel}=$ 0.2 T
followed by ZBCPs. The nonlocal conductances are strongly antisymmetric with
respect to bias, undergo a gap-reopening transition at $B_{\parallel}=$ 0.2 T,
and show a reopened gap that is devoid of sub-gap structure. Figure S8:
Conductance matrix measured as a function of out-of-plane magnetic field at
different values of in-plane magnetic field, at $V_{1}~{}$= 85 mV in Device 1.
(a)–(d) $B_{\parallel}=0$. All four conductance matrix elements exhibit
periodic modulation of the superconducting gap. (e)–(h) $B_{\parallel}=0.15$
T. Subgap states are lowered in energy. These states create a phase-asymmetric
modulation of the gap within each flux lobe. (i)–(l) $B_{\parallel}$=0.2 T.
The superconducting gap is closed at all values of $\Phi$. (m)–(p)
$B_{\parallel}=0.3$ T. The superconducting gap reopens with a phase
independent ZBCP visible in both $G_{\rm TT}$ and $G_{\rm BB}$. The nonlocal
conductance matrix elements show a phase-modulated superconducting gap, with
no sub-gap structure. Figure S9:  Conductance matrix measured as a function of
top gate voltage $V_{1}$ at different values of in-plane magnetic field in
Device 1 (a)–(d) $B_{\parallel}=0$. All four conductance matrix elements show
a superconducting gap of $\sim 150~{}\mu$V, with some dependence on $V_{1}$.
(e)–(h) $B_{\parallel}=0.15$ T. Sub-gap states are lowered in energy, but the
superconducting gap remains finite as a function of $V_{1}$. (i)–(l)
$B_{\parallel}=0.3$ T. Finite superconducting gap measured in the reopened
state, particularly visible in the nonlocal conductance matrix elements. Both
$G_{\rm TT}$ and $G_{\rm BB}$ show ZBCPs existing in voltage range $V_{1}=$
0.08 – 0.1 V. (m)–(p) $B_{\parallel}=0.5$ T. Final closure of the gap, visible
in all four conductance matrix elements. $\Phi=0$ for all plots. Figure S10:
Conductance matrix in Device 2 (a)–(d) In-plane magnetic field dependence of
the conductance matrix at $V_{1}~{}=185$ mV. Both local conductances show a
gap-reopening transition at $B_{\parallel}=$ 0.2 T, with a ZBCP appearing only
at the top end. The nonlocal conductances are strongly antisymmetric with
respect to bias, undergo a gap-reopening transition at $B_{\parallel}=0.2$ T.
(e)–(h) Conductance matrix at $B_{\parallel}=0.3$ T, $\Phi=0$ measured as a
function of $V_{1}$. $G_{\rm TT}$ shows a stable ZBCP, whereas $G_{\rm BB}$
remains pinched off for $V_{1}<0.19$ V. Both $G_{\rm TB}$ and $G_{\rm BT}$
exhibit finite sub-gap nonlocal conductance for $V_{1}=0.19-0.21$ V Figure
S11: Phase dependence of the conductance matrix in Device 2. Conductance
matrix measured as a function of out-of-plane magnetic field at different
values of in-plane magnetic field, at $V_{1}~{}=0.185$ V. (a)–(d)
$B_{\parallel}=0$. All four conductance matrix elements exhibit periodic
modulation of the superconducting gap. (e)–(h) $B_{\parallel}=0.18$ T. Subgap
states are lowered in energy. These states create a phase-asymmetric
modulation of the gap within each flux lobe, particularly visible in $G_{\rm
TT}$. (i)–(l) $B_{\parallel}=0.3$ T. The superconducting gap reopens with a
phase independent ZBCP visible in $G_{\rm TT}$. The nonlocal conductance
matrix elements show a phase-modulated superconducting gap. (m)–(p)
$B_{\parallel}=0.5$ T. The superconducting gap is closed at all values of
$\Phi$. Figure S12:  Voltage-divider effects in Device 1. (a) Conductance
matrix, uncorrected for voltage-divider effects, as a function of in-plane
magnetic field at $V_{1}~{}=~{}85$ mV (same data as Fig. 2). (b) Conductance
matrix corrected for voltage-divider effects, as a function of in-plane
magnetic field at $V_{1}~{}=~{}85$ mV. $U_{\rm T}$ and $U_{\rm B}$ represent
the estimated voltage biases appearing at the top and bottom terminals of the
device, after accounting for voltage drops across line resistances $R_{\rm
T}\simeq 2~{}{\rm k}\Omega$ connected to the top terminal, $R_{\rm B}\simeq
2~{}{\rm k}\Omega$ connected to the bottom terminal, and $R_{\rm G}\simeq
0.9~{}{\rm k}\Omega$ connected to the ground terminal. The conductance
corrections arising from voltage-divider effects do not result in any
qualitative changes.