††thanks: These authors contributed equally to this work.††thanks: These
authors contributed equally to this work.††thanks: These authors contributed
equally to this work.

# Supplementary Materials:
Exploring Hilbert-Space Fragmentation on a Superconducting Processor

Yong-Yi Wang  Institute of Physics, Chinese Academy of Sciences, Beijing
100190, China  School of Physical Sciences, University of Chinese Academy of
Sciences, Beijing 100049, China  Yun-Hao Shi  Institute of Physics, Chinese
Academy of Sciences, Beijing 100190, China  School of Physical Sciences,
University of Chinese Academy of Sciences, Beijing 100049, China  Beijing
Academy of Quantum Information Sciences, Beijing 100193, China  Zheng-Hang Sun
Institute of Physics, Chinese Academy of Sciences, Beijing 100190, China
Theoretical Physics III, Center for Electronic Correlations and Magnetism,
Institute of Physics, University of Augsburg, D-86135 Augsburg, Germany  Chi-
Tong Chen  Institute of Physics, Chinese Academy of Sciences, Beijing 100190,
China  School of Physical Sciences, University of Chinese Academy of Sciences,
Beijing 100049, China  Zheng-An Wang  Beijing Academy of Quantum Information
Sciences, Beijing 100193, China  Hefei National Laboratory, Hefei 230088,
China  Kui Zhao  Beijing Academy of Quantum Information Sciences, Beijing
100193, China  Hao-Tian Liu  Institute of Physics, Chinese Academy of
Sciences, Beijing 100190, China  School of Physical Sciences, University of
Chinese Academy of Sciences, Beijing 100049, China  Wei-Guo Ma  Institute of
Physics, Chinese Academy of Sciences, Beijing 100190, China  School of
Physical Sciences, University of Chinese Academy of Sciences, Beijing 100049,
China  Ziting Wang  Beijing Academy of Quantum Information Sciences, Beijing
100193, China  Hao Li  Beijing Academy of Quantum Information Sciences,
Beijing 100193, China  Jia-Chi Zhang  Institute of Physics, Chinese Academy of
Sciences, Beijing 100190, China  School of Physical Sciences, University of
Chinese Academy of Sciences, Beijing 100049, China  Yu Liu  Institute of
Physics, Chinese Academy of Sciences, Beijing 100190, China  School of
Physical Sciences, University of Chinese Academy of Sciences, Beijing 100049,
China  Cheng-Lin Deng  Institute of Physics, Chinese Academy of Sciences,
Beijing 100190, China  School of Physical Sciences, University of Chinese
Academy of Sciences, Beijing 100049, China  Tian-Ming Li  Institute of
Physics, Chinese Academy of Sciences, Beijing 100190, China  School of
Physical Sciences, University of Chinese Academy of Sciences, Beijing 100049,
China  Yang He  Institute of Physics, Chinese Academy of Sciences, Beijing
100190, China  School of Physical Sciences, University of Chinese Academy of
Sciences, Beijing 100049, China  Zheng-He Liu  Institute of Physics, Chinese
Academy of Sciences, Beijing 100190, China  School of Physical Sciences,
University of Chinese Academy of Sciences, Beijing 100049, China  Zhen-Yu Peng
Institute of Physics, Chinese Academy of Sciences, Beijing 100190, China
School of Physical Sciences, University of Chinese Academy of Sciences,
Beijing 100049, China  Xiaohui Song  Institute of Physics, Chinese Academy of
Sciences, Beijing 100190, China  School of Physical Sciences, University of
Chinese Academy of Sciences, Beijing 100049, China  Guangming Xue  Beijing
Academy of Quantum Information Sciences, Beijing 100193, China  Haifeng Yu
Beijing Academy of Quantum Information Sciences, Beijing 100193, China
Kaixuan Huang<EMAIL_ADDRESS>Institute of Physics, Chinese Academy of
Sciences, Beijing 100190, China  Beijing Academy of Quantum Information
Sciences, Beijing 100193, China  Zhongcheng Xiang<EMAIL_ADDRESS>Institute of Physics, Chinese Academy of Sciences, Beijing 100190, China
School of Physical Sciences, University of Chinese Academy of Sciences,
Beijing 100049, China  Dongning Zheng  Institute of Physics, Chinese Academy
of Sciences, Beijing 100190, China  School of Physical Sciences, University of
Chinese Academy of Sciences, Beijing 100049, China  Songshan Lake Materials
Laboratory, Dongguan, Guangdong 523808, China  CAS Center for Excellence in
Topological Quantum Computation, UCAS, Beijing 100190, China, and Mozi
Laboratory, Zhengzhou 450001, China  Kai Xu<EMAIL_ADDRESS>Institute of
Physics, Chinese Academy of Sciences, Beijing 100190, China  School of
Physical Sciences, University of Chinese Academy of Sciences, Beijing 100049,
China  Beijing Academy of Quantum Information Sciences, Beijing 100193, China
Songshan Lake Materials Laboratory, Dongguan, Guangdong 523808, China  CAS
Center for Excellence in Topological Quantum Computation, UCAS, Beijing
100190, China, and Mozi Laboratory, Zhengzhou 450001, China  Heng Fan
<EMAIL_ADDRESS>Institute of Physics, Chinese Academy of Sciences, Beijing
100190, China  School of Physical Sciences, University of Chinese Academy of
Sciences, Beijing 100049, China  Beijing Academy of Quantum Information
Sciences, Beijing 100193, China  Songshan Lake Materials Laboratory, Dongguan,
Guangdong 523808, China  CAS Center for Excellence in Topological Quantum
Computation, UCAS, Beijing 100190, China, and Mozi Laboratory, Zhengzhou
450001, China

###### Contents

  1. I Frequency Arrangement
  2. II Multi-Qubit Calibration
    1. II.1 Z Crosstalk Correction
    2. II.2 Timing Calibration
    3. II.3 On-Site Potential
  3. III Model and Hamiltonian
    1. III.1 Jordan-Wigner Transformation for the $XX$-Ladder Model
    2. III.2 Level Statistics
    3. III.3 Participation entropy for Gaussian orthogonal ensemble and Gaussian unitary ensemble
      1. III.3.1 Derivation
      2. III.3.2 Numerical verification
    4. III.4 Statistical Properties of Eigenstates in Stark Systems
      1. III.4.1 Participation entropy
      2. III.4.2 Entanglement entropy
      3. III.4.3 Domain wall number
  4. IV Numerical Krylov Subspace
  5. V Time Evolution with Krylov Space Method
  6. VI Post-Selection of Measurement Results

Figure S1: Coupling strengths of the 24 qubits used in the ladder-type chip.

##  I Frequency Arrangement

Before multi-qubit experiments, we also need to arrange the idle frequency of
each single qubit to ensure a long coherence time and a low level of XY
crosstalk. We optimize the frequency arrangement by considering the two-level
system (TLS), the overlap of energy levels, and the residual coupling between
qubits.

  * •

To avoid TLS, we can measure the energy relaxation time $T_{1}$ near the idle
point. If the TLS is coupled to the qubit near the idle point, then the curve
of $T_{1}$ will drop sharply or oscillate. Besides, the dephasing time $T_{2}$
also needs to be considered.

  * •

If the qubit at the idle point has the transition frequency of
$\omega_{\text{10}}$ from $\ket{0}$ to $\ket{1}$ and the transition frequency
of $\omega_{\text{21}}$ from $\ket{1}$ to $\ket{2}$, we should avoid
overlapping $\omega_{\text{10}}$, $\omega_{\text{21}}$, and the two-photons
frequency $(\omega_{\text{10}}+\omega_{\text{21}})/2$. In addition, the
leakage frequencies when we excite the qubit cannot overlap with these
frequencies. Specifically, the leakage frequencies include $\omega_{\text{c}}$
(intrinsic leakage) and $2\omega_{\text{c}}-\omega_{\text{10}}$ (mirror
leakage), where $\omega_{\text{c}}$ denotes the carrier frequency.

  * •

The XY crosstalk includes classical microwave crosstalk and quantum crosstalk.
The classical microwave crosstalk can be effectively suppressed by
appropriately increasing the duration of the excitation microwave pulse
(reducing the broadening in the frequency domain). In our experiment, we set
the duration of $\pi$ pulse to $120~{}\mathrm{ns}$. However, the quantum
crosstalk resulting from the residual coupling between qubits can only be
reduced by changing the frequencies and increasing the qubit frequency
difference.

Parameter | Median | Mean | Stdev. | Units
---|---|---|---|---
Qubit maximum frequency | $5.025$ | $5.032$ | $0.240$ | GHz
Qubit idle frequency | $4.723$ | $4.728$ | $0.346$ | GHz
Qubit anharmonicity | $-0.222$ | $-0.222$ | $0.022$ | GHz
Readout frequency | $6.715$ | $6.714$ | $0.061$ | GHz
Mean energy relaxation time $\overline{T}_{1}$ | $33.2$ | $32.1$ | $7.5$ |  $\mu$s
Pure dephasing time at idle frequency $T_{2}^{\ast}$ | $1.0$ | $2.4$ | $4.2$ |  $\mu$s

Table S1: List of device parameters. Reproduced from Ref. [1]. Figure S2:
Typical experimental data of measuring Z crosstalk. (a) Pulse sequences for
measuring the Z crosstalk coefficient of $Q_{j}$ to $Q_{i}$, i.e., $M_{i,j}$.
Here we fix the delay at $1.5~{}\mu\mathrm{s}$ to suppress the broadening of
data. (b) Heatmap showing the probability of $Q_{i}$ in $\ket{1}$ before and
after the Z crosstalk correction. The red dash line is the result of linear
fitting to obtain $M_{i,j}$. Figure S3: Off-diagonal Z crosstalk matrix. The
heatmap shows the Z crosstalk coefficient $M_{i,j}$ ($i\neq{j}$). The size of
the colored bubble represents $|M_{i,j}|$. Here we show the crosstalk
coefficients of the 24 qubits used in the ladder-type chip.

##  II Multi-Qubit Calibration

After single-qubit calibration, we implement multi-qubit calibration to
correct the crosstalk, align the timing, and find the control parameters of
qubits to achieve the target Hamiltonian with high precision. The reference
frequency is fixed at $4.534~{}\mathrm{GHz}$, and the on-site potential
manipulated in the experiment is the difference between the qubit frequency
and this reference frequency.

###  II.1 Z Crosstalk Correction

The Z crosstalk results from the fact that the low-frequency Z bias signal is
not completely localized to a single qubit. The individual Z bias signal of
each qubit has a spatial distribution across the chip, but the strength decays
with the distance between qubits. Assuming that the Z pulse amplitude (Zpa) of
the $j$-th qubit $Q_{j}$ is $z_{j}$, and the vertical distance between its Z
control line and the $i$-th qubit $Q_{i}$ is $r_{i,j}$, then the magnetic
induction strength felt by $Q_{i}$ from the Z line of the $Q_{j}$ can be
expressed as $B_{{i}\leftarrow{j}}\propto{z_{j}/r_{i,j}}$. The corresponding
crosstalk flux thus is
$\Phi_{{i}\leftarrow{j}}=B_{{i}\leftarrow{j}}S_{i}=c_{{i}\leftarrow{j}}z_{j}$,
where $S_{i}$ represents the loop area of $Q_{i}$’s SQUID and
$c_{{i}\leftarrow{j}}\propto{S_{i}/r_{i,j}}$ denotes the flux crosstalk per
unit Zpa. To compensate for the crosstalk $\Phi_{{i}\leftarrow{j}}$, we apply
$\Phi_{{i}\leftarrow{i}}=c_{{i}\leftarrow{i}}z_{i}$ on the Z line of $Q_{i}$
to satisfy

$\Phi_{{i}\leftarrow{i}}=-\Phi_{{i}\leftarrow{j}},$ (S1)

thus the Z crosstalk coefficient of $Q_{j}$ to $Q_{i}$ can be calculated as

$M_{i,j}=\frac{c_{{i}\leftarrow{j}}}{c_{{i}\leftarrow{i}}}=-\frac{z_{i}}{z_{j}}.$
(S2)

For each $z_{j}$, we scan the Zpa of $Q_{i}$ to find the corresponding $z_{i}$
for compensation. Hence, the Z crosstalk coefficient $M_{i,j}$ can be
determined by linear fitting. As shown in Fig. S2(a), we apply an XY drive
pulse to the target qubit $Q_{i}$ while scanning its Zpa to compensate for the
crosstalk from the Zpa of the source $Q_{j}$. If the crosstalk is compensated
exactly, then we will measure a high probability of $Q_{i}$ in $\ket{1}$.
Typical experimental data are shown in Fig. S2(c) and S2(d). Note that the
above Z crosstalk is linear and does not take into account the anti-cross
effect of energy levels when qubits are close to each other. When measuring
such linear crosstalk, it is preferable to adjust the step of the scanned
source Zpa so that the frequency change corresponding to each scanning step is
greater than twice the coupling strength. Thus, the impact of the anti-cross
effect on the measured data will be reduced.

After measuring the Z crosstalk coefficient for each pair of qubits, we obtain
the whole crosstalk matrix $\mathbf{M}$ (Fig. S3). In the subsequent multi-
qubit calibration and experiments, we will implement the following correction
to the Z pulses of all the qubits:

$\mathbf{z}_{\text{corrected}}=\mathbf{M}^{-1}\mathbf{z},$ (S3)

where $\mathbf{M}^{-1}$ is the inverse of $\mathbf{M}$,
$\mathbf{z}=[z_{1},z_{2},\dots]^{T}$ denotes the vector of all the Zpas, and
$\mathbf{z}_{\text{corrected}}$ represents the corrected vector.

Figure S4: Typical experimental data of timing calibration. (a) Pulse
sequences of two-qubit swap experiment for timing calibration. (b) Joint
probabilities $P_{\ket{01}}$ and $P_{\ket{10}}$ as the function of the Zpa of
$Q_{j}$. By fitting this data with the cosine function, we obtain the swap
period $T$. Then we determine the delay $t_{\text{d}}$ slightly less than
$T/4$. (c) Joint probabilities as the function of delay. (d) Joint
probabilities as the function of time shift. The optimal time shift
$t_{\text{s}}$ applied to $Q_{j}$ is calculated as the average of the Gaussian
fitting results.

###  II.2 Timing Calibration

In single-qubit timing calibration, we have calibrated the timing between each
qubit’s own XY and Z signals. However, different qubits have different control
line lengths and this causes a timing shift between qubits. Here we introduce
an efficient method to automatically calibrate the timing between qubits (Fig.
S4):

  * •

Find the Zpa corresponding to the two-qubit resonance: We initially prepare
$\ket{10}$ by applying a $\pi$ pulse to $Q_{i}$, then vary the Zpa of $Q_{j}$
(the Zpa of $Q_{i}$ is fixed near the resonance point) and measure the joint
probabilities $P_{\ket{01}}$ and $P_{\ket{10}}$, see Fig. S4(b). For these
probabilities to change significantly as the Zpa of $Q_{j}$ approaches the
resonance point, we choose the delay to be approximately half of the swap
period, i.e., $t_{\text{d}}\approx{T/2}$, where $T=\pi/g$ is the swap period
and $g$ denotes the coupling strength of qubits. Here the coupling strength
only needs to be given a rough value and does not need to be finely measured.

  * •

Determine a sensitive delay time of the swap: With the Zpa of $Q_{j}$
corresponding to the resonance, we measure $P_{\ket{01}}$ and $P_{\ket{10}}$
as the function of delay, see Fig. S4(c). By fitting this data with the cosine
function, we can obtain the exact swap period $T$. The sensitive delay time is
near where the derivative is at its maximum ($\approx T/4$). Here we choose a
delay slightly less than $T/4$ for the next step.

  * •

Optimize the relative time shift between two qubits: We measure $P_{\ket{01}}$
and $P_{\ket{10}}$ as the function of $Q_{j}$’s time shift and fit them with
Gaussian functions, see Fig. S4(d). The optimal time shift is at the peak (or
dip) of the Gaussian function.

###  II.3 On-Site Potential

Our experiment requires a high-precision control of on-site disorder and Stark
potential. Even though we have corrected the linear Z crosstalk and aligned
the timing, the qubit (dressed-)frequency will still deviate from the target
frequency when other qubits are simultaneously biased. This effect mainly
comes from the interaction with other qubits (if the residual Z crosstalk is
relatively negligible). Here we introduce a method to suppress this effect and
calibrate qubit frequencies with a high precision.

As shown in Fig. S5, there are two configurations for the calibration, namely
A and B. Configuration A requires that the nearest-neighbor (NN) qubits, the
next-nearest-neighbor (NNN) qubits, and the rest qubits biased from the target
frequency of $+\Delta_{\text{NN}}$, $+\Delta_{\text{NNN}}$, and
$+\Delta_{\text{rest}}$, respectively. Thus, the dressed frequency of the
target qubit $Q_{i}$ can be approximately expressed as

$\tilde{\omega}_{i}^{\text{A}}\approx\omega_{i}+\sum_{j\neq{i}}\frac{g^{2}_{i,j}}{\Delta_{j}},$
(S4)

where $\omega_{i}$ denotes the target frequency, $g_{i,j}$ is the coupling
strength, and $\Delta_{j}=\omega_{i}-\omega_{j}$ is the frequency detuning.
Similarly, for configuration B with the opposite detuning, we have

$\tilde{\omega}_{i}^{\text{B}}\approx\omega_{i}-\sum_{j\neq{i}}\frac{g^{2}_{i,j}}{\Delta_{j}}.$
(S5)

We perform the Rabi oscillation experiments under these two configurations,
and thus determine the corresponding Zpas of $\tilde{\omega}_{i}^{\text{A}}$
and $\tilde{\omega}_{i}^{\text{B}}$, i.e., $z_{i}^{\text{A}}$ and
$z_{i}^{\text{B}}$ (Fig. S6).

To ensure a large adjustment range of Zpa, we usually set the target frequency
below the maximum frequency of the qubit with the smallest maximum frequency.
The relationship between the qubit’s Zpa and frequency near such target
frequencies is monotonic. Hence, the corresponding Zpa of the target
$\omega_{i}$ can be approximately calculated as the average of
$z_{i}^{\text{A}}$ and $z_{i}^{\text{B}}$:

$z_{i}=f(\omega_{i})\approx\frac{1}{2}\left[f(\tilde{\omega}_{i}^{\text{A}})+f(\tilde{\omega}_{i}^{\text{B}})\right]=\frac{1}{2}(z_{i}^{\text{A}}+z_{i}^{\text{B}}),$
(S6)

where $f$ is a monotonic function that represents the relationship between
qubit’s Zpa and frequency near the target frequency. For each qubit, we repeat
the above calibration procedure and finally obtain the Zpas of all qubits
corresponding to the target on-site potentials. The schematic of experimental
pulse sequences is shown in Fig. S7.

Figure S5: Configurations for multi-qubit calibration. As an example, here we
take $Q_{7,\uparrow}$ as the target qubit. The frequency differences are
$\Delta_{\text{NN}}/2\pi\approx 160~{}\mathrm{MHz}$,
$\Delta_{\text{NN}}/2\pi\approx 50~{}\mathrm{MHz}$ and
$\Delta_{\text{rest}}/2\pi\approx 25~{}\mathrm{MHz}$. All unused qubits
($Q_{13,\uparrow}$, $Q_{13,\downarrow}$, $Q_{14,\uparrow}$,
$Q_{14,\downarrow}$, $Q_{15,\uparrow}$, $Q_{15,\downarrow}$) are biased away
from the reference frequency. Figure S6: Typical experimental data of
calibrating on-site potential. The results of Rabi oscillation experiments
under configurations A (a) and B (b). The heatmap shows the probability of the
target qubit in $\ket{1}$. (c) and (d) correspond to (a) and (b),
respectively, representing the Rabi frequency as a function of Zpa. Figure S7:
Schematic of experimental pulse sequences for different initial states. Here
we take the case of $L=8$ as an example, where the initial states
$\ket{\psi_{n_{\text{DW}}=2}}=\ket{\mathbb{11000011}}$ and
$\ket{\psi_{n_{\text{DW}}=6}}=\ket{\mathbb{10100101}}$ correspond to (a) and
(b), respectively. The on-site potential $W_{jm}$ is denoted as the difference
between the qubit frequency $\omega_{jm}$ and the reference frequency
$\omega_{\text{ref}}$. By controlling the Zpas of qubits, we manipulate
$\omega_{jm}$ to realize the Stark (linear) potentials
$W_{jm}/\bar{J}_{\text{NN}}=-j\gamma$, or random potentials
$W_{jm}/\bar{J}_{\text{NN}}$ with a uniform distribution in $[-W,W]$.

##  III Model and Hamiltonian

###  III.1 Jordan-Wigner Transformation for the $XX$-Ladder Model

The main interacting term of the Hamiltonian describing our superconducting
processor can be approximated as the $XX$-ladder model. The approximating
ladder-XX model can be mapped to an interacting spinless fermion model,
utilizing a generalized form of Jordan-Wigner transformation, described as

$\hat{S}_{j,1}^{-}=\hat{c}_{j,1}\exp\left[i\pi\sum_{l=0}^{j-1}\left(\hat{n}_{l,1}+\hat{n}_{l,2}\right)\right],$
(S7)

for the first row of chain, and

$\hat{S}_{j,2}^{-}=\hat{c}_{j,2}\exp\left[i\pi\left[\sum_{l=0}^{j}\hat{n}_{l,1}+\sum_{l=0}^{j-1}\hat{n}_{l,2}\right]\right],$
(S8)

for the second row of chain, where
$\hat{S}_{j,m}^{-}=(\hat{\sigma}_{j,m}^{x}-i\hat{\sigma}_{j,m}^{y})/2$
($j=1,2,\ldots,L$ and $m=1,2$) and
$\hat{n}_{j,m}=\hat{c}^{\dagger}_{j,m}\hat{c}_{j,m}$ is the spinless fermion
number operator for site $(j,m)$. Therefore, the Hamiltonian in the spinless
fermion representation reads

$\displaystyle\hat{H}=-2J\sum_{j,\delta}[$
$\displaystyle\hat{c}_{j,1}\exp\left(-i\pi\hat{\phi}_{j,j+\delta}^{(1)}\right)\hat{c}_{j+\delta,1}^{\dagger}$
(S9) $\displaystyle+$
$\displaystyle\hat{c}_{j,1}\exp\left(-i\pi\hat{\phi}_{j,j+\delta}^{(2)}\right)\hat{c}_{j+\delta,2}^{\dagger}]$
$\displaystyle-$ $\displaystyle
2J\sum_{j}\hat{c}_{j,1}\hat{c}_{j,2}^{\dagger},$

where $\delta$ denotes the nearest neighbor, and
$\hat{\phi}_{j,j+1}^{(1)}=\hat{n}_{j,2},\hat{\phi}_{j,j-1}^{(1)}=-\hat{n}_{j-1,2},\hat{\phi}_{j,j+1}^{(2)}=\hat{n}_{j,1},\hat{\phi}_{j,j-1}^{(2)}=-\hat{n}_{j-1,1}$.
The many-body interaction in the Hamiltonian leads to the nonintegrability and
many-body quantum chaos.

###  III.2 Level Statistics

In this work, we study systems in both Stark potentials
$W_{jm}/\bar{J}_{NN}=-j\gamma$, and random potentials $W_{jm}/\bar{J}_{NN}$
drawn from a uniform distribution in $[-W,W]$, as illustrated in the upper and
lower panel of Fig. S8(a), respectively. With the increasing of the gradient
$\gamma$ or the disorder strength $W$, a transition from ergodicity to
localization is expected in this typical nonintegrable model.

A common numerical way to characterize the ergodicity and MBL is the level
statistics based on random matrix theory, quantified by the ratio of adjacent
level spacings
$r^{(n)}=\min\left(\delta^{(n)},\delta^{(n+1)}\right)/\max\left(\delta^{(n)},\delta^{(n+1)}\right)$
with $\delta^{(n)}=E_{n}-E_{n-1}>0$. As the gradient $\gamma$ (or the disorder
strength $W$) increases, the level statistics changes from the Wigner-Dyson
distribution with the average $\langle r\rangle\simeq 0.531$ in the ergodic
phase, to the Poisson distribution with $\langle r\rangle\simeq 0.386$ in the
many-body localized phase (see Fig. S8(b)).

It is worth noting that, across all parameters depicted in Fig. LABEL:fig3 and
Fig. LABEL:fig4, the average ratio of adjacent level spacings $\langle
r\rangle$ consistently remains close to $0.531$, signifying thermalization for
at least the majority of spectrum.

Figure S8: Ergodicity to localization transition in Stark and random
potentials. (a) Schematic representation of the Stark potentials and random
potentials sampled from a uniform distribution to study the transition from
ergodicity to localization. All qubits are adjusted to corresponding
frequencies deviated from the center frequency
$\approx{4.534}~{}\mathrm{GHz}$. (b) The ratio of adjacent level spacings as
the function of gradient $\gamma$ or the disorder strength $W$, with level
statistics from the ergodic Wigner-Dyson distribution ($\langle r\rangle\simeq
0.531$) to the many-body localized Poisson distribution ($\langle
r\rangle\simeq 0.386$).

###  III.3 Participation entropy for Gaussian orthogonal ensemble and Gaussian
unitary ensemble

####  III.3.1 Derivation

Here, we give a direct derivation of the participation entropy predicted by
the Gaussian orthogonal ensemble (GOE) and the Gaussian unitary ensemble
(GUE), which corresponds to the eigenstates and the generic long-time-evolved
states, respectively, for a typical nonintegrable model [2, 3].

For the GOE, the $N$ components of an eigenvector is real number, and its norm
must be one. Therefore,

$P_{\text{GOE}}(\mathbf{c})=C_{N}\delta\left(1-\sum_{n=1}^{N}c_{n}^{2}\right).$
(S10)

For the function of random variants $f(\mathbf{x})$ subjects to distribution
$P(f)$, we have
$P(f)=\int{\mathbf{dx}\delta(P-f({\mathbf{x}}))p(\mathbf{x})}$. So the
marginal distribution of a single component by

$P_{\text{GOE}}(p)=\int dc_{1}\cdots
dc_{N}\delta\left(p-c_{i}^{2}\right)P_{\text{GOE}}(\mathbf{c}).$ (S11)

Defining

$P_{\text{GOE}}(p;t)=\int dc_{1}\cdots
dc_{N}\delta\left(p-c_{i}^{2}\right)C_{N}\delta\left(t-\sum_{n=1}^{N}c_{n}^{2}\right),$
(S12)

and taking the Laplace transform, we have:

$\displaystyle\int_{0}^{\infty}dte^{-st}P_{\text{GOE}}(p;t)$ (S13)
$\displaystyle=C_{N}\int
dc_{i}\delta\left(p-c_{i}^{2}\right)e^{-sc_{i}^{2}}\left(\int
dce^{-sc^{2}}\right)^{N-1}$
$\displaystyle\propto\frac{e^{-sp}}{\sqrt{p}s^{(N-1)/2}}.$

Then inverting the Laplace transform, we finally get

$\displaystyle P_{\text{GOE}}(p;t)\propto(t-p)^{(N-3)/2}\theta(t-p)/\sqrt{p}$
(S14) $\displaystyle\Rightarrow$ $\displaystyle
P_{\text{GOE}}(p)=P_{\text{GOE}}(p;1)=\frac{1}{\sqrt{\pi}}\frac{\Gamma(N/2)}{\Gamma((N-1)/2)}\frac{(1-p)^{(N-3)/2}}{\sqrt{p}}.$

Therefore, $P_{\text{GOE}}(p)\propto\exp{(-Np/2)}/\sqrt{p}$ for large $N\gg
1$, and

$\displaystyle S^{\text{PE}}_{\text{GOE}}=-\sum_{i}^{N}{{p_{i}}\log{p_{i}}}$
(S15) $\displaystyle=-N\int_{0}^{\infty}P_{\text{GOE}}(p)p\log(p)dp$
$\displaystyle\propto-N\int_{0}^{\infty}\exp{(-Np/2)}\sqrt{p}\log(p)dp.$

Finally, we get $S^{\text{PE}}_{\text{GOE}}=\log
N-2+\log(2)+\gamma_{e}\approx\log{0.482N}$, with $\gamma_{e}$ being Euler’s
constant.

Similarly, since the only constraint on $N$ components of an eigenvector for
the GUE is that its norm must be one, we have

$P_{\text{GUE}}(\mathbf{c})=C_{N}\delta\left(1-\sum_{n=1}^{N}\left|c_{n}\right|^{2}\right).$
(S16)

Defining

$P_{\text{GUE}}(p;t)=\int d^{2}c_{1}\cdots
d^{2}c_{N}\delta\left(p-\left|c_{i}\right|^{2}\right)C_{N}\delta\left(t-\sum_{n=1}^{N}\left|c_{n}\right|^{2}\right),$
(S17)

and taking the Laplace transform:

$\displaystyle\int_{0}^{\infty}dte^{-st}P_{\text{GUE}}(p;t)$ (S18)
$\displaystyle=C_{N}\int
d^{2}c_{i}\delta\left(p-\left|c_{i}\right|^{2}\right)e^{-s\left|c_{i}\right|^{2}}\left(\int
d^{2}ce^{-s|c|^{2}}\right)^{N-1}$
$\displaystyle\propto\frac{e^{-sp}}{s^{N-1}},$

then inverting the Laplace transform, we finally get

$\displaystyle P_{\text{GUE}}(p;t)\propto(t-p)^{N-2}\theta(t-p)$ (S19)
$\displaystyle\Rightarrow$ $\displaystyle
P_{\text{GUE}}(p)=P_{\text{GUE}}(p;1)=(N-1)(1-p)^{N-2}.$

Therefore, $P_{\text{GUE}}(p)=N\exp{(-Np)}$ for large $N\gg 1$, and

$\displaystyle S^{\text{PE}}_{\text{GUE}}=-\sum_{i}^{N}{{p_{i}}\log{p_{i}}}$
(S20) $\displaystyle=-N\int_{0}^{\infty}pNe^{-Np}\log pdp$ $\displaystyle=\log
N-1+\gamma_{e}$

Figure S9: Participation entropy for the $XX$-ladder model. (a)Participation
entropy of the eigenstate against its normalized energy. The dashed line
represent $S^{\text{PE}}_{\text{GOE}}\approx 8.73$. (b) Time-evolution of
participation entropy for the initial states
$\ket{\psi_{n_{\text{DW}}=6}}=\ket{\mathbb{10100101}}$ and
$\ket{\psi_{n_{\text{DW}}=2}}=\ket{\mathbb{11000011}}$. The dashed line
represent $S^{\text{PE}}_{\text{GUE}}\approx 9.04$.

####  III.3.2 Numerical verification

As stated before, the Hamiltonian [Eq. (LABEL:Ham1)] in the main text without
the on-site potential is a typical nonintegrable model, for which the PE of
highly excited eigenstates will be close to
$S^{\text{PE}}_{\text{GOE}}=\log\mathcal{N}-2+\log(2)+\gamma_{e}$. At the same
time, the time-evolved states for generic initial states correspond to
eigenvectors of a random matrix in the GUE, of which the PE will be
$S^{\text{PE}}_{\text{GUE}}=\log\mathcal{N}-1+\gamma_{e}$.

Here, we numerically calculate the PE of all eigenstates for the $XX$-ladder
model with the system size $8\times 2$. The PE is plotted against the
normalized energy
$\epsilon=\frac{E-E_{\text{min}}}{E_{\text{max}}-E_{\text{min}}}$ for each
eigenstate in Fig. S9(a). The PE of highly excited eigenstates
($\epsilon\approx 0.5$) is very close to $S^{\text{PE}}_{\text{GOE}}\approx
8.73$. We also plot the time-evolution of PE for initial states
$\ket{\psi_{0}}=\ket{\psi_{n_{\text{DW}}=6}}$ and
$\ket{\psi_{n_{\text{DW}}=2}}$ in Fig. S9(b). We can see that in a few hundred
nanoseconds, the PE closely approximates $S^{\text{PE}}_{\text{GUE}}\approx
9.04$ for both initial states.

###  III.4 Statistical Properties of Eigenstates in Stark Systems

Figure S10: Participation entropy of Eigenstates in Stark Systems. (a–c)
Participation entropy of eigenstates for a Stark potential with the gradient
$\gamma=0,\ 2$ and $4$. The colorbar indicates the logarithm of the eigenstate
occupation numbers
$\log_{10}{\left|{\langle\psi_{0}|E_{n}\rangle}\right|^{2}}$ for the initial
states $\ket{\psi_{n_{\text{DW}}=2}}$ and $\ket{\psi_{n_{\text{DW}}=6}}$.
(d–f) Histogram of the participation entropy at $\gamma=0,\ 2$ and $4$ in the
diagonal ensembles corresponding to the initial states
$\ket{\psi_{n_{\text{DW}}=2}}$ and $\ket{\psi_{n_{\text{DW}}=6}}$. Figure S11:
Entanglement entropy of eigenstates in Stark Systems. (a) Half-chain
entanglement cut (region S) of the ladder. (b–d) Entanglement entropy of
eigenstates for a Stark potential with the gradient $\gamma=0,\ 2$ and $4$.
The colorbar indicates the logarithm of the eigenstate occupation numbers
$\log_{10}{\left|{\langle\psi_{0}|E_{n}\rangle}\right|^{2}}$ for the initial
states $\ket{\psi_{n_{\text{DW}}=2}}$ and $\ket{\psi_{n_{\text{DW}}=6}}$.
Figure S12: Domain wall number of eigenstates in Stark Systems. (a–c) Domain
wall number of eigenstates for a Stark potential with the gradient $\gamma=0,\
2$ and $4$. The colorbar indicates the logarithm of the eigenstate occupation
numbers $\log_{10}{\left|{\langle\psi_{0}|E_{n}\rangle}\right|^{2}}$ for the
initial states $\ket{\psi_{n_{\text{DW}}=2}}$ and
$\ket{\psi_{n_{\text{DW}}=6}}$. (d) Participation entropy of eigenstates at
$\gamma=4$ with the colorbar representing the domain wall number
$n_{\text{dw}}$. (e) Participation entropy as a function of $n_{\text{DW}}$
for eigenstates with $\epsilon\in[0.4,0.6]$, corresponding to the region
between the dashed lines in (d).

####  III.4.1 Participation entropy

To quantify the localization in Hilbert space of spin configurations, we can
calculate the PE of all eigenstates with exact diagonalization, defined as

$S^{\text{PE}}(n)=-\sum_{i}^{\mathcal{N}}{{p_{i}(n)}\log{p_{i}(n)}},$ (S21)

where $p_{i}(n)={\left|{\langle i|E_{n}\rangle}\right|^{2}}$, with
$\\{\ket{i}\\}$ being spin configuration basis, and $\mathcal{N}$ is the
dimension of Hilbert space ($\mathcal{N}={16\choose 8}=12870$ for the $Q=L$
sector with $L=8$).

In Fig. S10(a–c), the PE of eigenstates is displayed for the Hamiltonian [Eq.
(LABEL:Ham1)] with a Stark potential at $\gamma=0,2$ and $4$. The colorbars
represent the eigenstate occupation numbers
$\left|{\langle\psi_{0}|E_{n}\rangle}\right|^{2}$ (i.e., the weights that
eigenstates $\ket{E_{n}}$ have in the diagonal ensemble corresponding to
$\ket{\psi_{0}}$). At $\gamma=0$, the PE $S^{\text{PE}}\approx
S^{\text{PE}}_{\text{GOE}}$ for eigenstates with $\epsilon\approx 0.5$ (same
as in Fig. S9(a)). As $\gamma$ increases, we can observe that the PE of
eigenstates decreases overall, with the range of fluctuations increasing.

We also displayed the histogram of the PE in the diagonal ensembles for
$\ket{\psi_{0}}=\ket{\psi_{n_{\text{DW}}=2}}$ and
$\ket{\psi_{n_{\text{DW}}=6}}$ in Fig. S10(d–f). For nonzero $\gamma$, the
histogram of the PE in the diagonal ensemble for
$\ket{\psi_{n_{\text{DW}}=2}}$ tends to distribute more at lower values of PE,
compared to the diagonal ensemble for $\ket{\psi_{n_{\text{DW}}=6}}$. From the
distribution of the PE in the diagonal ensemble, we can infer that in quench
dynamics, the PE at long times also differs for different initial states (as
illustrated in Fig. LABEL:fig4 of the main text).

####  III.4.2 Entanglement entropy

Another quantity that is often used to characterize MBL is the entanglement
entropy (EE) [4, 5], defined as

$S^{\text{EE}}(n)=-\text{Tr}[\hat{\rho}_{S}(n)\log{\hat{\rho}_{S}(n)}],$ (S22)

where $\hat{\rho}_{S}(n)=\text{Tr}_{\bar{S}}[\ket{E_{n}}\bra{E_{n}}]$ is the
reduced density matrix of the subsystem S (with $\bar{S}$ being the complement
of $S$) for the eigenstate $\ket{E_{n}}$. As shown in Fig. S11(a), the
subsystem S considered here is the left half of the ladder. In Fig. S11(b–d),
we display the EE of eigenstates for the Hamiltonian [Eq. (LABEL:Ham1)] with a
Stark potential at $\gamma=0,2$ and $4$. The colorbars represent the
eigenstate occupation numbers
$\left|{\langle\psi_{0}|E_{n}\rangle}\right|^{2}$ in the diagonal ensemble.
Similar to the PE, the EE of eigenstates decreases overall as $\gamma$
increases, and substantial fluctuations in the EE are observed across the
eigenstates for nonzero $\gamma$, indicating violations of ETH. At the same
time, compared to $\ket{\psi_{0}}=\ket{\psi_{n_{\text{DW}}=6}}$, eigenstates
with lower EE exhibit larger eigenstate occupation numbers in the diagonal
ensemble for $\ket{\psi_{n_{\text{DW}}=2}}$, which implies that the system
also exhibits initial-state dependent dynamics in terms of the EE.

####  III.4.3 Domain wall number

Besides the PE and the EE, violations of ETH can also be diagnosed from
observables like domain wall number, defined as

$n_{\text{DW}}=\sum_{j=1}^{L-1}(1-\langle\hat{\bar{\sigma}}_{j}^{z}\hat{\bar{\sigma}}_{j+1}^{z}\rangle)/2,$
(S23)

with
$\hat{\bar{\sigma}}_{j}^{z}=(\hat{\sigma}_{j,1}^{z}+\hat{\sigma}_{j,2}^{z})/2$.
In Fig. S12(a–c), we display the domain wall number of the eigenstates, with
colorbars being the logarithm of the eigenstate occupation numbers. At
$\gamma=0$, according to the ETH, the domain wall numbers $n_{\text{DW}}$ for
the majority of eigenstates are close to the ETH prediction of
$n_{\text{DW}}^{\text{ETH}}\approx 3.73$ for the ladder model. As $\gamma$
increases, the $n_{\text{DW}}$ of eigenstates fluctuates, with a number of
eigenstates significantly deviating from the ETH prediction. For nonzero
$\gamma$, eigenstate occupation numbers predominantly distribute on
eigenstates with relatively smaller $n_{\text{DW}}$ in the diagonal ensemble
for $\ket{\psi_{n_{\text{DW}}=2}}$, but they mainly distribute on eigenstates
with larger $n_{\text{DW}}$ for $\ket{\psi_{n_{\text{DW}}=6}}$.

In Fig. S12(d), we plot the PE of eigenstates at $\gamma=4$, where the
colorbar represents $n_{\text{dw}}$. It is evident that eigenstates with
higher PE tend to correspond to larger values of $n_{\text{DW}}$. To elucidate
this relationship more clearly, in Fig. S12(e), the PE is plotted against
$n_{\text{DW}}$ for eigenstates in the middle of the spectrum,
$\epsilon\in[0.4,0.6]$. This illustrates a discernible positive correlation
between the PE and $n_{\text{DW}}$.

Figure S13: Numerical Krylov subspace in Stark systems. (a), (b) The
cumulative distribution function of the time-evolved states for
$\ket{\psi_{0}}=\ket{\psi_{n_{\text{DW}}=2}}$ and
$\ket{\psi_{n_{\text{DW}}=6}}$ at $\gamma=0$, and $4$. The cumulative
distribution function for the state with equal weight on all the bases, and
for a random pure state corresponding to the GUE are also provided for
comparison. The shaded region in (a) corresponds to the bases which construct
the numerical Krylov subspace. (c) The dimension of the numerical Krylov
subspace for the initial states $\ket{\psi_{n_{\text{DW}}=2}}$ and
$\ket{\psi_{n_{\text{DW}}=L-2}}$ with $L=8$ and $12$ at $\gamma=2$. The inset
shows the time evolution of the dimension of $\mathcal{K}_{\delta}(t)$ from
$t=0$ to $40\mu$s. (d) The overlap between numerical Krylov subspaces
$\eta\equiv\frac{{\text{dim}(\mathcal{K}_{\delta}(\psi_{0})\cap\mathcal{K}_{\delta}(\overline{\psi_{0}}))}}{{\text{dim}(\mathcal{K}_{\delta}(\psi_{0})\cup\mathcal{K}_{\delta}(\overline{\psi_{0}}))}}$
for $\ket{\psi_{0}}$ and $\ket{\overline{\psi_{0}}}$, with
$\ket{\psi_{0}}=\ket{\psi_{n_{\text{DW}}=2}}$ and
$\ket{\psi_{n_{\text{DW}}=L-2}}$, respectively, and
$\ket{\overline{\psi_{0}}}=\prod_{j,m}\hat{\sigma}_{jm}^{x}\ket{\psi_{0}}$.
The dashed line denotes the typical value of the overlap between the subspaces
corresponding to two random pure states.

##  IV Numerical Krylov Subspace

To illustrate the construction of the numerical Krylov subspace, we firstly
introduce the cumulative distribution function (CDF) for a time-evolved state
$\ket{\psi(t)}$. The calculation of CDF involves the following steps:

  * i)

Calculate the probabilities of spin configuration basis
$p_{i}=\left|{\langle\psi(t)|i\rangle}\right|^{2}$ (i.e., the moduli squared
of the wave function coefficients expressed in spin configuration basis
$\\{\ket{i}\\}$);

  * ii)

Sort the probabilities in non-decreasing order, denoted as $p_{i^{\prime}}$,
with $\\{i^{\prime}\\}$ representing the sorted indices;

  * iii)

Compute the CDF as the sum of the sorted probabilities to the index $n$, given
by $\text{CDF}(n)=\sum_{i^{\prime}=1}^{n}p_{i^{\prime}}$.

The CDF defined here illustrates the bases with the probabilities arranged
from low to high, contributing to the state, with the curvature representing
the degree of unevenness in the distribution. One can expect the CDF would be
a straight line for the state with equal weight on all the bases
$\\{\ket{i}\\}$, while if the state is a single configuration, the CDF would
be $0$ for indices $n\leq\mathcal{N}-1$ of the sorted basis, and $1$ at index
$n=\mathcal{N}$.

In Fig. S13(a) and (b), we displays the CDF of the time-evolved states for
$\ket{\psi_{0}}=\ket{\psi_{n_{\text{DW}}=2}}$ and
$\ket{\psi_{n_{\text{DW}}=6}}$ at $\gamma=0$, and $4$. For both initial
states, the curve of CDF for the time-evolved states lies slightly beneath
that for the state with equal weight on all the bases, nearly coinciding with
the curve for a random pure state corresponding to the GUE at $\gamma=0$. At
$\gamma=4$, the CDF further curve downward for both initial states. Notably,
for $\ket{\psi_{0}}=\ket{\psi_{n_{\text{DW}}=2}}$, the curve is very close to
that of a single configuration state, suggesting that the time-evolved state
predominantly exhibits high probabilities on only a limited number of bases.

Based on the CDF, we can define the subset $\mathcal{K}_{\delta}(t)$ for an
initial state as the complement of the cumulative portion of the basis
corresponding to the cumulative probabilities $\delta$. This subset satisfies
$\sum_{i\in\mathcal{K}_{\delta}(t)}\left|{\langle\psi(t)|i\rangle}\right|^{2}=1-\delta$,
as depicted by the shaded region in Fig. S13(a). Therefore, the numerical
Krylov subspace $\mathcal{K}_{\delta}$ can be further defined as largest
subset satisfying $\mathop{\arg\max}\limits_{\mathcal{K}_{\delta}(t),\
t\leq\tau}\\{\text{dim}(\mathcal{K}_{\delta}(t))\\}$ within the evolution time
$\tau$ [6].

Here, we calculate $\mathcal{K}_{\delta}$ for $\delta=0.01$ and $\tau=40\mu$s,
at $\gamma=2$. The dimension of the numerical Krylov subspace is plotted in
Fig. S13(c) for the initial states $\ket{\psi_{n_{\text{DW}}=2}}$ and
$\ket{\psi_{n_{\text{DW}}=L-2}}$ with $L=8$ and $12$, with the inset showing
the dimension of $\mathcal{K}_{\delta}(t)$ from $t=0$ to $40\mu$s. We observe
that the dimension $\text{dim}\mathcal{K}_{\delta}(t)$ quickly saturates for
$\ket{\psi_{n_{\text{DW}}=L-2}}$, regardless of the system size; while for
$\ket{\psi_{n_{\text{DW}}=2}}$, the growth of
$\text{dim}(\mathcal{K}_{\delta}(t))$ slows considerably with increasing
system length $L$. Correspondingly, the dimension of the numerical Krylov
Subspace for $\ket{\psi_{n_{\text{DW}}=2}}$ decreases much more dramatically
with increasing $L$, compared to $\ket{\psi_{n_{\text{DW}}=6}}$. Thus, the
numerical Krylov subspaces $\mathcal{K}_{\delta}$ exhibit different scaling
behaviors for the two initial states, which leads to the initial-state
dependent dynamics of the PE as shown in Fig. LABEL:fig4 in the main text.

We have shown that initial states with different $n_{\text{DW}}$ belong to
different fragments, therefore exhibiting distinct dynamics. However, it is
worth noting that the reverse is not necessarily valid: initial states sharing
the same $n_{\text{DW}}$ may also each reside within nearly disconnected
fragments. Here, we consider a pair of initial states
$\\{\ket{\psi_{0}},\ket{\overline{\psi_{0}}\\}}$, which have the maximum
Hamming distance $L\times 2$ between each other [7], with
$\ket{\overline{\psi_{0}}}=\prod_{j,m}\hat{\sigma}_{jm}^{x}\ket{\psi_{0}}$.
For example, $\\{\ket{\psi_{n_{\text{DW}}=2}}=\ket{\mathbb{11000011}},\
\ket{\overline{\psi}_{n_{\text{DW}}=2}}=\ket{\mathbb{00111100}}\\}$ and
$\\{\ket{\psi_{n_{\text{DW}}=L-2}}=\ket{\mathbb{101000101}},\
\ket{\overline{\psi}_{n_{\text{DW}}=L-2}}=\ket{\mathbb{01011010}}\\}$ for
$L=8$. From symmetry considerations, one can expect that the initial states
$\\{\ket{\psi_{0}},\ket{\overline{\psi_{0}}\\}}$ exhibit the same dynamics in
terms of PE, EE, imbalance and other quantities. However, this does not
necessarily mean these initial states belong to the same fragment. To explore
this aspect, we introduce the following quantity:

$\eta\equiv\frac{{\text{dim}(\mathcal{K}_{\delta}(\psi_{0})\cap\mathcal{K}_{\delta}(\overline{\psi_{0}}))}}{{\text{dim}(\mathcal{K}_{\delta}(\psi_{0})\cup\mathcal{K}_{\delta}(\overline{\psi_{0}}))}},$
(S24)

which measures the overlap between numerical Krylov subspaces for
$\ket{\psi_{0}}$ and $\ket{\overline{\psi_{0}}}$.

The results are presented in Fig. S13(d). At $\gamma=2$, long-time $\eta$
between $\ket{\psi_{n_{\text{DW}}=L-2}}$ and
$\ket{\overline{\psi}_{n_{\text{DW}}=L-2}}$ is close to $\eta$ between two
random pure states (indicated by dashed lines in Fig. S13(d)) regardless of
the system size, which implies that they belong to the same Krylov subspace.
In contrast, $\eta$ between $\ket{\psi_{n_{\text{DW}}=2}}$ and
$\ket{\overline{\psi}_{n_{\text{DW}}=2}}$ fails to reach the value between two
random pure states even after a long time, and this discrepancy becomes more
pronounced with increasing system size. Specially, for $L=12$, $\eta$ remains
almost zero for a long time, and then slowly increase when $t\gtrsim 20\mu$s
($\bar{J}_{\text{NN}}t\gtrsim 10^{3}$), which provides evidence that
$\ket{\psi_{n_{\text{DW}}=2}}$ and $\ket{\overline{\psi}_{n_{\text{DW}}=2}}$
belong to different fragments, dynamically disconnected from each other. Note
that the increase of $\eta$ after $t\gtrsim 20\mu$s stems from the coupling of
different fragments by higher-order terms in the perturbative expansion of
Hamiltonian [Eq. (LABEL:Ham1)], which break the conservation of the emergent
dipole moment and eventually leads to the destruction of fragmentation.

These findings raise intriguing questions of whether the subspace
corresponding to $\ket{\psi_{n_{\text{DW}}=L-2}}$ and
$\ket{\overline{\psi}_{n_{\text{DW}}=L-2}}$ will further fragment at higher
$\gamma$, and whether there is a transition from weak to strong fragmentation
[8, 9] with increasing $\gamma$ in Stark systems. These questions requires a
detailed investigation involving scaling analysis with respect to evolution
time $t$ and system length $L$, which is left for future exploration.

##  V Time Evolution with Krylov Space Method

In dealing with large systems, the full diagonalization becomes unfeasible. In
our study involving up to 24 qubits, we employ the Krylov space method, which
offers a more efficient approach for computing the time evolution of the
system.

The unitary evolution of a isolated system governed by the Schrödinger
equation can be expressed as

$\ket{\psi(t+\tau)}=\exp{(-i\tau\hat{H}/\hbar)}\ket{\psi(t)},$ (S25)

where $\ket{\psi(t)}$ denotes the state at time $t$. The Krylov subspace of
dimension $m$ is spanned by the successive powers of $H$ applied to
$\ket{\psi(t)}$

$\mathcal{K}_{m}(\mathbf{H},\mathbf{v}_{0})=\text{Span}\\{\mathbf{v}_{0},\mathbf{H}\mathbf{v}_{0},\mathbf{H}^{2}\mathbf{v}_{0},\ldots,\mathbf{H}^{m-1}\mathbf{v}_{0}\\},$
(S26)

with $\mathbf{H}$ and $\mathbf{v}_{0}$ being the matrix form of the
Hamiltonian and the state, respectively. For Hermitian matrix $\mathbf{H}$,
via the Lanczos algorithm, an orthogonal transformation matrix
$\mathbf{K}_{m}=(\mathbf{v}_{0}\ \mathbf{v}_{1}\ \mathbf{v}_{2}\ \ldots\
\mathbf{v}_{m-1})$ can be constructed, which turns $\mathbf{H}$ into a
tridiagonal $m\times m$ matrix form,

$\mathbf{H}_{m}=\mathbf{K}_{m}^{\dagger}\mathbf{HK}_{m}=\left(\begin{array}[]{ccccc}\alpha_{0}&\beta_{1}&&&0\\\
\beta_{1}&\alpha_{1}&\beta_{2}&&\\\ &\beta_{2}&\alpha_{2}&\ddots&\\\
&&\ddots&\ddots&\beta_{m-1}\\\
0&&&\beta_{m-1}&\alpha_{m-1}\end{array}\right).$ (S27)

Here, $\alpha_{j}$, $\beta_{j}$ and $\mathbf{v_{j}}$ are determined by the
Lanczos algorithm, which is an iterative procedure

$\displaystyle\alpha_{i}=\mathbf{v}_{i}\cdot\mathbf{H}\mathbf{v}_{i},$ (S28)
$\displaystyle\beta_{i}\mathbf{v}_{i+1}=\mathbf{H}\mathbf{v}_{i}-\alpha_{i}\mathbf{v}_{i}-\beta_{i-1}\mathbf{v}_{i-1}.$

For a short time step $\tau$, the time evolution in S25 can be approximately
calculated as

$\ket{\psi(t+\tau)}\approx\mathbf{K}_{m}\exp{(-i\tau\mathbf{H}_{m}/\hbar)}\mathbf{K}_{m}^{\dagger}\mathbf{v}_{0}.$
(S29)

Given that the eigenvalues of $\mathbf{H}_{m}$ are the most crucial $m$
eigenvalues governing the dynamics at the current time $t$, the Krylov space
method can yields highly accurate results. In this work, we adaptively select
the time steps while ensuring convergence with $m=15$ Krylov vectors.

##  VI Post-Selection of Measurement Results

Figure S14: Effect of Post-Selection on Experimental Results. The experimental
data for the participation entropy (PE) with and without (grey data point)
post-selection (PS) within the $Q=L$ sector, in comparison with the numerical
results (solid lines), for specific initial states in both Stark ($\gamma=2$)
and random ($W=4$) potentials.

The Hamiltonian [Eq. (LABEL:Ham1)] in the main text conserves the total $U(1)$
charge $\hat{Q}\equiv\sum_{j,m}{\hat{\sigma}_{jm}^{+}\hat{\sigma}_{jm}^{-}}$
during the time evolution. Despite the experimental evolution time being up to
$800$ ns, significantly less than the qubit energy relaxation times $T_{1}\sim
33.2$ $\mu$s listed in Supplementary Tab. S1, the leakage out of the $Q=L$
sector could still lead to an overestimation of the participation entropy
(PE). To address this, we perform post-selection on the measured probabilities
within the $Q=L$ sector for the data presented in Fig. LABEL:fig4 in the main
text.

Comparisons between the experimental results with and without post-selection,
for $\psi_{n_{\text{DW}}=2}$ and $\psi_{n_{\text{DW}}=L-2}$ in a Stark and
random potential, are depicted in Fig. S14. Notably, the inclusion of post-
selection demonstrates a significantly improved agreement between the
experimental and numerical results.

## References

  * [1] Y.-H. Shi, Z.-H. Sun, Y.-Y. Wang, Z.-A. Wang, Y.-R. Zhang, W.-G. Ma, H.-T. Liu, K. Zhao, J.-C. Song, G.-H. Liang, Z.-Y. Mei, J.-C. Zhang, H. Li, C.-T. Chen, X. Song, J. Wang, G. Xue, H. Yu, K. Huang, Z. Xiang, K. Xu, D. Zheng, and H. Fan, Probing spin hydrodynamics on a superconducting quantum simulator,  arXiv:2310.06565 .
  * Torres-Herrera  _et al._ [2016] E. Torres-Herrera, J. Karp, M. Távora, and L. Santos, Realistic Many-Body Quantum Systems vs. Full Random Matrices: Static and Dynamical Properties, Entropy 18, 359 (2016).
  * Boixo  _et al._ [2018] S. Boixo, S. V. Isakov, V. N. Smelyanskiy, R. Babbush, N. Ding, Z. Jiang, M. J. Bremner, J. M. Martinis, and H. Neven, Characterizing quantum supremacy in near-term devices, Nat. Phys. 14, 595 (2018).
  * Nandkishore and Huse [2015] R. Nandkishore and D. A. Huse, Many-Body Localization and Thermalization in Quantum Statistical Mechanics, Annu. Rev. Condens. Matter Phys. 6, 15 (2015).
  * Abanin  _et al._ [2019] D. A. Abanin, E. Altman, I. Bloch, and M. Serbyn, Colloquium : Many-body localization, thermalization, and entanglement, Rev. Mod. Phys. 91, 021001 (2019).
  * Scherg  _et al._ [2021] S. Scherg, T. Kohlert, P. Sala, F. Pollmann, B. Hebbe Madhusudhana, I. Bloch, and M. Aidelsburger, Observing non-ergodicity due to kinetic constraints in tilted Fermi-Hubbard chains, Nat. Commun. 12, 4490 (2021).
  * Yao  _et al._ [2023] Y. Yao, L. Xiang, Z. Guo, Z. Bao, Y.-f. Yang, Z. Song, H. Shi, X. Zhu, F. Jin, J. Chen, S. Xu, Z. Zhu, F. Shen, N. Wang, C. Zhang, Y. Wu, Y. Zou, P. Zhang, H. Li, Z. Wang, C. Song, C. Cheng, R. Mondaini, H. Wang, J. Q. You, S.-y. Zhu, L. Ying, and Q. Guo, Observation of many-body Fock space dynamics in two dimensions, Nat. Phys. 19, 1459 (2023).
  * Taylor  _et al._ [2020] S. R. Taylor, M. Schulz, F. Pollmann, and R. Moessner, Experimental probes of Stark many-body localization, Phys. Rev. B 102, 054206 (2020).
  * Sala  _et al._ [2020] P. Sala, T. Rakovszky, R. Verresen, M. Knap, and F. Pollmann, Ergodicity Breaking Arising from Hilbert Space Fragmentation in Dipole-Conserving Hamiltonians, Phys. Rev. X 10, 011047 (2020).


# Entanglement Distribution and Quantum Teleportation in Higher Dimension over
the Superposition of Causal Orders of Quantum Channels ††thanks: This material
is based upon work supported by Science Foundation Ireland (SFI) and is co-
funded under the European Regional Development Fund under Grant Numbers
13/RC/2077 and 13/RC/2077-P2.
©2023 IEEE. Personal use of this material is permitted. Permission from IEEE
must be obtained for all other uses, in any current or future media, including
reprinting/republishing this material for advertising or promotional purposes,
creating new collective works, for resale or redistribution to servers or
lists, or reuse of any copyrighted component of this work in other works.

I. DEY and N. MARCHETTI

###### Abstract

Multiple photonic degrees of freedom can be explored to generate high-
dimensional quantum states; commonly referred to as ‘qudits’. Qudits offer
several advantages for quantum communications, including higher information
capacity, noise resilience and data throughput, and lower information loss
over different propagation mediums (free space, optical fibre, underwater) as
compared to conventional qubits based communication system. However, qudits
have been little exploited in literature, owing to their difficulty in
transmission and detection. In this paper, for the first time, we develop and
formulate the theoretical framework for transmission of classical information
through entanglement distribution of qudits over two quantum channels in
superposition of alternative causal order. For the first time we i) engineer
quantum switch operation for 2-qudit systems and ii) formulate theoretical
system model for entanglement distribution of qudits via quantum switch.
Results show that entanglement distribution of a qudit provides a considerable
gain in fidelity even with increase in noise.

_K_ eywords $d$-dimensional qudits, entanglement distribution, generalized
quantum gates, quantum communication, quantum switch

##  1 Introduction

Communication networks are gradually emerging as critical infrastructure for
all walks of life. Presently, the choicest mode of communication for any
network is combination of radio frequency (RF) wireless links for fronthaul
and optical fibre links for backhaul connection. However, accommodating
increasing number of users and applications in the RF band will be impossible
in terms of carrier frequency and bandwidth. On the other hand, optical fibre
links impose high latency over long-distance connections. Quantum
communication networks, in future, will be the most promising alternative that
can solve problems with wireless and fibre links [1].

Quantum communication between any two points is viable through a property
called ‘shared entanglement’. Shared entanglement can be accomplished by
distributing an Einstein-Podolsky-Rosen (EPR) pair between two classically or
quantumly operated nodes communicating through a quantum channel [2].
Entangled qubits are distributed between two nodes with high fidelity. Since
quantum states cannot be copied or amplified, success of quantum communication
over long distances strongly depends on the use of one or more quantum
repeaters that forward information by entanglement distribution [3].

Quantum particles experiencing two quantum processes in superposition of
alternative orders (i.e. which channel is traversed before the other is not
imperative) are capable of conveying more information than when travelling
through each process of arbitrary order separately and consecutively [4]. For
example, in normal condition, no quantum information can be transmitted over
either of the pair of dephasing channels [5]. However, a qubit travelling
through the superposition of two dephasing channels of alternative orders can
reach the receiver error-free with 25% probability. The success of an error-
free transmission can only be validated by the outcome of a measurement on the
order of the qubit [6], where the order depends on the state of the quantum
degree of freedom. Similarly, for partially dephasing channels, quantum
capacity of two superposed channels is larger than that achievable over
individual channels.

The travelling quantum particles can arrive intact at the receiver with a
finite probability of success. This advantage can be capitalized on through
the implementation of a quantum ‘switch’ in a quantum teleportation chain for
entanglement distribution, a phenomenon by which a pair of maximally entangled
qubits must be generated and shared between the transmitter and the receiver
[7]. In a quantum switch, the transmitted quantum particle experiences two
quantum channels $\mathcal{A}$ and $\mathcal{B}$ in a superposition of
different orders, $\mathcal{A}\to\mathcal{B}$ and $\mathcal{B}\to\mathcal{A}$.
The order depends on the state of a quantum degree of freedom, which can be
two-dimensional for qubits.

This phenomenon can be leveraged to create a perfect noiseless communication
channel through quantum superposition of two noisy channels of alternative
orders. As a result, fidelity of a teleported qubit via an ordinary channel
can be enhanced by using a maximally entangled state as a quantum resource in
a superposition of alternative orders [8]. Moreover, quantum entanglement is
easily affected by noise resulting in degeneration of the transmitted
information, an effect that can be significantly ameliorated by the
application of quantum switch enabling superposition of noisy quantum channels
of different causal orders [9].

Qubits can only exploit two photonic degrees of freedom and are therefore
limited in channel capacity, data throughput and noise resilience, thereby
lowering the possibility of breaking the Shannon’s limit on achievable
capacity. This has spurred a recent interest in teleporting qudits as units
for information exchange with arbitrary dimensions exploiting $d$-degrees of
freedom [10]. Teleportation of photonic qudits increases the quantum
information sent per carrier photon. Teleportation of qudits can be achieved
either by preparing $d$ additional photons in a highly entangled state [11] or
by transcribing qudit encoded on a single photon to $d$ qubits carried by
light modes which propagate along different optical paths [12]. Moreover,
qudits can carry information over a plethora of media ranging from free space
optical links to multi-core and multi-mode optical fibres, and peer-to-peer
underwater acoustic communications [13].

Complementary to the preliminary efforts of the quantum physics [14] and
quantum communication engineering community [15] in exploiting quantum switch
for entanglement distribution of 2-qubits system, we investigate the
possibility of employing quantum switch for entanglement distribution of
qudits teleportation system. To achieve this, we formulate the theoretical
system model for quantum switch for entanglement distribution of 2-qudit
systems. We also derive closed-form expressions for connecting the teleported
qudit received by Bob to the degradation suffered by the entangled pair during
the distribution process.

The rest of the paper is organized as below. Section 2 introduces the concept
and mathematical expressions for the functionality of a quantum switch for
qudits. Section 3 develops the basic framework for quantum teleportation of
qudits. Section 4 follows up with the introduction to different circuit
configurations possible for realizing a qudit switch and the theoretical
framework for entanglement distribution of qudits. Finally, we develop the
mathematical formulation for teleportation of qudits through the
implementation of entanglement distribution of a 2-qudit system by means of a
qudit-based quantum switch in Section 5. The concluding remarks are provided
in Section 6.

##  2 Quantum Switch - Concept for Qudits

Let us consider an arbitrary qudit $|\psi\rangle=\\{|k\rangle\\}_{k=0}^{d-1}$
of dimension $d$ and unknown state travelling through two noisy quantum
channels $\mathcal{A}$ and $\mathcal{B}$, and let us assume that $\mathcal{A}$
is the _cyclic shift_ channel and $\mathcal{B}$ is the _cyclic clock_ channel.
The cyclic shift channel $\mathcal{A}$ performs a cyclic permutation of the
basis states and maps the set, $\\{|0\rangle,|1\rangle,\dotso,|d-1\rangle\\}$
to $\\{|1\rangle,|2\rangle,\dotso,|0\rangle\\}$ and vice-versa with a
probability $p$, leaving the qudit unaltered with probability $1-p$, such
that,

$\displaystyle\mathcal{A}(|\psi\rangle)=(1-p)|\psi\rangle+pX_{d}|\psi\rangle$
(1)

where $X_{d}$ denotes the generalized Pauli-X gate given in Table 1. Here (1)
is the qudit generalization of the bit-flip channel that introduces a flip of
the state of a qubit from $|0\rangle$ to $|1\rangle$ with a certain
probability.

The cyclic clock channel $\mathcal{B}$ introduces relative phase shift of
$\omega=e^{2\pi i/d}$ between complex set of amplitudes with probability $q$,
leaving the qudit unaltered with probability $1-q$, such that,

$\displaystyle\mathcal{B}(|\psi\rangle)=(1-q)|\psi\rangle+qZ_{d}|\psi\rangle$
(2)

where $Z_{d}$ denotes the generalized Pauli-Z gate given in Table 1. Here, (2)
is the qudit generalization of the phase-flip channel that introduces a
relative phase-shift between complex amplitudes. The complex set of amplitudes
are given by $\alpha_{k}$ such that,
$|\psi\rangle=\sum_{k=0}^{d-1}\alpha_{k}|k\rangle$, where
$\alpha_{k}\in\mathbb{C}^{d},\sum_{k=0}^{d-1}|\alpha_{k}|^{2}=1$,
$\mathbb{C}^{d}$ is the computational basis of $d$-dimensional Hilbert space,
$\mathcal{H}^{d}$. A $n$-qudit state, in the tensor product Hilbert space is
given by, $\mathcal{H}^{\otimes n}=(\mathbb{C}^{d})^{\otimes n}$.

The standard basis of $\mathcal{H}^{\otimes n}$ is the orthonormal basis given
by the $d^{n}$ classical $n$-qudits, such that

$|k_{1}\dotso k_{n}\rangle=|k_{1}\rangle\otimes\dotso|k_{n}\rangle$

where $0\leq k_{i}\leq d-1$ (for $i=1,\dotso,n$). Since a qudit itself is an
arbitrary superposition of orthonormal basis states $|k\rangle$ for
$d$-dimensional quantum system, quantum switch for qudits will operate in a
slightly different way with respect to quantum switch for qubits. For example,
quantum circuit architectures completely described by instances of the CNOT
gate cannot implement a transposition of a pair of qudits for dimension $d>2$.
Therefore, in order to design quantum switch for qudits, we need to consider
the problem of generalizing the design beyond the qubit settings.

Next we introduce the concept of quantum switch for qudits where the
traversing particles experience superposition of different orders, depending
on the $d$ degrees of freedom. Whenever the control qudit is initialized to
one of the basis states $|\xi\rangle$, such that $|\xi\rangle$ is a prime
$|\xi\rangle\geq 2$, the quantum switch enables the message $m$ to experience
the alternative classical trajectory $\mathcal{B}\to\mathcal{A}$ representing
channel $\mathcal{B}$ being traversed before channel $\mathcal{A}$. Whenever
the control qudit is initialized to a non-prime basis state, the quantum
switch enables the message $m$ to experience the classical trajectory
$\mathcal{A}\to\mathcal{B}$ representing channel $\mathcal{A}$ being traversed
before channel $\mathcal{B}$. Whenever the control qudit is initialized to a
superposition of the basis states, such as,

$\displaystyle|\varphi\rangle=\frac{1}{\sqrt{d}}\sum_{j=0}^{d-1}e^{2\pi
i\varphi j/d}|j\rangle$ (3)

where $\varphi$ is an integer from the set $\\{0,\dotso,d-1\\}$ and
$|\varphi\rangle$ is an eigenstate of the cyclic shift operator
$X_{d}|\cdot\rangle$. Now, since a quantum superposition of two alternative
orders of noisy channels, specifically $\mathcal{A}\to\mathcal{B}$ and
$\mathcal{B}\to\mathcal{A}$, can behave as a perfect quantum communication
channel [4], we will use it to conceptualize quantum teleportation of qudits
through entanglement generation and distribution via the application of
quantum switch for qudits.

Table 1: Higher dimensional Quantum Gates Gate |  Symbol |  align |  Description
---|---|---|---
Generalized Pauli-$X$ |  $X_{d}$ |  $X_{d}^{q}|k\rangle=|k\oplus q\rangle$ |  $q$ is the shift change in standard basis
Generalized Pauli-$Z$ |  $Z_{d}$ |  $Z_{d}^{l}|k\rangle=\omega^{lk}|k\rangle$ |  $l$ is the phase change in standard basis
Generalized Hadamard |  $H_{d}$ |  $H_{d}|k\rangle=\frac{1}{\sqrt{d}}\sum_{m=0}^{d-1}x_{m}$, $x_{m}=(-1)^{k\cdot m}|m\rangle$ |  $k\cdot m=\oplus_{\rho=0}^{n-1}k_{\rho}m_{\rho}$ is the bitwise inner product
CNOT Right-Shift |  $R_{c_{d}}$ |  ${R_{c}}_{d}(|k\rangle\otimes|x\rangle)=|k\rangle\otimes|x\oplus k\rangle$ |  $\oplus$ denotes modulo-$d$ addition and $|k\rangle\otimes|x\rangle$ is the standard basis state of $\mathcal{H}_{\mathcal{A}}\otimes\mathcal{H}_{\mathcal{B}}$
CNOT
Left-Shift |  $L_{c_{d}}$ |  ${L_{c}}_{d}(|k\rangle\otimes|x\rangle)=|k\rangle\otimes|x\ominus k\rangle$ |  $\ominus$ denotes modulo-$d$ subtraction
Generalized XOR (GXOR) |  $\text{GXOR}_{d}$ |  GXOR${}_{d}|k,l\rangle=|k,k\ominus l\rangle$ |  GXOR is the generalized quantum XOR-gate

##  3 Quantum Teleportation using Qudits

Our model of qudit teleportation protocol consists of two processes: Alice and
Bob. The sender Alice possesses a qudit of unknown state

$\displaystyle|\psi\rangle_{\alpha}=\sum_{k=0}^{d-1}\alpha_{k}|k\rangle_{d}$
(4)

that is to be teleported to Bob. This qudit is teleported using maximally
entangled Einstein, Podolsky, and Rosen (EPR) pair of $|\Phi^{+}\rangle$
generated by entangling qudits $|\Phi\rangle_{1}$ and $|\Phi\rangle_{2}$ such
that,
$|\Phi\rangle_{12}=\frac{1}{\sqrt{d}}\sum_{k=0}^{d-1}|k\rangle_{1}|k\rangle_{2}$,
through local quantum operations and classical communications [16], as shown
in Fig. 1.

Figure 1: Quantum Teleportation Circuit for Qudits, where $|\Psi\rangle$ is
the qudit of unknown state that is to be teleported from the sender Alice to
the receiver Bob, $|\Phi^{+}\rangle$ is the maximally entangled EPR pair to be
generated by entangling qudits $|\Phi\rangle_{1}$ and $|\Phi\rangle_{2}$ where
$|\Phi\rangle_{1}$ is stored at the sender and $|\Phi\rangle_{2}$ is stored at
the receiver. The double line represents the transmission of a classical
$d$-ary digit from the sender to the receiver.

Firstly, CNOT Left shift (defined in Table 1) operation is applied to qudits
$|\psi\rangle_{\alpha}$ and $|\Phi^{+}\rangle$ followed by Hadamard operation
(also defined in Table 1) applied to $|\psi\rangle_{\alpha}$. Finally, the
qudits $|\psi\rangle_{\alpha}$ and $|\Phi^{+}\rangle$ are measured resulting
in classical values ranging between 0 and $d-1$. Using the classical values
Bob can perform necessary unitary operations to recover the original state
$|\psi\rangle_{\alpha}$. To begin with, Alice will make a joint von Neumann
measurement of the qudit $|\psi\rangle_{\alpha}$ and the qudit
$|\Phi\rangle_{1}$. To make this measurement, Alice will use the Bell basis of
the qudits,

$\displaystyle|\Psi_{yz}\rangle=\frac{1}{\sqrt{d}}\sum_{k=0}^{d-1}e^{2\pi
iyk/d}|k\oplus z\rangle|k\rangle$ (5)

where $k\oplus z$ means sum of $k$ and $z$ modulo $d$, $y,z=0,1,\dotso,d-1$,
$y$ and $z$ can take integer values between 0 and $d-1$ and form the Bell
basis for the qudits such that $|\Psi_{yz}\rangle$ form a set of basis vectors
for a two-qudit systems. $|\Psi_{yz}\rangle$ can be inverted to obtain,

$\displaystyle|st\rangle=\frac{1}{\sqrt{d}}\sum_{y,w=0}^{d-1}e^{-2\pi
iyk/d}\delta_{t,s\oplus w}|\Phi^{yw}\rangle.$ (6)

The combined state of the system $|\psi\rangle_{\alpha}|\Phi\rangle_{12}$ in
terms of the above Bell basis vectors for the system
$|\psi\rangle_{\alpha}|\Phi\rangle_{1}$ is as follows,

$\displaystyle|\psi\rangle_{\alpha}|\Phi\rangle_{12}=\frac{1}{d}\sum_{y,v=0}^{d-1}|\Psi_{yv}\rangle_{\alpha
1}U^{\dagger}_{yv}|\chi\rangle_{2}$ (7)

where these unitary operators $U_{vw}$ are given by,

$\displaystyle U_{vw}=\sum_{k=0}^{d-1}e^{2\pi ivk/d}|k\rangle\langle k\oplus
w|$ (8)

and obey the following orthogonality conditions, $U^{\dagger}$ is the adjoint
of $U$ and $U^{\dagger}$ is the conjugate transpose of $U$ in basis states,
$\text{Tr}(U^{\dagger}_{vw}U_{yz})=d~{}\delta_{vy}\delta_{wz}$, where Tr
denotes trace of a square matrix and $\delta_{vy}$ and $\delta_{wz}$ are
Kronecker delta functions. We introduced the indices $v$ and $w$, such that
they can take any integer between 0 and $d-1$. We introduced these indices to
differentiate between the unitary operators used by Alice ($U_{yv}$) and Bob
($U_{vw}$). Alice uses the unitary operator to control the qudit state for
teleportation and Bob uses a different set of operators to convert its qudit
state back to the input state.

Figure 2: Different circuit realizations for distributing and teleporting
qudit entanglement where the $U$-gates control the trajectory path of the
entanglement carrier from $\mathcal{A}\to\mathcal{B}$ and from
$\mathcal{B}\to\mathcal{A}$ depending on the control qudit $\rho_{c}$. The
$U^{\dagger}$-gates recombine the paths $\mathcal{A}\to\mathcal{B}$
trajectories of the entanglement carriers to generate the final density
$\rho_{e}^{QS}$ of the EPR pair distributed between the sender and the
receiver through the generalized SWAP gate for qudits.

In summary, Alice makes the von Neumann measurement in the Bell basis to
obtain one of the possible $d^{2}$ results. She conveys the results of her
measurement to Bob by sending $2\log_{2}d$ classical bits of information.
After receiving this information Bob uses appropriate unitary operator
$U_{vw}$ on his qudit to convert its state to that of the input state, thereby
completing the standard teleportation of a qudit $|\psi\rangle_{\alpha}$ of
arbitrary state over quantum channels. Though quantum teleportation with the
aid of quantum entanglement looks promising, it is a very fragile resource and
is easily degraded by noise [22], resulting in loss of teleported information.
However, the situation can be ameliorated by exploiting quantum superposition
of different causal orders realized through a quantum switch.

##  4 Entanglement Distribution of Qudits

We aim to engineer entanglement distribution of qudits and thereby develop the
theoretical framework for entanglement distribution of qudits via quantum
switch from a communications engineer’s point of view, and different circuit
realizations to achieve this are given in Fig. 2.

The operator $U$ routes the entanglement carrier $|\Phi^{+}\rangle_{2}$
through either channel $\mathcal{A}$ or $\mathcal{B}$ depending on the state
of $|\varphi_{c}\rangle$. The combination of $X_{d}$ and GXOR gates performs
the function of routing the entanglement carrier through the other portion of
the circuit through the trajectories $\mathcal{A}\to\mathcal{B}$ and
$\mathcal{B}\to\mathcal{A}$ respectively (refer to Fig. 2-Configuration 1).
The $U^{\dagger}$-gate is responsible for recombining the
$\mathcal{A}\to\mathcal{B}$ trajectories and generating the final set of EPR
pairs to be distributed between the transmitter and the receiver. It is worth-
mentioning here that the GXOR is the generalized quantum XOR-gate (refer to
Table 1) working according to

$\displaystyle\text{GXOR}_{d}|k,l\rangle=|k,k\ominus l\rangle.$ (9)

A quantum switch for qudits can be implemented through the generalized SWAP
gate operation. The generalized SWAP gate is not a single basic quantum gate
for qudits; it has to be engineered using combinations of basic qudit gates.
The main circuit idea for qudit SWAP, in this paper, is taken from [17] where
three GXOR and two $X_{d}$ gates are implemented for interchanging control and
target qudits between two traversing quantum channels.

Generalized SWAP gate can also be realized through several other combinations
of basic higher dimensional quantum gates. A few common combinations involve
CNOT Right-Shift ($R_{c_{d}}$), CNOT Left-Shift ($L_{c_{d}}$) and generalized
Pauli-X ($X_{d}$) gates (refer to Table 1). Using these combinations, we
extend our design for quantum switch in Fig. 2-Configuration 1 to different
implementational ways depicted in Figs. 2-Configuration 2, 2-Configuration 3
and 2-Configuration 4. The SWAP Gate for qudits can be implemented in many
different ways depending on the requirement. The design can involve from at
least 4 gates to 6 gates or more. The qudit GXOR gate can be combined with
qudit Pauli $X$-gates in a cyclic fashion to generate the SWAP operation [18,
19, 20]. On the other hand, the qudit CNOT gates can also be combined with
Pauli-$X$ gates in a cyclic fashion without the involvement of any GXOR gates.
In that case, they are extensions of SWAP operation for qubits to qudits. This
is because, a combination of a single CNOT and a single Pauli $X$-gate is used
to implement SWAP operation for qubits only. It is worth-mentioning here, that
three CNOT gates can be concatenated to implement two-qubit SWAP operation,
where the control and the target qubits are exchanged in the second gate [21].
However, the choice of the design of the SWAP gate will have an impact on the
overall operation of the quantum switch. By comparing the quantum fidelity and
channel capacity achievable through qudit teleportation over each of the
proposed design for quantum SWAP gates (refer to Fig. 2), we can choose the
best possible combination for the design of a quantum switch. We leave this
interesting design detail for our future work.

To analytically model a quantum switch, we extend the concept of switch for
one-qubit system to that for a one-qudit system
$\mathcal{P}(\mathcal{A},\mathcal{B},\rho_{c})(\rho)$ as the output function
of the two channels $\mathcal{A}$ and $\mathcal{B}$ along with the control
state $\rho_{c}=|\varphi_{c}\rangle\langle\varphi_{c}|$ of the control qudit,

$\displaystyle\mathcal{P}(\mathcal{A},\mathcal{B},\rho_{c})(\rho)=\sum_{s,t}W_{st}(\rho\otimes\rho_{c})W^{\dagger}_{st}$
(10)

where $\\{W_{st}\\}$ denotes the set of Kraus operators associated with the
superposed channel trajectories given by,

$\displaystyle{}{W_{st}}={}{\mathcal{A}_{s}\mathcal{B}_{t}\otimes|k\rangle\langle
k|}$ (11)

with $\\{\mathcal{A}_{s}\\}$ and $\\{\mathcal{B}_{t}\\}$ denoting the Kraus
operators associated with the channels $\mathcal{A}$ and $\mathcal{B}$
respectively. If the control quantum state is a qubit then,

$\displaystyle{}{W_{st}}$
$\displaystyle={}{\mathcal{A}_{s}\mathcal{B}_{t}\otimes|0\rangle\langle
0|_{c}+\mathcal{B}_{t}\mathcal{A}_{s}\otimes|1\rangle\langle 1|_{c}}.$ (12)

If the control quantum state is a qudit then
$\rho_{c}=|\psi_{c}\rangle\langle\psi_{c}|$, and,

$\displaystyle{}{W_{st}}$
$\displaystyle={}{\sum_{k,k^{\prime}=0}^{d-1}\big{[}\mathcal{A}_{s}\mathcal{B}_{t}\otimes|k\rangle\langle
k|+\mathcal{B}_{t}\mathcal{A}_{s}\otimes|k^{\prime}\rangle\langle
k^{\prime}|\big{]}}.$ (13)

The Kraus operators $W_{st}$ and $W_{st}^{+}$ satisfy the completeness
property,

$\displaystyle{}{\sum_{\\{s,t\\}|_{s,t=0}^{d-1}}W_{st}W_{st}^{\dagger}}$
$\displaystyle={}{I_{T}\otimes I_{c}}.$ (14)

where $I_{T}$ and $I_{c}$ are the identity operators in the target and control
systems respectively. This check of completeness suggests how the indices $s$
and $t$ allow for systematic reordering of the sums by isolating or
regrouping. It is worth-mentioning here that the entanglement-pair member
$|\Phi^{+}\rangle_{1}$ is already at Alice’s side and therefore does not need
to go through the communication channel, while the second qudit
$|\Phi^{+}\rangle_{2}$ is distributed to Bob. Therefore, the first qudit
basically travels virtually through an ideal channel represented by identity
matrix $I_{k\times k}$.

Figure 3: Quantum Teleportation Process of qudits in terms of density matrices
where $\rho_{\Psi}$, $\tilde{\rho}_{e}$ and $\rho_{\tau}$ denote the density
matrices of the actual qudit sent by the sender, EPR pair of qudits shared
between the sender and the receiver, and the actual qudit received by the
receiver respectively, and the double line represents the classical
trajectories.

Next we extend quantum switch for one qudit system to the case of a two-qudit
system, represented by the density matrix $\rho_{e}$ that is a $d^{2}\times
d^{2}$ matrix taking into account the cyclic shift and the cyclic clock
channels. In order to formulate $\rho_{e}$, we start by defining the Hermitian
generators in the orthonormal basis $|j\rangle$ which takes on the values
$j=1,\dots,d$, given by,

$\displaystyle\Gamma_{j}$
$\displaystyle\equiv\frac{1}{\sqrt{j(j-1)}}\Big{(}\sum_{u=1}^{j-1}|u\rangle\langle
u|-(j-1)|j\rangle\langle j|\Big{)}$ $\displaystyle\Gamma_{ju}^{(+)}$
$\displaystyle\equiv\frac{1}{\sqrt{2}}\Big{(}|j\rangle\langle
u|+|u\rangle\langle j|\Big{)}$ $\displaystyle\Gamma_{ju}^{(-)}$
$\displaystyle\equiv\frac{-i}{\sqrt{2}}\Big{(}|j\rangle\langle
u|-|u\rangle\langle j|\Big{)}$ (15)

for $1\leq j<u\leq d$. It is noteworthy that the generators are unitary
matrices. Using orthonormality relation, we can invert the generators to give,

$\displaystyle|j\rangle\langle j|$
$\displaystyle=\frac{I}{d}+\frac{1}{\sqrt{j(j-1)}}\Big{(}-(j-1)\Gamma_{j}+\sum_{u=j+1}^{d}\Gamma_{u}\Big{)}$
$\displaystyle|j\rangle\langle u|$
$\displaystyle=\frac{1}{\sqrt{2}}\Big{(}\Gamma_{ju}^{(+)}+i\Gamma_{ju}^{(-)}\Big{)}$
$\displaystyle|u\rangle\langle j|$
$\displaystyle=\frac{1}{\sqrt{2}}\Big{(}\Gamma_{ju}^{(+)}-i\Gamma_{ju}^{(-)}\Big{)}$
(16)

We consider a class of two qudit states $|\Phi\rangle_{1}$ and
$|\Phi\rangle_{2}$ with a maximally entangled state, which can be chosen to be
$|\Phi^{+}\rangle=\frac{1}{\sqrt{d}}\sum_{j=1}^{d}|j\rangle\otimes|j\rangle$.
By assuming without any loss of generality $\rho_{e}$ being the $d^{2}\times
d^{2}$ density matrix associated with the EPR pair $|\Phi^{+}\rangle$, we can
calculate it from the inverted generators to give,

$\displaystyle\rho_{e}=|\Phi^{+}\rangle\langle\Phi^{+}|=\frac{1}{d^{2}}(I\otimes
I)$
$\displaystyle+\frac{1}{d}\sum_{j}\Big{[}\Gamma_{j}\otimes\Gamma_{j}+\Gamma_{ju}^{(+)}\otimes\Gamma_{ju}^{(+)}+\Gamma_{ju}^{(-)}\otimes\Gamma_{ju}^{(-)}\Big{]}.$
(17)

$\displaystyle{}{\mathcal{P}(\mathcal{A},\mathcal{B},\rho_{c})(\rho)=\bigg{\\{}(1-p)(1-q)\rho_{e}+(1-p)q\bigg{[}\sum_{s}\sum_{k=0}^{d-1}|\lambda_{sk}\rangle\langle\lambda_{sk}|\bigg{]}\rho_{e}\bigg{[}\sum_{s}\sum_{k=0}^{d-1}|\lambda_{sk}\rangle\langle\lambda_{sk}|\bigg{]}^{\dagger}+p(1-q)\bigg{[}\sum_{s,t}|st\rangle\langle\lambda_{st}|\bigg{]}}$
$\displaystyle{}{\times\rho_{e}\bigg{[}\sum_{s,t}|st\rangle\langle\lambda_{st}|\bigg{]}^{\dagger}\bigg{\\}}\otimes\bigg{[}\frac{1}{\sqrt{d}}\sum_{j=0}^{d-1}\omega^{\varphi_{c}j}|j\rangle\langle
j|\bigg{]}+pq\bigg{[}\sum_{s,t}\iota_{t}|st\rangle\langle\lambda_{st}|\bigg{]}\rho_{e}\bigg{[}\sum_{s,t}\iota_{t}|st\rangle\langle\lambda_{st}|\bigg{]}^{\dagger}\otimes\bigg{[}\frac{1}{\sqrt{d}}\sum_{j=0}^{d-1}\omega^{\varphi_{c}j}|j\rangle\langle
j|\bigg{]}.}$ (18)

Following the aligns for $\rho_{e}$, the global quantum state at the output of
the switch can be calculated by exploiting the tensor product properties after
some algebraic manipulation to obtain (4), where the control qudit
$\varphi_{c}$ is transformed into a mixed state of basis states,
$\lambda_{sk}=s_{\oplus(d-1)}k$, $\lambda_{st}=st_{\oplus(d-1)}$,
$\omega^{\varphi_{c}j}=e^{i2\pi\varphi_{c}j/d}$ and
$\iota_{t}=(-1)^{t\oplus(d-1)}$, where $\oplus(d-1)$ in subscript represents
the direct sum of $(d-1)$ copies of the fundamental variable. The global state
$\mathcal{P}(\mathcal{A},\mathcal{B},\rho_{c})(\rho_{e})$ at the output of the
quantum switch is a mixture of the pure states (i.e.
$\langle\psi_{c}||\psi_{c}\rangle=1$ in a complex Hilbert space) of the
control qudit. Based on the Hadamard basis, if the measurement outcome is
$|-\rangle$, the global state collapses with probability $pq$ to,

$\displaystyle\rho^{QS}_{e_{|-\rangle}}=\Big{[}\sum_{s,t}\iota_{t}|st\rangle\langle\lambda_{st}|\Big{]}\rho_{e}\Big{[}\sum_{s,t}\iota_{t}|st\rangle\langle\lambda_{st}|\Big{]}^{\dagger}.$
(19)

The entanglement distribution in this case is a noiseless process. Bob
receives the particle $|\Phi\rangle_{2}$ without any error if
$\rho_{e}^{QS}=\rho_{e}$, thereby yielding to a noiseless entanglement
distribution through the $d$-dimensional quantum switch with probability $pq$.
Alternatively, the distribution can be a noisy process degraded by a noisy
communication channel with a probability $(1-pq)$. In that case, the
measurement outcome is $|+\rangle$ and the global state collapses to,

$\displaystyle\rho^{QS}_{e_{|+\rangle}}=\frac{(1-p)(1-q)\rho_{e}}{1-pq}+\frac{p(1-q)\Big{[}\sum_{s,t}|st\rangle\langle\lambda_{st}|\Big{]}\rho_{e}\Big{[}\sum_{s,t}|st\rangle\langle\lambda_{st}|\Big{]}^{\dagger}}{1-pq}$
$\displaystyle\quad+\frac{(1-p)q\Big{[}\sum_{s,k}|\lambda_{sk}\rangle\langle\lambda_{sk}|\Big{]}\rho_{e}\Big{[}\sum_{s,k}|\lambda_{sk}\rangle\langle\lambda_{sk}|\Big{]}^{\dagger}}{1-pq}.$
(20)

Subsequently we can derive the expression for the density matrix
$\rho_{e}^{QS}$ of the EPR pair distributed between Alice and Bob at the
output of the quantum switch which will be equal to $\rho_{e}$ with
probability $pq$ and otherwise as,

$\displaystyle\rho^{QS}_{e}=\Big{(}(1-p)(1-q)\rho_{e}+p(1-q)\Big{[}\sum_{s,t}|st\rangle\langle\lambda_{st}|\Big{]}\times\rho_{e}\Big{[}\sum_{s,t}|st\rangle\langle\lambda_{st}|\Big{]}^{\dagger}+(1-p)q\Big{[}\sum_{s}\sum_{k=0}^{d-1}|\lambda_{sk}\rangle\langle\lambda_{sk}|\Big{]}$
$\displaystyle\quad\times\rho_{e}\Big{[}\sum_{s}\sum_{k=0}^{d-1}|\lambda_{sk}\rangle\langle\lambda_{sk}|\Big{]}^{\dagger}\Big{)}/(1-pq).$
(21)

With probability $pq$ for the measurement of the control qudit corresponding
to the basis states, the entanglement distribution is noiseless. Bob receives
the particle $|\Phi\rangle_{2}$ without any error if $\rho_{e}^{QS}=\rho_{e}$,
while it is degraded by noise with the probability $(1-pq)$. Since
$\mathcal{A}$ and $\mathcal{B}$ denote the bit-flip and the phase-flip
channels respectively, traversed in a well-defined order
$\mathcal{A}\to\mathcal{B}$, the density matrix of the entangled pair at Bob’s
side is given by,

$\displaystyle\rho^{CT}_{e}$
$\displaystyle=\mathcal{B}\big{[}\mathcal{A}(\rho_{e})\big{]}=\mathcal{B}\big{[}\sum_{s}\mathcal{A}_{s}\rho_{e}\mathcal{A}_{s}^{\dagger}\big{]}$
$\displaystyle=\sum_{t}\mathcal{B}_{t}\big{[}\sum_{s}\mathcal{A}_{s}\rho_{e}\mathcal{A}_{s}^{\dagger}\big{]}\mathcal{B}_{t}^{\dagger}.$
(22)

$\displaystyle\rho^{CT}_{e}$
$\displaystyle=(1-p)(1-q)\rho_{e}+p(1-q)\Big{[}\sum_{s,t}|st\rangle\langle\lambda_{st}|\Big{]}\rho_{e}\Big{[}\sum_{s,t}|st\rangle\langle\lambda_{st}|\Big{]}^{\dagger}+(1-p)q\Big{[}\sum_{s}\sum_{k=0}^{d-1}|\lambda_{sk}\rangle\langle\lambda_{sk}|\Big{]}\rho_{e}$
$\displaystyle\times\Big{[}\sum_{s}\sum_{k=0}^{d-1}|\lambda_{sk}\rangle\langle\lambda_{sk}|\Big{]}^{\dagger}+pq\Big{[}\sum_{s,t}\iota_{t}|st\rangle\langle\lambda_{st}|\Big{]}\rho_{e}\Big{[}\sum_{s,t}\iota_{t}|st\rangle\langle\lambda_{st}|\Big{]}^{\dagger}.$
(23)

Consequently, the expression for the density matrix $\rho_{e}^{CT}$ of the EPR
pair distributed between Alice and Bob for classical trajectories is given by
(4). Using (4, 4, 4), it is possible to calculate the performance gain
achievable through the teleportation of qudit over that achievable through the
teleportation of qubit. Also it is possible to compute the performance gain
achievable through entanglement distribution via quantum switch (for qudits)
over that achievable without any entanglement distribution.

##  5 Qudit Teleportation through Quantum Switch

Using the expression for the density matrix of the teleported qudit at Bob’s
side, we can quantify the effect of noise on EPR distribution on the
teleported qudit. Towards that end, let us assume
$\rho_{\psi}\text{Tr}iangleq|\psi\rangle_{a}\langle\psi|_{a}$ be the $d\times
d$ density matrix of the unknown pure qudit state $|\psi\rangle_{a}$ that
Alice sends to Bob, such that,

$\displaystyle\rho_{\psi}=\begin{bmatrix}\rho^{11}_{\psi}&\cdots&\rho^{1d}_{\psi}\\\
\vdots&\ddots&\vdots\\\ \rho^{d1}_{\psi}&\cdots&\rho^{dd}_{\psi}\end{bmatrix}$
(24)

Let us assume that the expression for the density matrix of the distributed
EPR pair at the output of the quantum switch corrupted by noise be equal to
$\tilde{\rho}_{e}$. If the entanglement distribution is perfect,
$\tilde{\rho}_{e}=\rho_{e}$. Here $\tilde{\rho}_{e}$ is a $d^{2}\times d^{2}$
matrix consisting of $d\times d$ block matrices to yield,

$\displaystyle\tilde{\rho}_{e}=\begin{bmatrix}\tilde{\rho}_{e_{11}}&\cdots&\tilde{\rho}_{e_{1d}}\\\
\vdots&\ddots&\vdots\\\
\tilde{\rho}_{e_{d1}}&\cdots&\tilde{\rho}_{e_{dd}}\end{bmatrix}$ (25)

We also assume that the density matrix of the teleported qudit at Bob’s side
is equal to $\rho_{\tau}$. Using these notations, we represent the quantum
teleportation process of a qudit in terms of density matrices in Fig. 3.

Figure 4: Conditional fidelity of the teleported qudit as a function of the
error probabilities $p$ and $q$ of the noisy bit flip channel where
measurement of the control qudit $|\psi_{c}\rangle$ yields $|+\rangle$. Figure
5: Average fidelity of the teleported qudit as a function of both the error
probabilities $p$ and $q$.

Following the procedure in [15], we first calculate $\rho_{\tau}$ if the
measurement outcome is $|kk\rangle$ (for $k=0,\dotso,d-1$) to obtain,

$\displaystyle\rho^{|kk\rangle}_{\tau}=d\Big{[}\sum_{k=1}^{d}\rho^{kk}_{\psi}\tilde{\rho}_{e_{kk}}+\sum_{v=1}^{j}\sum_{u=j+1}^{d}\rho^{vu}_{\psi}\tilde{\rho}_{e_{vu}}\Big{]}.$
(26)

If the measurement outcome is $|vu\rangle$ (for $v=0,\dotso,j-1$,
$u=j,\dotso,d-1$) $\rho_{\tau}$ can be given by,

$\displaystyle\rho^{|vu\rangle}_{\tau}=d\Big{[}\sum_{k=1}^{d}\rho^{kk}_{\psi}\tilde{\rho}_{e_{kk}}-\sum_{v=1}^{j}\sum_{u=j+1}^{d}\rho^{vu}_{\psi}\tilde{\rho}_{e_{vu}}\Big{]}.$
(27)

Using $\rho_{\tau}$ we can calculate the performance gain achievable through
the superposition of causal orders via the quantum switch in terms of quantum
fidelity $F$. The conditional fidelity $F^{QS}_{|+\rangle}$ associated with
the teleported qudit at Bob’s side can be calculated in terms of $\rho_{\tau}$
to obtain [15],

$\displaystyle
F^{QS}_{|+\rangle}=\text{Tr}\big{[}\rho_{\tau}\rho_{\psi}\big{]}=\frac{d-dp-
dq+pq+1}{(1-pq)(d+1)}$ (28)

while the conditional fidelity $F^{QS}_{|-\rangle}=1$ with probability $pq$,
when the measurement of the control qudit corresponds to the state $|+\rangle$
and $|-\rangle$ respectively. In this case, the control qudit
$|\psi_{c}\rangle=|+\rangle\langle+|$ where
$|+\rangle\equiv\frac{1}{\sqrt{d}}(|1\rangle+|2\rangle+\dotso+|d\rangle)$ or
$|\psi_{c}\rangle=|-\rangle\langle-|$ where
$|-\rangle\equiv\frac{1}{\sqrt{d}}(|1\rangle-|2\rangle-\dotso-|d\rangle)$. The
average fidelity $\bar{F}^{QS}$ of the teleported qudit on Bob’s side using
quantum switch can be given by,

$\displaystyle\bar{F}^{QS}$
$\displaystyle=pqF^{QS}_{|-\rangle}+(1-pq)F^{QS}_{|+\rangle}$
$\displaystyle=\text{Tr}\Bigg{[}d\Big{(}\sum_{k=1}^{d}\rho^{kk}_{\psi}\tilde{\rho}_{e_{kk}}+\sum_{v=1}^{j}\sum_{u=j+1}^{d}\rho^{vu}_{\psi}\tilde{\rho}_{e_{vu}}\Big{)}$
$\displaystyle\times\begin{bmatrix}\rho^{11}_{\psi}&\cdots&\rho^{1d}_{\psi}\\\
\vdots&\ddots&\vdots\\\
\rho^{d1}_{\psi}&\cdots&\rho^{dd}_{\psi}\end{bmatrix}\Bigg{]}$ (29)

which can be manipulated using mathematical operations, like matrix
multiplication, integration and trace of square matrices, to obtain,

$\displaystyle\bar{F}^{QS}=\frac{d-dp-dq+d^{2}pq+1}{d+1}$ (30)

with $p$ and $q$ being the error probabilities of the two considered noisy
channels $\mathcal{A}$ and $\mathcal{B}$. It is worth-mentioning that a
fidelity of 1 corresponds to the case where the received qudit on Bob’s side
is identical to the original qudit and the measurement of the control qudit is
equal to $|-\rangle$. The gain in performance on using quantum switch for
qudit teleportation can be calculated in terms of quantum fidelity and
achievable channel capacity by applying quantum switch with classical
trajectories, a detailed treatise of which we leave for our future work.

In Fig. 4, we plot the conditional fidelity $|F^{QS}_{|+\rangle}$ with respect
to $p$ and $q$, where $|F^{QS}_{|+\rangle}$ is the fidelity corresponding to
the control qudit $|+\rangle$. Fig. 4 demonstrates that the conditional
fidelity deteriorates irreversibly with the increase in error probability $p$
of the bit flip channel and $q$ of the phase flip channel, if the control
qudit $|\psi_{c}\rangle$ is measured to be in state $|+\rangle$ for any
$p,q>0$. This degradation in fidelity results from the degradation in
entanglement and thereby in the quality of the received information.

Fig. 5 portrays the average fidelity as a function of the error probabilities
$p$ and $q$ of the noisy channel, as long as the control qudit
$|\psi_{c}\rangle$ is measured to be in state $|-\rangle$. The noise in the
channel causes a degradation in fidelity up to a certain threshold and then
the fidelity starts increasing with increase in $p$ and $q$. This happens due
to the advantage offered by the qudit switch, which makes it possible to send
information even when the noise increases; a task impossible to achieve within
the standard quantum Shannon theory. This improvement in average fidelity
after a certain point is utilized in our design for quantum switch to yield an
overall improvement in fidelity even in a very noisy quantum channel. Highest
possible fidelity is achieved for the case when $p=0$ and $q=0$; that is the
case where the channels are technically noiseless. Here the threshold value
refers to the values of the error probabilities $p,q$. With further increase
in $p$ and $q$ after the threshold value, average fidelity increases. This
threshold is actually the limiting condition of the zero-capacity quantum
channel, i.e. $p=q=0.5$. For the values of $p$ and $q$, lower than the
threshold value, fidelity decreases with increase in error probabilities. This
means that the advantage in terms of fidelity obtained through employing a
quantum switch is only attainable after a certain threshold for $p$ and $q$ is
reached.

Figure 6: Average fidelity of the teleported qudit as a function of $p$ when
$q=p$ and entanglement distribution is achieved through a quantum switch.

In Fig. 6, we compare the average fidelity achievable with a quantum switch
where a $d$-dimensional qudit is teleported by varying $d$, $d=2,5,10,15,20$
as a function of the error probabilities $p$ and $q$ of the bit-flip and the
phase-flip channels respectively for the condition when $q=p$. When $p=q=0.5$
(like a 50% erasure channel), no quantum information can be sent through any
classical trajectory traversing the channels $\mathcal{A}$ and $\mathcal{B}$
[9]. Hence, achievable fidelity increases with the dimension of the quantum
particle used for teleportation. Increased dimension refers to an enlarged
Hilbert space and a larger Hilbert space offers the advantage of communicating
with larger information capacity and higher noise resilience than traditional
qubit-based quantum systems. It is evident from Fig. 6, average fidelity
${\bar{F}}^{QS}$ decreases as $p$ and $q$ increase (i.e. increase in error
probabilities), up to a certain threshold like $p=q=1/2$ for the case of
$d=2$. Fidelity starts increasing after the threshold even as noise affecting
the channel increases, thereby corroborating the results obtained in [15].

For the case when $p=q\to 0.5$ (the limiting condition of zero-capacity
quantum channel), teleportation of qubit over a quantum switch offers a
fidelity of approaching 0.5, i.e. $\bar{F}^{QS}_{d=2}\to 0.5$. While as we
increase $d$ to 5, fidelity approaches 0.65. If $d=10$,
$\bar{F}^{QS}_{d=10}\to 0.78$, $d=15$, if $\bar{F}^{QS}_{d=15}\to 0.89$,
$d=20$, and if $\bar{F}^{QS}_{d=20}\to 0.94$. Therefore, entanglement
distribution of a qudit over a quantum switch provides a remarkable gain in
fidelity of the teleported qudit at Bob’s side with increase in both noise
plaguing the communication channel and dimensionality of the quantum particle.

For any scenario of quantum teleportation, the EPR pair while being
distributed between transmitter and receiver may lose its coherence, and
become a mixed state due to interaction with the environment and noise
introduced by the quantum channel. This in turn reduces the fidelity of
teleportation, thereby reducing the range of states that can be accurately
teleported. However, the results obtained here demonstrate that fidelity can
be improved even in presence of noise, if a quantum switch is used. Towards
this end, if entanglement is distributed between sender and receiver using a
quantum switch, fidelity of teleportation can be enhanced even in presence of
noise, thereby offering a wide range of states that can be accurately
teleported.

##  6 Conclusion

In this paper, we formulated the first ever theoretical framework for
implementing quantum switch for $d$-dimensional quantum teleportation system
from a communication engineering point of view. With the help of generalized
quantum gates, we have designed the basic scheme for teleportation of qudits
and then applied it to quantum switch for entanglement distribution. In
future, we will use the results obtained in this paper directly to evaluate
the performance gain achievable by teleporting and distributing entanglement
of qudits over that of qubits in noisy quantum channels, in terms of
information carrying capacity and resilience to noise. We will also analyze
the impact of designing the generalized SWAP gate (different combinations of
higher dimensional quantum gates as demonstrated in Fig. 2) on the overall
performance of qudit teleportation in our future work.

## References

  * [1] L. Hanzo, H. Haas, S. Imre, D. O’Brien, M. Rupp, and L. Gyongyosi, “Wireless myths, realities, and futures: from 3G/4G to optical and quantum wireless, " _Proceedings of the IEEE_ , vol. 100, no. 13, pp. 1853–1888, 2012.
  * [2] C. H. Bennett, G. Brassard, C. Crpeau, R. Josza, A. Peres, and W. K. Wootters, “Teleporting an unknown quantum state via dual classical and Einstein-Podolsky-Rosen channels," _Phys. Rev. Lett._ , vol. 70, pp. 1895–1899, 1993.
  * [3] R. Matsumoto, “Improved Gilbert–Varshamov Bound for Entanglement-Assisted Asymmetric Quantum Error Correction by Symplectic Orthogonality," in _IEEE Trans. Quantum Engineering_ , vol. 1, Apr. 2020, 4 pages, Early Access, DOI: 10.1109/TQE.2020.2978454.
  * [4] S. Salek, D. Ebler, and G. Chiribella, “Quantum communication in a superposition of causal orders," Sep. 2018, _arXiv:1809.06655_. [Online]. Available: https://arxiv.org/abs/1809.06655.
  * [5] G. Chiribella and H. Kristjánsson, “Quantum Shannon theory with superpositions of trajectories," _Proc. Roy. Soc. A_ , vol. 475, no. 2225, May 2019, Art. no. 20180903.
  * [6] A. A. Abbott, J. Wechs, D. Horsman, M. Mhalla, and C. Branciard, “Communication through coherent control of quantum channels," Oct. 2018, _arXiv:1810.09826_. [Online]. Available: https://arxiv.org/abs/1810.09826.
  * [7] G. Chiribella, “Perfect discrimination of no-signalling channels via quantum superposition of causal structures," _Phys. Rev. A, Gen. Phys._ , vol. 86, Oct. 2012, Art. no. 040301(R).
  * [8] K. Chakraborty, D. Elkouss, B Rijsman and S. Wehner, “Entanglement Distribution in a Quantum Network: A Multicommodity Flow-Based Approach," in _IEEE Trans. Quantum Engineering_ , vol. 1, Oct. 2020, Early Access, DOI: 10.1109/TQE.2020.3028172.
  * [9] A. S. Cacciapuoti, M. Caleffi, R. Van Meter, and L. Hanzo, “When entanglement meets classical communications: Quantum teleportation for the quantum Internet (invited paper)," Jul. 2019, _arXiv:1907.06197._ [Online]. Available: https://arxiv.org/abs/1907.06197.
  * [10] D. T. Pegg, L. S. Phillips and S. M. Barnett, “Optical state truncation by projection synthesis." _Phys. Rev. Lett._ , vol. 81, pp. 1604–1606, 1998).
  * [11] S. K. Goyal, P. E. Boukama-Dzoussi, S. Ghosh, F. S. Roux, and T. Konrad, “Qudit-teleportation for photons with linear optics." _arXiv:1212:5115_ , 2012.
  * [12] K. Ozdemir, A. Miranowicz, M. Koashi, and N. Imoto, “Pulse-mode quantum projection synthesis:‘Effects of mode mismatch on optical state truncation and preparation." _Phys. Rev. A_ , vol. 66, 2002.
  * [13] M. Lanzagorta, _Underwater Communications_ , Synthesis Lectures on Communications vol. 5, Morgan & Claypool Publishers, San Rafael, CA 2012.
  * [14] C. H. Bennett, G. Brassard, C. Crépeau, R. Jozsa, A. Peres, and W. K. Wootters, “Teleporting an unknown quantum state via dual classical and Einstein-Podolsky-Rosen channels," _Phys. Rev. Lett._ , vol. 70, no. 13, pp. 1895–1899, Jul. 2002.
  * [15] M. Caleffi and A. S. Cacciapuoti,“Quantum Switch for the Quantum Internet: Noiseless Communications Through Noisy Channels," _IEEE Journal on Selected Areas in Comms._ , vol. 38, no. 3, pp. 575–588, Mar. 2020.
  * [16] C. Zhang, J. F. Chen, C. Cui, J. P. Dowling, Z. Y. Ou, and T. Byrnes, “Quantum teleportation of photonic qudits using linear optics", _Phys. Rev. A_ , vol. 100, no. 3, 6 pages, Sep. 2019.
  * [17] S. Balakrishnan, “Various Constructions of Qudit SWAP Gate," _Hindawi Physics International Research_ , vol. 2014, 5 pages, Available: https://doi.org/10.1155/2014/479320.
  * [18] K. Fujii, “Exchange gate on the qudit space and Fock space," _Journal of Optics B_ , vol. 5, no. 6, pp. S613–S618, 2003.
  * [19] G. Paz-Silva, S. Rebic, J. Twamley, and T. Duty, “Perfect mirror transport protocol with higher dimensional quantum chains," _Physical Review Letters_ , vol. 102, Article ID 020503, 2009.
  * [20] X. Wang, “Continuous-variable and hybrid quantum gates," _Journal of Physics A: Mathematical and General_ , vol. 34, no. 44, pp. 9577–9584, 2001.
  * [21] N. D. Mermin, “From classical state swapping to quantum teleportation," _Physical Review A_ , vol. 65, no. 1, Article ID 012320, 2001.
  * [22] A. A.-Saki, M. Alam and S. Ghosh, “Experimental Characterization, Modeling, and Analysis of Crosstalk in a Quantum Computer," in _IEEE Trans. Quantum Engineering_ , vol. 1, Sept. 2020, Early Access, DOI: 10.1109/TQE.2020.3023338.