# Multi-purpose robot for rehabilitation of small diameter water pipes

Julien Feiguel1, Mouhamed NDiaye1, Pascal Chambaud1, Adrien Chambellan1,
Pierre Blanc1, Steve Bourgeois1,
Lucas Labarussiat1, Clémence Dubois1, Audrey Vigneron1, Thomas Desrez1, Alain
Riwan1, and Caroline Vienne1 1Université Paris Saclay, CEA, LIST, Palaiseau,
France<EMAIL_ADDRESS>
###### Abstract

Rehabilitating cast iron pipes through lining offers several advantages,
including increased durability, reduced water leaks, and minimal disruption.
This approach presents a cost-effective and environmentally-friendly solution
by sealing cracks and joints, extending the pipeline’s lifespan, and reducing
water wastage, all while avoiding the need for trench excavation. However, due
to the relining process, branch connections are sealed and need to be
reestablished. To address the issue of rehabilitating small-diameter water
pipes, we have designed a modular robot capable of traversing and working
within 200-meter-long, 100mm-diameter cast iron pipes. This robot is equipped
with perception functions to detect, locate, and characterize the branch
connections in cast iron pipes and relocate them after lining, as well as
machining functions. A first prototype of this system has been developed and
validated on an 8-meter-long section, in a laboratory environment.

##  I INTRODUCTION

The integrity and efficiency of drinking water distribution systems are
critical to environmental sustainability and economic stability. A significant
challenge facing these infrastructures globally is the deterioration leading
to leaks, with some networks reporting as much as 20% water loss. This not
only results in substantial resource wastage but also increases the financial
burden on utilities and consumers alike. Traditional pipeline rehabilitation
methods, involving extensive trench-digging, disrupt local communities and
landscapes, and incur high costs due to the excavation and displacement of
large volumes of soil. In response, trenchless rehabilitation technologies
have gained prominence due to their less invasive nature and reduced
environmental impact.

One effective trenchless technique involves the insertion of a high-density
polyethylene (HDPE) liner inside existing cast iron pipes, which seals leaks
and restores structural integrity. This method, typically applied to straight
pipeline sections up to 200 meters in length, requires only access pits at
either end of the segment for liner insertion. However, this retubing process
seals off branch connections to customer properties, necessitating subsequent
interventions to re-establish these essential service connections.

To address this issue, we propose a novel robotic solution in collaboration
with SADE, aiming to realize the restoration of branch connections post-
retubing from inside the HDPE tube. The adopted strategy involves a first pass
through the bare cast iron to locate and map the branch connections through
visual perception as well as characterize their state. During this pass, the
branch connections will be re-bored to ensure a shape and surface condition
suitable for the future insertion of a sealing component. A second pass is
made, after relining, to pierce the PE precisely at the branch connection
location, and finally, a third pass will be dedicated to the restoration of
the branch connections. This last operation is not handled by the current
robotic system but will be addressed in future developments.

The large majority of the literature of in-pipe robot development adresses
inspection task [1] or cleaning tasks in large diameters pipe [2]. To our
knowledge, there is no other robotic system capable of carrying out detection,
characterization and machining operations in such narrow pipes. The
contributions that we highlight therefore concern the general design of the
solution and the mechanical and electrical realization of the system as well
as the integration of advanced perception functions by vision and eddy
currents.

##  II Prototype design

###  II-A Rehabilitation concept and design requirements

The water pipes we aim to rehabilitate have a nominal diameter of 100 mm, but
this diameter is reduced to 80 mm after the installation of the HDPE pipe. In
order to incorporate all the necessary movement, perception, and action
functionalities required for the rehabilitation of these water pipes from
within, while adhering to the constraint of an 80 mm diameter, the
cylindrical-shaped robotic system reaches a length of approximately 2 meters.
However, due to the slight angle between different segments of the pipeline,
it is not feasible to create a monolithic robot exceeding 2 meters in length
with a diameter close to the 80 mm limit. To ensure passage through the
junctions, a design with multiple wagon-like units articulated together is
necessary.

The robot is thus composed of various modules, all having an outside diameter
of 75 mm to fit both within the cast-iron and HDPE pipes:

  * •

A traction module equipped with two motors for propulsion and clamping within
the tube, along with two variators.

  * •

Electronic assembly made up of several sub-modules

    * –

an high-speed optical network module for long-distance communication inside
the tube

    * –

an embedded IPC to manage the robot’s actuators and sensors

    * –

an embedded computer to manage the video stream from the cameras

    * –

Eddy current acquisition electronics

    * –

voltage conversion sub-modules for the various systems

  * •

A controller and conversion card module for controlling the motors of the
modules.

  * •

A final interchangeable operational module, which is respectively used to bore
the cast iron (first pass) and to drill the HDPE (second pass).

The mechatronic system is powered by an external power source of 230Vac - 16A.
The control station is located outside the pipeline. An umbilical tether
provides power to the robot and fiber optic communication facilitates digital
data exchange. We ensure the effective coordination of the perception,
locomotion, and action functionalities through ROS2 middleware [3]. This
middleware facilitates the communication by managing data streams while
adhering to a common clock for temporal synchronization of data. Additionally,
it has the ability to abstract the underlying hardware architecture, thus
enabling the system to run on one or multiple computers operating on
heterogeneous operating systems.

###  II-B Locomotion

The traction module is necessary to move the entire system through the 200 m
section of the pipe. It can be separated into four distinct parts (see Fig.
1):

  * •

The vision system, which will be detailed subsequently.

  * •

A pantograph spacing mechanism, which is powered by a motor combined with a
3x120° wedge screw mechanism. This mechanism pushes the parallelograms
upwards, enabling the traction module to remain concentric to the pipes within
a wide range of diameters, starting from 75 mm to 120 mm. It exerts a force of
up to 300 N.

  * •

The traction mechanism itself utilizes tracks to navigate the robot through
the 200 m section. It is driven by a 100 W motoreducer connected to a 45°
helix angle worm screw, which in turn drives three worm wheels at 120° each to
power individual tracks. This configuration allows the robot to generate a
traction force of 370 N at a speed of 10 cm/s or 510 N at 5 cm/s by adjusting
the screw gear ratio.

  * •

The motor controllers that are positioned at the back end of the traction
module.

Figure 1: Locomotion function of the robot (CAD model and picture).

###  II-C Machining

Since the branch connections of the pipe can be situated on the sides or on
the top, the operational module must be adaptable to each scenario.
Consequently, it is linked to a 400° rotating and locking mechanism, ensuring
concentricity within the pipe. This enables scanning of the various holes to
determine their relative center and axis before initiating machining (refer to
Fig. 2). Three translations are necessary to perform the machining operation.
We selected a delta-type robotic kinematics, which provides high stiffness and
workspace within the harsh geometrical constraints.

Figure 2: Machining function of the robot (CAD model and picture). On the
right, the delta mechanism integrates a pen to visualize and validate the
trajectory.

The machining module is composed of three distinct parts:

  * •

The delta mechanism powered by three 8W motors connected to ballscrews and
linear guides to maintain a sufficient rigidity while machining provides up to
350N/ screw in continuous mode. The mechanism is totally backlash free to
maintain the machining precision as high as possible.

  * •

The machining head employs an 80 W motor connected to a spur gear to maximize
torque while maintaining a suitable machining speed. Tools can be interchanged
using a 6 mm lathe collet. For PE, a larger diameter tool is utilized to
expand the range of hole diameters the machine can handle. Depending on the
material, the machining head can sustain speeds of up to 8000 rpm and torque
of 0.5 Nm, hence cutting depths are kept shallow to minimize radial loads.

  * •

Sensors and electronics to detect precisely the center of the holes using
vision, laser and Eddy current probe.

###  II-D Sensors

While navigating in the cast pipe, visual perception is used to determine the
location of the robot through a simultaneous localization and mapping (SLAM)
approach [4] and to build a map of the branch connections. To carry out these
two tasks, a front vision system is made up of three elements (see Fig. 5):

  * •

A laser allowing detection of branch connections

  * •

A pair of cameras with a wide field of view and depth of field adapted to the
use case, along with cable connectors, allowing both the electronic board to
be separated from the cameras and providing freedom to choose the parallax.

  * •

Lighting for illuminating the scene so that the cameras can observe the
texture of the pipes, allowing the SLAM algorithm to utilize this information
to estimate the 3D movement of the robot

Figure 3: Integration of the front vision system: CAD model (left) and picture
(right).

In addition to this front vision system, an additional back vision system is
integrated in order to characterize the state of the branch connections. Its
objective is twofold: first, to identify non-compliance cases that do not
allow for a branch repair from the inside (significant offset between the
branch axis and the valve), and second, to accurately locate the center of
this branch to enable precise positioning of the drill for machining. The
proposed solution leverages the robot’s ability to rotate around its axis to
perform a 3D reconstruction through profilometry, i.e., by projecting a laser
stripe whose deformation observed by a camera allows deduction of the 3D
shape.

Once the internal lining step has been completed, the robot can no longer use
the vision modality to detect the holes. A non-destructive testing modality
allowing detection through HDPE has therefore been proposed. This measurement
system, based on Eddy Currents, includes two axial coils dedicated to locate
the connections in the axis of the pipeline through differential measurement
and an additional coil dedicated to identifying the center of the connection.
More details are given in the following section.

##  III Methodology and results

In this section, we specifically describe the methodologies of perception
through front vision and through eddy currents implemented to locate the
branches during the advancement of the robot.

###  III-A Detection and localization of branch connections through vision

As mentioned earlier, the front vision system is used both for robot
localization using SLAM and for branch detection and location based on laser
input. To prevent interference between these tasks, we propose projecting
light in different spectral domains so that each light projection does not
affect the same color channels of the camera. To determine if this option of
spectral separation is functional, we observed the spectral response of a cast
iron pipe for the different spectral bands of a color camera. To do this, we
illuminated the pipe with a perfectly white light source and observed the
light intensities for each of the three channels (red, green, and blue) of the
camera. As illustrated in Figure 4, it appears that the spectral response of
the cast iron pipe mainly covers the red and green domains, whereas the
response in the blue domain is significantly weaker. The lack of contrast in
the blue channel compared to the red and green channels indicates that it will
have a negligible contribution to the detection and characterization of visual
landmarks used by the SLAM. The similarity of the red and green channels
indicates that they provide similar information, and the use of either should
be equivalent for the SLAM algorithm.

|  |  |
---|---|---|---
color | red | green | blue
Figure 4: Observation of a sample of cast iron pipe across all channels and
each channel separately.

Based on this analysis, the chosen approach involves using a red laser for
branch detection and lighting projecting in the green spectral domain for
SLAM.

An initial experiment is conducted to assess the system’s capability to
estimate the distance between consecutive branches. The results presented (see
Tab. I) compare the actual distances (measured by a human operator) with those
derived from the 3D point cloud reconstructed by the SLAM algorithm.

TABLE I: Measured minimal and maximal errors in position estimation through vision on six different pipe sections. |  Min error (mm) |  Min error (%) |  Max error (mm) |  Max error (%)
---|---|---|---|---
S1 (0.23m) |  1 |  0.5% |  4 |  0.19%
S2 (0.308m) |  2 |  0.6% |  5 |  0.16%
S3 (0.265m) |  6 |  2.5% |  6 |  2.5%
S4 (0.63m) |  5 |  0.8% |  14 |  2.2%
S5 (1.635m) |  0 |  0.0% |  33 |  2%
S6 (1.50m) |  8 |  0.5% |  8 |  0.5%

The experiments conducted revealed that the branch connections aligned with
the laser were all successfully detected (see Fig. 5). The maximum observed
error reached a peak of 2.5% of the distance between two branches and reached
an absolute maximum of 3.3 cm on the longest section. Due to the methodology
employed, the detected position primarily corresponds to the intersection of
the rear wall of the branch connection with the laser, rather than its center.
Nonetheless, this location is generally sufficient as it allows for indicating
the position of the branch connection within a few centimeters, although
further refinement of the position with the rear system is required in all
cases.

Figure 5: Example of branch detection. Left: the two images of the lasers (red
channel of each video stream) with the detection of the laser lines (colored
lines) and the detection of the branch connection (green circle). Right: image
of the green channel of the left video stream with features detection using
SLAM approach.

###  III-B Characterization of branch connections through profilometry

After initially positioning the robot roughly using its front vision
capability, the characterization of the branch connection is performed through
a profilometry approach embedded in the final module. The objective is to
reconstruct the point cloud of the branch to confirm its suitability for
robotic intervention and to precisely locate its center. To assess the
performance of this perception task, several reconstructions of the same
branch are conducted by varying the position of the valve relative to the
branch connection and the position of the machining head relative to the
branch. As there is no ground truth regarding the orientation of the branches,
this evaluation relies on assessing the diameter of the branch connection and
the valve. In the current configuration, the time required to complete the
reconstruction is approximately 10 minutes for a 360° scan with an angular
step of 1°. However, it should be noted that a 360° scan is not necessary
since the position of the branch is approximately known during the
characterization step. Therefore, this time can be significantly reduced. The
results obtained are illustrated in Figure 6. The measurements taken indicate
an average reconstruction error of approximately 0.3 mm.

Figure 6: Example of a combination of shifted 3D reconstructions (on the
left), along with a zoom on the part corresponding to the branch (on the
right). An error of 0.35 mm is observed between the most extreme
reconstructions.

###  III-C Eddy current perception

To relocate the branch connections across the HDPE tube, the robot is equipped
with coils and sensors allowing for Eddy Current measurement. The principle of
this method, used notably in Non-Destructive Testing, involves supplying a
coil with a sinusoidal electrical current to create a magnetic field. In the
presence of a conductive material, so-called Eddy currents are induced, and a
reactive magnetic field is created in turn. In our case, the presence of holes
in the cast iron at the branch connections induces a variation in the
electrical impedance of the coil, which is automatically detected and used to
position the robot for machining the HDPE.

The localization process consists of two steps using two distinct probes. The
first probe is an axial probe, constituted of two coils used in a differential
mode, which provides an axial localization of the branch connection center.
Then, a point coil mounted on the last operational module of the robot, is
positioned at the axial position provided by the first probe. A scan is
performed along the axis and circumference of the pipeline to determine the
radial position of the branch connection as well as the precise axial position
since the initial axial detection may suffer from robot positioning errors.

The acquisition of data from the measurement system and the position of the
robot is carried out by a system synchronized to the microsecond. The function
is executed within the expected tolerances on two consecutive branch
connections without impact of the presence or absence of the connection
itself. In Fig. 7, we see in particular on the same hole the passage of the
coils triggering the stop (left) and the stabilization of a marker
representing the position of the radial coil in the center of the tap (right).

Figure 7: View of the tests validating the performance of the branch detection
through Eddy current.

The Figure 8 illustrates two positions a) and b) of the radial coil during
this phase, associating each of the captured photos with a diagram indicating
the corresponding position of the coil. In the diagram, the outer blue circle
represents the cast iron pipe, the interrupted area corresponding to the
branch connection, and the inner black circle represents the HDPE tube. The
white and orange rectangle represents the coil, offset from the hole of the
branch connection in the left configuration and aligned with it in the right
configuration. Finally, the blue arrow represents the rotational movement
performed by the probe.

Figure 8: Photo and diagram of the coil offset (a) and aligned (b) with the
hole in the pipe.

Positions a) and b) are also identified in Figure 9, which represents the
measured signal. This signal (blue curve) corresponds to the modulus
calculated from the real and imaginary parts of the complex impedance of the
coil (in green), filtered by a moving average over 15 values. These real and
imaginary parts are obtained from the demodulation of the digitized signal.
The demodulation frequency corresponds to the excitation frequency. The signal
amplitude (y-axis) is displayed against time (x-axis), expressed in seconds.
The hole position is finally determined from the information about the robot’s
speed of movement within the pipeline.

Figure 9: Display of the raw (green) and filtered (blue) measured EC signal as
a function of the time.

The signal processing from the radial coil ensures the radial centering of
this coil on the hole. For radial measurement, a balancing phase then
amplification of the signal before digitization is carried out in order to
maximize the measurement response.

###  III-D Prototype validation

The different functionalities of the robotic system are tested sequentially in
two 8-meter pipes installed in the laboratory. For the system’s first pass, we
have validated the ability of the system to detect, locate, and characterize a
branch of the cast-iron pipe as well as the ability of the operational module
to bore a preformed hole. The total time for the machining operation is around
30 mn to enlarge a hole from a diameter of 20 mm to 24.4 mm.

The robot is then inserted into the second pipe, which has been relined with a
HDPE tube with an internal diameter of 80 mm. The tested sequence then
involves detecting the branches using the axial probe for coarse pre-
positioning of the robot, and then conducting a fine mapping to precisely
position the drill at the center of the branch. Finally, a PE machining test
is performed (see Fig. 10).

Figure 10: Machining tests are realized inside the pipe on cast iron (left)
and PE (right).

The machining of the PE is performed using a Ø 20mm milling cutter with Z=2.
However, first tests demonstrated that this tool is suboptimal due to the low
number of teeth, which results in a high helix angle. Consequently, there is a
tendency for the cutter to screw and jam into the material if the rotational
speed is not high enough. To mitigate the occurrence of cutter jamming during
drilling, we vary the spindle speed between 2800 and 3500 rpm, instead of the
recommended 1000-1300 rpm for plastics. This higher speed helps counteract
screwing effects caused by tool inertia but may lead to local material
melting, thereby increasing viscosity. However, this phenomenon primarily
occurs during drilling, the operation that demands the most effort. Drilling
remains the most delicate operation. By adjusting the spindle rotation speed
parameters, we can achieve drilling in approximately 9 minutes, followed by
reaming up to a diameter of 23 mm, with a machining time of about 4 minutes.

##  IV Conclusion and Future works

Up until now, we have tested the robot’s ability to detect, characterize and
bore branches. This involves the collaboration of vision functionalities
(branch detection, motion estimation, characterization), locomotion (motor
control for positioning in characterization/machining position, control of
machining head movement for characterization) and the machining itself. After
PE retubing, we have assessed the robot’s ability to locate a branch through
the PE and to machine the PE. This involves the collaboration of Eddy current
functionalities (axial and radial) with locomotion (motor control) and PE
machining. In future development, we will have to assess quantitatively the
performances in terms of precision of positioning, duration of the different
robotic operations, and quality of the machining. In parallel, further
validation is planned through field tests to validate the performance and
robustness of the solution on a 100-meter section of pipeline before and after
retubing.

## ACKNOWLEDGMENT

Part of this work is supported by SADE, our industrial partner in this
project.

## APPENDIX

In this section we display the general overview of the robotic system (see
Fig. 11).

Figure 11: General overview of the proposed system for rehabilitation of small
diameters water pipes.

## References

  * [1] L. A. Mateos and M. Vincze, In-pipe Cleaning Mechanical System for DeWaLoP Robot - Developing Water Loss Prevention in Plastics, Proceedings of the Austrian Robotics Workshop 2013.
  * [2] A. Verma, A. Kaiwart, N. D. Dubey, F. Naseer, S. Pradhan, A review on various types of in-pipe inspection robot, in Proceedings of Materials Today, Vol. 50, Pages 1425-1434, 2022.
  * [3] S. Macenski, A. Soragna, M. Carroll, Z. Ge, Impact of ROS 2 Node Composition in Robotic Systems, IEEE Robotics and Autonomous Letters (RA-L), 2023.
  * [4] K. Yousif, A. Bab-Hadiashar, and R. Hoseinnezhad, An overview to visual odometry and visual slam: Applications to mobile robotics, Intelligent Industrial Systems, vol. 1, no. 4, pp. 289–311, 2015


# Topological pumping induced by spatiotemporal modulation of interaction

Boning Huang 1,2, Yongguan Ke 1,2,4,∗, Wenjie Liu 5, Chaohong Lee 2,3,5,†
1Laboratory of Quantum Engineering and Quantum Metrology, School of Physics
and Astronomy, Sun Yat-Sen University (Zhuhai Campus), Zhuhai 519082, China
2Institute of Quantum Precision Measurement, State Key Laboratory of Radio
Frequency Heterogeneous Integration, Shenzhen University, Shenzhen 518060,
China  3College of Physics and Optoelectronic Engineering, Shenzhen
University, Shenzhen 518060, China  4International Quantum Academy, Shenzhen
518048, China  5Quantum Science Center of Guangdong-Hongkong-Macao Greater Bay
Area (Guangdong), Shenzhen 518045, China<EMAIL_ADDRESS><EMAIL_ADDRESS>
###### Abstract

Particle-particle interaction provides a new degree of freedom to induce novel
topological phenomena. Here, we propose to use spatiotemporal modulation of
interaction to realize topological pumping without single-particle
counterpart. Because the modulation breaks time-reversal symmetry, the
multiparticle energy bands of bound states have none-zero Chern number, and
support topological bound edge states. In a Thouless pump, a bound state that
uniformly occupies a topological energy band can be shifted by integer unit
cells per cycle, consistent with the corresponding Chern number. We can also
realize topological pumping of bound edge state from one end to another. The
entanglement entropy between particles rapidly increases at transition points,
which is related to the spatial spread of a bounded pair. In addition, we
propose to realize hybridized pumping with fractional displacement per cycle
by adding an extra tilt potential to separate topological pumping of the bound
state and Bloch oscillations of single particle. Our work could trigger
further studies of correlated topological phenomena that do not have a single-
particle counterpart.

Keywords: Interaction-induced topological pumping, spatiotemporal modulated
interaction, bound states, Bose-Hubbard model

##  1 Introduction

Topological pumping, topologically protected transport of particles assisted
by the topological structure of spatiotemporal modulation, has attracted
intense interest and attention in recent years [1]. Topological pumping can be
mainly classified as quantized Thouless pumping [2], nonquantized geometrical
pumping [3], and topological pumping of edge states [4, 5]. In Thouless
pumping, the displacement of an initial state uniformly filling a topological
band is related to the nontrivial topological invariant known as the Chern
number [6, 7]. However, in geometrical pumping only a single-momentum state of
the topological band is involved, and the displacement is determined by the
nonquantized Berry phase [3]. According to bulk-edge correspondence,
topological edge states exist in the band gap if the topological invariant of
the below bands is nontrivial, providing a topological pumping channel from
one edge to the other edge. Due to the rapid development of quantum
technologies [8, 9, 10, 11], topological pumping has been experimentally
realized in a variety of quantum simulators [4, 12, 13, 14, 3, 15, 16, 17, 18,
5, 19, 20, 21].

Interactions between particles cause novel correlation and entanglement, which
may introduce rich correlated topological phenomena [22, 23, 24, 25, 26, 27,
28, 29, 30, 31, 32, 33, 34, 35, 36]. Understanding the interaction effect in
topological pumping could open a window to the largely unexplored field.
Although interactions break translational symmetry of individual particles,
interacting systems may still preserve cotranslation sysmetry, that is, the
energy is invariant if all particles as a whole are shifted by multiple unit
cells [23, 24, 37]. With the help of cotranslation symmetry, single-particle
Thouless pumping has been generalized to the multiparticle case [23, 38, 39].
Recently, Thouless pumping of bound state and topologically resonant tunneling
have been experimentally realized in cold atomic systems [40] and
superconducting qubits [41]. In particular, topologically resonant tunneling
describes particles transported one-by-one by unit cells per cycle, which has
no counterpart of single-particle topological pumping. Up to now, interaction-
induced topological pumping relies on spatiotemporal modulation of two more
parameters, such as hopping and onsite potential, or hopping and interaction.
It is still unclear whether multiparticle topological pumping can be induced
by only modulating the interaction. Another fascinating topological phenomenon
without single-particle counterpart is fractional topological pumping [42, 43,
44], in which long-range interaction is generally needed to realize fractional
displacement per cycle. It is highly appealing to search for a novel and
simple way to realize fractional topological pumping.

In this paper, we predict multiparticle topological pumping of the bound bulk
state and bound edge state in a Bose-Hubbard model with modulated interaction.
With an initial state given by a maximally localized multiparticle Wannier
state, the mean position shift in a pumping cycle will be determined by the
Chern number of the filled band. Under open boundary conditions, there are
edge modes and transport from one side to another can be performed. The
evolution of the entanglement entropy between particles is also explored.
Without uniform occupation, the multiparticle Bloch states occupying only a
single momentum correspond to a geomerical pumping with nonquantized charge
pump. In addition, we propose a hybridized pumping with fractional mean
position shift per cycle by combining the multiparticle Thouless pumping and
Bloch oscillations of a single particle. Our work provides a new manner to
realize multiparticle Thouless pumping and hybridized fractional pumping, and
brings a new sight to the interacting topological systems.

##  2 Bose-Hubbard model with modulated interction

Figure 1: Schematics for two bosons in a lattice with spatiotemporal
modulation of interaction. If the two bosons are separated, the interaction
plays no role. If two bosons are at the same site, they can feel
spatiotemporal modulation of onsite potential.

We consider $N$ bosons in a one dimensional lattice with $M$ sites, where the
onsite interactions between bosons are periodically modulated, as depicted in
Fig. 1. The modulation of onsite interaction works on bosons at the same site
but plays no role for a single boson or bosons at different sites. The system
is governed by the following Hamiltonian,

$\hat{H}_{M}=J\sum_{j}(\hat{a}_{j+1}^{{\dagger}}\hat{a}_{j}+h.c.)+\frac{1}{2}\sum_{j}(U_{0}+\delta\cos(2\pi\beta
j+\phi))\hat{n}_{j}(\hat{n}_{j}-1).$ (1)

Here, $\hat{a}_{j}^{{\dagger}}$ ($\hat{a}_{j}$) are the bosonic creation
(annihilation) operators at the $j$th site.
$\hat{n}_{j}=\hat{a}_{j}^{{\dagger}}\hat{a}_{j}$ is the particle number
operator at the $j$th site. $J$ is the nearest neighboring hopping strength,
and $U_{j}=U_{0}+\delta\cos(2\pi\beta j+\phi)$ is the periodically modulated
onsite interaction strength between the particles, where $\beta=p/q$ is the
spatial modulation frequency and $\phi$ can be tuned temporally as
$\phi=\omega t$ with modulation frequency $\omega$ and period $T=2\pi/\omega$.
$\delta$ is the modulation strength of the onsite interaction. The
spatiotemporal modulation of interaction can be readily realized in some
quantum simulators. In optical lattices of cold atoms, it is possible to tune
the onsite interaction independently and dynamically through the Feshbach
resonance [40, 45]. In arrays of superconducting circuits, modulation of the
onsite interaction can be achieved by controlling the anharmonicity of
individual qubits [41].

Since the modulation of interaction only takes effect on bound states, we
focus on the subspace of two-particle bound state and derive an effective
Hamiltonian of bounded pair through perturbation theory [46]. The hopping term

$\hat{H}^{\prime}=J\sum_{j=1}^{M}(\hat{a}^{{\dagger}}_{j+1}\hat{a}_{j}+h.c.)$
(2)

can be treated as a perturbation to the interaction term

$\hat{H}_{0}=\frac{1}{2}\sum_{j=1}^{M}(U_{0}+\delta\cos(2\pi\beta
j+\phi))\hat{n}_{j}(\hat{n}_{j}-1).$ (3)

For convenience, considering two bosons, there are eigenstates $|2\rangle_{j}$
of $\hat{H}_{0}$ with eigenvalues $E_{j}=U_{0}+\delta\cos(2\pi\beta i+\phi)$
forming a subspace $\mathcal{U}$, and eigenstates
$|1\rangle_{j}|1\rangle_{k}(j\neq k)$ with eigenvalues $E_{j,k}=0$ forming a
complement subspace $\mathcal{V}$. Here, $|{2}\rangle_{j}$ are short for Fock
states $|0,...,n_{j}=2,...,0\rangle$ and $|1\rangle_{j}|1\rangle_{k}$ are
short for $|0,...,n_{j}=1,...,n_{k}=1,...,0\rangle$. The projection operators
on $\mathcal{U}$ and $\mathcal{V}$ are respectively defined as

$\displaystyle\hat{P}$ $\displaystyle=\sum_{j}|2\rangle_{j}\langle 2|_{j},$
(4) $\displaystyle\hat{S}$ $\displaystyle=\frac{1}{2}\sum_{j\neq
k}(\frac{1}{E_{j}-E_{j,k}}+\frac{1}{E_{k}-E_{j,k}})|1\rangle_{j}|1\rangle_{k}\langle
1|_{k}\langle 1|_{j}.$

Applying the perturbation theory up to the second order [47, 48], the
effective Hamiltonian of subspace $\mathcal{U}$ is given by

$\hat{H}_{\rm
eff}=\hat{P}\hat{H}\hat{P}+\hat{P}\hat{H}^{\prime}\hat{S}\hat{H}^{\prime}\hat{P}.$
(5)

After a series simplifications, we can obtain the effective Hamiltonian as

$\displaystyle\hat{H}_{\rm eff}$
$\displaystyle=\sum_{j=1}^{M}[U_{0}+\delta\cos(2\pi\beta
j+\phi)]\hat{b}_{j}^{{\dagger}}\hat{b}_{j}$ (6)
$\displaystyle+J^{2}\left(\frac{2}{U_{0}+\delta\cos(2\pi\beta
j+\phi)}+\frac{1}{U_{0}+\delta\cos(2\pi\beta(j+1)+\phi)}+\frac{1}{U_{0}+\delta\cos(2\pi\beta(j-1)+\phi)}\right)\hat{b}_{j}^{{\dagger}}\hat{b}_{j}$
$\displaystyle+\sum_{j=1}^{M}J^{2}\left(\frac{1}{U_{0}+\delta\cos(2\pi\beta
j+\phi)}+\frac{1}{U_{0}+\delta\cos(2\pi\beta(j+1)+\phi)}\right)(\hat{b}_{j+1}^{{\dagger}}\hat{b}_{j}+h.c.).$

where
$\hat{b}_{j}^{{\dagger}}=\frac{1}{\sqrt{2}}\hat{a}_{j}^{{\dagger}}\hat{a}_{j}^{{\dagger}}$
simultaneously creates two bosons at the $j$th site. For
$|U_{0}|>>(|J|,|\delta|)$, small variables can be neglected, the effective
Hamiltonian becomes

$\hat{H}_{\rm eff}=\sum_{j=1}^{M}[U_{0}+\delta\cos(2\pi\beta j+\omega
t)]\hat{b}_{j}^{{\dagger}}\hat{b}_{j}+\frac{2J^{2}}{U_{0}}\sum_{j=1}^{M}(\hat{b}_{j+1}^{{\dagger}}\hat{b}_{j}+h.c.).$
(7)

It turns out to be an AAH model which may support nontrivial topological
invariants.

##  3 Topological pumping in two particle system

Figure 2: Energy spectrum and edge modes. (a) Multiparticle Bloch bands in 3D
view. (b) Energy spectrum as a function of time under the open boundary
condition, and the red dots mark the edge modes. (c) and (d): Second-order
correlation of two edge states with energies indicated by the blue asterisk
and cross in (b), respectively. Parameters are chosen as
$N=2,M=87,U_{0}=30,\delta=5,J=3,\beta=1/3$.

Motivated by the effective Hamiltonians (7), we may realize interaction-
induced topological pumping in two-particle system. To this end, we first
briefly review the theory of multiparticle Thouless pumping [23]. If an
interacting system has cotranslation symmetry, that is, the system is
invariant if all the particles are shifted by unit cell, we can separate the
center-of-mass motion and the relative motion, and find that the center-of-
mass (c.m.) momentum $\kappa$ is a good quantum number. The energies as
functions of c.m. momentum form multiparticle Bloch bands, with the
corresponding multiparticle Bloch states $|\psi_{m}(\kappa,t)\rangle$. If the
system is periodically modulated by breaking the time-reversal symmetry, the
Bloch bands may have nontrivial Chern numbers defined by
$|u_{m}(\kappa,t)\rangle$ in the momentum-time space,

$C_{m}=\frac{1}{2\pi}\int_{0}^{\frac{2\pi}{q}}d\kappa\int_{0}^{T}dti(\langle\partial_{t}u_{m}(\kappa,t)|\partial_{\kappa}u_{m}(\kappa,t)\rangle-\langle\partial_{\kappa}u_{m}(\kappa,t)|\partial_{t}u_{m}(\kappa,t)\rangle),$
(8)

where $T$ is the modulation period and $|u_{m}(\kappa,t)\rangle$ is the cell-
periodic part of the Bloch state $|\psi_{m}(\kappa,t)\rangle$. By preparing a
multiparticle Wannier state centered in the $R$th cell,

$|W_{m}(R,t=0)\rangle=\frac{1}{\sqrt{\beta M}}\sum_{\kappa}e^{i\kappa
qR}|\psi_{m}(\kappa,0)\rangle,$ (9)

the state will adiabatically sweep the occupied multiparticle Bloch band. In
one pumping cycle, the c.m. shift $\Delta P$ of particles is related to the
Chern number of the occupied band, $\Delta P=qC_{m}$. In the following, we
will first study the topological properties of multiparticle Bloch band, and
then show different kinds of topological pumping.

###  3.1 Multiparticle Bloch band and topological invariant

In individual subspaces of c.m. momentum $\kappa$ under periodic boundary
condition, we can obtain eigenvalues by diagonalizing
$\hat{H}(\kappa)|u_{m}(\kappa)\rangle=E_{m}(\kappa)|u_{m}(\kappa)\rangle$,
where $\hat{H}(\kappa)$ can be constructed in the basis with translation
symmetry and momentum $\kappa$. The eigenvalues changing with c.m. momentum
form Bloch bands. In Fig. 2(a), we show multiparticle Bloch bands as time
varies in one period. The parameters are chosen as
$N=2,M=87,U_{0}=30,\delta=5,J=3,\beta=1/3$. The energy bands can be classified
as two types, that is, the lower continuum bands correspond to scattering
states, and the higher three separated bands correspond to bound states. The
two particles at the same sites form a bound state, due to strong repulsive
interaction $U_{0}>>(|J|,|\delta|)$ [49, 50, 48, 51]. Because the scattering-
state bands are highly degenerate, the definition of their Chern numbers is
vague. However, different bound-state bands are isolated, and we can calculate
the Chern numbers of the three bound-state bands, which turn out to be
$C=\\{-1,2,-1\\}$ from top to bottom.

The topological properties of the bound-state band can also be related to
topological edge states [52]. According to bulk-edge correspondence, if the
bulk bands are topologically nontrivial, there should exist topological edge
states in the band gaps under open boundary condition. Fig. 2(b) shows the
multiparticle energy spectrum under open boundary conditions in the $t-E$ view
with the same parameters, and there are edge modes in band gaps marked with
red dots. These edge states are dubbed as topological bound edge states,
because the two particles are correlated and localized at left or right ends;
see the second-order correlation functions $C_{ij}^{(2)}$ in Fig. 2(c) and (d)
respectively, where

$C_{ij}^{(2)}=\langle\psi|\hat{a}_{i}^{{\dagger}}\hat{a}_{j}^{{\dagger}}\hat{a}_{j}\hat{a}_{i}|\psi\rangle.$
(10)

###  3.2 Thouless pumping of bound state

Figure 3: Topological pumping of bound bulk state. (a) The trajectory in the
parameter space ($U_{21}\equiv U_{3i+2}-U_{3i+1},U_{32}\equiv
U_{3i+3}-U_{3i+2}$) in one pumping cycle. (b) Density distribution as a
function of time in a pumping cycle with a bound bulk state $|2\rangle_{42}$
as the initial state. (c) C.m. shift (blue) and instantaneous entanglement
entropy (tangerine) between particles as functions of time in one pumping
cycle. (d) Second-order correlation of instantaneous state at $t=T/2$. (e)
Density distribution as a function of time in a pumping cycle with initial
state $|1\rangle_{42}|1\rangle_{82}$, where the two particles are separated.
The parameters are set as $U_{0}=30$, $J=1$, $\delta=2$, $\omega=0.001$ and
$M=87$.

We numerically simulate the Thouless pumping of the two-particle bound state.
Without loss of generality, we consider the bound state $|2\rangle_{42}$ as an
initial state, which is almost the maximally localized multiparticle Wannier
states [53, 54] of the highest multiparticle Bloch bands with Chern number
$-1$. Note that the middle and lowest bound-state bands with Chern numbers $2$
and $-1$ can also work; see A. The parameters are set as $U_{0}=30$, $J=1$,
$\delta=2$, $\omega=0.001$, and the system size is $87$ sites. Dynamics is
governed by
$|\psi(t)\rangle=\mathcal{T}\exp[-i\int\hat{H}(t)dt]|\psi(0)\rangle$, where
$\mathcal{T}$ is a time ordering operator. Fig. 3(a) shows the trajectory in a
parameter space per pumping cycle, where $U_{21}=U_{3i+2}-U_{3i+1}$,
$U_{32}=U_{3i+3}-U_{3i+2}$, and $i$ can be an integer from $0$ to $M-1$.
Similarly to the previous pumping scheme governed by the Rice-Mele model [14,
13], the trajectory also forms a closed loop that encircles the singular point
$(0,0)$, indicating the nontrivial topology of this pumping scheme. Fig. 3(b)
shows the time evolution of density distribution
$\langle\hat{n}_{j}(t)\rangle=\langle\psi(t)|\hat{n}_{j}|\psi(t)\rangle$ for
the bound state. Despite diffusion, the bound state is unidirectionally
shifted by $0.9897$ unit cells during a pumping cycle from $t=0$ to $T$,
consistent with the corresponding Chern number. The blue solid line in Fig.
3(c) shows the c.m. position shift $\Delta P$ in one cycle, which is defined
as

$\Delta
P=\frac{1}{q}(\langle\psi(t)|\hat{x}|\psi(t)\rangle-\langle\psi(0)|\hat{x}|\psi(0)\rangle),$
(11)

where $\hat{x}=\frac{1}{N}\sum_{j}j\hat{n}_{j}$. During the pumping process,
the two bosons are shifted as a whole because of the strong interaction
between them. Fig. 3(d) shows the second-order correlation function of state
at time $t=T/2$. We can find that the diagonal line of correlation function is
dominated, indicating that the two particles are bound states. To further
illustrate the role of interaction in Thouless pumping, we show the time
evolution of density distribution for the initial state
$|1\rangle_{42}|1\rangle_{82}$ with the same parameters; see Fig. 3(e).
Obviously, there is no charge pumping in the pumping process and the particles
spread over the whole system. This happens because the interaction and its
modulation do not work on the particles that are initially separated.

Another interesting phenomenon is the time evolution of the entanglement
entropy [55] between particles, which is defined as

$S=-\sum_{\nu}\lambda_{\nu}ln\lambda_{\nu}.$ (12)

$\lambda_{\nu}$ are the coefficients when performing singular value
decomposition (SVD) on a two-particle state
$|\psi\rangle=\psi_{i,j}|\psi_{i,j}\rangle$, where $i,j$ represent one of the
two particles at the $i$th site and the other at the $j$th site. Obviously,
due to the exchange symmetry of bosonic particles, swapping positions between
any two particles does not change the amplitude, that is,
$\psi_{i,j}=\psi_{j,i}$. $\lambda_{\nu}$ can be obtained by

$\psi_{ij}=\sum_{\nu}\sqrt{\lambda_{\nu}}\psi_{i}^{\nu}\psi_{j}^{\nu},$ (13)

with $\sum_{\nu}\lambda_{\nu}=1$. The tangerine solid line in Fig. 3(c) shows
the entanglement entropy of instantaneous states during the pumping process.
We find that the overall trend of the entanglement entropy is increasing along
with the diffusion of instantaneous state, and there are peaks near the
transition points. This is because the evolved state at the transition time is
bound pair in superposition of two adjacent lattice sites, which can be viewed
as a joint effect of both particle and spatial entanglements.

###  3.3 Geometrical pumping of bound bloch state

Geometrical pumping is known as nonquantized transport of a single-momentum
state, which is determined by the local geometrical properties of the band
structure. This is in stark contrast to quantization transport in Thouless
pumping which is determined by the global topological properties of the filled
band [3]. The pumped charge $Q$ can be calculated using a local current
operator $\hat{J}_{j}(t)$ [56], where the density of particles flown through a
site $j$ is measured by

$Q^{(j)}(\kappa)=\int_{0}^{T}\langle\psi_{m}(\kappa,t)|\hat{J}_{j}(t)|\psi_{m}(\kappa,t)\rangle
dt$ (14)

in a pumping cycle $\hat{J}_{j}(t)$ is obtained via the continuity equation
$\hat{J}_{j}(t)-\hat{J}_{j-1}(t)=i[\hat{n}_{j},\hat{H}_{M}(t)]$, with the form

$\hat{J}_{j}(t)=iJ(\hat{a}^{{\dagger}}_{j}\hat{a}_{j+1}-h.c.).$ (15)

Here, we apply the period boundary condition. Because a unit cell contains
three sites, the pumped charge is averaged over one unit by

$Q(\kappa)=\frac{1}{q}\left[Q^{(j)}(\kappa)+Q^{(j+1)}(\kappa)+Q^{(j+2)}(\kappa)\right].$
(16)

Here, $Q(\kappa)$ is independent of the lattice site due to the cotranslation
symmetry. $Q(\kappa)$ is generally nonquantized, but the sum of all
$Q(\kappa)$ in the Brillouim zone of a certain band gives the quantized pumped
charge,

$\displaystyle\frac{1}{q}\sum_{j}\int_{0}^{T}\langle
W_{m}(R,t)|\hat{J}_{j}|W_{m}(R,t)\rangle dt$ $\displaystyle=$
$\displaystyle\frac{1}{q\beta
M}\sum_{j,\kappa}\int_{0}^{T}\langle\psi_{m}(\kappa,t)|\hat{J}_{j}|\psi_{m}(\kappa,t)\rangle
dt$ (17) $\displaystyle=$ $\displaystyle\sum_{\kappa}Q(\kappa)=NC_{m}.$

Figure 4: Geometrical pumping of single bound Bloch state. Left panel: Average
pumped charge per cycle as a function of c.m. momentum (red cross) for the
highest (a), middle (c) and lowest (e) bound-state bands. Right panel:
Displacement per cycle as a function of c.m. momentum (red asterisk) for the
highest (b), middle (d) and lowest (f) bound-state bands. The parameters are
set as $U_{0}=30,\delta=2,J=1,\omega=0.001,M=87$ for (a) and (b) and
$U_{0}=30,\delta=5,J=1,\omega=0.0005,M=87$ for (c)-(f).

When the particle number $N$ is equal to $1$, our result returns to the
quantized pumping of non-interacting charge [2]. The quantized pumped charge
is also consistent with the quantized displacement per cycle. Similarly to the
group velocity of a single-particle momentum state [7, 36], the group velocity
of a single multiparticle Bloch state is given by

$v_{m}(\kappa,t)=\frac{\partial E_{m}(\kappa,t)}{\hbar\partial
k}+i(\langle\partial_{t}u_{m}(\kappa,t)|\partial_{\kappa}u_{m}(\kappa,t)\rangle-\langle\partial_{\kappa}u_{m}(\kappa,t)|\partial_{t}u_{m}(\kappa,t)\rangle).$
(18)

Then the displacement of a multiparticle Bloch state per cycle is given by

$\Delta X(\kappa)=\int_{0}^{T}v_{m}(\kappa,t)dt.$ (19)

Average displacement over the Brillouim zone of the occupied band gives the
c.m. shift of the corresponding multiparticle Wannier state

$\frac{q}{2\pi}\int_{0}^{2\pi/q}\Delta X(\kappa)d\kappa=qC_{m}.$ (20)

Here, the first term of $v_{m}(\kappa,t)$ vanishes under the integral of
$\kappa$ because of the periodicity of the energy band structure. Eq. (20) is
exactly the same as the formula for multiparticle Thouless pumping [23].
Compared Eq. (20) with Eq. (17), we can find that both the pumped charge and
the displacement per cycle are given by the Chern number of the filled band.

As shown in Fig. 4(a), under periodic boundary condition, we consider the
initial state as the mutiparticle Bloch state with different $\kappa$ in the
highest bound-state band, the pumped charge $Q$ changes with $\kappa$ in the
Brillouin zone. The sum of pumped charge in the Brillouin zone gives
$\sum_{\kappa}Q(\kappa)=-1.9898$, consistent with Eq. (17). The parameters are
set as $U_{0}=30,\delta=2,J=1,\omega=0.001,M=87$. Fig. 4(b) shows $\Delta X$
as a function of $\kappa$ with the same condition as Fig. 4(a). The average of
$\Delta X$ over the Brillouin zone gives $(1/L)\sum_{\kappa}\Delta
X(\kappa)=-3$, which is consistent with Thouless pumping of the initial state
as a maximally localized Wannier state uniformly filling the highest bound-
state band, where the bound pair is shifted by one unit cell per cycle. For
completeness, Figs. 4 (c) and (d) respectively show the average pumped charge
and displacement per cycle as a function of $\kappa$ with initial states as
multiparticle Bloch states of the middle bound band at $t=0$, and the lowest
bound band for (e) and (f). The parameters are
$U_{0}=30,\delta=5,J=1,\omega=0.0005,M=87$, and there is an initial phase
$\phi_{0}=\pi/5$, which is the same as the topological pumping of the bulk
states in A. In a pumping cycle, we obtain $\sum_{\kappa}Q(\kappa)=3.9483$ and
$(1/L)\sum_{\kappa}\Delta X(\kappa)=6$ for (c) and (d) respectively;
$\sum_{\kappa}Q(\kappa)=-1.9909$ and $(1/L)\sum_{\kappa}\Delta X(\kappa)=-3$
for (e) and (f), respectively. We find that the profiles of pumped charge and
the displacement in geometrical pumping are almost the same, and only differ
by a factor related to the particle number and size of a unit cell.

To save computing resource, the numerical simulation of geometrical pumping is
performed in a subspace in which the Hamiltonian is expanded only by the basis
$\mathcal{V}=\\{|2\rangle_{j}\\}$ and
$\mathcal{U}=\\{|1\rangle_{j}|1\rangle_{j+1}\\}$. The reason is that
considering the strong onsite interaction, the instantaneous states during the
adiabatic pumping with a two-particle bound state as initial state dominantly
occupy the subspace $\mathcal{V}$, besides $|2\rangle_{j}$ and
$|2\rangle_{j+1}$ are mainly coupled through $|1\rangle_{j}|1\rangle_{j+1}$.
The validity of this method is numerically verified; see B.

###  3.4 Topological pumping of bound edge state

Figure 5: Topological pumping of bound edge state. (a) Energy spectrum of
bound states as a function of time under open boundary condition, where colors
mark the projection probability of instantaneous state onto eigenstates
$P=|\langle\psi_{m}(t)|\psi(t)\rangle|^{2}$. The initial state is the edge
mode below the lowest bound-state band at $t=0$. (b) Density distribution as a
function of time. (c) C.m. shift (blue) and instantaneous entanglement entropy
(tangerine) between particles as functions of time. (d) Second-order
correlation of instantaneous state at $t=T/3$. The parameters are $U_{0}=30$,
$J=3$, $\delta=2$, $\omega=0.00001$, $M=27$.

We have already shown the topological bound edge states in the band gaps of
bound-state bands in previous section. The energies of the bound edge states
connect different bulk bands, and the left bound edge states can be
continuously transferred to the right bound edge states through bulk states.
In this section, we numerically simulate topological pumping of the bound edge
state [57, 27, 41]. Without loss of generality, we consider the edge-to-edge
transport channel near the lowest bound-state band from $t=0$ to $2T/3$, as
shown in Figs. 5(a) and (b). The initial state is the eigenstate below the
lowest bound-state band, which can be approximately two bosons at the leftmost
site. The parameters are set as $U_{0}=30$, $J=3$, $\delta=2$,
$\omega=0.00001$, and $M=27$. During the pumping process, the initial left
bound edge state is adiabatically transferred to the right side, and the
transition occurs at time $t=1/3T$. Besides, the entanglement entropy of
instantaneous state peaks at the transition time [Fig. 5(c)], where the
diagonal correlation function has broadest spatial distribution [Fig. 5(d)].
The entanglement entropy is nearly the same at the initial and final time
because of the same divergence of spatial distribution.

##  4 Hybridized pumping in three particle system

Figure 6: Hybridized fractional pumping and its breakdown. Top: time evolution
of (a) density distribution and (b) displacement for the initial state
$|1\rangle_{13}|2\rangle_{30}$. Bottom: similar to the top, but the initial
state is changed to $|1\rangle_{15}|2\rangle_{21}$. The white dashed line
denotes the maximal spread width of the single particle with Bloch
oscillations. Parameters are set as
$U_{0}=30,J=1,\omega=0.001,\delta=2,F=0.5,M=33$.

In the previous section, we have shown interaction-induced Thouless pumping of
two-particle bound states. Note that a single particle will only undergo
quantum walks with zero displacement in our system. If we consider that a
bound pair and a single particle cannot meet each other, we can predict that
the overall displacement per particle is $(1\times 2+0)/3=2/3$ in one cycle,
which naturally is a kind of hybridized pumping in few-body systems with
fractional c.m. shift per cycle. Following this similar idea, we can realize
hybridized fractional pumping with any designed fractional values by combining
more spatially-separated bound pairs and individual particles. However, the
only trouble is that a single particle will spread over the entire lattices
without constraint, and the interaction between a single particle and a bound
pair is unavoidable. To overcome this problem, we add an extra term of tilt
potential,

$\hat{H}_{F}=\sum_{j}Fj\hat{n}_{j}.$ (21)

Here, $F=\hbar\omega_{F}$ with $\omega_{F}$ being the frequency of Bloch
oscillations and $\hbar$ is the Planck constant. The frequency of Bloch
oscillations needs to be a multiple of the driven frequency in Thouless
pumping. Due to the tilted potential, a single particle initially at a single
site will undergo Bloch oscillations of breathing modes [58, 59, 60]. The
single particle will return to its initial position in multiple periods of
Thouless pumping. The maximal spatial spread of the single-particle wavepacket
is inversely proportional to the frequency of Bloch oscillations, which is
given by [61]

$L=4J/F.$ (22)

In addition, the weak tilt potential does not affect the quantization
displacement of the bound pair in the Thouless pump. We can choose a proper
frequency of Bloch oscillations, so that both the adiabatic condition and
spatial separation can be satisfied at the same time.

In our simulation, the parameters are set as
$U_{0}=30,J=1,\omega=0.001,\delta=2$, and $F=0.5$, and the size of the lattice
system is $M=33$. The initial state is $|1\rangle_{13}|2\rangle_{30}$. With
adiabatic dynamics
$|\psi(t)\rangle=\mathcal{T}\exp[-i\int\hat{H}(t)dt]|\psi(0)\rangle$, Figs. 6
(a) and (b) show the time evolution of the density distribution
$\langle\hat{n}_{j}\rangle$ and the c.m. shift $\Delta P$ in a pumping cycle,
respectively. As expected, the response of the two-boson bound state is still
a topological pumping with less diffusion, and the behavior of the single one
is Bloch oscillations. The c.m. position shifts $-0.6454$ unit cells overall,
consistent with a $2/3$ fractional shift of displacement, which is attributed
to the combination of one unit cell shift of two-boson bound state and that
the single particle goes back to its initial position in multiple periods of
Bloch oscillations. During the pumping process, there is little overlap
between the density distribution of the two-boson bound state and the other
single particle.

However, if there is spatial overlap between the bound pair and single
particle, the bound pair and single particle will interact with each other.
The interaction will break both the Thouless pumping of bound pair and the
Bloch oscillations of the single particle. Figs. 6(c) and (d) correspond to
the initial state $|1\rangle_{15}|2\rangle_{21}$, where the bound pair and the
single particle are quite close. Because the distance is within the maximal
spread width of the single particle, The overlap between the bound pair and
single particle is large. The parameters are the same as those in Figs. 6(a)
and (b). With large density overlap, the phenomenon becomes complex and
disordered, and the fractioanl value deviates from $2/3$. That is why we have
to add the tilted potential to separate bound pair and single particle.

##  5 Summary and discussion

We have suggested a form of multiparticle Thouless pumping induced by
periodically modulated interaction. With spatiotenporal modulated onsite
interaction in a Bose-Hubbard model, topological pumping of the bound bulk
state given by the maximally localized multiparticle Wannier state and bound
edge state under open boundary condition can be performed. The entanglement
entropy between particles also evolves during the pumping process, with a peak
value at the transition point related to the dispersion of the single-particle
wavepacket. The topological pumping of bound bulk states can be revealed by a
summation of geometrical pumping for different multiparticle Bloch states with
different c.m. momenta. In the next step, it is interesting to explore
unexpected richer topological states and phenomena induced by interactions in
higher dimensions [37, 62, 63, 64, 65]. Besides, we propose a hybridized
pumping with fractional c.m. shift per cycle realized through the combination
of Thouless pumping of multiparticle bound state and Bloch oscillations of a
single particle. It is worth studying nontrivial factional topological states
in few-body systems [66, 67, 68, 69, 70].

##  6 Acknowledgements

The authors acknowledge useful discussions with Ling Lin, Li Zhang, Na Zhang
and Zhoutao Lei. This work was supported by the National Key Research and
Development Program of China (Grant No. 2022YFA1404104), the National Natural
Science Foundation of China (Grants Nos. 12025509, 12275365), the Key-Area
Research and Development Program of Guangdong Province (Grant No.
2019B030330001), and the Natural Science Foundation of Guangdong (Grant No.
2023A1515012099).

##  Appendix A Topological pumping of bound states in the middle and lowest
bound-state bands

Figure 7: Topological pumping of bound bulk states assisted with the middle
and lowest bound-state bands. Density distribution as a function of time in
one pumping cycle with an initial phase $\phi_{0}=\pi/5$. The initial states
are set as $|2\rangle_{41}$ for (a) and $|2\rangle_{40}$ for (c); c.m. shift
and entanglement entropy as a function of time in the pumping cycle for the
initial states $|2\rangle_{41}$ for (b) and $|2\rangle_{40}$ for (d).
Parameters are set as $U_{0}=30,\delta=5,J=1,\omega=0.0005$ and $M=87$.

In Figs. 7(a) and (b), we show the time evolution of density distribution,
displacement, and entanglement entropy of an initial state $|2\rangle_{41}$.
In one pumping cycle, the c.m. position is shifted by $1.97$ unit cells, and
the peaks of entanglement entropy appear at the transition time. Similar
results for the initial state $|2\rangle_{40}$ are shown in Figs. 7(c) and
(d), where the c.m. position of particles is shifted by $-0.9728$ unit cells
in one cycle. The displacements per cycle are consistent with the
corresponding Chern numbers, accompanied by an overall increase of the
entanglement entropy with peaks emerging at the transition points.

##  Appendix B Validity of calculation in subspace

Figure 8: Comparison of energy spectral and density evolution in full Hilbert
space and bound state dominated subspace. The bound-state bands shown in $t-E$
view in full Hilbert space for (a) and the subspace for (b); Density
distribution as a function of time with initial state as $|2\rangle_{42}$
simulated in the full space for (c) and subspace for (d). The parameters are
chosen as $U_{0}=30,\delta=2,J=1,\omega=0.001,M=87$.

We numerically verify the equivalence of simulating the adiabatic evolution of
bound states in the total Hilbert space or in a dominated subspace expanded by
$\mathcal{V}=\\{|2\rangle_{j}\\}$ and
$\mathcal{U}=\\{|1\rangle_{j}|1\rangle_{j+1}\\}$. Figs. 8(a) and (b) show the
multiparticle bound-state band in $t-E$ view for total space and the dominated
subspaces, respectively. Figs. 8(c) and (d) show the density distribution as a
function of time in the topological pumping process with initial state
$|2\rangle_{42}$, which almost uniformly fills the highest bound-state band.
The parameters are the same as the corresponding simulation in the main text
that $U_{0}=30,\delta=2,J=1,\omega=0.001,M=87$. Through the comparison of the
two cases, it is obvious that the method of projecting to a reduced dominated
subspace is valid and is beneficial for numerical calculation in the many-body
region.

## References

## References

  * [1] Citro R and Aidelsburger M 2023 Nature Reviews Physics 5 87–101 URL https://www.nature.com/articles/s42254-022-00545-0
  * [2] Thouless D J 1983 Phys. Rev. B 27(10) 6083–6087 URL https://link.aps.org/doi/10.1103/PhysRevB.27.6083
  * [3] Lu H I, Schemmer M, Aycock L M, Genkina D, Sugawa S and Spielman I B 2016 Phys. Rev. Lett. 116(20) 200402 URL https://link.aps.org/doi/10.1103/PhysRevLett.116.200402
  * [4] Kraus Y E, Lahini Y, Ringel Z, Verbin M and Zilberberg O 2012 Phys. Rev. Lett. 109(10) 106402 URL https://link.aps.org/doi/10.1103/PhysRevLett.109.106402
  * [5] Cheng Q, Wang H, Ke Y, Chen T, Yu Y, Kivshar Y S, Lee C and Pan Y 2022 Nature Communications 13 249 URL https://www.nature.com/articles/s41467-021-27773-9
  * [6] King-Smith R D and Vanderbilt D 1993 Phys. Rev. B 47(3) 1651–1654 URL https://link.aps.org/doi/10.1103/PhysRevB.47.1651
  * [7] Xiao D, Chang M C and Niu Q 2010 Rev. Mod. Phys. 82(3) 1959–2007 URL https://link.aps.org/doi/10.1103/RevModPhys.82.1959
  * [8] Wang L, Troyer M and Dai X 2013 Phys. Rev. Lett. 111(2) 026802 URL https://link.aps.org/doi/10.1103/PhysRevLett.111.026802
  * [9] Mei F, You J B, Zhang D W, Yang X C, Fazio R, Zhu S L and Kwek L C 2014 Phys. Rev. A 90(6) 063638 URL https://link.aps.org/doi/10.1103/PhysRevA.90.063638
  * [10] Wei R and Mueller E J 2015 Phys. Rev. A 92(1) 013609 URL https://link.aps.org/doi/10.1103/PhysRevA.92.013609
  * [11] Ke Y, Qin X, Mei F, Zhong H, Kivshar Y S and Lee C 2016 Laser & Photonics Reviews 10 1064–1064 (Preprint https://onlinelibrary.wiley.com/doi/pdf/10.1002/lpor.201670069) URL https://onlinelibrary.wiley.com/doi/abs/10.1002/lpor.201670069
  * [12] Verbin M, Zilberberg O, Lahini Y, Kraus Y E and Silberberg Y 2015 Phys. Rev. B 91(6) 064201 URL https://link.aps.org/doi/10.1103/PhysRevB.91.064201
  * [13] Lohse M, Schweizer C, Zilberberg O, Aidelsburger M and Bloch I 2016 Nature Physics 12 350–354 URL https://www.nature.com/articles/nphys3584
  * [14] Nakajima S, Tomita T, Taie S, Ichinose T, Ozawa H, Wang L, Troyer M and Takahashi Y 2016 Nature Physics 12 296–300 URL https://www.nature.com/articles/nphys3622
  * [15] Cerjan A, Wang M, Huang S, Chen K P and Rechtsman M C 2020 Light: Science & Applications 9 178 URL https://www.nature.com/articles/s41377-020-00408-2
  * [16] Grinberg I H, Lin M, Harris C, Benalcazar W A, Peterson C W, Hughes T L and Bahl G 2020 Nature communications 11 974 URL https://www.nature.com/articles/s41467-020-14804-0
  * [17] Xia Y, Riva E, Rosa M I N, Cazzulani G, Erturk A, Braghin F and Ruzzene M 2021 Phys. Rev. Lett. 126(9) 095501 URL https://link.aps.org/doi/10.1103/PhysRevLett.126.095501
  * [18] Jürgensen M, Mukherjee S and Rechtsman M C 2021 Nature 596 63–67 URL https://www.nature.com/articles/s41586-021-03688-9
  * [19] Sun Y K, Zhang X L, Yu F, Tian Z N, Chen Q D and Sun H B 2022 Nature Physics 18 1080–1085 URL https://www.nature.com/articles/s41567-022-01669-x
  * [20] You O, Liang S, Xie B, Gao W, Ye W, Zhu J and Zhang S 2022 Phys. Rev. Lett. 128(24) 244302 URL https://link.aps.org/doi/10.1103/PhysRevLett.128.244302
  * [21] Jürgensen M, Mukherjee S, Jörg C and Rechtsman M C 2023 Nature Physics 19 420–426 URL https://www.nature.com/articles/s41567-022-01871-x
  * [22] Gorlach M A and Poddubny A N 2017 Phys. Rev. A 95(5) 053866 URL https://link.aps.org/doi/10.1103/PhysRevA.95.053866
  * [23] Ke Y, Qin X, Kivshar Y S and Lee C 2017 Phys. Rev. A 95(6) 063630 URL https://link.aps.org/doi/10.1103/PhysRevA.95.063630
  * [24] Qin X, Mei F, Ke Y, Zhang L and Lee C 2017 Phys. Rev. B 96(19) 195134 URL https://link.aps.org/doi/10.1103/PhysRevB.96.195134
  * [25] Hayward A, Schweizer C, Lohse M, Aidelsburger M and Heidrich-Meisner F 2018 Phys. Rev. B 98(24) 245148 URL https://link.aps.org/doi/10.1103/PhysRevB.98.245148
  * [26] van Voorden B A and Schoutens K 2019 New Journal of Physics 21 013026 URL https://dx.doi.org/10.1088/1367-2630/aaf748
  * [27] Olekhno N A, Kretov E I, Stepanenko A A, Ivanova P A, Yaroshenko V V, Puhtina E M, Filonov D S, Cappello B, Matekovits L and Gorlach M A 2020 Nature communications 11 1436 URL https://www.nature.com/articles/s41467-020-14994-7
  * [28] Stepanenko A A and Gorlach M A 2020 Phys. Rev. A 102(1) 013510 URL https://link.aps.org/doi/10.1103/PhysRevA.102.013510
  * [29] Kuno Y and Hatsugai Y 2020 Phys. Rev. Res. 2(4) 042024 URL https://link.aps.org/doi/10.1103/PhysRevResearch.2.042024
  * [30] Pelegrí G, Marques A M, Ahufinger V, Mompart J and Dias R G 2020 Phys. Rev. Res. 2(3) 033267 URL https://link.aps.org/doi/10.1103/PhysRevResearch.2.033267
  * [31] Ke Y, Zhong J, Poshakinskiy A V, Kivshar Y S, Poddubny A N and Lee C 2020 Phys. Rev. Res. 2(3) 033190 URL https://link.aps.org/doi/10.1103/PhysRevResearch.2.033190
  * [32] Ke Y, Hu S, Zhu B, Gong J, Kivshar Y and Lee C 2020 Phys. Rev. Res. 2(3) 033143 URL https://link.aps.org/doi/10.1103/PhysRevResearch.2.033143
  * [33] Azcona P M and Downing C 2021 Scientific Reports 11 12540 URL https://www.nature.com/articles/s41598-021-91778-z
  * [34] Zheng Y and Yang S J 2023 Journal of Physics B: Atomic, Molecular and Optical Physics 56 125301 URL https://dx.doi.org/10.1088/1361-6455/acd66f
  * [35] Lin L, Ke Y and Lee C 2023 Phys. Rev. B 107(12) 125161 URL https://link.aps.org/doi/10.1103/PhysRevB.107.125161
  * [36] Liu W, Hu S, Zhang L, Ke Y and Lee C 2023 Phys. Rev. Res. 5(1) 013020 URL https://link.aps.org/doi/10.1103/PhysRevResearch.5.013020
  * [37] Qin X, Mei F, Ke Y, Zhang L and Lee C 2018 New Journal of Physics 20 013003 URL https://dx.doi.org/10.1088/1367-2630/aa9556
  * [38] Lin L, Ke Y and Lee C 2020 Phys. Rev. A 101(2) 023620 URL https://link.aps.org/doi/10.1103/PhysRevA.101.023620
  * [39] Ke Y and Lee C 2023 Nature Physics 19 1387–1388 URL https://www.nature.com/articles/s41567-023-02169-2
  * [40] Walter A S, Zhu Z, Gächter M, Minguzzi J, Roschinski S, Sandholzer K, Viebahn K and Esslinger T 2023 Nature Physics 1–5 URL https://www.nature.com/articles/s41567-023-02145-w
  * [41] Tao Z, Huang W, Niu J, Zhang L, Ke Y, Gu X, Lin L, Qiu J, Sun X, Yang X, Zhang J, Zhang J, Zhao S, Zhou Y, Deng X, Hu C, Hu L, Li J, Liu Y, Tan D, Xu Y, Yan T, Chen Y, Lee C, Zhong Y, Liu S and Yu D 2023 Interaction-induced topological pumping in a solid-state quantum system (Preprint 2303.04582) URL https://arxiv.org/abs/2303.04582
  * [42] Zeng T S, Wang C and Zhai H 2015 Phys. Rev. Lett. 115(9) 095302 URL https://link.aps.org/doi/10.1103/PhysRevLett.115.095302
  * [43] Zeng T S, Zhu W and Sheng D N 2016 Phys. Rev. B 94(23) 235139 URL https://link.aps.org/doi/10.1103/PhysRevB.94.235139
  * [44] Taddia L, Cornfeld E, Rossini D, Mazza L, Sela E and Fazio R 2017 Phys. Rev. Lett. 118(23) 230402 URL https://link.aps.org/doi/10.1103/PhysRevLett.118.230402
  * [45] Meinert F, Mark M J, Lauber K, Daley A J and Nägerl H C 2016 Phys. Rev. Lett. 116(20) 205301 URL https://link.aps.org/doi/10.1103/PhysRevLett.116.205301
  * [46] Bravyi S, DiVincenzo D P and Loss D 2011 Annals of physics 326 2793–2826 URL https://www.sciencedirect.com/science/article/abs/pii/S0003491611001059
  * [47] Takahashi M 1977 Journal of Physics C: Solid State Physics 10 1289 URL https://dx.doi.org/10.1088/0022-3719/10/8/031
  * [48] Qin X, Ke Y, Guan X, Li Z, Andrei N and Lee C 2014 Phys. Rev. A 90(6) 062301 URL https://link.aps.org/doi/10.1103/PhysRevA.90.062301
  * [49] Winkler K, Thalhammer G, Lang F, Grimm R, Hecker Denschlag J, Daley A, Kantian A, Büchler H and Zoller P 2006 Nature 441 853–856 URL https://www.nature.com/articles/nature04918
  * [50] Valiente M and Petrosyan D 2008 Journal of Physics B: Atomic, Molecular and Optical Physics 41 161002 URL https://iopscience.iop.org/article/10.1088/0953-4075/41/16/161002/meta
  * [51] Fukuhara T, Schauß P, Endres M, Hild S, Cheneau M, Bloch I and Gross C 2013 Nature 502 76–79 URL https://www.nature.com/articles/nature12541
  * [52] Gorlach M A and Poddubny A N 2017 Phys. Rev. A 95(5) 053866 URL https://link.aps.org/doi/10.1103/PhysRevA.95.053866
  * [53] Marzari N and Vanderbilt D 1997 Phys. Rev. B 56(20) 12847–12865 URL https://link.aps.org/doi/10.1103/PhysRevB.56.12847
  * [54] Marzari N, Mostofi A A, Yates J R, Souza I and Vanderbilt D 2012 Rev. Mod. Phys. 84(4) 1419–1475 URL https://link.aps.org/doi/10.1103/RevModPhys.84.1419
  * [55] Poshakinskiy A V, Zhong J, Ke Y, Olekhno N A, Lee C, Kivshar Y S and Poddubny A N 2021 npj Quantum Information 7 34 URL https://www.nature.com/articles/s41534-021-00372-8
  * [56] Rossini D, Gibertini M, Giovannetti V and Fazio R 2013 Phys. Rev. B 87(8) 085131 URL https://link.aps.org/doi/10.1103/PhysRevB.87.085131
  * [57] Qi L, Han N, Hu S and He A L 2023 Phys. Rev. A 108(3) 032402 URL https://link.aps.org/doi/10.1103/PhysRevA.108.032402
  * [58] Bloch F 1929 Zeitschrift für physik 52 555–600
  * [59] Feldmann J, Leo K, Shah J, Miller D A B, Cunningham J E, Meier T, von Plessen G, Schulze A, Thomas P and Schmitt-Rink S 1992 Phys. Rev. B 46(11) 7252–7255 URL https://link.aps.org/doi/10.1103/PhysRevB.46.7252
  * [60] Leo K, Bolivar P H, Brüggemann F, Schwedler R and Köhler K 1992 Solid State Communications 84 943–946 URL https://www.sciencedirect.com/science/article/abs/pii/003810989290798E
  * [61] Hartmann T, Keck F, Korsch H J and Mossmann S 2004 New Journal of Physics 6 2 URL https://dx.doi.org/10.1088/1367-2630/6/1/002
  * [62] Rachel S 2018 Reports on Progress in Physics 81 116501 URL https://iopscience.iop.org/article/10.1088/1361-6633/aad6a6/meta
  * [63] Salerno G, Palumbo G, Goldman N and Di Liberto M 2020 Phys. Rev. Res. 2(1) 013348 URL https://link.aps.org/doi/10.1103/PhysRevResearch.2.013348
  * [64] Berti A and Carusotto I 2022 Phys. Rev. A 105(2) 023329 URL https://link.aps.org/doi/10.1103/PhysRevA.105.023329
  * [65] Iskin M 2023 Phys. Rev. A 107(5) 053323 URL https://link.aps.org/doi/10.1103/PhysRevA.107.053323
  * [66] Gemelke N, Sarajlic E and Chu S 2010 Rotating few-body atomic systems in the fractional quantum hall regime (Preprint 1007.2677) URL https://arxiv.org/abs/1007.2677
  * [67] Repellin C, Léonard J and Goldman N 2020 Phys. Rev. A 102(6) 063316 URL https://link.aps.org/doi/10.1103/PhysRevA.102.063316
  * [68] Fabre A, Bouhiron J B, Satoor T, Lopes R and Nascimbene S 2022 Phys. Rev. Lett. 128(17) 173202 URL https://link.aps.org/doi/10.1103/PhysRevLett.128.173202
  * [69] Okuma N and Mizoguchi T 2023 Phys. Rev. Res. 5(1) 013112 URL https://link.aps.org/doi/10.1103/PhysRevResearch.5.013112
  * [70] Léonard J, Kim S, Kwan J, Segura P, Grusdt F, Repellin C, Goldman N and Greiner M 2023 Nature 1–5 URL https://www.nature.com/articles/s41586-023-06122-4#citeas


Determinantal Expressions in Multi-Species TASEP

Determinantal Expressions in Multi-Species TASEP

Jeffrey KUAN

J. Kuan

Texas A&M University, Department of Mathematics,
Mailstop 3368, College Station, TX 77843-3368, USA<EMAIL_ADDRESS>http://math.tamu.edu/~jkuan/

Received July 14, 2020, in final form December 03, 2020; Published online
December 11, 2020

Consider an inhomogeneous multi-species TASEP with drift to the left, and
define a height function which equals the maximum species number to the left
of a lattice site. For each fixed time, the multi-point distributions of these
height functions have a determinantal structure. In the homogeneous case and
for certain initial conditions, the fluctuations of the height function
converge to Gaussian random variables in the large-time limit. The proof
utilizes a coupling between the multi-species TASEP and a coalescing random
walk, and previously known results for coalescing random walks.

determinantal; multi-species; TASEP; coalescing

60J27

##  1 Introduction

Consider a multi-species TASEP (totally asymmetric simple exclusion process)
on the infinite lattice $\mathbb{Z}$. The species numbers are indexed by the
integers, so any particle configuration can be expressed as
$\eta(t)=(\eta_{x}(t))_{x\in\mathbb{Z}}$, where $\eta_{x}(t)$ denotes the
species number of the particle at site $x$. By convention, let
$\eta_{x}(t)=-\infty$ denote a hole (no particle) at lattice site $x$. The
particle at site $x$ jumps one step to the left with exponential jump rates of
rate $\lambda_{x}$, and particles with higher species have priority to jump.
In other words, if a particle of species $j$ attempts to jump to a site of
species $i$, then the particles switch places if $i<j$, and the jump is
blocked if $i>j$. All jumps occur independently of each other. For
$\mu\in\mathbb{Z}\cup\\{\pm\infty\\}$, let $\xi_{\mu}(t)$ denote the location
of a single random walker which starts at $\mu$, and jumps one step to the
left with exponential rate $\lambda_{x}$ when at site $x$. If
$\mu=\\{\pm\infty\\}$ then the particle remains at $\\{\pm\infty\\}$ for all
times $t\geq 0$. In other words, $\xi_{\mu}(t)$ describes the evolution of the
single-species inhomogeneous TASEP with only a single particle starting at
lattice site $\mu$.

Define

$g_{y}(t)=\max\\{\eta_{x}(t)\colon x<y\\}.$

In words, this is the maximum species number among all particles to the left
of $y$. The main result of this paper is the following theorem.

######  Theorem 1.1.

Consider any initial conditions $\eta_{x}(0)$. Fix $k$ lattice points
$x_{1}<\dots<x_{k}$ and $k$ species numbers $m_{1}<\dots<m_{k}$, and let

$\mu_{i}=\inf\\{x\colon\eta_{x}(0)\geq m_{i}\\},$

where $\inf\varnothing=+\infty$ by convention. Then

$\mathbb{P}(g_{x_{1}}(t)\geq m_{1},\dots,g_{x_{k}}(t)\geq
m_{k})=\det[G_{ij}]_{1\leq i,j\leq k},$

where $G$ is the $k\times k$ matrix with entries

$G_{ij}=\mathbb{P}(\xi_{\mu_{i}}(t)<x_{j})-1_{\\{{i<j}\\}}.$

In the homogeneous case when all $\lambda_{x}\equiv 1$, it is readily seen (by
the central limit theorem) that the fluctuations are Gaussian.

######  Corollary 1.2.

Suppose that all $\lambda_{x}$ are equal to $1$. Assume that
$x_{1},\dots,x_{k}$ and $m_{1},\dots,m_{k}$ depend on $t$ in such a way that

$\frac{\mu_{i}-\nu_{i}t}{t^{1/2}}\rightarrow
w_{i},\qquad\frac{x_{j}-\nu_{j}^{\prime}t}{t^{1/2}}\rightarrow
w_{j}^{\prime}.$

for some $\nu_{i},w_{i},\nu_{j}^{\prime},w_{j}^{\prime}\in\mathbb{R}$. Then

$\lim_{t\rightarrow\infty}\mathbb{P}(g_{x_{1}}(t)\geq
m_{1},\dots,g_{x_{k}}(t)\geq m_{k})=\det[M_{ij}]_{1\leq i,j\leq k},$

where $M$ has entries

$M_{ij}=\Phi(w_{i}-w_{j}^{\prime})-1_{\\{i<j\\}}.$

Here, $\Phi$ is the cumulative distribution of a standard Gaussian

$\Phi(z)=\int_{-\infty}^{z}\frac{1}{\sqrt{2\pi}}{\rm e}^{-y^{2}/(2)}{\rm d}y.$

######  Remark 1.3.

If $m_{i}\geq m_{j}$ for some $i<j$ in Theorem 1.1, then

$\\{g_{x_{i}}(t)\geq m_{i}\\}\subseteq\\{g_{x_{j}}(t)\geq m_{j}\\},$

so there is no loss of generality in assuming that $m_{1}<\dots<m_{k}$. In
other words, if $m_{i}\geq m_{j}$ and $i<j$, there is the equality

$\mathbb{P}\bigg{(}\bigcap_{a\in\\{1,\dots,k\\}}\\{g_{x_{a}}(t)\geq
m_{a}\\}\bigg{)}=\mathbb{P}\bigg{(}\bigcap_{a\in\\{1,\dots,k\\}-\\{j\\}}\\{g_{x_{a}}(t)\geq
m_{a}\\}\bigg{)},$

so one can successively remove indices until what remains is an increasing
sequence

$m_{a_{1}}<\dots<m_{a_{l}}.$

######  Remark 1.4.

Since $m_{1}<\dots<m_{k}$, then

$\\{x\colon\eta_{x}(0)\geq
m_{1}\\}\supseteq\cdots\supseteq\\{x\colon\eta_{x}(0)\geq m_{k}\\},$

so $\mu_{1}\leq\dots\leq\mu_{k}$. If $\mu_{i}=\infty$ for some $i$, then
$\mu_{k}=\infty$, so the $k$th row of $G$ consists entirely of zeroes. This is
consistent with $\mathbb{P}(h_{x_{k}}(t)\geq\infty)=0$. Similarly, if
$\mu_{i}=-\infty$ for some $i$, then $\mu_{l}=-\infty$ for $1\leq l\leq i$, so
the $l$th row of $G$ equals $(\underbrace{1,\dots,1}_{l\ \text{times
}},0,0,\dots,0)$ for $1\leq l\leq i$. Then the matrix $G$ is a block matrix of
the form

$\left(\begin{array}[]{@{}c|c@{}}\begin{matrix}1&0&\cdots&0\\\ 1&1&\cdots&0\\\
\vdots&\vdots&\ddots&0\\\ 1&1&\cdots&1\end{matrix}&\mbox{\Large*}\\\
\hline\cr\mbox{\Large*}&\begin{matrix}\mbox{\Large*}\end{matrix}\end{array}\right),$

so its determinant equals the determinant of the lower right
$(k-l)\times(k-l)$ block. This is consistent with
$\mathbb{P}(h_{x_{1}}(t)\geq-\infty,\dots,h_{x_{i}}(t)\geq-\infty)=1$.

######  Remark 1.5.

Two initial conditions which satisfy the assumptions of Corollary 1.2 are the
step initial conditions

$\eta_{x}(t)=x\qquad\text{for}\quad x\in\mathbb{Z},$

and the flat initial conditions

$\eta_{x}(t)=\begin{cases}x/2&\text{for }x\text{ even},\\\ -\infty&\text{for
}x\text{ odd}.\end{cases}$

In these two cases, respectively, $\mu_{i}=m_{i}$ and $\mu_{i}=2m_{i}$. Then

$\\{g_{x_{i}}(t)\geq
m_{i}\\}=\left\\{\frac{g_{x_{i}}(t)-\nu_{i}t}{t^{1/2}}\geq
w_{i}\right\\}\qquad\text{or}\qquad\left\\{\frac{2g_{x_{i}}(t)-\nu_{i}t}{t^{1/2}}\geq
w_{i}\right\\}.$

######  Remark 1.6.

For the (single-species) TASEP with various initial conditions (step, flat,
stationary, and their mixtures), it has been shown in [2, 4, 5, 6, 7, 10] that
the fluctuations in the usual height function converges to various Airy
processes; see also [12]. Intuitively, these fluctuations are asymptotically
independent of the fluctuations of $h_{x_{i}}$ in Corollary 1.2. For example,
letting $h_{x}(t)$ denote the number of particles to the left of $x$ at time
$t$, then with step initial conditions

$\lim_{t\rightarrow\infty}\mathbb{P}\left(\frac{h_{x}(t)-a_{1}t}{a_{2}t^{1/3}}\geq-s,\frac{g_{x}(t)-b_{1}t}{b_{2}t^{1/2}}\leq
s^{\prime}\right)=F_{2}(s)\Phi(s^{\prime}),$

where $F_{2}$ is the Tracy–Widom distribution. The position $x$ depends on $t$
as $x=-\nu t$ for $\nu\in(0,1)$, and the constants are given by
$a_{1}=\tfrac{1}{4}(1-\nu)^{2}$,
$a_{2}=2^{-4/3}\big{(}1-\nu^{2}\big{)}^{2/3},$ and $b_{1}=1-\nu$,
$b_{2}=(1-\nu)^{1/2}$. To see independence, simply note that the fluctuations
of $h_{x}(t)$ depend on the exponential clocks of the particles starting with
distance less than $a_{1}t+O\big{(}t^{1/3}\big{)}$ of $x$, whereas the
fluctuations of $g_{x}(t)$ depend on the exponential clocks of the particles
starting with distance $b_{1}t+O\big{(}t^{1/2}\big{)}$ from $x$. (A rigorous
proof requires a quantification of the statement that far-away clocks have
negligible contribution to the fluctuations). These sets of clocks are
independent of each other. Note that an independent Tracy–Widom and Gaussian
also appear in the asymptotic crossing probability of the AHR model [9].

######  Remark 1.7.

Determinantal expressions in the multi-species TASEP have appeared before in
different contexts: [8, 11], see also [13].

##  2 Proof of Theorem 1.1

The theorem is proved by coupling the multi-species TASEP to a coalescing
random walk, and using a known determinantal expression for the distributions
of coalescing random walk from Proposition 2.5 of [1] (see also Proposition 9
of [14] for a similar expression in the context of coalescing Brownian
motions). This coupling does not seem to have appeared in the literature so
far; but note the upcoming preprint [3], in which a similar coupling is being
used.

Before describing the coupling, we first need to modify the initial conditions
$\eta_{x}(0)$. There may be values of $x$ and $y$ where $x<y$ and
$\eta_{x}(0)>\eta_{y}(0)$. In other words, this means that the particle
starting at $x$ has higher species number and is to the left of the particle
starting at $y$. Thus, the particle starting at $y$ can never contribute to
any value of $g_{z}(t)$ for $z\in\mathbb{Z}$ and $t\geq 0$, because any
contribution made by the particle starting at $y$ has already been made by the
particle starting at $x$. Thus, removing the particle at $y$ in the initial
conditions will not change the joint distribution of $g_{x_{i}}(t)$. Therefore
we may assume that the particles in the initial conditions are ordered by
species number.

In the coalescing random walk, if a particle tries to jump to a site that is
already occupied, then the two particles coalesce into a single particle. The
particle at site $x$ has an exponential clock of rate $\lambda_{x}$, and jumps
one step to the left when the clock rings, just as in the multi-species TASEP.
The coupling is constructed as follows. The coalescing (multi-species) random
walk and the multi-species TASEP begin with the same initial conditions. We
couple the exponential clocks of each particle in the multi-species TASEP with
the exponential clocks of the corresponding particle in the coalescing random
walk. In this way, each particle in the multi-species TASEP is coupled with a
particle in the coalescing random walk at time $t=0$. The evolution occurs
with the following rules:

  1. 1.

If a particle in the multi-species TASEP attempts to jump but the jump is
blocked, then no change occurs in either system.

  2. 2.

If a particle in the multi-species TASEP jumps one step to the left and that
site is unoccupied, then the corresponding coupled particle (if there is one)
in the coalescing random walk makes the same jump.

  3. 3.

If a particle in the multi-species TASEP jumps one step to the left, and that
site is occupied by a lower species particle, then the corresponding coupled
particle (if there is one) makes the same jump in the coalescing random walk.
Additionally, the lower species particle is then absorbed into the higher
species particle in the coalescing random walk (if it had not already been
absorbed previously), and is no longer coupled with the corresponding particle
in the multi-species TASEP.

  4. 4.

If a particle in the multi-species TASEP attempts a jump but there is no
coupled particle in the coalescing random walk, then no jump occurs in the
coalescing random walk.

One can see that two particles are coupled if and only if they occupy the same
location. Below is an example of the particle updates illustrating each of the
four rules. The color red is used to indicate that particle attempts to jump
to the left:

$\displaystyle\underline{\ }\ \underline{1}\
\underline{{\color[rgb]{1,0,0}2}}\ \underline{3}\ \underline{4}\
\stackrel{{\scriptstyle\text{rule 3}}}{{\longrightarrow}}\underline{\ }\
\underline{{\color[rgb]{0,0,0}2}}\ \underline{{\color[rgb]{1,0,0}1}}\
\underline{3}\ \underline{4}\ \stackrel{{\scriptstyle\text{rule
1}}}{{\longrightarrow}}\underline{\ }\ \underline{{\color[rgb]{1,0,0}2}}\
\underline{1}\ \underline{3}\ \underline{4}\ \stackrel{{\scriptstyle\text{rule
2}}}{{\longrightarrow}}\underline{2}\ \underline{\ }\
\underline{{\color[rgb]{1,0,0}1}}\ \underline{3}\
\underline{4}\stackrel{{\scriptstyle\text{rule
4}}}{{\longrightarrow}}\underline{2}\ \underline{1}\ \underline{\ }\
\underline{3}\ \underline{4}\,,$ $\displaystyle\underline{\ }\ \underline{1}\
\underline{{\color[rgb]{1,0,0}2}}\ \underline{3}\ \underline{4}\
\stackrel{{\scriptstyle\text{rule 3}}}{{\longrightarrow}}\underline{\ }\
\underline{2}\ \underline{{\color[rgb]{1,1,1}1}}\ \underline{3}\
\underline{4}\ \stackrel{{\scriptstyle\text{rule
1}}}{{\longrightarrow}}\underline{\ }\ \underline{{\color[rgb]{1,0,0}2}}\
\underline{{\color[rgb]{1,1,1}1}}\ \underline{3}\ \underline{4}\
\stackrel{{\scriptstyle\text{rule 2}}}{{\longrightarrow}}\underline{2}\
\underline{\ }\ \underline{{\color[rgb]{1,1,1}1}}\ \underline{3}\
\underline{4}\stackrel{{\scriptstyle\text{rule
4}}}{{\longrightarrow}}\underline{2}\ \underline{\ }\
\underline{{\color[rgb]{1,1,1}1}}\ \underline{3}\ \underline{4}\,.$

Also note that there is no rule to couple the system after updates of the form

$\displaystyle\underline{2}\
\underline{{\color[rgb]{1,0,0}1}}\longrightarrow\underline{2}\
\underline{1}\,,$ $\displaystyle\underline{\circ}\
\underline{{\color[rgb]{1,0,0}{\circ}}}\longrightarrow\underline{\circ}\
\underline{\ }\,.$

However, the initial conditions preclude these configurations from occurring,
because the particles are ordered by species number.

Let $\tilde{\eta}_{x}(t)$ denote the species of the particle located at site
$x$ at time $t$ in the coalescing random walk. As before, let
$\tilde{\eta}_{x}(t)=-\infty$ by convention if there is no particle at site
$x$ and time $t$. Now let

$\displaystyle\tilde{g}_{y}(t)=\max\\{\tilde{\eta}_{x}(t)\colon x<y\\}.$

######  Lemma 2.1.

For $x_{1}<\dots<x_{k}$ and any fixed time $t$, the joint distribution of
$(g_{x_{1}}(t),\dots,g_{x_{k}}(t))$ equals the joint distribution of
$(\tilde{g}_{x_{1}}(t),\dots,\tilde{g}_{x_{k}}(t))$. In other words, for any
$m_{1}<\dots<m_{k}$,

$\mathbb{P}(g_{x_{1}}(t)\geq m_{1},\dots,g_{x_{k}}(t)\geq
m_{k})=\mathbb{P}(\tilde{g}_{x_{1}}(t)\geq
m_{1},\dots,\tilde{g}_{x_{k}}(t)\geq m_{k}).$ (2.1)

###### Proof.

At time $t=0$, the multi-species TASEP and the coalescing random walk are in
the same configuration, so (2.1) holds. The only potential way that the two
sides of (2.1) could be different is when a particle jump causes the two
processes to be updated in a different manner – that is, update rules 3 or 4.

Suppose that an update under rule 4 occurs. This means that a particle has
jumped in multi-species TASEP when there is no coupled particle in the
coalescing random walk. This means that the particle in the coalescing random
walk has already been absorbed into a higher species particle.
Correspondingly, this means that in the multi-species TASEP, there is a higher
species particle to the left of that particle. This means that the evolution
of the uncoupled/absorbed particle no longer affects the values of
$g_{x_{i}}(t)$ or $\tilde{g}_{x_{i}}(t)$, since the particle cannot cross
$x_{i}$ unless the higher species particle has also crossed $x_{i}$.
Similarly, if an update under rule 3 occurs, then the multi-species TASEP and
the coalescing random walk differ in the uncoupled/absorbed particle. For
identical reasons, the evolution of this particle does not affect the values
of $g_{x_{i}}(t)$ or $\tilde{g}_{x_{i}}(t)$. ∎

By Lemma 2.1, it suffices to find the right-hand side of (2.1). As mentioned
in Remark 1.4, $\mu_{1}\leq\dots\leq\mu_{k}$, where $\mu_{i}$ is defined as in
Theorem 1.1. Let $\tilde{x}_{\mu}(t)$ denote the position of the particle in
the coalescing random walk that started at lattice site $\mu$ at time $0$.
Then it is straightforward to see that there is the inclusion of events

$\\{\tilde{g}_{x_{i}}(t)\geq m_{i}\\}\supseteq\\{\tilde{x}_{\mu_{i}}(t)\leq
x_{i}\\}.$

To see the reverse inclusion, note that $\tilde{x}_{\mu_{i}}(t)$ is the left-
most particle which has species higher than $m_{i}$, and remains the left-most
such particle after coalescence. Therefore,

$\mathbb{P}(\tilde{g}_{x_{1}}(t)\geq m_{1},\dots,\tilde{g}_{x_{k}}(t)\geq
m_{k})=\mathbb{P}(\tilde{x}_{\mu_{1}}(t)\leq
x_{1},\dots,\tilde{x}_{\mu_{k}}(t)\leq x_{k}).$

Note that the random variables $\tilde{x}_{\mu}(t)$ on the right-hand side do
not depend on the particles’ species. The probability on the right-hand side
can be found from Proposition 2.5 of [1]. The cited proposition applies more
generally (it allows for right jumps as well), and includes a reflecting
boundary. We can take the limit of this boundary point to $-\infty$ and apply
the proposition. The result is that the probability equals $\det[G_{ij}]$, as
needed.

### Acknowledgements

The author thanks Alexei Borodin and Alexey Bufetov for helpful conversations.

## References

  * [1] Assiotis T., Random surface growth and Karlin–McGregor polynomials, Electron. J. Probab. 23 (2018), 106, 81 pages, arXiv:1709.10444.
  * [2] Baik J., Ferrari P.L., Péché S., Limit process of stationary TASEP near the characteristic line, Comm. Pure Appl. Math. 63 (2010), 1017–1070, arXiv:0907.0226.
  * [3] Borodin A., Bufetov A., Leading particles in multi-type TASEP, in preparation.
  * [4] Borodin A., Ferrari P.L., Prähofer M., Fluctuations in the discrete TASEP with periodic initial configurations and the ${\rm Airy}_{1}$ process, Int. Math. Res. Pap. 2007 (2007), rpm002, 47 pages, arXiv:math-ph/0611071.
  * [5] Borodin A., Ferrari P.L., Prähofer M., Sasamoto T., Fluctuation properties of the TASEP with periodic initial configuration, J. Stat. Phys. 129 (2007), 1055–1080, arXiv:math-ph/0608056.
  * [6] Borodin A., Ferrari P.L., Sasamoto T., Large time asymptotics of growth models on space-like paths. II. PNG and parallel TASEP, Comm. Math. Phys. 283 (2008), 417–449, arXiv:0707.4207.
  * [7] Borodin A., Ferrari P.L., Sasamoto T., Transition between ${\rm Airy}_{1}$ and ${\rm Airy}_{2}$ processes and TASEP fluctuations, Comm. Pure Appl. Math. 61 (2008), 1603–1629, arXiv:math-ph/0703023.
  * [8] Chatterjee S., Schütz G.M., Determinant representation for some transition probabilities in the TASEP with second class particles, J. Stat. Phys. 140 (2010), 900–916, arXiv:1003.5815.
  * [9] Chen Z., de Gier J., Hiki I., Sasamoto T., Exact confirmation of 1d nonlinear fluctuating hydrodynamics for a two-species exclusion process, Phys. Rev. Lett. 120 (2018), 240601, 6 pages, arXiv:1803.06829.
  * [10] Johansson K., Discrete polynuclear growth and determinantal processes, Comm. Math. Phys. 242 (2003), 277–329, arXiv:math.PR/0206208.
  * [11] Lee E., On the TASEP with second class particles, SIGMA 14 (2018), 006, 17 pages, arXiv:1705.10544.
  * [12] Quastel J., Remenik D., KP governs random growth of a one dimensional substrate, arXiv:1908.10353.
  * [13] Schütz G.M., Exact solution of the master equation for the asymmetric exclusion process, J. Stat. Phys. 88 (1997), 427–445, arXiv:cond-mat/9701019.
  * [14] Warren J., Dyson’s Brownian motions, intertwining and interlacing, Electron. J. Probab. 12 (2007), no. 19, 573–590, arXiv:math.PR/0509720.