Solutions of a Two-Particle Interacting Quantum Walk

We study the solutions of an interacting Fermionic cellular automaton which is the analogue of the Thirring model with both space and time discrete. We present a derivation of the two-particle solutions of the automaton recently in the literature, which exploits the symmetries of the evolution operator. In the two-particle sector, the evolution operator is given by the sequence of two steps, the first one corresponding to a unitary interaction activated by two-particle excitation at the same site, and the second one to two independent one-dimensional Dirac quantum walks. The interaction step can be regarded as the discrete-time version of the interacting term of some Hamiltonian integrable system, such as the Hubbard or the Thirring model. The present automaton exhibits scattering solutions with nontrivial momentum transfer, jumping between different regions of the Brillouin zone that can be interpreted as Fermion-doubled particles, in stark contrast with the customary momentum-exchange of the one-dimensional Hamiltonian systems. A further difference compared to the Hamiltonian model is that there exist bound states for every value of the total momentum and of the coupling constant. Even in the special case of vanishing coupling, the walk manifests bound states, for finitely many isolated values of the total momentum. As a complement to the analytical derivations we show numerical simulations of the interacting evolution.


Introduction
Quantum walks (QWs) describe the evolution of one-particle quantum states on a lattice, or more generally, on a graph.The quantum walk evolution is linear in the quantum state and the quantum aspect of the evolution occurs in the interference between the different paths available to the walker.There are two kinds of quantum walks: continuous time QWs, where the evolution operator of the system given in terms of an Hamiltonian can be applied at any time (see Farhi et al. [2]), and discrete-time QWs, where the evolution operator is applied in discrete unitary time-steps.The discrete-time model, which appeared already in the Feynman discretization of the Dirac equation [3], was later rediscovered in quantum information [4][5][6][7][8], and proved to be a versatile platform for various scopes.For example, QWs have been used for empowering quantum algorithms, such as database search [9,10], or graph isomorphism [11,12].Moreover, quantum walks have been studied as a simulation tool for relativistic quantum fields [13][14][15][16][17][18][19][20][21][22][23][24][25][26][27][28], and they have been used as discrete models of spacetime [29][30][31][32].
New research perspectives are unfolding in the scenario of multi-particle interacting quantum walks where two or more walking particles are coupled via non-linear (in the field) unitary operators.The properties of these systems are still largely unexplored.Both continuous-time [40] and discrete-time [41] quantum walks on sparse unweighted graphs are equivalent in power to the quantum circuit model.However, it is highly non trivial to design a suitable architecture for universal quantum computation based on quantum walks.Within this perspective a possible route has been suggested in Ref. [42] based on interacting multi-particle quantum walks with indistinguishable particles (Bosons or Fermions), proving that "almost any interaction" is universal.Among the universal interacting many-body systems are the models with coupling term of the form χδ x 1 ,x 2 n(x 1 ) n(x 2 ), with n(x) the number operator at site x.The latter two-body interaction lies at the basis of notable integrable quantum systems in one space dimension such as the Hubbard and the Thirring Hamiltonian models.
The first attempt at the analysis of interacting quantum walks was carried out in Ref. [43].More recently, in Ref. [1], the authors proposed a discrete-time analogue of the Thirring model which is indeed a Fermionic quantum cellular automaton, whose dynamics in the two-particles sector reduces to an interacting two-particle quantum walk.As for its Hamiltonian counterpart, the discrete-time interacting walk has been solved analytically in the case of two Fermions.Analogously to any Hamiltonian integrable system, also in the discrete-time case the solution is based on the Bethe Ansatz technique.However, discreteness of the evolution prevents the application of the usual Ansatz, and a new Ansatz has been introduced successfully [1].
In this paper we present an original simplified derivation of the solution of Ref. [1] which exploits the symmetries of the interacting walk.We present the diagonalization of the evolution operator and the characterization of its spectrum.We explicitly write the two particle states corresponding to the scattering solutions of the system, having eigenvalues in the continuous spectrum of the evolution operator.We then show how the present model predicts the formation of bound states, which are eigenstates of the interacting walk corresponding to the discrete spectrum.We provide also in this case the analytic expression of such molecular states.
We remark the phenomenological differences between the Hamiltonian model and the discrete-time one.First we see that the set of possible scattering solutions is larger in the discrete-time case: for a fixed value total momentum, a non trivial transfer of relative momentum can occur besides the simple exchange of momentum between the two particles, differently from the Hamiltonian case.Also the family of bound states appearing in the discrete-time scenario is larger than the corresponding Hamiltonian one.Indeed, for any fixed value of the coupling constant, a bound state exists with any possible value of the total momentum, while for Hamiltonian systems bound states cannot have arbitrary total momentum.
Finally we show that in the set of solutions for the interacting walk there are perfectly localized states (namely states which lie on a finite number of lattice sites) and, differently from the Hamiltonian systems, bound states exist also for vanishing coupling constant.In addition to the exact analytical solution of the dynamics we show the simulation of some significant initial state.

The Dirac Quantum Walk
In this section, we review the Dirac walk on the line describing the free evolution of a two-component Fermionic field.The walk evolution is provided by the unitary operator W acting on the single particle Hilbert space H := C 2 ⊗ 2 (Z) for which we employ the factorized basis |a |x , with a ∈ {↑, ↓} and x ∈ Z. Being the evolution of a quantum walk linear in the field operators, the single-step evolution is expressed by the following equation: where T x denotes the translation operator on 2 (Z), defined by Since the walk W is translation invariant (it commutes with the translation operator), it can be diagonalized in momentum space.In the momentum representation, defining |p := (2π) −1/2 ∑ x∈Z e −ipx |x , with p ∈ B := (−π, π], the walk operator can be written as where |ν| 2 + |µ| 2 = 1.The spectrum of the walk is given by {e −iω(p) , e iω(p) }, where the dispersion relation ω(p) is given by where Arccos denotes the principal value of the arccosine function.The single-particle eigenstates, solving the eigenvalue problem can be conveniently written as

The Thirring Quantum Walk
In this section we present a Fermionic cellular automaton in one spatial dimension with an on-site interaction, namely two particles interact only when they lie at the same lattice site.The linear part corresponds to the Dirac QW [44] and the interaction term is the most general number-preserving coupling in one dimension [45].The same kind of interaction characterizes also the most studied integrable quantum systems, such as the Thirring [46] and the Hubbard [47] models.
The linear part of the N-particle walk is described by the operator W N := W ⊗N , acting on the Hilbert space H N = H ⊗N and describing the free evolution of the particles.In order to introduce an interaction, we modify the update rule of the walk with an extra step V int : U N := W N V int .In the present case the term V int has the form where n a (x), a ∈ {↑, ↓}, represents the particle number at site x, namely n a (x) = ψ † a (x)ψ a (x), and χ is a real coupling constant.Since the interaction term preserves the total number operator we can study the walk dynamics for a fixed number of particles.In this work we focus on the two-particle sector whose solutions has been derived in Ref. [1].As we will see, the Thirring walk features molecule states besides scattering solutions.This features is shared also by the Hadamard walk with the same on-site interaction [48].
Since we focus on the solutions involving the interaction of two particles, it is convenient to write the walk in the centre of mass basis |a 1 , a 2 |y |w , with a 1 , a 2 ∈ {↑, ↓}, y = x 1 − x 2 and w = x 1 + x 2 .Therefore in this basis the generic Fermionic state is |ψ = ∑ a 1 ,a 2 ,y,w c(a 1 , a 2 , y, w) |a 1 , a 2 |y |w with c(a 2 , a 1 , y, w) = −c(a 1 , a 2 , −y, w).Notice that only the pairs y, w with y and w both even or odd correspond to physical points in the original basis x 1 , x 2 .
We define the two-particle walk with both y and w in Z, so that the linear part of walk can be written as where T y represents the translation operator in the relative coordinate y, and T w the translation operator in the centre of mass coordinate w, whereas the interacting term reads This definition gives a walk U 2 = W 2 V 2 (χ) that can be decomposed in two identical copies of the original walk.Indeed, defining as C the projector on the physical center of mass coordinates, one has where CU 2 C and (I − C)U 2 (I − C) are unitarily equivalent.We will then diagonalize the operator U 2 , reminding that the physical solutions will be given by projecting the eigenvectors with C.
Introducing the (half) relative momentum k = 1 2 (p 1 − p 2 ) and the (half) total momentum p = 1 2 (p 1 + p 2 ), the free evolution of the two particles is written in the momentum representation as where the matrix W 2 (p, k) is defined as Furthermore, we introduce the vectors v sr k := v s p+k ⊗ v r p−k , with s, r = ±, such that where is the dispersion relation of the two-particle walk.Explicitly, the vectors v sr k are given by We focus in this work on Fermionic solutions satisfying the eigenvalue equation with |ψ(y) ∈ C 4 .In the centre of mass basis the antisymmetry condition reads .

Symmetries of the Thirring Quantum Walk
The Thirring walk manifests some symmetries that allow to simplify the derivation and the study of the solutions.First of all, as we already mentioned, one can show that the interaction V(χ) commutes with the total number operator.This means that one can study the walk dynamics separately for each fixed number of particles.We focus here on the two-particle walk Since the interacting walk U 2 commutes with the translations in the centre of mass coordinate w, the total momentum is a conserved quantity, so it is convenient to study the walk parameterized by the total momentum p.To this end we consider the basis |a 1 , a 2 |y |p , so that for fixed values of p the interacting walk of two particles can be expressed in terms of a one-dimensional QW U 2 (χ, p) = W 2 (p)V(χ) with a four dimensional coin: , Although the range of the variable p is the interval (−π, π], it is possible to show that one can restrict the study of the walk to the interval [0, π/2].On the one hand, the two-particle walk transforms unitarily under a parity transformation in the momentum space.Starting from the single particle walk, W(p) transforms under a parity transformation as so that for the two-particle walk we have the relation On the other hand, a translation of π of the total momentum p entails that while the interaction term remains unaffected in both cases.The Thirring walk features also another symmetry that can be exploited to simplify the derivation of the solutions.It is easy to check that the walk operator U 2 (p, χ) = W 2 (p)V(χ) commutes with the projector defined by where P e and P o are the projectors on the even and the odd subspaces, respectively: The projector P induces a splitting of the total Hilbert space H in two subspaces PH and (I − P)H , with the interaction term acting non trivially only in the subspace PH .In the complementary subspace (I − P)H the evolution is free for Fermionic particles.This means that solutions of the free theory are also solutions of the interacting one, as opposed to the Bosonic case for which the interaction is non-trivial also in (I − P)H .

Review of the solutions
We focus in this section on the antisymmetric solutions of the Thirring walk which actually feel the interaction.From the remarks that we have made in the previous section, such solutions can only be found in the subspace PH .Formally, we have to solve the eigenvalue equation PU 2 (χ, p) |ψ = e −iω |ψ , with |ψ ∈ PH .Conveniently, we write a vector |ψ ∈ PH in the form and the antisymmetry condition becomes: The restriction of the walk to the subspace PH entails that the eigenvalue problem is equivalent to the following system of equations: The most general solution of Eq. ( 10) for p ∈ {0, π/2} has two forms: and Continuous spectrum of the two-particle walk as a function of the total momentum p ∈ [0, π/2] with mass parameter m = 0.7.The continuous spectrum is the same as in the free case.
The solid blue curves are described by the functions ω = ±2ω(p), and the red ones by ω = ±(π − 2 Arccos(n sin p)).As one can notice, the light-red lines ω = ±2p lie entirely in the gaps between the solid curves, highlighting the fact that e ±i2p is not in the range of e −iω sr (p,k) for p = 0, π/2 (see text). with }, and g sr ω satisfying the condition Solving Eq. ( 10) corresponds now to find the function g sr ω .
Let us now study the equation Since e −iω sr (p,k) has to be an eigenvalue of U 2 (χ, p), ω sr (p, k) must be real and thus k ∈ Γ f or k ∈ Γ l with l = 0, ±1, 2, so we conveniently define the sets: It is easy to see that Ω sr f ∩ Ω sr l = ∅ for all s, r and l, and the range of the function e −iω sr (p,k) covers the entire unit circle except for the points e ±i2p .Therefore, we can discuss separately the case e −iω ∈ Ω sr f and the case e −iω ∈ Ω sr l .A solution with e −iω = e ±i2p actually exists, corresponding to the function of Eq. (11), and it will be discussed in Section 5.3.
Let us start with the case e −iω ∈ Ω sr f which will lead to the characterization of the continuous spectrum of the Thirring walk U 2 (χ, p) and of the scattering solutions.

Scattering solutions
In this section we assume p ∈ {0, π/2} with e −iω ∈ Ω sr f .This implies that e −iω = e ±i2p : indeed, as one can notice from Fig. 1, the lines ω = ±2p lie entirely in the gaps between the curves ω = ±2ω(p) and ω = ±(π − 2 Arccos(n sin p)).The solution is thus the one given in Eq. (12).One can prove that . Furthermore, as one can notice from Fig. 2, there are four values of the triple (s, r, k) such that e −iω sr (p,k) = e −iω for a given value of e −iω : if the triple (+, +, k) is a solution, so are (+, +, π − k), (−, −, −k) and (−, −, k − π); and if (+, −, k) is a solution, then also (+, −, π − k), (−, +, −k) and (−, +, k − π) are solutions.This result greatly simplifies Eq. (12).Indeed the sum over s, r and the integral over k reduces to the sum of four terms: As we will see, the original problem can be simplified in this way to an algebraic problem with a finite set of equations.We remark that the fact that the equation e −iω sr (p,k) = e −iω has a finite number of solutions is a consequence of the fact that we are considering a model in one spatial dimension.However, in analogous one-dimensional Hamiltonian models (e.g. the Hubbard model) the degeneracy of the eigenvalues is two.
Let us consider for the sake of simplicity the solution of the kind ψ +,j k (z), since the other one can be analysed in a similar way.Using the notation of Appendix A, Eq. ( 13) reduces to the expressions (dropping the + superscript) We notice that now the number of unknown parameters is further reduced to three, namely λ, ρ, and ξ.
Clearly, one of the parameters can be fixed by choosing arbitrarily the normalization.From now on we fix λ = 1 and define T + := ρ.Eq. ( 14) has to satisfy the recurrence relations of Eq. ( 10) for z = 0 and z = 1, while for z > 1 it is automatically satisfied.For z = 0, Eq. ( 10) becomes Starting from Eq. (15a), we can notice that ν 2 e i2p a − iµνe ip e ik b − iµνe ip e −ik c − µ 2 d = e −iω a, where we employed the notation of Appendix A, so that we obtain ξ = e −iχ (b − T + c).We can then substitute this expression in Eq. (15c) and use the relations  For these values of ξ and T + one can verify that Eq. ( 10) is satisfied also for z = 1, thus concluding the derivation.For the solution of the kind ψ −,j k (z) we can follow a similar reasoning, obtaining the analogous quantity T − : It is worth noticing that T ± is of unit modulus for k ∈ (−π, π].
The final form of the solution results to be which in terms of the relative coordinate y can be written as We can interpret such a solution as a scattering of plane waves for which the coefficient T ± plays the role of the transmission coefficient.Being the total momentum a conserved quantity, the two particles can only exchange their momenta, as expected from a theory in one-dimension.Furthermore, for each value k of the relative momentum, the two particles can also acquire an additional phase of π.As the interaction is a compact perturbation of the free evolution, the continuous spectrum is the same as that of the free walk.Eq. ( 18) provides the generalized eigenvector if U 2 (χ, p) corresponding to the continuous spectrum

Bound states
In the previous section, we derived the solutions in the continuous spectrum, which can be interpreted as scattering plane waves in one spatial dimension.We seek now the solutions corresponding to the discrete spectrum, namely solutions with eigenvalue in any one of the sets highlighting the appearance of bound states components along the diagonal, namely when the two particles are at the same site (i.e.x 1 = x 2 ), where x 1 and x 2 denote the positions of the two particles.The plots show the probability distribution p(x 1 , x 2 ) in position space after t = 32 time-steps.The chosen value of the mass parameter is m = 0.6 and the coupling constant is χ = π/2.The two particles are initially prepared in a singlet state located at the origin.
Ω sr l .The derivation of the solution follows similar steps as for the scattering solutions.In particular, the degeneracy in k is the same: there are four solutions to the equation e −iω sr (p,k) = e −iω even in this case, as proved in Ref. [1].Therefore the general form of the solution in this case can be written again as in Eq. ( 13) and, following the same reasoning, one obtains the same set of solutions as in Eq. (18).At this stage we did not imposed that the solution is a proper eigenvector in the Hilbert space H .To this end, we have to set T ± = 0 to eliminate the exponentially-divergent terms in Eq. (18).As one can prove, the equation T ± = 0 has only one solution for fixed values of χ and p.More precisely, there is a unique k ∈ Γ 0 ∪ Γ −1 ∪ Γ 1 ∪ Γ 2 , with k I < 0 and e iχ ∈ {1, −1}, such that either T + = 0 or T − = 0.
In other words, for each pair of values (χ, p) the walk U 2 (p) has one and only eigenvector corresponding to an eigenvalue in the point spectrum.Such eigenvector can be written as where k is the solution of T + = 0 or T − = 0 and ± chosen accordingly.More compactly, in the y coordinate, the solution can be written as Referring to Fig. 3, we show the evolution of two particles initially prepared in a singlet state localized at the origin.From the figure one can appreciate the appearance of the bound state component which has non-vanishing overlapping with the initial state.The bound state, being exponentially decaying in the relative coordinate y, is localized on the diagonal of the plot, that is when the two particles lie at the same point.In Fig. 5 is depicted the probability distribution of the bound state corresponding to choice of the parameters χ = 0.2π and p = 0.035π.The plot highlights the exponential decay of the tails, which is the characterizing feature of the bound state.

Solution for e −iω = e ±i2p
So far we have studied proper eigenvectors which decay exponentially as the two particles are further apart.However, the previous analysis failed to cover the particular case when e −iω = e ±i2p , since the range of e −iω sr (p,k) does not include the two points of the unit circle e ±i2p .
We study now the solutions with e −iω = e ±i2p having the form given in Eq. (11).One can prove that such solutions are non-vanishing only for z = 0 on PH , namely we look for a solution of the form Subtracting the first and the last equations of (10) using (20), we obtain the following equation: If both ζ and ζ are non-zero, one can prove that a solution does not exist and thus we have to consider the two cases ζ = 0 and ζ = 0 separately.Starting from ζ = 0, Eq. ( 21) imposes that e −iω = e i2p , meaning that if a solution exists in this case, it is an eigenvector corresponding to the eigenvalue e i2p .From the second equation of (10) we obtain the relation and, using the first equation of (10), it turns out that a solution exists only if e iχ = e i2p , as expected since otherwise we would have been in the case of Section 5.2 would hold.The other case, namely e −iω = e −i2p , can be studied analogously.Let us, then, denote as |ψ ±∞ such proper eigenvectors with eigenvalue e ±i2p for χ = e ±i2p and, choosing η = µ ν as the value for the free parameter η, we obtain the following expression for |ψ ±∞ : Such solutions provide a special case of molecule states (namely, proper eigenvectors of U 2 (χ, p)), being localized on few sites, and differ from the previous solutions showing an exponential decay in the relative coordinate.

Solutions for p ∈ {0, π/2}
The solutions that we presented in the previous discussion do not cover the extreme values p = 0, π/2 (see Ref. [1] for a reference).Let us consider for definiteness the case p = 0 since the other case is obtained in a similar way.For e −iω = 1 the previous analysis still holds.Indeed, noticing that ω ±± (0, k) = ±2ω(k), we have ω(k) ∈ R and ω(k) = 0 if and only if k ∈ Γ f ∪ Γ 0 ∪ Γ 2 , whereas ω ±∓ (0, k) = 0 for all k ∈ C.This means that the solutions |ψ + k of Eq. ( 18) are actually eigenvectors of U 2 (χ, 0).Thus, the spectrum is made by a continuous part, given by the arc of the unit circle containing −1 and having e ±2iω(0) as extremes, and a point spectrum with two points: e −2iω( k) , where k is the solution of T + = 0 for p = 0, and 1.As shown in Ref. [1], 1 is a separated part of the spectrum of U 2 (χ, 0) and the corresponding eigenspace is a separable Hilbert space of stationary bound states.This fact underlines an important feature of the Thirring walk not shared by analogous Hamiltonian models.It is remarkable that this behaviour occurs also for the free walk with χ = 0.In Fig. 6 we show the probability distribution of two states having he properties hereby discussed.It is worth noticing that all the states v +− k with k ∈ (−π, π] are eigenvectors relative to the eigenvalue 1, and thus they generate a subspace on which the walk acts identically.We remark that this behaviour relies on the fact that the dispersion relation in one dimension is an even function of k.
In Fig. 4 is depicted the discrete spectrum of the interacting walk together with the continuous spectrum as a function of the total momentum p.The solid curves in the gaps between the continuous bands denote the discrete spectrum for different values of the coupling constant χ = 2π/3, 3π/7, −3π/7, −2π/3.Molecule states appear also in the Hadamard walk with the same on-site interaction [48].

Conclusions
In this work we reviewed the Thirring quantum walk providing a simplified derivation of its solutions for Fermionic particles.The simplified derivation relies on the symmetric properties of the walk evolution operator allowing to separate the subspace of solutions affected by the interaction from the subspace where the interaction step acts trivially.The interaction term is the most general number-preserving interaction in one dimension, whereas the free evolution is provided by the Dirac QW [44].
We showed the explicit derivation of the scattering solutions (solutions for the continuous spectrum) as well as for the bound-state solutions.The Thirring walk features also localized bound states (namely, states whose support is finite on the lattice) when e −iω = e ±i2p .Such solutions exist only when the coupling constant is χ = 2p.In Fig. 3 is depicted the evolution of a perfectly localized state showing the overlapping with bound states components.In Fig. 5 we reported the evolution of a bound state of the two particles peaked around a certain value of the total momentum: one can appreciate that the probability distribution remains localized on the main diagonal during the evolution.
Finally, we discussed also the class of proper eigenvectors arising in the free theory highlighting another difference between the discrete model of the present work with analogous Hamiltonian models.
Explicitly, for the even case we have: Analogously for the odd case the eigenstates are In order to simplify the derivation of the solution, we adopt the following notation: .

Figure 2 .
Figure 2. Spectrum of the walk for m = 0.6 and p = π/6 as a function of k.The colours highlight the different ranges of eigenvalues corresponding to the dispersion ω sr (p, k).The range of ω sr (p, k) is understood to be computed mod(2π).One can notice that there are four values of the relative momentum k having the same value of the dispersion relation (ω = 2 in the figure).This is in contrast to the Hamiltonian model for which there are only two solutions.

Figure 3 .
Figure 3.We show for comparison the free evolution (left panel) and the interacting one (right panel) highlighting the appearance of bound states components along the diagonal, namely when the two particles are at the same site (i.e.x 1 = x 2 ), where x 1 and x 2 denote the positions of the two particles.The plots show the probability distribution p(x 1 , x 2 ) in position space after t = 32 time-steps.The chosen value of the mass parameter is m = 0.6 and the coupling constant is χ = π/2.The two particles are initially prepared in a singlet state located at the origin.

7 χ=-3 π 7 χ=-2 π 3 Figure 4 .
Figure 4. Complete spectrum of the two-particle Thirring walk as a function of the total momentum p with mass parameter m = 0.7.The continuous spectrum is as in Fig. 1.The solid lines in the gaps show the point spectrum for different values of the coupling constant: from top to bottom, χ = 2π/3, 3π/7, −3π/7, −2π/3.It is worth noticing that for each pair (χ, p) there is only one value in the discrete spectrum.The light-red lines ω = ±2p intersect the curves of the discrete spectrum for p = χ/2.

Figure 5 .Figure 6 .
Figure 5.We show the evolution of a bound state of the two particles peaked around the value of the total momentum p = 0.035π.The mass paramater is m = 0.6 and the coupling constant χ = 0.2π.On the left is depicted the probability distribution of the initial state and on the right the that of the evolved state after t = 128 time-steps.One can notice that in the relative coordinate x 1 − x 2 the probability distribution remains concentrated on the diagonal, highlighting that the two particles are in a bound state.The diffusion of the state happens only in the centre of mass coordinate.