Nonclassical trajectories in head-on collisions

Rutherford scattering is usually described by treating the projectile either classically or as quantum mechanical plane waves. Here we treat them as wave packets and study their head-on collisions with the stationary target nuclei. We simulate the quantum dynamics of this one-dimensional system and study deviations of the average quantum solution from the classical one. These deviations are traced back to the convexity properties of Coulomb potential. Finally, we sketch how these theoretical findings could be tested in experiments looking for the onset of nuclear reactions.


Introduction
The seminal theoretical discussion of scattering angles in the Rutherford experiment is based on the Coulomb interaction between point charges modeling the alpha projectiles and the stationary gold nuclei [1,2]. In the conventional quantum approach, the projectiles are described as incident plane waves, and the collision is studied in asymptotic limits under suitable approximations [3]. A fuller approach that we pursue here is to compute the time dependence of the quantum evolution. With this in mind we study the simplest case of Rutherford scattering, i.e., the head-on collision.
As in the original discussion we model the nuclei as stationary sources of the Coulomb potential, but in contradistinction we describe the alpha particles by incident Gaussian wave packets. This leads to several predictions that differ from their classical counterparts. We focus on the average quantum behavior and show that it does not recover the classical solutions. In particular, the quantum projectile does not approach the target as close as its classical counterpart, the quantum dynamics is not symmetric about the time of collision, and the expected quantum trajectories do not overlap with the classical ones even in the asymptotic limits after the collision.
We trace back these discrepancies to the Ehrenfest theorem and emphasize that the average quantum dynamics recovers classical solutions only for potentials that are at most quadratic in the position operator [4,5,6]. We derive inequalities between average quantum and classical forces that hold for any cubic potential as well as for potentials with fixed convexity properties. This is directly applicable to Coulomb or gravitational forces and clearly shows that average quantum dynamics are different from their classical counterparts in a plethora of physically interesting scenarios.
Finally, we briefly discuss tunneling through the Coulomb barrier to infer the distance of closest approach and its dependence on the initial spread of the incident Gaussian wave packet. The tunneling probability is an experimentally accessible parameter as the particles that have crossed the barrier give rise to nuclear reactions [7,8]. In addition to computing the probability in the dynamical quantum model we also show the conditions under which the Wentzel-Kramers-Brillouin formula accurately approximates it [9]. In particular, we note that the latter may be orders of magnitude off for low energy projectiles even for a negligible momentum dispersion (see also Refs. [10,11,12,13]).

Ehrenfest theorem and classical limit
Before moving to the collisions we would like to present a general discussion on a comparison between the average Schrödinger dynamics and the classical one, especially that some textbooks give an incorrect statement that the average quantum trajectory recovers the classical motion [14]. This comes in relation to the Ehrenfest theorem showing that the quantum equations of motion for average position and momentum: have the same general form as the Hamilton's equations of classical mechanics: where x is the position, p is the momentum, m is the mass of the particle and V stands for the gradient of potential to which we will refer as (negative) force. In what follows, we only consider the one-dimensional motion. While the general form of (1) and (2) is the same, they become identical only if the average of the force equals the force for the average position [5]: For arbitrary wave functions, this condition is satisfied only for potentials that are at most quadratic in x. Already a cubic potential shows that there is a consistent difference between the quantum and the classical forces. Namely, the derivative of the potential . Similar inequalities follow for derivatives of convex or concave potentials using Jensen's inequalities [15]. Many potentials of natural interactions have welldefined convexity properties, e.g., the repulsive Coulomb potential we study here. In this case, the wave packet's average quantum force can be very different from what a classical particle experiences in the same potential. This leads to consistent differences between the average quantum trajectory and the classical one.

Classical trajectories
For comparison, we first show how to compute the classical trajectories. Consider a projectile, initially at x = x 0 (negative), is shot towards the target, located at x = 0. We follow the tradition that projectiles are propagating from left to right. If the particle's initial kinetic energy is T 0 , the law of conservation of energy implies T (x) + V (x) = T 0 + V (x 0 ) ≡ E 0 , with the Coulomb potential given by where Z P and Z T are the atomic numbers of the projectile and the target, respectively, and α is the fine structure constant. Accordingly, the classical equation of motion is where d cl = Z P Z T α c/E 0 is the classical distance of closest approach. The above equation can be integrated and gives the following trajectories: where τ cl denotes the time of classical collision: Since Eq. (6) is transcendental, we solve for x(t) numerically.

Quantum trajectories
The quantum setup is shown schematically in Fig. 1. We describe here the methods used to obtain the average quantum dynamics of the projectile. In non-relativistic quantum theory, the properties of a particle are described by a wave function ψ(x, t), which evolves in accordance with the time-dependent Schrödinger equation. For time-independent interactions, one can alternatively use the evolution operator U to write where is the Hamiltonian. Any truncation in a series expansion of U leads to loss of unitarity, and consequently a change in the norm of wave function over time. We circumvent this problem and calculate the time-development by implementing the Cayley's form of evolution operator [16,17,18,19]. Since the wave function is expanding continuously during the collision event, the simulation becomes numerically demanding as time progresses. To circumvent this we have developed a method to dynamically allocate the size of discretization of space, details of which are given in Appendix A.
Finally, in the classical picture the projectile particles are modelled as point charges shot towards the target from a distance L with a kinetic energy T 0 . In our quantum simulations, we model the projectiles with Gaussian wave packets centered at the same distance L, with a width σ and an average momentum P = √ 2mT 0 . Therefore, the wave function describing the quantum state of the projectile at t = 0 is given by (x 0 = −L, We now describe several parameters that show differences between the classical and the average quantum dynamics.

Comparison of trajectories
We follow the original Rutherford scattering experiment and study the collision of alpha particles with gold nuclei. Hence we set Z P = 2 and Z T = 79 throughout this work. Typical classical and quantum trajectories are presented in Fig. 2a. In both cases L = 10 pm and T 0 = 5 MeV, i.e., the typical energy of alpha particles in the original Rutherford scattering experiment. We now systematically discuss the following quantities: distance of closest approach, the effects of a finite spread of the projectile's initial wave function, the time of collision, and the asymptotic behavior of trajectories.

Distance of closest approach
As seen in Fig. 2a the quantum projectiles are, on average, reflected from a bigger distance to the target. This can be intuitively understood by invoking either convexity of the Coulomb interaction or the Heisenberg uncertainty principle. In the latter case, note that a quantum particle cannot be stopped completely. There is always some momentum dispersion which leads to a non-zero kinetic energy. It follows that the maximal potential energy cannot be as big as in the classical case, and hence the distance of the closest approach increases. In the former case, Jensen's inequality at time t = 0 allows us to write: Accordingly, the average quantum mechanical force is stronger (more repulsive) than the force experienced by the classical particle. As a result, the quantum projectile moves slower than its classical counterpart at the very beginning of its journey. We emphasize that this inequality holds at the beginning of the evolution, and at later times it may reverse (as we will show later). Nevertheless, the quantum projectile never gets as close to the target as the classical one. In the next section we derive conditions for the quantum distance of closest approach to differ from its classical counterpart minimally. Fig.2b shows the wave function at different times for the case of σ = 100 fm. As the wave packet approaches closer to the target, the leading edge of the wave packet is reflected (a) Classical and quantum trajectories. The right panel shows a magnified view of the scattering event. The quantum trajectories are the results from the dynamical simulations.
Snapshots of the wave packet at different times for the case of σ = 100 fm. Figure 2: Time-dependent properties of the alpha particle wave packets shot in the Coulomb potential of a gold nucleus fixed at the origin. The projectiles are shot from the distance L = 10 pm with the initial width of Gaussian wave packet σ and average momentum P = √ 2mT 0 , where T 0 = 5 MeV, and m is the mass of an alpha particle. x cl represents the classical path, and x is the expected position of the wave packet. Numerical details are given in Appendix B.
back. This interferes with the trailing incident part to create rapid oscillations, which eventually die out as the wave packet travels back attaining a near-Gaussian shape [20,21,22]. As in the case of collision with a hard wall, the qualitative similarity between x and x cl is not evident when one looks at the wave packet itself [23].

Quantum distance of closest approach
We define the quantum distance of closest approach as the smallest average position to the target: Figure 3: Difference between the classical and quantum distances of closest approach. The numerically found minima of the difference are marked with dots for various launching distances. In the main text we argue that they are achieved for the initial width of the projectile given by σ 0 = L/2P , where P is the momentum of a classical alpha particle with kinetic energy T 0 . Numerical details are given in Appendix B.
As seen in Fig. 2a this quantity depends on the initial spread of the wave function of the projectile. We now give physical arguments which determine the optimal initial spread σ 0 , for which the difference d qm − d cl is the smallest.
Since the projectile is launched from a large distance L d cl , the time of collision satisfies τ cl ≈ mL/P . We approximate the position spread before the collision by the value obtained for a free quantum evolution: We established in Eq. (10) that the quantum wave packet feels a stronger force compared to its classical counterpart. As an implication the wave packet, on average, is reflected from a farther distance, i.e., d qm > d cl . For a given potential, the quantum averages at a given time are closer to their classical values for wave functions that are more and more concentrated in space. Therefore the difference d qm − d cl is smaller for narrow wave functions, and it is minimal when σ(τ cl ) is minimal. Eq. (12) suggests that this happens for the initial spread given by and corresponds to the spread at the collision time σ(τ cl ) = √ 2 σ 0 . Fig. 3 shows that this is in an excellent agreement with the numerically obtained minimas in d qm .
Furthermore, we establish the range accessible to the quantum distance of closest approach for the energies considered in this work. The left panel of Fig. 4 shows the position spread of wave packets launched from various distances. The collisions correspond to the sudden drops of the position spread, i.e., the wave packet gets compressed when the projectile approaches the target. The distance of closest approach matches the minima in these curves, which are all smaller than the initial spread of Eq. (13), i.e, the uncertainty in d qm is less than σ 0 . Next, we show in the right panel of Fig. 4 that the classical distance of closest approach is within the range of one standard deviation from the quantum distance of closest approach. Finally, this translates to

Origin of the optimal spread
Mathematically, Eq. (3) encoding the equivalence between the classical and Ehrenfest dynamics, is satisfied for any potential if the wave function is given by the Dirac delta. The more the wave function is concentrated in space, the more similar classical and average quantum positions are. Physically, however, due to the Heisenberg uncertainty principle, there is no evolution which preserves the Dirac delta, and we show that this leads to the optimal initial spread. For the values presented in Fig. 2a, the optimal spread is about 50 fm (see green dotted line). If the initial spread is larger (compare with the red dashed line), the wave function spreads further and the deviation from the classical trajectory is higher. But for smaller initial spreads (compare with the blue dashed-dotted line), the free quantum evolution dominating initially predicts faster spreading due to large initial momentum dispersion. Therefore, the standard deviation becomes larger by the time the wave packet reaches close to the target, and the deviation from the classical trajectory is again higher.

Time asymmetry of quantum collision
As already seen in Fig. 2a  A distinguishing feature of the quantum trajectory is its asymmetry with respect to the collision time. Recall that the classical trajectory is symmetric in time, i.e., x cl (τ cl − ∆t) = x cl (τ cl + ∆t), and similarly for momentum p cl (τ cl − ∆t) = −p cl (τ cl + ∆t). In particular, the projectile returns back to its original launch distance exactly at time T = 2 τ cl . Our calculations show that the quantum projectile makes a collision at a different time τ qm (> τ cl ), and the evolution for t > τ qm is not a mirror image of evolution for t < τ qm . The return journey of a quantum particle takes a longer time than the onward one, and the wave packet returns at its launch position at time T > 2 τ qm (> 2 τ cl ). For the case in Fig. 5 (σ = 159.85 fm), the time taken for the return journey is longer by ≈ 1.3 × 10 −22 s. Such effects are more prominent for wave functions with a larger momentum variance, e.g., the time difference is ≈ 1.1×10 −21 s for the case of σ = 20 fm in Fig. 2. The rate of increase of position spread itself increases with variance in the initial momentum, and given that the potential has a well-defined convexity, this tends to amplify the asymmetries over time.
The discussed asymmetry in quantum trajectories is yet another manifestation of the convexity in Coulomb potential. Consider a point at distance X from the target. Irrespective of whether it is traveling towards or away from the target, a classical particle feels the exact same force when it passes through that point, i.e. −V (X). Accordingly, the Hamilton's solution of Eq. (2) has to be symmetric about τ cl . The quantum case is much more interesting; Fig. 4 and 5 show that the position spread immediately after the collision is larger than what it was before. This means that the width of the wave packet, when it passes through point X, is larger during its return journey. As a consequence, the average force experienced by the quantum projectile at that point, i.e. − V (X) , is different in the onward and the return journeys. Accordingly the Ehrenfest's solution of Eq. (1) has to be asymmetric about τ qm . For the parameters used in this work, the classical trajectory is not even statistically within the quantum prediction, i.e. x cl and p cl are outside x ± ∆x and p ± ∆p (see Fig. 5).

Classical and quantum force
We started our discussion by showing that initially the quantum projectile is repelled more than the corresponding classical counterpart. For example, alpha particle of T 0 = 5 MeV launched from L = 50 pm feels an average force F ≈ 1.000031F cl at t = 0. Since for later times x < x cl (negative, see Fig. 2a) we obtain: The second inequality starts dominating closer to the target, and during the collision the classical force far exceeds the quantum one as seen in Fig. 5.

Proposed experiments
Having established the differences between the classical and average quantum dynamics, we now move to a possibility of their experimental verification. The basic idea is to look for the onset of nuclear reactions. We model the potential in the vicinity of the nucleus by the Coulomb barrier given in Fig. 6 and defined by: In our simulations we varied the launching distance about 50 pm, and observe that the results are practically independent of L. Under the classical model of the collision, if the projectile can cross the barrier it can be observed as a nuclear reaction. In such a case, the classical distance of closest approach is smaller than the size where the barrier is truncated. We take as a wellestablished fact that alpha particles cannot be prepared with only one momentum value and repeated momentum measurement on the projectile before the collision returns a normal distribution. For this reason, within the classical model, we deal with an ensemble of projectiles with Gaussian momentum probability density. Accordingly, they cross the Coulomb barrier only probabilistically with the corresponding probability P cl plotted by the dashed line in Fig. 7. Conversely, from the measured probability of nuclear reaction, one can estimate the variance in position and momentum of the classical ensemble and compute the mean distance of closest approach.
Within the quantum model, the initial randomness is encoded by the Gaussian wave function. We compute its evolution in the presence of the Coulomb barrier and determine the probability that the particle tunnels through: The result is plotted as a solid curve in the left panel of Fig. 7 and it is considerably different from the classical result for initial wave functions with position spread bigger than about 10 fm. For smaller position spreads the quantum and classical curves are approaching each other. Finally, for Dirac delta position distribution, in both cases we have flat momentum distribution, and hence half of the particles tunnel through. The range of very small position spreads is not accessible in our numerical calculations because the wave function spreads very fast such that the dynamics of approaching the barrier already consumes all of the computational resources. We again emphasize that the nuclear reaction cross-section directly translates to the spread of the projectile's initial wave function. This provides an interesting way of experimentally determining this spread. Note also that simultaneous measurement of the tunneling probability and the initial position or momentum spread is capable of excluding the classical model. For example, at σ = 10 fm the classical probability is more than three orders of magnitude smaller than the quantum one, which is at the measurable level of 10 −3 .

Validity of the WKB approximation
Quantum tunneling is of course a well-studied phenomenon, with many interesting applications even in astrophysics [7,8]. Traditionally, the projectile is assumed as an incident plane wave and the tunneling probability is approximated within the time-independent theory. One such treatment is the Wentzel-Kramers-Brillouin (WKB) approximation [24,25].
We now provide conditions under which this approximation matches the time-dependent approach presented here. Some advancements have already been made in the limiting case of low-energy projectiles with a negligible momentum variance [13], and it is expected that many peculiar effects arise with a large momentum spread [9,10,11,12,13]. In the left panel of Fig. 7 we plotted the tunneling probability in the WKB approximation by the horizontal dashed line at the bottom. The WKB result tends to the dynamical one in the limit of small initial momentum spread (large position spread). However, it does not exactly approach the tunneling probability obtained in the dynamical simulation. The right panel of Fig. 7 compares the two results for different kinetic energies. One can draw a boundary between the range where the two results are comparable and where they are very different. It turns out that it is given by the optimal σ 0 we have derived above.

Prospects and limitations
An interesting regime of Rutherford scattering is when the projectiles are moving at relativistic speeds, e.g., the photons moving through a medium with refractive index ∼ 1/x [26]. While the non-relativistic wave packets shown in Fig. 3 of the present work do not recover the classical solution, in the Rutherford scattering of photons the classical limit is achieved [26]. This is a consequence of the non-dispersive behaviour of photonic systems. For non-relativistic particles the variance in momentum directly translates to a variance in the velocity, which leads to a growing position spread in time. Unlike alpha particles, photonic systems are non-dispersive as the variance in momentum does not correspond to a variance in the velocity [27]. Accordingly, apart from the negligible impact of refractive index, there is no change in the position spread of a photonic wave function. One can therefore have an extremely narrow wave packet, and hence particle-like dynamics, all along the scattering event. This explains the recovery of classical solutions in the photonic Rutherford scattering [26].
The time-dependent scattering problems, like the one discussed in this work, involve the notions of positions and trajectories. While one can talk about the expected positions for a wave packet, no such analogous counterpart exists for individual plane waves [28]. Hence a consistent theory can only be built on the dynamics of wave functions with finite position spreads [26,29].
In this work we represented the projectile as an incident Gaussian wave packet and the target with a simplistic Coulomb potential. Even though this works well in the considered energy range, this is an approximation on a few fronts. For one thing the target should also be considered as a wave packet and the scattering must be understood from the dynamics of two colliding wave packets [30]. Secondly the internal structures of the projectile and the target will give rise to various exchange effects [31,32,33,34]. These effects will be more prominent at higher energies and must be accounted for by an inclusion of exchange terms in the Hamiltonian [35,36,37]. An even fuller approach would be to replace the Coulomb potential with more realistic interactions [38,39,40].
Finally, if we consider a two or three-dimensional setting, different parts of the wave packet will acquire different phases as time passes. This will give rise to a diffraction pattern in the angular distribution of the total amplitude [41]. An extension of our work in two-dimensions promises to calculate the possible dependence of this pattern on the initial spread of the projectile wave function.

Conclusions
Average quantum dynamics in potentials with well-defined convexity properties does not approach its classical counterpart. We have shown this for the case of Coulomb interaction between the alpha particles and the stationary gold nuclei in a head-on collision. The differences include the distance of closest approach, time of collision, and time symmetry of the dynamics. We sketched an experiment aimed at verifying these predictions. It could be rather challenging as we focused on head-on collisions in this work. It would be interesting to work out possible differences in the 2D setting and compute how the differential cross-sections depend on the spread of the initial wave function of the alpha particle. Such predictions might be easier to verify in a laboratory. Finally, the model could be extended to study fully quantum mechanically astrophysical nuclear reactions. elapsed in the collision events. This makes it challenging to study the scattering in a fully quantum mechanical time-dependent theory [42,43,44]. Within the scope of this work, the wave function remains confined near a single point; this allows us to perform a dynamic allocation of the grid and recast the problem as a sequential chain of small simulations. We start by defining a box around the initial launch distance such that the wave packet is highly concentrated at its center. The simulation is initiated, and the wave packet evolves in the Coulomb potential. We keep track of the peak and the width of the probability density function to make sure there is no reflection from the boundaries of the box. Once the tail of the probability distribution starts approaching either of the boundaries, we define a new box that re-confines the wave packet at its center. The final wave function from the old grid (box) is interpolated and used as an initial condition to start the calculation onto the new grid. The grid allocation and interpolation keep repeating until the wave packet returns at its original launch distance after the collision. This chain of small simulations is verified to precisely simulate the whole of the scattering event.

B Numerical details
• Numerical calculations are performed in natural units of c = 1. Accordingly, the conversion constant c = 197.3269804 MeV fm. We follow these units within this section.
• We have employed two different convergence tests to perform the error analysis: (i) the convergence of the quantum distance of closest approach, and (ii) the constancy of the expected Hamiltonian all along the scattering event. In Fig. 8 we present the error analysis assuming L = 10 pm and σ = 100 fm (the typical width of Gaussian wave packets considered in this work).    8a shows that there is a very good convergence in d qm for grid size 0.3 fm and time step 2 fm. The first panel of Fig. 8b shows the variation of the relative error in expected Hamiltonian, i.e. 1 − H(t) / H(0) , as a function of time for a fixed grid size and time step. For most times this error is negligible, except at t = τ qm where it attains its peak. The second panel shows this peak error is negligible for grid size 0.3 fm. The last panel shows that an increased time step has no significant impact on the error in energy. Since this was not the case for the error in d qm (see the second panel in Fig. 8a), we use two tests of convergence for reliable calculations.
• A grid size of 0.2 fm, with a time step of 1 fm, ensures a good precision in all of our calculations. We have accordingly set the same throughout this work.