# Scaling rule in critical non-Hermitian skin effect

Kazuki Yokomizo  Condensed Matter Theory Laboratory, RIKEN, 2-1 Hirosawa,
Wako, Saitama, 351-0198, Japan  Shuichi Murakami  Department of Physics, Tokyo
Institute of Technology, 2-12-1 Ookayama, Meguro-ku, Tokyo, 152-8551, Japan
TIES, Tokyo Institute of Technology, 2-12-1 Ookayama, Meguro-ku, Tokyo,
152-8551, Japan

###### Abstract

Non-Hermitian systems show a non-Hermitian skin effect, where the bulk states
are localized at a boundary of the systems with open boundary conditions. In
this paper, we study the dependence of the localization length of the
eigenstates on a system size in a specific non-Hermitian model with a critical
non-Hermitian skin effect, where the energy spectrum undergoes discontinuous
transition in the thermodynamic limit. We analytically show that the
eigenstates exhibit remarkable localization, known as scale-free localization,
where the localization length is proportional to a system size. Our result
gives theoretical support for the scale-free localization, which has been
proposed only numerically in previous works.

###### pacs:

Valid PACS appear here

††preprint: APS/123-QED

##  I Introduction

Nonequilibrium and open systems are effectively described in terms of non-
Hermitian Hamiltonians. Thus far, non-Hermitian systems have been attracting
much attention to study remarkable phenomena, such as unidirectional
transmission [1, 2] and retroreflection [3]. Among theoretical and
experimental studies on non-Hermitian systems, a non-Hermitian skin effect
plays an important role in condensed matter physics. In the non-Hermitian skin
effect, eigenstates of the bulk are localized at the end of a system with an
open boundary condition [4, 5, 6, 7, 8, 9, 10, 11]. Furthermore, this effect
gives rise to the difference between energy spectra in an open chain and those
in a periodic chain. These phenomena have been experimentally observed in
various physical systems [12, 13, 14, 15, 16, 17, 18]. It is significant to
investigate some phenomena induced by the non-Hermitian skin effect because
this effect can lead to non-Hermitian phenomena [19, 20, 21, 22, 23, 24, 25,
26, 27].

Additionally, in recent years, non-Hermitian systems have been much
investigated from various perspectives. For example, Ref. [28] theoretically
proposed a class of criticality. The authors showed a critical non-Hermitian
skin effect, where an energy spectrum and localization of eigenstates
discontinuously jump across a critical point. We note that, since the critical
non-Hermitian skin effect occurs in the thermodynamic limit, this effect can
be systematically understood in terms of a non-Bloch band theory [4, 29, 30,
31, 32, 24, 33, 34]. This is because, in non-Hermitian systems with open
boundary conditions, the non-Bloch band theory can calculate continuum energy
bands in the limit of a large system size. On the other hand, Ref. [28] also
proposed that when a system with the critical non-Hermitian skin effect has a
finite system size, the eigenstates exhibit remarkable localization behavior.
It is called scale-free localization, where the localization length is
proportional to the system size. Importantly, the non-Bloch band theory cannot
study the scale-free localization because this localization occurs only in a
finite open chain. Therefore, the nature of the scale-free localization has
not been revealed yet.

Thus far, finite non-Hermitian systems have been investigated only by a
numerical calculation. Hence, it is necessary to study non-Hermitian systems
with a finite system size in more detail to reveal the nature of non-Hermitian
phenomena unique to a finite system, such as the scale-free localization.
Furthermore analytical studies of such systems will give us useful insights
into experimental studies on the non-Hermitian skin effect.

In this paper, we analyze a scaling rule of the localization length of the
eigenstates in a ladder model. As a result, we find remarkable phenomena
unique to a finite open non-Hermitian system as follows. When the system size
is sufficiently small, the energy eigenvalues and the eigenstates are
insensitive to the change of the system size. On the other hand, when the
system size is larger than a critical value, the localization length changes
in proportion to the system size, and the energy eigenvalues also change as
the system size changes. Then in our analysis, it is shown that these
eigenstates exhibit the scale-free localization. Furthermore, we can derive a
formula of the critical value of system size.

This paper is organized as follows. In Sec. II, we introduce the ladder model
and describe the critical non-Hermitian skin effect. Furthermore we explain
the concept of our study on a finite open chain. In Sec. III, we derive an
analytical expression of the localization length and the critical system size
in the ladder model. Then we confirm that our result can precisely predict the
critical value compared with the result of Ref. [35]. Finally, we summarize
the results of this paper and comment on extension to general non-Hermitian
systems.

##  II Critical non-Hermitian skin effect

Figure 1: (a) Ladder model in an infinite chain. The red and blue circles
represent the sites in sublattices ${\rm A}$ and ${\rm B}$, respectively. We
set the hopping amplitude to be $t_{L}>t_{R}>0$. The parameter $\Delta$
represents the coupling between two Hatano-Nelson models. (b) Energy spectra
of the ladder model in an infinite chain with open (red) and periodic (blue)
boundary conditions. The energy spectra in an open chain discontinuously jump
across the critical point $\Delta=0$.

In this section, we briefly review the critical non-Hermitian skin effect.
First, we introduce a tight-binding model on a one-dimensional ladder, as
shown in Fig. 1(a). It consists of two one-dimensional chains with lattice
sites ${\rm A}$ and ${\rm B}$. The real-space Hamiltonian of this system is
given by

$\displaystyle H$ $\displaystyle=$ $\displaystyle\sum_{n}\left(t_{L}c_{n,{\rm
A}}^{\dagger}c_{n+1,{\rm A}}+t_{R}c_{n+1,{\rm A}}^{\dagger}c_{n,{\rm
A}}+\Delta c_{n,{\rm A}}^{\dagger}c_{n,{\rm B}}\right.$ (1)
$\displaystyle\left.+t_{R}c_{n,{\rm B}}^{\dagger}c_{n+1,{\rm
B}}+t_{L}c_{n+1,{\rm B}}^{\dagger}c_{n,{\rm B}}+\Delta c_{n,{\rm
B}}^{\dagger}c_{n,{\rm A}}\right),$

where $c_{n,\mu}~{}(\mu={\rm A},{\rm B})$ represents a fermionic annihilation
operator at the $n$th unit cell in the $\mu$ sublattice. For simplicity, the
parameters $t_{L},t_{R}$, and $\Delta$ are set to be real, and we assume
$t_{L}>t_{R}>0$. Below, we explain the behavior of the ladder model in the
limit of a large system size.

When $\Delta=0$, the system is decoupled into two Hatano-Nelson models [36].
In this case, the two Hatano-Nelson models exhibit the non-Hermitian skin
effect, and the bulk eigenstates in one chain and those in the other chain are
localized at the opposite ends of the chain. Furthermore, the energy spectrum
in an open chain takes different values from that in a periodic chain, as
shown in Fig. 1(b2). On the other hand, by coupling the two chains with
$\Delta\neq 0$, the non-Hermitian skin effect disappears because the two skin
modes are coupled to each other. Then in contrast to the case of $\Delta=0$,
the energy spectrum in an open chain coincides with that in a periodic chain,
as shown in Figs. 1(b1) and 1(b3). Thus, the energy spectrum in the ladder
model with open boundary conditions discontinuously jumps across the critical
point $\Delta=0$. In Ref. [28], this discontinuous phase transition was
explained in terms of the critical non-Hermitian skin effect.

In contrast, when the system size $L$ is finite, the energy eigenvalues behave
continuously as a function of $\Delta$. Depending on the order of the limits
$L\rightarrow\infty$ and $\Delta\rightarrow 0$, the behavior of the system is
critically different [35]. Namely, after taking the thermodynamic limit, the
infinitesimal coupling discontinuously changes the energy spectrum [Fig.
2(a)]. On the other hand, under $0<\Delta\ll{\cal O}\left(t_{L},t_{R}\right)$,
an increase of the system size amplifies the coupling strength. As a result,
the energy eigenvalues in an open chain gradually change as the system size
increases, and finally, they approach those in a periodic chain in the limit
of $L\rightarrow\infty$ [Fig. 2(b)]. In this paper, we theoretically
investigate the energy eigenvalues and the corresponding eigenstates of the
ladder model as the system size changes.

Figure 2: Energy eigenvalues of the ladder model (a) in an infinite open chain
and (b) in a finite open chain. (a) In the limit of a large system size, by
adding the coupling parameter $\Delta$, the energy spectrum discontinuously
jumps. (b) When $\Delta\neq 0$, the energy eigenvalues in a finite open chain
(orange) gradually approach those in a periodic chain (gray) as the system
size increases. This change is represented by the cyan arrows.

##  III Scaling rule

In this section, first, we show a way to get the energy eigenvalues in the
ladder model. Next, we obtain system size dependence of the solutions of the
characteristic equation for $\beta={\rm e}^{ik}$, where $k$ corresponds to the
Bloch wave number. Accordingly, it is shown that the eigenstates exhibit the
scale-free localization in the ladder model with a finite system size. We also
reveal dependence of the corresponding energy eigenvalues on the system size.
In the following, we assume $0<\Delta\ll{\cal O}\left(t_{L},t_{R}\right)$ for
simplicity.

###  III.1 Periodic and open chains

In the ladder model with a periodic boundary condition, the energy eigenvalues
$E_{\rm PBC}$ are obtained from the Bloch Hamiltonian:

$\displaystyle H\left(k\right)=\left(\begin{array}[]{cc}t_{L}{\rm
e}^{ik}+t_{R}{\rm e}^{-ik}&\Delta\vspace{5pt}\\\ \Delta&t_{R}{\rm
e}^{ik}+t_{L}{\rm e}^{-ik}\end{array}\right),$ (4)

as

$\displaystyle E_{\rm PBC}^{\pm}\left(k\right)=\left(t_{L}+t_{R}\right)\cos
k\pm i\sqrt{\left(t_{L}-t_{R}\right)^{2}\sin^{2}k-\Delta^{2}},$ (5)

where $k$ is the real Bloch wave number.

On the other hand, in the ladder model with an open boundary condition, the
energy eigenvalues are calculated from a real-space Schrödinger equation. In
the following, let $E_{\rm OBC}$ denote the energy eigenvalues in the open
chain. Then we can explicitly write the equation $H\ket{\psi}=E_{\rm
OBC}\ket{\psi}$, where $\ket{\psi}=\left(\dots,\psi_{1,{\rm A}},\psi_{1,{\rm
B}},\dots,\psi_{n,{\rm A}},\psi_{n,{\rm B}},\dots\right)^{\rm T}$, as

$\displaystyle\left\\{\begin{array}[]{l}t_{R}\psi_{n-1,{\rm
A}}+\Delta\psi_{n,{\rm B}}+t_{L}\psi_{n+1,{\rm A}}=E_{\rm OBC}\psi_{n,{\rm
A}},\vspace{3pt}\\\ t_{L}\psi_{n-1,{\rm B}}+\Delta\psi_{n,{\rm
A}}+t_{R}\psi_{n+1,{\rm B}}=E_{\rm OBC}\psi_{n,{\rm B}}.\end{array}\right.$
(8)

According to the theory of linear difference equations, we can take an ansatz
for the eigenstates as a linear combination:

$\displaystyle\left(\begin{array}[]{c}\psi_{n,{\rm A}}\vspace{5pt}\\\
\psi_{n,{\rm
B}}\end{array}\right)=\sum_{j=1}^{4}\left(\beta_{j}\right)^{n}\left(\begin{array}[]{c}\phi_{\rm
A}^{\left(j\right)}\vspace{5pt}\\\ \phi_{\rm
B}^{\left(j\right)}\end{array}\right).$ (13)

Hence Eq. (8) can be rewritten as

$\displaystyle\left(\begin{array}[]{cc}t_{L}\beta+t_{R}\beta^{-1}&\Delta\vspace{5pt}\\\
\Delta&t_{R}\beta+t_{L}\beta^{-1}\end{array}\right)\left(\begin{array}[]{c}\phi_{\rm
A}\vspace{5pt}\\\ \phi_{\rm B}\end{array}\right)=E_{\rm
OBC}\left(\begin{array}[]{c}\phi_{\rm A}\vspace{5pt}\\\ \phi_{\rm
B}\end{array}\right),$ (20) (21)

with $\beta=\beta_{j}$, and
$\phi_{\alpha}=\phi_{\alpha}^{\left(j\right)}~{}(\alpha={\rm A},{\rm B})$.
Here, we can obtain the characteristic equation as the condition that the
coefficients $\phi_{\rm A}$ and $\phi_{\rm B}$ take nonzero values, explicitly
written as

$\displaystyle\beta^{2}$ $\displaystyle-$
$\displaystyle\left(\frac{1}{t_{L}}+\frac{1}{t_{R}}\right)E_{\rm
OBC}\beta+\frac{1}{t_{L}t_{R}}\left(t_{L}^{2}+t_{R}^{2}+E_{\rm
OBC}^{2}-\Delta^{2}\right)$ (22) $\displaystyle-$
$\displaystyle\left(\frac{1}{t_{L}}+\frac{1}{t_{R}}\right)E_{\rm
OBC}\beta^{-1}+\beta^{-2}=0.$

Now we assume that the system has $L$ unit cells, and we impose the open
boundary condition $\psi_{0,\mu}=\psi_{L+1,\mu}=0~{}(\mu={\rm A},{\rm B})$ to
this system. Then from Eq. (8), we can derive the boundary equation for the
solutions of the characteristic equation (22) as

$\displaystyle
X_{1,4}X_{2,3}\left[\left(\beta_{1}\beta_{4}\right)^{L+1}+\left(\beta_{2}\beta_{3}\right)^{L+1}\right]$
$\displaystyle-
X_{1,3}X_{2,4}\left[\left(\beta_{1}\beta_{3}\right)^{L+1}+\left(\beta_{2}\beta_{4}\right)^{L+1}\right]$
$\displaystyle+X_{1,2}X_{3,4}\left[\left(\beta_{1}\beta_{2}\right)^{L+1}+\left(\beta_{3}\beta_{4}\right)^{L+1}\right]=0,$
(23)

where $\beta_{j}~{}(j=1,\dots,4)$ satisfy
$\left|\beta_{1}\right|\leq\dots\leq\left|\beta_{4}\right|$, and
$X_{i,j}~{}(i,j=1,\dots,4)$ are defined as

$X_{i,j}=t_{L}\left(\beta_{j}-\beta_{i}\right)+t_{R}\left(\beta_{j}^{-1}-\beta_{i}^{-1}\right),~{}\left(i,j=1,\dots,4\right).$
(24)

We note that the solutions satisfy

$\beta_{1}=\frac{1}{\beta_{4}},~{}\beta_{2}=\frac{1}{\beta_{3}}$ (25)

because Eq. (22) is a reciprocal equation for $\beta$. In Appendix A, we give
the derivation of Eq. (23) in the ladder model. Since the characteristic
equation (22) and the boundary equation (23) determine a set of solutions for
$\left(\beta,E\right)$, the energy eigenvalues in the finite open chain can be
calculated by combining these equations.

We note that, in the ladder model with $\Delta\neq 0$, the energy spectrum in
an open chain is equal to that in a periodic chain in the thermodynamic limit.
This can be understood in terms of the non-Bloch band theory [29, 31]. In
Appendix B, we show a way to get the energy spectrum in an infinite open
chain.

###  III.2 Scale-free localization

Now we focus on the boundary equation (23). In the following, we assume that
the system size is large. In this case, we can approximate Eq. (23) by
neglecting the terms other than the two dominant terms
$-X_{1,3}X_{2,4}\left(\beta_{2}\beta_{4}\right)^{L+1}$ and
$-X_{1,2}X_{3,4}\left(\beta_{3}\beta_{4}\right)^{L+1}$ on the left-hand side.
By substituting the solutions of the characteristic equation (22) into this
approximated boundary equation, we can rewrite Eq. (23) as

$\left(\beta_{2}\right)^{2L+2}\simeq\frac{\left(t_{L}+t_{R}\right)^{2}}{\left(t_{L}-t_{R}\right)^{2}}\frac{E_{\rm
OBC}^{2}-4t_{L}t_{R}}{\left[\left(t_{L}+t_{R}\right)^{2}-E_{\rm
OBC}^{2}\right]^{2}}\Delta^{2}$ (26)

up to the second order of the coupling parameter $\Delta$. We explain a way to
derive Eq. (26) in Appendix C.

Next, we estimate the right-hand side in Eq. (26). As mentioned in Sec. II, in
the thermodynamic limit, the energy spectrum in an open chain coincides with
that in a periodic chain. Accordingly, in a system with a large system size,
the energy eigenvalues in an open chain asymptotically approach those in a
periodic chain, as shown in Fig. 2(b). Hence, we can approximate the energy
eigenvalues $E_{\rm OBC}$ by the energy eigenvalues $E_{\rm
PBC}^{\pm}\left(k\right)$ [Eq. (5)] in a periodic chain. Here, we take $E_{\rm
PBC}^{+}=\left(t_{L}+t_{R}\right)\cos
k+i\sqrt{\left(t_{L}-t_{R}\right)^{2}\sin^{2}k-\Delta^{2}}~{}(0\leq
k\leq\pi/2)$ as $E_{\rm OBC}$. As a result, we can get an analytical
expression of the absolute value of $\beta_{2}$ as

$\left|\beta_{2}\right|\simeq\left[\frac{t_{L}+t_{R}}{t_{L}-t_{R}}f\left(k\right)\Delta\right]^{\frac{1}{L+1}},$
(27)

where

$f\left(k\right)=\frac{\sqrt{\left|t_{L}^{2}{\rm e}^{2ik}+t_{R}^{2}{\rm
e}^{-2ik}-2t_{L}t_{R}\right|}}{\left|2\left(t_{L}^{2}+t_{R}^{2}\right)\sin^{2}k-i\left(t_{L}^{2}-t_{R}^{2}\right)\sin
2k\right|}.$ (28)

Equation (27) gives the system size dependence of the solutions of the
characteristic equation (22), which determines localization behavior of the
eigenstates. This is because, in a large system, $\beta_{2}(=1/\beta_{3})$
determines the localization length of the eigenstates given by Eq. (13). Hence
we can obtain the localization length of the eigenstates from Eq. (27). For
example, when $k=\pi/2$, which corresponds to the energy eigenvalue with a
maximum imaginary part [Fig. 3(a)], Eq. (27) can be explicitly written as

$\left|\beta_{2}\right|\simeq\left[\frac{\left(t_{L}+t_{R}\right)^{2}}{2\left(t_{L}-t_{R}\right)\left(t_{L}^{2}+t_{R}^{2}\right)}\Delta\right]^{\frac{1}{L+1}}.$
(29)

Therefore, its localization length, given by $-1/\log\left|\beta_{2}\right|$,
is proportional to the system size, and this localization length becomes
larger as the system size increases. This behavior of the eigenstates is
called the scale-free localization [28, 37]. In fact, our numerical result of
exact diagonalization shows that the localization lengths of the eigenstates
divided by the system sizes are the same for various system sizes [Fig. 3(b)].
Similar scale-free localization appears when $k\neq 0$. On the other hand,
when $k=0$ in Eq. (27), $\beta_{2}$ approximately becomes $1$ to the same
order in Eq. (29). This indicates that the corresponding eigenstate is
delocalized, and it does not exhibit the scale-free localization.

Figure 3: Energy eigenvalues and eigenstates in the ladder model with
$t_{L}=1,t_{R}=0.5$, and $\Delta=10^{-3}$. (a) Energy eigenvalues in an open
chain with the system size $L=30$ ($E_{\rm OBC}$, orange dots) and those in a
periodic chain ($E_{\rm PBC}$, black line). The red dot on the black line
expresses the maximum value of the imaginary part of $E_{\rm PBC}$ given by
$k=\pi/2$ in Eq. (5). (b) Spatial distribution of the density of states
corresponding to $\max\left[\Im\left(E_{\rm OBC}\right)\right]$ for various
system sizes. Their localization lengths divided by the system sizes are
equal, which means that the eigenstates have the scale-free localization.

###  III.3 Energy eigenvalues

We can also get the system size dependence of the energy eigenvalues from Eq.
(27). For simplicity, we focus only on the case with $k=\pi/2$ corresponding
to the maximum value of an imaginary part of the energy eigenvalues
$\max\left[\Im\left(E_{\rm OBC}\right)\right]$ in the open chain [Fig. 3(a)].
Now by substituting Eq. (29) into Eq. (22), $\max\left[\Im\left(E_{\rm
OBC}\right)\right]$ can be written as

$\max\left[\Im\left(E_{\rm OBC}\right)\right]\simeq
t_{L}\left[f\left(\frac{\pi}{2}\right)\Delta\right]^{\frac{1}{L+1}}-t_{R}\left[f\left(\frac{\pi}{2}\right)\Delta\right]^{-\frac{1}{L+1}}.$
(30)

As mentioned in Sec. II, in the limit of $L\rightarrow\infty$, it is expected
to approach the maximum value of an imaginary part of the energy eigenvalues
in a periodic chain given by

$E_{0}=\sqrt{\left(t_{L}-t_{R}\right)^{2}-\Delta^{2}}.$ (31)

In Fig. 4, we show the value of $E_{0}-\max\left[\Im\left(E_{\rm
OBC}\right)\right]$ both from the asymptotic formula in Eq. (30) and from
numerical diagonalization. We find that, when the system size is sufficiently
large, Eq. (30) agrees with the result from numerical diagonalization. On the
other hand, when the system size is smaller than a critical value, these
results deviate from each other, and in fact, the energy eigenvalues in a
finite open chain are independent of the system size. In Sec. III.4, we
investigate this evolution and the critical system size.

Figure 4: System size dependence of the energy eigenvalues in the ladder model
with the open boundary condition. The parameters are set to be
$t_{L}=1,t_{R}=0.5$, and $\Delta=10^{-3}$. We show the result from the
asymptotic formula in Eq. (30) in a red line and that from numerical
diagonalization in blue dots. $E_{0}$ is given in Eq. (31). $L_{c}$ is the
critical value of our result given in Eq. (32), and $L_{c}^{\prime}$ is the
critical value given in Eq. (33) from Ref. [35]. When the system size is
smaller than $L_{c}$, Eq. (30) deviates from the result from numerical
diagonalization.

###  III.4 Critical value of system size

Now we derive an expression of the critical system size as mentioned in Sec.
III.3. When the system size is sufficiently small, the skin mode of one
Hatano-Nelson model in the ladder model is decoupled to the skin mode of the
other. We note that the decay of one skin mode is described by the value of
$\bar{\beta}=\sqrt{t_{R}/t_{L}}$. In fact, its decay length is equal to
$-1/\log\bar{\beta}$. Then the energy eigenvalues and the eigenstates are
insensitive to the change of the system size. On the other hand, when the
system size is sufficiently large, the two skin modes are coupled to each
other, and the eigenstates in the ladder model exhibit the scale-free
localization. As mentioned in Sec. II, the coupling strength between the two
skin modes localized at the opposite ends is enhanced exponentially for a
large system size. According to Ref. [28], we can estimate the scaled coupling
as $\bar{\Delta}\sim\bar{\beta}^{-L/2}$. It may look counterintuitive because
the skin modes at the left and right ends become faraway from each other as
the system size becomes larger.

To understand this counterintuitive behavior, we note that the skin modes form
the energy spectra and correspond to not edge modes but bulk modes, although
they are localized at either end of the system. Namely, in an open chain, the
energy spectra of the skin modes are insensitive to boundary conditions, which
indeed means that the skin mode reflects the bulk physics. Thus, while the
skin modes of each Hatano-Nelson model are localized at either end, they are
sensitive to the bulk physics, and the coupling $\Delta$ makes the non-
Hermitian skin effect disappear. Based on these observations, it is natural
that the increase of the system size amplifies the effect of the coupling.

Importantly, the scale-free localization competes with the localization of the
original skin mode described by $\bar{\beta}$. For example, in the case of
$k=\pi/2$, the degree of the delocalization of the eigenstates is obtained
from Eq. (29) because the more $\left|\beta_{2}\right|$ approaches $1$, the
broader the eigenstates extend over the whole system. In conclusion, the
scale-free localization occurs when Eq. (29) exceeds
$\bar{\beta}=\sqrt{t_{R}/t_{L}}$. From the above discussion, we can get the
critical system size as

$L_{c}=\frac{2}{\log\left(t_{R}/t_{L}\right)}\log\left[\frac{\left(t_{L}+t_{R}\right)^{2}\Delta}{2\left(t_{L}-t_{R}\right)\left(t_{L}^{2}+t_{R}^{2}\right)}\right]-1.$
(32)

In Ref. [35], the critical value was evaluated in a qualitative manner, given
by

$L_{c}^{\prime}=\frac{2\log\Delta}{\log\left(t_{R}/t_{L}\right)}.$ (33)

In Fig. 4, we show the two critical values for Eqs. (32) and (33) in green and
black dashed lines, respectively. We find that our result in Eq. (33) can more
precisely predict the change from the localization of the original skin mode
to the scale-free localization in the ladder model than the result of the
previous work in Eq. (32).

##  IV Summary

In this paper, we study the localization phenomenon in the ladder model. As a
manifestation of the critical non-Hermitian skin effect, the eigenstates
exhibit the scale-free localization in a finite open chain. Accordingly, the
energy eigenvalues gradually change as the system size increases. Furthermore,
our analysis can predict the evolution from the localization stemming from the
non-Hermitian skin effect to the scale-free localization as the system size
becomes larger. It is worth noting that the non-Bloch band theory cannot
describe such behavior in the ladder model because the scale-free localization
does not occur in the thermodynamic limit.

We discuss universality of the scale-free localization. We focus on non-
Hermitian multiband tight-binding systems with open boundary conditions. In
this case, we have the characteristic equation $\det\left[{\cal
H}\left(\beta\right)-E\right]=0$, and it is an algebraic equation for $\beta$
with an even degree $2M$ in general. Then the behavior of the systems is
determined by $\beta_{M}$ and $\beta_{M+1}$, which are the solutions of the
characteristic equation with
$\left|\beta_{1}\right|\leq\dots\leq\left|\beta_{2M}\right|$. Here, we assume
that the systems have some symmetries which reverse the direction of a chain,
such as inversion symmetry, as is the case in our ladder model. In this case,
the characteristic equation is a reciprocal equation, and it follows that
$\beta_{j}=1/\beta_{2M+1-j}~{}(j=1,\dots,M)$. Hence, we obtain
$\beta_{M}=1/\beta_{M+1}$, and from open boundary conditions, we can get the
system size dependence of $\beta_{M}$ as

$\left|\beta_{M}\right|=\left|G\left(\beta_{1},\dots,\beta_{2M},E_{\rm
OBC}\right)\right|^{\frac{1}{2L}},$ (34)

In Appendix D, we explain a way to derive this equation. Equation (34)
indicates that the eigenstates are localized with the localization length
$-1/\log\beta_{M}$ proportional to a system size. Therefore, the eigenstates
exhibit the scale-free localization. In terms of the non-Bloch band theory,
this condition is equivalent to that the trajectories of $\beta_{M}$ and
$\beta_{M+1}$, called the generalized Brillouin zone, become a unit circle. We
note that the method of the self-consistent calculation discussed in Sec.
III.2 is useful to analyze behavior of localization in eigenstates.

In this model, the generalized Brillouin zone is given by
$\left|\beta_{M}\right|=\left|\beta_{M+1}\right|=1$ because of the inversion
symmetry. From our derivation of the scale-free localization, our theory also
applies to a wide class of non-Hermitian systems where
$\left|\beta_{M}\right|=\left|\beta_{M+1}\right|=1$ in the limit of a large
system size is guaranteed by some symmetries. Namely, in such systems, if we
include a correction due to a finite system size, $\left|\beta_{M}\right|$ and
$\left|\beta_{M+1}\right|$ deviate from 1, and it necessarily leads to the
scale-free localization in the same way as in the present case.

Recently, various non-Hermitian systems with the non-Hermitian skin effect
were experimentally realized. Therefore, by combining two non-Hermitian
systems, we expect the realization of a non-Hermitian system which has
eigenstates with the scale-free localization in various physical systems,
e.g., an electrical circuit system.

###### Acknowledgements.

This paper was supported by JSPS KAKENHI (Grant No. JP18H03678) and the MEXT
Elements Strategy Initiative to Form Core Research Center (TIES) (Grant No.
JPMXP0112101001). K.Y. was also supported by JSPS KAKENHI Grants No.
JP18J22113 and No. JP21J01409.

##  Appendix A Derivation of Eq. (23)

In this appendix, we describe the derivation of the boundary equation (23) in
the ladder model with the system size $L$. From the real-space eigenequation
in Eq. (8) and the open boundary condition
$\psi_{0,\mu}=\psi_{L+1,\mu}=0~{}(\mu={\rm A},{\rm B})$, we can get the
equations for the eigenstates in real space as

$\displaystyle\left\\{\begin{array}[]{l}\Delta\psi_{1,{\rm
B}}+t_{L}\psi_{2,{\rm A}}=E_{\rm OBC}\psi_{1,{\rm A}},\vspace{5pt}\\\
\Delta\psi_{1,{\rm A}}+t_{R}\psi_{2,{\rm B}}=E_{\rm OBC}\psi_{1,{\rm
B}},\vspace{5pt}\\\ t_{R}\psi_{L-1,{\rm A}}+\Delta\psi_{L,{\rm B}}=E_{\rm
OBC}\psi_{L,{\rm A}},\vspace{5pt}\\\ t_{L}\psi_{L-1,{\rm
B}}+\Delta\psi_{L,{\rm A}}=E_{\rm OBC}\psi_{L,{\rm B}}.\end{array}\right.$
(39)

Now Eq. (39) can be rewritten into coupled equations for the coefficients
$\phi_{\mu}^{\left(j\right)}~{}(\mu={\rm A},{\rm B},j=1,\dots,4)$ by
substituting the general solution in Eq. (13). Furthermore, by using the bulk
eigenequation in Eq. (21), we can get the coupled equations including only
$\phi_{\rm A}^{\left(j\right)}~{}(j=1,\dots,4)$. As a result, we can get the
condition that $\phi_{\rm A}^{\left(j\right)}~{}(j=1,\dots,4)$ have nonzero
values, written as

$\displaystyle\left|\begin{array}[]{cccc}1&1&1&1\vspace{5pt}\\\
X_{1}&X_{2}&X_{3}&X_{4}\vspace{5pt}\\\
\left(\beta_{1}\right)^{L+1}&\left(\beta_{2}\right)^{L+1}&\left(\beta_{3}\right)^{L+1}&\left(\beta_{4}\right)^{L+1}\vspace{5pt}\\\
X_{1}\left(\beta_{1}\right)^{L+1}&X_{2}\left(\beta_{2}\right)^{L+1}&X_{3}\left(\beta_{3}\right)^{L+1}&X_{4}\left(\beta_{4}\right)^{L+1}\end{array}\right|=0$
(44) (45)

with $\left|\beta_{1}\right|\leq\dots\leq\left|\beta_{4}\right|$. Here,
$X_{j}~{}(j=1,\dots,4)$ are defined as

$X_{j}=E_{\rm
OBC}-t_{L}\beta_{j}-t_{R}\beta_{j}^{-1},~{}\left(j=1,\dots,4\right).$ (46)

Finally we can obtain the boundary equation (23) from Eq. (45).

##  Appendix B Non-Bloch band theory

In this appendix, we briefly review the non-Bloch band theory [29]. First, we
give the condition for continuum energy bands in non-Hermitian systems with
open boundary conditions. Next, we show that, in the ladder model, the
continuum energy band in an open chain is equal to that in a periodic chain.

###  B.1 Condition for continuum energy bands

In a non-Hermitian system, the real-space Hamiltonian is given as

$H=\sum_{n}\sum_{i=-N}^{N}\sum_{\mu,\nu=1}^{q}t_{i,\mu\nu}c_{n+i,\mu}^{\dagger}c_{n,\nu},$
(47)

where $t_{i,\mu\nu}$ is the hopping amplitude to the $i$th nearest unit cells,
and it is not equal to $t_{-i,\nu\mu}^{\ast}$. In this Hamiltonian,
$c_{n,\mu}$ is an annihilation operator of an electron with sublattice
$\mu~{}(\mu=1,\dots,q)$ in the $n$th unit cell, and the electrons hop up to
the $N$th nearest unit cells. Here, the solutions of the real-space
eigenequation $H\ket{\psi}=E_{\rm OBC}\ket{\psi}$, where
$\ket{\psi}=\left(\dots,\psi_{1,1},\dots,\psi_{1,q},\dots,\psi_{n,1},\dots,\psi_{n,q},\dots\right)^{\rm
T}$, can be written as

$\psi_{n,\mu}=\sum_{j=1}^{2M}\left(\beta_{j}\right)^{n}\phi_{\mu}^{\left(j\right)},~{}\left(\mu=1,\dots,q\right),$
(48)

where $\beta=\beta_{j}$ is the solution of the characteristic equation

$\det\left[{\cal H}\left(\beta\right)-E_{\rm OBC}\right]=0,$ (49)

with the non-Bloch matrix

$\left[{\cal
H}\left(\beta\right)\right]_{\mu\nu}=\sum_{i=-N}^{N}t_{i,\mu\nu}\beta^{i},~{}\left(\mu,\nu=1,\dots,q\right).$
(50)

We note that the characteristic equation is an algebraic equation for $\beta$
with an even degree $2M=2qN$ in general. Furthermore, we number the solutions
of the characteristic equation so as to satisfy

$\left|\beta_{1}\right|\leq\dots\leq\left|\beta_{2M}\right|.$ (51)

Then the continuum energy bands in an infinite open chain can be obtained by
the condition to the solutions of the characteristic equation, given by [29]

$\left|\beta_{M}\right|=\left|\beta_{M+1}\right|.$ (52)

This is the main conclusion of the non-Bloch band theory. The trajectories of
$\beta_{M}$ and $\beta_{M+1}$ give the generalized Brillouin zone, from which
we can obtain the continuum energy bands in the limit of a large system size.

###  B.2 Ladder model

The ladder model in Eq. (1) with $\Delta\neq 0$ is the case of $q=2$ and
$N=1$, i.e., $M=2$ in Eq. (52). Then because of Eq. (25), the condition for
continuum energy bands can be explicitly written as

$\left|\beta_{2}\right|=\left|\beta_{3}\right|=1.$ (53)

Therefore, the generalized Brillouin zone is a unit circle, and the energy
spectrum is described by the real Bloch wave number, i.e., $\beta={\rm
e}^{ik},~{}k\in{\mathbb{R}}$. By substituting $\beta={\rm e}^{ik}$ into Eq.
(22), we can get the energy spectrum in an open chain in the form of Eq. (5).
Therefore, in the ladder model, the energy spectrum in an open chain matches
with that in a periodic chain in the limit of a large system size.

##  Appendix C Derivation of Eq. (26)

In this appendix, we derive Eq. (26) from the characteristic equation (22) and
the boundary equation (23). In the following, we focus on the case with the
sufficiently large system size. In this case, the factors
$\left(\beta_{2}\beta_{4}\right)^{L+1}$ and
$\left(\beta_{3}\beta_{4}\right)^{L+1}$ are dominant on the left-hand side of
Eq. (23). Hence, we can approximate Eq. (23) by

$\left(\beta_{2}\right)^{2L+2}\simeq\frac{X_{1,2}X_{3,4}}{X_{1,3}X_{2,4}}$
(54)

by using Eq. (25).

Next, we give perturbative solutions of Eq. (22) up to the second order of the
coupling parameter $\Delta$. They can be explicitly written as

$\displaystyle\left\\{\begin{array}[]{l}\beta_{1}\simeq
x_{-}^{\left(1\right)}+y_{-}^{\left(1\right)}\Delta^{2},\vspace{5pt}\\\
\beta_{2}\simeq
x_{+}^{\left(1\right)}+y_{+}^{\left(1\right)}\Delta^{2},\vspace{5pt}\\\
\beta_{3}\simeq
x_{-}^{\left(2\right)}+y_{-}^{\left(2\right)}\Delta^{2},\vspace{5pt}\\\
\beta_{4}\simeq
x_{+}^{\left(2\right)}+y_{+}^{\left(2\right)}\Delta^{2},\end{array}\right.$
(59)

where

$\displaystyle\left\\{\begin{array}[]{l}\displaystyle
x_{\pm}^{\left(1\right)}=\frac{1}{2t_{L}}\left(E_{\rm
OBC}\pm\delta\right),\vspace{5pt}\\\ \displaystyle
x_{\pm}^{\left(2\right)}=\frac{1}{2t_{R}}\left(E_{\rm
OBC}\pm\delta\right),\vspace{5pt}\\\ \displaystyle
y_{\pm}^{\left(1\right)}=\frac{\pm\left(E_{\rm OBC}^{2}-2t_{L}t_{R}\pm\delta
E_{\rm
OBC}\right)}{\left(t_{L}-t_{R}\right)\delta\left[2t_{L}\left(t_{L}+t_{R}\right)-E_{\rm
OBC}^{2}\mp\delta E_{\rm OBC}\right]},\vspace{5pt}\\\ \displaystyle
y_{\pm}^{\left(2\right)}=\frac{\pm\left(E_{\rm OBC}^{2}-2t_{L}t_{R}\pm\delta
E_{\rm
OBC}\right)}{\left(t_{L}-t_{R}\right)\delta\left[-2t_{R}\left(t_{L}+t_{R}\right)+E_{\rm
OBC}^{2}\pm\delta E_{\rm OBC}\right]},\end{array}\right.$ (64) (65)

and

$\delta=\sqrt{E_{\rm OBC}^{2}-4t_{L}t_{R}}.$ (66)

From Eq. (59), we can calculate the right-hand side of Eq. (54). In fact, the
factors included in the right-hand side of this equation can be given by

$\displaystyle\left\\{\begin{array}[]{l}\displaystyle
X_{1,2}=\frac{t_{L}+t_{R}}{t_{L}-t_{R}}\frac{\delta}{\left(t_{L}+t_{R}\right)^{2}-E_{\rm
OBC}^{2}}\Delta^{2}+{\cal O}\left(\Delta^{4}\right),\vspace{5pt}\\\
\displaystyle
X_{3,4}=\left(\frac{t_{L}}{t_{R}}-\frac{t_{R}}{t_{L}}\right)\delta+{\cal
O}\left(\Delta^{2}\right),\vspace{5pt}\\\ \displaystyle
X_{1,3}=\left(t_{L}-t_{R}\right)\left(\frac{E_{\rm
OBC}-\delta}{2t_{R}}-\frac{2t_{R}}{E_{\rm OBC}-\delta}\right)+{\cal
O}\left(\Delta^{2}\right),\vspace{5pt}\\\ \displaystyle
X_{2,4}=\left(t_{L}-t_{R}\right)\left(\frac{E_{\rm
OBC}+\delta}{2t_{R}}-\frac{2t_{R}}{E_{\rm OBC}+\delta}\right)+{\cal
O}\left(\Delta^{2}\right).\end{array}\right.$ (71) (72)

As a result, we can get Eq. (26) by substituting Eq. (72) into Eq. (54) and
approximating the equation up to the second order of $\Delta$.

##  Appendix D Derivation of Eq. (34)

In this appendix, we briefly explain the derivation of Eq. (34). In the
following, we focus on a non-Hermitian system in Eq. (47) with an open
boundary condition with the system size $L$. In the following, we assume that
the system has some symmetries which reverse the direction of the chain, such
as inversion symmetry. We note that, in this case, the characteristic equation
(49) is a reciprocal equation, and its solutions satisfy
$\beta_{j}=1/\beta_{2M+1-j}~{}(j=1,\dots,M)$. Now we have $2M$ boundary
conditions to the eigenstates in Eq. (48). Although Eq. (48) has the $2qM$
unknown combination coefficients, from the real-space eigenequation
$H\ket{\psi}=E_{\rm OBC}\ket{\psi}$, we can reduce the $2qM$ coefficients to
the $2M$ coefficients, e.g., $\phi_{1}^{\left(j\right)}~{}(j=1,\dots,2M)$. As
a result, the $2M$ boundary conditions to the $2M$ coefficients can be
obtained from the left and right ends of an open chain. Finally, we have the
condition so that the coefficients
$\phi_{1}^{\left(j\right)}~{}(j=1,\dots,2M)$ take nonzero values, written as

$\sum_{P,Q}F\left(\beta_{i\in P},\beta_{j\in Q},E_{\rm OBC}\right)\prod_{k\in
P}\left(\beta_{k}\right)^{L}=0,$ (73)

where the sets $P$ and $Q$ are two disjoint subsets of the set
$\\{1,\dots,2M\\}$ with $M$ elements, and the sets of the solutions of Eq.
(49) determine the values of the function $F$. In the limit of a large system
size, we can reduce Eq. (73) to

$\left(\frac{\beta_{M}}{\beta_{M+1}}\right)^{L}=-\frac{F\left(\beta_{i\in
P_{0}},\beta_{j\in Q_{0}},E_{\rm OBC}\right)}{F\left(\beta_{i\in
P_{1}},\beta_{j\in Q_{1}},E_{\rm OBC}\right)},$ (74)

with
$P_{0}=\\{M+1,\dots,2M\\},Q_{0}=\\{1,\dots,M\\},P_{1}=\\{M,M+2,\dots,2M\\}$,
and $Q_{1}=\\{1,\dots,M-1,M+1\\}$. This is because Eq. (73) has two leading
terms proportional to $\left(\beta_{M}\beta_{M+2}\cdots\beta_{2M}\right)^{L}$
and $\left(\beta_{M+1}\beta_{M+2}\cdots\beta_{2M}\right)^{L}$. Then, because
of $\beta_{M}=1/\beta_{M+1}$, we can write Eq. (74) in the form of Eq. (34).

## References

  * Guo  _et al._ [2009] A. Guo, G. J. Salamo, D. Duchesne, R. Morandotti, M. Volatier-Ravat, V. Aimez, G. A. Siviloglou, and D. N. Christodoulides, Phys. Rev. Lett. 103, 093902 (2009).
  * Feng  _et al._ [2013] L. Feng, Y.-L. Xu, W. S. Fegadolli, M.-H. Lu, J. E. Oliveira, V. R. Almeida, Y.-F. Chen, and A. Scherer, Nat. Mater. 12, 108 (2013).
  * Wang  _et al._ [2019] X. Wang, X. Fang, D. Mao, Y. Jing, and Y. Li, Phys. Rev. Lett. 123, 214302 (2019).
  * Yao and Wang [2018] S. Yao and Z. Wang, Phys. Rev. Lett. 121, 086803 (2018).
  * Yoshida  _et al._ [2020] T. Yoshida, T. Mizoguchi, and Y. Hatsugai, Phys. Rev. Research 2, 022062(R) (2020).
  * Okuma  _et al._ [2020] N. Okuma, K. Kawabata, K. Shiozaki, and M. Sato, Phys. Rev. Lett. 124, 086801 (2020).
  * Zhang  _et al._ [2020] K. Zhang, Z. Yang, and C. Fang, Phys. Rev. Lett. 125, 126402 (2020).
  * Kawabata  _et al._ [2020a] K. Kawabata, M. Sato, and K. Shiozaki, Phys. Rev. B 102, 205118 (2020a).
  * Okugawa  _et al._ [2020] R. Okugawa, R. Takahashi, and K. Yokomizo, Phys. Rev. B 102, 241202(R) (2020).
  * Fu  _et al._ [2021] Y. Fu, J. Hu, and S. Wan, Phys. Rev. B 103, 045420 (2021).
  * Okugawa  _et al._ [2021] R. Okugawa, R. Takahashi, and K. Yokomizo, Phys. Rev. B 103, 205205 (2021).
  * Brandenbourger  _et al._ [2019] M. Brandenbourger, X. Locsin, E. Lerner, and C. Coulais, Nat. Commun. 10, 4608 (2019).
  * Gou  _et al._ [2020] W. Gou, T. Chen, D. Xie, T. Xiao, T.-S. Deng, B. Gadway, W. Yi, and B. Yan, Phys. Rev. Lett. 124, 070402 (2020).
  * Xiao  _et al._ [2020] L. Xiao, T. Deng, K. Wang, G. Zhu, Z. Wang, W. Yi, and P. Xue, Nat. Phys. 16, 761 (2020).
  * Weidemann  _et al._ [2020] S. Weidemann, M. Kremer, T. Helbig, T. Hofmann, A. Stegmaier, M. Greiter, R. Thomale, and A. Szameit, Science 368, 311 (2020).
  * Hofmann  _et al._ [2020] T. Hofmann, T. Helbig, F. Schindler, N. Salgo, M. Brzezińska, M. Greiter, T. Kiessling, D. Wolf, A. Vollhardt, A. Kabaši, C. H. Lee, A. Bilušić, R. Thomale, and T. Neupert, Phys. Rev. Research 2, 023265 (2020).
  * Helbig  _et al._ [2020] T. Helbig, T. Hofmann, S. Imhof, M. Abdelghany, T. Kiessling, L. Molenkamp, C. Lee, A. Szameit, M. Greiter, and R. Thomale, Nat. Phys. 16, 747 (2020).
  * Ghatak  _et al._ [2020] A. Ghatak, M. Brandenbourger, J. van Wezel, and C. Coulais, Proc. Nat. Ac. Sc. USA 117, 29561 (2020).
  * Song  _et al._ [2019] F. Song, S. Yao, and Z. Wang, Phys. Rev. Lett. 123, 170401 (2019).
  * Longhi [2020] S. Longhi, Phys. Rev. Lett. 124, 066602 (2020).
  * Yu  _et al._ [2020a] T. Yu, Y.-X. Zhang, S. Sharma, X. Zhang, Y. M. Blanter, and G. E. W. Bauer, Phys. Rev. Lett. 124, 107202 (2020a).
  * Yu  _et al._ [2020b] T. Yu, X. Zhang, S. Sharma, Y. M. Blanter, and G. E. W. Bauer, Phys. Rev. B 101, 094414 (2020b).
  * Lee and Longhi [2020] C. H. Lee and S. Longhi, Commun. Phys. 3, 147 (2020).
  * Yi and Yang [2020] Y. Yi and Z. Yang, Phys. Rev. Lett. 125, 186802 (2020).
  * Liu  _et al._ [2020] C.-H. Liu, K. Zhang, Z. Yang, and S. Chen, Phys. Rev. Research 2, 043167 (2020).
  * Flebus  _et al._ [2020] B. Flebus, R. A. Duine, and H. M. Hurst, Phys. Rev. B 102, 180408(R) (2020).
  * Xue  _et al._ [2021] W.-T. Xue, M.-R. Li, Y.-M. Hu, F. Song, and Z. Wang, Phys. Rev. B 103, L241408 (2021).
  * Li  _et al._ [2020] L. Li, C. H. Lee, S. Mu, and J. Gong, Nat. Commun. 11, 5491 (2020).
  * Yokomizo and Murakami [2019] K. Yokomizo and S. Murakami, Phys. Rev. Lett. 123, 066404 (2019).
  * Kawabata  _et al._ [2020b] K. Kawabata, N. Okuma, and M. Sato, Phys. Rev. B 101, 195147 (2020b).
  * Yokomizo and Murakami [2020a] K. Yokomizo and S. Murakami, Prog. Theor. Exp. Phys. 2020, 12A102 (2020a).
  * Yokomizo and Murakami [2020b] K. Yokomizo and S. Murakami, Phys. Rev. Research 2, 043045 (2020b).
  * Yang  _et al._ [2020] Z. Yang, K. Zhang, C. Fang, and J. Hu, Phys. Rev. Lett. 125, 226402 (2020).
  * Yokomizo and Murakami [2021] K. Yokomizo and S. Murakami, Phys. Rev. B 103, 165123 (2021).
  * Okuma and Sato [2019] N. Okuma and M. Sato, Phys. Rev. Lett. 123, 097701 (2019).
  * Hatano and Nelson [1996] N. Hatano and D. R. Nelson, Phys. Rev. Lett. 77, 570 (1996).
  * Li  _et al._ [2021] L. Li, C. H. Lee, and J. Gong, Commun. Phys. 4, 42 (2021).


# Optimal Fidelity-Aware Entanglement Distribution
in Linear Quantum Networks

Iordanis Koutsopoulos
Department of Informatics
Athens University of Economics and Business
Athens, Greece

###### Abstract

We study the problem of entanglement distribution in terms of maximizing a
utility function that captures the total fidelity of end-to-end entanglements
in a two-link linear quantum network with a source, a repeater, and a
destination. The nodes have several quantum memories, and the problem is how
to coordinate entanglement purification in each of the links, and entanglement
swapping across links, so as to achieve the goal above. We show that
entanglement swapping (i.e, deciding on the pair of qubits from each link to
perform swapping on) is equivalent to finding a max-weight matching on a
bipartite graph. Further, entanglement purification (i.e, deciding which pairs
of qubits in a link will undergo purification) is equivalent to finding a max-
weight matching on a non-bipartite graph. We propose two polynomial
algorithms, the Purify-then-Swap (PtS) and the Swap-then-Purify (StP) ones,
where the decisions about purification and swapping are taken with different
order. Numerical results show that PtS performs better than StP, and also that
the omission of purification in StP gives substantial benefits.

##  I Introduction

The vision of the quantum internet [1] is currently underway, and in the next
years we will witness ground-breaking changes on a variety of areas applicable
in tactical networks, such as secure information exchange using Quantum Key
Distribution (QKD), improved sensing, and quantum processor networking for
quantum computations. The quantum internet can be realized through terrestrial
optical-fiber or wireless links, or satellite-to-ground-station links [2]. The
latter seem a good option, due to low attenuation of quantum signals through
the atmosphere, and long communication ranges.

Regardless of the communication link, quantum communication is achieved
through pairs of quantum bits (qubits), namely pairs of entangled photons,
i.e, photons with a coupled state. Photon pairs are generated by a photon
source and are distributed to two nodes at distinct locations through wireless
or optical-fiber links. By measuring the information carried by one photon, it
is possible to obtain information on the entangled one, based on the Einstein-
Podolsky-Rosen (EPR) paradox of quantum mechanics. If two neighboring nodes,
say Alice (A) and Bob (B) share a pair of entangled qubits, A can secretly
send one information bit to B with the help of a quantum measurement, using
the classical internet channel. This process is called teleportation. Secret
communication is achieved due to the laws of quantum mechanics for quantum
entanglement that make it impossible to eavesdrop or measure the communicated
bit. By repeating the process above, A and B can exchange a secret key with
many information bits.

Figure 1: (a) Quantum Teleportation, with which one bit of information is
transmitted from A to B, if A,B share an entangled pair of qubits; (b)
Entanglement swapping that leads to end-to-end entanglement between A and B,
through BSM measurement at the qubits in the memories of repeater R. After
end-to-end entanglement is established, A can send an information bit to B
through teleportation.

Since transmitted signals over optical fibers or over the air attenuate with
distance, the probability of successful entanglement between two nodes
decreases at higher distance. Thus, repeater nodes are employed to maintain
high probability of successful entanglement, while extending the distance
range for the established entanglement. If nodes A and B need to share a pair
of entangled qubits and they are not directly connected, a repeater node R
that holds an entanglement to A through physical link $(A,R)$ and an
entanglement to B through physical link $(R,B)$, can convert the two one-hop
entanglements into one direct, longer-distance entanglement between A and B
through entanglement swapping, with a so-called Bell State Measurement (BSM)
at node R. After swapping, the entangled pairs at each of the two links
$(A,R)$ and $(R,B)$ do not exist, and an entangled pair between nodes A and B
is created. Thus, a quantum repeater facilitates the transfer of quantum
states between nodes attached to it through entanglement swapping. Swapping
can be repeated in more intermediate repeaters to build entanglements over
long distances. Teleportation and swapping are shown in Fig. 1.

Entangled qubits reside in quantum memories in repeaters. The quality of an
entangled qubit pair is measured through its fidelity, which captures at a
given time the purity of quantum state, namely the similarity of the state of
the entangled pair to its state when it was created. Stored qubits in quantum
memories have limited lifetime since they are subject to decoherence, which
means that fidelity decreases with time, and after some time, they are
rendered non-functional for quantum applications. Finally, two or more pairs
of entangled qubits can be sacrificed to create one entangled pair of higher
fidelity, and this process is known as entanglement purification.

Entangled pairs of qubits are the most important resources in the quantum
internet and are necessary for quantum protocols and procedures. A prime goal
is to efficiently distribute a large number of high-fidelity entanglements
across the network, from source to destination. In a network in which nodes
have several quantum memories, the main question is how to tune entanglement
purification of stored pairs of qubits in quantum memories of the two links,
and entanglement swapping across the two links so as to achieve the goal
above.

We study the problem of efficient entanglement distribution in terms of
maximizing a utility function that captures total fidelity of end-to-end
entanglements in a two-link linear quantum network, where each node has
multiple stored qubits in its quantum memories. At each time slot, a
controller observes the system state, that is, the set of entangled pairs of
qubits in each link, and their fidelity. We ask the question: given a
continuous backlog of requests for end-to-end entanglement, how can we
increase the derived utility through appropriate management of stored qubits
in each link with entanglement swapping and purification? Our two-hop model
may represent a terrestrial network with optical-fiber or wireless links, or
one with satellite-to-ground links, and it is the simplest non-trivial model
on which we show the joint impact of entanglement swapping and purification on
entanglement distribution without involving entanglement routing. Also, this
model stands for entanglement distribution over a general network topology and
several source-destination pairs. If the route for each source-destination
pair is given and routes do not intersect, then each route corresponds to the
repeater chain model of this paper. The contributions of our work are as
follows:

  * •

We study the problem of increasing utility of end-to-end entanglement through
purification and swapping operations in the stored qubits of quantum memories.
We take into account essential assumptions such as decoherence of stored
qubits with time, uncertain swapping and purification operations, and the
impact of swapping and purification on the fidelity of entangled pairs.

  * •

We show that entanglement swapping in the two-link network, i.e. deciding
which pairs of qubits from each link to swap, is equivalent to finding a max-
weight matching on a bipartite graph. Further, entanglement purification, i.e,
deciding which pairs of qubits of a link will undergo purification, is
equivalent to finding a max-weight matching on a non-bipartite graph.

  * •

We propose two polynomial algorithms, the Purify-then-Swap (PtS) and the Swap-
then-Purify (StP) ones, where the decisions about entanglement purification
and swapping are taken with different order. Numerical results show that PtS
performs better than StP, and that the omission of purification in StP gives
substantial benefits.

Our work is the first one that studies the problem of entanglement
distribution through nodes with multiple quantum memories, through the
interplay of entanglement purification and swapping. The paper is organized as
follows. In section II, we provide a literature overview, while in section III
we present the model and problem statement, and in section IV, we present two
algorithms. Sections V and VI include numerical results and the conclusion
respectively.

##  II Related work

Quantum entanglement routing. Quantum entanglement routing studies the
distribution of entanglement over multiple hops. The most important metric is
entanglement generation rate which measures the number of entangled pairs
shared between source and destination per unit of time. In [3], the authors
propose entanglement distribution protocols under realistic assumptions on
entanglement swapping, number and quality of quantum memories, link losses,
and local link state information (e.g. in terms of entanglement success or
failure).

A survey of entanglement routing protocols is [4], with a view towards
multiple commodities, where a commodity is a source-destination pair with an
entanglement request. The work [5] is among the first that builds a quantum
entanglement routing protocol by taking into account different topologies,
link state exchanges, qubit capacity limitations, and commodities that compete
for resources. The work [6] also studies a multi-commodity flow optimization
problem to find the set of paths for each commodity such that the total end-
to-end entanglement generation rate is maximized, subject to a minimum end-to-
end fidelity for each commodity, reflected as a maximum path length
constraint. In [7], the authors propose entanglement routing algorithms with
purification decision to provide fidelity guarantees for many commodities.

In [8], the authors use Markov Decision Processes to find the policy that
minimizes delivery time of end-to-end entanglement in a repeater chain, where
repeater memories have cut-off times due to decoherence. The optimal policy
decides when to perform entanglement swapping in each repeater after observing
the state, which is the age of each entangled qubit pair. The work [9] also
considers cut-off times. An entangled pair is kept for some time in the
memory, and then it is discarded so that a new generated pair with higher
fidelity is entered. The authors optimize cut-offs to maximize end-to-end
entanglement generation rate. In [10], the authors design a quantum data plane
protocol for end-to-end communication and provide mechanisms for coordinating
entanglement swapping, handling decoherence of stored qubits, and ensuring
high end-to-end fidelity through distillation. They also introduce the
NetSquid simulator for quantum networks.

The work [11] defines three logarithmic utility functions that involve
entanglement generation rate and fidelity. Two of them respectively are
interpreted as the utility from entanglement rate and distillation, the latter
being a function of fidelity, and as the utility from the amount of bits in
the secret key, which also involves the entanglement rate and entropy of
Werner states. The third utility function is concave in entanglement rate and
fidelity, thus mathematically tractable.

Entanglement purification and routing. The works [12] and [13] study the
integration of purification protocols into entanglement distribution with
realistic system constraints such as imperfect fidelities and limited storage
time. The works highlight the importance of tuning purification protocols for
better end-to-end entanglement rate. In [14], the authors show that under
perfect swapping and purification operations, the per-memory entanglement
distribution rate increases as the number of memories increases even under
memory decoherence; however this is not the case for imperfect operations.

Quantum switches. The authors in [15] study a quantum switch that serves
through its memories incoming bipartite end-to-end entanglement requests made
by users users in a star topology. The authors use Continuous-Time Markov
Chains to characterize the switch capacity, namely the maximum possible number
of end-to-end entanglements served by the switch per time unit, as well as the
probability distribution and expected number of qubits in the memories.
Finally, the authors in [16] propose a policy for allocating the limited
quantum memories of a switch to link-level entanglements of clients under
dynamic arrivals of entanglement requests. They characterize the switch
capacity region, namely the set of request arrival vectors for which the
queues of serving requests are stable and propose a memory allocation policy
that is based on the max-weight algorithm [17] and is throughput optimal.

##  III System Model and Problem Statement

###  III-A System model

We consider a two-hop linear quantum network with a source node $s$, a
repeater $r$, and a destination node $d$, and two physical links, $(s,r)$ and
$(r,d)$. The network may represent a terrestrial system with optical-fiber or
wireless physical links, or a hybrid terrestrial-satellite system, where nodes
$s$, $d$ are ground stations, and the repeater $r$ is a satellite. Nodes $s$
and $d$ engage in applications such as QKD, thus entangled qubit pairs need to
be generated and distributed end-to-end from $s$ to $d$ so that they are used
in applications. Each of $s$ and $d$ has $M$ quantum memories, while the
repeater has $2M$ quantum memories; $M$ of those are used to attempt link-
level entanglements (LLEs) with the $M$ memories of $s$, while the other $M$
memories are used to attempt LLEs with the $M$ memories of $d$. Each memory
can store one qubit.

####  III-A1 Entanglement generation

Time is slotted, and the slot duration $\delta$ is assumed to be larger than
the time it takes for purification and swapping together to take place, and it
can be of the order of tens of msecs. At each time slot, in link $(s,r)$,
link-level entanglement (LLEs) are attempted as follows. Node $s$ has a photon
generator that generates an entangled qubit (photon) pair. One photon is kept
locally, and the other is sent to $r$ via an optical-fiber or free-space link.
The probability of successful photon transmission through link $(s,r)$,
$p_{sr}$, decays exponentially with distance due to photon loss.

If the entangled qubit is received correctly, $r$ sends an acknowledgment
message to $s$ through a classical communication channel, while one qubit of
the entangled qubit pair is stored in a memory of $s$ and the other qubit in a
memory of $r$. The entanglement generation process is repeated for $M$ times
in parallel, and each successfully generated entangled qubit is stored at a
memory $i$ of $s$, and the corresponding memory $i$ of $r$, for
$i=1,\ldots,M$. Similarly, $M$ LLEs are attempted on link $(r,d)$ between
memory $M+i$ of $r$ and memory $i$ of $d$, $i=1,\ldots,M$, with success
probability is $p_{rd}$.

At each slot, the entanglement generation attempts in link $(s,r)$ are
continual as long as there are unoccupied pairs of memories at node $s$ and
$r$, and the indices of successful EPR pairs between $s$ and $r$ are known.
The same holds for link $(r,d)$. If a pair of entangled qubits are shared by
$s$ and $r$, then $s$ can send one bit of secret information to $r$ with
teleportation.

####  III-A2 Fidelity and decoherence

The quality of an LLE is measured by fidelity $F$, which relates to the
coherence of the two entangled qubits. When an LLE is generated, it has
initial fidelity $F_{0}$. The quality of the LLE decays with time due to
interactions of stored qubits with the environment, and this is known as
decoherence. If the LLE has been generated at time slot $t_{0}$, its fidelity
at time slot $t>t_{0}$ is [4], [8]:

$F(t)=\frac{1}{4}+\left(F_{0}-\frac{1}{4}\right)e^{-\delta(t-t_{0})/\tau}$ (1)

where $\tau$ characterizes the exponential decay. If fidelity is smaller than
a threshold, the qubit is useless for applications.

####  III-A3 Entanglement Swapping

We assume bipartite entanglements. An end-to-end entanglement between nodes
$s$ and $d$ can be created through entanglement swapping of one LLE of link
$(s,r)$ and one LLE of link $(r,d)$. The swapping is realized with BSM of the
two qubits at the corresponding memories of node $r$. The two LLEs disappear,
and after swapping, there exists an entangled qubit pair between nodes $s$ and
$d$. Entanglement swapping succeeds with probability $p_{s}$. If an end-to-end
entanglement is successful, it comes with a fidelity that is a deterministic
function of the fidelity values of the involved LLEs. We assume for simplicity
that swapping operations are perfect, that is, Bell states preserve their form
after swapping. If $F_{1},F_{2}$ are the fidelity values of the two LLEs
before swapping, the resulting end-to-end entanglement after swapping has
fidelity:

$F_{s}(F_{1},F_{2})=F_{1}F_{2}+(1-F_{1})(1-F_{2})\,.$ (2)

####  III-A4 Entanglement Purification

The fidelity of entangled qubits (either LLEs or end-to-end) decreases with
time. In order to increase fidelity, we use entanglement purification. We
assume $2$-to-$1$ purification, where two entangled qubit pairs in the same
link are sacrificed to generate a single pair over that link, with higher
fidelity than the sacrificed ones. In this case, the new purified qubit pair
substitutes one of the two existing pairs in memories, while the other qubit
pair is discarded and the memories are freed. We assume perfect purification
operation which succeeds with probability
$p_{pur}=F_{1}F_{2}+(1-F_{1})(1-F_{2})$ [18]. The fidelity of an entanglement
resulting from purification of two entangled qubit pairs with fidelities
$F_{1}$, $F_{2}$ is [18]

$F_{p}(F_{1},F_{2})=\frac{F_{1}F_{2}}{F_{1}F_{2}+(1-F_{1})(1-F_{2})}\,.$ (3)

After purification, the link nodes receive the fidelity value. We do not
assume multiple levels of purification, namely once an LLE is formed from
purification, it is not purified again.

####  III-A5 Sequence of events

At the end of time slot $t-1$, entanglements are attempted between unoccupied
memories of nodes $s$ and $r$, and between memories of $r$ and $d$ and succeed
with probability $p_{s}$. At the beginning of time slot $t$, a controller at
node $s$ observes the system state, namely the successful LLEs of links
$(s,r)$ and $(r,d)$, and their fidelity values. Then, it decides on the order
of purification and swapping, and on how they will be performed. Swapping and
purification operations are realized during slot $t$, and $L$ attempts for
purification may be made. After swapping and after purification, the resulting
entanglement obtains the age of the oldest LLE involved in the operation. In
the end of slot $t$, entanglements are again attempted between unoccupied
memories of $s$, $r$ and those of $r$, $d$, and the process continues. At the
end of each slot, the fidelity of each entangled qubit decays according to
(1). When a number of time slots elapses and fidelity drops below a threshold,
the qubit pair is discarded, and the corresponding memories are emptied.

###  III-B Problem statement

Since nodes have $M$ quantum memories, we can create at most $M$ entanglements
from $s$ to $d$. For quantum applications, it is important to have a large
number of entangled qubits so that we securely transmit more bits. On the
other hand, the fidelity of entangled qubits needs to be high, so that they
can better support applications, and for longer time before they decohere. Let
$\mathbf{x}$ be the $M\times 1$ binary vector of end-to-end entanglements and
$\mathbf{F}$ the vector of fidelity values, where $x_{i}=1$ means that there
exists an end-to-end entanglement that originates from the $i$-th memory pair
of $(s,r)$, and $0$ if not, and $F_{i}$ is the fidelity of this entanglement.
The importance of the number and fidelity values of entanglements is captured
by a utility function

$u(\mathbf{x},\mathbf{F})=\displaystyle\log(\sum_{i=1}^{M}x_{i}g(F_{i}))\,$
(4)

Figure 2: Upper part (a): Equivalence of the optimal purification problem to
max-weight matching in a non-bipartite graph. We show a hypothetical best
solution: LLEs b and c are purified, while LLE a is not. This corresponds to a
max-weight matching in the non-bipartite graph on the right. Bottom part (b):
Equivalence of the optimal swapping problem to max-weight matching in a
bipartite graph. We show a hypothetical best solution: LLEs a,b are swapped,
and so are LLEs b,h, and LLEs c,e. This corresponds to the max-weight matching
in the bipartite graph on the right.

where function $g(\cdot)$ depends on the application [11] and denotes utility
out of fidelity. On the one hand, it is important to perform the swapping of
qubit pairs from links $(s,r)$ and $(r,d)$ so as to have high total fidelity.
On the other hand, purification leads to a smaller number of entangled pairs,
albeit of higher fidelity. In order to achieve good utility, we need to
balance the trade-off above by managing the swapping and purification
operations. It is not clear whether we should strive for a small number of
entangled qubit pairs of large fidelity or for a higher number of qubit pairs
of lower fidelity. The order of swapping and purification plays an important
role.

##  IV Entanglement Swapping and Purification

We propose two algorithms for high utility of end-to-end entanglements. The
algorithms differ in the order they perform swapping and purification. First,
we show how to do each of these operations efficiently.

###  IV-A Optimal Entanglement Swapping

First, consider entanglement swapping between LLEs of link $(s,r)$ and LLEs of
link $(r,d)$ in order to form end-to-end entanglements with high utility. At
each time slot, the controller sees the successful entangled qubit pairs and
their fidelity values in each of the two links. We define a bipartite graph
$\mathcal{G}=(\mathcal{U}\cup\mathcal{V},\mathcal{E})$. Each node
$i\in\mathcal{U}$ stands for an LLE in link $(s,r)$, namely a pair of memories
$(i,i)$ with a successful entangled qubit pair (one memory in node $s$ and
another in $r$). Each node $j\in\mathcal{V}$ stands for an LLE in link
$(r,d)$, namely a pair of memories $(M+i,i)$ with a successful entangled qubit
pair (one memory in node $r$ and the another in $d$). Let $F_{i},F_{j}$ be the
fidelity values of those LLEs. We define the weight of link $(i,j)$ as
$w_{ij}=g(F_{s}(F_{i},F_{j}))$, which is the utility of fidelity
$F_{s}(F_{i},F_{j})$ after swapping between LLEs $i$ and $j$.

The problem of finding the LLE swapping operation that maximizes total utility
of end-to-end entanglements is equivalent to finding a max-weight matching
$\mathcal{M}^{*}$ in bipartite graph $\mathcal{G}$. The edges of the matching
denote the LLEs of links $(s,r)$ and $r,d)$ that will be swapped (Fig. 2(b)).
The max-weight matching can be found with the Hungarian algorithm [19].

###  IV-B Optimal Entanglement Purification

Next, we consider $2$-to $1$ purification of LLEs in link $(s,r)$ (and
similarly for $(r,d)$), namely finding the pairs of LLEs to be merged so that
the resulting LLEs have high utility. Focus on link $(s,r)$. We define the
non-bipartite graph $\mathcal{G}_{sr}=(\mathcal{V}_{sr},\mathcal{E}_{sr})$.
The node set $\mathcal{V}_{sr}$ consists of one node $i$ for each successful
LLE of $(s,r)$, and one replica node $i^{\prime}$ of each node $i$. The set of
edges $\mathcal{E}_{sr}$ is formed as follows. We connect nodes $i$ and $j$,
$i\neq j$, that stand for successful LLEs. We assign to link $(i,j)$ weight
$w_{ij}=F_{p}(F_{i},F_{j})$, which is the utility after purification of LLEs
corresponding to $i$ and $j$, with initial fidelity $F_{i}$ and $F_{j}$. We
also connect each node $i$ only with its own replica node $i^{\prime}$, and we
assign to link $(i,i^{\prime})$ weight $w_{ii^{\prime}}=F_{i}$. This link
denotes that the LLE corresponding to $i$ will not be purified, thus its
fidelity remains the same. Thus, links $(i,j)$ with $i\neq j$ stand for the
tentative pairing of LLEs $i$ and $j$ for purification, while links
$(i,i^{\prime})$ stand for tentative non-pairing of LLE $i$ with another LLE.

The problem of purification on link $(s,r)$ to maximize utility is equivalent
to finding a max-weight matching $\mathcal{M}^{*}_{sr}$ on $\mathcal{G}_{sr}$.
The edges in $\mathcal{M}^{*}_{sr}$ dictate the way to perform the
purification: edges $(i,j)$, $i\neq j$, in $\mathcal{M}^{*}_{sr}$ mean that
LLE $i$ is merged with LLE $j$ and purified; edges $(i,i^{\prime})$ in
$\mathcal{M}^{*}_{sr}$ mean that LLE $i$ is not merged with another LLE and
maintains its fidelity (Fig. 2(a)). The max-weight matching can be found in
polynomial time with the Edmonds blossom algorithm [20].

###  IV-C Purify-then-Swap algorithm

The Purify-then-Swap (PtS) algorithm performs first purification and then
swapping. For the purification of LLEs on link $(s,r)$, we find a max-weight
matching $\mathcal{M}^{*}_{sr}$ in the non-bipartite graph $\mathcal{G}_{sr}$
defined above. Similarly, we purify LLEs on link $(r,d)$ by finding a max-
weight matching $\mathcal{M}^{*}_{rd}$ in the corresponding graph
$\mathcal{G}_{rd}$. The set of edges in each matching denotes for each link
the resulting purified LLEs after merging, and the LLEs that will not undergo
purification.

Next, we define bipartite graph
$\mathcal{G}_{s}=\mathcal{U}_{s}\cup\mathcal{V}_{s},\mathcal{E}_{s}$ as
follows, for the swapping. For each link in $\mathcal{M}^{*}_{sr}$, we create
a node $k\in\mathcal{U}_{s}$, and for each link in $\mathcal{M}^{*}_{rd}$, we
create a node $\ell\in\mathcal{V}_{s}$. We assign a weight to link $(k,\ell)$
as follows. If node $k$ corresponds to a matched edge of the form $(i,j)$,
$i\neq j$ for link $(s,r)$, and node $\ell$ corresponds to a matched link of
the form $(i_{0},j_{0})$, $i_{0}\neq j_{0}$ for link $(r,d)$, the weight is
$w_{k\ell}=g(F_{s}(F_{p}(F_{i},F_{j}),F_{p}(F_{i_{0}},F_{j_{0}})))$, which
corresponds to swapping of LLEs $k$ and $\ell$ that have emerged through
purification. If node $k$ corresponds to a matched edge of the form
$(i,i^{\prime})$, for link $(s,r)$, and node $\ell$ corresponds to a matched
edge of the form $(i_{0},i^{\prime}_{0})$ for link $(r,d)$, the weight is
$w_{k\ell}=g(F_{s}(F_{i},F_{i_{0}}))$. The weights for the other two cases of
links are assigned similarly. The LLEs that will be swapped result from a max-
weight matching on $\mathcal{G}_{s}$. An example of the PtS algorithm is shown
in Fig. 3.

Figure 3: Example of the Purify-then-Swap algorithm. Purification: In link
$(s,r)$, LLEs a and b are purified, while c is not. The new purified LLE is
stored in memory 1 of s and 1 of r. In link $(r,d)$, LLEs e and h are
purified, while q is not. The new purified LLE is stored in memory 4 or r and
2 of d. Swapping:LLEs a and q are swapped, and also LLEs c and e are swapped.
The outcome is two end-to-end entanglements.

###  IV-D Swap-then-Purify algorithm

In the Swap-then-Purify (StP) algorithm, swapping is performed first by
finding a max-weight matching $\mathcal{M}^{*}$ on bipartite graph
$\mathcal{G}$ described in subsection IV-A. The edges in the matching denote
the pairs of LLEs, one from link $(s,r)$ and one from link $(r,d)$ that will
be swapped. Let $\mathcal{N}$ be the set of end-to-end entanglements created
after swapping. We index these entanglements as $n=1,\ldots,|\mathcal{N}|$ and
denote by $F_{n}$ the fidelity of the $n$-th end-to-end entanglement after
swapping.

Next, we form the non-bipartite graph $\mathcal{G}_{p}$ with one node $n$ for
each end-to-end entanglement, and one replica node $n^{\prime}$ for each $n$.
We connect node $n$ with each other node $m\neq n$, and we connect $n$ with
$n^{\prime}$. To each edge $(n,m)$, $m\neq n$, we assign weight
$w_{mn}=g(F_{p}(F_{n},F_{m}))$, which is the fidelity after purification of
the $n$-th and $m$-th end-to-end entanglements. To each edge $(n,n^{\prime})$,
we assign weight $w_{nn^{\prime}}=g(F_{n})$, as this denotes entanglements
that are not purified. The max-weight matching $\mathcal{M}^{*}_{p}$ on
$\mathcal{G}_{p}$ denotes the way that pairs of end-to-end entanglements will
be purified or not.

##  V Numerical results

We consider a two-hop linear network $s\rightarrow r\rightarrow d$, where $s$,
$d$ have $M=3$ memories each, and $r$ has $6$ memories. We define two utility
functions for the utility of the $n\leq 3$ end-to-end entanglements [11]:

  * •

$u_{A}(\mathbf{F})=\log(\sum_{i=1}^{n}F_{i})$. The utility is the logarithm of
the total fidelity of entanglements between $s$ and $d$.

  * •

$u_{B}(\mathbf{F})=\log(\sum_{i=1}^{n}D(F_{i}))$, where $D(\cdot)$ measures
the rate at which perfect Bell states are distilled from bipartite states and
is given by

$D(x)=1+x\log_{2}x+(1-x)\log_{2}\frac{1-x}{3}\,.$ (5)

As discussed in [11], $u_{A}(\cdot)$ and $u_{B}(\cdot)$ are appropriate for
applications with emphasis on many entanglements and high fidelity
respectively. We compare the following algorithms: (i Purify-then-Swap (PtS);
(ii) Swap-then-Purify (StP); (iii) Swap-only, with no purification after
swapping.

We carried 500 experiments; in each experiment, we created 6 LLES: 3 on link
$(s,r)$ and 3 on link $(r,d)$, with fidelity values uniformly distributed in
$[0.8,1]$. Figures 4 and 5 show the comparative results of the 3 algorithms
for utility functions $u_{A}(\cdot)$ and $u_{B}(\cdot)$ respectively.
Interestingly, the Swap-only approach achieves better perfomance than PtS and
StP. Also, for $u_{A}(\cdot)$, PtS achieves $4-5\%$ better performance than
StP, suggesting that it is better to do purification first. For
$u_{B}(\cdot)$, PtS is $20\%$ better than StP, and Swap-only is again best of
all. The additive nature of utility functions favors more entanglements, and
thus it favors PtS over StP. The superiority of Swap-only is again due to the
additive utility functions.

##  VI Conclusion

We studied the interplay of swapping and purification towards achieving high
utility of end-to-end entanglements in a two-hop linear network. First, we
showed that the optimal swapping and purification problems can be viewed as
max-weight matching problems in bipartite and non-bipartite graphs
respectively. Then, we proposed two algorithms that execute these operations
in different order. This simple setting offers some first insights which are
seemingly due to the nature of the utility function that is additive in the
fidelity values of entanglements. These insights are the superiority of Swap-
only poilcy, and the better performance of PtS over StP.

The methods can be extended to larger repeater chains and can be viewed
jointly with routing i.e, determination of the repeater path that forms the
end-to-end entanglement. Further, new utility functions should quantify the
priority between the number of entanglements and their fidelity, depending on
the application. For example, if the application requires handling large
amounts of data or can tolerate higher error rates with good error correction
methods, more entanglements of lower fidelity might make more sense, while in
other cases priority could be placed on higher fidelity. Imperfections in
purification and swapping can also be modeled, as well as generalizations of
purification so as to substitute $m>2$ LLEs with $n<m$ entanglements of higher
fidelity, where $n$ is to be decided so that we free different number of
memories. Another question whether deciding the purification and swapping
jointly rather than sequentially is better. Finally, in the general case, with
dynamic requests for end-to-end entanglements between different source-
destination pairs, the problem is to manage repeater memories, purify LLEs of
links, find the route and form end-to-end entanglements so as to maximize the
utility, while keeping the queue of requests stable.

## References

  * [1] A. Singh, K. Dev, H. Siljak, H. D. Joshi and M. Magarini, “Quantum Internet-Applications, Functionalities, Enabling Technologies, Challenges, and Research Directions,” IEEE Commun. Surveys & Tut., vol. 23, no. 4, pp. 2218-2247, 4th quarter 2021.
  * [2] S. Khatri, A.J.Brady, R. A.Desporte, M.P. Bart, and J.P.Dowling, “Spooky action at a global distance: analysis of space-based entanglement distribution for the quantum internet”, Quantum inform., vol.7, no.4, Jan. 2021.
  * [3] M. Pant, H. Krovi, D. Towsley, L. Tassiulas, L. Jiang, P. Basu, D. Englund, and S. Guha, “Routing entanglement for the quantum internet”, Quantum inform., vol.5, no.25, 2019.
  * [4] F. Dupuy, C. Goursaud, and F. Guillemin, “A Survey of Quantum Entanglement Routing Protocols-Challenges for Wide-Area Networks”, Wiley Advanced Quantum Tech., vol.6, no.5, 2023.
  * [5] S. Shi and C. Qian, “Concurrent Entanglement Routing for Quantum Networks: Model and Designs”, in Proc. ACM SIGCOMM, 2020\.

Figure 4: Total utility of end-to-end entanglements vs. number of purification
attempts for utility function $\log(\sum_{i}F_{i})$. Figure 5: Total utility
of end-to-end entanglements vs. number of purification attempts for utility
function $\log(\sum_{i}D(F_{i}))$.

* [6] K. Chakraborty, D. Elkouss, B. Rijsman, and S. Wehner, “Entanglement Distribution in a Quantum Network: A Multicommodity Flow-Based Approach”, IEEE Trans. on Quantum Eng., vol.1, pp.1-21, Nov. 2020.
* [7] J. Li, M. Wang, K. Xue, R. Li, N. Yu, Q. Sun, and J. Lu, “Fidelity-Guaranteed Entanglement Routing in Quantum Networks”, IEEE Trans. on Commun., vol. 70, no.1, pp.6748-6763, Oct. 2022.
* [8] A. G. Inesta, G. Vardoyan, L. Scavuzzo, and S. Wehner, “Optimal entanglement distribution policies in homogeneous repeater chains with cutoffs”, Quantum inform., vol.9, no.46, 2023.
* [9] B. Li, T. Coopmans and D. Elkouss,“Efficient Optimization of Cutoffs in Quantum Repeater Chains,” IEEE Trans. Quantum Eng., vol.2, pp. 1-15, 2021.
* [10] W. Kozlowski, A. Dahlberg, and S. Wehner, “Designing a Quantum Network Protocol”, in Proc. ACM CoNEXT, 2020\.
* [11] G. Vardoyan and S. Wehner, “Quantum network utility maximization,” in Proc. IEEE Int. Conf. on Quantum Computing and Eng. (QCE), 2023.
* [12] M. Victora, S. Tserkis, S. Krastanov, A. Sanchez de la Cerda, S. Willis, and P. Narang, “Entanglement purification on quantum networks”, Phys. Review Research, vol.5, no.3, pp.033171, Sept. 2023.
* [13] M. Victora, S. Krastanov, A. Sanchez de la Cerda, S. Willis, and P. Narang, “Purification and Entanglement Routing on Quantum Networks”, arXiv, http://www.arxiv.org/abs/2011.11644.
* [14] A. Zang, X. Chen, A. Kolar, J. Chung, M. Suchara, T. Zhong, R. Kettimuthu, “Entanglement Distribution in Quantum Repeater with Purification and Optimized Buffer Time”, arXiv, https://arxiv.org/pdf/2305.14573.
* [15] G. Vardoyan, S. Guha, P. Nain, and D. Towsley, “On the Stochastic Analysis of a Quantum Entanglement Distribution Switch”, IEEE Trans. on Quantum Eng., vol. 2, pp. 1-16, 2021.
* [16] P. Promponas, V. Valls, S. Guha, and L. Tassiulas, “Maximizing entanglement rates via efficient memory management in flexible quantum switches”, IEEE J. Select. Areas Commun., 2024, to appear.
* [17] L. Tassiulas and A. Ephremides, “Stability properties of constrained queueing systems and scheduling policies for maximum throughput in multihop radio networks,’ IEEE Trans. on Aut. Control, vol. 37, no. 12, pp. 1936–1948, 1992.
* [18] R. van Meter, Quantum Networking, Wiley, 2014.
* [19] H.W. Kuhn, “The Hungarian Method for the assignment problem”, Naval Research Log. Quarterly, vol.2, pp.83-97, 1955.
* [20] J. Edmonds, “Maximum matching and a polyhedron with 0-1 vertices” J. of Research Nat. Bureau of Standards, vol.69B, pp. 125–130, 1965.