# Multiterminal transport spectroscopy of subgap states in Coulomb-blockaded
superconductors

Rubén Seoane Souto  Center for Quantum Devices, Niels Bohr Institute,
University of Copenhagen, 2100 Copenhagen, Denmark  Division of Solid State
Physics and NanoLund, Lund University, S-22100 Lund, Sweden  Matteo M. Wauters
Center for Quantum Devices, Niels Bohr Institute, University of Copenhagen,
2100 Copenhagen, Denmark  Niels Bohr International Academy, Niels Bohr
Institute, University of Copenhagen, 2100 Copenhagen, Denmark  Karsten
Flensberg  Center for Quantum Devices, Niels Bohr Institute, University of
Copenhagen, 2100 Copenhagen, Denmark  Martin Leijnse  Center for Quantum
Devices, Niels Bohr Institute, University of Copenhagen, 2100 Copenhagen,
Denmark  Division of Solid State Physics and NanoLund, Lund University,
S-22100 Lund, Sweden  Michele Burrello  Center for Quantum Devices, Niels Bohr
Institute, University of Copenhagen, 2100 Copenhagen, Denmark  Niels Bohr
International Academy, Niels Bohr Institute, University of Copenhagen, 2100
Copenhagen, Denmark

###### Abstract

Subgap states are responsible for the low-bias transport features of hybrid
superconducting–semiconducting devices. Here, we analyze the local and
nonlocal differential conductance of Coulomb-blockaded multiterminal
superconducting islands that host subgap states with different spatial
structures. The emerging patterns of their transport spectroscopy are used to
characterize the possible topological nature of these devices and offer the
possibility of controlling their transport properties. We develop a next-to-
leading order master equation to describe the multiterminal transport in
superconductors with both strong Coulomb interactions and multiple subgap
states, coupled with metallic leads. We show that the nonlocal differential
conductance characterizes the spatial extension of the subgap states and
signals the presence of degenerate bound states with a finite support on
different parts of the device. Additionally, it displays sharp sign changes as
a function of the induced charge of the superconductor, signaling energy
crossings among its lowest excited states.

##  I Introduction

Hybrid systems fabricated with superconducting and semiconducting materials
are considered a key ingredient for the development of quantum technologies
and have demonstrated a potential for quantum information storage and
processing [1, 2, 3, 4]. In particular, the advent of topological
superconductors has raised the hope of realizing protected quantum devices,
robust to decoherence induced by local perturbations [5].

In one dimension, topological superconductors host Majorana subgap modes [6],
predicted to be nonabelian anyons and to give rise to several nonlocal
properties [7, 8, 9, 10]. These Majorana modes are exponentially localized at
the system edges but, importantly, each pair of them defines a nonlocal
quasiparticle state, with vanishing energy in the ideal case. In this context,
distinguishing between local and nonlocal subgap states becomes relevant for
identifying devices that are potentially in the topological regime.

Subgap states dominate the transport properties of superconducting devices at
low voltage bias and they can be experimentally detected with tunneling
spectroscopy. In particular, zero-energy Majorana modes in grounded
topological superconductors reveal themselves through a quantized conductance
peak at zero bias [11, 12, 13, 14]. However, also trivial states can mimic
similar transport features in local spectroscopy, see for example Refs. [15,
16, 17, 18, 19, 20, 21, 22, 23, 24]. This motivated further analyses of the
nonlocal conductance in multiterminal setups, where the device couples to a
grounded superconductor [25, 26, 27, 28, 29, 30, 31, 32], with experiments
performed in nanowires [33, 34, 35, 36, 37, 38], and devices based on two-
dimensional electron gases [39, 40, 41, 42].

When the superconductor is not grounded but exhibits a strong charging energy,
the intricate interplay between superconductivity and many-body interactions
yields a rich phenomenology. In the topological regime, Coulomb blockade
effects provide effective tools to initialize, read out, and measure coherence
times in topological qubits, thus constituting a key element to design
platforms for quantum information devices [43, 44, 45, 46, 47, 48, 49, 50, 51,
52, 53].

Transport properties of floating superconducting islands have been extensively
studied experimentally in semiconducting nanowires, with the goal of
extracting quasiparticle poisoning times [54, 55], showing the doubling of
Coulomb peaks in the high-field regime [56, 57, 58, 59], and investigating the
exponential protection [56] and spin polarization [60] of subgap states. These
systems have been explored using interferometry [61] and reflectometry
techniques [62, 63]. However, most of the studies so far, both experimental
and theoretical, focused on local conductance measurements of two-terminal
devices; nonlocal transport properties of multiterminal Coulomb-blockaded
devices remain largely unexplored.

A more complex typology of devices is based on two nanowires coupled via a
floating superconductor, sharing a common charging energy. In the case where
the device couples to more than two leads, the topological Kondo effect is
predicted to appear, and its realization would be an unequivocal signature of
the Majorana nonlocality [64, 65, 66]. As such, systems with Coulomb-blockaded
double nanowires are at the focus of recent experiments [67, 68] and
theoretical investigations [69].

In this work, we propose multiterminal quantum transport measurements as a way
to distinguish between local and nonlocal subgap states in Coulomb-blockaded
superconducting devices. We focus on two geometries, namely a single and a
double superconducting nanowire device. We employ a second-order master
equation approach that allows for studying both sequential and cotunneling
transport in superconducting devices with strong Coulomb interaction, multiple
subgap states, and the coupling to two or more metallic leads. In particular,
we show that sequential transport at low bias voltages is blocked whenever the
lowest-energy subgap state does not couple two of the leads. Indeed, even
though the Coulomb interaction is highly nonlocal, transport is dominated by
the local projection of subgap states, providing indications about their
spatial structure. Moreover, a three-terminal Coulomb-blockaded device
displays abrupt changes in the current direction at fixed bias potential as a
function of the charge $n_{g}$ induced on the superconducting island, due to
energy crossings in the excitation spectrum, similarly to non-superconducting
Coulomb-blockaded quantum dots [70]. Concerning the double nanowire geometry,
we investigate the local and nonlocal conductance associated to different
spatial structures of the subgap states, aiming to facilitate the
interpretation of future experiments.

The rest of the paper is organized as follows. In Sec. II we present the
effective zero bandwidth model we use to describe the low-energy behavior of
the devices and the formalism used to compute the sequential and cotunneling
current signals. In Sec. III we investigate both the two-terminal and three-
terminal conductance of a single nanowire with one or multiple subgap states.
In Sec. IV we analyze the transport signatures associated to different
structures of the subgap states in the double nanowire geometry. We draw our
conclusions and discuss future developments in Sec. V. The appendices contain
more supporting results and a detailed derivation of the tunneling rates
calculation.

Figure 1: (a) Sketch of a Coulomb-blockaded hybrid superconducting nanowire
device, hosting a subgap state (red). In this example, the nanowire couples to
two leads (indicated with $L$ and $R$). (b) Single-particle energies of the
island, containing a subgap state at energy $\varepsilon$ (black line),
quasiparticle excitations above the gap $\Delta$, and the lead states
displaying Fermi-Dirac distributions (yellow) at different chemical
potentials. (c) Lowest eigenergies in the presence of interactions as a
function of the offset charge, $n_{g}=\alpha V_{g}$ with $\alpha$ being the
lever arm. Blue (orange) lines correspond to states with odd (even) number of
electrons in the island. We show an example where the nanowire hosts a
spinless subgap state, where the solid and dashed lines correspond to the
ground state and the excited state in each charge sector.

##  II Model and methods

We consider a superconducting island coupled to normal metallic electrodes,
given by the Hamiltonian

$H=H_{\rm L}+H_{\rm I}+H_{\rm T}\,,$ (1)

where the leads are described by

$H_{\rm L}=\sum_{\rm\nu,k,\sigma}\left(\xi_{\nu
k\sigma}-\mu_{\nu}\right)c^{\dagger}_{\nu k\sigma}c_{\nu k\sigma}\,,$ (2)

with the electron energy, $\xi_{\nu k\sigma}$, and the annihilation operator,
$c_{\nu k\sigma}$, referring to an electron in lead $\nu$ with momentum $k$
and spin $\sigma\in\\{\uparrow,\downarrow\\}$. We assume that each lead
remains in internal equilibrium, described by a Fermi-Dirac distribution
$n_{\rm{F}}$ with chemical potential $\mu_{\nu}$ and temperature $T$.

The superconducting island is described by

$H_{\rm
I}=\sum_{j,\sigma}\varepsilon_{j\sigma}\gamma_{j\sigma}^{\dagger}\gamma_{j\sigma}+E_{\rm
el}(N)\,,$ (3)

where $j$ labels a state with energy $\varepsilon_{j\sigma}$, and the
Bogoliubov-de Gennes (BdG) operators read

$\gamma_{j\sigma}=u_{j}d_{j\sigma}-\rho_{\sigma}\,v_{j}e^{-i\phi}d_{j\bar{\sigma}}^{\dagger}\,.$
(4)

Here, $d$ and $d^{\dagger}$ are annihilation and creation operators of
electrons in the island, $u_{j}$ and $v_{j}$ are the BdG coefficients for the
subgap states, $\phi$ is the superconducting phase operator, such that
$e^{-i\phi}$ annihilates a Cooper pair in the island, and $\rho_{\sigma}=\pm$
for $\sigma=\uparrow,\downarrow$.

The electrostatic repulsion term is given by

$E_{\rm el}(N)=E_{\rm C}(N-n_{\rm g})^{2}\,,$ (5)

where $n_{\rm g}$ is the dimensionless gate-induced charge offset and $N$ the
(excess) electron number operator, accounting for the fermionic occupation of
the states and the number of Cooper pairs on the island. Figure 1(a) sketches
the device described by this model, in the case of only two leads and a single
nanowire. Panels (b) and (c) represent the superconducting quasiparticle
spectrum and the many-body energy levels, respectively.

Finally, the tunneling between leads and superconducting island is described
by the Hamiltonian

$\begin{split}H_{T}&=\sum_{\nu,k,j,\sigma}\left(t_{\nu k\sigma}c_{\nu
k\sigma}^{\dagger}d_{j\sigma}+{\rm H.c.}\right)\end{split}$ (6)

where $t$ is the tunneling amplitude. The corresponding tunneling rates are
defined as $\Gamma_{\nu,j,\sigma}=2\pi\rho_{F\nu}|t_{\nu,j,\sigma}|^{2}$, and
we consider them $k$-independent (wideband limit), where $\rho_{F\nu}$ is the
density of states at the Fermi level of the lead $\nu$.

In this work, we employ a zero-bandwidth approximation for the superconducting
island, meaning that we consider only two relevant energy levels
$\varepsilon_{0,\uparrow}=\varepsilon$ and
$\varepsilon_{1,\sigma}=\Delta>\varepsilon$. The first represents a low-energy
spin-polarized subgap state whose spatial structure and degeneracy will be the
main focus of our work. The state at energy $\Delta$ is a spin-degenerate
level and represents an effective description of the proximitized
superconducting continuum which extends through the whole hybrid device. This
approximation models a short nanowire where the quasiparticle spectrum above
the proximity induced superconducting gap displays a quantization determined
by the finite system size. The corresponding lowest energy excitation can
mediate elastic cotunelling processes through the device. A possible
alternative is to consider a number of degenerate states at energy $\Delta$,
each coupled to a single lead. Each of these states models a metallic
quasiparticle continuum with a strong relaxation towards the bottom of the
continuum band. With this choice, the states at energy $\Delta$ can not
mediate elastic cotunneling transport. Throughout the rest of the paper, we
focus on the first scenario, where each nanowire hosts a single extended
superconducting continuum, to avoid an excessive increase of the computational
run time. However, this choice does not qualitatively influence the transport
features at low bias, where the conductance is dominated by the subgap states.
See Fig. 6 in App. B for a comparison of the two cases.

###  II.1 Formalism

In this work, we focus on the regime where $\Gamma_{\nu,j,\sigma}\ll T$. In
this weak coupling limit, a perturbation theory in the tunneling amplitude
provides accurate results. Since we are not interested in the thermal
excitations of the superconducting device, we additionally consider
temperatures much smaller than the remaining energy scales. To compute the
tunneling rates, we use the $\mathcal{T}$-matrix formalism, describing the
transition probability between two states

$\Gamma_{i\to f}=2\pi\left|\left\langle
f\left|\mathcal{T}\right|i\right\rangle\right|^{2}W_{if}\delta(E_{i}-E_{f})\,,$
(7)

where $W_{if}$ weights the rate through thermal distributions of the electrons
in the leads, $\delta$ is the Dirac delta functional ensuring energy
conservation between the initial and final states ($E_{i}=E_{f}$) and

$\mathcal{T}=H_{\rm T}+H_{\rm T}\frac{1}{E_{i}-H_{\rm L}-H_{\rm
I}+i0^{+}}\mathcal{T}\,,$ (8)

which can be truncated at the desired order. In Eq. (7) and in the rest of the
paper, we set $\hbar=1$. The linear term in Eq. (8) describes the sequential
tunneling; the higher order terms define, among other processes, the
cotunneling contributions, which become progressively more important when the
tunnel amplitudes $t$ increase. In Sec. A of the appendix, we provide
expressions for the sequential and the (second-order) cotunneling rates, used
in this work to evaluate the transport through the superconducting islands.

The quantum state of the island is described by
$\left|a\right\rangle$=$\left|N,N_{\rm C},\textbf{n}\right\rangle$, where
$N_{C}$ is the number of Cooper pairs in the island and n is a vector
representing the occupation of the subgap and continuum excited states. The
time derivative of the occupation probability of a given state is given by the
master equation

$\dot{P}_{a}=\sum_{b}\left[-\Gamma_{a\to b}P_{a}+\Gamma_{b\to
a}P_{b}\right]\,.$ (9)

Here, we assume that coherences between different states in the
superconducting island are negligible; indeed, we consider either models
without near degeneracies where coherences cannot develop or with no coupling
between degenerate states, due to different parity or spatial support. In the
stationary limit, the system does not evolve in time and
$\dot{P}^{\rm{stat}}_{b}=0$. These conditions, together with the normalization
$\sum_{b}P^{\rm{stat}}_{b}=1$, form a linear system of equations for the
occupation probabilities of the island states. Using the resulting stationary
distribution, the current flowing from lead $\nu$ to the device is determined
by the island transition rate to other states as

$\displaystyle I_{\nu}=\sum_{a,b}\left[s\Gamma^{\rm{seq}}_{\nu,b\to
a}+\sum_{\nu^{\prime}\neq\nu}\delta
n(\nu)\,\Gamma^{\rm{cot}}_{\nu,\nu^{\prime},b\to
a}\vphantom{\frac{1}{2}}\right]P^{\rm{stat}}_{b}\,,$ (10)

where we take $s=+1$ ($s=-1$) for electrons tunneling in (out) of the island
and we set the electron charge to unity. Here, $\delta n(\nu)$ is the net
charge transferred from lead $\nu$ in a cotunneling process, and the sum runs
over all the possible rates connecting the island state $\left|b\right\rangle$
with any state $\left|a\right\rangle$ (see Appendix A.2 for more detail on the
cotunneling rates $\Gamma^{\rm{cot}}$).

In Eqs. (9) and (10) we are considering processes that change the island
occupation by $\pm$1 electron, while neglecting local and crossed Andreev
reflection processes [71]. This approximation is justified when the charging
energy is the dominant energy scale, in particular $E_{c}\gg\Delta$, meaning
that there is a large energy penalty for changing the charge of the device by
$2e$. Throughout the rest of the paper, we will indeed consider this strongly
Coulomb-blockaded regime.

##  III Single nanowire geometry

In this section, we focus on resolving the spatial structure of subgap states
in a single nanowire device coupled to two or three normal metallic leads.
Regarding the spatial profile of the subgap state, we consider three
situations: (i) a single subgap state delocalized at both ends of the nanowire
but vanishing in its bulk, which we refer to as a “Majorana-like scenario”,
(ii) two degenerate subgap states each localized close to one end of the
device, (iii) an extended subgap state that couples to all leads independently
from their position. In all these cases, we consider that the subgap states
are strongly spin-split, behaving as effectively spinless. This situation can
be achieved by a strong magnetic field or by proximity with ferromagnetic
materials, and corresponds to standard models to achieve the topological
regime [72]. Our goal is to distinguish between these possible scenarios based
on the structure of the cotunneling steps and the nonlocal conductance.

###  III.1 Two-terminal device

We first focus on a two-terminal device, with leads coupled to the nanowire
ends, which can host either a single subgap state or two degenerate ones
located at the ends, see Figs. 2 (a) and (b). For simplicity, we choose a
symmetric voltage drop at the contacts between the device and $L$ and $R$
leads by setting their chemical to $\mu_{L}=+V/2$ and $\mu_{R}=-V/2$, where
$V$ is the voltage bias.

Figure 2: Sketches of a nanowire coupled to two leads, hosting an extended (a)
and two local states (b). Panels (c) and (d) show the differential conductance
(in logscale) corresponding to cases (a) and (b), respectively, for subgap
states with energies $\varepsilon=0.3\Delta$. We show results for
$\Delta=0.2E_{C}$ and $\Gamma=5\cdot 10^{-5}E_{C}$. The color tics denote the
dominant inelastic cotunneling processes, appearing at $V=2\Delta$ (blue),
$V=\Delta\pm\varepsilon$ for the even/odd Coulomb valley (red), and
$V=2\varepsilon$, $\Delta+\varepsilon$ (green). Panels (e-j) show examples of
the different cotunneling processes.

If the subgap levels have zero energy $\varepsilon=0$, which is the case for
ideal Majorana zero energy modes, the Coulomb blockade structure is
qualitatively similar in the two cases (data not shown).Instead, when
$\varepsilon\neq 0$, the most evident difference is the strong reduction of
the conductance peaks at $V=0$ in the case of local states. This can be seen
by comparing the conductance at $V\sim 0$ and half-integer $n_{g}$ in panels
(c) and (d) of Fig. 2. Indeed, when the subgap states do not couple to both
leads, transport is suppressed at small biases because the subgap states
cannot support resonant tunneling of electrons. This results in a much weaker
conductance for $|V|<2\varepsilon,\ \Delta-\varepsilon$ with respect to the
situation where the subgap state connects both ends of the nanowire. However,
zero-bias tunneling is still possible through elastic cotunneling, with an
amplitude $\sim\Gamma_{L}\Gamma_{R}$ [71].

We also observe a clear difference in the structure of the cotunneling steps,
highlighted by the colored tics in Figs. 2 (c) and (d). In the case of an
extended state, two inelastic cotunneling steps appear at
$V=\Delta\pm\varepsilon$ for the even and the odd valleys, corresponding to
transitions to the lowest excited states (red tics) [60]. In the even valley,
it corresponds to the splitting of a Cooper pair, whose electrons end up in
the subgap and the continuum of states. In the odd valley, the lowest
excitation corresponds to promoting the electron in the subgap state to the
continuum. There are two additional steps due to the breaking of a Cooper pair
into two electrons that enter the quasiparticle continuum (blue tics).

In the case with two states with local support, each coupled to a single lead,
the number of cotunneling steps increases and, in particular, the threshold
for inelastic cotunneling in the even-parity valleys decreases substantially.
Indeed, the lowest excited state with even parity has an energy given by the
sum of the energies of the two subgap states [see panel (j) in Fig. 2]. In the
case considered in the figure with the two states having the same energy, this
threshold is $2\varepsilon$. We note that this process does not require the
creation of a quasiparticle excitation above the superconducting gap, reducing
the conductance step onset with respect to the nondegenerate situation.

###  III.2 Three-terminal device

Next, we analyse a three terminal device, where two leads ($L$ and $R$) couple
to the nanowire ends and the third ($M$) probes the center of the device.
Therefore, the $M$ lead acts as a probe of the subgap support at the middle of
the wire, distinguishing local and nonlocal states. The $L$ and $R$ leads are
biased, with a chemical potential $\mu_{L}=+V/2$ and $\mu_{R}=-V/2$,
respectively, while we consider that $M$ is grounded ($\mu_{M}=0$). The goal
is to use the nonlocal conductance $G_{M}=\frac{\mathrm{d}I_{M}}{\mathrm{d}V}$
to detect whether the subgap state extends inside the device or is localized
only on the edges. These two scenarios are sketched in Figs. 3 (a) and (b).
Furthermore, we explore the situations where the subgap states’ energy
$\varepsilon$ is zero, Figs. 3 (c)-(d), or finite, Figs. 3(e)-(f).

Figure 3: Nonlocal differential conductance
$G_{M}=\frac{\mathrm{d}I_{M}}{\mathrm{d}V}$ through a grounded lead coupled to
the middle of the wire. (a) and (b) respectively depict situations where the
middle lead does not and does couple to the subgap state (red blobs). Panels
(c) and (d) show the nonlocal conductance for a case where the subgap state
energy is $\varepsilon=0$, while panels (e) and (f) show cases for
$\varepsilon=0.3\Delta$. Panels (c) and (e) refer to the subgap state
structure shown in (a) while panels (d) and (f) correspond to (b). The
remaining parameters are the same as in Fig. 2 and all the tunnel couplings
are considered to be equal.

While the local conductance (through $L$ and $R$) shows little dependence on
the presence of a third lead and how it couples with the device (data not
shown), the nonlocal conductance displays interesting features. First, notice
that $G_{M}$ is an odd function of the voltage bias $V$, thanks to the
symmetry $R\leftrightarrow L$ we impose in our model: an inversion of the bias
only changes the direction of the current flow between the $L$ and $R$ lead,
while the current $I_{M}$ remains insensitive to the sign of $V$ since the
lead is always grounded. Hence, the differential conductance
($G_{M}=\frac{\mathrm{d}I_{M}}{\mathrm{d}V}$) is an odd function of the
voltage bias. When the lead does not couple to the subgap state, the
conductance at $V=0$ is suppressed, as clearly seen by comparing left and
right panels in Fig. 3. Even though this effect is somehow trivial because low
bias transport only involves the lowest energy level in a gapped system, its
consequence is of great importance: the middle lead probes the density of
states inside the wire. Hence, the $M$ lead can resolve the spatial profile of
subgap states and detect whether they have a nonvanishing projection on a
specific portion of the device.

The non-local conductance also shows a peculiar sign dependence on $n_{g}$,
illustrated by the sharp $G_{M}$ jumps in Fig. 3. These sign changes appear
whenever there is a crossing between ground states or excited states with
different parities, corresponding to the crossing of the parabolas in Fig. 1
(c) (see Appendix B for more details). Their dependence on $V$, instead,
arises from the necessity of having a chemical potential on the $L$ and $R$
leads large enough to populate the excited states coupled also with the $M$
lead. This behavior can be better understood by considering the energy
dependence of the superconducting island states on $n_{g}$ and the sequential
tunneling rates (see App. A.1 for their derivation). For instance, the rates
describing the tunneling of an electron from the $M$ lead to an empty
quasiparticle state $\gamma_{j}$ is $\Gamma_{M,i\to
f}=\Gamma_{M}\left|u_{j}\right|^{2}n_{\rm F}(E_{f}-E_{i})$, where $i$ and $f$
label the initial and final states with energy $E_{i}$ and $E_{f}$.
$\Gamma_{M}\left|u_{j}\right|^{2}$ is the effective tunneling rate in the
wideband limit multiplied by the local projector on the particle-like
component of $\gamma_{j}$. Since the temperature is small, the Fermi factor
$n_{\rm F}$ is a step-like function, activating the tunneling process mainly
when $E_{f}-E_{i}<0$. The opposite is true when we examine the process for
extracting one electron from the occupied quasiparticle state $\gamma_{j}$,
meaning that it is activated when $E_{f}-E_{i}>0$. Hence, the sign of the
current in the $M$ lead changes abruptly when $n_{g}$ tunes the energy levels
of the device across the resonance $E_{i}=E_{f}$, as adding/removing electrons
from M becomes more favorable. We show more details in Fig. 5.

In conclusion, a three terminal Coulomb-blockaded device has a highly tunable
differential conductance alongside the standard transport suppression inside
the Coulomb diamonds. As shown in panels (c)-(f) of Fig. 3, this is a robust
feature that requires neither the presence of a zero energy subgap state nor a
third lead coupled only to the superconducting continuum, even though the
latter is useful to have an extended bias window where the $M$ lead is
effectively decoupled from the device.

##  IV Double nanowire

Figure 4: Double wire system coupled via a trivial superconductor allowing the
exchange of Cooper pairs, so the system has a common charging energy.
\raisebox{-.9pt} {$1$}⃝ One subgap state connects to $L$ and $R$, while
another couples to $M$ (a). Panels (b) and (c) are the local conductance
through the lead $L$ and the nonlocal conductance through the grounded $M$
lead, respectively. \raisebox{-.9pt} {$2$}⃝ Situation for three local states
coupled to the different leads, \raisebox{-.9pt} {$3$}⃝ an extended state
coupled to all the leads, and \raisebox{-.9pt} {$4$}⃝ a local state that
couples to $R$, while $L$ and $R$ couple to a different extended state. The
parameters are the same as in Fig. 2.

In this section, we extend the analysis of the nonlocal conductance to double
nanowire setups; two parallel semiconductor nanowires are covered by the same
floating superconducting island. The common superconductor allows the exchange
of Cooper pairs between the wires, making them share a common charging energy,
but it does not allow the tunneling of excitations above the gap from one wire
to the other. In this context, we therefore introduce two separate states at
energy $\Delta$ to model the quasiparticle continua in the two wires. In
contrast, we assume that subgap states can delocalize between the two wires.
This setup for the device is inspired by recent experimental achievements [73,
74, 75, 68] and the scope of exploring the topological Kondo effect when the
system couples to more than two leads [64, 65, 66]. We consider that the
nanowire ends couple to three different leads, as sketched in panels (a), (d),
(g), and (j) of Fig. 4. As before, we consider that the $L/R$ leads are
symmetrically biased, while $M$ is grounded. A different biasing condition is
shown in Fig. 7 in the appendix.

Regarding the spatial structure of subagap states, we identify four possible
situations: \raisebox{-.9pt} {$1$}⃝ two subgap states that extend along a
single nanowire, \raisebox{-.9pt} {$2$}⃝ a subgap state localized at each
lead-device interface, \raisebox{-.9pt} {$3$}⃝ a common subgap state that
couples with all three terminals, \raisebox{-.9pt} {$4$}⃝ two subgap states
localized respectively at the left and right ends of the double nanowire
device. These scenarios are sketched in panels (a), (d), (g), and (j) of Fig.
4, while the panels below each sketch show the related local $G_{L}$ and
nonlocal $G_{M}$ differential conductances. In all cases, we consider
degenerate subgap states with energy $\varepsilon=0.3\Delta$. Subgap states
with different energies would be easily detected by the different sizes of the
Coulomb diamonds. The differences between the four cases are summarized in
Table 1.

\raisebox{-.9pt} {$1$}⃝ |  \raisebox{-.9pt} {$2$}⃝
---|---
Sequential tunneling absent in $G_{M}$ |  Sequential tunneling absent in all terminals
2 cotunneling steps |  3 cotunneling steps
[1.2pt,gray] \raisebox{-.9pt} {$3$}⃝ |  \raisebox{-.9pt} {$4$}⃝
Sequential tunneling present in all terminals |  Sequential tunneling absent in $G_{R}$
2 cotunneling steps |  3 cotunneling steps
|  Negative differential conductance in $G_{L}$
Table 1: Summary of the main transport features for the four situations shown
in Fig. 4.

Let us start from case \raisebox{-.9pt} {$1$}⃝, where each nanowire hosts a
single subgap state with support on both ends. In this case, there is no probe
at the middle of any of the wires. Therefore, the transport features cannot
distinguish between trivial extended states or nonlocal Majorana-like subgap
states. The associated local conductance,
$G_{L}=\frac{\mathrm{d}I_{L}}{\mathrm{d}V}$, and the non local one,
$G_{M}=\frac{\mathrm{d}I_{M}}{\mathrm{d}V}$, are reported in Figs. 4 (b) and
(c), respectively. Sequential tunneling processes contribute to the current
between $L$ and $R$ at low bias, while they are suppressed for the $M$ lead. A
single sequential tunneling line is observed in $G_{L}$ at low bias voltages;
in the odd valleys, the stationary distribution $P^{\rm stat}$ has indeed a
significant contribution from the configuration where a quasiparticle is
frozen in the subgap state in the lower wire. This is due to the imbalance
between the $\Gamma_{M}$ rates for incoming and outgoing particles at small
temperatures. When the lower wire subgap state is occupied, the
superconducting island cannot be excited by removing this quasiparticle via
sequential tunneling through leads $L$ or $R$, thus suppressing the leading-
order transport mechanism in $G_{L}$. This is evident when comparing the odd
valleys in Fig. 2(c) and Fig. 4(b): in proximity of the charge degeneracy
points the sequential tunneling lines moving towards the odd valleys disappear
and cotunneling becomes the dominating process at low bias. Also the weak
conductance through the $M$ lead at low bias voltage is due to cotunneling
processes, where the device exchanges one electron with either $L$ or $R$, and
$M$, keeping the total charge on the superconducting island constant.

In case \raisebox{-.9pt} {$2$}⃝, when a trivial subgap state localizes at the
interface with each lead, the conductance is suppressed for $|V|<\varepsilon$.
This is due to the local support of the subgap states, which cannot directly
mediate sequential transport between different leads. Therefore, transport is
dominated by elastic cotunneling at low bias voltages. The nonlocal
conductance $G_{M}$, Fig. 4(f), displays feint resonances in the even valleys.
The reason behind is that tunneling processes between the superconducting
island and the $M$ lead involve (virtual) changes of the number of Cooper
pairs in the island. The sign of the current depends on the charge difference
between the ground and the lowest excited state.

Case \raisebox{-.9pt} {$3$}⃝ corresponds to a single subgap state coupled with
all three leads. Therefore, the same Coulomb structure is present in the local
and the nonlocal conductance, although the latter is an odd function of the
bias $V$ and exhibits sign changes when there are level crossings in the many-
body spectrum of the superconducting device. The sequential tunneling of
electrons dominates transport, leading to strong conductance features at the
charge degeneracy points for small bias voltages. Qualitatively, the signal of
the $M$ lead is equivalent to that shown in Fig. 3(f) for a single nanowire
with an extended bound state.

The last possibility we analyse, case \raisebox{-.9pt} {$4$}⃝, corresponds to
two degenerate subgap states, localized one in the left part and one in the
right part of the device. This situation is of particular interest, as it
might correspond to a double nanowire geometry where Majorana modes localized
at each device ends strongly hybridize due to their small spatial separation.
The local conductance $G_{L}$ displays a parity dependent negative
differential conductance (NDC) region at small bias, see Fig. 4(k). This is
caused by a quasiparticle being trapped in the subgap state coupled to the $R$
lead, blocking the current between the $L$ and $M$ leads. The appearance of
the NDC region for either positive or negative bias is due to the electron- or
hole-like excitation trapped in the right subgap state. In this regime,
sequential transport through the right lead is suppressed, leading to a
current reduction (not shown). For $|V|>2\varepsilon$, the sequential
transport channel to the right lead is open again, therefore making
$|G_{L}|\neq|G_{M}|$. It also provides a relatively fast relaxation mechanism
for the trapped quasiparticles in the right lead and the NDC region
disappears.

The four cases show conductance steps inside the Coulomb valleys at finite
voltage values. These steps appear when $V$ matches the system excitation
energy, opening the inelastic cotunneling channel: exchange of two electrons
between the leads and the island, leaving it in an excited state (although
keeping the island charge invariant). Depending on whether a bound state
couples simultaneously to the $L$ and $R$ leads, cases \raisebox{-.9pt} {$1$}⃝
and \raisebox{-.9pt} {$3$}⃝, or not, cases \raisebox{-.9pt} {$2$}⃝ and
\raisebox{-.9pt} {$4$}⃝, the number of cotunneling steps vary from 2 to 3
(colored tics in Fig. 4). A similar behavior has been described for the single
wire situation, see discussion around Fig. 2.

Finally, we show in Appendix B another biasing situation, where $L$ and $M$
are symmetrically biased and $R$ is grounded. We find that local and nonlocal
transport can distinguish between the considered four situations also in that
case.

##  V conclusions

In this paper, we analysed how the transport features of a multiterminal
superconducting device with strong charging energy depend on the number and
spatial structure of subgap bound states. In particular, we investigated the
role of the spatial extent of the subgap states, which might not couple to all
the leads attached to the device; we showed that the nonlocal differential
conductance in multiterminal devices allows for a qualitative characterization
of their spatial profiles for biases below the superconducting gap. We
considered a zero bandwidth model to describe the low-energy features of a
superconducting floating island consisting of a proximitized single
semiconducting nanowire or a pair of parallel nanowires. These systems are
indeed known to host subgap states that determine the transport properties of
the device. We focused on a situation where the charging energy $E_{\rm C}$ is
the dominant energy scale, and we adopted a second-order master equations
approach which allows us to characterize both the sequential tunneling of
single electrons and the inelastic cotunneling events.

The sequential tunneling signal, dominating for small leads-device tunneling
strength, gives information on the energy and spatial structure of the subgap
states mediating transport. When two leads are not coupled by the lowest-
energy states, the zero bias conductance is strongly suppressed, leaving only
a faint cotunneling feature. In this way, transport can determine whether a
state has support on the two ends of the wire. In the same way, additional
leads can be added to gain spatial resolution inside the wire. Therefore, the
absence of sequential tunneling conductance at low bias is a way to
discriminate between trivial extended states and possible topological states
with only support at the ends of the wires. Moreover, the number and voltage
of inelastic cotunneling steps allow us to determine the number and the
energies of the subgap states. The cotunneling signal is therefore another
indicator that can be used to characterize the spatial structure of the subgap
states. The described features hold for an arbitrary number of discrete subgap
states with energy $|\varepsilon|>T,\ \Gamma$. For $\varepsilon<T,\ \Gamma$,
instead, different subgap state configurations may result in the same
qualitative transport features, thus hindering their spatial characterization.

Finally, electron transport in multiterminal Coulomb-blockaded devices results
in an interesting pattern of peaks with positive and negative differential
nonlocal conductance, depending on the induced charge $n_{g}$ and on the
voltage bias $V$. This allows to switch the direction of the current flowing
in the grounded lead, or suppress it, without changing the potential
difference between the source and the drain but only by tuning the induced
charge on the whole superconducting island.

We presented results based on a perturbative analysis of the transport
properties, valid when the temperature is larger than the coupling between the
leads and the device. Complementary methods are needed to describe the low-
temperature and strong-coupling regimes, where electron correlations effects
are important, in the non-equilibrium situation [76, 77, 78]. The master
equation approach we presented, however, is less computationally intensive and
provides a clear picture of the transport mechanisms as long as non-
perturbative effects can be neglected. .

## aknowledgements

We acknowledge support from the Danish National Research Foundation, the
Danish Council for Independent Research $|$ Natural Sciences, the European
Research Council (Grant Agreement No. 856526), the Swedish Research Council,
and NanoLund. M.W. and M.B. are supported by the Villum Foundation (Research
Grant No. 25310). This project has received funding from the European Union’s
Horizon 2020 research and innovation program under the Marie Sklodowska-Curie
grant agreement No. 847523 “INTERACTIONS”.

##  Appendix A Tunneling rates

###  A.1 Sequential tunneling rates

In this article, we focus on the limit where the charging energy is the
largest energy scale ($E_{\rm C}\gg\max\left(|\Delta|,\,|V|\equiv|\mu_{\rm
L}-\mu_{\rm R}|\right)$). In this regime, transport is dominated by processes
where one electron is exchanged between the superconducting island and the
leads. For $\Gamma$ much smaller than any other energy scale, sequential
tunneling processes yield the most prominent features in the conductance,
providing the usual Coulomb diamond structure. The sequential rates are given
by

$\displaystyle\begin{split}\Gamma^{\rm{seq}}_{\nu,\left|N,N_{\rm
C},\textbf{n}\right\rangle\to\left|N+1,N_{\rm
C},\textbf{n'}\right\rangle}&=\Gamma_{\nu}\left|u_{j}\right|^{2}n_{\rm
F}(E_{f}-E_{i}-\mu_{\nu})\,\\\ \Gamma^{\rm{seq}}_{\nu,\left|N,N_{\rm
C},\textbf{n}\right\rangle\to\left|N-1,N_{\rm
C},\textbf{n'}\right\rangle}&=\Gamma_{\nu}\left|u_{j}\right|^{2}n_{\rm
F}(\mu_{\nu}+E_{f}-E_{i})\,\\\ \Gamma^{\rm{seq}}_{\nu,\left|N,N_{\rm
C},\textbf{n}\right\rangle\to\left|N+1,N_{\rm
C}+1,\textbf{n'}\right\rangle}&=\Gamma_{\nu}\left|v_{j}\right|^{2}n_{\rm
F}(E_{f}-E_{i}-\mu_{\nu})\,\\\ \Gamma^{\rm{seq}}_{\nu,\left|N,N_{\rm
C},\textbf{n}\right\rangle\to\left|N-1,N_{\rm
C}-1,\textbf{n'}\right\rangle}&=\Gamma_{\nu}\left|v_{j}\right|^{2}n_{\rm
F}(\mu_{\nu}+E_{f}-E_{i})\,,\end{split}$ (11)

where $n_{\rm F}$ is the Fermi-Dirac distribution function and $\mu_{\nu}$ is
the chemical potential of the lead $\nu=L,R,M$. Here, we have used the
wideband approximation, where the tunneling rates are energy independent and
$\Gamma_{\nu}=2\pi\rho_{F}|t_{\nu}|^{2}$, with the lead density of states at
the Fermi level $\rho_{F}$. These rates induce transitions between charge
states differing by one electron, denoted through the vectors n and
$\textbf{n}^{\prime}$ accounting for the fermionic occupation of the island
states.

###  A.2 Cotunneling rates

Inside the Coulomb-blockaded region, inelastic cotunneling creates a series of
steps in the differential conductance, where the voltage bias corresponds to
energy differences between different states of the superconducting island.
These are processes where one electron is transferred between two leads,
leaving the island in an excited state. We also consider the elastic
cotunneling that, instead, describes processes in which the island energy is
conserved, and gives rise to the conductance background inside the diamond.
However, we disregard local and crossed Andreev processes, where the island
charge changes by 2e, as they are suppressed by the strong charging energy in
the system. The corresponding rates are given by Eq. (7), where

$\begin{split}F_{\nu,\nu^{\prime}}(\omega_{1})&\equiv\left|\left\langle
f\left|\mathcal{T}\right|i\right\rangle\right|^{2}=\\\
&\left|\sum_{m_{1},\nu}\frac{\left\langle f\left|H_{\rm
T}(\nu)\right|m_{1}\right\rangle\left\langle m_{1}\left|H_{\rm
T}({\nu^{\prime}})\right|i\right\rangle}{E_{m_{1}}-E_{i}-\omega_{1}}\right.\\\
&\left.+\sum_{m_{2},\nu}\frac{\left\langle f\left|H_{\rm
T}({\nu^{\prime}})\right|m_{2}\right\rangle\left\langle m_{2}\left|H_{\rm
T}(\nu)\right|i\right\rangle}{E_{m_{2}}-E_{f}+\omega_{1}}\right|^{2}\,.\end{split}$
(12)

Here, $H_{\rm T}(\nu)$ describes the tunneling between the lead $\nu$ and the
island, $\nu^{\prime}$ denotes a lead different from $\nu$, and $m_{1,2}$ are
virtual intermediate states. To derive this expression, we have imposed energy
conservation, which leads to a function dependent on the energy of the
tunneling electron from/to one of the leads, $\omega_{1}$. The cotunneling
rate can be written as

$\begin{split}\Gamma^{\rm{cot}}_{\nu,\nu^{\prime},i\to f}=2\pi\int
d\omega_{1}\,&F_{\nu,\nu^{\prime}}(\omega_{1})\,n_{\rm F}(\omega_{1}-\mu_{\rm
L})\,\\\ \times\,&n_{\rm F}(\mu_{\rm R}+E_{f}-E_{i}-\omega_{1})\,.\end{split}$
(13)

This expression for the cotunneling rates is divergent when the denominator in
Eq. (12) vanishes. To avoid the divergent behaviour, we regularize the
divergences as explained in Ref. [79]. The resulting integral can be formally
solved analytically [80], which leads to a complicated expression involving
special functions. We find that for $T/E_{\rm C}\gtrsim 10^{-3}$, it is
computationally more efficient to expand the Fermi distribution function into
a sum of complex Matsubara-Ozaki frequencies [81]

$n_{\rm
F}(\omega-\mu)=\sum_{\alpha}r_{\alpha}\frac{1}{\omega-\mu+i\beta_{\alpha}}\,,$
(14)

where $\beta_{\alpha}$ and $r_{\alpha}$ are the approximated Matsubara
frequencies and residues, respectively. Finally, Eq. (13) can be evaluated
using the residue theorem, yielding to a rather compact expression

$\begin{split}&\Gamma^{\rm{cot}}_{\nu,\nu^{\prime},i\to
f}=\,8\pi\,\mbox{Im}\sum_{\alpha=0}^{\infty}r_{\alpha}\left\\{n_{\rm
F}(\mu_{\rm R}+E_{f}-E_{i}-\mu_{\rm L}+i\beta_{\alpha})\right.\\\
&\left.\times\left[F_{\nu,\nu^{\prime}}(\mu_{\rm
L}-i\beta_{\alpha})-F_{\nu,\nu^{\prime}}(\mu_{\rm
R}+E_{f}-E_{i}+i\beta_{\alpha})\right]\right\\}\,.\end{split}$ (15)

This sum can be truncated at $\alpha\approx 100$ Matsubara-Ozaki frequencies
for the parameters used in the calculations.

##  Appendix B Supporting results

In this appendix, we present additional results to support the chosen
modelling of the considered devices and describe some details of the transport
features presented in the main text. In particular we first focus on the
processes describing the transport of the three-terminal single-nanowire
devices, which determine the transport across the lead $M$; then we address
the comparison between different modelling of the quasiparticle states above
the superconducting gap; finally, we consider a different choice of the
voltage drop across the double-nanowire devices.

Figure 5: (a) Many-body eigenergies of the island for a spinless state at zero
energy ($\Delta=0.2E_{\rm C}$). (b) Many-body energies for a subgap state at
finite energy, $\varepsilon=0.3\Delta$. (c) Differential conductance through a
middle lead that does not couple to the subgap state [see the sketch in Fig. 3
(a)]. Here we consider only the sequential tunneling contribution to the
current. Panel (d) shows the corresponding results for
$\varepsilon=0.3\Delta$. The lines in the lower panels correspond to the
excitation threshold denoted by the arrows in the upper ones. In panels (a)
and (b), sequential tunneling via the $M$ lead connects only states with
different parity and charge occupation in the states at $\Delta$. These
processes always involve a relaxation from a higher energy state to a lower
energy one, thus explaining the changes in sign of the differential
conductance when two parabolas cross. The remaining parameters are the same as
in Fig. 2

###  B.1 Features of the current across the lead $M$ in the three-terminal
nanowire

To understand the three-terminal transport features in a Coulomb-blockaded
superconducting device, it is useful to compare the conductance dependence on
the voltage bias and on the induced charge with the low-energy many-body
spectrum. In Fig. 5, we report the eigenvalue structure and the nonlocal
conductance $G_{M}=\frac{\mathrm{d}I_{M}}{\mathrm{d}V}$, close to the charge
degeneracy point, for the device sketched in Fig. 3(a): a single Coulomb-
blockaded nanowire hosting a subgap state that does not couple with the
grounded lead $M$. The lowest energies of states with odd and even parity are
represented by the continuous blue and orange lines, respectively, in Fig.
5(a)-(b). These states are connected via sequential tunneling involving a
particle transfer to or from the subgap state, which does not contribute to
the current flowing through the $M$ lead. Hence $G_{M}$ is suppressed at low
bias, as can be seen in Figs. 5(c) and (d). The current in the $M$ lead is
activated only when $V$ is large enough to populate excited states (via the
$L$ or the $R$ leads) that can then relax to a lower energy state through a
tunneling event between the superconducting continuum and the lead $M$. These
processes correspond to transitions between a dashed line (i.e., states with
excited quasiparticles in the superconducting continuum) and a continuous line
with different parity (colors) in Figs. 5 (a) and (b). These events are
highlighted by the oblique lines in panels (c) and (d), representing energy
thresholds for populating states that can contribute to the current through
$M$. The sign of the current, and hence of the conductance, changes at the
crossing between energy levels with different parity. This phenomenology is
independent from the energy $\varepsilon$ of the subgap states, whether it is
zero, as in panels (a) and (c), or finite, as in panels (b) and (d). The
latter situation only has a slightly richer structure of the Coulomb diamonds
due to the energy difference between the first excited states in the even and
odd sectors.

Figure 6: Two-terminal conductance for different subgap and above gap states
situations. Panels (a) and (c) show schematically the situation where a bound
state extends from left to right (red), while the state at higher energy can
be either extended (a) or local (c). Panels (b) and (d) show the corresponding
local conductance. Panels (e) and (g) show the situation where the island
hosts local subgap states, whose corresponding conductance is shown in (f) and
(h). Notice the NDC region appearing in panel (d): they correspond to a
particle blocked in one of the two state at energy $\Delta$ and unable to
tunnel out from the device on the other side, thus suppressing transport.
Below the energy threshold required to populate the superconducting continuum
($|V|=2\Delta$), transport is dominated by the properties of the subgap state
alone: the low-bias peaks of panels (b) and (d) are almost identical and so
are those of panels (f) and (h).

Another interesting effect shown in Fig. 3 is the sign mismatch between the
sequential tunneling and the cotunneling contribution close to the lower edge
of the Coulomb diamonds. This can be seen at $V=\Delta$, where cotunneling
(blurred signal) gives rise to negative conductance while sequential tunneling
(sharp lines) contributes with positive conductance when $n_{g}$ approaches
$5.5$ from below. Both can be understood by considering the processes mediated
by tunneling through the $M$ lead for specific values of $n_{g}$ and $V$. When
$V>0$ and $n_{g}\to 5.5$, the lowest energy excitation that can relax through
the central lead is the odd parity state with a quasiparticle in the
superconducting continuum that couples with an incoming electron to create a
Cooper pair. This process is thus associated with a positive (ingoing) current
from lead $M$. Instead, the lowest energy inelastic cotunneling step in the
odd valley corresponds to the excitation of a quasiparticle from the subgap
state to the continuum, mediated by the virtual occupation of the even parity
state with no quasiparticle present. In this second-order process, the only
tunneling event through the $M$ lead is the destruction of a Cooper pair into
the high-energy quasiparticle and an outgoing electron, which carries a
negative particle current and, thus, negative conductance. Similar arguments
explain the sign change between sequential tunneling and cotunneling in other
regions of the Coulomb diamonds.

###  B.2 One- vs two-state approximation to model the Bogoliubov quasiparticle
continuum

Concerning the modelling of the quasiparticle continuum above the gap, in Fig.
6 we compare the two-terminal conductance for a nanowire with a single
quasiparticle state at energy $\Delta$ coupled with both leads and with two
degenerate quasiparticle states at the same energy with finite supports on
each end, as sketched in the upper panels. The nanowire hosts one or two
subgap states at energy $\varepsilon$. As long as the voltage bias is smaller
than $\Delta-\varepsilon$, transport is dominated by the properties of the
subgap state and the spatial structure of the superconducting continuum does
not affect the qualitative features in conductance of the device. When the
bias is larger and the states at energy $\Delta$ become accessible via
sequential tunneling events, the “broken” continuum can trap a quasiparticle
in one of the two edges of the device, suppressing transport and inducing a
negative differential conductance region in the Coulomb diamonds. Notice,
however, how this NDC is qualitatively different from that appearing in all
data presented in the main text, as it is a sequential tunneling feature
appearing in the local conductance of two-terminal devices.

Figure 7: Upper left panels: sketch of the double wire system.
\raisebox{-.9pt} {$1$}⃝ one bound state connects $L$ and $R$, while $M$
connects to another one (a). Panels (b) and (c) are the conductance through
the grounded $M$ lead. \raisebox{-.9pt} {$2$}⃝ three local states coupled to
the different leads, \raisebox{-.9pt} {$3$}⃝ an extended state coupled to all
the leads, and \raisebox{-.9pt} {$4$}⃝ a local state that couples to $R$,
while $L$ and $R$ couple to another extended state. The parameters are the
same as in Fig. 2. The conductance we show for lead $M$ is
$G_{M}=\frac{\mathrm{d}I_{M}}{\mathrm{d}V_{R}}$, hence it is always negative,
with the exception of some small regions.

###  B.3 Example of the multiterminal transport with different voltage drops

Finally, in Fig. 7 we present the conductance results for the double nanowire
geometry with a different bias choice with respect to the data shown in Fig. 4
in the main text: here we bias symmetrically leads $L$ and $M$, while lead $R$
is left grounded. The main difference with respect to Fig. 4 is that now the
biased leads are not connected by a single superconducting continuum at energy
$\Delta$. In particular, notice the different natures of the NDC regions
appearing in panel (b) and (h): in the former, it is due to a quasiparticle
trapped in the subgap state connected to lead $M$ and, indeed, it is not
reflected in $G_{M}$. In the latter, it corresponds to a quasiparticle trapped
in one of the two superconducting continua, in the upper or in the lower
nanowire. This second case is similar to the data presented in Fig.6(d), where
the NDC is due to the broken superconducting continuum on a single nanowire,
although the conductance in Fig. 7(h) is not symmetric in $V$ because of the
presence of the third terminal.

## References

  * Janvier  _et al._ [2015] C. Janvier, L. Tosi, L. Bretheau, Ç. Ö. Girit, M. Stern, P. Bertet, P. Joyez, D. Vion, D. Esteve, M. F. Goffman, H. Pothier, and C. Urbina, Coherent manipulation of Andreev states in superconducting atomic contacts, Science 349, 1199 (2015).
  * Hays  _et al._ [2018] M. Hays, G. de Lange, K. Serniak, D. J. van Woerkom, D. Bouman, P. Krogstrup, J. Nygård, A. Geresdi, and M. H. Devoret, Direct microwave measurement of Andreev-bound-state dynamics in a semiconductor-nanowire Josephson junction, Phys. Rev. Lett. 121, 047001 (2018).
  * Hays  _et al._ [2021] M. Hays, V. Fatemi, D. Bouman, J. Cerrillo, S. Diamond, K. Serniak, T. Connolly, P. Krogstrup, J. Nygård, A. L. Yeyati, A. Geresdi, and M. H. Devoret, Coherent manipulation of an Andreev spin qubit, Science 373, 430 (2021).
  * Prada  _et al._ [2020] E. Prada, P. San-Jose, M. W. A. de Moor, A. Geresdi, E. J. H. Lee, J. Klinovaja, D. Loss, J. Nygård, R. Aguado, and L. P. Kouwenhoven, From Andreev to Majorana bound states in hybrid superconductor–semiconductor nanowires, Nature Reviews Physics 2, 575 (2020).
  * Nayak  _et al._ [2008] C. Nayak, S. H. Simon, A. Stern, M. Freedman, and S. Das Sarma, Non-Abelian anyons and topological quantum computation, Rev. Mod. Phys. 80, 1083 (2008).
  * Kitaev [2001] A. Y. Kitaev, Unpaired Majorana fermions in quantum wires, Physics-Uspekhi 44, 131–136 (2001).
  * Alicea [2012] J. Alicea, New directions in the pursuit of Majorana fermions in solid state systems, Reports on Progress in Physics 75, 076501 (2012).
  * Leijnse and Flensberg [2012] M. Leijnse and K. Flensberg, Introduction to topological superconductivity and Majorana fermions,  27, 124003 (2012).
  * Beenakker [2020] C. W. J. Beenakker, Search for non-Abelian Majorana braiding statistics in superconductors, SciPost Phys. Lect. Notes , 15 (2020).
  * Aguado [2017] R. Aguado, Majorana quasiparticles in condensed matter, La Rivista del Nuovo Cimento 40, 523 (2017).
  * Sengupta  _et al._ [2001] K. Sengupta, I. Žutić, H.-J. Kwon, V. M. Yakovenko, and S. Das Sarma, Midgap edge states and pairing symmetry of quasi-one-dimensional organic superconductors, Phys. Rev. B 63, 144531 (2001).
  * Law  _et al._ [2009] K. T. Law, P. A. Lee, and T. K. Ng, Majorana Fermion Induced Resonant Andreev Reflection, Phys. Rev. Lett. 103, 237001 (2009).
  * Flensberg [2010] K. Flensberg, Tunneling characteristics of a chain of Majorana bound states, Phys. Rev. B 82, 180516 (2010).
  * Wimmer  _et al._ [2011] M. Wimmer, A. R. Akhmerov, J. P. Dahlhaus, and C. W. J. Beenakker, Quantum point contact as a probe of a topological superconductor, New Journal of Physics 13, 053016 (2011).
  * Pikulin  _et al._ [2012] D. I. Pikulin, J. P. Dahlhaus, M. Wimmer, H. Schomerus, and C. W. J. Beenakker, A zero-voltage conductance peak from weak antilocalization in a Majorana nanowire, New Jour. Phys. 14, 125011 (2012).
  * Prada  _et al._ [2012] E. Prada, P. San-Jose, and R. Aguado, Transport spectroscopy of $NS$ nanowire junctions with Majorana fermions, Phys. Rev. B 86, 180503 (2012).
  * Kells  _et al._ [2012] G. Kells, D. Meidan, and P. W. Brouwer, Near-zero-energy end states in topologically trivial spin-orbit coupled superconducting nanowires with a smooth confinement, Phys. Rev. B 86, 100503 (2012).
  * Liu  _et al._ [2017] C.-X. Liu, J. D. Sau, T. D. Stanescu, and S. Das Sarma, Andreev bound states versus Majorana bound states in quantum dot-nanowire-superconductor hybrid structures: Trivial versus topological zero-bias conductance peaks, Phys. Rev. B 96, 075161 (2017).
  * Moore  _et al._ [2018] C. Moore, C. Zeng, T. D. Stanescu, and S. Tewari, Quantized zero-bias conductance plateau in semiconductor-superconductor heterostructures without topological Majorana zero modes, Phys. Rev. B 98, 155314 (2018).
  * Vuik  _et al._ [2019] A. Vuik, B. Nijholt, A. R. Akhmerov, and M. Wimmer, Reproducing topological properties with quasi-Majorana states, SciPost Phys. 7, 61 (2019).
  * Awoga  _et al._ [2019] O. A. Awoga, J. Cayao, and A. M. Black-Schaffer, Supercurrent detection of topologically trivial zero-energy states in nanowire junctions, Phys. Rev. Lett. 123, 117001 (2019).
  * Pan and Das Sarma [2020] H. Pan and S. Das Sarma, Physical mechanisms for zero-bias conductance peaks in Majorana nanowires, Phys. Rev. Research 2, 013377 (2020).
  * Avila  _et al._ [2019] J. Avila, F. Peñaranda, E. Prada, P. San-Jose, and R. Aguado, Non-hermitian topology as a unifying framework for the Andreev versus Majorana states controversy, Communications Physics 2, 133 (2019).
  * Cayao and Burset [2021] J. Cayao and P. Burset, Confinement-induced zero-bias peaks in conventional superconductor hybrids, Phys. Rev. B 104, 134507 (2021).
  * Rosdahl  _et al._ [2018] T. O. Rosdahl, A. Vuik, M. Kjaergaard, and A. R. Akhmerov, Andreev rectifier: A nonlocal conductance signature of topological phase transitions, Phys. Rev. B 97, 045421 (2018).
  * Danon  _et al._ [2020] J. Danon, A. B. Hellenes, E. B. Hansen, L. Casparis, A. P. Higginbotham, and K. Flensberg, Nonlocal conductance spectroscopy of Andreev bound states: Symmetry relations and bcs charges, Phys. Rev. Lett. 124, 036801 (2020).
  * Hess  _et al._ [2021] R. Hess, H. F. Legg, D. Loss, and J. Klinovaja, Local and nonlocal quantum transport due to Andreev bound states in finite rashba nanowires with superconducting and normal sections, Phys. Rev. B 104, 075405 (2021).
  * Melo  _et al._ [2021] A. Melo, C.-X. Liu, P. Rożek, T. Örn Rosdahl, and M. Wimmer, Conductance asymmetries in mesoscopic superconducting devices due to finite bias, SciPost Phys. 10, 37 (2021).
  * Pikulin  _et al._ [2021] D. I. Pikulin, B. van Heck, T. Karzig, E. A. Martinez, B. Nijholt, T. Laeven, G. W. Winkler, J. D. Watson, S. Heedt, M. Temurhan, V. Svidenko, R. M. Lutchyn, M. Thomas, G. de Lange, L. Casparis, and C. Nayak, Protocol to identify a topological superconducting phase in a three-terminal device, arXiv:2103.12217 (2021).
  * Singh and Muralidharan [2022] R. Singh and B. Muralidharan, Conductance spectroscopy of Majorana zero modes in superconductor-magnetic insulating nanowire hybrid systems, arXiv:2203.08413 (2022).
  * Maiani  _et al._ [2022] A. Maiani, M. Geier, and K. Flensberg, Conductance-matrix symmetries of multiterminal semiconductor-superconductor devices, arXiv:2205.11193 (2022).
  * [32] A. Tsintzis, R. Seoane Souto, and M. Leijnse, Creating and detecting poor man’s Majorana bound states in interacting quantum dots, To be soon submitted .
  * Gramich  _et al._ [2017] J. Gramich, A. Baumgartner, and C. Schönenberger, Andreev bound states probed in three-terminal quantum dots, Phys. Rev. B 96, 195418 (2017).
  * 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, Conductance-matrix symmetries of a three-terminal hybrid device, Phys. Rev. Lett. 124, 036802 (2020).
  * 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, Closing of the induced gap in a hybrid superconductor-semiconductor nanowire, Phys. Rev. B 103, 235201 (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, Measurement circuit effects in three-terminal electrical transport measurements, arXiv:2104.02671 (2021).
  * Wang  _et al._ [2022] G. Wang, T. Dvir, N. van Loo, G. P. Mazur, S. Gazibegovic, G. Badawy, E. P. A. M. Bakkers, L. P. Kouwenhoven, and G. de Lange, Non-local measurement of quasiparticle distribution in proximitized semiconductor nanowires using quantum dots, arXiv:2110.05373 (2022).
  * Dvir  _et al._ [2022] T. Dvir, G. Wang, N. van Loo, C.-X. Liu, G. P. Mazur, A. Bordin, S. L. D. ten Haaf, J.-Y. Wang, D. van Driel, F. Zatelli, X. Li, F. K. Malinowski, S. Gazibegovic, G. Badawy, E. P. A. M. Bakkers, M. Wimmer, and L. P. Kouwenhoven, Realization of a minimal Kitaev chain in coupled quantum dots, arXiv:2206.08045 (2022).
  * Pöschl  _et al._ [2022a] A. Pöschl, A. Danilenko, D. Sabonis, K. Kristjuhan, T. Lindemann, C. Thomas, M. J. Manfra, and C. M. Marcus, Nonlocal signatures of hybridization between quantum dot and Andreev bound states, arXiv:2201.03687 (2022a).
  * Pöschl  _et al._ [2022b] A. Pöschl, A. Danilenko, D. Sabonis, K. Kristjuhan, T. Lindemann, C. Thomas, M. J. Manfra, and C. M. Marcus, Nonlocal conductance spectroscopy of Andreev bound states in gate-defined InAs/Al nanowires, arXiv:2204.02430 (2022b).
  * Banerjee  _et al._ [2022] A. Banerjee, O. Lesser, M. A. Rahman, C. Thomas, T. Wang, M. J. Manfra, E. Berg, Y. Oreg, A. Stern, and C. M. Marcus, Local and nonlocal transport spectroscopy in planar Josephson junctions, arXiv:2205.09419 (2022).
  * Aghaee  _et al._ [2022] M. Aghaee _et al._ , InAs-Al Hybrid Devices Passing the Topological Gap Protocol, arXiv:2207.02472 (2022).
  * van Heck  _et al._ [2012] B. van Heck, A. R. Akhmerov, F. Hassler, M. Burrello, and C. W. J. Beenakker, Coulomb-assisted braiding of Majorana fermions in a Josephson junction array, New Journal of Physics 14, 035019 (2012).
  * Terhal  _et al._ [2012] B. M. Terhal, F. Hassler, and D. P. DiVincenzo, From Majorana fermions to topological order, Phys. Rev. Lett. 108, 260504 (2012).
  * Hyart  _et al._ [2013] T. Hyart, B. van Heck, I. C. Fulga, M. Burrello, A. R. Akhmerov, and C. W. J. Beenakker, Flux-controlled quantum computation with Majorana fermions, Phys. Rev. B 88, 035121 (2013).
  * Plugge  _et al._ [2016] S. Plugge, L. A. Landau, E. Sela, A. Altland, K. Flensberg, and R. Egger, Roadmap to Majorana surface codes, Phys. Rev. B 94, 174514 (2016).
  * Aasen  _et al._ [2016] D. Aasen, M. Hell, R. V. Mishmash, A. Higginbotham, J. Danon, M. Leijnse, T. S. Jespersen, J. A. Folk, C. M. Marcus, K. Flensberg, and J. Alicea, Milestones Toward Majorana-Based Quantum Computing, Phys. Rev. X 6, 031016 (2016).
  * Landau  _et al._ [2016] L. A. Landau, S. Plugge, E. Sela, A. Altland, S. M. Albrecht, and R. Egger, Towards realistic implementations of a Majorana surface code, Phys. Rev. Lett. 116, 050501 (2016).
  * Plugge  _et al._ [2017] S. Plugge, A. Rasmussen, R. Egger, and K. Flensberg, Majorana box qubits, New Journal of Physics 19, 012001 (2017).
  * Karzig  _et al._ [2017] T. Karzig, C. Knapp, R. M. Lutchyn, P. Bonderson, M. B. Hastings, C. Nayak, J. Alicea, K. Flensberg, S. Plugge, Y. Oreg, C. M. Marcus, and M. H. Freedman, Scalable designs for quasiparticle-poisoning-protected topological quantum computation with Majorana zero modes, Phys. Rev. B 95, 235305 (2017).
  * Manousakis  _et al._ [2020] J. Manousakis, C. Wille, A. Altland, R. Egger, K. Flensberg, and F. Hassler, Weak measurement protocols for Majorana bound state identification, Phys. Rev. Lett. 124, 096801 (2020).
  * Nitsch  _et al._ [2022] M. Nitsch, R. Seoane Souto, and M. Leijnse, Interference and parity blockade in transport through a Majorana box, arXiv:2205.10002 10.48550/arXiv.2205.10002 (2022).
  * Souto and Leijnse [2022] R. S. Souto and M. Leijnse, Fusion rules in a Majorana single-charge transistor, SciPost Phys. 12, 161 (2022).
  * Higginbotham  _et al._ [2015] A. P. Higginbotham, S. M. Albrecht, G. Kiršanskas, W. Chang, F. Kuemmeth, P. Krogstrup, T. S. Jespersen, J. Nygård, K. Flensberg, and C. M. Marcus, Parity lifetime of bound states in a proximitized semiconductor nanowire, Nature Physics 11, 1017 (2015).
  * Albrecht  _et al._ [2017] S. M. Albrecht, E. B. Hansen, A. P. Higginbotham, F. Kuemmeth, T. S. Jespersen, J. Nygård, P. Krogstrup, J. Danon, K. Flensberg, and C. M. Marcus, Transport signatures of quasiparticle poisoning in a Majorana island, Phys. Rev. Lett. 118, 137701 (2017).
  * Albrecht  _et al._ [2016] S. M. Albrecht, A. P. Higginbotham, M. Madsen, F. Kuemmeth, T. S. Jespersen, J. Nygård, P. Krogstrup, and C. M. Marcus, Exponential protection of zero modes in Majorana islands, Nature 531, 206 (2016).
  * O’Farrell  _et al._ [2018] E. C. T. O’Farrell, A. C. C. Drachmann, M. Hell, A. Fornieri, A. M. Whiticar, E. B. Hansen, S. Gronin, G. C. Gardner, C. Thomas, M. J. Manfra, K. Flensberg, C. M. Marcus, and F. Nichele, Hybridization of subgap states in one-dimensional superconductor-semiconductor Coulomb islands, Phys. Rev. Lett. 121, 256803 (2018).
  * Hansen  _et al._ [2018] E. B. Hansen, J. Danon, and K. Flensberg, Probing electron-hole components of subgap states in Coulomb blockaded Majorana islands, Phys. Rev. B 97, 041411 (2018).
  * Vaitiekėnas  _et al._ [2020] S. Vaitiekėnas, G. W. Winkler, B. van Heck, T. Karzig, M.-T. Deng, K. Flensberg, L. I. Glazman, C. Nayak, P. Krogstrup, R. M. Lutchyn, and C. M. Marcus, Flux-induced topological superconductivity in full-shell nanowires, Science 367, eaav3392 (2020).
  * Vaitiekėnas  _et al._ [2022] S. Vaitiekėnas, R. S. Souto, Y. Liu, P. Krogstrup, K. Flensberg, M. Leijnse, and C. M. Marcus, Evidence for spin-polarized bound states in semiconductor–superconductor–ferromagnetic-insulator islands, Phys. Rev. B 105, L041304 (2022).
  * Whiticar  _et al._ [2020] A. M. Whiticar, A. Fornieri, E. C. T. O’Farrell, A. C. C. Drachmann, T. Wang, C. Thomas, S. Gronin, R. Kallaher, G. C. Gardner, M. J. Manfra, C. M. Marcus, and F. Nichele, Coherent transport through a Majorana island in an Aharonov-Bohm interferometer, Nature Communications 11, 3212 (2020).
  * Razmadze  _et al._ [2019] D. Razmadze, D. Sabonis, F. K. Malinowski, G. C. Ménard, S. Pauka, H. Nguyen, D. M. van Zanten, E. C. O′Farrell, J. Suter, P. Krogstrup, F. Kuemmeth, and C. M. Marcus, Radio-frequency methods for Majorana-based quantum devices: Fast charge sensing and phase-diagram mapping, Phys. Rev. Applied 11, 064011 (2019).
  * Sabonis  _et al._ [2019] D. Sabonis, E. C. T. O’Farrell, D. Razmadze, D. M. T. van Zanten, J. Suter, P. Krogstrup, and C. M. Marcus, Dispersive sensing in hybrid InAs/Al nanowires, Applied Physics Letters 115, 102601 (2019).
  * Béri and Cooper [2012] B. Béri and N. R. Cooper, Topological Kondo effect with Majorana fermions, Phys. Rev. Lett. 109, 156803 (2012).
  * Altland  _et al._ [2014] A. Altland, B. Béri, R. Egger, and A. M. Tsvelik, Multichannel Kondo impurity dynamics in a Majorana device, Phys. Rev. Lett. 113, 076401 (2014).
  * Zazunov  _et al._ [2014] A. Zazunov, A. Altland, and R. Egger, Transport properties of the Coulomb–Majorana junction, New Journal of Physics 16, 015010 (2014).
  * Kanne  _et al._ [2022] T. Kanne, D. Olsteins, M. Marnauza, A. Vekris, J. C. Estrada Saldaña, S. Loric, R. D. Schlosser, D. Ross, S. Csonka, K. Grove-Rasmussen, and J. Nygård, Double nanowires for hybrid quantum devices, Advanced Functional Materials 32, 2107926 (2022).
  * Vekris  _et al._ [2022] A. Vekris, J. C. Estrada Saldaña, T. Kanne, T. Hvid-Olsen, M. Marnauza, D. Olsteins, M. M. Wauters, M. Burrello, J. Nygård, and K. Grove-Rasmussen, Electronic transport in double-nanowire superconducting islands with multiple terminals, Nano Letters 22, 5765 (2022).
  * Ekström  _et al._ [2020] J. Ekström, P. Recher, and T. L. Schmidt, Transport properties of coupled Majorana bound states in the Coulomb blockade regime, Phys. Rev. B 101, 195420 (2020).
  * Żak and Flensberg [2008] R. A. Żak and K. Flensberg, Coulomb blockade of a three-terminal quantum dot, Phys. Rev. B 77, 045329 (2008).
  * van Heck  _et al._ [2016] B. van Heck, R. M. Lutchyn, and L. I. Glazman, Conductance of a proximitized nanowire in the Coulomb blockade regime, Phys. Rev. B 93, 235431 (2016).
  * Lutchyn  _et al._ [2018] R. M. Lutchyn, E. P. A. M. Bakkers, L. P. Kouwenhoven, P. Krogstrup, C. M. Marcus, and Y. Oreg, Majorana zero modes in superconductor–semiconductor heterostructures, Nature Reviews Materials 3, 52 (2018).
  * Kürtössy  _et al._ [2021] O. Kürtössy, Z. Scherübl, G. Fülöp, I. E. Lukács, T. Kanne, J. Nygård, P. Makk, and S. Csonka, Andreev Molecule in Parallel InAs Nanowires, Nano Lett. 21, 7929 (2021).
  * Vekris  _et al._ [2021a] A. Vekris, J. C. E. Saldaña, T. Kanne, M. Marnauza, D. Olsteins, F. Fan, X. Li, T. Hvid-Olsen, X. Qiu, H. Xu, J. Nygård, and K. Grove-Rasmussen, Josephson junctions in double nanowires bridged by in-situ deposited superconductors, Phys. Rev. Research 3, 033240 (2021a).
  * Vekris  _et al._ [2021b] A. Vekris, J. C. Estrada Saldaña, J. de Bruijckere, S. Lorić, T. Kanne, M. Marnauza, D. Olsteins, J. Nygård, and K. Grove-Rasmussen, Asymmetric Little–Parks oscillations in full shell double nanowires - Scientific Reports, Sci. Rep. 11, 1 (2021b).
  * Schmitteckert [2004] P. Schmitteckert, Nonequilibrium electron transport using the density matrix renormalization group method, Phys. Rev. B 70, 121302 (2004).
  * Seoane Souto  _et al._ [2021] R. Seoane Souto, A. E. Feiguin, A. Martín-Rodero, and A. L. Yeyati, Transient dynamics of a magnetic impurity coupled to superconducting electrodes: Exact numerics versus perturbation theory, Phys. Rev. B 104, 214506 (2021).
  * Chung  _et al._ [2022] C.-M. Chung, M. M. Wauters, and M. Burrello, Matrix product state simulations of quantum quenches and transport in Coulomb blockaded superconducting devices (2022).
  * Koller  _et al._ [2010] S. Koller, M. Grifoni, M. Leijnse, and M. R. Wegewijs, Density-operator approaches to transport through interacting quantum dots: Simplifications in fourth-order perturbation theory, Phys. Rev. B 82, 235307 (2010).
  * Koch  _et al._ [2004] J. Koch, F. von Oppen, Y. Oreg, and E. Sela, Thermopower of single-molecule devices, Phys. Rev. B 70, 195107 (2004).
  * Ozaki [2007] T. Ozaki, Continued fraction representation of the Fermi-Dirac function for large-scale electronic structure calculations, Phys. Rev. B 75, 035123 (2007).