# Stabilization of Hubbard-Thouless pumps through nonlocal fermionic repulsion

Javier Argüello-Luengo  ICFO - Institut de Ciencies Fotoniques, The Barcelona
Institute of Science and Technology, Av. Carl Friedrich Gauss 3, 08860
Castelldefels (Barcelona), Spain  Manfred J. Mark  Institut für Quantenoptik
und Quanteninformation, Österreichische Akademie der Wissenschaften,
Technikerstraße 21a, 6020 Innsbruck, Austria  Institut für Experimentalphysik,
Universität Innsbruck, Technikerstraße 25, 6020 Innsbruck, Austria  Francesca
Ferlaino  Institut für Quantenoptik und Quanteninformation, Österreichische
Akademie der Wissenschaften, Technikerstraße 21a, 6020 Innsbruck, Austria
Institut für Experimentalphysik, Universität Innsbruck, Technikerstraße 25,
6020 Innsbruck, Austria  Maciej Lewenstein  ICFO - Institut de Ciencies
Fotoniques, The Barcelona Institute of Science and Technology, Av. Carl
Friedrich Gauss 3, 08860 Castelldefels (Barcelona), Spain  ICREA, Pg. Lluís
Companys 23, 08010 Barcelona, Spain  Luca Barbiero  Institute for Condensed
Matter Physics and Complex Systems, DISAT, Politecnico di Torino, I-10129
Torino, Italy  Sergi Julià-Farré<EMAIL_ADDRESS>ICFO - Institut de
Ciencies Fotoniques, The Barcelona Institute of Science and Technology, Av.
Carl Friedrich Gauss 3, 08860 Castelldefels (Barcelona), Spain

###### Abstract

Thouless pumping represents a powerful concept to probe quantized topological
invariants in quantum systems. We explore this mechanism in a generalized
Rice-Mele Fermi-Hubbard model characterized by the presence of competing
onsite and intersite interactions. Contrary to recent experimental and
theoretical results, showing a breakdown of quantized pumping induced by the
onsite repulsion, we prove that sufficiently large intersite interactions
allow for an interaction-induced recovery of Thouless pumps. Our analysis
further reveals that the occurrence of stable topological transport at large
interactions is connected to the presence of a spontaneous bond-order-wave in
the ground-state phase diagram of the model. Finally, we discuss a concrete
experimental setup based on ultracold magnetic atoms in an optical lattice to
realize the newly introduced Thouless pump. Our results provide a new
mechanism to stabilize Thouless pumps in interacting quantum systems.

##  1 Introduction

Since the discovery of the two-dimensional integer quantum Hall effect (IQHE)
[1, 2], the classification of phases of matter through global topological
invariants has become an intense research area [3, 4]. Topological phases are
of great interest both from a fundamental viewpoint, as they are not captured
by the standard Ginzburg-Landau-Wilson theory [5, 6] of spontaneously
symmetry-broken phases and, potentially, at the technological level, as they
represent promising tools for quantum metrology [7] and quantum computation
[8]. While noninteracting topological matter can be accurately described in
terms of the topological invariant and symmetries of single-particle bands,
such classification becomes challenging in the presence of interactions [9].

In this context, Thouless pumps [10, 11] have emerged as a powerful tool to
characterize the topology of interacting one-dimensional systems through real-
time dynamics, i.e., without relying on static groundstate properties.
Specifically, these pumps represent a one-dimensional reduction of the IQHE,
as they exploit the direct correspondence between quantized particle pumping
and the topological invariant during cyclic adiabatic dynamics around a
topological singularity. This scheme provides a natural theoretical and
experimental framework to study the effect of interactions in topological
systems undergoing such adiabatic transport dynamics. From a theoretical
perspective, recent works showed that topological quantized transport can
survive, and even be induced by many-body correlations in interacting bosonic
systems [12, 13, 14, 15, 16, 17, 18, 19]. In contrast, a repulsive onsite
Hubbard $U$ interaction [20, 21] or electron-phonon coupling [22, 23] lead to
a breakdown of the quantized Thouless pump in fermionic systems. At the
experimental level, while initially conceived within the realm of solid-state
systems [10, 11], the advent of analog quantum simulators [24, 25, 26, 27,
28], offering unprecedented levels of control and tunability, has
revolutionized the study of Thouless pumps [29, 30, 31]. These advances have
enabled the first experimental detection of topological quantized transport in
photonic platforms [32, 33, 34, 31, 35], and in systems where ultracold
bosonic [36] and fermionic [37, 38, 39, 40] atoms are trapped in optical
lattices [41, 42].

A crucial advantage offered by ultracold atomic systems, in particular, is
their impressive flexibility in engineering highly tunable interactions, as
evidenced by the recent realization of longly sought interacting topological
states of matter [43, 44]. In the case of Thouless pumps, although the first
experiments focused on single-particle transport [36, 37, 38], the pioneering
experiment [39], consisting of a fermionic spin mixture of neutral atoms
trapped in an optical lattice, reported the breakdown of Thouless pumps at
large repulsive Hubbard $U$. This observation was in agreement with the
mechanism described in Ref. [20], which explained this breakdown by the
appearance of a gapless Mott insulator (MI) state, in which the transport of
fermionic pairs is suppressed, and lead to the conclusion that for a dominant
repulsive $U$ quantized charge pump is absent in the thermodynamic limit.
Remarkably, the subsequent study [21] provided an alternative interpretation
of this breakdown in terms of the splitting of the topological singularity,
opening up the possibility to observe half-quantized transport by encircling
the displaced singularities at large $U$ in small systems exhibiting a finite-
size gap. This latter exotic Thouless pump has been recently studied in the
experiment [40], which showed that under realistic conditions half-quantized
transport can be observed for single or few Thouless cycles, despite this
effect being absent in the thermodynamic limit. While the abovementioned
results allow for a deep understanding of the efficiency of Thouless pumps in
systems characterized by local fermionic interactions, configurations where
fermions interact through different and competing processes are unexplored
and, therefore, novel phenomena might appear.

Figure 1: (a) Rice-Mele extended Fermi-Hubbard model. (b) Thouless pump
encircling the singularity $(\delta,\Delta)=(0,0)$. (c) Sketch of the
adiabatic transport properties in the interacting $U-V$ plane, for fixed
$\delta_{0}=0.9J,\ \Delta_{0}=J/2,\ TJ=50$. For shorter times satisfying
$T<(2\Delta_{0})^{-1}$, Thouless pumping would become nonadiabatic in the
noninteracting limit (red region). In the dark blue region, the ground state
of the system is the BOW phase for $\Delta=\delta=0$ [45]. (d) Expected phases
during a Thouless pump. Note that $U$ ($V$) can prevent the interaction-
induced CDW (MI) order that can break the pump at $\theta=0$ ($\theta=\pi/2$).

In order to shed light on this regime, in this work we investigate fermionic
Thouless pumps in the presence of both onsite and intersite (nonlocal)
interactions. This scenario results particularly important because, on the one
hand, electrons in solid-state systems are subject to interactions with
different range and, on the other, nonlocal fermionic interactions can be
efficiently engineered in quantum simulators made of magnetic atoms [46, 47].
To this end, we study Thouless pumps [Fig. 1(b)] in a generalized Rice-Mele
model with onsite $U$ and nearest-neighbor (NN) $V$ interactions, i.e., the
Rice-Mele extended Fermi-Hubbard model sketched in Fig. 1(a) using tensor
network methods [48, 49] and exact diagonalization. As main results [see Fig.
1(c)], we find that (i) $V$ naturally allows one to preserve the adiabaticity
of Thouless pumps in finite-time protocols, which is highly desirable from an
experimental perspective, (ii) in contrast with previous experiments and
numerical works [39, 20, 21], we show that robust Thouless pumps can be
realized even for large repulsive $U$, when realistic values of $V$ are
considered. That is, we find a revival of Thouless pumps caused by the
competition between $V$ and $U$, which is not limited to finite-size systems.
And (iii), the quantization of Thouless pumps in this competition regime is
connected to the presence of a largely debated spontaneously dimerized phase,
the bond-order-wave (BOW) [50, 51, 52, 53, 54, 55, 56, 57, 45, 58, 59], in the
$U-V$ groundstate phase diagram of the extended Fermi-Hubbard model. Finally,
we also discuss a concrete experimental setup and protocol based on magnetic
atoms in an optical lattice, to highlight that the newly introduced Thouless
pumps in the regime of competing interactions can be observed in state-of-the-
art dipolar quantum simulators.

This article is organized as follows. In Sections 2 and 3, we describe,
respectively, the Hamiltonian model and numerical methods used to analyze
interacting Thouless pumps. In Section 4, we present the results of Thouless
pumps in the presence of solely intersite interactions ($U=0$). In Section 5,
we study the competition between intersite and onsite interactions. Section 6
contains the analysis of the connection between the quantization of Thouless
pumps at moderately large interactions and the presence of a BOW phase in the
phase diagram of a parent model. In Section 7, we discuss an experimental
proposal with magnetic atoms, and we conclude in Section 8.

##  2 The Rice-Mele extended Fermi-Hubbard model

We consider an extended version of the Rice-Mele Fermi-Hubbard model, sketched
in Fig. 1(a), describing a chain of length $L$ with $N$ spinful fermions,
labeled by $\sigma=\uparrow,\downarrow$. Here we restrict to the half-filling
case where both $N$ and the total magnetization
$\hat{S}_{z}\equiv\sum_{i}(\hat{n}_{i,\uparrow}-\hat{n}_{i,\downarrow})/2$ are
conserved with $N=L$ and $S_{z}=0$. The full Hamiltonian reads

$\begin{split}&\hat{H}=-\sum_{j=1,\sigma=\uparrow,\downarrow}^{L}\left[J-\delta(-1)^{j}\right](\hat{c}^{\dagger}_{j,\sigma}\hat{c}_{j+1,\sigma}+\textrm{H.c.})+\\\
&+\Delta\sum_{j}(-1)^{j}\hat{n}_{j}+U\sum_{i}\hat{n}_{j,\uparrow}\hat{n}_{j,\downarrow}+V\sum_{j}\hat{n}_{j}\hat{n}_{j+1}.\end{split}$
(1)

Here $J=1$ parametrizes the NN hopping with staggered dimerization strength
$\delta$, and $\Delta$ is the strength of the staggered onsite energy.
Concerning the interactions, $U>0$ accounts for an onsite repulsive Hubbard
term, while $V>0$ describes the repulsion between fermions in NN sites. Note
that $\hat{H}$ faithfully describes a fermionic dipolar mixture trapped in a
one-dimensional optical lattice, in the presence of a superlattice giving rise
to the staggered $\delta$ and $\Delta$ terms.

The noninteracting limit of $\hat{H}$ corresponds to the Rice-Mele model [60],
or to the static Su-Schrieffer-Heegger (SSH) model [61] if chiral and
inversion symmetries are further imposed, i.e., with $\Delta=0$. The SSH model
exhibits symmetry-protected topological and trivial insulators for $\delta<0$
and $\delta>0$, respectively [62]. These phases are characterized by a global
topological invariant, the winding number $\varphi$, which has an integer-
quantized value $\varphi_{\text{topo}}=2$ ($\varphi_{\text{triv}}=0$) in the
topological (trivial) sector [63]. In the presence of inversion symmetry,
$\varphi$ can only change at the singularity $\delta=0$, where a gap closing
associated with a topological phase transition $|\Delta\varphi|=2$ takes
place.

The effect of interactions in the ground-state topological properties of
$\hat{H}$ in the inversion-symmetric case ($\Delta=0$) has attracted great
interest recently [64, 63, 65, 66, 67, 68, 69, 70, 71, 72]. In a nutshell, the
ground state of the system remains topologically protected at $\delta<0$ for
sufficiently small values of $U$ and $V$, as these interactions preserve the
protecting symmetries of the static SSH model. However, any finite $U$ ($V$)
has an impact on the bulk-edge correspondence of the topological phase [65],
and only the edge modes in the charge (spin) sector preserve their topological
degeneracy. Moreover, large interactions are in general detrimental to the
topology. For sufficiently strong $V$, the system features a transition to a
charge-density-wave (CDW), characterized by alternating empty and doubly
occupied sites. Such a phase breaks the inversion and chiral symmetries that
protect the topological phase, leading to a trivial band insulator. Concerning
onsite interaction, for increasing $U$ the system evolves continuously towards
a MI. While such a state does preserve the protecting symmetries, the bulk
spin gap vanishes in the limit $U\to\infty$, and the degenerate spin edge
states merge with the bulk modes. Finally, notice that $\hat{H}$ can also
exhibit a spontaneously dimerized insulator at $\delta=0$, the BOW phase,
arising at the transition between competing orders. For instance, in the
extended $J-U-V$ Fermi-Hubbard limit of $\hat{H}$, which is believed to
capture the physics of several solid-state chains, a BOW appears in the phase
diagram for moderate interactions [50, 51, 52, 53, 54, 55, 56, 57, 45, 58,
59], in a finite region between the CDW and the MI, preventing a direct phase
transition between these two phases 111It is worth mentioning that another
spontaneous BOW phase is also present in the ionic Hubbard limit of $\hat{H}$,
namely at $\delta=V=0$, due to the competition between the onsite repulsion
$U$ and the staggered potential $\Delta$, which we do not consider in this
work. However, we note that this phase is not topologically protected, as the
ionic Hubbard model breaks explicitly the protecting inversion and chiral
symmetries. The latter explains the absence of quantized Thouless pumps in the
regime of the spontaneous BOW phase of the ionic Hubbard model, described in
Ref. [20], and its appearance in the regime of the chiral-symmetric
spontaneous BOW phase of the extended Fermi-Hubbard model discussed here..
Interestingly, such a BOW phase has the same properties as the trivial and
topological ground states of the static SSH model, and therefore it can be
viewed as an interaction-induced topological phase [59]. In this work, we show
that the abovementioned richness of the ground state phase diagram of
$\hat{H}$ in the SSH symmetric case ($\Delta=0$) also leads to a plethora of
dynamical effects in the Rice-Mele regime ($\Delta\neq 0$), in which the
adiabatic breaking of the protecting inversion symmetry can lead to quantized
transport related to the groundstate topology. In particular, we study the
resulting interacting Thouless pumps, with a focus on the unexplored effect of
$V$ and its competition with the onsite term $U$.

##  3 Numerical methods for interacting Thouless pumps

To study Thouless pumps in systems described by $\hat{H}$, we consider the
simultaneous periodic modulation in time of the bond superlattice $\delta$,
and the staggered site potential $\Delta$, which breaks inversion and chiral
symmetries. As sketched in Fig. 1(b), this protocol allows for adiabatic
transport by slowly encircling the singularity $\delta=0$. We parametrize the
Hamiltonian $\hat{H}(\theta)$ for a pump period as

$\begin{split}\delta(\theta)=\delta_{0}\cos(\theta),\
\Delta(\theta)=\Delta_{0}\sin(\theta),\ \theta&\in[0,2\pi].\end{split}$ (2)

Here $\theta$ is related to the time $t$ by $\theta=2\pi t/T$, where $T$ is
the time of a single cycle, and we consider that at $t=0$ the system is in the
ground state of $H(0)$. We are particularly interested in the charge that is
pumped from one edge to the other of the chain during the cycle, which can be
expressed as the integral over time of the instantaneous bulk current, that
is:

$\Delta
Q_{\theta}=\int_{0}^{t(\theta)}\mathcal{J}(t^{\prime})\,\text{dt}^{\prime},$
(3)

where

$\mathcal{J}(t)=\frac{i}{2}\sum_{j\in(0,1),\sigma}J_{j}\left\langle{\hat{c}^{\dagger}_{L/2+j,\sigma}\hat{c}_{L/2+j+1,\sigma}-\text{H.c.}}\right\rangle,$
(4)

is the current operator averaged over a two-site unit cell in the bulk, and
$J_{j}=(J\pm\delta)$. The connection between charge transport and the presence
of a topological invariant comes from the fact that the transferred charge at
the end of a cycle, $\Delta Q_{2\pi}$, can be explicitly written in terms of
the change in the topological invariant across the topological singularity
[73]

$|\Delta Q_{2\pi}|=|\Delta\varphi|\in\mathbb{Z}.$ (5)

We also define the regime of slow time dynamics through the condition
$\left[T\min_{\theta}E_{g}(\theta)\right]^{-1}\ll 1$, where $E_{g}(\theta)$ is
the instantaneous bulk gap of $\hat{H}(\theta)$ during the pump. That is, the
pump cycle frequency has to be much smaller than the lowest excitation of the
system during the pump, to prevent excitations from higher bands. For
instance, in the noninteracting Rice-Mele case, the adiabatic timescale is
given by $\left[2T\min(\Delta_{0},\delta_{0})\right]^{-1}\ll 1$. In the more
complex spinful interacting scenario, we consider the energy gaps at half
filling ($N=L$) with respect to the ground state energy $E_{0}\equiv
E_{0}(S_{z}=0)$. We define both the internal energy gaps $\Delta E_{i}$ in the
unpolarized subspace, as well as the spin gaps $\Delta E^{s}_{i}$
corresponding to flipping one spin as

$\begin{split}\Delta E_{i}\equiv E_{i}(S_{z}=0)-E_{0},\ \Delta E^{s}_{i}\equiv
E_{i}(S_{z}=1)-E_{0}.\end{split}$ (6)

To study the adiabatic dynamics generated by $\hat{H}(\theta)$, we fix the
maximal bond dimerization to $\delta_{0}=0.9J$, which in the absence of
interactions induces a large gap $E_{g}(0)=2\delta_{0}=1.8J$. We also consider
relatively fast cycles of duration $TJ=50$, which are convenient for cold atom
experiments with limited lifetimes. For the numerical calculations, we combine
tensor network simulations in infinite systems using the infinite density-
matrix-renormalization group (DMRG) [74, 48] and time-evolving-block-
decimation (TEBD) algorithms [75], with a maximum bond dimension of
$\chi_{\text{max}}=400$, and exact diagonalization in small periodic systems
to estimate the bulk energy gaps.

##  4 Vanishing onsite interaction $U$

Let us first explore the effect of solely nonlocal $V$ interactions, i.e., at
$U=0$, in a Thouless pump. As depicted in Fig. 1(d), the main effect of a
moderate $V$ term in our model is to suppress the uniform charge distribution
in the regions $\theta=\pi/2$ and $\theta=3\pi/2$. This enhances the CDW
order, and its associated insulating gap, in cooperation with the $\Delta$
term. A direct consequence is that finite values of $V$ can allow one to
reduce the time period $T$ required for having adiabatic (and thus quantized)
Thouless pumps. This behavior can be readily seen in Fig. 2(a), which shows
the evolution of the transferred charge $\Delta Q_{\theta}$ for different
choices of $\Delta_{0}$ and $V$. We observe perfectly quantized pumps for the
case $(\Delta_{0},V)=(J,0)$, as expected for such an adiabatic noninteracting
pump, and also for a finite interaction $(\Delta_{0},V)=(J/4,J)$ (blue lines).
In both cases, the final value of $\Delta Q_{2\pi}=2$ signals the presence of
a topological singularity encircled by the closed adiabatic path, following
Eq. (5). In contrast, note that, in the absence of interactions, the pump with
$(\Delta_{0},V)=(J/4,0)$ is not adiabatic and $\Delta Q_{2\pi}$ reaches a
nonquantized value, as we observe strong oscillations of $\Delta Q_{\theta}$
in Fig. 2(a) (pink line), caused by the coupling to excited states. In this
regime of relatively small $\Delta_{0}=J/4$, we can understand the advantage
of a moderate $V$ in preserving the adiabatic nature of the pump from the
instantaneous gap behavior of $\hat{H}(\theta)$, shown in Fig. 2(b)-(c). In
the noninteracting case of Fig. 2(b), the smallest gaps approach a zero value.
In contrast, in Fig. 2(c) one can see a clear enhancement of the gaps due to
the finite $V$.

Figure 2:  Effect of long-range interactions in Thouless pumps at $U=0$. (a)
Evolution of the pumped charged for different values of $V$ and $\Delta_{0}$,
simulated with the iTEBD algorithm. (b)-(c) First excitation energies of
$\hat{H}(\theta)$ during the pumps with $\Delta_{0}=J/4$ for (b) $\ V=0$, (c)
$V=J$, and (d) $V=4J$, computed with exact diagonalization and $L=8$.

Importantly, note that while $V$ can be helpful to reduce the time cycles of
the pump close to the noninteracting Rice-Mele limit, as sketched in the red-
to-quantized transition of Fig. 1(c), this is not true for arbitrarily large
values of $V$. As sketched in Fig. 1(c), at large $V$ the system exhibits a
spontaneous CDW, and no transport is observed, as can be seen in Fig. 2(a) for
$V=4J$ (green line). In this case, the initially dimerized state at
$\theta=0$, depicted in Fig. 1(d), is replaced by the twofold degenerate
ground state consisting of two states with the form $\ket{02\dots
02},\ket{20\dots 20}$, which have a gap from the rest of the spectrum during
the pump, as shown in Fig. 2(d). We note that, in fact, any finite CDW order
in $\hat{H}(0)$ leads to the breakdown of the topological pump, as Eq. (5) is
only strictly valid when the closed cycle encounters inversion-symmetric
points 222For the similar model considered in Ref. [76], the authors argued
that, even in such a regime dominated by a local spontaneous charge order, one
can observe quantized adiabatic pumps. We note that their scheme relies on the
system being able to jump discontinuously from the low-entangled CDW ground
states, with broken translational symmetry, to the highly entangled symmetric
superposition of the two local CDW ground states. Therefore, while this can
provide valuable insights of the underlying topology in numerical ground state
simulations, such effects cannot be observed in real-time dynamics.. This can
also be understood from the point of view of the noninteracting singularity,
which becomes a twofold degenerate CDW line in the $\delta$ axis, and thus
ceases being encircled by the closed adiabatic path, as depicted in Fig. 4(c).

Finally, we also note another limitation of the gap enhancement driven by the
$V$ term when $\Delta_{0}\to 0$. One could think that, even for
infinitesimally small $\Delta_{0}$, $V$ can amplify the CDW insulating order
and drive an interaction-induced Thouless pump, as indeed suggested in Ref.
[76]. However, we note that such an adiabatic scheme is not possible. In this
limit, $\hat{H}$ preserves inversion symmetry and thus, it cannot induce net
charge currents unless this symmetry is spontaneously broken in a
thermodynamic phase transition. Such a transition to the spontaneous CDW
implies an unavoidable gap closing leading to a dramatic breaking of the
adiabatic condition.

##  5 Interplay between onsite and nonlocal interactions

Let us now consider the combined effect of $U$ and $V$, which constitutes the
primary motivation of this work, due to the possibility of simulating this
regime of competing interactions with magnetic atoms in optical lattices. To
this aim, here we focus on the regime of $U\gg\Delta_{0}$ by fixing $U=8J$ and
$\Delta_{0}=J/2$. We note that, in the absence of $V$, this limit is already
quite well understood [20, 21, 39]: $U$ brings the system into a MI which has
a dramatic effect in Thouless pumps, as the suppression of double occupancies
prevents the breaking of inversion symmetry generated by the $\Delta$ term
during the pump at $\theta=\pi/2$ and $\theta=3\pi/2$ [see Figs. 1(c)-(d)].

Figure 3: Effect of intersite interactions in Thouless pumps at $U=8J$. Here
we also fix $\Delta_{0}=J/2$. (a) Evolution of the pumped charged for
different values of $V$, simulated with the iTEBD algorithm. (b)-(c) First
excitation energies of $\hat{H}(\theta)$ during the pumps for (b) $\ V=2J$,
(c) $V=4J$, and (d) $V=8J$, computed with exact diagonalization and $L=8$.

Figure 3 proves how this scenario is modified when considering also the
presence of $V$. In particular, Fig. 3(a) and Figs. 3(b)-(d) show,
respectively, the pumped charge and instantaneous gaps for fixed $U=8J$ and
$\Delta_{0}=J/2$ and different values of $V$. For a relatively small $V=2J$,
we still observe in Fig. 3(a) the absence of transport (yellow dashed line),
showing that $U$ is still the dominant interaction in the system. The presence
of a MI insulator is further confirmed by the strong reduction of the internal
and spin gap observed in Fig. 3(b). However, when $V$ becomes larger, the
transport is progressively enhanced and there is a revival of the quantized
Thouless pump in a region around $V\sim U/2$ [blue lines of Fig. 3(a)]. By
observing Fig. 1(d), we can explain such interaction-induced Thouless pumps
from the simultaneous presence of large $U$ and $V$. This leads to a
competition between the breakdown caused by the MI [20, 21, 39], and the CDW
breakdown discussed in the previous section. That is, at such a large $U$, the
$V$ term promotes double occupancies at $\theta=\pi/2$ and $3\pi/2$, and
transport is recovered. At the same time, in this regime of strong NN
repulsion, the $U$ term prevents any spontaneous CDW in the initial state at
$\theta=0$. The absence of a MI or spontaneous CDW is also confirmed in Fig.
3(c), which shows a robust gap during the pump. For an even larger $V\gtrsim
4.5J$, $U$ is not able to prevent the emergence of the spontaneous CDW order,
and the quantized pump progressively disappears again, as shown in Fig. 3(a).
This effect is due to the appearance of a phase transition to the twofold
degenerate CDW, as observed in Fig. 3(d). Note that this CDW breakdown is
analogous to the one occurring at $U=0$ described in Sec. 4, but here it takes
place at a larger $V$ due to the competition induced by $U$. We also note that
the behavior of Thouless pumps at $U=8J$ shown in Fig. 3 is found for other
values of $U/J$ until a certain critical value $U_{c}\sim 10J$. Beyond this
critical point, we find either a MI or a CDW breakdown of the pump, without an
intermediate region of quantized transport. The latter is expected, as the
tunneling dynamics do not play any role in the large interaction limit. This
leads to a phase diagram for Thouless pumps in the $U/J-V/J$ plane of the
Hamiltonian $\hat{H}$ in Eq. (1), for fixed values of $\delta_{0}/J=0.9$,
$\Delta_{0}/J=0.5$, and $TJ=50$, sketched in Fig. 1(c).

Figure 4: Sketches of the ground-state singularities and gapless lines that
appear in the $(\Delta,\delta)$ plane for different interaction regimes. The
closed path used for Thouless pumps is also depicted. The orange lines in (b)
and (d) represent a gapless MI. The green lines in (c) and (f) represent a
twofold degenerate CDW. The blue short line in (e) corresponds to a twofold
degenerate BOW. In panels (d)-(f) we considered the regime of finite hopping
$J\sim U,V$.

Interestingly, one can also understand the interacting Thouless pumps shown in
Fig. 3 and Fig. 1(c) from the behavior of the topological singularity. To this
end, in Fig. 4 we sketch the evolution of the topological singularities in the
$(\Delta,\delta)$ plane in the presence of interactions. For dominant $U$, the
noninteracting singular point becomes the gapless MI vertical line [21] of
Fig. 4(b), which is not enclosed by the closed path. Alternatively, a dominant
$V$ leads to the horizontal twofold degenerate line of Fig. 4(c),
corresponding to the spontaneous CDW order. For a very large interaction scale
$U,V\gg J$, increasing $V$ in the MI regime [Fig. 4(b)] leads to a first-order
transition [45] to the CDW regime of Fig. 4(c) around $U\approx 2V$. On the
contrary, when interactions are of the order of the hopping, $U,V\approx J$,
the enhanced quantum fluctuations lead to a competing regime in which such a
first-order transition is replaced by a progressive reduction of the gapless
lines, and an intermediate BOW phase [45], as sketched in Figs. 4(d)-(f). In
this case, the closed path encircles the shorter gapless lines corresponding
to the MI, BOW or CDW regimes, and quantized transport is recovered.

Remarkably, the above discussion shows that the behavior of Thouless pumps in
the $U/J-V/J$ phase diagram sketched in Fig. 1(c) is very reminiscent of the
ground-state phase diagram of the extended Fermi-Hubbard ($J-U-V$) model [50,
51, 52, 53, 54, 55, 56, 57, 45, 58, 59], where the spontaneous BOW phase
appears in a similar region as the quantized Thouless pumps. In what follows,
we explore this connection and prospects to use it as a tool for detecting the
spontaneous BOW phase.

##  6 Thouless pumps and the bond-order-wave phase

The fact that, in the regime given by $V\sim U/2,\ U<10J$, we can recover
quantized Thouless pumps by increasing the strength of the $V$ interaction at
fixed $U$ cannot be understood from a simple noninteracting picture. However,
we can explain this phenomenon by noting that, in the absence of the
$\Delta,\delta$ terms used during the pump, the ground state of the
Hamiltonian (1) is a spontaneous BOW phase [50, 51, 52, 53, 54, 55, 56, 57,
45, 58] also in the regime $V\sim U/2,\ U<10J$, as sketched in Fig 1(c). In
such a phase, which is also induced by the competition between the MI and CDW
orders, there is a ground-state bulk degeneracy between a topological and a
trivial state [59]. From the perspective of Thouless pumps, in the BOW region,
the main effect of $\delta$ is to break the degeneracy between the two
spontaneous BOW ground states and to control/guide its spontaneous symmetry
breaking during the pump, such that the system is able to explore the trivial
and topological dimerized sectors.

Figure 5: (a) Accuracy of Thouless pumps (iTEBD) as a function of the nonlocal
interaction $V/J$ and at different cycle numbers $C$, for fixed $U=8J$ and
$\Delta_{0}=J/2$. We only observe quantized pumps independent of $C$ around
the narrow region where the groundstate of the $J-U-V$ model (with
$\delta=\Delta=0$) is the spontaneous BOW phase (gray shaded area, extracted
from Ref. [45]). The latter can be more clearly seen with the inset
logarithmic scale. The nonomonotonic behavior of $|\delta Q|$ in the window
$V/J\in(4.5,5.5)$ is related to the spontaneous CDW phase transition, as
detailed in the main text.

Interestingly, one can therefore interpret the quantized Thouless pumping as
an indirect proof of the presence of a BOW phase, whose detection and
engineering remain elusive in cold atoms experiments, despite recent
theoretical proposals [77].

To illustrate such a relation more explicitly, in Fig. 5 we plot the breaking
of the quantization in Thouless pumps consisting of $C$ consecutive periodic
cycles, defined as $\delta Q\equiv|\Delta Q_{2\pi C}|-2C$, for a fixed $U=8J$
and as a function of $V$, which corresponds to the dashed line in Fig. 1(c).
Remarkably, we observe perfectly quantized and cycle-independent transport,
i.e., $|\delta Q|=0$, for Thouless pumps in the region $V/J\in(3.9,4.3)$
around the narrow window of the BOW phase $V/J\in(4,4.15)$. As sketched in
Fig. 1(c), the same behavior is also observed for other values of $U/J$. Two
important comments are in order at this point, concerning the quantitative
behavior of the curves in Fig. 5. In first place, note that the fact that the
BOW region does not coincide exactly with the region of quantized pumps is
expected. This is because, in the adiabatic dynamics considered here, the
staggered field $\delta$ ($\Delta$) used to perform the pump prevents the CDW
(MI) in a slightly wider range of $V$ values than in the ground-phase diagram
of the model in the absence of such guiding fields [45]. In second place, the
nonmonotonic behavior of $|\delta Q|$ after the quantized plateau in Fig. 5 is
explained by the fact that there is a phase transition to a local CDW order at
$V/J\approx 4.75$. At this critical point, $|\delta Q|$ exhibits a local
maximum due to the strong deviation from the adiabatic regime in our finite-
time protocol. The latter also leads to dynamical results that slighlty depend
on the bond dimension $\chi_{\text{max}}$. For $V$ larger than the critical
point, the system reenters the adiabatic regime due to the insulating nature
of the CDW phase. Inside this CDW phase, increasing $V$ enhances the gap by
moving away from the critical point, thus favoring adiabatic transport, but it
also enhances the amplitude of the CDW spontaneous order. While the
competition between these two phenomena initially leads to the
counterintuitive decrease of $|\delta Q|$ for increasing $V$, we emphasize
that in this CDW regime Thouless pumps do not encircle a topological
singularity and cannot be associated to any topological invariant. That is,
$|\delta Q|$ is cycle-dependent and not quantized. Further increasing $V/J$
progressively leads to a maximally polarized CDW phase, suppressing any
hopping (transport) process during the pump.

As a general remark, note that the scheme depicted in Fig. 5 overcomes the
main problems that make it difficult to detect the BOW phase of the extended
Fermi-Hubbard model from direct measurements. First, with Thouless pumps the
topology of the interacting model (1) can be probed by measuring the average
global charge transport in the system, a topological property robust to small
imperfections such as local disorder or low finite temperature. In contrast, a
direct measurement scheme of the topological nature of the BOW phase would
require, e.g., measuring nonlocal correlators [78, 59] with single-site
resolution through a quantum gas microscope [79, 80, 81]. Similarly, even a
direct measurement of the BOW local order parameter requires involved probes
capable of measuring the bond density in-situ, or superlattice modulation
spectroscopy [77]. Second, the spontaneous nature of the BOW phase, and its
very small intrinsic spin gap, require very small temperatures and complex
state-preparation protocols to avoid the presence of solitonic excitations. On
the contrary, in Thouless pumps the spin gap is enhanced due to the presence
of the bond guiding field $\delta$, which also lifts the groundstate
degeneracy of the BOW phase, thus preventing the appearance of solitons.
Third, as shown in Fig. 1(c) and Fig. 5, within experimental errors quantized
pumps are observed in a much wider parameter region than the BOW phase.

##  7 Experimental proposal with magnetic atoms

A prime candidate for an experimental realization of the Rice-Mele extended
Fermi-Hubbard Hamiltonian (1) are quantum simulators using ultracold
lanthanide atoms in optical lattices [82, 83]. For erbium, the preparation of
a two-component Fermi gas in a lattice with long lifetimes and broad interspin
Feshbach resonances for interaction control has been already demonstrated
[84]. To manipulate the control parameters $\Delta$ and $\delta$ of the
Thouless pump, a superlattice similar to experiments with neutral atoms can be
used, see e.g., Refs. [37, 39]. Such a superlattice can also be used to
prepare the initial state e.g., a dimerized chain along the superlattice
direction. As the proposed scenario is based on single one-dimensional
systems, experimental setups need sufficient spacing between individual tubes
to avoid strong coupling between the tubes due to the long-range character of
the dipole-dipole interaction. This can be realized either by choosing a
larger lattice spacing perpendicular to the superlattice, or by isolating
single stripes in a two-dimensional lattice plane via the removal of atoms in
between using single-site addressing of a quantum gas microscope.

As already discussed in Ref. [59], taking realistic values currently
achievable in experiments should allow one to reach the quantized pump regime
with e.g., $U=4J$ and $V=2J$ [see Fig. 1(c)] using erbium in a lattice with
$532\,$nm unit cell ($266\,$nm short lattice spacing), giving a maximum
interaction strength of $V/h\approx 30\,$Hz. This can be further increased by
using dysprosium and shorter wavelength lattices as planned for future
experiments [85, 86], reaching interaction strengths of up to $V/h\approx
200\,$Hz. To be able to observe the revival it would be preferential to probe
the dynamics with and without $V$. To probe the dynamics without $V$, one can
either rotate the quantization axis to the magic angle, fully eliminating $V$,
or scale $U$ and $J$ up considerably $(V\ll J,U$) to minimize the influence.
The regime $V\gg U$ can also be achieved, since $U$ can be tuned independently
by means of Feshbach resonances [87]. The induced transport due to the
Thouless pump mechanism can be directly measured using the center-of-mass
displacement of the atoms using either high-resolution in-situ imaging (see
e.g. [37]) or a quantum gas microscope. For $J=2\pi\times 100\text{Hz}$, the
adiabatic condition $TJ=50$ used in our numerical simulations would be
achieved for realistic pump periods of half a second. We also expect that the
averaged center-of-mass displacement is robust against defects appearing
during the dynamics, making the lifetime requirements less stringent.

##  8 Conclusions and outlook

We investigated the appearance of quantized topological transport, i.e.,
Thouless pumping, in a Rice-Mele extended Fermi-Hubbard model where fermions
are subject to competing interactions of different ranges. Contrary to
previous analysis showing a detrimental effect of onsite repulsive
interactions to quantized Thouless pumps [20, 21, 39], we unveiled that
including nonlocal repulsions leads to novel phenomena.

Specifically, in configurations with only intersite repulsive interactions
$V$, we showed that a moderate fermionic repulsion favors the presence of
quantized Thouless pumps. We proved this phenomenon to be totally induced by
this specific interaction which tends to increase the instantaneous energy
gaps and thus stabilize the system topology. In contrast, at large nonlocal
interactions, we find an expected breakdown of the transport due to the
appearance of a charge-density-wave.

When further considering a repulsive onsite Hubbard $U$, which is known to
destroy the topological transport [20, 21, 39], we still identified a sizeable
region where topological transport persists. Here the nonlocal repulsion
counteracts the destructive role of a moderately large local repulsion and
vice versa, and therefore the quantized transport remains stable. Noticeably,
this competing mechanism allows for a robust quantized transport at
interactions up to one order of magnitude larger than the noninteracting
energy, i.e., at $U/J\sim 10$.

Interestingly, this stabilization of topological Thouless pumps through
competing interactions is reminiscent of the ground state properties of the
parent extended Fermi-Hubbard model. There, the presence of moderate and
competing contact and intersite interactions favors the appearance of a
spontaneous bond-order-wave phase, whose topological properties have been
recently revealed [59]. Based on this fact, we thus expect that quantized
topological transport can take place in different strongly correlated models
where competing interactions give rise to spontaneously generated bond-order-
wave regimes [17]. Beyond our numerical analysis, we have been further able to
propose a realistic experimental scheme where our findings can be probed. The
latter relies on trapping an ultracold mixture of magnetic atoms in an optical
lattice. The advantage of using such particles comes from the fact that they
allow for a fine-tuning of both contact and intersite interaction, thus
providing access to simulate the Hamiltonian of Eq. (1) in a wide regime of
parameters. Importantly, note that these configurations where fermions are
subject to local and nonlocal repulsions are the most common in solid-state
platforms. Moreover, whilst we focused on repulsive interactions and charge
transport, our numerical analysis and experimental proposal could be adapted
to study the role of attractive terms and spin transport in the same model.
We, therefore, believe that our results can have an important impact on a
large variety of physical systems, thus representing an essential step towards
a complete understanding of many-body phases of matter characterized by
nontrivial topological properties.

Acknowledgments. The authors acknowledge fruitful discussions with T.
Esslinger, Z. Zhu, A. Dauphin, F. Heidrich-Meisner, A. A. Aligia, and E.
Bertok. The DMRG and TEBD calculations were performed using the TeNPy Library
[49]. ICFO group acknowledges support from: ERC AdG NOQIA; MICIN/AEI
(PGC2018-0910.13039/501100011033, CEX2019-000910-S/10.13039/501100011033, Plan
National FIDEUA PID2019-106901GB-I00, FPI; MICIIN with funding from European
Union NextGenerationEU (PRTR-C17.I1): QUANTERA MAQS PCI2019-111828-2);
MCIN/AEI/ 10.13039/501100011033 and by the “European Union NextGeneration
EU/PRTR” QUANTERA DYNAMITE PCI2022-132919 within the QuantERA II Programme
that has received funding from the European Union’s Horizon 2020 research and
innovation programme under Grant Agreement No 101017733Proyectos de I+D+I
“Retos Colaboración” QUSPIN RTC2019-007196-7); Fundació Cellex; Fundació Mir-
Puig; Generalitat de Catalunya (European Social Fund FEDER and CERCA program,
AGAUR Grant No. 2021 SGR 01452, QuantumCAT U16-011424, co-funded by ERDF
Operational Program of Catalonia 2014-2020); Barcelona Supercomputing Center
MareNostrum (FI-2023-1-0013); EU (PASQuanS2.1, 101113690); EU Horizon 2020
FET-OPEN OPTOlogic (Grant No 899794); EU Horizon Europe Program (Grant
Agreement 101080086 — NeQST), National Science Centre, Poland (Symfonia Grant
No. 2016/20/W/ST4/00314); ICFO Internal “QuantumGaudi” project; European
Union’s Horizon 2020 research and innovation program under the Marie-
Skłodowska-Curie grant agreement No 101029393 (STREDCH) and No 847648 (“La
Caixa” Junior Leaders fellowships ID100010434: LCF/BQ/PI19/11690013,
LCF/BQ/PI20/11760031,LCF/BQ/PR20/ 11770012, LCF/BQ/PR21/11840013). Views and
opinions expressed are, however, those of the author(s) only and do not
necessarily reflect those of the European Union, European Commission, European
Climate, Infrastructure and Environment Executive Agency (CINEA), nor any
other granting authority. Neither the European Union nor any granting
authority can be held responsible for them. The Innsbruck group acknowledges
support from a QuantERA grant MAQS from the Austrian Science Fund (FWF No.
I4391-N). L. B. acknowledges financial support within the DiQut Grant No.
2022523NA7 funded by European Union – Next Generation EU, PRIN 2022 program
(D.D. 104 - 02/02/2022 Ministero dell’Università e della Ricerca).

## References

  * Klitzing  _et al._ [1980] K. v. Klitzing, G. Dorda, and M. Pepper, Phys. Rev. Lett. 45, 494 (1980).
  * Thouless  _et al._ [1982a] D. J. Thouless, M. Kohmoto, M. P. Nightingale, and M. den Nijs, Phys. Rev. Lett. 49, 405 (1982a).
  * Hasan and Kane [2010] M. Z. Hasan and C. L. Kane, Rev. Mod. Phys. 82, 3045 (2010).
  * Chiu  _et al._ [2016] C.-K. Chiu, J. C. Y. Teo, A. P. Schnyder, and S. Ryu, Rev. Mod. Phys. 88, 035005 (2016).
  * Landau  _et al._ [1999] L. D. Landau, E. M. Lifshitz, and M. Pitaevskii,  _Statistical Physics_ (Butterworth-Heinemann, New York, 1999).
  * Wilson and Kogut [1974] K. G. Wilson and J. Kogut, Phys. Rep. 12, 75 (1974).
  * von Klitzing [2017] K. von Klitzing, Nat. Phys. 13, 198 (2017).
  * Nayak  _et al._ [2008] C. Nayak, S. H. Simon, A. Stern, M. Freedman, and S. Das Sarma, Rev. Mod. Phys. 80, 1083 (2008).
  * Rachel [2018] S. Rachel, Rep. Prog. Phys. 81, 116501 (2018).
  * Thouless [1983] D. J. Thouless, Phys. Rev. B 27, 6083 (1983).
  * Niu and Thouless [1984] Q. Niu and D. J. Thouless, Journal of Physics A: Mathematical and General 17, 2453 (1984).
  * Berg  _et al._ [2011] E. Berg, M. Levin, and E. Altman, Phys. Rev. Lett. 106, 110405 (2011).
  * Greschner  _et al._ [2020] S. Greschner, S. Mondal, and T. Mishra, Phys. Rev. A 101, 053630 (2020).
  * Hayward  _et al._ [2018] A. Hayward, C. Schweizer, M. Lohse, M. Aidelsburger, and F. Heidrich-Meisner, Phys. Rev. B 98, 245148 (2018).
  * Mondal  _et al._ [2021] S. Mondal, S. Greschner, L. Santos, and T. Mishra, Phys. Rev. A 104, 013315 (2021).
  * Lin  _et al._ [2020a] L. Lin, Y. Ke, and C. Lee, Phys. Rev. A 101, 023620 (2020a).
  * Mondal  _et al._ [2022a] S. Mondal, A. Padhan, and T. Mishra, Phys. Rev. B 106, L201106 (2022a).
  * Kuno and Hatsugai [2020] Y. Kuno and Y. Hatsugai, Phys. Rev. Res. 2, 042024 (2020).
  * Padhan  _et al._ [2023] A. Padhan, S. Mondal, S. Vishveshwara, and T. Mishra, “Interacting bosons on a Su-Schrieffer-Heeger ladder: Topological phases and Thouless pumping,” (2023), arXiv:2306.09325 [cond-mat.quant-gas] .
  * Nakagawa  _et al._ [2018] M. Nakagawa, T. Yoshida, R. Peters, and N. Kawakami, Phys. Rev. B 98, 115147 (2018).
  * Bertok  _et al._ [2022] E. Bertok, F. Heidrich-Meisner, and A. A. Aligia, Phys. Rev. B 106, 045141 (2022).
  * Mondal  _et al._ [2022b] S. Mondal, E. Bertok, and F. Heidrich-Meisner, Phys. Rev. B 106, 235118 (2022b).
  * Mondal  _et al._ [2023] S. Mondal, E. Bertok, and F. Heidrich-Meisner, Phys. Rev. B 107, 239903 (2023).
  * Feynman [1982] R. P. Feynman, Int. J. Theor. Phys. 21, 467 (1982).
  * Cirac and Zoller [2012] J. I. Cirac and P. Zoller, Nat. Phys. 8, 264 (2012).
  * Georgescu  _et al._ [2014] I. M. Georgescu, S. Ashhab, and F. Nori, Rev. Mod. Phys. 86, 153 (2014).
  * Daley  _et al._ [2022] A. J. Daley, I. Bloch, C. Kokail, S. Flannigan, N. Pearson, M. Troyer, and P. Zoller, Nature 607, 667 (2022).
  * Altman  _et al._ [2021] E. Altman, K. R. Brown, G. Carleo, L. D. Carr, E. Demler, C. Chin, B. DeMarco, S. E. Economou, M. A. Eriksson, K.-M. C. Fu, M. Greiner, K. R. Hazzard, R. G. Hulet, A. J. Kollár, B. L. Lev, M. D. Lukin, R. Ma, X. Mi, S. Misra, C. Monroe, K. Murch, Z. Nazario, K.-K. Ni, A. C. Potter, P. Roushan, M. Saffman, M. Schleier-Smith, I. Siddiqi, R. Simmonds, M. Singh, I. Spielman, K. Temme, D. S. Weiss, J. Vučković, V. Vuletić, J. Ye, and M. Zwierlein, PRX Quantum 2, 017003 (2021).
  * Cooper  _et al._ [2019] N. R. Cooper, J. Dalibard, and I. B. Spielman, Rev. Mod. Phys. 91, 015005 (2019).
  * Citro and Aidelsburger [2023] R. Citro and M. Aidelsburger, Nat. Rev. Phys. 5, 87 (2023).
  * Ozawa  _et al._ [2019] T. Ozawa, H. M. Price, A. Amo, N. Goldman, M. Hafezi, L. Lu, M. C. Rechtsman, D. Schuster, J. Simon, O. Zilberberg, and I. Carusotto, Rev. Mod. Phys. 91, 015006 (2019).
  * Kraus  _et al._ [2012] Y. E. Kraus, Y. Lahini, Z. Ringel, M. Verbin, and O. Zilberberg, Phys. Rev. Lett. 109, 106402 (2012).
  * Cerjan  _et al._ [2020] A. Cerjan, M. Wang, S. Huang, K. P. Chen, and M. C. Rechtsman, Light: Science & Applications 9, 178 (2020).
  * Jürgensen  _et al._ [2021] M. Jürgensen, S. Mukherjee, and M. C. Rechtsman, Nature 596, 63 (2021).
  * Jürgensen  _et al._ [2023] M. Jürgensen, S. Mukherjee, C. Jörg, and M. C. Rechtsman, Nat. Phys. 19, 420 (2023).
  * Lohse  _et al._ [2016] M. Lohse, C. Schweizer, O. Zilberberg, M. Aidelsburger, and I. Bloch, Nat. Phys. 12, 350 (2016).
  * Nakajima  _et al._ [2016] S. Nakajima, T. Tomita, S. Taie, T. Ichinose, H. Ozawa, L. Wang, M. Troyer, and Y. Takahashi, Nat. Phys. 12, 296 (2016).
  * Minguzzi  _et al._ [2022] J. Minguzzi, Z. Zhu, K. Sandholzer, A.-S. Walter, K. Viebahn, and T. Esslinger, Phys. Rev. Lett. 129, 053201 (2022).
  * Walter  _et al._ [2023] A.-S. Walter, Z. Zhu, M. Gächter, J. Minguzzi, S. Roschinski, K. Sandholzer, K. Viebahn, and T. Esslinger, Nat. Phys. 19, 1471 (2023).
  * Viebahn  _et al._ [2023] K. Viebahn, A.-S. Walter, E. Bertok, Z. Zhu, M. Gächter, A. A. Aligia, F. Heidrich-Meisner, and T. Esslinger, “Interaction-induced charge pumping in a topological many-body system,” (2023), arXiv:2308.03756 [cond-mat.quant-gas] .
  * Lewenstein  _et al._ [2012] M. Lewenstein, A. Sanpera, and V. Ahufinger,  _Ultracold Atoms in Optical Lattices: Simulating Quantum many-body systems_, Vol. 54 (Oxford University Press, Oxford, 2012).
  * Bloch  _et al._ [2008] I. Bloch, J. Dalibard, and W. Zwerger, Rev. Mod. Phys. 80, 885 (2008).
  * Sompet  _et al._ [2022] P. Sompet, S. Hirthe, D. Bourgund, T. Chalopin, J. Bibo, J. Koepsell, P. Bojović, R. Verresen, F. Pollmann, G. Salomon, C. Gross, T. A. Hilker, and I. Bloch, Nature 606, 484 (2022).
  * Léonard  _et al._ [2023] J. Léonard, S. Kim, J. Kwan, P. Segura, F. Grusdt, C. Repellin, N. Goldman, and M. Greiner, Nature 619, 495 (2023).
  * Ejima and Nishimoto [2007] S. Ejima and S. Nishimoto, Phys. Rev. Lett. 99, 216403 (2007).
  * Lahaye  _et al._ [2009] T. Lahaye, C. Menotti, L. Santos, M. Lewenstein, and T. Pfau, Rep. Prog. Phys. 72, 126401 (2009).
  * Chomaz  _et al._ [2022] L. Chomaz, I. Ferrier-Barbut, F. Ferlaino, B. Laburthe-Tolra, B. L. Lev, and T. Pfau, Reports on Progress in Physics 86, 026401 (2022).
  * Schollwöck [2011] U. Schollwöck, Ann. Phys. 326, 96 (2011).
  * Hauschild and Pollmann [2018] J. Hauschild and F. Pollmann, SciPost Phys. Lect. Notes , 5 (2018).
  * Nakamura [1999] M. Nakamura, J. Phys. Soc. Japan 68, 3123 (1999).
  * Nakamura [2000] M. Nakamura, Phys. Rev. B 61, 16377 (2000).
  * Jeckelmann [2002] E. Jeckelmann, Phys. Rev. Lett. 89, 236401 (2002).
  * Sengupta  _et al._ [2002] P. Sengupta, A. W. Sandvik, and D. K. Campbell, Phys. Rev. B 65, 155113 (2002).
  * Sandvik  _et al._ [2004] A. W. Sandvik, L. Balents, and D. K. Campbell, Phys. Rev. Lett. 92, 236401 (2004).
  * Zhang [2004] Y. Z. Zhang, Phys. Rev. Lett. 92, 246404 (2004).
  * Tam  _et al._ [2006] K.-M. Tam, S.-W. Tsai, and D. K. Campbell, Phys. Rev. Lett. 96, 036408 (2006).
  * Glocke  _et al._ [2007] S. Glocke, A. Klümper, and J. Sirker, Phys. Rev. B 76, 155121 (2007).
  * Di Dio  _et al._ [2014] M. Di Dio, L. Barbiero, A. Recati, and M. Dalmonte, Phys. Rev. A 90, 063608 (2014).
  * Julià-Farré  _et al._ [2022] S. Julià-Farré, D. González-Cuadra, A. Patscheider, M. J. Mark, F. Ferlaino, M. Lewenstein, L. Barbiero, and A. Dauphin, Phys. Rev. Res. 4, L032005 (2022).
  * Rice and Mele [1982] M. J. Rice and E. J. Mele, Phys. Rev. Lett. 49, 1455 (1982).
  * Su  _et al._ [1979] W. P. Su, J. R. Schrieffer, and A. J. Heeger, Phys. Rev. Lett. 42, 1698 (1979).
  * Ryu  _et al._ [2010] S. Ryu, A. P. Schnyder, A. Furusaki, and A. W. W. Ludwig, New J. Phys. 12, 065010 (2010).
  * Manmana  _et al._ [2012] S. R. Manmana, A. M. Essin, R. M. Noack, and V. Gurarie, Phys. Rev. B 86, 205119 (2012).
  * Gurarie [2011] V. Gurarie, Phys. Rev. B 83, 085426 (2011).
  * Yoshida  _et al._ [2014] T. Yoshida, R. Peters, S. Fujimoto, and N. Kawakami, Phys. Rev. Lett. 112, 196404 (2014).
  * Wang  _et al._ [2015] D. Wang, S. Xu, Y. Wang, and C. Wu, Phys. Rev. B 91, 115118 (2015).
  * Ye  _et al._ [2016] B.-T. Ye, L.-Z. Mu, and H. Fan, Phys. Rev. B 94, 165167 (2016).
  * Sbierski and Karrasch [2018] B. Sbierski and C. Karrasch, Phys. Rev. B 98, 165101 (2018).
  * Barbiero  _et al._ [2018] L. Barbiero, L. Santos, and N. Goldman, Phys. Rev. B 97, 201115 (2018).
  * Le  _et al._ [2020] N. H. Le, A. J. Fisher, N. J. Curson, and E. Ginossar, npj Quantum Inf. 6, 24 (2020).
  * Lin  _et al._ [2020b] Y.-T. Lin, D. M. Kennes, M. Pletyukhov, C. S. Weber, H. Schoeller, and V. Meden, Phys. Rev. B 102, 085122 (2020b).
  * Montorsi  _et al._ [2022] A. Montorsi, U. Bhattacharya, D. González-Cuadra, M. Lewenstein, G. Palumbo, and L. Barbiero, Phys. Rev. B 106, L241115 (2022).
  * Thouless  _et al._ [1982b] D. J. Thouless, M. Kohmoto, M. P. Nightingale, and M. den Nijs, Phys. Rev. Lett. 49, 405 (1982b).
  * White [1992] S. R. White, Phys. Rev. Lett. 69, 2863 (1992).
  * Orús and Vidal [2008] R. Orús and G. Vidal, Phys. Rev. B 78, 155117 (2008).
  * Marks  _et al._ [2021] J. A. Marks, M. Schüler, J. C. Budich, and T. P. Devereaux, Phys. Rev. B 103, 035112 (2021).
  * Loida  _et al._ [2017] K. Loida, J.-S. Bernier, R. Citro, E. Orignac, and C. Kollath, Phys. Rev. Lett. 119, 230403 (2017).
  * Barbiero  _et al._ [2013] L. Barbiero, A. Montorsi, and M. Roncaglia, Phys. Rev. B 88, 035109 (2013).
  * Bakr  _et al._ [2009] W. S. Bakr, J. I. Gillen, A. Peng, S. Fölling, and M. Greiner, Nature 462, 74 (2009).
  * Endres  _et al._ [2011] M. Endres, M. Cheneau, T. Fukuhara, C. Weitenberg, P. Schauß, C. Gross, L. Mazza, M. C. Bañuls, L. Pollet, I. Bloch, and S. Kuhr, Science 334, 200 (2011).
  * Hilker  _et al._ [2017] T. A. Hilker, G. Salomon, F. Grusdt, A. Omran, M. Boll, E. Demler, I. Bloch, and C. Gross, Science 357, 484 (2017).
  * Patscheider  _et al._ [2020] A. Patscheider, B. Zhu, L. Chomaz, D. Petter, S. Baier, A.-M. Rey, F. Ferlaino, and M. J. Mark, Phys. Rev. Research 2, 023050 (2020).
  * Su  _et al._ [2023] L. Su, A. Douglas, M. Szurek, R. Groth, S. F. Ozturk, A. Krahn, A. H. Hébert, G. A. Phelps, S. Ebadi, S. Dickerson, F. Ferlaino, O. Marković, and M. Greiner, Nature 622, 724 (2023).
  * Baier  _et al._ [2018] S. Baier, D. Petter, J. H. Becher, A. Patscheider, G. Natale, L. Chomaz, M. J. Mark, and F. Ferlaino, Phys. Rev. Lett. 121, 093602 (2018).
  * Fraxanet  _et al._ [2022] J. Fraxanet, D. González-Cuadra, T. Pfau, M. Lewenstein, T. Langen, and L. Barbiero, Phys. Rev. Lett. 128, 043402 (2022).
  * Sohmen  _et al._ [2023] M. Sohmen, M. J. Mark, M. Greiner, and F. Ferlaino, SciPost Phys. 15, 182 (2023).
  * Lange  _et al._ [2009] A. D. Lange, K. Pilch, A. Prantner, F. Ferlaino, B. Engeser, H.-C. Nägerl, R. Grimm, and C. Chin, Phys. Rev. A 79, 013622 (2009).


∎

11institutetext:  D. Wingard 22institutetext: Department of Physics and
Astronomy
California State University Long Beach
Long Beach, California, USA 33institutetext: B. Kónya 44institutetext:
Department of Physics
Lund University, Lund, Sweden 55institutetext: Z. Papp∗ 66institutetext:
Department of Physics and Astronomy
California State University Long Beach
Long Beach, California, USA
66email<EMAIL_ADDRESS>
# Calculation of Relativistic Single-Particle States

D. Wingard  B. Kónya  Z. Papp∗

(Received: date / Accepted: date)

###### Abstract

A computational method is proposed to calculate bound and resonant states by
solving the Klein-Gordon and Dirac equations for real and complex energies,
respectively. The method is an extension of a non-relativistic one, where the
potential is represented in a Coulomb-Sturmian basis. This basis facilitates
the exact analytic evaluation of the Coulomb Green’s operator in terms of a
continued fraction. In the extension to relativistic problems, we cast the
Klein-Gordon and Dirac equations into an effective Schrödinger form. Then the
solution method is basically an analytic continuation of non-relativistic
quantities like the angular momentum, charge, energy and potential into the
effective relativistic counterparts.

###### Keywords:

Relativistic quantum mechanics Klein-Gordon equation Dirac equation Resonances
Integral equation Separable interactions Analytic continuation Continued
fraction

##  1 Introduction

The most often used theoretical tool for atomic and nuclear physics is quantum
mechanics, mostly its non-relativistic version. Usually, the effects of
relativity are taken into account as the non-relativistic limit of the
relativistic equations. While there are a large number of methods for non-
relativistic calculations, methods for relativistic calculations are rather
scarce.

The aim of this work is to generalize an approximation method, that has been
rather successful in non-relativistic quantum mechanics, to relativistic
calculations. We present a method that can equally solve the Schrödinger, the
Klein-Gordon and the Dirac equations for bound and resonant states and for
Coulomb plus short range potentials.

The computational method has been developed a while ago papp1987bound ;
papp1988potential ; PhysRevC.61.034302  and has been applied for solving
various problems in nuclear and atomic physics like the Faddeev equation with
Coulomb-like interactions PhysRevC.54.50 ; PhysRevC.55.1080 ;
PhysRevA.63.062721 .

In this work, we want to develop a computational method that will allow us to
incorporate relativistic quantum mechanics in our studies. The design of this
work is as follows. In Section 2 we review the method applied to the
Schrödinger equation. Then, in Sections 3 and 4 we show how to extend it for
solving the Klein-Gordon and the Dirac equations, respectively. In Section 5
we present some numerical illustrations and in Section 6 we summarize our
findings.

##  2 Solution of the Schrödinger equation

We consider a Hamiltonian with a Coulomb plus short-range potential in angular
momentum channel $l$

$h_{l}=h_{l}^{0}+Z/r+v_{l}^{(s)},$ (1)

where $h_{l}^{0}$ is the non-relativistic kinetic energy operator, $Z$ is the
charge number, and $v^{(s)}$ is a short-range potential. This Hamiltonian
gives us the Schrödinger eigenvalue equation

$(h_{l}^{0}+Z/r+v_{l}^{(s)})|\psi_{l}\rangle=E|\psi_{l}\rangle.$ (2)

If we represent the momentum operator by a derivative in terms of the spatial
variable, the eigenvalue problem becomes a differential equation, which can be
solved with the appropriate boundary conditions.

We can also cast the Schrödinger eigenvalue equation into a Lippmann-Schwinger
form. If we are concerned about bound and resonant states, we need to solve
the homogenous Lippmann-Schwinger equation

$|\psi_{l}\rangle=g_{l}^{C}(Z,E)v_{l}^{(s)}|\psi_{l}\rangle~{},$ (3)

with negative real and positive complex energies, respectively. Here

$g_{l}^{C}(Z,E)=(E-h_{l}^{0}-Z/r)^{-1}$ (4)

is the Coulomb Green’s operator.

We solve the Lippmann-Schwinger equation by approximating the short-range
potential $v_{l}^{(s)}$ on a Hilbert-space basis. For that purpose we take the
Coulomb-Sturmian (CS) basis. The CS functions in angular momentum $l$ are
defined by

$\langle
r|nl\rangle=\left(\frac{{\Gamma(n+1)}}{{\Gamma(n+2l+2)}}\right)^{1/2}\exp(-br)(2br)^{l+1}L_{n}^{2l+1}(2br)~{},$
(5)

where $n=0,1,2,\dots$, $L$ is the Laguerre polynomial and $b$ is a parameter.
Together with $\langle r|\widetilde{nl}\rangle=\langle r|{nl}\rangle/r$ these
functions are orthonormal $\langle
nl|\widetilde{n^{\prime}l}\rangle=\delta_{nn^{\prime}}$ and form a complete
set
$\lim_{N\to\infty}\sum_{n=0}^{N}|nl\rangle\langle\widetilde{n^{\prime}l}|=1$.

The finite dimensional representation of the short-range potential is given by

$v_{l}^{(s)}\approx\sum_{nn^{\prime}}^{N}|\widetilde{nl}\rangle\underline{\tilde{v}}^{(s)}_{l,nn^{\prime}}\langle\widetilde{n^{\prime}l}|~{}.$
(6)

To construct the matrix $\underline{\tilde{v}}^{(s)}_{l,nn^{\prime}}$ we
calculate the matrix elements
$\underline{v}^{(N^{\prime})}=\langle{nl}|v_{l}^{(s)}|n^{\prime}l\rangle$ up
to $N^{\prime}\geq N$, numerically in general, invert
$\underline{v}^{(N^{\prime})}$, then truncate the $N^{\prime}\times
N^{\prime}$ inverse matrix to a $N\times N$ matrix and invert again to obtain
the $N\times N$ matrix $\underline{\tilde{v}}^{(s)}_{l,nn^{\prime}}$
brown2013approximations .

With this approximation, the Lippmann-Schwinger equation (3) becomes

$|\psi_{l}\rangle=\sum_{nn^{\prime}}^{N}g_{l}^{C}(Z,E)|\widetilde{nl}\rangle\underline{\tilde{v}}^{(s)}_{l,nn^{\prime}}\langle\widetilde{n^{\prime}l}|\psi_{l}\rangle.$
(7)

We can see that $|\psi_{l}\rangle$ is determined by the coefficients
$\langle\widetilde{n^{\prime}l}|\psi_{l}\rangle$ where $n^{\prime}$ goes only
up to $N$. Therefore, to determine these coefficients we multiply from the
left by $\bra{\widetilde{n^{\prime\prime}l}}$ with $n^{\prime\prime}$ up to
$N$ as well. This results in a matrix equation for the CS coefficients of the
wave function $\underline{\psi}=\langle\widetilde{nl}|\psi_{l}\rangle$

$\underline{\psi}_{l}=\underline{g}_{l}^{C}(Z,E)\underline{\tilde{v}}_{l}^{(s)}\underline{\psi}_{l},$
(8)

where

$\underline{g}_{l}^{C}(Z,E)=\langle\widetilde{nl}|g_{l}^{C}(Z,E)|\widetilde{n^{\prime}l}\rangle.$
(9)

The equation is solvable if the determinant is zero

$|(\underline{g}_{l}^{C}(Z,E))^{-1}-\underline{\tilde{v}}_{l}^{(s)}|=0.$ (10)

The calculation of the matrix $(\underline{g}_{l}^{C})^{-1}$ is based on the
infinite symmetric tridiagonal representation

$\matrixelement{nl;b}{(z-\hat{h}_{l}^{(C)})}{n^{\prime}l;b}=J_{nn^{\prime}}=\begin{cases}\displaystyle{\frac{k^{2}-b^{2}}{2m/\hbar^{2}\,b}}(n+l+1)-Z&\text{for}\quad
n^{\prime}=n\\\
\displaystyle{-\frac{k^{2}+b^{2}}{4m/\hbar^{2}\,b}}\sqrt{(n+1)(n+2l+2)}&\text{for}\quad
n^{\prime}=n+1\\\
\displaystyle{-\frac{k^{2}+b^{2}}{4m/\hbar^{2}\,b}}\sqrt{n(n+2l+1)}&\text{for}\quad
n^{\prime}=n-1\\\ 0&\text{otherwise},\end{cases}$ (11)

where $k=\sqrt{2m/\hbar^{2}\,E}$. It has been shown in Refs. Konya:1997JMP ;
PRADemir2006  that the $N\times N$ matrix $(\underline{g}_{l}^{C})^{-1}$ is
identical to the $N\times N$ matrix $\underline{J}$ plus a correction term in
the bottom-right matrix element

$[\underline{g}_{l}^{C}(Z,E)]^{-1}=\underline{J}_{l}^{C}-\delta_{i,N}\delta_{j,N}J_{N,N+1}C_{N+1}J_{N+1,N}.$
(12)

This correction is given in terms of ${}_{2}F_{1}$ hypergeometric functions

$C_{N}=-\frac{4m/\hbar^{2}\,b}{(b-ik)^{2}(N+l+i\gamma)}\frac{\mbox{}_{2}F_{1}(-l+i\gamma,N+1;N+l+2+i\gamma;(b+ik)^{2}/(b-ik)^{2})}{\mbox{}_{2}F_{1}(-l+i\gamma,N;N+l+1+i\gamma;(b+ik)^{2}/(b-ik)^{2})},$
(13)

where $\gamma=Zm/k$. The ratio of hypergeometric functions with this
combination of indexes can be evaluated by a continued fraction (see eq. Ref.
baker1975essentials ).

In this approach the only approximation is the finite-basis representation of
the short-range potential. As $N\to\infty$ the convergence is guaranteed,
although the method is not variational. Only the short-range potential is
approximated, not the whole Hamiltonian. As a result, the convergence to the
energy is not from above like in an usual Hilbert-space basis approximation
scheme. If the parameter $b$ matches the size of the potential, the
convergence is very fast.

The evaluation of the energy dependent $\underline{g}_{l}^{C}(Z,E)$ is exact
and analytic, so the method can readily be extended to complex energies to
calculate resonances, as has been shown in Refs. papp1987bound ; Konya:1997JMP
.

In this method the solution is not a linear combination of basis functions.
Rather, as Eq. (7) shows

$\innerproduct{r}{\psi_{l}}=\sum_{n}^{N}c_{n}\,\matrixelement*{r}{g_{l}^{C}(Z,E)}{\widetilde{nl}},$
(14)

where
$c_{n}=\sum_{n^{\prime}}\underline{\tilde{v}}^{(s)}_{l,nn^{\prime}}\innerproduct*{\widetilde{n^{\prime}l}}{\psi_{l}}$.
The Green’s function in configuration space is given by the regular and
irregular Coulomb functions $\varphi_{l}^{C}$ and $f_{l}^{C}$, respectively.
Therefore

$\matrixelement*{r}{g_{l}^{C}(Z,E)}{\widetilde{nl}}\sim\int_{0}^{\infty}\differential{r^{\prime}}\varphi_{l}^{C}(r_{<})f_{l}^{C}(r_{>})\innerproduct*{r^{\prime}}{\widetilde{nl}},$
(15)

where $r_{<}=\min(r,r^{\prime})$ and $r_{<}=\max(r,r^{\prime})$. This integral
behaves like $\phi_{l}^{C}(r)$ as $r\to 0$ and like $f_{l}^{C}(r)$ as
$r\to\infty$, which is the correct Coulomb-like asymptotic behavior
PhysRevA.46.4437 , irrespectively of the basis parameter $b$ and the energy
$E$. It should be noted however, that usually we don’t need the wave function,
we need matrix elements of operators representing physical quantities. If an
observable is represented by operator $O$, we can approximate it on Hilbert
space basis like in Eq. (6)

$O\approx\sum_{nn^{\prime}}^{N}|\widetilde{nl}\rangle\underline{\tilde{O}}^{(s)}_{nn^{\prime}}\langle\widetilde{n^{\prime}l}|.$
(16)

Then, the expectation value between eigenstates reads

$\expectationvalue{O}{\psi}\approx\sum^{N}_{nmm^{\prime}n^{\prime}}c_{n}c_{n^{\prime}}\underline{g}^{C}_{nm}\underline{\tilde{O}}_{mm^{\prime}}\underline{g}^{C}_{m^{\prime}n^{\prime}}.$
(17)

##  3 Extension to the Klein-Gordon equation

The relativistic spin-$0$ Klein-Gordon equation with potential term associated
with the energy is given by

$\left(E-V\right)^{2}\psi={{p}}^{2}c^{2}\psi+{m^{2}c^{4}}\psi.$ (18)

If we divide by $2mc^{2}$ we obtain

$\left[\frac{1}{2m}{{p}}^{2}+\frac{mc^{2}}{2}-\frac{1}{2mc^{2}}(E^{2}-2EV+V^{2})\right]\psi=0.$
(19)

By introducing the effective energy

$\varepsilon=\frac{1}{2mc^{2}}(E^{2}-m^{2}c^{4})$ (20)

and effective potential

$\tilde{V}=V\left(\frac{E}{mc^{2}}-\frac{V}{2mc^{2}}\right),$ (21)

we can write the Klein-Gordon equation in a more familiar form

$\tilde{H}\psi=\varepsilon\psi,$ (22)

where

$\tilde{H}=\frac{1}{2m}{{p}}^{2}+\tilde{V}.$ (23)

If we separate off the rest energy, $E=mc^{2}+E^{\prime}$, we find

$\varepsilon=E^{\prime}\left(1+\frac{E^{\prime}}{2mc^{2}}\right)$ (24)

and

$\tilde{V}=V\left(1+\frac{E^{\prime}}{mc^{2}}-\frac{V}{2mc^{2}}\right).$ (25)

We can see that in the non-relativistic limit $\epsilon\sim E^{\prime}$ and
$\tilde{V}\sim V$.

For spherical potentials, the effective Hamiltonian commutes with the angular
momentum operators. Thus, $\\{\tilde{H},{L}^{2},{L}_{z}\\}$ form a complete
set of commuting observables. Assuming a Coulomb plus short range potential
again

$V(r)=Z/r+v_{l}^{(s)}(r),$ (26)

the effective Hamiltonian becomes

$\tilde{H}=-\frac{\hbar^{2}}{2m}\partialderivative[2]{r}+\frac{\hbar^{2}l(l+1)}{2mr^{2}}+\left(\frac{Z}{r}+v_{l}^{(s)}\right)\left(1+\frac{E^{\prime}}{mc^{2}}-\frac{1}{2mc^{2}}\left(\frac{Z}{r}+v_{l}^{(s)}\right)\right).$
(27)

Reorganizing, we find

$\tilde{H}=-\frac{\hbar^{2}}{2m}\partialderivative[2]{r}+\frac{\hbar^{2}}{2m}\frac{l(l+1)-Z^{2}\alpha^{2}}{r^{2}}+\frac{Z^{\prime}}{r}+w_{l}^{(s)},$
(28)

where $\alpha=e^{2}/(\hbar c)$ is the fine structure constant with $e$
electric charge,

$Z^{\prime}=Z\left(1+\frac{E^{\prime}}{mc^{2}}\right)$ (29)

and

$w_{l}^{(s)}=v_{l}^{(s)}\left(1+\frac{E^{\prime}}{mc^{2}}-\frac{1}{mc^{2}}\frac{Z}{r}-\frac{1}{2mc^{2}}v_{l}^{(s)}\right).$
(30)

If we equate

$l(l+1)-Z^{2}\alpha^{2}=\lambda(\lambda+1),$ (31)

we find that

$\lambda=-1/2+\sqrt{(l+1/2)^{2}-Z^{2}\alpha^{2}}.$ (32)

Then Eq. (28) becomes

$\left[-\frac{\hbar^{2}}{2m}\partialderivative[2]{r}+\frac{\hbar^{2}\lambda(\lambda+1)}{2mr^{2}}+\frac{Z^{\prime}}{r}+w_{l}^{(s)}\right]\psi_{l}(r)=\varepsilon\psi_{l}(r).$
(33)

This equation looks just like an ordinary radial Schrödinger equation with
correspondences $l\to\lambda$, $Z\to Z^{\prime}$ and $v_{l}^{(s)}\to
w_{l}^{(s)}$. Then, the corresponding Lippmann-Schwinger equation reads

$|\psi_{l}\rangle=g_{\lambda}^{C}(Z^{\prime},\varepsilon)w_{l}^{(s)}|\psi_{l}\rangle~{}.$
(34)

To determine the solution we perform an analytic continuation in the CS basis
$|nl\rangle\to|n\lambda\rangle$ and the determinant equation becomes

$|(\underline{g}_{\lambda}^{C}(Z^{\prime},\varepsilon))^{-1}-\underline{\tilde{w}}_{l}^{(s)}|=0,$
(35)

where

$\underline{g}_{\lambda}^{C}(Z^{\prime},\varepsilon)=\langle\widetilde{n\lambda}|g_{\lambda}^{C}(Z^{\prime},\varepsilon)|\widetilde{n^{\prime}\lambda}\rangle$
(36)

and

$\underline{\tilde{w}}^{(s)}=\langle{n\lambda}|w_{l}^{(s)}|n^{\prime}\lambda\rangle.$
(37)

We can see in Eq. (5) that the analytic continuation $l\to\lambda$ does not
pose any technical problem. The situation is the same with the Green’s matrix.
Both the matrix elements of $J$ and the ${}_{2}F_{1}$ are analytic in terms of
variables, so the analytic continuation amounts of straightforward
substitutions $l\to\lambda$ and $Z\to Z^{\prime}$ konya1999coulomb .

##  4 Extension to the Dirac equation

The Dirac equation is a first order differential equation for the four-
component wave function. Feynman and Gell-Mann ”squared” it and obtained a
second order differential equation for a two-component wave function
feynman1958theory . If we assume that the potential term is associated only
with the energy we have

$\left(\nabla^{2}+\frac{1}{\hbar^{2}c^{2}}(E-V)^{2}-\frac{m^{2}c^{2}}{\hbar^{2}}\right)\psi+\frac{i}{\hbar
c}\gradient V\cdot\mathrm{\sigma}\psi=0,$ (38)

where $\mathrm{\sigma}$ denotes the Pauli matrices. If we separate off the
rest energy, $E=E^{\prime}+mc^{2}$, we get

$\left[-\frac{\hbar^{2}}{2m}\nabla^{2}+\tilde{V}\right]\psi-\frac{i\hbar
c}{2mc^{2}}\gradient V\cdot\mathrm{\sigma}\psi=\varepsilon\psi.$ (39)

Assume that the potential is spherical and it is a sum of a Coulomb plus
short-range terms. Now the Hamiltonian forms a complete set of commuting
observables with the total angular momentum operators $J^{2}$ and $J_{z}$. We
follow the method of Ref. holstein2013topics , with the difference that our
formulae are valid even if the particle is not charged. The total wave
function is a product radial and angular terms

$\psi=\frac{1}{r}\psi_{j}^{(\pm)}(r)\Phi^{(\pm)}_{j,m}(\theta,\phi)~{},$ (40)

where

$\Phi^{(\pm)}_{j,m}(\theta,\phi)=\sum_{m_{l},m_{s}}\langle
l_{\pm},1/2;m_{l},m_{s}|j,m\rangle Y_{l_{\pm}}(\theta,\phi)\chi_{1/2,m_{s}}$
(41)

are constructed by coupling orbital angular momentum $l_{\pm}$ to the spin
such that $j=l_{+}+1/2=l_{-}-1/2$. Then, from Eq. (39) we obtain the
Hamiltonian

$\tilde{H}=-\frac{\hbar^{2}}{2m}\partialderivative[2]{r}+\frac{\hbar^{2}}{2m}\frac{l_{\pm}(l_{\pm}+1)-Z^{2}\alpha^{2}}{r^{2}}-\frac{\hbar^{2}}{2m}\frac{iZ\alpha}{r^{2}}(\hat{\mathrm{r}}\cdot\mathrm{\sigma})+\frac{Z^{\prime}}{r}+w_{l}^{(s)}+w^{\prime}_{l}\left(\hat{\mathrm{r}}\cdot\mathrm{\sigma}\right)$
(42)

where

$w^{\prime}_{l}(r)=-\frac{i\hbar c}{2mc^{2}}\derivative{v_{l}^{(s)}(r)}{r}.$
(43)

We should recall that the parity operator ${\cal P}$, the mirroring of the
coordinates, in polar coordinates, entails the transformation
$\theta\to\pi-\theta$ and $\phi\to\phi+\pi$. The spherical harmonics transform
as ${\cal P}\>Y_{lm}=(-)^{l}Y_{lm}$ and the electron has positive intrinsic
parity. Consequently

${\cal P}\>\Phi^{(\pm)}_{jm}=(-)^{l_{\pm}}\Phi^{(\pm)}_{jm},$ (44)

i.e. the states $\Phi^{(+)}_{jm}$ and $\Phi^{(-)}_{jm}$ have opposite
parities. We can also see that

$\hat{\mathrm{r}}\cdot\mathrm{\sigma}=\sigma_{x}\sin\theta\cos\phi+\sigma_{y}\sin\theta\sin\phi+\sigma_{z}\cos\theta=\begin{pmatrix}\cos\theta&e^{-i\phi}\sin\theta\\\
e^{i\phi}\sin\theta&-\cos\theta\end{pmatrix}$ (45)

is an odd operator under parity, i.e. ${\cal
P}\hat{\mathrm{r}}\cdot\mathrm{\sigma}=-\hat{\mathrm{r}}\cdot\mathrm{\sigma}$,
and also $(\hat{\mathrm{r}}\cdot\mathrm{\sigma})^{2}=1$. Additionally, we can
easily verify by explicitly calculating the commutator that it commutes with
${\mathrm{J}}$ angular momentum operator

$\commutator{{\mathrm{J}}}{\hat{\mathrm{r}}\cdot\mathrm{\sigma}}=0.$ (46)

So, $\hat{\mathrm{r}}\cdot\mathrm{\sigma}$ acting on
$\ket{\Phi_{j,m}^{(\pm)}}$ does not change the eigenvalue $j$ . On the other
hand, $\hat{\mathrm{r}}\cdot\mathrm{\sigma}$ is an odd operator whose square
is a unit operator. Its action on $\ket{\phi_{j,m}^{(\pm)}}$ should result in
a state with opposite parity, i.e. $\hat{\mathrm{r}}\cdot\mathrm{\sigma}$
should transform the states $\Phi^{(\pm)}_{j,m}$ into each other. So,
$(\hat{\mathrm{r}}\cdot\mathrm{\sigma})\,\ket{\Phi_{j,m}^{(\pm)}}=\ket{\Phi_{j,m}^{(\mp)}}$,
consequently
$\matrixelement{\Phi^{(\pm)}_{j,m}}{\hat{\mathrm{r}}\cdot\mathrm{\sigma}}{\Phi^{(\pm)}_{j,m}}=0$
and
$\matrixelement{\Phi^{(\mp)}_{j,m}}{\hat{\mathrm{r}}\cdot\mathrm{\sigma}}{\Phi^{(\pm)}_{j,m}}=1$.

Eq. (42) is a set of two-component coupled equations. The terms that are
proportional to $\hbar^{2}/(2mr^{2})$ are given by

$\begin{split}&\matrixelement{\Phi_{j,m}^{(\pm)}}{l_{\pm}(l_{\pm}+1)-Z^{2}\alpha^{2}-iZ\alpha\,{\mathrm{\sigma}}\cdot\hat{r}}{\Phi_{j,m}^{(\pm)}}\\\
&=\begin{pmatrix}(j-1/2)(j+1/2)-Z^{2}\alpha^{2}&-iZ\alpha\\\
-iZ\alpha&(j+1/2)(j+3/2)-Z^{2}\alpha^{2}\end{pmatrix}.\end{split}$ (47)

By solving the matrix eigenvalue problem, we find the eigenstates
$\ket{\eta^{(\pm)}}$, which are linear combinations of
$\ket{\Phi_{j,m}^{(\pm)}}$

$\matrixquantity(\eta^{(+)}\\\
\eta^{(-)})=\frac{1}{2}\sqrt{\frac{j+1/2+s}{j+1/2}}\matrixquantity(1&iZ\alpha/(j+1/2+s)\\\
-iZ\alpha/(j+1/2+s)&1)\matrixquantity(\Phi^{(+)}\\\ \Phi^{(-)}),$ (48)

where

$s=\sqrt{(j+1/2)^{2}-Z^{2}\alpha^{2}}.$ (49)

We can equate the eigenvalues by $\lambda_{\pm}(\lambda_{\pm}+1)$ and find
that

$\lambda_{\pm}=s-1/2\mp 1/2.$ (50)

So, Eq. (47) in the $\ket{\eta^{(\pm)}}$ basis becomes diagonal

$\matrixelement{\chi_{j,m}^{(\pm)}}{\hat{L}^{2}/\hbar^{2}-Z^{2}\alpha^{2}-iZ\alpha\,{\mathrm{\sigma}}\cdot\hat{r}}{\chi_{j,m}^{(\pm)}}=\begin{pmatrix}\lambda_{+}(\lambda_{+}+1)&0\\\
0&\lambda_{-}(\lambda_{-}+1)\end{pmatrix}.$ (51)

Consequently, for Eq. (42) we obtain

$\left[-\frac{\hbar^{2}}{2m}\partialderivative[2]{r}+\frac{\hbar^{2}}{2m}\frac{\lambda_{\pm}(\lambda_{\pm}+1)}{r^{2}}+\frac{Z^{\prime}}{r}+w_{l}^{(s)}+w^{\prime}_{l}\matrixquantity(0&1\\\
1&0)\right]\matrixquantity(\psi_{j}^{(+)}\\\
\psi_{j}^{(-)})=\varepsilon\matrixquantity(\psi_{j}^{(+)}\\\ \psi_{j}^{(-)}).$
(52)

We can turn this differential equation into a Lippmann-Schwinger form

$\matrixquantity(\psi_{j}^{(+)}\\\
\psi_{j}^{(-)})=\matrixquantity(g^{C}_{\lambda_{+}}(Z^{\prime},\varepsilon)&0\\\
0&g^{C}_{\lambda_{-}}(Z^{\prime},\varepsilon))\matrixquantity(\underline{\tilde{w}}_{\lambda_{+},\lambda_{+}}&\underline{\tilde{w}}_{\lambda_{+},\lambda_{-}}\\\
\underline{\tilde{w}}_{\lambda_{-},\lambda_{+}}&\underline{\tilde{w}}_{\lambda_{-},\lambda_{-}})\matrixquantity(\psi_{j}^{(+)}\\\
\psi_{j}^{(-)}),$ (53)

where
$\underline{\tilde{w}}_{\lambda_{\pm},\lambda_{\pm}}=\langle{n\lambda_{\pm}}|w_{l}^{(s)}|n^{\prime}\lambda_{\pm}\rangle$
and
$\underline{\tilde{w}}_{\lambda_{\pm},\lambda_{\mp}}=\langle{n\lambda_{\pm}}|w^{\prime}_{l}|n^{\prime}\lambda_{\mp}\rangle$.
Then, the energy $\varepsilon$ can be determined by the zeros of the
determinant

$\matrixquantity|(g^{C}_{\lambda_{+}}(Z^{\prime},\varepsilon))^{-1}-\underline{\tilde{w}}_{\lambda_{+},\lambda_{+}}&-\underline{\tilde{w}}_{\lambda_{+},\lambda_{-}}\\\
-\underline{\tilde{w}}_{\lambda_{-},\lambda_{+}}&(g^{C}_{\lambda_{-}}(Z^{\prime},\varepsilon))^{-1}-\underline{\tilde{w}}_{\lambda_{-},\lambda_{-}}|=0.$
(54)

The corresponding matrix elements can be calculated the same way as before, by
performing an analytic continuation in the non-relativistic formulae
$l\to\lambda_{\pm}$ and $Z\to Z^{\prime}$.

##  5 Numerical Illustrations

As numerical illustrations we consider the model with $m=1$, $\hbar=1$,
$e^{2}=1$ and $\alpha=e^{2}/\hbar c=1/137.03604$. We take $Z=50$ and

$v_{l}^{(s)}(r)=-240\exp(-r)/r+320\exp(-4r)/r.$ (55)

The Schrödinger, the Klein-Gordon and Dirac bound and resonant state results
for $l=0,1,2$ are given in Table 1. The complex energies are given by
$E=E_{r}-i\Gamma/2$, where $E_{r}$ is the resonance energy and $\Gamma/2$ is
the lifetime of the resonant state.

Table 1: Non-relativistic energies.  Non-relativistic energies
---
l=0 | l=1 | l=2
-92.264199 | -86.36494 | -75.76312
-54.224609 | -49.69048 | -41.60855
-26.210528 | -22.84595 | -16.91107
-6.5302229 | -4.175213 | -0.091564
6.139886 | 7.600636 | 10.019283
-0.00000002 i |  -0.0000000003 i | -0.000017 i
Klein-Gordon energies
-92.27553 | -86.37913 | -75.78641
-54.25825 | -49.72854 | -41.65660
-26.26132 | -22.89969 | -16.97094
-6.583528 | -4.229067 | -0.146437
6.098560 | 7.561164 | 10.019284
-0.000000003 i |  -0.0000000002 i |  -0.000014 i
Dirac energies
j=1/2 | j=3/2 | j=5/2
-91.73292 | -86.20067 | -75.66290
-86.80452 | -75.92108 | -62.28595
-54.03927 | -49.62861 | -41.59556
-49.89311 | -41.72801 | -31.45510
-26.13929 | -22.84399 | -16.94087
-22.98923 | -17.01021 | -9.587586
-6.503713 | -4.195174 | -0.130783
-4.285842 | -0.169263 | 4.786557
6.1480973 | 7.580558 |  -0.0000000000006 i
-0.00000008 i |  -0.0000000004 i | 9.991357
| 9.972041 | -0.000015 i
| -0.000014 i |

##  6 Summary and Conclusions

In this work we have extended a quantum mechanical approximation method that
has been rather successful in non-relativistic calculations to calculate bound
and resonant states of the relativistic Klein-Gordon and Dirac equations. We
brought the relativistic equations in a form similar to the non-relativistic
Schrödinger equation. We accomplish this by redefining the angular momentum
$l\to\lambda$, the charge $Z\to Z^{\prime}$, the energy $E\to\varepsilon$, the
short-range potential $v^{(s)}\to w^{(s)}$ and the Green’s operator
$g_{l}^{C}(Z,E)\to g_{l}^{C}(Z^{\prime},\varepsilon)$. This way all the
advantages of the method have been retained and transferred to study
relativistic problems.

## References

  * (1) Z. Papp, Journal of Physics A: Mathematical and General 20(1), 153 (1987)
  * (2) Z. Papp, Physical Review C 38(5), 2457 (1988)
  * (3) B. Kónya, G. Lévai, Z. Papp, Phys. Rev. C 61, 034302 (2000)
  * (4) Z. Papp, W. Plessas, Phys. Rev. C 54(1), 50 (1996)
  * (5) Z. Papp, Phys. Rev. C 55, 1080 (1997)
  * (6) Z. Papp, C.Y. Hu, Z.T. Hlousek, B. Kónya, S. Yakovlev, Phys. Rev. A 63(6), 062721 (2001)
  * (7) N.C. Brown, S.E. Grefe, Z. Papp, Physical Review C 88(4), 047001 (2013)
  * (8) B. Kónya, G. Lévai, Z. Papp, J.Math.Phys. 38, 4832 (1997)
  * (9) F. Demir, Z.T. Hlousek, Z. Papp, Phys. Rev. A 74, 014701 (2006)
  * (10) G. Baker Jr, _Essentials of Pade approximants_ (Academic Press, New York, 1975)
  * (11) Z. Papp, Phys. Rev. A 46, 4437 (1992)
  * (12) B. Kónya, Z. Papp, Journal of Mathematical Physics 40(5), 2307 (1999)
  * (13) R.P. Feynman, M. Gell-Mann, Physical Review 109(1), 193 (1958)
  * (14) B.R. Holstein, _Topics in advanced quantum mechanics_ (Courier Corporation, 2013)