# Space-time evolution of Dirac wave packets

###### Abstract

In this work we study the dynamics of free D relativistic Gaussian wave packets with different spin polarizations. We analyze the connection between the symmetry of initial state and the dynamical characteristics of moving particle. The corresponding solutions of Dirac equation having different types of symmetry were evaluated analytically and numerically and after that the electron probability densities, as well as, the spin densities were visualized. The average values of velocity of the packet center and the average spin were calculated analytically, and the parameters of transient Zitterbewegung in different directions were obtained. These results can be useful for the interpretation of future experiments with trapped ions.

###### pacs:

73.22.-f, 73.63.Fg, 78.67.Ch, 03.65.Pm## I Introduction

The Dirac equation belongs to the most important equations in modern physics. It predicts the existence of electron spin and magnetic moment, gives the natural description of the positron states, describes fine structure of the energy spectrum of hydrogen-like atoms. The quantized solutions of Dirac equation are considered to be a natural transition to quantum field theory. Moreover, the one-particle relativistic quantum mechanics describes the unexpected electron dynamics including Schrödinger’s Zitterbewegung (ZB)SchBarut () and Klein paradox.Klein () The trembling motion of relativistic particles or ZB is caused by the interference between positive and negative energy states which form the electron wave packet. The frequency of ZB is determined by the gap between two energy bands, and the amplitude of oscillation of the wave packet center is of the order of the Compton wavelength.

The results of the first experimental observation of ZB phenomena were published recently in the paper by Gerrisma et. al.GKZ () For the ZB simulation the experimentalists used a linear Paul trap where ion motion can be described by one-dimensional Dirac equation.Lamata () The authors of Ref.[3] study the motion of ion and determined its position as a function of time for different initial conditions. As was shown in Ref.[4] the solution of the Dirac equation can also be simulated using a single trapped ion with four ionic internal states. In this case the ion position and momentum are associated with respective characteristics of D Dirac particle.

The dynamics of relativistic one-dimensional wave packet for the first time was investigated numerically by Thaller.Thal () He visualized the initially localized solutions of single particle Dirac equation and observed the trembling motion of the wave packet centers as well as some other phenomena which are caused by the interference of positive- and negative-energy states. J. Lock found that ZB oscillations of localized initial states have transient rather than sustained character.Lock () Some interesting examples of the relativistic dynamics of D electron wave packets were presented in Ref.[7]. It should be stressed that ZB oscillations exist only in the one-particle relativistic theory. Krekora, Su, and Grobe demonstrated analytically and numericallyKrekora () that quantum field theory does not permit any Zitterbewegung of real electrons.

The oscillatory ZB motion of D electron wave packets in crystalline solids for the first time was predicted in Ref. [9]. This phenomenon has been considered for 2D electron gas with Rashba spin-orbit coupling in Ref. [10, 11], in narrow gap semiconductors in Ref. [12], in carbon nanotubesZaw72 (), in single and bilayer grapheneKat (); RusZaw (); MDF () and also in superconductorsLur ().

In the present work we investigate the relativistic dynamics of D wave packets. In Sec. II the general properties of symmetry of Dirac equation solutions which determine the dynamics of wave packet are analyzed. After that, in Sec.III we consider the evolution of initial states with different initial symmetry. The electron probability densities evaluated analytically and numerically are visualized. The directions of average velocity of wave packets as well as the phenomenon of high-frequency Zitterbewegung for different initial polarization are determined. In Sec IV the symmetry and structure of spin densities for relativistic packets are analyzed. Finally, in Sec. V, we conclude with the discussion of results. Some auxiliary results are found in Appendixes A and B.

## Ii The symmetry of the relativistic wave packets

In this section we shall first study some general symmetry properties of the electron wave function which determine the dynamics of wave packets. So, we start with the famous Dirac equation for a four-component wave function

where

is the light velocity, is the electron mass and

The four independent free-particle solutions of Eq.(1) for given momentum and energy , can be written in the form

where , , , and are the free Dirac spinorsSak ()

Here we use the normalization condition

which means

Let us consider now the symmetry of the solutions of Dirac equation with respect to the space reflection. Suppose that initial wave packet is symmetric: . As it known this symmetry is conserved, i.e. only if the full function including the spinor part is the eigenfunction of the parity operator , where the matrix is determined by Eq.(3) and is the inversion operator for : . Also the Dirac equation is invariant under the replacement

where , () is the inversion operator for ()- component, and is the corresponding component of spin operator . So, if the wave function satisfies (in particular at ) the relation

then the parity in -plane is conserved,

Besides, we may determine the reflection transform in Hilbert space as the operator which, just as the operators and , commutes with Dirac Hamiltonian (2). If the initial wave packet has the certain party

then satisfies this equation for all times that leads to the conservation of the symmetry

Similarly we may introduce the operators and connected with the reflection transforms and correspondingly.

Thus the symmetry of the solution with respect to full () or partial space symmetry depends not only on the space symmetry of the initial wave function, but also on the ratio between its components. This statement is valid with respect to other types of symmetry. Below we are interested in the dynamics of the initial wave packet of the form

where are the complex numbers and the space part satisfies the normalization condition

Specifically we suppose that the probability density at is spherically or axially symmetric. So let has the form

where , , are cylindrical coordinates and is an integer. The considered wave packet remains an axial symmetrical (in -plane) for all the times if its spinor part is one of the eigenfunctions of operator

Indeed it is readily to show that the function

is the eigenstate of -component of total angular momentum operator

where .

For time the general expression for the wave function is

But is a conserved quantity so that the obeys Eq.(20) too. Solving this equation we find the -dependence of components

which immediately leads to conclusion that the probability density is axially symmetric in -plane.

Note that the initial state is the eigenfunction of too with an eigenvalue . Notice also that as would be shown below the spherical symmetry of the initial wave packet, Eq. (19) (with and , ) turns into axial one. The above statements concerning the wave packet symmetry are justified by our analytical and numerical calculations for different initial conditions.

## Iii Dirac wave packets dynamics

i) The propagation of cylindrically symmetric wave packet

Now we describe some peculiarities of the striking kinematics of three-dimensional relativistic wave packets. As a first example let us consider the time evolution of the initially localized wave function of the form

As was shown above (see Eq.(18)) for such polarization of the wave packet the probability density conserves its axial symmetry if the space part of the wave function (24) has the form (17).

The wave function can be expanded in plane waves

In the momentum space the components of the bispinor wave function

are given as a linear superposition of the positive- and negative-energy solutions (5) and (6)

where and . Here is to be determined by the Fourier expansion of . Straightforward calculation using (5), (6) and (24) gives

Here is the Fourier transform of the function . To obtain the components of wave function we substitute Eqs. (28)-(30) into Eq. (25). Then for the axially symmetric initial wave packet , where , i.e. , , after integrating over the angular variable we shall have

where , are Bessel functions. Using these expressions one may find the dependence of the electron probability density on coordinates and time for the wave packet (24). To illustrate our results in this Section and everywhere below we choose the function in the Gaussian form

where and determine the width of wave packet in the -plane and in direction correspondingly, and is an average momentum of the wave packet.

In Fig.1(a) we plot the electron probability density at for Gaussian wave packet with initial momentum and , at time . Here and below all distances are measured in units of the Compton wavelength , the time is in units of . This result was obtained by using the numerical method (”leap-frog” algorithm) described in Appendix A. Just as we expected Fig.1(a) demonstrates the axial symmetry of the considered distribution in -plane; the probability density has a form of cylindrical wave propagating from the point of origin and having some maxima.

To analyze the motion of the packet we have to find the average value of velocity of the packet center. In the momentum representation

Substituting Eqs.(28)-(30) into Eq.(35) and using the expressions for matrices (Eq.(3)) we obtain

Since the Fourier transform of Gaussian wave packet (34) has an axial symmetry in - plane

the components of velocity as it follows from Eqs.(36),(38). Otherwise, owing to the axial symmetry of spacial distribution of the electron density, the average coordinates of packet . As a result, mean components of velocity in -plane are equal to zero.

The motion of the packet center in -direction experiences rapid oscillations commonly known as Zitterbewegung (the second term in square brackets in Eq.(38)). Besides, the wave packet displaces slowly with constant velocity (the first term in square brackets in Eq.(38); see also Eq.(B.8)) even if its momentum .

The existence of the constant component of average velocity in this case can be understood from the relation between velocity and momentum depending on the sign of energy: for the wave packet consisting of the states with positive (negative) energy a positive momentum corresponds to a positive (negative) velocity . Let us represent the time-independent probability density in the momentum space as a superposition of the positive energy part and the negative energy part

where , . Using Eqs.(28)-(30) one may readily show that

Fig.2 shows the dependence (in arbitrary units) on (in units of ) for the Gaussian initial wave packet, Eqs. (24), (39) for (a) and for (in units of ) (b).

We see that for the case the momentum distribution for positive- and negative- energy components is shifted towards positive and negative momentum , respectively. But in both cases such dependence leads to the motion of each parts (and the whole wave packet) with a positive velocity along -axis. We also see in Fig. 2(a) that the function has essential overlap with the function that is a necessary condition for existence of Zitterbewegung of the packet center (see Eq.(38)). The value of constant component of velocity depends on the ratio between the initial width of the wave packet and the Compton wavelength. In the limiting case is much less than the light velocity: . In particular, for the symmetrical wave packet () of width .

Let now the initial wave function describes the wave packet moving along -axis with average momentum . Then the distribution of the full electron density in -plane is similar to one shown in Fig.1(a). However the essential difference appears in the character of evolution of the wave packet in (or )-plane (Fig. 1(b)). The dependencies for the states with positive- and negative- energy parts for this case is shown in Fig. 2(b). We see that both components consist of positive momentum . For the smaller negative-energy parts, this corresponds to the motion with negative velocity along -axis. Thus in the position space the initial wave packet splits into two packets propagating in the opposite directions along axis.

ii) Asymmetrical wave packet evolution

We next consider another example of the initial spin polarization of the electron

where as before the function is determined by Eq.(34).

Note that the example under review is invariant with respect to the reflection transformation , i.e. the expression (42) satisfies Eq.(13). As was shown above this means that the probability density is an even function of at all times. Performing the same kind of calculations as for symmetrical wave packet we find the components of the initial wave function (42) at in the momentum space.

The components of can be obtained directly by the Fourier transform of Eqs.(43)-(45).

Using these expressions one may find the full electron density . Figure 3(a) shows the corresponding probability density distribution in the plane at time . The parameters of wave packet are: and , . Comparing Figs. 1(a) and 3(a) we see that the change of initial spin polarization leads to the fact that the wave packet being axially symmetric in -plane at loses its symmetry at .

The kinematics of the wave packet can be characterized by the average velocity of its center with components

As in the previous case the wave packet center drifts with constant velocity (the first term in square brackets in Eq.(49)), but now it is directed along axis. Besides, one can see from Eqs. (49) and (50) that the packet center performs damped oscillations (Fig. 4) along and directions. We also see from Eq.(51) that that is the result of the symmetry under the replacement , Eq.(14).

We now consider the behavior of the initial wave packet with nonzero momentum . Obviously such initial state is not an eigenfunction of the parity operator . Nevertheless, the probability density remain to be a symmetrical function relatively to the reflection transform . Indeed, one may check that the Dirac equation (1) is invariant under the transformation

So that if the initial wave function satisfies the equation

then, as one can check using Eqs.(46)-(48), this relationship is valid at and consequently Eq. (14) holds.

The character of the motion of the packet center in -plane is similar to the case . Namely, the center of the wave packet drifts along -direction with constant component of velocity and oscillates along and axis (Zitterbewegung). As it known the ZB is significant if the subpackets with positive and negative energy have overlap in the position space. But if the initial average momentum of the wave packet is nonzero, both subpackets move with the opposite velocities along axis that leads to their spatial separation. So, the amplitude of the ZB decreases more rapidly than for the case (compare Fig. 4 and Fig. 5). Notice that this result is valid for the narrow enough wave packets with width and (for very large the ZB oscillations are almost undamped). The constant component of packet center velocity also depends on the initial momentum . In fact, as it follows from Eq. (49) it decreases as increases.

## Iv Spin dynamics

At present we consider the average value of the spin operator for Dirac particle

or in the momentum representation

One can verify that the spin densities for the arbitrary wave function with components are

i) Cylindrically symmetric wave packet

The form of Eqs.(55)-(57) and the expressions (31)-(33) for the components of wave function show that for an axially symmetric wave packet, Eq.(24) the component of spin density is an axially symmetric too both for and (Fig. 6(c)). One can also verify that and components can be written in the form

where the expression for the function is enough cumbrous and will not be presented here. Thus the -plane projection of vector is directed along the unit vector of cylindrical system . The spin densities and are represented in Fig. 6(a),(b) for the Gaussian wave packet with the average momentum . We see that () is the antisymmetrical function of () and the spin density conserves its axial symmetry.

The average values of spin operators for this polarization can be found by using Eqs.(28)-(30) for wave function in the momentum representation and previous definition, Eq.(54b).

As it follows from these relations only is not equal to zero for the considered wave packet. The first term in square brackets of the last formula corresponds to the constant component of and the second one describes the typical transient Zitterbewegung.

ii) Asymmetrical wave packet

One can easy verify that the second example, Eq.(42) considered in our work corresponds to the initial spin density . Really, as it follows from Eqs.(46)-(48) the components of wave function for the packet with at obey the relations

that together with Eqs.(55)-(57) leads to the result .

The analysis of expressions (46)-(48) and (57) shows that in -plane the -component of spin density is an antisymmetric function of that is connected with discussed above symmetry of probability density with respect to the replacement .

Fig. 7 illustrates the distributions of and for Gaussian wave packet with and , at .

By inserting Eqs.(43)-(45) into Eq.(54b) we find for the average components of spin

Obviously for a symmetric wave function (that means ) these expressions give . Note that remains to be equal to zero also for . It may be shown that in this case Eqs.(62),(63) describe the spin ”precession” (in -plane about vector ) which has a transient character. Such phenomenon for a hole system, described by Luttinger model was discussed in the recent work of authors.DMF ()

## V Summary

In this work we have studied the quantum dynamics of relativistic particles represented by three-dimensional Gaussian wave packets with different initial spin polarizations, described by the Dirac equation. The analysis of the general symmetry properties of solutions of one-particle Dirac equation allows to predict the direction of average electron velocity as well as the direction of trembling motion. In particular, the evolution of spherically and cylindrically symmetric initial Gaussian wave packet demonstrates that the wave packet with initial polarization, which is determined by bispinor has cylindrical symmetry at all times , but the wave packet with initial polarization loses its cylindrical symmetry at time . The influence of the symmetry of initial wave packet on the distribution of spin densities is analyzed.

## Appendix A

The ”leap-frog” algorithmS () is applied in a spatial grid of bin-sizes , , and with time step :

The spatial derivatives are computed symmetrically. Reflecting boundary conditions are applied on a very large grid (running stops before reflections occur if necessary). We use the norm as the stability measure of the algorithm (1). In accordance with von Neumann stability analysisGKO () (for large component plane waves) the stability region (-spatial grid bin, -time step) is:

Thus, for a single precision calculation the loss of norm can be kept within - in a time step run. It should be noted that in the case of Zitterbewegung, i.e. of the spatial oscillation of the wave packet, one more condition have to be imposed to the lattice sizes:

The conditions (A.2) and (A.3) were fulfilled in all our calculations.

## Appendix B. Drift velocity for the arbitrary initial wave function

As was shown in previous investigations the motion of the Dirac wave packet center does not obey classical relativistic kinematics. In particular, besides the rapid oscillations (ZB) the wave packet can drift with constant velocity although its average momentum is zero. One can show that in this case the direction of such motion coincides with the direction of initial average velocity. In fact, for the second example considered in this work (Eq. (42)) only and . Therefore such initial spin polarization leads to the motion of the wave packet with constant velocity along axis (see Eqs. (21),(26)). It is easy to see that in the first example, Eq. (10), the wave packet at has the velocity directed along axis, so the motion at occurs in this direction.

Let now find the drift velocity of particle for the case of arbitrary initial polarization

where are the complex coefficients, , and is to be determined from the Fourier expansion of coordinate wave function . (We do not suppose that and have any symmetry).

At the wave function in momentum representation is

Using the Eqs (5),(6) and (8) for free Dirac spinors we find from Eq. (B.2)

The density of -component of velocity () in the momentum space is

Obviously the time-independent part of is defined as

where and corresponds to the contribution of positive and negative energy into Eq.(B.2). One may check that

So, using Eqs. (B.2), (B.5) and (B.6) we obtain

Substituting the expression for the coefficients , Eq.(B.3) into the Eq.(B.7), we find the constant velocity of wave packet center

It is convenient to represent this expression using the initial components of velocity , so that

(Note that in this equation there is no summation over ). Let the initial Gaussian wave packet be spherically symmetric, i.e. is determined by Eq.(39) with . Then if the average momentum the second term in Eq.(B.9) equals to zero and the value of integral in the first term does not depend on . Thus the direction of the constant velocity of wave packet coincides with the initial one. In the case as it follows from Eq.(B.9) , and the asymptotic direction of the average velocity is along -axis, i.e. along the average momentum . A similar result was obtained in our previous workMDF (), concerning the propagation of the wave packet in graphene.

Note that Eq.(A.8) is valid also for the most general initial wave packet of the form

if we rename in Eq.(B.8) . It is not difficult to check that the expressions for the constant components of velocity obtained earlier for the examples Eq.(24) and Eq.(42) follow also from the general equation (B.8).

## Acknowledgments

This work was supported by the Program of the Russian Ministry of Education and Science ”Development of scientific potential of High education” (Project No. 2.1.1.2686) and Grant of Russian Foundation for Basic Research (No. 09-02-01241-a), and by the President of RF Grant for Young Researchers MK-1652.2009.2.

## References

- (1) E. Schrödinger, Sitzungsber. Peuss. Akad. Wiss., Phys. Math. Kl. 24, 418 (1930). See also A. O. Barut and A. J. Bracken, where Schrödinger’s work on Zitterbewegung of the electrons is reexamined, Phys. Rev. D 23, 2454 (1981).
- (2) O. Klein, Z.Phys. 53, 157 (1929).
- (3) R. Gerritsma, G. Kirchmair, F. Zähringer, E. Solano, R. Blatt , C. F. Roos, Nature 463, 68 (2010).
- (4) L. Lamata, J. Leon, T. Schätz, and E. Solana, Phys. Rev. Lett. 98, 253005 (2007).
- (5) B. Thaller, arXiv:quant-ph/0409079 (2004, unpublished).
- (6) J. A. Lock, Am. J. Phys. 47, 797 (1979).
- (7) N. Simicevic, arXiv: 0812.1807v1 [physics.comp-ph] (2008); N. Simicevic, arXiv: 0901.3765v1 [quant-ph] (2009).
- (8) P. Krekora, Q. Su, and R. Grobe, Phys. Rev. Lett. 93, 043004 (2004).
- (9) L. Ferrari and G. Russo, Phys. Rev. B 42, 7454 (1990).
- (10) J. Schliemann, D. Loss, and R. M. Westervelt, Phys. Rev. Lett. 94, 206801 (2005).
- (11) V. Ya. Demikhovskii, G. M. Maksimova, and E. V. Frolova, Phys. Rev. B 78, 115401 (2008).
- (12) W. Zawadzki, Phys. Rev. B 72, 085217 (2005).
- (13) M. Katsnelson, Eur. Phys. J. B 51, 157 (2006).
- (14) T. M. Rusin and W. Zawadzki, Phys. Rev. B 76, 195439 (2007).
- (15) G. M. Maksimova, V. Ya. Demikhovskii, and E. V. Frolova, Phys. Rev. B 78, 235321 (2008).
- (16) D. Lurie and S. Cremer, Physica (Amsterdam) 50, 224 (1970).
- (17) J. J. Sakurai, Advanced Quantum Mechanics, Addison-Wesly, Reading, MA, 92, (1967).
- (18) V. Ya. Demikhovskii, G. M. Maksimova, and E. V. Frolova, Phys. Rev. B 81, 115206 (2010).
- (19) U. Schumann, J. Comput. Phys. 18, N. 4, 465 (1975).
- (20) B. Gustafsson , H.-O. Kreiss , J. Oliger, Time Dependent Problems and Difference methods, Wiley, N.Y. (1995).