Quantum backflow for a free-particle hermite wavepacket

Quantum backflow is the unexpected effect that wavepackets consisting of only positive momentum components can apparently move in the negative direction. This is usually described in terms of the backflow constant, which is a dimensionless quantity describing least upper bound on the amount of probability that can flow backwards during a given time interval. Backflow is usually calculated for wavepackets that can be written as a sum of positive momentum plane waves. Here we present a calculation of the backflow constant using the localised free particle hermite wavefunctions where equal weights of positive and negative momentum eigenfunctions occur. The resulting backflow constant is substantially smaller than the accepted value. The reasons for this are discussed and finally we draw conclusions about the calculation of backflow more generally.


Introduction
Quantum backflow is a remarkable, and yet relatively unknown, phenomenon that occurs in quantum mechanics.It is the initially unexpected effect that for a free particle described by a wavefunction containing only positive momenta, the probability density can flow in the negative direction.i.e. probability can flow in the opposite direction to the momentum in certain cases.Allcock [1] was the first to identify backflow in his study of the arrival time problem in quantum mechanics.However it was Bracken and Melloy [2] who first studied nonrelativistic quantum backflow systematically and showed that, although a period of backflow can be arbitrarily long, the increase in probability flowing backwards cannot exceed a limited amount given by a dimensionless number which they calculated to be approximately c bf = 0.04.The most precise estimates of this quantity have been computed numerically by Penz and co-workers [3] to be c bf = 0.0384517.A number of authors have provided further insight into this topic.More recently Berry [4] has derived a theory of the spatial-temporal extent of the backflow wavefunctions.Although he found large quantum backflow the wavefunctions used were not normalizable and so that work does not really address the problem of the magnitude of the quantum backflow constant.While most investigations of quantum backflow previously have considered free particles, Bostelmann and co-workers have extended the theory recently to systems in which there is a scattering potential [5].Halliwell [6] has considered a novel approach to the arrival time problem and shown that backflow provides an interesting example in which to test their formalism.They have also discussed backflow in terms of the Leggett-Garg inequalities.Strange [7] considered backflow in angular momentum, Goussev has examined backflow in the case when the particle is forced to move on a circular path [8] which has some advantages over the linear case, the backflow constant becomes the highest eigenvalue of a matrix in this case.Following up work by Melloy and Bracken [9] on backflow in the presence of a potential, Goussev [10] introduced the quantum reentry problem and discussed its connection to quantum backflow in the presence of a constant force.Here he was able to show an equivalence between quantum backflow and quantum reentry, a property we make use of in this paper.He has looked at backflow for correlated quantum states [11].In a very interesting recent paper Trillo et al [12] have shown that the backflow constant can also be seen as a measure of maximum quantum advantage at least under some circumstances.Van Dijk and Toyama [13] have also made a very comprehensive study of quantum reentry in the decay of quasi-stable quantum systems.
A number of authors have attempted to find analytic expressions for the maximum backflow eigenvector by fitting to the numerical values.This has proved difficult, but some progress has been made.Yearsley et al [14,15] and Halliwell et al [16] have found a number of expressions, including one based on Fresnel integrals which yield up to about 70% of the maximum possible backflow.O'Mullane [17] has attempted to fit Bessel functions to the numerical eigenvector and found an expression that yielded of order 20% of the maximum.A relativistic theory of quantum backflow was written down by Melloy and Bracken [18], while Ashfaque et al [19] and Su and Chen [20] have searched for analytic expressions for the eigenvectors in the relativistic case.
There are two main challenges to be addressed in the theory of quantum backflow.The first is to demonstrate the effect experimentally with particles.This would be an observation of a new non-classical effect.To this end an 'experimentally friendly' description of backflow has been provided by Miller et al [21].In this work the authors attempted to look at backflow when both positive and negative momentum eigenstates were present.This approach was commented on by Barbier and Goussev [22], who showed that this formulation and the standard one are not necessarily equivalent.Backflow has been observed theoretically in a PT symmetric ring which implies that it could be observed experimentally in such systems [23].The main theoretical challenge is to find an analytic expression for the quantum backflow constant c bf .Quantum backflow has been known for over 25 years and no such expression has been found.It is clear that a new approach is required.Bracken [24] has thought about this and tried to extend our view of backflow.The project reported here is also in this category.Previously backflow has been defined using wavepackets written as a sum of plane waves.As plane waves are eigenfunctions of the momentum operator this has the advantage of allowing us to define the momentum unambiguously and plane waves are mathematically easy to deal with, however they have led to an intractable integral equation, which has been the key barrier preventing us from discovering an analytic formula for the backflow constant.In this work we set up a model that is written in terms of the free particle hermite wavefunctions [25][26][27][28][29], rather than plane waves, to examine quantum backflow.In section II the model is described and the mathematical theory underlying the model is written down in section III.We then present the results of the calculations based upon this model and go on to discuss how they fit in with the bigger picture of quantum backflow in section V. Throughout this paper we retain constants in equations, but figures are presented in units where ÿ = 1 and m = 1/2.

The model
Here we describe the calculation and what we expect to achieve with it and in the following section we discuss the mathematics of the model.All calculations are one-dimensional.We set up a region of space defined by -∞ < x < L where L > 0. In this region we set up an initially localised wavepacket and allow it to evolve according to the Schrodinger equation.The wavepacket will be written in terms of the Hermite solutions to the free particle Schrödinger equation.Each of these solutions can be written as a sum of plane waves and it is only plane wave components with positive values of momentum that can leave this region of space by passing through the point x = L. Components with negative momentum travel towards x = -∞ and, hence, stay within the defined region forever.We will then consider the probability of finding the particle in the defined region as a function of time.
The quantum reentry problem has been discussed in detail by Goussev [10,11] and by Trillo et al [12] and co-workers and quantum reentry has been shown to have an equivalence with quantum backflow.The model described in this paper is similar to that model, but has two crucial differences.Firstly, in quantum reentry the initial wavepacket is completely delocalised while in the present model we do have a localised initial wavepacket.Secondly in this work the initial state has a symmetry which means that its net momentum is zero.
The current model is also similar to the original quantum backflow model.This is unsurprising because Goussev showed a direct equivalence between the quantum backflow and the quantum reentry problems.Our model has three distinct differences from the standard backflow model: (i) The initial wavefunction can, and does, contain plane wave components with negative momentum.These cannot leave the defined region, but can interfere with themselves and with the positive momentum components; (ii) The original work defined the initial wavepacket as the integral over momentum of a plane wave multiplied by a momentum dependent envelope function.Here we define the wavepacket as a sum over discrete functions; (iii) The difference defined in (ii) means that we end up with a set of simultaneous equations that can be solved numerically.Bracken and Melloy ended up with an integral equation that cannot be solved analytically, but which can also be solved numerically.We expect to find that, even though part of the wavepacket, and hence the probability density, leaves the defined region, the probability of finding the particle within that region can increase.This is an interference effect and can be regarded as being due to the wavepacket spreading backwards faster than it moves forwards.This is the same as occurs in standard backflow calculations, but here we do not exclude negative momentum contributions.As this is a different physical situation to the standard one we might expect to find a modified value for the backflow constant for this case.

Theory
The time-dependent free particle Schrödinger equation for a wavefunction function f is and a solution of this, first discovered by Miller [25] and discussed by Guerrero and co-workers [26,27] and Strange [28,29] is Here the symbols have their conventional meaning and τ is a positive constant with the dimensions of time which simply sets a time scale.Throughout this work we have chosen τ = 1.We can create a wavepacket composed of a linear combination of such solutions where the c n are coefficients to be determined and we want to evaluate the current density The derivatives here can be evaluated straightforwardly using simple properties of the Hermite polynomials or mathematical software such as Mathematica.Now, we want to calculate the quantity where P T is the probability of finding the particle in the region − ∞ < x < L at time T, and T T max min

>
. The wavefunctions described by equation (2) are localised around x = 0 at t = 0 for all values of the quantum number n.This means that the probability density ρ(x, t) exists in a finite region of space at all times.We choose L such that at time t = 0 the wavepacket is exponentially small above x = L.So where we have used the conservation of probability in one dimension and the fact that the current density at x = -∞ is always zero.To find the maximum value of η we have to maximise the time integral of the current density at the point L subject to the normalisation which has the property It is straightforward to show that λ = η and the task now is to identify the largest eigenvalue of C for differing values of !T and the corresponding eigenfunctions, the c n coefficients.

Implementation
Solution of the theory above has been accomplished successfully making use of the mathematical software packages Mathematica and Matlab.A first problem is that the size of the matrices K and C is infinite.In practice this cannot be realised of course and so we limit the size by setting the upper limit on the sums in equation ( 9) to a maximum value n max .This limits the accuracy, but we have done these calculations for increasing values of n max and have been able to extrapolate the results to give an estimate of the maximum value of the backflow constant.
By setting m = 1/2, ÿ = 1 and τ = 1 we have defined our units.The next problem was how to define L accurately.The wavefunctions of equation (2) are never zero, but become exponentially small beyond |x| = L.The choice of L contains some arbitrariness, but we selected L correct to four decimal places such that at t = 0 x dx , 0 0.99995 13 is the probability density associated with the wavefunction of equation ( 2) with n n max = .This definition sets a limit on the accuracy of our calculations, but we can be more restrictive later if required.Once this criterion had been established we found that the dimensionless quantity is constant.If we change our units, i.e. changes the values of m, τ or ÿ then L changes to compensate and y remains constant for a given value of n max .The reason for the square root in equation ( 14) will become apparent later.As we increase n max we increase L and decrease the dimensionless quantity y.For a given n max , y is constant so increasing L has an equivalence to decreasing τ, the unit of time.This suggests that the wavepackets created with larger values of n max will take longer (in units of τ) to leave the designated region than those with lower n max .
The main problem in the implementation of the theory above is that we don't know T min and T max .These were searched for systematically.A value of T min was chosen and we tried increasing values of T max until the current at L had decayed to near zero.T min was then increased and the process was repeated until T min was a time at which the current at x = L was very small.Initially we changed both T min and T max in units of 0.1τ and found the largest eigenvalue of the matrix C. We then refined the values of T min and T max so that, to an accuracy of four decimal places, they gave the maximum eigenvalue.
The time integrals in equation (10) were done using both Simpson's rule and the trapezoidal rule which gave identical answers to a greater accuracy than the other uncertainties in the calculation.
As n max increases the calculations take longer.This was the final limitation on the calculation.For the larger values of n max the matrix contains thousands of elements each of which requires a time integral and we then have to find the eigenvalues and eigenvectors of an n n max max ´matrix and this has to be done thousands of times.
The calculation for n 100 max = took about one week to perform on a standard IMac.

Results
The results of the calculations are summarised in table 1.We can see that, as expected, the backflow increases with increasing n max , or equivalently, decreasing y, but it increases fairly slowly.The value of backflow shown is the maximum value found for any value of T T T max min  = for the given values of y.We see that as n max increases the value of L increases slowly and hence y decreases, there is a general trend that !T increases as the maximum backflow increases.For y < 0.1 the maximum eigenvalue increases smoothly and nearly linearly with decreasing y.This is shown in explicitly figure 1.The straight line is the reason for the square root in equation (14).If we fit a linear polynomial to it the curve crosses the vertical axis (corresponding to n max = ¥) at 0.0241 max l = .If, instead, we fit a quadratic polynomial the value we extrapolate to at n max = ¥ is 0.0244 max l = .These values do assume linear behaviour of the backflow constant as a function of y.While this is the best we can do, it is possible that η could be a nonlinear function of y.
Figure 2 is a plot of the probability of finding the particle described by the wavepacket within the region -∞ < x < L as a function of time.We have plotted this for two different values of n max .These graphs are similar and seem to be characteristic of the states with high backflow.There is a rapid decrease in P(t) followed by a long slow increase, representing backflow, and then the probability density drops towards a value of P(t → ∞ ) = 1/ 2. These lines start at P(t) = 1 at t = 0 showing that the particle is certainly at a position x < L at t = 0.

Discussion
Here we put the calculation we have done within the context of the theory of quantum backflow and quantum reentry more widely.Firstly we note that the reentry depends on the quantity y and, despite being a purely quantum phenomenon does not depend on Planck's constant or the mass of the particle or the constant τ.
In figure 2 we plot P(t) for the cases n 50 max = and n 70 max = .Clearly the curves are similar.We see that there is a substantial period of time in which the probability of finding the particle within the region increases, even though probability density can only flow out of the region.When n 50 max = we have attained the basic shape that has maximum backflow.When n 50 max > the wavepacket has more degrees of freedom and  backflow can be supported for a longer time, but the basic shape of the P(t) curve remains the same.This is illustrated with the n 70 max = curve in figure 2. Figure 2 also illustrates the fact that wavepackets created with a higher value of n max take longer to pass the point x = L.This was discussed in relation to equation (14).
In this work we have performed a calculation which is, in principle, a hybrid of the original work of Bracken and Melloy and the quantum reentry work of Goussev [10].The wavefunctions described by equation (2) can be written as sums over plane waves.If this is done the symmetry of the wavefunction shows that for every positive momentum plane wave there is an equivalent negative momentum plane wave.i.e. 50% of the wavepacket has positive momentum and 50% has negative momentum.So, as the wave evolves 50% moves towards x = -∞ and, thus, stays in the region -∞ < x < L forever.However 50% moves in the positive x direction and will eventually pass the point x = L and leave the region completely.Recall that quantum backflow/reentry are interference effects, so one might expect that the inclusion of the negative momentum states gives us an extra degree of freedom which would allow greater backflow than in the case where only positive momentum states are allowed.However the symmetry of the wavefunctions imposes a strong limitation on the present work.We cannot produce wavepackets with more momentum in one direction than the other.This restriction reduces our capacity to maximise the backflow and our final result, assuming linear or quadratic behaviour, was 0.02425 0.0002 max l =  .While it may well not be significant, we note that this is very close to 2c bf /π.This work implies the following.The Bracken-Melloy constant is a measure of the maximum backflow that can occur when the states involved are restricted to being composed of plane waves with positive momentum only.That is the original and probably the most well-defined backflow.In this work we have attempted to calculate the maximum backflow with different basis functions.These contain equal amounts of both positive and negative momentum eigenstates and this yields a backflow constant of lower magnitude.It seems likely that there is a spectrum of backflow constants depending on the relationship between the positive and negative momentum components of the wavepacket.The Bracken-Melloy constant is the appropriate measure when only positive, or only negative, components of momentum are present.The value found in this work is the value that is relevant when the nett momentum is zero.

3 )F*
and (4) enable us to develop an explicit expression for the current at point L, however the time integral in equation(7) can only be performed numerically.The maximum of equation (7) can be found using the method of Lagrange multipliers.To this end we define where λ is the Lagrange multiplier and

( 10 )
obeys time reversal symmetry.Implementing the Lagrange multiplier method results in a matrix C whose eigenvalues are λ.The elements of C are

Figure 1 .
Figure 1.Maximum eigenvalue of the backflow matrix C plotted against y, the crosses are the data points and the line is a linear fit to the points.

Table 1 .
Table of results.The first column measure the number of wavefunctions of equation (2) in our sum.The second column is the value of L used.Column 3 is the value of y corresponding to n max .Columns 4 and 5 give T min and T max respectively and the final column is the value of the backflow constant for these values.