# Realization of two-qutrit quantum algorithms on a programmable
superconducting processor

Tanay Roy† Present address: Superconducting Quantum Materials and Systems
Center, Fermi National Accelerator Laboratory (FNAL), Batavia, IL 60510, USA
James Franck Institute, University of Chicago, Chicago, Illinois 60637, USA
Department of Physics, University of Chicago, Chicago, Illinois 60637, USA
These authors contributed equally to this work  Ziqian Li† James Franck
Institute, University of Chicago, Chicago, Illinois 60637, USA  Department of
Physics, University of Chicago, Chicago, Illinois 60637, USA  These authors
contributed equally to this work  Eliot Kapit  Department of Physics, Colorado
School of Mines, Golden, CO 80401  David I. Schuster  James Franck Institute,
University of Chicago, Chicago, Illinois 60637, USA  Department of Physics,
University of Chicago, Chicago, Illinois 60637, USA  Pritzker School of
Molecular Engineering, University of Chicago, Chicago, Illinois 60637, USA

###### Abstract

Processing quantum information using quantum three-level systems or qutrits as
the fundamental unit is an alternative to contemporary qubit-based
architectures with the potential to provide significant computational
advantages. We demonstrate a fully programmable two-qutrit quantum processor
by utilizing the third energy eigenstates of two transmons. We develop a
parametric coupler to achieve excellent connectivity in the nine-dimensional
Hilbert space enabling efficient implementations of two-qutrit gates. We
characterize our processor by realizing several algorithms like Deutsch-Jozsa,
Bernstein-Vazirani, and Grover’s search. Our efficient ancilla-free protocols
allow us to show that two stages of Grover’s amplification can improve the
success rates of an unstructured search with quantum advantage. Our results
pave the way for building fully programmable ternary quantum processors using
transmons as building blocks for a universal quantum computer.

##  I Introduction

Noisy intermediate-scale quantum (NISQ) computers [1] are making rapid
progress toward practical applications with demonstrated supremacy over
classical algorithms for specific problems [2, 3]. However, a majority of
those processors utilize quantum two-level systems. Processing information
using $d$-dimensional quantum systems or qudits can boost performance through
access to a larger computational space and with fewer entangling gates for
certain algorithms [4, 5, 6, 7, 8]. In the family of qudits, quantum three-
level systems or qutrits are the closest members to qubits and provide the
immediate opportunity to explore possibilities beyond two levels as basic
units.

Theoretical studies show that qutrit-based processors can be beneficial for
quantum error correction by providing compact logical encoding to protect
against erasure [9, 10] and ternary errors [11], enhance fault tolerance [12],
and enable magic state distillation [13]. Other proposals claim improved
implementations of quantum algorithms [14], simulations [15], cryptography
[16, 17], and communication [18] using qutrits. Due to these advantages,
several platforms including photonic circuits [19, 20], trapped ions [21, 22],
and superconducting circuits [23, 24, 25] have started to explore qutrits as
computational units. In the superconducting circuits platform [26], the
transmon [27] is a natural choice to be utilized as a qutrit due to its fast
gates [23, 25], long coherence times [28, 29], weak anharmonicity, and easy
measurement [30]. While a significant effort is being made in realizing
efficient two-qutrit gates [23, 24, 25] with demonstrations of solving
specific problems [31], implementing multiple generic quantum algorithms on a
single processor has remained a challenge in the circuit QED architecture [26]
due to limitations arising from decoherence, slow, static inter-qutrit
interaction or native entangling operations restricted to a smaller subspace.

In this work, we present a versatile two-qutrit processor with excellent
connectivity between different states of the full Hilbert space facilitating a
rich set of native entangling gates. The connectivity is achieved through
multiple beam-splitting and two-photon squeezing-like primitive operations
enabled by a linear parametric coupler [32, 33]. We demonstrate two-qutrit
versions of Deutsch-Jozsa [34, 35], Bernstein-Vazirani [36, 35], and Grover’s
search [37] algorithms without using any ancilla [38]. Deutsch-Jozsa and
Bernstein-Vazirani algorithms provide exponential and linear speed-ups
respectively over corresponding classical algorithms, whereas Grover’s search
provides a quadratic improvement. We perform two stages of Grover’s
amplification with success probabilities significantly larger than classically
achievable values. To our knowledge, ours is the first successful
demonstration of a qutrit-based Grover’s search across any quantum computing
platform.

##  II Device description

Our two-qutrit processor, shown in Fig. 1(a), is comprised of two transmons
[27] ($Q_{1}$ and $Q_{2}$), each having its readout resonator. Parametric
coupling [32] between the two transmons is realized by grounding the transmon
junctions [33] through a superconducting quantum interference device (SQUID).
The SQUID acts as a tunable inductor when an external DC magnetic field
$\Phi_{\rm ext}$ is threaded into the loop. The transmon pads are also
designed to provide a static capacitive coupling with a strength larger than
the inductive coupling at zero external flux. By applying a finite $\Phi_{\rm
ext}$ one can, thus, nullify the competing capacitive and inductive energies
to minimize the cross-Kerr coupling. The fast-flux line is used to enable
various inter-qutrit interactions by modulating the SQUID at appropriate radio
frequencies, whereas the charge lines can be utilized to drive qutrit-
resonator sidebands. Individual manipulation of the qutrits are done by
sending appropriate microwave pulses to the resonators, and the states of the
qutrits are determined through transmission measurements. We keep large
Josephson-to-charging energy ratios ($\gtrsim$ 90) to maintain long dephasing
times for the second excited states $(\left|2\right\rangle)$ for both
transmons [27]. Figure 1(b) shows the equivalent circuit diagram.

Figure 1: The two-qutrit device. (a) False colored optical micrograph of the
device. Two transmons (red) are inductively connected through a SQUID loop
(purple, inset is an SEM image) that enables inter-qutrit parametric coupling.
The transmons are capacitively coupled to linear resonators (green) for
individual readout. Single-qutrit and readout pulses are sent through the
input lines. (b) Simplified circuit diagram of the device. (c) Energy level
diagram with various sideband interaction strengths (in MHz). The
computational subspace is highlighted with a gray background and the other
states are used as ancillary levels to implement controlled-phase gates.

##  III Hamiltonian and operating point

The system can be described by the following Hamiltonian when no parametric
drives are present

$\dfrac{H}{\hbar}=\sum_{k=1}^{2}\omega_{k}\hat{n}_{k}+\dfrac{\alpha_{k}}{2}\hat{n}_{k}(\hat{n}_{k}-1)+\sum_{j,k\geq
1}J_{jk}(\hat{n}_{1})^{j}(\hat{n}_{2})^{k},$ (1)

where $\omega_{k}$ is
$\left|0\right\rangle\leftrightarrow\left|1\right\rangle$ transition
frequency, $\alpha_{k}$ represents anharmonicity, and
$\hat{n}_{k}=a^{\dagger}_{k}a_{k}$ is the photon number operator with
$a^{\dagger}_{k}(a_{k})$ being the creation (annihilation) operator for the
$k$-th qutrit. At $\Phi_{\rm ext}=0$, only second order parametric processes
are allowed, whereas biasing the SQUID loop sufficiently away from the sweet-
spot enables first order processes. The cross-Kerr couplings $J_{jk}$ also
tune with $\Phi_{\rm ext}$ and we choose a bias flux $\Phi_{\rm min}$ so that
the leading order dispersive shift $J_{11}$ is minimized. We treat $\Phi_{\rm
min}$ as the coupler’s ‘off’ point as the higher order cross-Kerr shifts
$J_{12}$ and $J_{21}$ are very close to minimum at this bias (see Appendix C
for details) ensuring reduced crosstalk during idling periods or gate
operations. The value of $\Phi_{\rm min}$ is experimentally determined by
performing Ramsey-fringe experiments on $Q_{1}$ conditioned on $Q_{2}$ being
in the ground or first excited state and minimizing the difference in the
oscillation frequencies with respect to $\Phi_{\rm ext}$.

We activate qutrit-qutrit (QQ) red and blue sidebands corresponding to photon-
exchange (beam splitting)
$\left|m,n\right\rangle\leftrightarrow\left|m-1,n+1\right\rangle$ and two-
photon pumping (squeezing)
$\left|m,n\right\rangle\leftrightarrow\left|m+1,n+1\right\rangle$ interactions
respectively by modulating the SQUID at corresponding transition frequencies
(see Fig. 1(c)). Here the state of the system $\left|m,n\right\rangle$
indicates the qutrits being in the energy eigenstates $\left|m\right\rangle$
and $\left|n\right\rangle$ respectively. The SQUID (coupler) mode is designed
to have a frequency much larger than any relevant radio-frequency (RF) drives
during the experiment at the operating point and is never populated. The blue
sideband frequencies (around 6.5 - 6.9 GHz) are usually one order of magnitude
higher than the red sidebands (around 0.2 - 0.7 GHz). For the current
geometry, stray capacitive coupling of the flux line to the SQUID [33] limits
the blue sideband rates due to frequency-enhanced unwanted interactions.

##  IV Gate operations

Single-qutrit operations are performed by applying separate microwave pulses
at $\left|0\right\rangle\leftrightarrow\left|1\right\rangle$ and
$\left|1\right\rangle\leftrightarrow\left|2\right\rangle$ transition
frequencies of both qutrits through the readout resonators. We use Gaussian-
edge rectangular pulses having $2\sigma$ edge lengths. Single-qutrit drives
use roughly 5 MHz of Rabi rates, whereas phase gates are realized by simply
advancing the phases of the appropriate subsequent pulses. Consequently, the
virtual phase gates take no time and are near-perfect. For the quantum
algorithms we use qutrit-Hadamard $H$, bit-shift $X$ (extension of the qubit’s
NOT gate), and qutrit-phase $Z$ gates having the following matrix
representations

$H=\dfrac{1}{\sqrt{3}}\begin{bmatrix}1&1&1\\\ 1&\omega&\omega^{2}\\\
1&\omega^{2}&\omega\end{bmatrix},X=\begin{bmatrix}0&0&1\\\ 1&0&0\\\
0&1&0\end{bmatrix},Z=\begin{bmatrix}1&0&0\\\ 0&\omega&0\\\
0&0&\omega^{2}\end{bmatrix},$ (2)

where $\omega=e^{2\pi i/3}$, the cube-root of unity. Implementation of $H,X,$
and $Z$ gates require the application of three, two, and zero physical pulses
respectively (see Appendix A).

Figure 2: Deutsch-Jozsa and Bernstein-Vazirani algorithms. (a) Quantum circuit
for the algorithms. In DJ algorithm, gates $W_{1}$,
$W_{2}\in{\\{I,X,X^{2},Z,Z^{2}\\}}$ are applied to implement a constant or a
balanced oracle. The final output state being in $\left|00\right\rangle$ or
non $\left|00\right\rangle$ distinguishes the two cases. (b) Experimental
results for DJ algorithm. The rows and columns represent gates applied to
$Q_{1}$ and $Q_{2}$ respectively. The average SPs are 75.5(3)% and 98.5(1)%
for the constant (hatched boxes) and balanced (plain boxes) oracles
respectively beating the classical rate of 50%. (c) Experimental results for
BV algorithm. Each row corresponds to a specific oracle with the mapping
$\\{I,Z,Z^{2}\\}\rightarrow\\{0,1,2\\}$. The diagonal terms show the SPs for
all nine strings mapped to the basis states with an average of 78.3(3)%, which
is much larger than the classical SP of 33.3%.

We can implement arbitrary two-qutrit unitaries by combining single-qutrit
rotations along with generalized controlled-phase (CPhase) gates
$C_{p}(\theta,\left|mn\right\rangle)=\mathcal{I}-(1-e^{i\theta})\left|mn\right\rangle$
representing an accumulation of $\theta$ phase on the eigenstate
$\left|mn\right\rangle$, where $\mathcal{I}$ is the nine-dimensional identity
matrix. We realize $C_{p}(\pi,\left|22\right\rangle)$ gate by applying a 94 ns
long $2\pi$-pulse at the red sideband
$\left|22\right\rangle\leftrightarrow\left|31\right\rangle$. Similarly, a 56
ns $2\pi$-pulse at the red sideband
$\left|21\right\rangle\leftrightarrow\left|30\right\rangle$ enables a
$C_{p}(\pi,\left|21\right\rangle)$ gate (see Fig. 1(c)). We sandwich these
$C_{\theta}$ gates between single-qutrit $\pi$ gates to achieve phase flip on
any computational basis (see Appendix D). While we have access to a native
$C_{p}(\pi,\left|12\right\rangle)$ gate, we avoid that due to lower fidelity.
Note that the bosonic enhancement of the coupling to an external field makes
$C_{p}(\pi,\left|22\right\rangle)$ twice as fast compared to
$C_{p}(\pi,\left|11\right\rangle)$ when restricting to the qubit subspace.
Next, we implement three two-qutrit quantum algorithms on our processor.

##  V Deutsch-Jozsa algorithm

The Deutsch-Jozsa (DJ) algorithm[34] is one of the earliest quantum algorithms
showing an exponential advantage over any classical algorithm. The original DJ
algorithm applies to qubits, and here we extend it to qutrits [35] for our
quantum processor. For an $n$-qutrit system, the task of the DJ algorithm is
to distinguish a function $f:\\{0,1,2\\}^{n}\rightarrow\\{0,1,2\\}$, which
takes $n$-trits as an input and outputs one trit, between two cases, a
balanced or a constant function. The constant function always results in the
same output (0, 1, or 2) independent of the input, whereas the balanced
function outputs each of the three possibilities for exactly one-third of the
possible inputs. Note that implementing the different test functions is often
termed as the oracle. A deterministic classical algorithm needs $3^{n-1}+1$
queries (with at least two queries in the best case) to distinguish the two
cases, whereas DJ algorithm needs only one and hence provides the exponential
speed-up.

Figure. 2(a) depicts our (ancilla-free) circuit implementation of the DJ
algorithm. Two Hadamard gates are simultaneously applied to both qutrits
initialized to $\left|0\right\rangle$ (ground state) to prepare the state
$\frac{1}{3}\left(\left|0\right\rangle+\left|1\right\rangle+\left|2\right\rangle\right)^{\otimes
2}$. The oracles (gray gates in Fig. 2(a)) are implemented by applying gates
$(W)$ to the qutrits chosen from the set
$\mathcal{S}=\\{I,X,X^{2},Z,Z^{2}\\}$. For the constant case, $W_{1}$ and
$W_{2}$ are picked up from the subset $\mathcal{S}_{c}=\\{I,X,X^{2}\\}$, and
the total number of $X$ gates modulo $3$ specifies the constant output. For
example, the gate $X\otimes X^{2}$ would implement constant 0, whereas
$I\otimes X^{2}$ would be the case of constant 2. A balanced oracle can be
realized by choosing any combination of elements from the set $\mathcal{S}$
except for those cases that result in a constant function. 16 different kinds
of balanced functions (see Appendix F for mapping to various addition
functions) are implemented in our experiment as shown in Fig. 2(a). Finally,
two $H^{\dagger}$ gates are applied before simultaneous readout. A final
measured state of $\left|00\right\rangle$ indicates a constant function,
whereas any other output implies a balanced function. The theoretical success
probability (SP) for each case is 100%, and the experimental results are
summarized in Fig. 2(b) where the hatched (unhatched) boxes represent constant
(balanced) cases. The average SPs for the three constant cases with outputs
$\\{0,1,2\\}$ are separately $72.8(3)\%$, $76.5(3)\%$ and $77.2(3)\%$, marked
with horizontal, left and right hatching. The numbers in parentheses represent
standard error of mean obtained after 20,000 repetitions of each oracle. All
experimental data are corrected for measurement error (see Appendix E). For
the 16 balanced cases, the average SP is $98.5(1)\%$. The SPs for all cases
are well above the classical case, which would be 50% after a single query.

##  VI Bernstein-Vazirani algorithm

The Bernstein-Vazirani (BV) algorithm[36] for qutrits can be restated as
follows: given an oracle
$f\left(\boldsymbol{x,s}\right)\equiv\sum_{j=1}^{n}x_{j}s_{j}\pmod{3}$ that
performs inner product between two strings of ternaries followed by modulo 3,
the goal is to determine the unknown string
$\boldsymbol{s}=\\{s_{1},s_{2},\ldots,s_{n}\\}$ where the user has control
over the input string $\boldsymbol{x}$. The most efficient classical algorithm
will need $n$ oracle queries to find all digits of $\boldsymbol{s}$. BV
algorithm, on the contrary, needs only one query and the quantum circuit is
identical to the case of DJ as shown in Fig. 2(a). Oracles representing 9
different strings for our two-qutrit system are implemented by choosing gates
($W$) from the set $\\{I,Z,Z^{2}\\}^{\otimes 2}$ with the mapping
$\\{I,Z,Z^{2}\\}\rightarrow\\{0,1,2\\}$. The final state of the system after
measurement directly reveals the unknown string with 100% theoretical success
rate. Fig. 2(c) tabulates the experimental results where the vertical axis
represents gates applied to the qutrits corresponding to different unknown
strings, and the horizontal axis shows the measured probability for each
state. The diagonal entries indicate individual SPs for each input string
mapped to the final state. The average SP for all $9$ cases is $78.3(3)\%$,
which is far above the classical SP of 33.3% after one query.

Figure 3: Grover’s search algorithm for two-qutrits. (a) Quantum circuit. The
oracles are implemented by CPhase gates
$C_{\theta}\left(\pi,\left|jk\right\rangle\right)$. The diffusion operator
amplifies the detection probability of the marked state. (b) Experimental
results. Detection probabilities (corrected for measurement error) after one
(top panel) and two (bottom panel) rounds of amplitude amplification are
obtained with 20,000 averages. All individual success rates are far beyond the
corresponding classical SPs of 11.1% and 22.2%.

##  VII Grover’s Search

Grover’s algorithm [37] provides a quadratic speed-up for searching an
unstructured database. For a database of size $\mathcal{N}$, the algorithm can
find the unique input that satisfies a certain condition using
$O\left(\sqrt{\mathcal{N}}\right)$ search queries, while a classical algorithm
requires on average ${\mathcal{N}}/{2}$ repetitions. Several groups have
recently realized Grover’s search on qubit-based platforms [39, 38, 40, 41].
For the two-qutrit case with $\mathcal{N}=3^{2}=9$, the classical SPs with one
and two rounds of search are $\frac{1}{9}=11.1\%$ and
$\frac{1}{9}+\frac{8}{9}\cdot\frac{1}{8}=22.2\%$ respectively. The
corresponding theoretical SPs for the original Grover’s search are $72.6\%$
and $98.4\%$, and can also be modified to achieve determinism [42, 43].

Our ancilla-free version of the two-qutrit Grover’s search is illustrated in
Fig. 3(a). It has four stages: initialization, oracle implementation,
amplitude amplification, and measurement. Starting from the ground state, we
apply Hadamard gates on both qutrits to initialize the system to the equal
superposition state
$\frac{1}{3}\left(\left|0\right\rangle+\left|1\right\rangle+\left|2\right\rangle\right)^{\otimes
2}$. A controlled-phase gate
$C_{\theta}\left(\pi,\left|jk\right\rangle\right)$, that flips the phase of
the target state $\left|jk\right\rangle$, is used to realize the oracle. Due
to the structure of the Hilbert space (see Fig 1(c)), we have access to two
high-fidelity CPhase gates $C_{\theta}\left(\pi,\left|22\right\rangle\right)$
and $C_{\theta}\left(\pi,\left|21\right\rangle\right)$ and all other CPhase
gates are realized in conjunction with single-qutrit rotations (see Appendix
D). Amplitude amplification of the marked state happens through the Grover’s
diffusion or reflection unitary, which is constructed using a phase flip of
the $\left|22\right\rangle$ state sandwiched between Hadamard and $Z$ gates.
Here, we utilized the decomposition
$C_{\theta}\left(\pi,\left|00\right\rangle\right)=(ZZ)\otimes
C_{\theta}\left(\pi,\left|22\right\rangle\right)\otimes(Z^{\dagger}Z^{\dagger})$.
Simultaneous measurements on the qutrits are performed after one and two
iterations of Grover’s search for each target state. During a similar search
using three qubits ($\mathcal{N}=8$), each oracle (and amplification) step
requires eight CNOT gates (for a linear chain) [41], resulting in an eight-
fold rise in entangling operations compared to our efficient two-qutrit
implementation.

Figure 3(b) shows the experimentally obtained detection probabilities for the
9 different marked states after single (top panel) and double (bottom panel)
rounds of the Grover’s iteration. Each row represents a probability
distribution acquired with 20,000 repetitions and after correcting for
measurement error. The diagonal terms represent successful detection rates
with an average SP of $44.4(3)\%$ after the first round, which increases to
$49.6(3)\%$ with the second iteration. The performance degradation of target
states closer to $\left|00\right\rangle$ is caused by the less-efficient
implementations of the corresponding oracles, where more single-qutrit
rotations are required for the CPhase gate decomposition (see Appendix D). As
promised by the algorithm, experimental detection probabilities for the
individual correct states increased after the second iteration for all cases
(except for $\left|02\right\rangle$ and $\left|12\right\rangle$, which we
attribute to the lower lifetime of $Q_{2}$’s $\left|2\right\rangle$ level).
The experimental outcome also has a good agreement with our Master equation
simulation as described in Appendix G indicating dephasing limited
performance. The average experimental SPs beat the classical rates of 11.1%
and 22.2% for the two rounds by more than a factor of 2, with clear
improvement in performance after the second iteration.

##  VIII Conclusion

We demonstrate a fully programmable two-qutrit superconducting processor based
on transmons. We utilize a linear coupler to obtain excellent connectivity in
the underlying Hilbert space through parametrically activated fast sidebands.
The success probabilities of all the algorithms implemented, namely Deutsch-
Jozsa, Bernstein-Vazirani, and Grover’s search, are significantly higher than
corresponding classical rates. Notably, we achieve improved success rates in
finding the correct answer after performing the second iteration of Grover’s
search. To our knowledge, this is the first experimental demonstration of two-
qutrit Grover’s search across all quantum computing platforms. The performance
of the processor is currently limited by low dephasing times, and we expect
significant improvement with increased coherence [28, 29], pulse shaping [44],
geometrical optimization, and improved readout [26]. Our results indicate that
qutrit-based processors could be a promising candidate for building a large-
scale quantum computer and we believe, will foster investigations in the field
of ternary quantum information processing.

##  IX Acknowledgements

This work was supported by AFOSR Grant No. FA9550-19-1-0399 and ARO Grant No.
W911NF-17-S-0001. Devices are fabricated in the Pritzker Nanofabrication
Facility at the University of Chicago, which receives support from Soft and
Hybrid Nanotechnology Experimental (SHyNE) Resource (NSF ECCS-1542205), a node
of the National Science Foundation’s National Nanotechnology Coordinated
Infrastructure. This work also made use of the shared facilities at the
University of Chicago Materials Research Science and Engineering Center.

##  Appendix A Single-qutrit gates

The native operations accessible for a single qutrit are

$\displaystyle R_{01}\left(\phi,\theta\right)$
$\displaystyle=\begin{bmatrix}\cos{\frac{\theta}{2}}&-e^{-i\phi}\sin{\left(\theta/2\right)}&0\\\
e^{i\phi}\sin{\frac{\theta}{2}}&\cos{\frac{\theta}{2}}&0\\\
0&0&1\end{bmatrix},$ (3a) $\displaystyle R_{12}\left(\phi,\theta\right)$
$\displaystyle=\begin{bmatrix}1&0&0\\\
0&\cos{\frac{\theta}{2}}&-e^{-i\phi}\sin{\left(\theta/2\right)}\\\
0&e^{i\phi}\sin{\frac{\theta}{2}}&\cos{\frac{\theta}{2}}\end{bmatrix},$ (3b)
$\displaystyle\Theta\left(x,y\right)$ $\displaystyle=\begin{bmatrix}1&0&0\\\
0&e^{ix}&0\\\ 0&0&e^{i\left(x+y\right)}\end{bmatrix},$ (3c)

where, $R_{01}\left(\phi,\theta\right)$ and $R_{12}\left(\phi,\theta\right)$
are performed through sending drives at
$\left|0\right\rangle\leftrightarrow\left|1\right\rangle$ and
$\left|1\right\rangle\leftrightarrow\left|2\right\rangle$ transitions
respectively with appropriate lengths and phases. However, the generic phase
gate $\Theta\left(x,y\right)$ is implemented virtually as it requires no
physical pulses. For each qutrit we record two phase parameters $\theta_{01}$
and $\theta_{12}$ corresponding to the pulses applied to the
$\left|0\right\rangle\leftrightarrow\left|1\right\rangle$ and
$\left|1\right\rangle\leftrightarrow\left|2\right\rangle$ transitions
respectively. To apply the $\Theta\left(x,y\right)$ gate, we advance both
$\theta_{01}$ and $\theta_{12}$ by $x$ and $y$ for all the subsequent pulses.
The single-qutrit $Z$ gate becomes a special case of the generic phase gate:
$Z=\Theta\left(\frac{2\pi}{3},\frac{2\pi}{3}\right)$. Since, all phase updates
are performed in software, the $Z$ gate (or a $\Theta\left(x,y\right)$ gate)
has nearly $100\%$ fidelity.

An arbitrary single-qutrit unitary can be constructed using combinations of
Eqs. (3). For example, the $H$ gate is deconstructed as

$H=\frac{1}{\sqrt{3}}\begin{bmatrix}1&1&1\\\
1&e^{\frac{2\pi}{3}i}&e^{\frac{4\pi}{3}i}\\\
1&e^{\frac{4\pi}{3}i}&e^{\frac{2\pi}{3}i}\end{bmatrix}\\\
=R_{12}\left(0,\frac{\pi}{2}\right)\cdot
R_{01}\left(0,\beta\right)\cdot\Theta\left(\pi,\frac{\pi}{2}\right)\cdot
R_{12}\left(0,\frac{\pi}{2}\right)\cdot\Theta\left(0,\pi\right),$ (4)

with $\beta=2\tan^{-1}(\sqrt{2})$. Similarly, the bit-shift gate is decomposed
as

$X=\begin{bmatrix}0&0&1\\\ 1&0&0\\\ 0&1&0\end{bmatrix}=R_{01}(0,\pi)\cdot
R_{12}(0,\pi).$ (5)

We obtain process fidelities of 98.96% ($Q_{1}$) and 97.06% ($Q_{2}$) for the
$H$ gate as shown in Fig. 4. The same for the $Z$ gate are 97.48% ($Q_{1}$)
and 96.76% ($Q_{2}$). Even though the $Z$ gates should be nearly perfect, the
process fidelities are limited by state preparation and measurement (SPAM)
error.

Figure 4: Process tomography for Hadamard gates. The top and bottom panels show real and imaginary components of the process matrix ($\chi$) for (a) an ideal case, (b) $Q_{1}$, and (c) $Q_{2}$. Gate type | On $Q_{1}$ (ns) | On $Q_{2}$ (ns)
---|---|---
$\pi/2_{01}$ | 49.50 | 49.71
$\pi/2_{12}$ | 41.27 | 44.15
$\pi_{01}$ | 94.98 | 95.41
$\pi_{12}$ | 78.52 | 84.28
$Z$ | 0 | 0
$H$ | 141.88 | 147.90
Table 1: Total gate lengths for different single-qutrit operations. The $Z$
gates are implemented virtually.

All single-qutrit gates use rectangular pulses with Gaussian edges. The
envelope shape $h\left(t\right)$ is defined in Eq. 6 with $\sigma=2.5$ ns,

$h\left(t\right)=\begin{cases}A_{0}e^{\frac{-\left(t-t_{0}-2\sigma\right)^{2}}{2\sigma^{2}}}&\text{if
$t_{0}<t<t_{0}+2\sigma$,}\\\ A_{0}&\text{if
$t_{0}+2\sigma<t<t_{1}-2\sigma$,}\\\
A_{0}e^{\frac{-\left(t_{1}-2\sigma-t\right)^{2}}{2\sigma^{2}}}&\text{if
$t_{1}-2\sigma<t<t_{1}$,}\\\ 0&\text{Other cases},\\\ \end{cases}$ (6)

where, $2\sigma$ is the Gaussian tail length, $A_{0}$ is the amplitude, and
$t_{1}-t_{0}$ is the total pulse length.

##  Appendix B Device fabrication and parameters

The device was fabricated on a 430 $\mu$m thick C-plane sapphire substrate
annealed at $1200^{o}$C for 2 hours. A 200 nm thick film of Tantalum was
sputtered at $800^{o}$C to form the ground plane. The large features
(excluding the Josephson junctions) were made via optical lithography followed
by wet-etching (dipped for 20 seconds) at wafer-scale using Transene Tantalum
etchant 111. A 600 nm thick layer of AZ 1518 was used as the (positive)
photoresist, and a Heidelberg MLA 150 Direct Writer was used for the
photolithography. The junction mask was fabricated via electron-beam
lithography using a Raith EBPG5000 Plus E-Beam Writer on a bi-layer resist
(MMA-PMMA) comprising of MMA EL11 and 950PMMA A7. Both small (transmon) and
large (coupler) Josephson junctions were made with the Dolan bridge technique.
Those were subsequently evaporated in a Plassys electron-beam evaporator with
double angle evaporation ($\pm 23^{o}$). The wafer was then diced into
$7\times 7$ mm chips, lifted-off, mounted on a printed circuit board, and
subsequently wire-bonded. The coherence times and frequency parameters are
listed in Table. 2 and Table. 3.

Qutrit |  $T_{1}^{01}(\mu$s) |  $T_{1}^{12}(\mu$s) |  $T_{2R}^{01}(\mu$s) |  $T_{2R}^{12}(\mu$s)
---|---|---|---|---
$Q_{1}$ | $47.9$ | 21.7 | $4.5$ | $2.0$
$Q_{2}$ | $35.1$ | 3.9 | $3.2$ | $2.4$
Table 2: Device coherence parameters. $T^{jk}_{1}$ and $T^{jk}_{2R}$ respectively represent the relaxation and Ramsey dephasing time constants for the $\left|j\right\rangle\leftrightarrow\left|k\right\rangle$ transition. Parameter | Symbol | Value$/2\pi$
---|---|---
$Q_{1}$ $\left|0\right\rangle\leftrightarrow\left|1\right\rangle$ transition | $\omega_{1}$ |  $3.3494$ (GHz)
$Q_{2}$ $\left|0\right\rangle\leftrightarrow\left|1\right\rangle$ transition | $\omega_{2}$ |  $3.8310$ (GHz)
$Q_{1}$ anharmonicity | $\alpha_{1}$ |  $-115.2$ (MHz)
$Q_{2}$ anharmonicity | $\alpha_{2}$ |  $-159.8$ (MHz)
Readout1 frequency | $\omega_{r1}$ |  $4.9602$ (GHz)
Readout2 frequency | $\omega_{r2}$ |  $5.4225$ (GHz)
$\left(E_{\left|11\right\rangle}-E_{\left|01\right\rangle}\right)-\left(E_{\left|10\right\rangle}-E_{\left|00\right\rangle}\right)$ | $ZZ$ | -238 (kHz)
$\left(E_{\left|21\right\rangle}-E_{\left|11\right\rangle}\right)-\left(E_{\left|20\right\rangle}-E_{\left|10\right\rangle}\right)$ | $ZZ_{2110}$ | -148 (kHz)
$\left(E_{\left|12\right\rangle}-E_{\left|11\right\rangle}\right)-\left(E_{\left|02\right\rangle}-E_{\left|01\right\rangle}\right)$ | $ZZ_{1021}$ | -183 (kHz)
$\left(E_{\left|22\right\rangle}-E_{\left|12\right\rangle}\right)-\left(E_{\left|20\right\rangle}-E_{\left|10\right\rangle}\right)$ | $ZZ_{2120}$ | -211 (kHz)
$\left(E_{\left|22\right\rangle}-E_{\left|21\right\rangle}\right)-\left(E_{\left|02\right\rangle}-E_{\left|01\right\rangle}\right)$ | $ZZ_{2021}$ | -262 (kHz)
$\left(E_{\left|12\right\rangle}-E_{\left|02\right\rangle}\right)-\left(E_{\left|10\right\rangle}-E_{\left|00\right\rangle}\right)$ | $ZZ_{1020}$ | -402 (kHz)
$\left(E_{\left|21\right\rangle}-E_{\left|20\right\rangle}\right)-\left(E_{\left|01\right\rangle}-E_{\left|00\right\rangle}\right)$ | $ZZ_{2010}$ | -403 (kHz)
Coefficient of $(a^{\dagger}a)(b^{\dagger}b)$ | $J_{11}$ | -304.3 (kHz)
Coefficient of $(a^{\dagger}a)^{2}(b^{\dagger}b)$ | $J_{21}$ | 37.8 (kHz)
Coefficient of $(a^{\dagger}a)(b^{\dagger}b)^{2}$ | $J_{12}$ | 23.6 (kHz)
Coefficient of $(a^{\dagger}a)^{2}(b^{\dagger}b)^{2}$ | $J_{22}$ | 5.4 (kHz)
Table 3: Experimentally obtained various frequency parameters.

##  Appendix C Hamiltonian

Capacitance | (fF) | Josephson energy | (GHz)
---|---|---|---
$C_{q1}$ | $178.0$ | $E_{J1}$ | $13.6$
$C_{q2}$ | $131.0$ | $E_{J2}$ | $13.3$
$C_{c}$ | $193.6$ | $E_{Jc}$ | $1140.0$
$C_{q12}$ | $2.0$ |  |
Table 4: Estimated capacitances and Josephson energies. Capacitances are
extracted from Ansys Q3D simulations. Josephson energies are estimated using
Ambegaokar-Baratoff relation and room temperature resistances of on-chip
(nominally identical) test junctions.

We develop a simplified model to extract various parameters and explain the
sideband interactions. Omitting the readout resonators, the bare Hamiltonian
can be expressed as

$\displaystyle H=$
$\displaystyle\overrightarrow{q}^{\intercal}C_{L}^{-1}\overrightarrow{q}-E_{J1}\cos{\left(\phi_{c}-\phi_{1}\right)}-E_{J2}\cos{\left(\phi_{c}-\phi_{2}\right)},$
$\displaystyle-E_{Jc}\cos{\left(2\pi\frac{\Phi_{\rm
ext}}{\Phi_{0}}\right)}\cos{\left(\phi_{c}\right)},$ (7)

with the capacitance matrix

$C_{L}=\begin{bmatrix}C_{1}+C_{12}&-C_{12}&0\\\ -C_{12}&C_{2}+C_{12}&0\\\
0&0&C_{1}+C_{2}+C_{c}\end{bmatrix}.$ (8)

Here $\\{\phi_{1},\phi_{2},\phi_{c}\\}$ are the node fluxes for the two qubits
and the coupler respectively.
$\overrightarrow{q}^{\intercal}=\left(q_{1},q_{2},q_{c}\right)$ represent the
charge vector conjugate to the node fluxes, so that
$\left[\phi_{k},q_{k}\right]=i$. The coupler capacitance $C_{c}$ originates
from the self-capacitance of the large coupler junctions.

We first extract the linear part of the Hamiltonian,

$H_{\rm
lin}=\overrightarrow{q}^{\intercal}C_{L}^{-1}\overrightarrow{q}+\frac{E_{J1}}{2}\left(\phi_{c}-\phi_{1}\right)^{2}+\frac{E_{J2}}{2}\left(\phi_{c}-\phi_{2}\right)^{2}\\\
+\frac{E_{Jc}}{2}\cos{\left(2\pi\frac{\Phi_{\rm
ext}}{\Phi_{0}}\right)}\phi_{c}^{2}$ (9)

Table. 4 includes all coefficients used in the simulation. We rewrite the
charge and phase variables in the dressed basis such that $H_{\rm lin}$ is
simultaneously diagonalized for both $\overrightarrow{n}$ and
$\overrightarrow{\phi}$ with the transformation matrix $U$:

$\displaystyle H_{\rm
lin}=\sum_{k=1,2,c}\left(\tilde{C_{k}}\tilde{n}_{k}^{2}+\tilde{D}_{k}\tilde{\phi}_{k}^{2}\right)$
(10a)
$\displaystyle\overrightarrow{\tilde{n}}=\left(\tilde{n}_{1},\tilde{n}_{2},\tilde{n}_{c}\right)^{\intercal}=U^{-1}\overrightarrow{n}$
(10b)
$\displaystyle\overrightarrow{\tilde{\phi}}=\left(\tilde{\phi}_{1},\tilde{\phi}_{2},\tilde{\phi}_{c}\right)^{\intercal}=U^{-1}\overrightarrow{\phi}$
(10c) $\displaystyle U=\begin{bmatrix}U_{11}&U_{12}&U_{1c}\\\
U_{21}&U_{22}&U_{2c}\\\ U_{c1}&U_{c2}&U_{cc}\end{bmatrix}.$ (10d)

Then we quantize the circuit on the dressed basis,

$\displaystyle\tilde{n}_{k}=\frac{i}{\sqrt{2}}\sqrt{\frac{\tilde{D}_{k}}{\tilde{C}_{k}}}\left(a_{k}^{\dagger}-a_{k}\right),$
(11a)
$\displaystyle\tilde{\phi}_{k}=\frac{1}{\sqrt{2}}\sqrt{\frac{\tilde{C}_{k}}{\tilde{D}_{k}}}\left(a_{k}^{\dagger}+a_{k}\right),$
(11b)

and insert the nonlinear part back into the Hamiltonian,

$H=\sum_{k=1,2,c}\left(\tilde{C}_{k}\tilde{n_{k}}^{2}\right)-E_{J1}\cos{\left(\sum_{k=1,2,c}\left(U_{ck}\tilde{\phi}_{k}-U_{1k}\tilde{\phi}_{k}\right)\right)}\\\
-E_{J2}\cos{\left(\sum_{k=1,2,c}\left(U_{2k}\tilde{\phi}_{k}-U_{ck}\tilde{\phi_{k}}\right)\right)}\\\
-E_{Jc}\cos{\left(2\pi\frac{\Phi_{\rm
ext}}{\Phi_{0}}\right)}\cos{\left(\sum_{k=1,2,c}U_{ck}\tilde{\phi}_{k}\right)}.$
(12)

The resulting transition frequencies and various $ZZ$ interactions strengths
as a function of $\Phi_{\rm ext}$ obtained from the numerical diagonalization
along with the experimental data are shown in Fig. 5. We notice some deviation
when the DC flux is biased very close to $\frac{\pi}{2}$. This is caused by
the parasitic SQUID-loop inductance. However, around the DC flux position
($\Phi_{\rm ext}=$0.185$\Phi_{0}$) where we implement all quantum algorithms
(vertical dashed line in Fig. 5(b)), we have a good match between the
numerical and experimental data.

Figure 5: Circuit quantization results. (a) Comparison of the experimentally
obtained $\left|0\right\rangle\leftrightarrow\left|1\right\rangle$ transition
frequency $\omega_{01}$ (solid circle) and
$\left|1\right\rangle\leftrightarrow\left|2\right\rangle$ transition frequency
$\omega_{12}$ (hollow circle) as a function of $\Phi_{\rm ext}$ with
simulations (solid/dashed lines) for qutrits $Q_{1}$ (red) and $Q_{2}$ (blue).
(b) Different cross-Kerr parameters $J_{11}$ (magenta), $J_{21}$ (orange),
$J_{12}$ (green), $J_{22}$ (brown) as a function of $\Phi_{\rm ext}$. The
solid lines are obtained from the numerical calculation, and the experimental
values at the operating flux-bias (dashed line) are shown with plus symbols.

To understand how two-qutrit sidebands work, we follow Ref. 32 and apply
adiabatic approximation to the Hamiltonian. The frequency of the coupler mode
remains well above ($>15$ GHz) any other frequencies ($<4$ GHz) in the system
and thus can be assumed to stay in the ground state during any operation. The
static energy of the coupler mode is removed by minimizing the full circuit
Hamiltonian. Treating qutrits as Duffing-type oscillators, a toy model with
adiabatic approximation can be written as:

$\displaystyle H_{toy}=$
$\displaystyle\omega_{1}a_{1}^{\dagger}a_{1}+\omega_{2}a_{2}^{\dagger}a_{2}+\frac{\alpha_{1}}{2}a_{1}^{\dagger}a_{1}^{\dagger}a_{1}a_{1}+\frac{\alpha_{2}}{2}a_{2}^{\dagger}a_{2}^{\dagger}a_{2}a_{2}$
$\displaystyle+g_{1}\left(t\right)\left(a_{1}^{\dagger}+a_{1}\right)\left(a_{2}^{\dagger}+a_{2}\right)$
$\displaystyle+g_{2}\left(-a_{1}^{\dagger}+a_{1}\right)\left(-a_{2}^{\dagger}+a_{2}\right),$
(13a) $\displaystyle g_{1}\left(t\right)=$
$\displaystyle\frac{\sqrt{E_{J1}E_{J2}}}{2E_{Jc}\cos{\left(2\pi{\Phi_{\rm
ext}\left(t\right)}/\Phi_{0}\right)}}\sqrt{\omega_{1}\omega_{2}},$ (13b)
$\displaystyle g_{2}=$
$\displaystyle\frac{\sqrt{C_{1}C_{2}}}{2C_{12}}\sqrt{\omega_{1}\omega_{2}}.$
(13c)

Here $g_{1}\left(t\right)$ and $g_{2}$ are inductive and capacitive coupling
strengths respectively. $g_{1}\left(t\right)$ is flux tunable, and when the RF
flux modulation frequency matches with any two-photon transition frequency
$\left|jk\right\rangle\leftrightarrow\left|j+1,k-1\right\rangle$ or
$\left|jk\right\rangle\leftrightarrow\left|j+1,k+1\right\rangle$, the
corresponding sideband is activated in the system. This adiabatic
approximation requires sideband frequencies to be much smaller than qutrit
frequencies, which can capture the (red) sideband rates for
$\left|22\right\rangle\leftrightarrow\left|31\right\rangle$ and
$\left|21\right\rangle\leftrightarrow\left|30\right\rangle$ used in the
system. However, the (blue) sideband rate for
$\left|jk\right\rangle\leftrightarrow\left|j+1,k+1\right\rangle$ cannot be
explained with such an approximation and requires keeping the coupler mode in
the analysis [45].

##  Appendix D CPhase gates

Figure 6: CPhase gate construction. (a) Partitioned energy level diagram.
Flipping the phase of a specific target state in region $\left(0\right)$ is
performed through a $2\pi$ sideband rotation (shown in red arrow). Flipping
the other target states requires a decomposition. The application of a
$C_{p}(\pi,\left|mn\right\rangle)$ follows the path marked with brown arrows
starting from $\left|mn\right\rangle$. Here the numbers inside the yellow
circle indicate the total number of single-qutrit $\pi$ pulses required. (b)
An example of circuit decomposition for the
$C_{p}\left(\pi,\left|00\right\rangle\right)$ gate. (c) Single-qutrit phase
compensation calibration for a $C_{p}\left(\pi,\left|jk\right\rangle\right)$
gate.

As shown in Fig. 6(a), we utilize two native CPhase gates
$C_{p}\left(\theta,\left|21\right\rangle\right)$ and
$C_{p}\left(\theta,\left|22\right\rangle\right)$, which can be realized by
applying $2\pi$ rotations to
$\left|21\right\rangle\leftrightarrow\left|30\right\rangle$ and
$\left|22\right\rangle\leftrightarrow\left|31\right\rangle$ transitions with a
phase difference of $\pi-\theta$ between the first and the last $\pi$
rotations. The optimized gate lengths are $55.9$ ns and $94.0$ ns for
$C_{p}\left(\theta,\left|21\right\rangle\right)$ and
$C_{p}\left(\theta,\left|22\right\rangle\right)$ respectively.

The CPhase gates on the states in the region $\left(2\right)$,
$\left(4\right)$ and $\left(6\right)$ (highlighted with yellow) need to
decomposed where the number indicates the number single-qutrit rotations
required. The decomposition starts from the target state, followed by the
application of single-qutrit rotations to arrive at one of the states in the
region $\left(0\right)$ (following brown arrows). After applying the native
CPhase gate, reverse single-qutrit rotations are administered traversing the
same path back to the target state. Fig 6(b) shows the circuit decomposition
of $C_{p}\left(\pi,\left|00\right\rangle\right)$ as an example, which requires
the maximum number of pulses. Note that, even though we have access to
$C_{p}\left(\pi,\left|12\right\rangle\right)$, we do not use it due to poorer
fidelity.

Figure 7: Two-qutrit tomography of the state
$\frac{1}{3}\left(\left|0\right\rangle+\left|1\right\rangle+\left|2\right\rangle\right)^{\otimes
2}$ after applying the CPhase gate
$C_{p}\left(\frac{8\pi}{9},\left|22\right\rangle\right)$. The top and bottom
rows are experimental and ideal density matrices, with the real and imaginary
parts shown in the left and right columns.

The calibration of the CPhase gates involves compensation for the additional
phases acquired during the gate operation [46]. In our coupler, both the AC
stark shift and the rectification effect from SQUID flux modulation can cause
qutrit frequencies to change during the gate, resulting in extra phases
$\beta_{01}^{(j)}$ and $\beta_{12}^{(j)}$ for the
$\left|0\right\rangle\leftrightarrow\left|1\right\rangle$ and
$\left|1\right\rangle\leftrightarrow\left|2\right\rangle$ transitions of the
$j$-th qutrit. We use the circuit shown in Fig. 6(c) to extract these
additional phases. Starting from the ground state, we apply the gate
$G_{1}=R_{01}\left(0,\frac{\pi}{2}\right)$ to both qutrits if calibrating
$\beta_{01}^{(j)}$, followed by $C_{p}\left(0,\left|jk\right\rangle\right)$
and sweep the virtual phase $Z_{\theta}=\Theta(\theta,0)$ before the last gate
$G_{2}=R_{01}\left(\pi,\frac{\pi}{2}\right)$ on $Q_{j}$. By fitting the
readout on $Q_{j}$ to $C_{0}+C_{1}\sin{\left(\beta_{01}^{(j)}+\theta\right)}$,
one can extract the extra phase $\beta_{01}^{(j)}$ acquired. A similar
procedure is used to extract $\beta_{12}^{(j)}$ where we use
$G_{1}=R_{12}\left(0,\frac{\pi}{2}\right)\cdot R_{01}\left(0,\pi\right)$,
$G_{2}=R_{12}\left(\pi,\frac{\pi}{2}\right)$, and
$Z_{\theta}=\Theta(0,\theta)$. Hence, each CPhase gate has four virtual phases
corresponding to two native transitions for each qutrit. As a demonstration,
we apply the CPhase gate
$C_{p}\left(\frac{8\pi}{9},\left|22\right\rangle\right)$ on the state
$\frac{1}{3}\left(\left|0\right\rangle+\left|1\right\rangle+\left|2\right\rangle\right)^{\otimes
2}$ (with initial state fidelity of 87.7%). The final state (shown in Fig. 7)
obtained after two-qutrit tomography [47] shows a fidelity of 82.4%, which is
limited by the cross-Kerr coupling between qutrits and SPAM error.

##  Appendix E Measurement error mitigation

Figure 8: Heatmap of the confusion matrix. $3^{2}$ basis states are prepared
and measured for 20,000 times. The numbers represent average assignment
probabilities.

Figure 8 shows the single-shot assignment probability for the nine basis
states of our two-qutrit processor. We refer this 2d array as the confusion
matrix. Our readout fidelity is limited because there are no parametric
amplifiers on the output lines. To fairly demonstrate the performance of the
qutrit algorithms with single shot results, we apply the inverse of the
confusion matrix to the readout results to compensate for the measurement
error. After correcting for the measurement error, it is possible that some of
the readout counts (number of times the system found in a specific state)
become negative, which happens due to drifts in the calibration parameters.

This is corrected using Maximum-Likelihood-Estimation (MLE) with the
assumption that the minimum fluctuation of a measurement repeated $N$ times
should not be lower than $\sqrt{N}$ (assuming normal distribution). We define
the following cost function to avoid non-physical measurement counts:

$f_{c}(\overrightarrow{p},\overrightarrow{q})=\sum_{i=1}^{9}\left(\frac{p_{i}-q_{i}}{q_{i}}\right)^{2},$
(14)

with the restriction $p_{i}\geq{\sqrt{N}}$. Here $\overrightarrow{q}$ contains
the experiment counts (1d array of 9 elements corresponding to the basis
states) after applying the inverse of the confusion matrix, and
$\overrightarrow{p}$ is the extracted counts after MLE.

##  Appendix F Balanced oracles of DJ algorithm

We implement 16 different balanced functions whose equivalent classical
functions are tabulated in Table 5. We use $A$ and $B$ to represent the
classical ternary values (0,1, and 2) for the two qutrits. The operator
$\oplus$ and $\odot$ correspond to addition and multiplication modulo 3
respectively.

Oracle | Classical function | Oracle | Classical function
---|---|---|---
$Z\otimes I$ | A $\oplus$ 0 | $I\otimes Z$ | 0 $\oplus$ B
$Z\otimes X$ | A $\oplus$ 1 | $X\otimes Z$ | 1 $\oplus$ B
$Z\otimes X^{2}$ | A $\oplus$ 2 | $X^{2}\otimes Z$ | 2 $\oplus$ B
$Z^{2}\otimes I$ | (2 $\odot$ A) $\oplus$ 0 | $I\otimes Z^{2}$ | 0 $\oplus$ (2 $\odot$ B)
$Z^{2}\otimes X$ | (2 $\odot$ A) $\oplus$ 1 | $X\otimes Z^{2}$ | 1 $\oplus$ (2 $\odot$ B)
$Z^{2}\otimes X^{2}$ | (2 $\odot$ A) $\oplus$ 2 | $X^{2}\otimes Z^{2}$ | 2 $\oplus$ (2 $\odot$ B)
$Z\otimes Z$ | A $\oplus$ B | $Z\otimes Z^{2}$ | A $\oplus$ (2 $\odot$ B)
$Z^{2}\otimes Z$ | (2 $\odot$ A) $\oplus$ B | $Z^{2}\otimes Z^{2}$ | (2 $\odot$ A) $\oplus$ (2 $\odot$ B)
Table 5: Equivalent classical functions (ternary valued) for the 16 balanced
functions implemented in the DJ algorithm.

##  Appendix G Error analysis

Figure 9: Master equation simulation of two-qutrit Grover’s search with one
and two stages of amplitude amplification. These results are in a good
agreement with the experimental outcomes shown in Fig. 3(b).

Multiple error sources limiting us from approaching the theoretical success
rates for different quantum algorithms. The dominant one is the dephasing
noise. The deepest circuit (2 stages of Grover’s search) implemented includes
$17$ single-qutrit and $4$ two-qutrit gates, with a total execution time of
$2.11~{}\mu$s. This time becomes comparable to the qutrit’s Ramsey times (see
Table 2) and strongly affects the success rates. The second source is the
inter-qutrit dispersive coupling, including both the static-$ZZ$ values
measured in Table. 3 and the dynamic-$ZZ$ terms induced during (parametric)
gate operations. The $ZZ$ interactions introduce unwanted phase accumulations
that reduce the performance of both single-qutrit and CPhase gates. The third
source is the leakage to non-computational levels. However, our maximum gate
rates are much smaller than the energy gaps, and thus errors due to leakage
should not be significant. In order to verify that we have captured all the
main sources of error, we perform master-equation simulation using the
experimentally obtained device parameters. The simulated results for Grover’s
search are shown in Fig. 9, and those match very well with the averaged square
statistical overlap [38] of $93.27\%$ and $94.46\%$ respectively for one and
two stages of amplitude amplification.

##  Appendix H Scaling up

We envisage building a larger system using qutrit-based units, and there are
several choices. First, one can couple three transmons to the same parametric
coupler [33] to build a three-qutrit processor where inter-qutrit static ZZ
can be simultaneously minimized through an optimized layout ensuring almost
identical inter-transmon capacitances. Next, these units can be used as
building blocks to expand in a linear or planar geometry, as shown in Fig. 10.
Two units could be coupled via commonly used tunable couplers [48, 49, 50]
where a galvanic connection is not required. Any two qutrits between two
neighboring units can be made equally distant (with respect to hopping) when
the capacitive coupling for the inter-unit coupler is made to the central
superconducting islands. This way one can realize a moderately-sized qutrit
processor consisting of highly efficient two or three qutrit units with fast
gates and a slightly slower but switchable interaction between neighboring
pairs with almost identical rates.

Figure 10: Scaling-up schemes using qutrit-based units as building blocks.

##  Appendix I Setup

Figure 11: Detailed cryogenic and room-temperature measurement setup.

The room and cryogenic temperature wiring diagram is shown in Fig 11. The
device is placed inside a bilayer $\mu$-metal shield, mounted on the base
plate of a dilution fridge with 15 mK base temperature. A Tektronix 5014C
(1.2GSa/s) Arbitrary Waveform Generator (AWG) acts as the master trigger for
all other equipment. The readout pulses are generated through two CW tones
from RF sources (PSG-E8257D), modulated by the AWG 5014C. The qutrit pulses
are generated through another 4-channel AWG (Keysight M8195 65 GSa/s, 16 GSa/s
per channel). Before entering the fridge, the qutrit and readout pulses are
combined and sent through lines In1 and In2. One DC source (Yokogawa GS200) is
used to bias the coupler’s DC flux position. The RF flux modulation is
synthesized through the same 4-channel AWG. Two other direct charge drives,
also synthesized through the same AWG, are not used in this project. Inside
the fridge, at the 4K plate all input lines have 20 dB attenuation. At the
base plate, In1 and In2 lines have 10-dB attenuation each, followed by a
strong Eccosorb® providing 20 dB attenuation at 4 GHz; Charge1 and Charge2
lines have 20 dB attenuation, followed by a strong Eccosorb providing 20 dB
attenuation at 4 GHz, and a band pass filter with pass band $3.9-4.8$ GHz. The
DC Flux line has a low pass filter (DC $-~{}1.9$ MHz), followed by a weak
Eccosorb, and the RF Flux line passes through a weak Eccosorb first, followed
by a high-pass filter (with a cut-off at $200$ MHz). The DC and RF flux lines
are combined with a Bias Tee, and then the signal passes through a Step
Impedance Purcell Filter (SIPF), which blocks the frequency band $2.5-5.5$
GHz. The output signal goes through a weak Eccosorb, two circulators, a DC
block, and is amplified with one LNF HEMT amplifier at 4K. The output signals
are further amplified at room temperature and then demodulated with IQ mixers.
The demodulated signals pass through low-pass filters (DC $-~{}250$ MHz) and
are amplified again using the SRS Preamplifier. The final signals are
digitized with Alazar ATS 9870 (1 GSa/s) and analyzed in a computer.

## References

  * Preskill [2018] J. Preskill, Quantum computing in the nisq era and beyond, Quantum 2, 79 (2018).
  * Arute  _et al._ [2019] F. Arute, K. Arya, R. Babbush, D. Bacon, J. C. Bardin, R. Barends, R. Biswas, S. Boixo, F. G. Brandao, D. A. Buell, _et al._ , Quantum supremacy using a programmable superconducting processor, Nature 574, 505 (2019).
  * Wu  _et al._ [2021] Y. Wu, W.-S. Bao, S. Cao, F. Chen, M.-C. Chen, X. Chen, T.-H. Chung, H. Deng, Y. Du, D. Fan, M. Gong, C. Guo, C. Guo, S. Guo, L. Han, L. Hong, H.-L. Huang, Y.-H. Huo, L. Li, N. Li, S. Li, Y. Li, F. Liang, C. Lin, J. Lin, H. Qian, D. Qiao, H. Rong, H. Su, L. Sun, L. Wang, S. Wang, D. Wu, Y. Xu, K. Yan, W. Yang, Y. Yang, Y. Ye, J. Yin, C. Ying, J. Yu, C. Zha, C. Zhang, H. Zhang, K. Zhang, Y. Zhang, H. Zhao, Y. Zhao, L. Zhou, Q. Zhu, C.-Y. Lu, C.-Z. Peng, X. Zhu, and J.-W. Pan, Strong quantum computational advantage using a superconducting quantum processor, Phys. Rev. Lett. 127, 180501 (2021).
  * Bullock  _et al._ [2005] S. S. Bullock, D. P. O’Leary, and G. K. Brennen, Asymptotically optimal quantum circuits for $d$-level systems, Phys. Rev. Lett. 94, 230502 (2005).
  * Ralph  _et al._ [2007] T. C. Ralph, K. J. Resch, and A. Gilchrist, Efficient toffoli gates using qudits, Phys. Rev. A 75, 022313 (2007).
  * Lanyon  _et al._ [2009] B. P. Lanyon, M. Barbieri, M. P. Almeida, T. Jennewein, T. C. Ralph, K. J. Resch, G. J. Pryde, J. L. O’brien, A. Gilchrist, and A. G. White, Simplifying quantum logic using higher-dimensional hilbert spaces, Nature Physics 5, 134 (2009).
  * Gedik  _et al._ [2015] Z. Gedik, I. A. Silva, B. Çakmak, G. Karpat, E. L. G. Vidoto, D. O. Soares-Pinto, E. R. deAzevedo, and F. F. Fanchini, Computational speed-up with a single qudit, Scientific reports 5, 1 (2015).
  * Nikolaeva  _et al._ [2021] A. Nikolaeva, E. Kiktenko, and A. Fedorov, Efficient realization of quantum algorithms with qudits, arXiv preprint arXiv:2111.04384 (2021).
  * Kapit [2016] E. Kapit, Hardware-efficient and fully autonomous quantum error correction in superconducting circuits, Phys. Rev. Lett. 116, 150501 (2016).
  * Muralidharan  _et al._ [2017] S. Muralidharan, C.-L. Zou, L. Li, J. Wen, and L. Jiang, Overcoming erasure errors with multilevel systems, New Journal of Physics 19, 013026 (2017).
  * Majumdar  _et al._ [2018] R. Majumdar, S. Basu, S. Ghosh, and S. Sur-Kolay, Quantum error-correcting code for ternary logic, Phys. Rev. A 97, 052302 (2018).
  * Campbell [2014] E. T. Campbell, Enhanced fault-tolerant quantum computing in $d$-level systems, Phys. Rev. Lett. 113, 230501 (2014).
  * Campbell  _et al._ [2012] E. T. Campbell, H. Anwar, and D. E. Browne, Magic-state distillation in all prime dimensions using quantum reed-muller codes, Phys. Rev. X 2, 041021 (2012).
  * Gokhale  _et al._ [2019] P. Gokhale, J. M. Baker, C. Duckering, N. C. Brown, K. R. Brown, and F. T. Chong, Asymptotic improvements to quantum circuits via qutrits, in  _Proceedings of the 46th International Symposium on Computer Architecture_, ISCA ’19 (Association for Computing Machinery, New York, NY, USA, 2019) p. 554–566.
  * Gustafson [2022] E. Gustafson, Noise improvements in quantum simulations of sqed using qutrits, arXiv preprint arXiv:2201.04546 (2022).
  * Bechmann-Pasquinucci and Peres [2000] H. Bechmann-Pasquinucci and A. Peres, Quantum cryptography with 3-state systems, Phys. Rev. Lett. 85, 3313 (2000).
  * Bruß and Macchiavello [2002] D. Bruß and C. Macchiavello, Optimal eavesdropping in cryptography with three-dimensional quantum states, Phys. Rev. Lett. 88, 127901 (2002).
  * Vaziri  _et al._ [2002] A. Vaziri, G. Weihs, and A. Zeilinger, Experimental two-photon, three-dimensional entanglement for quantum communication, Phys. Rev. Lett. 89, 240401 (2002).
  * Luo  _et al._ [2019] Y.-H. Luo, H.-S. Zhong, M. Erhard, X.-L. Wang, L.-C. Peng, M. Krenn, X. Jiang, L. Li, N.-L. Liu, C.-Y. Lu, A. Zeilinger, and J.-W. Pan, Quantum teleportation in high dimensions, Phys. Rev. Lett. 123, 070505 (2019).
  * Hu  _et al._ [2020] X.-M. Hu, C. Zhang, B.-H. Liu, Y. Cai, X.-J. Ye, Y. Guo, W.-B. Xing, C.-X. Huang, Y.-F. Huang, C.-F. Li, and G.-C. Guo, Experimental high-dimensional quantum teleportation, Phys. Rev. Lett. 125, 230501 (2020).
  * Ringbauer  _et al._ [2021] M. Ringbauer, M. Meth, L. Postler, R. Stricker, R. Blatt, P. Schindler, and T. Monz, A universal qudit quantum processor with trapped ions, arXiv preprint arXiv:2109.06903 (2021).
  * Hrmo  _et al._ [2022] P. Hrmo, B. Wilhelm, L. Gerster, M. W. van Mourik, M. Huber, R. Blatt, P. Schindler, T. Monz, and M. Ringbauer, Native qudit entanglement in a trapped ion quantum processor, arXiv preprint arXiv:2206.04104 (2022).
  * Morvan  _et al._ [2021] A. Morvan, V. V. Ramasesh, M. S. Blok, J. M. Kreikebaum, K. O’Brien, L. Chen, B. K. Mitchell, R. K. Naik, D. I. Santiago, and I. Siddiqi, Qutrit randomized benchmarking, Phys. Rev. Lett. 126, 210504 (2021).
  * Goss  _et al._ [2022] N. Goss, A. Morvan, B. Marinelli, B. K. Mitchell, L. B. Nguyen, R. K. Nail, L. Chen, C. Jünger, J. M. Kreikebaum, D. I. Santiago, _et al._ , High-fidelity qutrit entangling gates for superconducting circuits, arXiv preprint arXiv:2206.07216 (2022).
  * Luo  _et al._ [2022] K. Luo, W. Huang, Z. Tao, L. Zhang, Y. Zhou, J. Chu, W. Liu, B. Wang, J. Cui, S. Liu, _et al._ , Experimental realization of two qutrits gate with tunable coupling in superconducting circuits, arXiv preprint arXiv:2206.11199 (2022).
  * Blais  _et al._ [2021] A. Blais, A. L. Grimsmo, S. M. Girvin, and A. Wallraff, Circuit quantum electrodynamics, Rev. Mod. Phys. 93, 025005 (2021).
  * Koch  _et al._ [2007] J. Koch, T. M. Yu, J. Gambetta, A. A. Houck, D. I. Schuster, J. Majer, A. Blais, M. H. Devoret, S. M. Girvin, and R. J. Schoelkopf, Charge-insensitive qubit design derived from the Cooper pair box, Physical Review A 76, 042319 (2007).
  * Place  _et al._ [2021] A. P. Place, L. V. Rodgers, P. Mundada, B. M. Smitham, M. Fitzpatrick, Z. Leng, A. Premkumar, J. Bryon, A. Vrajitoarea, S. Sussman, _et al._ , New material platform for superconducting transmon qubits with coherence times exceeding 0.3 milliseconds, Nature communications 12, 1 (2021).
  * Wang  _et al._ [2022] C. Wang, X. Li, H. Xu, Z. Li, J. Wang, Z. Yang, Z. Mi, X. Liang, T. Su, C. Yang, _et al._ , Towards practical quantum computers: transmon qubit with a lifetime approaching 0.5 milliseconds, npj Quantum Information 8, 1 (2022).
  * Bianchetti  _et al._ [2010] R. Bianchetti, S. Filipp, M. Baur, J. M. Fink, C. Lang, L. Steffen, M. Boissonneault, A. Blais, and A. Wallraff, Control and tomography of a three level superconducting artificial atom, Phys. Rev. Lett. 105, 223601 (2010).
  * Blok  _et al._ [2021] M. S. Blok, V. V. Ramasesh, T. Schuster, K. O’Brien, J. M. Kreikebaum, D. Dahlen, A. Morvan, B. Yoshida, N. Y. Yao, and I. Siddiqi, Quantum information scrambling on a superconducting qutrit processor, Phys. Rev. X 11, 021010 (2021).
  * Lu  _et al._ [2017] Y. Lu, S. Chakram, N. Leung, N. Earnest, R. K. Naik, Z. Huang, P. Groszkowski, E. Kapit, J. Koch, and D. I. Schuster, Universal stabilization of a parametrically coupled qubit, Phys. Rev. Lett. 119, 150502 (2017).
  * Brown  _et al._ [2022] T. Brown, E. Doucet, D. Ristè, G. Ribeill, K. Cicak, J. Aumentado, R. Simmonds, L. Govia, A. Kamal, and L. Ranzani, Trade off-free entanglement stabilization in a superconducting qutrit-qubit system, Nature communications 13, 1 (2022).
  * Deutsch and Jozsa [1992] D. Deutsch and R. Jozsa, Rapid solution of problems by quantum computation, Proceedings of the Royal Society of London A: Mathematical, Physical and Engineering Sciences 439, 553 (1992).
  * Wang  _et al._ [2020] Y. Wang, Z. Hu, B. C. Sanders, and S. Kais, Qudits and high-dimensional quantum computing, Frontiers in Physics 8, 589504 (2020).
  * Bernstein and Vazirani [1997] E. Bernstein and U. Vazirani, Quantum complexity theory, SIAM Journal on Computing 26, 1411 (1997).
  * Grover [1997] L. K. Grover, Quantum mechanics helps in searching for a needle in a haystack, Phys. Rev. Lett. 79, 325 (1997).
  * Roy  _et al._ [2020] T. Roy, S. Hazra, S. Kundu, M. Chand, M. P. Patankar, and R. Vijay, Programmable superconducting processor with native three-qubit gates, Phys. Rev. Applied 14, 014072 (2020).
  * Figgatt  _et al._ [2017] C. Figgatt, D. Maslov, K. A. Landsman, N. M. Linke, S. Debnath, and C. Monroe, Complete 3-qubit grover search on a programmable quantum computer, Nature Communications 8, 1918 (2017).
  * Mandviwalla  _et al._ [2018] A. Mandviwalla, K. Ohshiro, and B. Ji, Implementing grover’s algorithm on the ibm quantum computers, in  _2018 IEEE International Conference on Big Data (Big Data)_ (2018) pp. 2531–2537.
  * Zhang  _et al._ [2021] K. Zhang, P. Rao, K. Yu, H. Lim, and V. Korepin, Implementation of efficient quantum search algorithms on nisq computers, Quantum Information Processing 20, 1 (2021).
  * Roy  _et al._ [2022] T. Roy, L. Jiang, and D. I. Schuster, Deterministic grover search with a restricted oracle, Phys. Rev. Research 4, L022013 (2022).
  * Long [2001] G. L. Long, Grover algorithm with zero theoretical failure rate, Phys. Rev. A 64, 022307 (2001).
  * Werninghaus  _et al._ [2021] M. Werninghaus, D. J. Egger, F. Roy, S. Machnes, F. K. Wilhelm, and S. Filipp, Leakage reduction in fast superconducting qubit gates via optimal control, npj Quantum Information 7, 1 (2021).
  * Roth  _et al._ [2017] M. Roth, M. Ganzhorn, N. Moll, S. Filipp, G. Salis, and S. Schmidt, Analysis of a parametrically driven exchange-type gate and a two-photon excitation gate between superconducting qubits, Phys. Rev. A 96, 062323 (2017).
  * Ganzhorn  _et al._ [2020] M. Ganzhorn, G. Salis, D. J. Egger, A. Fuhrer, M. Mergenthaler, C. Müller, P. Müller, S. Paredes, M. Pechal, M. Werninghaus, and S. Filipp, Benchmarking the noise sensitivity of different parametric two-qubit gates in a single superconducting quantum computing platform, Phys. Rev. Research 2, 033447 (2020).
  * Roy  _et al._ [2021] T. Roy, Z. Li, E. Kapit, and D. I. Schuster, Tomography in the presence of stray inter-qubit coupling, arXiv preprint arXiv:2103.13611 (2021).
  * McKay  _et al._ [2016] D. C. McKay, S. Filipp, A. Mezzacapo, E. Magesan, J. M. Chow, and J. M. Gambetta, Universal gate for fixed-frequency qubits via a tunable bus, Phys. Rev. Applied 6, 064007 (2016).
  * Sung  _et al._ [2021] Y. Sung, L. Ding, J. Braumüller, A. Vepsäläinen, B. Kannan, M. Kjaergaard, A. Greene, G. O. Samach, C. McNally, D. Kim, A. Melville, B. M. Niedzielski, M. E. Schwartz, J. L. Yoder, T. P. Orlando, S. Gustavsson, and W. D. Oliver, Realization of high-fidelity cz and $zz$-free iswap gates with a tunable coupler, Phys. Rev. X 11, 021058 (2021).
  * Sete  _et al._ [2021] E. A. Sete, N. Didier, A. Q. Chen, S. Kulshreshtha, R. Manenti, and S. Poletto, Parametric-resonance entangling gates with a tunable coupler, Phys. Rev. Applied 16, 024050 (2021).