Simulation of Gaussian Wave Packets used to Illustrate Elementary Quantum Mechanics Scenarios

Moving free particles in Quantum Mechanics are some times represented by wave packets. Most of the time, the momentum of a particle is not known with complete certainty, it does not have a definite wavelength. Their wave functions are rather a sum of waves which cover a spectrum of wave numbers. A simple way to represent this is through a Gaussian wave packet, which is nothing more than a sum of plane waves with amplitudes that vary according to a Gaussian distribution. If the wave number of maximum amplitude corresponds to $k_{0}$ and the standard deviation of the distribution is $\sqrt{2}/a$, then the normalized wave function of such particle at $t=0$ is Cohen-Tannoudji et al. (1977):

The time evolution of such packet can be found by summing each plane wave $e^{ikx}$ multiplied by $e^{-iw(k)t}$, its time evolution. After doing the integral, the probability density of the particle happens to be Cohen-Tannoudji et al. (1977):

This spreading property of particles reminds us of the Quantum Safari, defined by George Gamow in his book Gamow (2012),Gamow and Edwards (1999). When Mr. Tompkins and the professor go on an expedition to this place, they soon face a mosquito. But because of its quantum properties, after being spotted, its position becomes less certain over time. After a while, they are covered by a uniform mosquito probability density. This is exactly what happens with the Gaussian packet, $\Delta x$ grows indefinitely until it covers the entire space.

Another thing that can be done with the evolution is illustrate the difference between group and phase velocities. It is still difficult to explain the difference between the group and phase velocities in a few snapshot, thus we have placed an animation of this effect at Gen .

The numerical parameters used for the numerical domain are $x\in[-1,3]$, $N_{x}=2000$, $CFL=0.125$, $N_{t}=40000$, whereas the initial conditions in (7) are $k_{0}=100$ and $a=0.1$.

We now proceed to evolve a moving Gaussian wave packet subject to the effects of a harmonic oscillator potential. The potential and the initial conditions are:

This result is independent of the value of $a$. Nevertheless, when $a$ takes on the special value $a=\sqrt{2/\omega}$, the total energy of the particle is minimum and the wave packet does not spread, a case called quasi-classical state Cohen-Tannoudji et al. (1977). In Fig. 2 we show snapshots of the evolution that illustrate how the pulse oscillates around the minimum. Two snapshots are taken at turning points and the other one at the minimum of the potential. Notice that at turning points the wave function has no nodes, whereas at the center the number of nodes is maximum, which is consistent with the fact that bigger the number of nodes the higher the kinetic energy. Also shown are the expectation values of position $\langle\hat{x}\rangle$ and momentum $\langle\hat{p}\rangle$ as functions of time, on top of them the points indicate the position and momentum of the classical system (11) at some times. Turning points can be identified as the maximae and minimae of $\langle\hat{x}\rangle$ at times $0.08k+0.04$ with $k$ integer.

Additional information in Fig. 2 are the uncertainties $\Delta x$ and $\Delta p$. As expected, since the wave packet does not spread, $\Delta x$ remains nearly constant whereas $\Delta p$ seems to change, however we verified that the variations are due to numerical errors, which is why we also show the uncertainties using a higher numerical resolution with dashed lines; in fact a Richardson extrapolation of the momentum uncertainty using the two resolutions reveals that $\Delta p$ is constant. The product $\Delta x\Delta p$ is also illustrative. Notice that, in agreement with the uncertainty principle, $\Delta x\Delta p$ is always $\geq 1/2$. Since the uncertainty in momentum is affected by numerical errors, then the product is also affected, reason why we also show in dashed line the product using high resolution. The extrapolation corresponds to a constant value in the continuum limit of $1/2$, which coincides with the theoretical value Cohen-Tannoudji et al. (1977).

Finally, in the figure we also show the conservation of the expectation value of the total energy, whereas the potential $\langle T\rangle$ and kinetic $\langle V\rangle$ energies oscillate. Dots indicate the values of these energies in the classical case (11) at various times.

In Fig. 3, we present snapshots of the wave packet evolution with $a=\sqrt{1/\omega}$, different from the quasi-classical profile. It can be seen that the Gaussian spreads as it approaches the turning points. The expectation values of the position and momentum are still the same as the classical analogs, however the uncertainties look very different to the quasi-classical wave packet, for example in this case $\Delta x$ is not constant in time anymore but oscillates. Unlike in the quasi-particle wave packet, the uncertainties genuinely oscillate and an extrapolation of their values using higher resolution does not lead to constant values anymore. The figure shows that the uncertainty principle $\Delta x\Delta p\geq 1/2$ holds, as well as the conservation of total energy. An animation of the evolution of this system can be fount in Gen .

The parameters that define the numerical domain are the same for the two wave packets, namely $x\in[-2,2]$, $N_{x}=2000$, $CFL=0.125$, $N_{t}=640000$, whereas the parameters for the initial conditions are $p_{0}=50$ and $\omega=12.5\pi$.

A more dynamical scenario that involves a time dependent potential is the forced harmonic oscillator. For this, we perturb the ground state of the harmonic oscillator with a sinusoidal wave, which can correspond to a harmonic trap perturbed with an electromagnetic wave, a case used to illustrate resonance Cohen-Tannoudji et al. (1977). The potential and initial condition used for this scenario are:

We solve Schrödinger equation for the evolution of this system with $\alpha=0.2$ and $\omega=12.5\pi$, on the numerical domain defined by $x\in[-2,2]$, $N_{x}=2000$, $CFL=0.125$ and $N_{t}=640000$. The results are shown in Fig. 4. The wave packet oscillates with increasing amplitude, so that turning points are each time further out from the center. The amplitude growth reflects the resonance caused by the sinusoidal perturbation with the appropriate frequency $\omega$. The right graph of Fig. 4 also shows that the Ehrenfest theorem continues to hold for this time dependent potential. Notice that the Gaussian packet preserves its shape during the evolution because the width used corresponds to the quasi-classical packet. For visual grasp we show the animation of this problem at Gen .

Finally, as far as we can tell, there is no exact solution to compare these results with.

Our first example is the quasi-classical analog of the one-dimensional case above, which does not show dispersion of the density of probability. The potential and wave function profile we use are the following:

Quasi-classical case. In Fig. 5 we illustrate the result of the evolution of this wave packet. The wave packet indeed moves on a circle of radius $0.2$ with a period of $0.16$ in the clockwise direction. The shape of the packet does not change because the width of the packet is set to the special value $a=\sqrt{2/\omega}$. This width is that of the ground state of the two-dimensional harmonic potential. For these reasons, we call this packet the quasi-classical solution. The circle in this Figure corresponds to the trajectory of a classical particle with the same initial conditions under the influence of the classical potential described by Eqs. (17) and (18), and also corresponds to the trayectory of the expectation value of the position $\langle\vec{x}\rangle(t)$. This solution of Schrödinger equation was calculated in the numerical domain defined by $x\in[-1,1]$, $y\in[-1,1]$, $Nx=Ny=200$, $CFL=0.125$ and $Nt=12800$, which corresponds to $t_{f}=T$

Another case. Likewise in the 1D scenario we also present the case with a non-educated Gaussian pulse with $a=\sqrt{1/\omega}$. The evolution is very interesting and the density appears in Fig. 6. The pulse evolves in clockwise direction. It starts with a concentrated probability distribution and like in the one dimensional case, it starts spreading as can be seen in the second snapshot, until it crosses the $y-$axis; afterwards, the packet begins to compress again until it acquires its original shape at the $x-$axis again in the third snapshot. The evolution repeats itself with the pulse compressing at the $x-$axis and most expanded at the $y-$axis.

The circle in this Figure corresponds to the expectation value $\langle\vec{x}\rangle(t)$ as well as the trajectory of the classical particle with the same initial conditions under the influence of the classical potential. One takeaway from this exercise is that the Ehrenfest theorem still holds, that is, the expectation value of the position keeps behaving like that of the equivalent classical particle when we jump to two dimensions.

An animation of the quasi-classical state and one with crazy initial conditions can be seen at Gen .

This problem is defined on a two space dimensions domain as follows. A Gaussian wave packet is launched towards a potential wall with a “hole”. The wave function partly bounces from the potential wall and partly passes through. What is illustrative is the well known diffraction pattern that has to be detected beyond the wall.

The initial condition for the wave function is $\psi(x,y,0)=\sqrt{\frac{2}{\pi a^{2}}}e^{{\rm i}p_{0}x}e^{-((x-x_{0})^{2}+y^{2})/a^{2}}$, which has a Gaussian profile with a momentum along the $x-$direction. The potential wall is centered at the $y-axis$ and technically has a finite thickness of five domain cells. The slit is implemented as a segment of the $y-$axis where the potential is zero. This potential can be expressed as follows:

In Fig. 7, one can see the evolution of the wave packet during and after passing through the slit. When the packet arrives at the slit, it starts interfering with itself and then a major part of the wave passes through. The reflected pulse looks interesting, with various fringes although not very well discussed in books. Concerning the pulse that passes through, the pattern forward is the one most analyzed. In order to picture this pattern, we implemented a detector screen placed at $x=0.35$, where we measure the probability density, shown on the fourth image of Fig. 7, measured at the time of the third snapshot. The pattern is composed of three intensity peaks and two minima located at $y=\pm 0.22$. This result can be compared with the Fraunhofer diffraction of a wave going through a slit, which is $\sin(u)/u$ function with minima located at the following positions Hecht (2012):

One of the most used examples of Quantum Mechanics is the experiment of the double slit. A single particle governed by Schrödinger equation is thrown against a barrier with two slits. After the barrier, a screen detector is placed to measure the position of the particle. The experiment is performed multiple times and a detecting plate records the distribution of arriving positions. If the particle were to behave classically, one should only see two spots in the screen, each one associated to the pass of the particle through each of the slits. Nevertheless, what is observed is a set of stripes that reassembles an interference pattern produced by a wave. This is what we will reproduce.

The particle state will be represented by the wave function $\psi$ on a two-dimensional space. The barrier with two slits will be represented by a potential barrier of height $V_{0}$ and $0.05$ units thick. The two slits are separated by the distance $d=0.4$ and $0.1$ units wide. The explicit potential is:

Figure 8 illustrates the behavior during the evolution of the wave packet. After the barrier, the expected interference pattern is formed. The location of the intensity peaks can be compared with the ones predicted by the Young interferometer, originally designed for electromagnetic waves Hecht (2012). According to Young, the intensity peaks are located at the angles given by the formula

The appealing interference pattern seen in the third snapshot of Fig. 8 reminds again of the Quantum Safari Gamow (2012),Gamow and Edwards (1999). In one of the scenes described, a flock of gazelles grazed peacefully when a lioness suddenly appeared. The gazelles ran scared against a row of trees with two small gaps. After crossing the trees, the gazelles divided in columns each one headed directly against a group of hungry lionesses. The gazelles were actually behaving like the quantum particle of Fig. 8, and the lionesses where position in the intensity peaks of the interference pattern.

Animations of the single and double slit problems are found at Gen .

In this example we illustrate what happens when a wave packet collides against a finite height step potential in two dimensions:

In Fig. 9 appears the density of a Gaussian wave packet initially located at $(-0.6,-1.2)$ with momentum $||\vec{p_{0}}||=15\pi$ heading towards the potential step with an incidence angle of $34{}^{o}$. In the second snapshot, the wave packet hits the interface and interferes with itself. In the third snapshot the packet splits into two as the result of the interaction with the step: one bump reflected and another one refracted. Notice that the black dashed line indicates the trajectory of $\langle\vec{x}\rangle(t)$ of the packet prior to the collision. Finally, in the last snapshot, we indicate with red and blue dashed lines the trajectory that $\langle\vec{x}\rangle(t)$ follows for each of the two wave packets. From these expectation values, we find that the reflection and transmission angles are approximately $\theta_{r}=38.93\textordmasculine$ and $\theta_{t}=52.13\textordmasculine$ respectively. A classical particle, as well as a monochromatic wave, would refract with an angle of $80{}^{o}$. We can see that the value obtained numerically doesn’t match the theoretical one very well. This is because the refraction angle depends on the wave number, not only on $\theta_{i}$. Since the wave packet is a sum of plane waves with different momentum, each of the plane waves refracts at a different angle, changing the shape of the packet and the refraction angle of $\langle\vec{x}\rangle$.

The domain parameters used for this simulations are $x\in[-2,2]$, $y\in[-2,2]$, $Nx=Ny=400$, $CFL=0.125$ and $Nt=3900$, and the initial conditions are set to $\vec{x}_{0}=-0.6\hat{\textbf{\i }}-1.2\hat{\textbf{\j }}$, $\vec{p}_{0}=15\pi(\cos(34{}^{o})\hat{\textbf{\i }}+\sin(34^{o})\hat{\textbf{\j }})$, packet width $a=0.15$ and $V_{0}=750$.

siddh

As a final example, we look at a case of dispersion by a central potential. In this experiment, a particle (or wave packet) is directed against a potential that is spherically symmetric. One not so obvious assumption that is made by Cohen-Tannoudji et al. (1977) for the treatment of this problem, is that after the plane wave interacts with the potential, the wave function becomes the sum of a plane wave that goes through and a spherical wave front with nonuniform amplitude. The objective of this section is to show numerically that the outgoing wave function truly has a spherical component.

We threat the two dimensional case, with potential and initial conditions as follows:

The domain parameters used for this simulations were: $x\in[-1.5,1.5]$, $y\in[-1.5,1.5]$, $Nx=Ny=300$, $CFL=0.125$ and $Nt=3900$, and initial conditions $x_{0}=-0.7$, $p_{0}=5\pi$, $a=0.15$, $V_{0}=120$ and $\sigma=0.3$.