# In Lieu of Privacy: Anonymous Contact Tracing

Rohit Bhat, Shranav Palakurthi, and Naman Tiwari
The Johns Hopkins Institute of Security Informatics

###### Abstract

We present Tracer Tokens, a hardware token of privacy-preserving contact
tracing utilizing Exposure Notification [4] protocol. Through subnetworks, we
show that any disease spread by proximity can be traced such as seasonal flu,
cold, regional strains of COVID-19, or Tuberculosis. Further, we show this
protocol to notify $n^{n}$ users in parallel, providing a speed of information
unmatched by current contact tracing methods.

Keywords: Contact Tracing, Perfect Forward Secrecy, Google-Apple Exposure
Notification

##  1 Let’s Call It An Introduction

The ongoing global pandemic of SARS-CoV-2 has shown reliance on the
overburdened American systems of healthcare and public health. We see the
stress reveal existing problems in our public health infrastructure that
compromises user privacy and sensitive information. Through this paper, we
seek to address the problem of social stigma in contact tracing. We will see
that this anonymous protocol can be generalized to any airborne pathogen that
is spread within an acceptable range of Bluetooth Low Energy (BLE) protocol.
Observing the history of public behavior with the HIV virus, we hope this
problem can be addressed through transparent anonymity and behavioral
economics. By removing the negative stimulus of interpersonal contact tracing,
we hope to encourage the prosocial behavior of reporting a positive viral
test.

While contact tracing has been a critical foundation of containment efforts,
we have seen many concerns around privacy of information storage ([13], [1]).
These concerns are only heightened upon inspection of our data collection
systems. The state of Maryland has implemented CovidLINK [7], a database for
contact tracing information such as phone numbers or levels of exposure. While
the Privacy Policy had originally linked a static PDF to be downloaded for
future reference, it has since been changed to the server-side Notice of
Privacy Practices on the MDH website [8]. This provides very limited
visibility into the growing collection and transmission of personal
information. The pdf is currently found under the phpa.health.maryland.gov
subdomain[5] via an internet search. The authors have provided a copy that was
downloaded in September 2020 [6]. We find their SHA-256 hash to be
$008a447af0bcc981d205bbdaff3b99354553046431ae269ab87385c5e1107b08$.

After confirming the Privacy Policy remains unchanged, we can make
observations of the logic that is permissible under this Policy. As security-
minded researchers, we will assume our data to be insecure unless we can show
it to be safe. Reading page 5 of the Privacy Policy, we see:

_”…In many cases, this [contact tracing] does not require sharing personally
identifiable data at all. In other cases, personally identifiable data (like
your phone number) may be required to deliver the Services to you.”_

Continuing on the same page and onto the next, we will look at the third
paragraph under ’Our Partners and How We Share or Disclose Your Information:’

_”When your data is collected on the Services, it may be shared with selected
third parties who assist us with our business operations…without limitation,
companies that provide data storage, support customer service, and facilitate
or deliver materials to you via e-mail, other electronic communications
platforms, or postal service…These other sites and services are not bound by
our Privacy Policy, and we are not responsible for their information
collection practices.” _

Observing this system, it is clear why an individual would feel they are
furthering their personal risk by informing a contact tracer - there is a
clear path for very sensitive data that must be handled with perfect forward
secrecy. Since beginning this research, we have since seen a data breach of
this very nature affect 72,000 residents of Pennsylvania ([2], [9]).

We look to the Google/Apple Exposure Notification (GAEN) Protocol [4].
Evaluating the system from a theory of information, we can see the protocol
relies on purely pseudorandom numbers and a separate channel of information
(pre-determined shared knowledge of possible exposure events) to inform the
individual of a possible exposure. We find this to be suitable from a
theoretic perspective, assuming appropriate pseudorandomness.

Unfortunately, while the protocol exists for the anonymous transfer of contact
tracing information, we see low adoption rates because the data is too
sensitive to risk against the number of attack surfaces presented by our cell
phones [14].

We present Tracer Tokens, a privacy-preserving contact tracing network that
seeks to retain provable forward secrecy. We remove our system of information
from the cell phone, and put it into a hardware token. While similar efforts
have been introduced through government mandate in Singapore [15], we seek an
open source protocol that can be accepted within the unique culture of
American individualism.

##  2 Behavior Section

###  2.1 Speed of Notification

Contact tracing is currently based on phone calls and careful tracing through
a graph of meaningful interactions. This requires a conversation that can take
anywhere from 5 to 20 minutes per person. Using a Tracer Token network,
notifications of possible exposures can occur in parallel at the speed of
network propagation and hash computation. This is orders of magnitude faster -
thus minimizing the risk of an asymptomatic carrier spreading disease.

Notification of the exposure is the Diagnosis Keys, a set of 14 $tek$ that
correspond to the most recent 14 days. These are sent from a single user to a
server network, which propagates through the server network. Servers then send
the Diagnosis Keys to up to $n$ Tokens. Each Tracer token will locally
calculate the $1440$ hash values for each $tek$, and comparing up to $20160$
hashes with their own list of collected hash values in the same time period.

So - for $n$ users/Tokens, a server network needs to distribute up to $n$ dk-
sets. The server network needs a throughput of only $n$. So $n$ dk-sets being
hashed on $n$ devices means a contact tracing network capable of size $n^{n}$
requires $n$ throughput.

###  2.2 Privacy

Given a notification from the token itself, it is impossible for the incident
of exposure to be known. This is important for participation by the end-user.
The individual will act according to their own interests, and a token
notification will simply mean they had been within infectious distance of
somebody that is reporting a diagnosis to the network.

Since each $tek$ creates a new hash every 10 minutes, it is possible to gain a
general understanding of when an exposure took place. However, sharing this
information is unwarranted - the user can keep such information to themselves,
and request a test for disease without any reporting of private information.

###  2.3 Trustless

By removing the human element of contact tracing, a Token holder can be
confident their anonymity will be preserved. No user information is collected,
or recorded. This is easily shown because no registration is required. By
itself, a token notification is meaningless - it could have come from
anywhere. Once discarded, even a forensic analyst would find the data
meaningless without information of location history. This is only known to the
holder of the individual Token.

To demonstrate this trustless network, Tokens would be just as meaningful if
swapped by two individuals. While this would ’reset’ the timer of meaningful
notification, it is clear that the utility of an exposure notification remains
useful without any further knowledge gained.

##  3 Device

The Tracer Token is designed from the ground up to be low-power, low-cost, and
extremely simple to manufacture and distribute. The Tracer project drew
inspiration from the commercially successful and widely deployed Tile
Bluetooth tracker. While Tile is functionally different from a Tracer Token,
the design constraints are similar. Both are low-power, portable, and have
Bluetooth Low Energy capabilities. However, while Tile’s intended purpose is
to be a beacon, the Tracers must act independently, scanning for peers and
transferring data between each other. However, most of this differentiation
happens in software. The Token hardware is relatively simple: a BLE
transceiver, microcontroller, battery, and supporting circuitry. Because the
design is straightforward and utilizes easily sourceable parts, the Token is
perfect for cheap and efficient mass production. The heart of the Token is the
Bluetooth transceiver and microcontroller. At the time of writing, there are
many microcontrollers which have integrated 2.4GHz modems that offer BLE
capabilities. These SoCs allow for a cheap system that is both low-power and
easy to program. The prototype uses the ESP32 chip from Espressif, which was
the only chip the author had access to when writing the initial code [12].
However, since the summer of 2020, newer chips have been released which offer
more effective solutions. Espressif released the ESP32-C3, which offers
massively increased battery life (a 70 decrease in sleep power consumption),
at a minor performance loss (one less CPU core). Additionally, the new
ESP32-C3 is cheaper, with assembled modules priced around $1.80, as opposed to
the ESP32’s module cost of approximately $2.50. We plan to use the ESP32-C3
microcontroller for the final product, resulting in a projected battery life
more than double that of the ESP32 prototype. We estimate that the
ESP32-C3-based system draws about 100mAh for every 24 hours of usage. However,
the power draw is approximately inversely proportional to the transmission
interval. Doubling the transmission interval from 5 seconds to 10 seconds
halves the power consumption to 50mAh per 24 hours. We plan on using a lithium
ion battery cell to power the system, enabling users to recharge and reuse
their Tokens. This makes the system easier to use and cheaper to operate on
the user’s end. A 500mAh cell costs around $2, which can power a Token for
approximately 5 days. This power source in combination with its supporting
circuitry which includes a battery charger, LEDs, and a buck/boost converter
to power the electronics, adds around $3 to the Token’s total cost. With all
hardware factors taken into consideration, the final production cost is around
$5 per Token. While we acknowledge this is a rough estimate, it serves to
emphasize that the Tokens themselves are cheap and cost-effective, especially
for large institutions who will enjoy the benefits of economies of scale.

We have shown a proof of concept with two Arduinos and a web-based enrollment
and key server [12].

Tokens will have at least one button, meant for intialization or re-
initialization in the event of changing owners.

###  3.1 Bill of Materials

Using a 400mAh LiFePO4 Battery, we find the essential hardware to be
$5.68/unit [11]. We hope this cost can be reduced through efficiencies of
scale, as well as alternative hardware yet to be determined.

###  3.2 Subnetworks

Utilizing hash salts, a GAEN-based protocol can be split into subnetworks
unique to each airborne-disease. We point to the hash function $HKDF()$ given
by the Exposure Notifications Internals [3]:

⬇

KeyDerivation.hkdfSha256(

mac,

temporaryExposureKey,

/* inputSalt =*/ null,

aemkHkdfInfoBytes,

associatedMetadataEncryptionKeySizeBytes);

We can see the $inputSalt$ is defaulted to $null$. By properties of
deterministic hash functions, we can change the $inputSalt$ to any value and
generate a unique hash. This gives us the ability to subnetwork our contact
tracing. The given $hkdf$ is defined by IETF RFC 5869 [10], allowing for non-
secret random value. By properties of a deterministic hash function, this
effectively creates a subnetwork - only values with a matching $Salt$ will
share a codomain.

⬇

SHA-256($tek_i$|01|UTF8("EN-AEMK"), 16) =
377d7b4053a85dcb47d7a7adc97c749271383216822b44ac4e841291a92fcec1

SHA-256($tek_i$|02|UTF8("EN-AEMK), 16) =
ebf7e504e179fdad6a6701c91c5f57b738741483af560e985a88325a6926fff6"

For a given $tek_{i}$, we can produce multiple hashes that are linked to
different diseases.

###  3.3 Distribution of Diagnosis Keys

Diagnosis Keys are sets of 14 $tek$ that correspond to the most recent 14
days. Upon receiving a Diagnosis Key, the Tracer Token will iterate through
each $tek$ and check every hash against its existing list of collected BLE
Payloads. If there is a match, that means that the Tracer Token was in
proximity of the Token reporting the Diagnosis Key. This can be indicated by
an LED emitter. The owner will then have the knowledge to be tested for the
given disease listed on the Tracer Token.

This information can be used to begin a process of isolation, or report to the
nearest health authorities as deemed necessary for the disease being traced.

##  4 Milestones achieved

Tracer Tokens are intended to provide a low-tech solution for immediate
accessibility. The Exposure Notification protocol that was designed for the
COVID-19 pandemic shows great potential for increased privacy in a digital
world. Due to the decentralization of Tracer Tokens, we can create a contact
tracing network of size $n^{n}$ while keeping local computational complexity
in polynomial time.

##  5 Forseeable Engineering Challenges

###  5.1 Key Server

By removing the Exposure Notification protocol from the cell phone, we also
remove the communications to a Key Server that is implied through cellular
networks. As such, we have a new problem of transmitting and receiving a
positive diagnosis. While this can be solved via a centralized server, we
believe that solutions can be engineered with further expertise.

###  5.2 Security Against Malicious Actors

The system is designed to assume honest actors. However, it is inevitable that
some individuals may act maliciously by falsely reporting a positive
diagnosis. This is particularly dangerous because of the anonymity provided by
design.

We propose that upon a positive diagnosis, the healthcare provider generates a
Public/Private Key Pair. The Public Key can be added to a centralized key
server, which a can be checked against the digital signature of the Diagnosis
Key. However, this creates additional burden on the individual who was
recently informed of an illness.

Alternatively, the Healthcare Provider could be the trusted party in charge of
distributing Diagnosis Keys. When an individual goes to be tested, they turn
in their Tracer Token to the Healthcare Provider. The Healthcare Provider then
be responsible for reporting the Diagnosis Keys from the Tracer Token using
the Public/Private Keypair described above.

###  5.3 Conciliation of Hardware vs Theoretical Systems

First - we will need a method by which to input the $SALT$. This creates
additional complexities than the desired one-button device.

Second - reporting of Positive Diagnosis and distribution of the associated
Diagnosis Keys.

##  6 Conclusion

Tracer Tokens are intended to provide a low-tech solution for immediate
accessibility. The Exposure Notification protocol that was designed for the
COVID-19 pandemic shows great potential for increased privacy in a digital
world. Due to the decentralization of Tracer Tokens, we can create a contact
tracing network of size $n^{n}$ while keeping local computational complexity
in polynomial time.

Each Tracer Token is designed to detect other tokens within a 5-10 foot radius
of itself. This means the system is agnostic to the disease for which it is
performing contact tracing - a Token will notify the owner when a hashed
positive diagnosis is matched. A Token labeled with its specific disease can
be left for 1-6 months in a bag and notify the owner via LED if they have been
exposed to the listed disease. It is then the decision of the individual to be
tested.

Further, we observe that nobody can collect meaningful data from these Tokens
without their greater context - they area easily swappable, so the value of an
exposure notification is only as long as an individual has been holding it.
Any data from before it is in the individual’s possession is meaningless.
Value of any data is also reset by throwing it in the trash.

Additionally, by adding a $SALT$, we are able to create multiple sub-networks
for each disease to trace. This allows network capabilities of tracing
regional strains of COVID-19, Tuberculosis, common colds, or seasonal flu on
the same $tek_{i}$.

## References

  * [1] Ronald Bayer Amy Lauren Fairchild  “Contact tracing’s long, turbulent history holds lessons for COVID-19”  URL: https://news.osu.edu/contact-tracings-long-turbulent-history-holds-lessons-for-covid-19/
  * [2] Insight Global  “Notice of Data Event Related to Pennsylvania Contact Tracing”  URL: https://web.archive.org/web/20210811211916/https:/insightglobal.com/notice-of-data-event
  * [3] Google  “AssociatedEncryptedMetadataHelper.java”  URL: https://github.com/google/exposure-notifications-internals/blob/aa75f5c834aacdcad2aa29d899ba882295b31d16/exposurenotification/src/main/java/com/google/samples/exposurenotification/data/generator/AssociatedEncryptedMetadataHelper.java
  * [4] Google/Apple  “Privacy Preserving Contact Tracing”
  * [5] Maryland Department Health  “COVID Link Privacy Policy”  URL: https://health.maryland.gov/phpa/Documents/COVID%5C%20Link%5C%20Privacy%5C%20Policy%5C%20(for%5C%20Short%5C%20Code)%5C%20-%5C%20FINAL.pdf
  * [6] Maryland Department Health  “COVID Link Privacy Policy”  URL: https://drive.google.com/file/d/1P43GDDABxSLFQ85Ex7gAEAf004-18atU/view?usp=sharing
  * [7] Maryland Department Health  “covidLINK”  URL: https://covidlink.maryland.gov/content/answer-the-call
  * [8] Maryland Department Health  “Maryland Department of Health and Your Health Information - Notice of Privacy Practices”  URL: https://health.maryland.gov/pages/privacy.aspx
  * [9] Jamie Martines  “Personal data from Pa. contact-tracing calls still online despite assurances it had been secured”  URL: https://www.inquirer.com/news/pennsylvania/spl/pa-contact-tracing-data-breach-compromised-insight-global-20210609.html
  * [10] NIST  “HMAC-based Extract-and-Expand Key Derivation Function (HKDF)”  URL: https://datatracker.ietf.org/doc/html/rfc5869
  * [11] Shranav Palakurthi  “Bill of Materials - Project Tracer”  URL: https://drive.google.com/file/d/10EaLi-9iit6kRzht1vCiOnCGFFKuMuvK/view?usp=sharing
  * [12] Shranav Palakurthi  “Project Tracer - Confidential Contact Tracing for the Masses!”  URL: https://create.arduino.cc/projecthub/epicface2304/project-tracer-confidential-contact-tracing-for-the-masses-a6e2dc
  * [13] Jessica Rich  “How our ourdated privacy laws doomed contact-tracing apps”  URL: https://www.brookings.edu/blog/techtank/2021/01/28/how-our-outdated-privacy-laws-doomed-contact-tracing-apps
  * [14] Matt Richtel  “Contact tracing could be much easier - but there are trade-offs”  URL: https://www.baltimoresun.com/coronavirus/sns-nyt-coronavirus-contact-tracing-apps-20200604-vpr5r7n5lbhy3i36pnlsu3u42q-story.html
  * [15] Blue Trace Singapore Government Digital Services  “Contact tracing could be much easier - but there are trade-offs”  URL: https://www.tracetogether.gov.sg/common/token/


# Minimal Kitaev–transmon qubit based on double quantum dots

D. Michel Pino  Instituto de Ciencia de Materiales de Madrid (ICMM), Consejo
Superior de Investigaciones Científicas (CSIC), Sor Juana Inés de la Cruz 3,
28049 Madrid, Spain  Rubén Seoane Souto  Instituto de Ciencia de Materiales de
Madrid (ICMM), Consejo Superior de Investigaciones Científicas (CSIC), Sor
Juana Inés de la Cruz 3, 28049 Madrid, Spain  Ramón Aguado
<EMAIL_ADDRESS>Instituto de Ciencia de Materiales de Madrid (ICMM),
Consejo Superior de Investigaciones Científicas (CSIC), Sor Juana Inés de la
Cruz 3, 28049 Madrid, Spain

###### Abstract

Minimal Kitaev chains composed of two semiconducting quantum dots coupled via
a grounded superconductor have emerged as a promising platform to realize and
study Majorana bound states (MBSs). We propose a hybrid qubit based on a
Josephson junction between two such double quantum dots (DQDs) embedded in a
superconducting qubit geometry. The qubit makes use of the $4\pi$-Josephson
effect in the Kitaev junction to create a subspace based on the even/odd
fermionic parities of the two DQD arrays hosting MBSs. Deep in the transmon
regime, we demonstrate that by performing circuit QED spectroscopy on such
hybrid Kitaev-Transmon "Kitmon" qubit one could observe distinct MBS features
in perfect agreement with precise analytical predictions in terms of DQD
parameters only. This agreement allows to extract the Majorana polarization in
the junction from the microwave response.

_Introduction_ – Majorana bound states (MBSs) appearing at the ends of one-
dimensional topological superconductors Leijnse_Review2012 ; Alicea_RPP2012 ;
beenakker2013search ; Aguado_Nuovo2017 ; BeenakkerReview_20 ;
flensberg2021engineered ; Marra_Review2022  feature non-abelian statistics
that can be exploited for robust quantum information processing Nayak_review .
Although early experiments showed signatures consistent with their presence,
other states mainly originated from disorder can mimic their behavior, making
it hard to distinguish between trivial and topological states Prada-Review .

Artificial Kitaev chains circumvent the inherent disorder issues that commonly
appear in other platforms. In their minimal version, two quantum dots (QDs)
couple via a narrow superconductor that allows for crossed Andreev reflection
(CAR) and single-electron elastic co–tunneling (ECT) Leijnse ; Sau_NatComm2012
; PhysRevLett.129.267701 ; PhysRevB.106.L201404 ; Souto_arXiv2023 ;
Bordin_PRX2023 . Minimal Kitaev chains can host localized MBSs when a so-
called sweet spot is reached with equal CAR and ECT amplitudes. Although the
states are not topologically protected, they share properties with their
topological counterparts, including non-abelian statistics tsintzis2023roadmap
; Boross_2023 . Recent experiments have shown measurements consistent with
predictions at the sweet spot regime Dvir-Nature2023 , breaking a new ground
for the investigation of MBSs and paving the way towards scaling a
topologically-protected long chain and Majorana qubits 10.1063/PT.3.4499  with
QDs.

Figure 1: Schematic illustration of the Kitaev-Transmon device. A
semiconductor (pink) can be gated (yellow) to create two minimal Kitaev chains
(labeled as $\alpha=L,R$) comprising two quantum dots (labeled as
$\beta=1,2$), connected via a middle superconductor (blue) and with chemical
potentials $\mu_{E}$ and $\mu_{I}$, external and internal, respectively. Each
quantum dot contains two Majorana states $\gamma_{\alpha,\beta}^{A}$ and
$\gamma_{\alpha,\beta}^{B}$. The two Kitaev chains are connected through a
weak link (hopping $t_{J}$, purple region) forming a minimal Majorana
Josephson junction. This minimal Kitaev junction is connected to a transmon
circuit, where the island, with charging energy $E_{C}$, is connected to
ground by a SQUID formed by the parallel combination of the Kitaev junction
and a reference Josephson junction $E_{J}$. The superconducting phase
difference $\phi$ across the Kitaev junction is fixed by an externally applied
magnetic flux $\Phi_{ext}$ applied through the SQUID loop.

Expanding on this idea, we here propose a qubit based on a minimal Kitaev
Josephson junction with four QDs and embedded in a superconducting qubit
geometry, Fig. 1. The Josephson potential of the QD array modifies the
superconducting qubit Hamiltonian and splits the microwave (MW) transitions
owing to the (nearly) degenerate fermionic parities of the Kitaev chains. Deep
in the transmon limit, the qubit frequency can be analytically written in
terms of QD parameters, Eq. (12), in perfect agreement with full numerics
(Fig. 4). This agreement allows to extract the Majorana polarization (MP) of
the QD chain, Eq. (10), a measure of the Majorana character of the ground
states wavefunction PhysRevB.106.L201404 ; Sedlmayr2015 ; Sedlmayr2016 ;
Aksenov2020 , from the microwave response.

_Model_ –The minimal realization of a DQD-based Kitaev chain can be written as

$H_{\mathrm{DQD}}=-\sum_{i}\mu_{i}c_{i}^{\dagger}c_{i}-tc_{1}^{\dagger}c_{2}+\Delta
c_{1}c_{2}+\mbox{H.c.}\,,$ (1)

where $c_{i}^{\dagger}$ ($c_{i}$) denote creation (annihilation) operators on
the $i\in 1,2$ quantum dot with a chemical potential $\mu_{i}$, while $t$ and
$\Delta$ are the coupling strengths mediated by CAR and ECT processes across a
middle superconducting segment, respectively 111For the sake of simplicity, we
assume in what follows that $t_{\alpha}$ and $\Delta_{\alpha}$ are parameters
of the model, but we note in passing that both can be obtained from a
microscopic description of the middle segments mediating the interdot
couplings PhysRevLett.129.267701 ; PhysRevB.106.L201404 . Using this idea, a
minimal Kitaev Josephson junction can be written as
$H^{JJ}_{\mathrm{DQD}}=H_{\mathrm{DQD}}^{L}+H_{\mathrm{DQD}}^{R}+H_{J}$, where
$H_{\mathrm{DQD}}^{L}$ and $H_{\mathrm{DQD}}^{R}$ are two left/right Kitaev
chains based on Eq. (1) and the Josephson coupling reads:

$H_{J}=-t_{J}e^{i\phi/2}c_{L,2}^{\dagger}c_{R,1}+\mbox{H.c.}\;,$ (2)

with $\phi=\phi_{R}-\phi_{L}$ being the superconducting phase difference and
$t_{J}$ the tunneling coupling between chains (see Fig. 1). The above model
can be written in Bogoliubov–de Gennes (BdG) form as
$H^{JJ}_{\mathrm{BdG}}=\frac{1}{2}\Psi^{\dagger}H^{JJ}_{\mathrm{DQD}}\Psi$, in
terms of an eight-Majorana Nambu spinor

$\Psi=\left(\begin{matrix}\gamma_{L,1}^{A}&\gamma_{L,1}^{B}&\gamma_{L,2}^{A}&\gamma_{L,2}^{B}&\gamma_{R,1}^{A}&\gamma_{R,1}^{B}&\gamma_{R,2}^{A}&\gamma_{R,2}^{B}\end{matrix}\right)^{T}\,.$
(3)

As we discuss below, the BdG model contains a standard Josephson coupling
$\sim\cos\phi$ involving the "bulk" fermions together with a Majorana-mediated
$4\pi$ Josephson effect of order $\sim\cos\frac{\phi}{2}$. The latter involves
coherent single-electron tunneling with a characteristic energy scale $E_{M}$.
From the perspective of circuit QED, previous papers have discussed how a
Majorana junction in a transmon circuit splits spectral lines corresponding to
different fermionic parities owing to $E_{M}\neq 0$ Ginossar ; Keselman ;
Yavilberg ; Li2018 ; Avila ; Avila2 ; Smith2020 ; Lupo2022 . In what follows,
we discuss this physics in the context of the DQD minimal Kitaev Josephson
junction and to analyse the novel aspects that arise when this promising new
platform is integrated into a superconducting circuit.

_Four Majoranas subspace_ –A convenient way of gaining physical intuition is
by projecting the above full model onto a low-energy subspace. The simplest
approach, widely used in previous literature PhysRevLett.108.257001 ;
PhysRevB.86.140504 ; PhysRevB.97.041415 ; Cayao2018 , is to use a subspace
spanned by just four MBSs: the two inner $\gamma_{L,2}^{B}$ and
$\gamma_{R,1}^{A}$, and the two external $\gamma_{L,1}^{A}$ and
$\gamma_{R,2}^{B}$. This results in an effective Josephson potential

$V^{JJ}_{\mathrm{DQD}}(\phi)=E_{M}\cos\frac{\phi}{2}\sigma_{x}+E_{M}^{S}\sin\frac{\phi}{2}\sigma_{y}+\lambda\sigma_{z},$
(4)

where $\sigma_{i}$ are Pauli matrices defined onto the pseudospin parity space
spanned by $|0\rangle\equiv|n_{L}=0,\,n_{R}=0\rangle$ and
$|1\rangle\equiv|n_{L}=1,\,n_{R}=1\rangle$, where $n_{L}=n_{L,1}+n_{L,2}$ and
$n_{R}=n_{R,1}+n_{R,2}$ are the fermion occupations in the left/right segments
of the junction. $E_{M}^{S}$ and $\lambda$ are due to additional inter and
intra Majorana couplings {$\gamma_{L,1}^{A}\leftrightarrow\gamma_{R,1}^{A}$,
$\gamma_{L,2}^{B}\leftrightarrow\gamma_{R,2}^{B}$} and
{$\gamma_{L,1}^{A}\leftrightarrow\gamma_{L,2}^{B}$,
$\gamma_{R,1}^{A}\leftrightarrow\gamma_{R,2}^{B}$}, respectively. In the
symmetric case $\mu_{L,1}=\mu_{R,2}=\mu_{E}$ and
$\mu_{L,2}=\mu_{R,1}=\mu_{I}$, $E_{M}^{S}=\lambda=0$, which gives

$V^{JJ}_{\mathrm{DQD}}(\phi)=\frac{t_{J}}{2}\left[1-\frac{\mu_{E}^{2}}{(t+\Delta)^{2}}\right]\cos\frac{\phi}{2}\sigma_{x}\,.$
(5)

While being able to capture some of the phenomenology, including the $E_{M}$
renormalization with the external gates, this four–Majorana projection has
important limitations. Most importantly, detuning the chemical potentials
$\mu_{E}$ and $\mu_{I}$ away from zero affects the localization of the MBSs
which acquire some weight in "bulk" sites removed from the projection (for
instance, a $\mu_{E}\neq 0$ induces weight of the order of $\sim{\mu_{E}\over
t}$ in the inner dots Leijnse ). This makes the four–Majorana projection
insufficient to describe the physics of the DQD junction (for a full
derivation of Eq. (5) and a detailed discussion about the limitations of this
projection, see Appendix I).

Figure 2: Majorana polarization and Majorana coupling. (a) $2E_{M}/t_{J}$ and
(b) $|\mathrm{MP}_{1}|$ as a function of $\mu_{E}$ and $\mu_{I}$.
$2E_{M}/t_{J}$, $|\mathrm{MP}_{1}|$ and $|\mathrm{MP}_{2}|$ as a function of
(c) $\mu_{E}$ with $\mu_{I}=0$ and (d) $\mu_{E}=\mu_{I}=\mu$ (blue and green
dotted lines in panel a, respectively). $\Delta/t=1$ for all panels. Figure 3:
MW spectroscopy in charging regime. Levels, parity texture
$\langle\tau_{z}\rangle$ and $S(\omega)$ from the solutions of Eq. (11)
against $n_{g}$ with $\mu_{E}=\mu_{I}=0$ for (a-f) $\Delta/t=0.5,\,1,\,1.5$
and $\phi_{ext}=0$; and (g-l) $\phi_{ext}=\pi/2,\pi,3\pi/2$ and $\Delta=t$
(from top to bottom). $E_{J}/E_{C}=1$ and $t_{J}/t=1$ in all panels.

_Beyond four Majoranas_ – To go beyond the previous projection and its
limitations, we choose the subspace spanned by the two lowest–energy many–body
eigenstates $\\{|O_{L}^{-},O_{R}^{-}\rangle,\,|E_{L}^{-},E_{R}^{-}\rangle\\}$
resulting from diagonalizing each isolated segment in the basis of occupation
states $\\{\ket{10},\,\ket{01},\,\ket{00},\,\ket{11}\\}$. The diagonal
Hamiltonian in the bipartite Hilbert space
$\mathcal{H}_{L}\otimes\mathcal{H}_{R}$ can be represented on the basis of
joint eigenstates
$\\{|i_{L},j_{R}\rangle=|i_{L}\rangle\otimes|j_{R}\rangle\\}$ with
$i,j=O^{\pm},E^{\pm}$ (see Appendix II):

$\tilde{H}_{L}+\tilde{H}_{R}=(P_{L}^{-1}H_{L}P_{L})\otimes\mathbb{I}_{R}+\mathbb{I}_{L}\otimes(P_{R}^{-1}H_{R}P_{R}),$
(6)

where $P_{\alpha}$ is the change–of–basis matrix onto the eigenbasis of each
chain with eigenenergies $\epsilon_{\alpha
O}^{\pm}=-\mu_{\alpha}\pm\sqrt{t_{\alpha}^{2}+\delta_{\alpha}^{2}}$ and
$\epsilon_{\alpha
E}^{\pm}=-\mu_{\alpha}\pm\sqrt{\Delta_{\alpha}^{2}+\mu_{\alpha}^{2}}$, where
we have defined $\mu_{\alpha}=(\mu_{\alpha,1}+\mu_{\alpha,2})/2$ and
$\delta_{\alpha}=(\mu_{\alpha,1}-\mu_{\alpha,2})/2$. The off–diagonal
Josephson term $\tilde{H}_{J}$ can be easily represented on the
joint–occupation basis
$\\{|n_{L,1},n_{L,2}\rangle\otimes|n_{R,1},n_{R,2}\rangle\\}_{n_{\alpha,i}=0,1}$
and then projected onto the eigenbasis by the change–of–basis matrix
$P_{LR}=P_{L}\otimes P_{R}$. Using this projection, the Josephson potential
can be obtained analytically (see Appendix II). Specifically, for the
mirror–symmetric case, $\mu_{L,1}=\mu_{R,2}=\mu_{E}$ and
$\mu_{L,2}=\mu_{R,1}=\mu_{I}$ (external vs. internal), such that
$\mu_{L}=\mu_{R}=(\mu_{E}+\mu_{I})/2=\mu$ and
$\delta_{L}=-\delta_{R}=(\mu_{E}-\mu_{I})/2=\delta$, and considering
$\Delta_{L}=\Delta_{R}$ and $t_{L}=t_{R}$, this Josephson potential reduces to
a very compact form

$V^{JJ}_{\mathrm{DQD}}(\phi)=\left(\begin{matrix}-2\mu-2\sqrt{t^{2}+\delta^{2}}&E_{M}\cos\frac{\phi}{2}\\\
E_{M}\cos\frac{\phi}{2}&-2\mu-2\sqrt{\Delta^{2}+\mu^{2}}\end{matrix}\right)$
(7)

with

$E_{M}=\frac{t_{J}\Delta
t}{2\sqrt{(t^{2}+\delta^{2})(\Delta^{2}+\mu^{2})}}\;.$ (8)

The diagonal terms in Eq. (7) originate from the MBSs overlapping within the
same chain. Taking a series expansion up to leading order of $\mu$ and
$\delta$, Eq. (8) reduces to $E_{M}$ in Eq. (5) for $t=\Delta$ and
$\mu=\delta$ ($\mu_{I}=0$).

_Majorana polarization_ –For $t_{J}=0$, the many body problem described above
can be separated into two independent blocks of even
($\\{|O_{L}^{\pm},O_{R}^{\pm}\rangle$, $|E_{L}^{\pm},E_{R}^{\pm}\rangle\\}$)
and odd
($\\{|E_{L}^{\pm},O_{R}^{\pm}\rangle,\,|E_{L}^{\pm},O_{R}^{\pm}\rangle\\}$)
total parity, which leads to a two–fold degenerate spectrum. To determine
whether these degeneracies are associated with MBSs, we use the Majorana
polarization (MP) defined on the left Kitaev chain as
$\mathrm{MP}_{i}(O,E)=\frac{w_{i}^{2}-z_{i}^{2}}{w_{i}^{2}+z_{i}^{2}}$, with
$w_{i}={\left\langle O\right|c_{i}+c_{i}^{\dagger}\left|E\right\rangle}$,
$z_{i}={\left\langle O\right|c_{i}-c_{i}^{\dagger}\left|E\right\rangle}$ and
$i\in 1,2$. For the left DQD, we take $\ket{E}=|O_{L}^{-},O_{R}^{-}\rangle$,
and $\ket{O}=|E_{L}^{-},O_{R}^{-}\rangle$, this gives

$\mathrm{MP}_{1/2}=\frac{t\Delta}{\pm\delta\mu-\sqrt{(t^{2}+\delta^{2})(\Delta^{2}+\mu^{2})}}\;,$
(9)

where we have omitted the left subscript for simplicity. A similar treatment
can be performed for the right chain.

For $t=\Delta$, $|\mathrm{MP}_{1}|$ ($|\mathrm{MP}_{2}|$) is maximum when
$\mu=\delta$ ($\mu=-\delta$), that is, when $\mu_{L,2}=0$ ($\mu_{L,1}=0$).
Interestingly, by comparison with Eq. (8), when $\mu_{L,1}=\mu_{R,2}=\mu_{E}$
and $\mu_{L,2}=\mu_{R,1}=\mu_{I}$ ($\mu_{L}=\mu_{R}=\mu$ and
$\delta_{L}=-\delta_{R}=\delta$), one can write:

$\mathrm{MP}_{I/E}=\frac{-E_{M}}{\frac{t_{J}}{2}\pm\frac{\delta\mu}{t\Delta}E_{M}}\;.$
(10)

Note that for $\delta=0$ or $\mu=0$ ($\mu_{E}=\mu_{I}$ or $\mu_{E}=-\mu_{I}$,
respectively), $\mathrm{MP}$ is equal on every QD and it is directly
proportional to $E_{M}$. Therefore, Eq. (10) directly relates the MP with
$E_{M}$, which allows its direct measurement via MW spectroscopy as we discuss
now.

_Hybrid superconducting qubit model_ –We now study a DQD-based Majorana
Josephson junction in a superconducting qubit geometry (namely a split
junction shunted by a capacitor, with charging energy $E_{C}$, see Fig. 1)
described by the Hamiltonian:

$H=4E_{C}(\hat{n}-n_{g})^{2}-E_{J}cos(\hat{\phi})+V^{JJ}_{\mathrm{DQD}}(\hat{\phi}-\phi_{ext})\;.$
(11)

Here, $\hat{n}=-i\frac{\partial}{\partial\hat{\phi}}$ is the Cooper-pair
number operator, conjugate to the junction superconducting phase difference
$\hat{\phi}$, and $n_{g}=Q_{g}/(2e)=V_{g}/(2eC_{g})$ the gate–induced offset
charge in the island (in units of Cooper pairs). The phase difference across
the DQD Josephson junction can be controlled by the magnetic flux through the
SQUID loop $\Phi_{ext}=\phi_{ext}\Phi_{0}/(2\pi)$, where $\Phi_{0}=h/2e$ is
the superconducting flux quantum. Using the solutions of (11) 222 In practice,
we solve the model as a tight–binding chain in charge space. Specifically, we
divide the phase interval $\phi\in[0,2\pi)$ in $N$ steps, constructing a
$2N\times 2N$ Hamiltonian matrix in tight–binding form. Then, we can move to
its dual space of charge states $\\{\ket{n},\ket{n+1/2}\\}_{n=-N}^{N}$ and
rewrite the tight–binding Hamiltonian in this basis by applying a Fourier
transformation of the quantum phase operators (see Appendix V)., the microwave
(MW) absorption spectrum333For graphical purposes, we have convolved this
quantity with a Cauchy–Lorentz distribution ($\gamma=0.008$), which yields a
finite–line broadening of the spectra. of the superconducting island can be
written in linear response as $S(\omega)=\sum_{k}\left|{\left\langle
k\right|\hat{n}\left|0\right\rangle}\right|^{2}\delta(\omega-\omega_{0k})\;$,
where the index $k$ orders the eigenstates of the system with increasing
energies. This response measures the energy transitions
$\omega_{0k}=\omega_{k}-\omega_{0}$ between the ground state
$E_{0}=\hbar\omega_{0}$ and the excited states $E_{k}=\hbar\omega_{k}$ and
with a probability weighted by the matrix elements of ${\hat{n}}$.

Single–electron tunneling processes mediated by the off–diagonal terms of the
DQD-based Josephson potential in Eq. (7) lead to very specific predictions in
the spectrum that should be easily detected using standard circuit QED
techniques. For example, crossing the sweet spot, while keeping
$\mu_{E}=\mu_{I}=0$, from the ECT-dominated regime ($t>\Delta$, Fig. 3a) to
the CAR-dominated regime ($t<\Delta$, Fig. 3c), changes the fermionic parity
of the GS. This is reflected as an _exact 1 $e$ shift in $n_{g}$ in the MW
spectra_ (compare Figs. 3d and f). At the sweet spot for $t=\Delta$, the
intraband coupling leads to maximally mixed parity states
$\langle\tau_{z}\rangle=0$ with avoided crossings around $n_{g}=0.25$ and
$n_{g}=0.75$, Fig. 3b. This results in an overall 1$e$-periodic MW spectrum
with a strong low-frequency response near these gates, Fig. 3e. Therefore, the
intraband transition $\omega_{01}$ is a direct indication of a $E_{M}\neq 0$
in the spectrum.

_Kitaev-Transmon regime_ –A way to check that the low-frequency MW transitions
$\omega_{01}$ near $n_{g}=0.25$ and $n_{g}=0.75$ are indeed due to parity
mixing mediated by MBSs in the DQD junction, instead of quasiparticle
poisoning SchreierPRB08 ; KringhojPRL20 ; BargerbosPRL20 , is to prove that
they can be tuned by $\phi_{ext}$, and reach a minimum value at
$\phi_{ext}=\pi$, Figs. 3j-l. Note however, that owing to quantum phase
fluctuations, the Josephson potential $V^{JJ}_{\mathrm{DQD}}$ in Eq. (11)
depends on a phase drop which deviates from the external phase imposed by the
circuit, hence resulting in a residual splitting at $n_{g}=0.25$ which does
not close completely at $\phi_{ext}=\pi$. This effect is shown in Fig. 4a,
where we plot the full $\phi_{ext}$ dependence corresponding to the MW spectra
of Figs. 3j-l at fixed $n_{g}=0.25$. Interestingly, parity changes due to
Majorana physics are already evident as a spectral hole near $\phi_{ext}=\pi$
in the transition $\omega_{02}$. By tracing such spectral hole in
$\omega_{02}$ (or, equivalently the appearance of the transition
$\omega_{03}$) we can identify when a true energy crossing occurs in the
system as a function of increasing $E_{J}/E_{C}$ ratios, Figs. 4b,c.

Figure 4: Kitaev-transmon qubit spectroscopy. (a) Full phase dependence of the
MW absorption spectrum of Fig. 3g-l at $n_{g}=0.25$. (b-c) Spectral weights
for transitions $\omega_{02}$ ($S_{2}$) and $\omega_{03}$ ($S_{3}$) as a
function of $\phi_{ext}$ and $E_{J}/E_{C}$ at the sweet spot ($\Delta=t$,
$\mu_{E}=\mu_{I}=0$). (d-g) MW absorption spectra as a function of (d)
$\phi_{ext}$ at the sweet spot; (e) $\mu_{E}$ with $\mu_{I}=0$ and $\Delta=t$
and $\phi_{ext}=0$; (f-g) $\mu_{E}=\mu_{I}=\mu$ with $\Delta=t$ and
$\phi_{ext}=0,\pi$; and (h-i) $\Delta/t$ with $\mu_{E}=\mu_{I}=0$ and
$\phi_{ext}=0,\pi$. Green dashed lines correspond to the analytical qubit
frequency $\omega_{KiT}$ in Eq. (12). For panels (d-g) we have fixed a ratio
$E_{J}/E_{C}=50$. $t_{J}/t=1$ for all panels.

While, generally, an analytical expression of the energy splitting at
$n_{g}=0.25$ would require knowing the explicit form of the qubit wave
functions, the deep transmon regime with $E_{J}/E_{C}\gg 1$ allows us to
approximate these eigenfunctions to two coupled (parity–defined)
harmonic–oscillator states sharpened around $\phi_{ext}$. In this regime, the
Kitmon qubit frequency $\omega_{KiT}\equiv\omega_{01}$ can be written as

$\omega_{KiT}\approx
2\sqrt{(\sqrt{t^{2}+\delta^{2}}-\sqrt{\Delta^{2}+\mu^{2}})^{2}+E_{M}^{2}\cos^{2}\frac{\phi_{ext}}{2}}$
(12)

(A detailed check of the validity of Eq. (12) for increasing values of
$E_{J}/E_{C}$ ratios can be found in Appendix IV). When $t=\Delta$ and
$\delta=\pm\mu$ ($\mu_{I}=0$ or $\mu_{E}=0$), the qubit frequency is directly
proportional to $E_{M}$,

$\omega_{KiT}\approx
2E_{M}\cos\frac{\phi_{ext}}{2}=\frac{t_{J}}{1+(\mu_{E}/\Delta)^{2}}\cos\frac{\phi_{ext}}{2}.$
(13)

A direct comparison between the full numerics and Eq. (12) against different
parameters of the junction, Figs. 4d-i, demonstrates an almost perfect
agreement. Therefore, MW measurements like the ones discussed here should
allow to check our predictions, e.g., the resonant behavior against $\mu_{E}$
in Eq. (13), see Fig. 4e. More importantly, a measurement like the one shown
in Figs. 4f and g (namely, $\omega_{KiT}$ versus $\mu=\mu_{E}=\mu_{I}$, hence
$\delta=0$) would allow to directly extract $E_{M}$ and hence determine the MP
polarization of the junction via Eq. (10).

In conclusion, we have proposed a minimal Kitaev-Transmon qubit based on a QD
Josephson junction array embedded in a superconducting circuit. Deep in the
transmon regime with $E_{J}/E_{C}\gg 1$ we have found an analytical expression
for the qubit frequency, Eq. (12), that allows to obtain very precise
predictions of its evolution against QD parameters, Fig. 4, and to extract the
Majorana polarization. The precise predictions in terms of analytics would
allow to experimentally distinguish the physics discussed here from either
quasiparticle poisoning or 4$\pi$ phase slips due to QD resonances Vakhtel23 .
This novel qubit architecture is a natural extension of the recent
experimental implementations of nanowire-based double island devices
Zanten_NatPhys2020 , gatemons gatemon1 ; gatemon2 ; Sabonis_PRL2020 ; Huo_2023
and Andreev spin qubits Hays2021 , although free from the uncertainties
originated from disorder. Most importantly, QD-based Josephson junctions
embedded in a transmon circuit have recently been implemented experimentally
KringhojPRL20 ; BargerbosPRL20 ; PRXQuantum.3.030311 . In the strong Coulomb
Blockade regime, they have been used to show spin-split MW transition lines
Bargerbos2022spectroscopy  forming a QD-based superconducting spin qubit
coherently coupled to a transmon Pita-Vidal-NaturePhys23 . In this context,
our DQD proposal could be seen as a minimal Majorana-based non-local parity
pseudospin, Eq. (7), coupled to a transmon. All this experimental progress,
together with the recent demonstration of poisoning times of the order of
milliseconds Hinderling_arXiv2023  and quasiparticle trapping engineering
Gerbold19 ; NguyenPRB23 ; uilhoorn2021quasiparticle , make the physics
discussed here within reach with superconductor-semiconductor hybrid devices
444Two-tone spectroscopy measurements used to detect the MW transitions
described here are typically integrated in time scales of the order of tens of
milliseconds, see e.g. BargerbosPRL20 .

###### Acknowledgements.

We acknowledge the support of the Spanish Ministry of Science through Grants
PID2021- 125343NB-I00 and TED2021-130292B-C43 funded by
MCIN/AEI/10.13039/501100011033, "ERDF A way of making Europe", the Spanish CM
“Talento Program” (project No. 2022-T1/IND-24070), and European Union
NextGenerationEU/PRTR. Support by the CSIC Interdisciplinary Thematic Platform
(PTI+) on Quantum Technologies (PTI-QTEP+) is also acknowledged.

## References

  * (1) M. Leijnse and K. Flensberg, “Introduction to topological superconductivity and majorana fermions,” Semicond. Sci. Technol., vol. 27, p. 124003, nov 2012.
  * (2) J. Alicea, “New directions in the pursuit of majorana fermions in solid state systems,” Rep. Prog. Phys., vol. 75, p. 076501, jun 2012.
  * (3) C. W. J. Beenakker, “Search for majorana fermions in superconductors,” Annu. Rev. Condens. Matter Phys., vol. 4, no. 1, pp. 113–136, 2013.
  * (4) R. Aguado, “Majorana quasiparticles in condensed matter,” Riv. Nuovo Cimento, vol. 40, p. 523, 2017.
  * (5) C. W. J. Beenakker, “Search for non-Abelian Majorana braiding statistics in superconductors,” SciPost Phys. Lect. Notes 15, 2020.
  * (6) K. Flensberg, F. von Oppen, and A. Stern, “Engineered platforms for topological superconductivity and majorana zero modes,” Nat. Rev. Mat., vol. 6, no. 10, pp. 944–958, 2021.
  * (7) P. Marra, “Majorana nanowires for topological quantum computation,” J. Appl. Phys., vol. 132, p. 231101, 12 2022.
  * (8) C. Nayak, S. H. Simon, A. Stern, M. Freedman, and S. Das Sarma, “Non-abelian anyons and topological quantum computation,” Rev. Mod. Phys., vol. 80, pp. 1083–1159, Sep 2008.
  * (9) 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, vol. 2, no. 10, pp. 575–594, 2020.
  * (10) M. Leijnse and K. Flensberg, “Parity qubits and poor man’s majorana bound states in double quantum dots,” Phys. Rev. B, vol. 86, p. 134528, 2012\.
  * (11) J. D. Sau and S. D. Sarma, “Realizing a robust practical majorana chain in a quantum-dot-superconductor linear array,” Nature Commun., vol. 3, p. 964, Jul 2012.
  * (12) C.-X. Liu, G. Wang, T. Dvir, and M. Wimmer, “Tunable superconducting coupling of quantum dots via andreev bound states in semiconductor-superconductor nanowires,” Phys. Rev. Lett., vol. 129, p. 267701, Dec 2022.
  * (13) A. Tsintzis, R. S. Souto, and M. Leijnse, “Creating and detecting poor man’s majorana bound states in interacting quantum dots,” Phys. Rev. B, vol. 106, p. L201404, Nov 2022.
  * (14) R. Seoane Souto, A. Tsintzis, M. Leijnse, and J. Danon, “Probing majorana localization in minimal kitaev chains through a quantum dot,” arXiv:2308.14751, 2023.
  * (15) A. Bordin, G. Wang, C.-X. Liu, S. L. D. ten Haaf, N. van Loo, G. P. Mazur, D. Xu, D. van Driel, F. Zatelli, S. Gazibegovic, G. Badawy, E. P. A. M. Bakkers, M. Wimmer, L. P. Kouwenhoven, and T. Dvir, “Tunable crossed andreev reflection and elastic cotunneling in hybrid nanowires,” Phys. Rev. X, vol. 13, p. 031031, Sep 2023.
  * (16) A. Tsintzis, R. Seoane Souto, K. Flensberg, J. Danon, and M. Leijnse, “Roadmap towards Majorana qubits and nonabelian physics in quantum dot-based minimal Kitaev chains,” arXiv:2306.16289, 2023.
  * (17) P. Boross and A. Pályi, “Braiding-based quantum control of a Majorana qubit built from quantum dots,” arXiv:2305.08464, 2023.
  * (18) 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,” Nature, vol. 614, no. 7948, pp. 445–450, 2023.
  * (19) R. Aguado and L. P. Kouwenhoven, “Majorana qubits for topological quantum computing,” Physics Today, vol. 73, pp. 44–50, 06 2020.
  * (20) N. Sedlmayr and C. Bena, “Visualizing majorana bound states in one and two dimensions using the generalized majorana polarization,” Phys. Rev. B, vol. 92, p. 115115, Sep 2015.
  * (21) N. Sedlmayr, J. M. Aguiar-Hualde, and C. Bena, “Majorana bound states in open quasi-one-dimensional and two-dimensional systems with transverse rashba coupling,” Phys. Rev. B, vol. 93, p. 155425, Apr 2016.
  * (22) S. V. Aksenov, A. O. Zlotnikov, and M. S. Shustin, “Strong coulomb interactions in the problem of majorana modes in a wire of the nontrivial topological class bdi,” Phys. Rev. B, vol. 101, p. 125431, Mar 2020.
  * (23) For the sake of simplicity, we assume in what follows that $t_{\alpha}$ and $\Delta_{\alpha}$ are parameters of the model, but we note in passing that both can be obtained from a microscopic description of the middle segments mediating the interdot couplings PhysRevLett.129.267701 ; PhysRevB.106.L201404 .
  * (24) E. Ginossar and E. Grosfeld, “Microwave transitions as a signature of coherent parity mixing effects in the majorana-transmon qubit,” Nature Communications, vol. 5, no. 1, p. 4772, 2014.
  * (25) A. Keselman, C. Murthy, B. van Heck, and B. Bauer, “Spectral response of Josephson junctions with low-energy quasiparticles,” SciPost Phys., vol. 7, p. 050, 2019.
  * (26) K. Yavilberg, E. Ginossar, and E. Grosfeld, “Fermion parity measurement and control in majorana circuit quantum electrodynamics,” Phys. Rev. B, vol. 92, p. 075143, Aug 2015.
  * (27) T. Li, W. A. Coish, M. Hell, K. Flensberg, and M. Leijnse, “Four-majorana qubit with charge readout: Dynamics and decoherence,” Phys. Rev. B, vol. 98, p. 205403, Nov 2018.
  * (28) J. Ávila, E. Prada, P. San-Jose, and R. Aguado, “Superconducting islands with topological josephson junctions based on semiconductor nanowires,” Phys. Rev. B, vol. 102, p. 094518, Sep 2020.
  * (29) J. Ávila, E. Prada, P. San-Jose, and R. Aguado, “Majorana oscillations and parity crossings in semiconductor nanowire-based transmon qubits,” Phys. Rev. Res., vol. 2, p. 033493, Sep 2020.
  * (30) T. B. Smith, M. C. Cassidy, D. J. Reilly, S. D. Bartlett, and A. L. Grimsmo, “Dispersive readout of majorana qubits,” PRX Quantum, vol. 1, p. 020313, Nov 2020.
  * (31) E. Lupo, E. Grosfeld, and E. Ginossar, “Implementation of single-qubit gates via parametric modulation in the majorana transmon,” PRX Quantum, vol. 3, p. 020340, May 2022.
  * (32) P. San-Jose, E. Prada, and R. Aguado, “ac josephson effect in finite-length nanowire junctions with majorana modes,” Phys. Rev. Lett., vol. 108, p. 257001, Jun 2012.
  * (33) D. I. Pikulin and Y. V. Nazarov, “Phenomenology and dynamics of a majorana josephson junction,” Phys. Rev. B, vol. 86, p. 140504, Oct 2012.
  * (34) M. Trif, O. Dmytruk, H. Bouchiat, R. Aguado, and P. Simon, “Dynamic current susceptibility as a probe of majorana bound states in nanowire-based josephson junctions,” Phys. Rev. B, vol. 97, p. 041415, Jan 2018.
  * (35) J. Cayao, A. M. Black-Schaffer, E. Prada, and R. Aguado, “Andreev spectrum and supercurrents in nanowire-based sns junctions containing majorana bound states,” Beilstein Journal of Nanotechnology, vol. 9, pp. 1339–1357, 2018\.
  * (36) In practice, we solve the model as a tight–binding chain in charge space. Specifically, we divide the phase interval $\phi\in[0,2\pi)$ in $N$ steps, constructing a $2N\times 2N$ Hamiltonian matrix in tight–binding form. Then, we can move to its dual space of charge states $\\{\ket{n},\ket{n+1/2}\\}_{n=-N}^{N}$ and rewrite the tight–binding Hamiltonian in this basis by applying a Fourier transformation of the quantum phase operators (see Appendix V).
  * (37) For graphical purposes, we have convolved this quantity with a Cauchy–Lorentz distribution ($\gamma=0.008$), which yields a finite–line broadening of the spectra.
  * (38) J. A. Schreier, A. A. Houck, J. Koch, D. I. Schuster, B. R. Johnson, J. M. Chow, J. M. Gambetta, J. Majer, L. Frunzio, M. H. Devoret, S. M. Girvin, and R. J. Schoelkopf, “Suppressing charge noise decoherence in superconducting charge qubits,” Phys. Rev. B, vol. 77, p. 180502, May 2008.
  * (39) A. Kringhøj, B. van Heck, T. W. Larsen, O. Erlandsson, D. Sabonis, P. Krogstrup, L. Casparis, K. D. Petersson, and C. M. Marcus, “Suppressed charge dispersion via resonant tunneling in a single-channel transmon,” Phys. Rev. Lett., vol. 124, p. 246803, Jun 2020.
  * (40) A. Bargerbos, W. Uilhoorn, C.-K. Yang, P. Krogstrup, L. P. Kouwenhoven, G. de Lange, B. van Heck, and A. Kou, “Observation of vanishing charge dispersion of a nearly open superconducting island,” Phys. Rev. Lett., vol. 124, p. 246802, Jun 2020.
  * (41) T. Vakhtel and B. van Heck, “Quantum phase slips in a resonant josephson junction,” Phys. Rev. B, vol. 107, p. 195405, May 2023.
  * (42) D. M. T. van Zanten, D. Sabonis, J. Suter, J. I. Väyrynen, T. Karzig, D. I. Pikulin, E. C. T. O’Farrell, D. Razmadze, K. D. Petersson, P. Krogstrup, and C. M. Marcus, “Photon-assisted tunnelling of zero modes in a majorana wire,” Nature Physics, vol. 16, pp. 663–668, Jun 2020.
  * (43) G. de Lange, B. van Heck, A. Bruno, D. J. van Woerkom, A. Geresdi, S. R. Plissard, E. P. A. M. Bakkers, A. R. Akhmerov, and L. DiCarlo, “Realization of microwave quantum circuits using hybrid superconducting-semiconducting nanowire josephson elements,” Phys. Rev. Lett., vol. 115, p. 127002, Sep 2015.
  * (44) T. W. Larsen, K. D. Petersson, F. Kuemmeth, T. S. Jespersen, P. Krogstrup, J. Nygård, and C. M. Marcus, “Semiconductor-nanowire-based superconducting qubit,” Phys. Rev. Lett., vol. 115, p. 127001, Sep 2015\.
  * (45) D. Sabonis, O. Erlandsson, A. Kringhøj, B. van Heck, T. W. Larsen, I. Petkovic, P. Krogstrup, K. D. Petersson, and C. M. Marcus, “Destructive little-parks effect in a full-shell nanowire-based transmon,” Phys. Rev. Lett., vol. 125, p. 156804, Oct 2020.
  * (46) J. Huo, Z. Xia, Z. Li, S. Zhang, Y. Wang, D. Pan, Q. Liu, Y. Liu, Z. Wang, Y. Gao, J. Zhao, T. Li, J. Ying, R. Shang, and H. Zhang, “Gatemon qubit based on a thin inas-al hybrid nanowire,” Chinese Physics Letters, vol. 40, p. 047302, mar 2023.
  * (47) 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, vol. 373, no. 6553, pp. 430–433, 2021.
  * (48) A. Bargerbos, M. Pita-Vidal, R. Žitko, J. Ávila, L. J. Splitthoff, L. Grünhaupt, J. J. Wesdorp, C. K. Andersen, Y. Liu, L. P. Kouwenhoven, R. Aguado, A. Kou, and B. van Heck, “Singlet-doublet transitions of a quantum dot josephson junction detected in a transmon circuit,” PRX Quantum, vol. 3, p. 030311, Jul 2022.
  * (49) A. Bargerbos, M. Pita-Vidal, R. Žitko, L. J. Splitthoff, L. Grünhaupt, J. J. Wesdorp, Y. Liu, L. P. Kouwenhoven, R. Aguado, C. K. Andersen, A. Kou, and B. van Heck, “Spectroscopy of spin-split andreev levels in a quantum dot with superconducting leads,” Phys. Rev. Lett., vol. 131, p. 097001, Aug 2023.
  * (50) M. Pita-Vidal, A. Bargerbos, R. Žitko, L. J. Splitthoff, L. Grünhaupt, J. J. Wesdorp, Y. Liu, L. P. Kouwenhoven, R. Aguado, B. van Heck, A. Kou, and C. K. Andersen, “Direct manipulation of a superconducting spin qubit strongly coupled to a transmon qubit,” Nature Physics, 2023.
  * (51) M. Hinderling, S. C. ten Kate, D. Z. Haxell, M. Coraiola, S. Paredes, F. K. E. Cheah, R. Schott, W. Wegscheider, D. Sabonis, and F. Nichele, “Flip-chip-based fast inductive parity readout of a planar superconducting island,” arXiv:2307.06718, 2023.
  * (52) G. C. Ménard, F. K. Malinowski, D. Puglia, D. I. Pikulin, T. Karzig, B. Bauer, P. Krogstrup, and C. M. Marcus, “Suppressing quasiparticle poisoning with a voltage-controlled filter,” Phys. Rev. B, vol. 100, p. 165307, Oct 2019.
  * (53) H. Q. Nguyen, D. Sabonis, D. Razmadze, E. T. Mannila, V. F. Maisi, D. M. T. van Zanten, E. C. T. O’Farrell, P. Krogstrup, F. Kuemmeth, J. P. Pekola, and C. M. Marcus, “Electrostatic control of quasiparticle poisoning in a hybrid semiconductor-superconductor island,” Phys. Rev. B, vol. 108, p. L041302, Jul 2023.
  * (54) W. Uilhoorn, J. G. Kroll, A. Bargerbos, S. D. Nabi, C.-K. Yang, P. Krogstrup, L. P. Kouwenhoven, A. Kou, and G. de Lange, “Quasiparticle trapping by orbital effect in a hybrid superconducting-semiconducting circuit,” arXiv.2105.11038, 2021.
  * (55) Two-tone spectroscopy measurements used to detect the MW transitions described here are typically integrated in time scales of the order of tens of milliseconds, see e.g. BargerbosPRL20 .

Supplemental Material

##  I Four Majoranas subspace

###  I.1 Effective low–energy projection

In order to derive a quantitative low–energy description of our system, we
project the full Hamiltonian $H^{JJ}_{\mathrm{DQD}}$ –Eqs. (1) and (2) of the
main text– onto the fermionic parity subspace that forms the superconducting
qubit. This procedure considers both standard Josephson events due to Cooper
pair tunneling, as well as anomalous Majorana–mediated events, where a single
electron is transferred across the junction. Hence, we can distinguish two
contributions of the Josephson potential,
$V_{J}=V_{J}^{\mathrm{bulk}}+V_{\mathrm{DQD}}^{JJ}$. The first one takes into
account the bulk contribution of the Bogoliubov–de Gennes (BdG) levels above
the gap to the ground state energy, which we just assume to be of the standard
form $V_{J}^{\mathrm{bulk}}(\phi)=-E_{J}\cos\phi$. The second contribution
corresponds to the subgap sector, and it can be expressed as the projection
onto a fermionic parity basis of an effective model of four Majorana
operators, $\gamma_{L,1}^{A},\gamma_{L,2}^{B}\in L$ and
$\gamma_{R,1}^{A},\gamma_{R,2}^{B}\in R$, corresponding to the end modes of
both chains. Its effective Hamiltonian takes the general BdG form

$H_{\gamma}=\frac{i}{2}\left(\begin{matrix}\gamma_{L,1}^{A}&\gamma_{L,2}^{B}&\gamma_{R,1}^{A}&\gamma_{R,2}^{B}\end{matrix}\right)\left(\begin{matrix}0&\lambda_{L1,L2}&\lambda_{L1,R1}&\lambda_{L1,R2}\\\
-\lambda_{L1,L2}&0&\lambda_{L2,R1}&\lambda_{L2,R2}\\\
-\lambda_{L1,R1}&-\lambda_{L2,R1}&0&\lambda_{R1,R2}\\\
-\lambda_{L1,R2}&-\lambda_{L2,R2}&-\lambda_{R1,R2}&0\end{matrix}\right)\left(\begin{matrix}\gamma_{L,1}^{A}\\\
\gamma_{L,2}^{B}\\\ \gamma_{R,1}^{A}\\\
\gamma_{R,2}^{B}\end{matrix}\right)\;.$ (S.1)

Our objective is now to relate $H^{JJ}_{\mathrm{DQD}}$ to this general
effective model of four Majoranas $H_{\gamma}$ to obtain an explicit
expression of its coefficients. Thus, we project the BdG form of the former,
$H^{JJ}_{\mathrm{BdG}}$ –Eqs. (1) and (2) of the main text using the Majorana
spinor in Eq. (3) of the main text–, onto the low–energy subspace of Majorana
operators. In order to do that, we define a basis of fermionic operators

$c_{\alpha}=\frac{1}{\sqrt{2}}(\gamma_{\alpha,1}^{A}+i\gamma_{\alpha,2}^{B})\quad,\quad
c_{\alpha}^{\dagger}=\frac{1}{\sqrt{2}}(\gamma_{\alpha,1}^{A}-i\gamma_{\alpha,2}^{B})\;,$
(S.2)

and we compute the matrix elements of the resolvent of
$H^{JJ}_{\mathrm{BdG}}$,

$G(\omega)=[(\omega+i\,\epsilon)\mathbb{I}-H^{JJ}_{\mathrm{BdG}}]^{-1}\;,\quad\epsilon\to
0^{+}\;,$ (S.3)

at $\omega=0$ on the
$\psi^{0}=(c_{L},c_{R},c_{L}^{\dagger},c_{R}^{\dagger})^{T}_{0}$ state basis.
The procedure is as follows: first of all, we calculate $G(\omega)$ by
inverting the matrix $(\omega+i\,\epsilon)\mathbb{I}-H^{JJ}_{\mathrm{BdG}}$
written on the state basis of the whole system,

$\Psi=\left(\begin{matrix}\gamma_{L,1}^{A}&\gamma_{L,1}^{B}&\gamma_{L,2}^{A}&\gamma_{L,2}^{B}&\gamma_{R,1}^{A}&\gamma_{R,1}^{B}&\gamma_{R,2}^{A}&\gamma_{R,2}^{B}&\end{matrix}\right)^{T}\quad,\quad\Psi^{\dagger}=\Psi^{T}\;.$
(S.4)

Then, we evaluate this resolvent matrix at $\omega=0$ and we project it onto
the $\psi^{0}$ basis, expressed in terms of $\Psi$ states as

$\displaystyle\left(\begin{matrix}1&0&0&0\end{matrix}\right)^{T}_{0}$
$\displaystyle\equiv\frac{1}{\sqrt{2}}\left(\begin{matrix}1&0&0&i&0&0&0&0\end{matrix}\right)^{T}\quad$
$\displaystyle,\quad\left(\begin{matrix}0&1&0&0\end{matrix}\right)^{T}_{0}$
$\displaystyle\equiv\frac{1}{\sqrt{2}}\left(\begin{matrix}0&0&0&0&1&0&0&i\end{matrix}\right)^{T}$
(S.5) $\displaystyle\left(\begin{matrix}0&0&1&0\end{matrix}\right)^{T}_{0}$
$\displaystyle\equiv\frac{1}{\sqrt{2}}\left(\begin{matrix}1&0&0&-i&0&0&0&0\end{matrix}\right)^{T}\quad$
$\displaystyle,\quad\left(\begin{matrix}0&0&0&1\end{matrix}\right)^{T}_{0}$
$\displaystyle\equiv\frac{1}{\sqrt{2}}\left(\begin{matrix}0&0&0&0&1&0&0&-i\end{matrix}\right)^{T}\;.$

This gives rise to a $4\times 4$ matrix

$(\mathcal{H}_{0}^{-1})_{ij}={\left\langle\psi^{0}_{i}\right|G(\omega=0)\left|\psi^{0}_{j}\right\rangle}\;,$
(S.6)

whose inverse

$H_{0}=\frac{1}{2}\sum_{i,j}\psi^{0\dagger}_{i}(\mathcal{H}_{0})_{ij}\psi^{0}_{j}\;,$
(S.7)

is the projection of $H^{JJ}_{\mathrm{DQD}}$ onto the subspace of low–energy
fermions. Finally, a simple change of basis
$\psi^{0}\to\psi^{\gamma}=(\gamma_{L,1}^{A},\gamma_{L,2}^{B},\gamma_{R,1}^{A},\gamma_{R,2}^{B})_{\gamma}^{T}$
will allow us to indentify this matrix $\mathcal{H}_{0}$ with the effective
sub–gap Hamiltonian (S.1). Indeed, writing the $\psi^{\gamma}$ basis states in
terms of $\psi^{0}$ components,

$\displaystyle\left(\begin{matrix}1&0&0&0\end{matrix}\right)^{T}_{\gamma}\equiv\frac{1}{\sqrt{2}}\left(\begin{matrix}1&0&1&0\end{matrix}\right)^{T}_{0}$
$\displaystyle,\quad\left(\begin{matrix}0&1&0&0\end{matrix}\right)^{T}_{\gamma}\equiv\frac{1}{\sqrt{2}}\left(\begin{matrix}-i&0&i&0\end{matrix}\right)^{T}_{0}\;,$
(S.8)
$\displaystyle\left(\begin{matrix}0&0&1&0\end{matrix}\right)^{T}_{\gamma}\equiv\frac{1}{\sqrt{2}}\left(\begin{matrix}0&1&0&1\end{matrix}\right)^{T}_{0}$
$\displaystyle,\quad\left(\begin{matrix}0&0&0&1\end{matrix}\right)^{T}_{\gamma}\equiv\frac{1}{\sqrt{2}}\left(\begin{matrix}0&-i&0&i\end{matrix}\right)^{T}_{0}\;,$

we can express the Hamiltonian $\mathcal{H}_{0}$ in this new basis as

$(\mathcal{H}_{\gamma})_{ij}=\langle\psi^{\gamma}_{i}|\mathcal{H}_{0}|\psi^{\gamma}_{j}\rangle\;,$
(S.9)

which yields

$\displaystyle
H_{\gamma}=\frac{1}{2}\sum_{ij}\psi^{\gamma\dagger}_{i}(\mathcal{H}_{\gamma})_{ij}\psi^{\gamma}_{j}=$
(S.10)
$\displaystyle\frac{i}{2}\psi^{\gamma\dagger}\left(\begin{matrix}0&\frac{\mu_{L,1}\mu_{L,2}-(t_{L}+\Delta_{L})(t_{L}-\Delta_{L})}{t_{L}+\Delta_{L}}&-\frac{t_{J}\mu_{L,1}\sin\frac{\phi}{2}}{t_{L}+\Delta_{L}}&-\frac{t_{J}\mu_{L,1}\mu_{R,2}\cos\frac{\phi}{2}}{(t_{L}+\Delta_{L})(t_{R}+\Delta_{R})}\\\
\frac{(t_{L}+\Delta_{L})(t_{L}-\Delta_{L})-\mu_{L,1}\mu_{L,2}}{t_{L}+\Delta_{L}}&0&t_{J}\cos\frac{\phi}{2}&-\frac{t_{J}\mu_{R,2}\sin\frac{\phi}{2}}{t_{R}+\Delta_{R}}\\\
\frac{t_{J}\mu_{L,1}\sin\frac{\phi}{2}}{t_{L}+\Delta_{L}}&-t_{J}\cos\frac{\phi}{2}&0&\frac{\mu_{R,1}\mu_{R,2}-(t_{R}+\Delta_{R})(t_{R}-\Delta_{R})}{t_{R}+\Delta_{R}}\\\
\frac{t_{J}\mu_{L,1}\mu_{R,2}\cos\frac{\phi}{2}}{(t_{L}+\Delta_{L})(t_{R}+\Delta_{R})}&\frac{t_{J}\mu_{R,2}\sin\frac{\phi}{2}}{t_{R}+\Delta_{R}}&\frac{(t_{R}+\Delta_{R})(t_{R}-\Delta_{R})-\mu_{R,1}\mu_{R,2}}{t_{R}+\Delta_{R}}&0\end{matrix}\right)\psi^{\gamma}\;.$

Therefore, we can identify each element of this matrix with one coefficient
$\lambda_{\alpha\beta}$ of Eq. (S.1). It should be noted that this
identification is an approximation; also, the separation between bulk and
subgap contributions is only well-defined if the subgap modes are well-
detached from the quasicontinuum.

###  I.2 Comparison between eight and four Majoranas

Since our main objective is to study the physics of a superconducting qubit
modified by the presence of the DQD Josephson junction, we first check the
limitations of the effective Josephson potential obtained previously. At this
level, it is enough to compare results from the projected potential in Eq.
(S.10) with the phase-dependent energy spectrum $E(\phi)$ of the BdG form of
the full Hamiltonian $H^{JJ}_{\mathrm{BdG}}$ before any projection, Fig. S.1.
At the sweet spot ($\Delta=t$, $\mu_{E}=\mu_{I}=0$, Fig. S.1a), the subgap
spectrum shows a $4\pi$ Josephson effect indicating the presence of Majorana
zero modes (thin grey lines). This spectrum originates from the fusion (energy
splitting) of the inner MBSs living in the junction $\gamma_{L,2}^{B}$ and
$\gamma_{R,1}^{A}$ (which is maximum at $\phi=0,2\pi$), but without breaking
the degeneracy point at $\phi=\pi$. Moreover, two states remain at zero energy
for all phases, corresponding to the Majorana states $\gamma_{L,1}^{A}$ and
$\gamma_{R,2}^{B}$ living in the outermost quantum dots. In this regime, both
the full solution (left panel) and the four MBSs projection (right panel)
coincide. Of course, the latter does not capture the bulk solutions that
disperse with phase near $2\Delta=2t$. Deviations from the sweet spot by
changing the internal chemical potential $\mu_{I}\neq 0$ do not affect the low
energy spectrum but open gaps in the bulk (colored lines). When moving away
from the sweet spot by tuning the external chemical potentials $\mu_{E}\neq
0$, while keeping $\mu_{I}=0$, the spectrum remains $4\pi$–periodic. In this
case, the low-energy states are lifted away from zero energy, Fig. S.1b
blue/green colored lines, resulting in a characteristic diamond-like shape.
The crossings forming the diamonds become avoided crossings for $\mu_{I}\neq
0$ and $\mu_{E}\neq 0$, Fig. S.1c, which also splits the crossings of the bulk
bands near $\phi=\pi$, giving an overall $2\pi$-periodic spectrum. In
contrast, a zero-energy state persists for $\mu_{E}=0$ and independently from
$\mu_{I}$, even at large values, Fig. S.1d, corresponding to the Majorana
states of the outermost dots having zero weight in the inner ones. In this
regime, the effect of detuning $\mu_{I}$ away from the sweet spot only affects
the localization of the inner Majorana state, decreasing the splitting between
the blue states, and resulting in a robust $4\pi$-periodic spectrum.

In all the cases described above, the approximation derived in Eq. (S.10)
using four Majorana states describes well the low-energy states of the system
close to the sweet spot. In contrast, this approximation largely deviates from
the results of the full Hamiltonian for sufficiently large
$\mu_{E}\gtrsim\Delta$ and irrespective of $\mu_{I}$, Figs. S.1e–f. In this
regime, the bulk solutions that appear at $E\sim 2\Delta$ at the sweet spot,
hybridize with the low-energy states, renormalizing their energy and strongly
affecting their dispersion against phase. Therefore, the low-energy states
cannot be described by only four Majorana states (one per dot).

Figure S.1: Evolution of the energy spectrum as a function of $\phi$ for the
parameter trajectory indicated in each panel. In each case, the leftmost
panels correspond to the BdG form of the full Hamiltonian –Eqs. (1) and (2)
using the Majorana spinor (3) of the main text– and the rightmost panels to
the four Majoranas projection –Eq. S.10–. Gray/colored levels denote the
beginning/end of each trajectory. We have fixed $t_{J}=t=\Delta$ for every
panel.

We demonstrate the importance of considering all the Majorana states in every
dot by calculating the real space–resolved distribution of the wave functions,
taken as the probability
$P_{j}(\gamma_{\alpha,i}^{A/B})=\langle\psi_{j}|\Psi_{\alpha,i}^{A/B}\rangle\langle\Psi_{\alpha,i}^{A/B}|\psi_{j}\rangle$
of the eigenstate $\ket{\psi_{j}}$ of $H^{JJ}_{\mathrm{BdG}}$ on each mode
$\gamma_{\alpha,i}^{A/B}$, represented in the Majorana basis (S.4). Here
indices $i=1,2$ and $\alpha=L,R$ denote the sites of each chain, whereas
$j=\text{green},\text{blue}$ labels the different levels that appear in Fig.
S.1. As we can see in Fig. S.2, at the sweet spot the outermost Majoranas are
pinned to zero energy (green states in Fig. S.1), whereas (oscillating) blue
states correspond to innermost Majoranas at $\phi=0$. Starting from this
point, varying $\phi$ causes the blue states to delocalize along the junction.
A similar behavior is found on the green states with variations of $\mu_{E}$
outside the sweet spot. Changing $t_{J}$, however, does not cause any change
in the wave functions of the sub–gap states.

The fact that the eigenstates of the system have non–negligible values outside
the low–energy subspace points to a limitation of the projection performed in
the previous section, which is only valid close to the sweet spot. As we
discuss in what follows, a low–energy subspace that is written in terms of
many–body occupations (even and odd) of the system is much more powerful.
Starting first from the four Majoranas projection written in the many–body
fermionic occupation basis (Appendix I.C), we obtain the corresponding subgap
Josephson potential (Eq. 5 in the main text). In Appendix I.D, we go beyond
this picture and describe the effective low–energy physics of the problem in
terms of total many–body occupations including contributions from the four QDs
(eight MBSs) forming the Josephson junction which allows us to obtain a subgap
Josephson potential that includes terms containing both $\mu_{E}$ and
$\mu_{I}$ on equal footing and to all orders (Eq. 8 of the main text).

(a)

(b)

(c)

(d)

Figure S.2: Evolution of the space distribution of sub–gap states as a
function of (a) $t_{J}$ with $\mu_{E}=\mu_{I}=0$ and $\phi=0$; (b) $\phi$ with
$\mu_{E}=\mu_{I}=0$; (c) $\mu_{E}$ with $\mu_{I}=0$ and $\phi=0$; and (d)
$\mu_{I}$ with $\mu_{E}=0$ and $\phi=0$. We have fixed $\Delta=t=t_{J}$ for
all panels, and subtitles refer to each eigenstate plotted in Fig. S.1.

###  I.3 Projection in the left/right fermionic parity basis

We can now write the matrix elements of $V_{\mathrm{DQD}}^{JJ}$ in the
fermionic parity basis $\ket{n_{L},n_{R}}$. For the total even parity state,
the effective Josephson coupling reads

$V_{\mathrm{DQD}}^{JJ}=\left(\begin{matrix}{\left\langle
00\right|H_{\gamma}\left|00\right\rangle}&{\left\langle
00\right|H_{\gamma}\left|11\right\rangle}\\\ {\left\langle
11\right|H_{\gamma}\left|00\right\rangle}&{\left\langle
11\right|H_{\gamma}\left|11\right\rangle}\end{matrix}\right)\;.$ (S.11)

Since the parity states are defined such that (similar for
$c_{R},c^{\dagger}_{R}$)

$\displaystyle
c_{L}^{\dagger}\ket{n_{L},n_{R}}=\sqrt{n_{L}+1}\ket{n_{L}+1,n_{R}}$
$\displaystyle,\quad
c_{L}\ket{n_{L},n_{R}}=\sqrt{n_{L}}\ket{n_{L}-1,n_{R}}\;,$ (S.12)
$\displaystyle\hat{n}_{L}\ket{n_{L},n_{R}}=c_{L}^{\dagger}c_{L}$
$\displaystyle\ket{n_{L},n_{R}}=n_{L}\ket{n_{L},n_{R}}\;,$

and, attending to the decomposition of these fermionic operators in Majorana
operators (S.2), we can write the following operations,

$\displaystyle
i\gamma_{\alpha,1}^{A}\gamma_{\alpha,2}^{B}\ket{00/11}=(2\hat{n}_{\alpha}-1)\ket{00/11}=-/+\ket{00/11}\;,$
(S.13)
$\displaystyle\gamma_{L,1}^{A}\gamma_{R,1}^{A}\ket{00/11}=(c_{L}c_{R}+c_{L}c_{R}^{\dagger}-c_{R}c_{L}^{\dagger}-c_{R}^{\dagger}c_{L}^{\dagger})\ket{00/11}=-/+\ket{11/00}\;,$
$\displaystyle
i\gamma_{L,1}^{A}\gamma_{R,2}^{B}\ket{00/11}=(c_{L}c_{R}-c_{L}c_{R}^{\dagger}-c_{R}c_{L}^{\dagger}+c_{R}^{\dagger}c_{L}^{\dagger})\ket{00/11}=\ket{11/00}\;,$
$\displaystyle
i\gamma_{L,2}^{B}\gamma_{R,1}^{A}\ket{00/11}=(c_{L}c_{R}+c_{L}c_{R}^{\dagger}-c_{R}c_{L}^{\dagger}+c_{R}^{\dagger}c_{L}^{\dagger})\ket{00/11}=\ket{11/00}\;,$
$\displaystyle\gamma_{L,2}^{B}\gamma_{R,2}^{B}\ket{00/11}=(-c_{L}c_{R}+c_{L}c_{R}^{\dagger}+c_{R}c_{L}^{\dagger}+c_{R}^{\dagger}c_{L}^{\dagger})\ket{00/11}=+/-\ket{11/00}\;.$

Therefore, the sub–gap contribution written in the even fermionic parity basis
is

$\displaystyle{\left\langle 00\right|H_{\gamma}\left|00\right\rangle}$
$\displaystyle=-(\lambda_{L1,L2}+\lambda_{R1,R2})\;,$ (S.14)
$\displaystyle{\left\langle 11\right|H_{\gamma}\left|11\right\rangle}$
$\displaystyle=\lambda_{L1,L2}+\lambda_{R1,R2}\;,$ $\displaystyle{\left\langle
00\right|H_{\gamma}\left|11\right\rangle}$
$\displaystyle=i\lambda_{L1,R1}+\lambda_{L1,R2}+\lambda_{L2,R1}-i\lambda_{L2,R2}\;,$
$\displaystyle{\left\langle 11\right|H_{\gamma}\left|00\right\rangle}$
$\displaystyle=-i\lambda_{L1,R1}+\lambda_{L1,R2}+\lambda_{L2,R1}+i\lambda_{L2,R2}\;,$

where $\lambda_{\alpha\beta}$ are the matrix elements of (S.10). Finally, the
sub–gap Josephson potential takes the form

$\displaystyle V_{\mathrm{DQD}}^{JJ}(\phi)=$ (S.15)
$\displaystyle\frac{1}{2}\left(\begin{matrix}\frac{2(t+\Delta)(t-\Delta)-(\mu_{L,1}\mu_{L,2}+\mu_{R,1}\mu_{R,2})}{t+\Delta}&t_{J}\left(1-\frac{\mu_{L,1}\mu_{R,2}}{(t+\Delta)^{2}}\right)\cos\frac{\phi}{2}-it_{J}\frac{\mu_{L,1}-\mu_{R,2}}{t+\Delta}\sin\frac{\phi}{2}\\\
t_{J}\left(1-\frac{\mu_{L,1}\mu_{R,2}}{(t+\Delta)^{2}}\right)\cos\frac{\phi}{2}+it_{J}\frac{\mu_{L,1}-\mu_{R,2}}{t+\Delta}\sin\frac{\phi}{2}&\frac{(\mu_{L,1}\mu_{L,2}+\mu_{R,1}\mu_{R,2})-2(t+\Delta)(t-\Delta)}{t+\Delta}\end{matrix}\right)\;.$

Therefore, we can split this sub–gap effective potential in three different
terms acting on a pseudospin parity space –Eq. (4) of the main text–,

$\displaystyle V_{\mathrm{DQD}}^{JJ}(\phi)$
$\displaystyle=E_{M}\cos\frac{\phi}{2}\sigma_{x}+E_{M}^{S}\sin\frac{\cos}{2}\sigma_{y}+\lambda\sigma_{z}\;,$
(S.16) $\displaystyle
E_{M}=\frac{t_{J}}{2}\left(1-\frac{\mu_{L,1}\mu_{R,2}}{(t+\Delta)^{2}}\right)\;,$
$\displaystyle E_{M}^{S}=t_{J}\frac{\mu_{L,1}-\mu_{R,2}}{2(t+\Delta)}\;,$
$\displaystyle\lambda=\frac{2(t+\Delta)(t-\Delta)-(\mu_{L,1}\mu_{L,2}+\mu_{R,1}\mu_{R,2})}{2(t+\Delta)}\;.$

It is straightforward to see that, when restricting ourselves to the symmetric
case $\mu_{L,1}=\mu_{R,2}=\mu_{E}$ and $\mu_{L,2}=\mu_{R,1}=\mu_{I}$, the
Josephson potential reduces to Eq. (5) of the main text.

##  II Beyond the four Majoranas projection: projection onto a full many–body
parity basis

A reasonable alternative treatment of the problem is to choose as our new
fermionic parity subspace the two lowest–energy many–body eigenstates
$\\{|O_{L}^{-},O_{R}^{-}\rangle,\,|E_{L}^{-},E_{R}^{-}\rangle\\}$ of both
chains isolated from each other ($t_{J}=0$), where

$H_{\alpha}=\left(\begin{matrix}0&0&0&\Delta_{\alpha}\\\
0&-\mu_{\alpha,1}&-t_{\alpha}&0\\\ 0&-t_{\alpha}&-\mu_{\alpha,2}&0\\\
\Delta_{\alpha}&0&0&-(\mu_{\alpha,1}+\mu_{\alpha,2})\end{matrix}\right)\;,$
(S.17)

is the many–body Hamiltonian of one chain in the basis of occupation states
$\\{\ket{00},\,\ket{10},\,\ket{01},\,\ket{11}\\}$. Defining
$\mu_{\alpha}=(\mu_{\alpha,1}+\mu_{\alpha,2})/2$ and
$\delta_{\alpha}=(\mu_{\alpha,1}-\mu_{\alpha,2})/2$, its eigenstates and
eigenenergies are

$\displaystyle\ket{O_{\alpha}^{-}}=\left(0,\,\Psi_{\alpha,1}^{A},\,\Psi_{\alpha,1}^{B},\,0\right)^{T}\propto\left(0,\,\frac{2\delta_{\alpha}+\epsilon_{\alpha
O}^{+}-\epsilon_{\alpha O}^{-}}{2t_{\alpha}},\,1,\,0\right)^{T}$
$\displaystyle,\quad\epsilon_{\alpha
O}^{-}=-\mu_{\alpha}-\sqrt{t_{\alpha}^{2}+\delta_{\alpha}^{2}}\;,$ (S.18)
$\displaystyle\ket{O_{\alpha}^{+}}=\left(0,\,\Psi_{\alpha,2}^{A},\,\Psi_{\alpha,2}^{B},\,0\right)^{T}\propto\left(0,\,\frac{2\delta_{\alpha}-\epsilon_{\alpha
O}^{+}+\epsilon_{\alpha O}^{-}}{2t_{\alpha}},\,1,\,0\right)^{T}$
$\displaystyle,\quad\epsilon_{\alpha
O}^{+}=-\mu_{\alpha}+\sqrt{t_{\alpha}^{2}+\delta_{\alpha}^{2}}\;,$
$\displaystyle\ket{E_{\alpha}^{-}}=\left(\Psi_{\alpha,3}^{A},\,0,\,0,\,\Psi_{\alpha,3}^{B}\right)^{T}\propto\left(\frac{-\epsilon_{\alpha
E}^{+}}{\Delta_{\alpha}},\,0,\,0,\,1\right)^{T}$
$\displaystyle,\quad\epsilon_{\alpha
E}^{-}=-\mu_{\alpha}-\sqrt{\Delta_{\alpha}^{2}+\mu_{\alpha}^{2}}\;,$
$\displaystyle\ket{E_{\alpha}^{+}}=\left(\Psi_{\alpha,4}^{A},\,0,\,0,\,\Psi_{\alpha,4}^{B}\right)^{T}\propto\left(\frac{-\epsilon_{\alpha
E}^{-}}{\Delta_{\alpha}},\,0,\,0,\,1\right)^{T}$
$\displaystyle,\quad\epsilon_{\alpha
E}^{+}=-\mu_{\alpha}+\sqrt{\Delta_{\alpha}^{2}+\mu_{\alpha}^{2}}\;.$

To construct the Hamiltonian of the junction living in the bipartite Hilbert
space $\mathcal{H}_{L}\otimes\mathcal{H}_{R}$, we represent it on the basis of
joint eigenstates
$\\{|i_{L},j_{R}\rangle=|i_{L}\rangle\otimes|j_{R}\rangle\\}$ with
$i,j=O^{\pm},E^{\pm}$. Thus, the Hamiltonian
$\tilde{H}^{JJ}_{\mathrm{DQD}}=\tilde{H}_{L}+\tilde{H}_{R}+\tilde{H}_{J}$ has
a diagonal term

$\displaystyle\tilde{H}_{L}+\tilde{H}_{R}$
$\displaystyle=(P_{L}^{-1}H_{L}P_{L})\otimes\mathbb{I}_{R}+\mathbb{I}_{L}\otimes(P_{R}^{-1}H_{R}P_{R})$
(S.19)
$\displaystyle=\mathrm{diag}\left(\epsilon_{LO}^{-},\,\epsilon_{LO}^{+},\,\epsilon_{LE}^{-},\,\epsilon_{LE}^{+}\right)\otimes\mathbb{I}_{R}+\mathbb{I}_{L}\otimes\mathrm{diag}\left(\epsilon_{RO}^{-},\,\epsilon_{RO}^{+},\,\epsilon_{RE}^{-},\,\epsilon_{RE}^{+}\right)\;,$

where $P_{\alpha}$ is the change–of–basis matrix onto the eigenbasis of each
chain. On the other hand, the off–diagonal term $\tilde{H}_{J}$ is due to the
Josephson tunneling between both chains, which can be easily represented on
the joint–occupation basis
$\\{|n_{L,1},n_{L,2}\rangle\otimes|n_{R,1},n_{R,2}\rangle\\}_{n_{\alpha,i}=0,1}$
and then projected onto the eigenbasis by the change–of–basis matrix
$P_{LR}=P_{L}\otimes P_{R}$.

Finally, the Josephson potential (ignoring higher–order contributions from the
rest of the eigenstates) can be written as

$V^{JJ}_{\mathrm{DQD}}=\left(\begin{matrix}\langle
O_{L}^{-},O_{R}^{-}|\tilde{H}^{JJ}_{\mathrm{DQD}}|O_{L}^{-},O_{R}^{-}\rangle&\langle
O_{L}^{-},O_{R}^{-}|\tilde{H}^{JJ}_{\mathrm{DQD}}|E_{L}^{-},E_{R}^{-}\rangle\\\
\langle
E_{L}^{-},E_{R}^{-}|\tilde{H}^{JJ}_{\mathrm{DQD}}|O_{L}^{-},O_{R}^{-}\rangle&\langle
E_{L}^{-},E_{R}^{-}|\tilde{H}^{JJ}_{\mathrm{DQD}}|E_{L}^{-},E_{R}^{-}\rangle\end{matrix}\right)\;,$
(S.20)

where

$\displaystyle\langle
O_{L}^{-},O_{R}^{-}|\tilde{H}^{JJ}_{\mathrm{DQD}}|O_{L}^{-},O_{R}^{-}\rangle$
$\displaystyle=\epsilon_{LO}^{-}+\epsilon_{RO}^{-}$ (S.21)
$\displaystyle\langle
E_{L}^{-},E_{R}^{-}|\tilde{H}^{JJ}_{\mathrm{DQD}}|E_{L}^{-},E_{R}^{-}\rangle$
$\displaystyle=\epsilon_{LE}^{-}+\epsilon_{RE}^{-}$ $\displaystyle\langle
E_{L}^{-},E_{R}^{-}|\tilde{H}^{JJ}_{\mathrm{DQD}}|O_{L}^{-},O_{R}^{-}\rangle$
$\displaystyle=t_{J}\left(4t^{2}\sqrt{\frac{\epsilon_{RE}^{+}}{\epsilon_{LE}^{+}}}e^{i\phi/2}+\sqrt{\frac{\epsilon_{LE}^{+}}{\epsilon_{RE}^{+}}}(2\delta_{L}+\epsilon_{LO}^{-}-\epsilon_{LO}^{+})(2\delta_{R}+\epsilon_{RO}^{-}-\epsilon_{RO}^{+})e^{-i\phi/2}\right)$
$\displaystyle\times\frac{\Delta\sqrt{(2\delta_{L}+\epsilon_{LO}^{+}-\epsilon_{LO}^{-})(2\delta_{R}+\epsilon_{RO}^{+}-\epsilon_{RO}^{-})}}{8t^{2}\sqrt{(\epsilon_{LO}^{+}-\epsilon_{LO}^{-})(\epsilon_{RO}^{+}-\epsilon_{RO}^{-})(\epsilon_{LE}^{+}-\epsilon_{LE}^{-})(\epsilon_{RE}^{+}-\epsilon_{RE}^{-})}}$
$\displaystyle=-t_{J}\Psi_{L,1}^{A}\Psi_{R,1}^{A}\left(\Psi_{L,3}^{B}\Psi_{R,3}^{A}e^{i\phi/2}-\Psi_{L,4}^{B}\Psi_{R,4}^{A}\frac{\Psi_{L,2}^{A}\Psi_{R,2}^{A}}{\Psi_{L,2}^{B}\Psi_{R,2}^{B}}e^{-i\phi/2}\right)\;.$

One can see that, if the chemical potentials are constrained to the special
symmetric choice $\mu_{L,1}=\mu_{R,2}=\mu_{E}$ and
$\mu_{L,2}=\mu_{R,1}=\mu_{E}$ (internal vs. external), such that
$\mu_{L}=\mu_{R}=\mu_{E}+\mu_{I}=\mu$ and
$\delta_{L}=-\delta_{R}=\mu_{E}-\mu_{I}=\delta$, and considering
$\Delta_{L}=\Delta_{R}$ and $t_{L}=t_{R}$, this Josephson potential reduces to
the simpler form –Eq. (7) of the main text–

$V^{JJ}_{\mathrm{DQD}}(\phi)=\left(\begin{matrix}-2\mu-2\sqrt{t^{2}+\delta^{2}}&\frac{t_{J}\Delta
t}{2\sqrt{(t^{2}+\delta^{2})(\Delta^{2}+\mu^{2})}}\cos(\phi/2)\\\
\frac{t_{J}\Delta
t}{2\sqrt{(t^{2}+\delta^{2})(\Delta^{2}+\mu^{2})}}\cos(\phi/2)&-2\mu-2\sqrt{\Delta^{2}+\mu^{2}}\end{matrix}\right)\;.$
(S.22)

##  III Majorana polarization

The Hamiltonian $H^{JJ}$ described above can be separated into two independent
blocks of even ($\\{|O_{L}^{\pm},O_{R}^{\pm}\rangle$,
$|E_{L}^{\pm},E_{R}^{\pm}\rangle\\}$) and odd
($\\{|E_{L}^{\pm},O_{R}^{\pm}\rangle,\,|E_{L}^{\pm},O_{R}^{\pm}\rangle\\}$)
total parity, which leads to a two–fold degenerate spectrum. To determine
whether these degeneracies are associated with MBSs, we use the Majorana
polarization (MP). This magnitude quantifies the MBS quality and is defined as
the degree that a Hermitian operator localized on one of the quantum dots can
switch between the lowest–energy states of even and odd blocks,

$\displaystyle\mathrm{MP}_{\alpha,i}(O,E)=\frac{w_{\alpha,i}^{2}-z_{\alpha,i}^{2}}{w_{\alpha,i}^{2}+z_{\alpha,i}^{2}}\;,$
(S.23) $\displaystyle w_{\alpha,i}={\left\langle
O\right|c_{\alpha,i}+c_{\alpha,i}^{\dagger}\left|E\right\rangle}\;,$
$\displaystyle z_{\alpha,i}={\left\langle
O\right|c_{\alpha,i}-c_{\alpha,i}^{\dagger}\left|E\right\rangle}\;.$

We can see that, for $t_{J}=0$, MP can be written as

$\mathrm{MP}_{\alpha,i}=\frac{t_{\alpha}\Delta_{\alpha}}{(-1)^{i+1}\delta_{\alpha}\mu_{\alpha}-\sqrt{(t_{\alpha}^{2}+\delta_{\alpha}^{2})(\Delta_{\alpha}^{2}+\mu_{\alpha}^{2})}}\;,$
(S.24)

where $\ket{E}=|O_{L}^{-},O_{R}^{-}\rangle$,
$\ket{O}_{\alpha=L}=|E_{L}^{-},O_{R}^{-}\rangle$,
$\ket{O}_{\alpha=R}=|O_{L}^{-},E_{R}^{-}\rangle$. Restricting ourselves to
$t_{\alpha}=\Delta_{\alpha}$, $|\mathrm{MP}_{\alpha,1}|$
($|\mathrm{MP}_{\alpha,2}|$) is maximum when $\mu_{\alpha}=\delta_{\alpha}$
($\mu_{\alpha}=-\delta_{\alpha}$), that is, when $\mu_{\alpha,2}=0$
($\mu_{\alpha,1}=0$).

Furthermore, from (S.22), the effective Majorana coupling $E_{M}$ is related
to this quantity such that

$E_{M}=\frac{-t_{J}\mathrm{MP}_{\alpha,i}/2}{1+(-1)^{i+\alpha}\frac{\delta\mu}{t\Delta}\mathrm{MP}_{\alpha,i}}\;,$
(S.25)

where $\alpha=\\{0\equiv L,\,1\equiv R\\}$. Thus, if $\mu_{E}=\mu_{I}$
($\mu_{E}=-\mu_{I}$), that is, $\delta=0$ ($\mu=0$), then $E_{M}$ is
proportional to MP: $E_{M}=-t_{J}\mathrm{MP}/2$.

##  IV Intraband splitting in transmon regime

At $n_{g}=0.25$, the energy splitting between the ground state and the first
excited state is merely due to the sub–gap Josephson potential since the rest
of terms on the qubit Hamiltonian give rise to a doubly degenerate state at
this point. Hence, it is reasonable to express the Kitmon qubit frequency
$\omega_{KiT}\equiv\omega_{01}$ as the difference between the two eigenvalues
of $V^{JJ}_{\mathrm{DQD}}(\phi)$,

$\Delta
E^{JJ}(\phi)=2\sqrt{(\sqrt{t^{2}+\delta^{2}}-\sqrt{\Delta^{2}+\mu^{2}})^{2}+E_{M}^{2}\cos^{2}\frac{\phi}{2}}\;.$
(S.26)

As we can see, this difference depends on $\phi$ and, hence, one should know
the explicit form of the qubit wave functions to relate this quantity to
$\omega_{01}$. Nevertheless, in the deep transmon regime ($E_{J}/E_{C}\gg 1$)
these eigenfunctions can be approximated to harmonic–oscilator states
sharpened around $\phi_{ext}$, so that the Kitmon frequency is
$\omega_{KiT}\approx\Delta E^{JJ}(\phi_{ext})$ –Eq. (12) of the main text.
Likewise, in transmon regime the qubit spectrum is insensitive to changes in
the charge offset $n_{g}$, being this approximation valid for every parametric
configuration of the system, even when diagonal terms of
$V^{JJ}_{\mathrm{DQD}}(\phi)$ are not equal and these avoided crossings do not
occur at $n_{g}=0.25$ in charging regime.

Fig. S.3 displays the transition frequency $\omega_{01}(n_{g}=0.25)$ as a
function of different parameters, showing their evolution with increasing
$E_{J}/E_{C}$ ratios. We show the convergence to $\Delta E^{JJ}(\phi_{ext})$
in the limit $E_{J}/E_{C}\gg 1$.

Figure S.3: Transition frequency $\omega_{01}$ for $E_{J}/E_{C}=2,4,10,50$,
compared to analytical result (S.26), black line, as a function of (a)
$\phi_{ext}$ at the sweet spot; (b) $\mu_{E}$ with $\mu_{I}=0$, $\Delta=t$ and
$\phi_{ext}=0$; (c,d) $\mu_{E}=\mu_{I}=\mu$ with $\Delta=t$ and
$\phi_{ext}=0,\pi$, respectively; and (e,f) $\Delta/t$ with
$\mu_{E}=\mu_{I}=0$ and $\phi_{ext}=0,\pi$, respectively. We have fixed
$t_{J}/t=1$ for all panels.

We can also check numerically this approximation by calculating the distance
between the curves that the analytical result (S.26) and $\omega_{01}$ trace
for increasing $E_{J}/E_{C}$ ratios. The distance between two curves described
by the functions $f(x)$ and $g(x)$ over a parametric trajectory
$x\in\mathcal{X}$ is written as

$d(f,g)=\left(\int_{\mathcal{X}}dx\,|f(x)-g(x)|^{2}\right)^{1/2}\;.$ (S.27)

As we can observe in Fig. S.4, increasing the ratio $E_{J}/E_{C}$ minimizes
the distance between numerical results and our analytical approximation, which
allows us to predict $\omega_{KiT}$ with great precision in the deep transmon
regimen.

Figure S.4: Distance between curves $\omega_{01}$ and $\Delta
E^{JJ}(\phi_{ext})$ as a function of $E_{J}/E_{C}$ for the same curves shown
in Fig. S.3 (see legend).

Finally, we include some additional results that show a full progression of
the energy spectrum and its MW response for increasing $E_{J}/E_{C}$ ratios.
In particular, we can see in Fig. S.5 an enhancement of the insensitivity to
the charge offset as the qubit enters in the transmon regime, with a dominant
transition $\omega_{02}$. Furthermore, Fig. S.6 shows how the spectral hole in
$\omega_{02}$ at $\phi_{ext}$ narrows until true energy crossing appears as
the $E_{J}/E_{C}$ ratio increases.

Figure S.5: Full evolution of the energy spectrum and its MW response as a
function of $n_{g}$ at the sweet spot ($\phi_{ext}=0$) for
$E_{J}/E_{C}=1.5,3,5,10$ (from left to right). Figure S.6: Full evolution of
the energy spectrum and its MW response as a function of $\phi_{ext}$ at the
sweet spot for $E_{J}/E_{C}=1.5,3,5,10$ (from left to right).

##  V Numerical methods for the Majorana–transmon qubit: tight–binding
treatment

###  V.1 Phase space

In phase space, the numerical solution of the qubit Hamiltonian

$H_{Q}=4E_{C}(\hat{n}-n_{g})^{2}+V_{J}(\phi)\;,$ (S.28)

is accomplished by discretizing the phase space as $\phi_{j}=2\pi j/l^{\phi}$,
with $j=1,\dots,l^{\phi}$, defining a set of sites arranged into a circular
chain. In so doing, the Hamiltonian acquires a tight–binding form and it
allows us to define a finite fermionic Hilbert space and operators
$b_{j}^{(\dagger)}$ such that their action on the ground state is
$b_{j}^{\dagger}\ket{0}=\Psi(\phi_{j})$, where $\Psi(\phi)$ is the eigenstate
at phase $\phi$.

Then, starting from the definition of the derivative

$\frac{df(x)}{dx}=\lim_{h\to 0}\frac{f(x+h)-f(x-h)}{2h}\;,$ (S.29)

we can express the operator $\hat{n}=-i\partial_{\phi}$ in the discretized
form

$-i\partial_{\phi}=-i\frac{(b_{i+1}^{\dagger}-b_{i-1}^{\dagger})b_{i}}{2a_{\phi}}\;,$
(S.30)

where $a_{\phi}=2\sin(\pi/l^{\phi})$ is a phase lattice constant. By
construction, the second derivative is defined as

$\frac{d^{2}f(x)}{dx^{2}}=\lim_{h\to 0}\frac{f(x+h)-2f(x)+f(x-h)}{h^{2}}\;,$
(S.31)

so we can write

$\partial^{2}_{\phi}=\frac{(b_{i+1}^{\dagger}-2b_{i}^{\dagger}+b_{i-1}^{\dagger})b_{i}}{a_{\phi}^{2}}\;.$
(S.32)

Hence, the Hamiltonian (S.28) reads

$\displaystyle H$
$\displaystyle=\sum_{j}b_{j}^{\dagger}h_{j}^{\phi}b_{j}+\sum_{\langle
j,k\rangle}b_{j}^{\dagger}v_{jk}^{\phi}b_{k}\;,$ (S.33) $\displaystyle
h_{j}^{\phi}=4E_{C}(2a_{\phi}^{-2}+n_{g}^{2})+V_{J}(\phi_{j})\;,$
$\displaystyle
v_{jk}^{\phi}=4E_{C}[\operatorname{sgn}(j-k)in_{g}a_{\phi}^{-1}-a_{\phi}^{-2}]\;,$

where each site element $h_{j}^{\phi},v_{jk}^{\phi}$ is a $2\times 2$ matrix,
owing to the pseudospin structure from even–odd projection.

Secondly, the eigenstates of the Hamiltonian (S.28) are defined as a
two–component spinor $\Psi_{k}=(f_{k}(\phi),g_{k}(\phi))^{T}$ with
periodic/antiperiodic boundary conditions in phase space,
$f(\phi+2\pi)=f(\phi)$ and $g(\phi+2\pi)=-g(\phi)$, due to their even/odd
fermionic parity. To make the Hamiltonian fully periodic, it is rotated
according to $H(\phi)\to UH(\phi)U^{\dagger}$, with
$U={\text{diag}\,}(1,e^{i\phi/2})$. Therefore, the final form of the
Hamiltonian (S.28) is

$H=\left(\begin{matrix}h(n_{g})+V_{J}^{11}&V_{J}^{12}e^{-i\frac{\phi}{2}}\\\
e^{i\frac{\phi}{2}}V_{J}^{21}&h\left(n_{g}+\frac{1}{2}\right)+e^{i\frac{\phi}{2}}V_{J}^{22}e^{-i\frac{\phi}{2}}\end{matrix}\right)\;,$
(S.34)

and hence the site elements $h_{j}^{\phi}$ and $v_{jk}^{\phi}$ change
according to this transformation.

###  V.2 Charge space

In charge representation, the set of states
$\\{\ket{n}\\}_{n=-\infty}^{\infty}$ form a orthonormal basis of such space.
Here, the number of Cooper pairs operator is defined as

$\hat{n}=\sum_{n=-\infty}^{\infty}n\ket{n}\bra{n}\;,$ (S.35)

whereas the action of its conjugate operator $\phi$ on each one of these
states is

$e^{ik\phi}\ket{n}=\ket{n+k}\;.$ (S.36)

Therefore, the Hamiltonian (S.28) can be expressed as

$H=\sum_{n=-\infty}^{\infty}(n-n_{g})^{2}\ket{n}\bra{n}+V_{J}(\phi)\;,$ (S.37)

where the form of the Josephson potential is conditioned by its
phase–dependent terms, being

$\displaystyle\cos(k\phi)$
$\displaystyle=\frac{1}{2}\sum_{n=-\infty}^{\infty}\left(\ket{n+k}\bra{n}+\operatorname{h.c.}\right)\;,$
(S.38) $\displaystyle\sin(k\phi)$
$\displaystyle=\frac{-i}{2}\sum_{n=-\infty}^{\infty}\left(\ket{n+k}\bra{n}-\operatorname{h.c.}\right)\;,$

the most usual of them. Indeed, for more complex potentials, we can perform a
Fourier transform which reduces it to a simple sum of these terms. This
representation gives rise to an identical spectrum to that calculated in phase
space. However, in this case, we require a smaller (truncated) number of sites
$N$ of the tight–binding Hamiltonian matrix, so this method needs less
computational power and time than the other one. Note that, in phase space,
$\dim=2N$ since each site is a spinor with two possible parities, whereas in
charge space we have a set of states $\\{\ket{n}\\}$
($n=-N,-N+1/2,\dots,0,1/2,\dots,N$, so that $\dim=2N+1$.

Indeed, Fig. S.7 shows the convergence of the first four states as a function
of $N$, defined as the maximum number of sites that discretize the
tight–binding space. This convergence is defined as the distance between the
curves that each eigenstate traces (as a function of $n_{g}$) with $N-1$ and
$N$ sites. It is straightforward to see that the tight–binding method
converges much faster in charge space than in phase space.

Figure S.7: Distance between curves $E_{i}^{N-1}(n_{g})$ and
$E_{i}^{N}(n_{g})$ (where $i=0,1,2,3$ labels eigenstates of increasing energy)
at the sweet spot as a function of a cutoff $N$. Numerical methods are
implemented in (a) charge space and (b) phase space.