Worked example · solved on paper

Quantum Tunnelling of a Wave Packet through a Rectangular Barrier (Analytical)

PHYS 3260, Mathematical Methods — worked example

The problem

A particle of mass \(m\) moves in one dimension and meets a rectangular barrier of height \(V_0\) and width \(a\),

\[ V(x) = V_0 \ \ \text{for } |x| < a/2, \qquad V(x) = 0 \ \ \text{otherwise,} \]

with \(\hbar = m = 1\), \(V_0 = 2\) and \(a = 1\). Its wavefunction obeys the time-dependent Schrödinger equation

\[ i\hbar\,\frac{\partial \psi}{\partial t} = -\frac{\hbar^{2}}{2m}\,\frac{\partial^{2}\psi}{\partial x^{2}} + V(x)\,\psi . \]

Abstract. The time-dependent Schrödinger equation is solved in closed form for a particle incident on a rectangular barrier. Separation of variables reduces it to a linear ordinary differential equation with constant coefficients in each of the three regions where the potential is constant; requiring \(\varphi\) and \(\varphi'\) to be continuous at the two walls gives four equations, which are assembled into a \(2\times2\) transfer matrix and solved for the reflection and transmission amplitudes. The result is a closed-form transmission coefficient \(T(E)\) valid at every energy, from which \(R+T=1\) follows as an algebraic identity and the above-barrier resonances \(T=1\) at \(qa=n\pi\) follow by analytic continuation. A Gaussian wave packet with mean energy \(0.64\,V_0\) is then expanded in these exact stationary states and evolved by superposition, giving \(T = 0.290\) and \(R = 0.710\) with \(R+T = 1.000\). The spectral decomposition shows the barrier acting as a high-pass filter: the transmitted packet's mean wavenumber is \(1.0\%\) higher than the incident packet's. The methods used are complex exponentials (Ch. 2), matrices (Ch. 3), the Fourier transform (Ch. 7), linear ordinary differential equations (Ch. 8) and separation of variables (Ch. 10).

1Separating the variables (Chapter 10)

The potential does not depend on time, so the equation admits solutions of the product form \(\psi(x,t) = \varphi(x)\,e^{-iEt/\hbar}\). Substituting this into the Schrödinger equation gives

\begin{equation} i\hbar\,\varphi(x)\,\frac{d}{dt}e^{-iEt/\hbar} = \left[-\frac{\hbar^{2}}{2m}\varphi''(x) + V(x)\varphi(x)\right]e^{-iEt/\hbar}, \end{equation}

and since \(i\hbar\,\frac{d}{dt}e^{-iEt/\hbar} = E\,e^{-iEt/\hbar}\), the exponential cancels from both sides. What is left is an ordinary differential equation for \(\varphi\) alone,

\begin{equation} \frac{d^{2}\varphi}{dx^{2}} + \frac{2m}{\hbar^{2}}\bigl(E - V(x)\bigr)\varphi = 0 . \label{eq:ode} \end{equation}

Equation \(\eqref{eq:ode}\) is the whole of the physics in this report. Everything that follows is Chapter 8 applied to it: a second-order linear equation with constant coefficients, solved region by region and joined at the boundaries.

2Solving in each region (Chapter 8)

The potential is piecewise constant, so \(\eqref{eq:ode}\) has constant coefficients in each of three regions, and its solutions there are exponentials.

Outside the barrier (\(|x| > a/2\)).

Here \(V = 0\), so the coefficient of \(\varphi\) is \(+2mE/\hbar^{2}\), which is positive. The characteristic equation \(\lambda^{2} + k^{2} = 0\) has imaginary roots, and the two solutions oscillate:

\begin{equation} \varphi = A\,e^{ikx} + B\,e^{-ikx}, \qquad k = \frac{\sqrt{2mE}}{\hbar} . \label{eq:outside} \end{equation}

Inside the barrier (\(|x| < a/2\)).

Here \(V = V_0 > E\), so the same coefficient is \(-2m(V_0-E)/\hbar^{2}\), which is negative. The characteristic equation \(\lambda^{2} - \kappa^{2} = 0\) has real roots, and the two solutions do not oscillate at all — one grows and one decays:

\begin{equation} \varphi = C\,e^{\kappa x} + D\,e^{-\kappa x}, \qquad \kappa = \frac{\sqrt{2m(V_0 - E)}}{\hbar} . \label{eq:inside} \end{equation}

That change of sign is the origin of tunnelling. A classical particle stops at the wall because no real motion is available to it at \(E < V_0\); the wave does not stop, it changes from a sinusoid into an exponential and continues, reduced. Both solutions in \(\eqref{eq:inside}\) are needed. The decaying one dominates near the first wall, but the growing one, small as its coefficient is, is what carries an amplitude to the second wall and out the far side.

Physical conditions now select which combination is wanted. The particle is incident from the left, so to the left of the barrier there is an incoming wave and a reflected one, and to the right there is only an outgoing wave — nothing returns from \(x = +\infty\). Normalising the incoming amplitude to \(1\),

\begin{equation} \varphi_k(x) = \begin{cases} e^{ikx} + r\,e^{-ikx}, & x < -a/2, \\[2pt] C\,e^{\kappa x} + D\,e^{-\kappa x}, & |x| < a/2, \\[2pt] t\,e^{ikx}, & x > +a/2 , \end{cases} \label{eq:ansatz} \end{equation}

which leaves four unknown constants: \(r\), \(t\), \(C\) and \(D\).

The three regions, the form the solution takes in each, and the two points at which the pieces must be joined. Four matching conditions determine the four unknown constants in eq:ansatz.
Figure 1. The three regions, the form the solution takes in each, and the two points at which the pieces must be joined. Four matching conditions determine the four unknown constants in \(\eqref{eq:ansatz}\).

3The matching conditions

At a step of finite height, both \(\varphi\) and \(\varphi'\) must be continuous. The reason is \(\eqref{eq:ode}\) itself: it expresses \(\varphi''\) in terms of \(\varphi\), and if \(\varphi\) and \(V\) are bounded then \(\varphi''\) is bounded, so \(\varphi'\) cannot jump and \(\varphi\) cannot have a kink. Continuity of both quantities at each of the two walls gives four equations:

\begin{align} e^{-ika/2} + r\,e^{ika/2} &= C\,e^{-\kappa a/2} + D\,e^{\kappa a/2}, \label{eq:m1}\\ ik\left(e^{-ika/2} - r\,e^{ika/2}\right) &= \kappa\left(C\,e^{-\kappa a/2} - D\,e^{\kappa a/2}\right), \label{eq:m2}\\ C\,e^{\kappa a/2} + D\,e^{-\kappa a/2} &= t\,e^{ika/2}, \label{eq:m3}\\ \kappa\left(C\,e^{\kappa a/2} - D\,e^{-\kappa a/2}\right) &= ik\,t\,e^{ika/2}. \label{eq:m4} \end{align}

Four linear equations in four unknowns: the problem is now solved in principle, and the rest is bookkeeping.

4The transfer matrix (Chapter 3)

Solving \(\eqref{eq:m1}\)\(\eqref{eq:m4}\) by substitution is possible but unpleasant. Written in matrix form the same elimination is two lines. For a solution built out of \(e^{\pm px}\), collect the values of \(\varphi\) and \(\varphi'\) at a point into a single matrix acting on the pair of amplitudes,

\begin{equation} M(p, x) = \begin{pmatrix} e^{px} & e^{-px} \\ p\,e^{px} & -p\,e^{-px} \end{pmatrix}, \qquad\text{so that}\qquad \begin{pmatrix} \varphi(x) \\ \varphi'(x) \end{pmatrix} = M(p, x) \begin{pmatrix} \text{amplitude of } e^{px} \\ \text{amplitude of } e^{-px} \end{pmatrix} . \label{eq:Mmatrix} \end{equation}

Continuity at a wall is then the statement that the two sides produce the same column vector there. Equations \(\eqref{eq:m1}\)\(\eqref{eq:m2}\) become

\begin{equation} M(ik, -a/2)\binom{1}{r} = M(\kappa, -a/2)\binom{C}{D}, \label{eq:matchleft} \end{equation}

and \(\eqref{eq:m3}\)\(\eqref{eq:m4}\) become

\begin{equation} M(\kappa, +a/2)\binom{C}{D} = M(ik, +a/2)\binom{t}{0}. \label{eq:matchright} \end{equation}

Solving \(\eqref{eq:matchright}\) for \((C,D)\) and substituting into \(\eqref{eq:matchleft}\) eliminates the interior amplitudes and leaves the whole barrier as a single matrix,

\begin{equation} \binom{1}{r} = \underbrace{M(ik,-\tfrac a2)^{-1} M(\kappa,-\tfrac a2)\, M(\kappa,+\tfrac a2)^{-1} M(ik,+\tfrac a2)}_{\textstyle M(E)} \binom{t}{0} . \label{eq:transfer} \end{equation}

The two rows of \(\eqref{eq:transfer}\) read \(1 = M_{11}t\) and \(r = M_{21}t\), so

\begin{equation} t = \frac{1}{M_{11}}, \qquad r = \frac{M_{21}}{M_{11}} . \label{eq:rt} \end{equation}

Carrying out the matrix product and collecting the exponentials into hyperbolic functions gives the transmission amplitude in closed form,

\begin{equation} t(k) = \frac{e^{-ika}} {\cosh(\kappa a) + \dfrac{i}{2}\left(\dfrac{\kappa}{k} - \dfrac{k}{\kappa}\right)\sinh(\kappa a)} . \label{eq:tamp} \end{equation}

The structure of \(\eqref{eq:transfer}\) is worth noting. A second barrier would insert two more factors into the same product and change nothing else in the derivation, which is what makes the transfer matrix the standard tool for layered structures.

5The transmission coefficient

Taking the squared modulus of \(\eqref{eq:tamp}\), and simplifying with \(k^{2} + \kappa^{2} = 2mV_0/\hbar^{2}\) and \(k^{2}\kappa^{2} = 4m^{2}E(V_0-E)/\hbar^{4}\), gives

\begin{equation} \boxed{\; T(E) = |t|^{2} = \left[\,1 + \frac{V_0^{2}\,\sinh^{2}(\kappa a)}{4E(V_0 - E)}\,\right]^{-1} , \qquad \kappa = \frac{\sqrt{2m(V_0-E)}}{\hbar} . \;} \label{eq:T} \end{equation}

\(T\) is non-zero for every \(E > 0\), however far below \(V_0\): that is tunnelling, and here it is a property of a formula rather than an observation about a particular run.

5.1\(R+T=1\) is an identity

The same calculation applied to \(r\) gives

\begin{equation} R = |r|^{2} = \left[\,1 + \frac{4E(V_0 - E)}{V_0^{2}\sinh^{2}(\kappa a)}\,\right]^{-1} . \label{eq:R} \end{equation}

Writing \(u = V_0^{2}\sinh^{2}(\kappa a)/4E(V_0-E)\), equations \(\eqref{eq:T}\) and \(\eqref{eq:R}\) are \(T = 1/(1+u)\) and \(R = 1/(1+u^{-1}) = u/(1+u)\), so

\begin{equation} R + T = \frac{u}{1+u} + \frac{1}{1+u} = 1 \end{equation}

for every \(u\) — that is, at every energy, for every barrier height and every width. Conservation of probability is therefore not a check to be performed on this solution but a property of it. Evaluating \(|r|^{2}+|t|^{2}-1\) numerically across the range of energies used below returns \(9\times10^{-16}\), which measures the arithmetic of the computer and not the correctness of the physics.

5.2Above the barrier

Equation \(\eqref{eq:T}\) was derived assuming \(E < V_0\), but nothing in the final expression requires it. For \(E > V_0\), \(\kappa\) becomes imaginary: writing \(\kappa = iq\) with \(q = \sqrt{2m(E-V_0)}/\hbar\) and using \(\sinh(iqa) = i\sin(qa)\) turns \(\sinh^{2}(\kappa a)\) into \(-\sin^{2}(qa)\), so

\begin{equation} T(E) = \left[\,1 - \frac{V_0^{2}\sin^{2}(qa)}{4E(V_0-E)}\,\right]^{-1} = \left[\,1 + \frac{V_0^{2}\sin^{2}(qa)}{4E(E-V_0)}\,\right]^{-1}, \qquad E > V_0 . \label{eq:Tabove} \end{equation}

This is a substitution, not a second derivation, and it is available only because the solutions were written as exponentials in the first place (Chapter 2). It predicts \(T = 1\) exactly whenever \(qa = n\pi\) — a barrier whose width is a whole number of half-wavelengths is perfectly transparent, the same interference condition that makes an anti-reflection coating work.

Equations eq:T and eq:Tabove, plotted with R = 1-T. Transmission is non-zero throughout the classically forbidden region E&lt;V0, and reaches exactly 1 at the resonances above the barrier. The marked point is the mean energy of the wave packet used in Step 6.
Figure 2. Equations \(\eqref{eq:T}\) and \(\eqref{eq:Tabove}\), plotted with \(R = 1-T\). Transmission is non-zero throughout the classically forbidden region \(E<V_0\), and reaches exactly \(1\) at the resonances above the barrier. The marked point is the mean energy of the wave packet used in Step 6.

At the energy used below, \(E = 1.28 = 0.64\,V_0\), equation \(\eqref{eq:T}\) gives \(T = 0.288\).

6The wave packet (Chapter 7)

A stationary state has a single energy and is spread over the whole line; a particle that is somewhere and moving is a superposition of them. Building that superposition is a Fourier transform, and it introduces no approximation.

The packet at \(t = 0\) is the Gaussian

\begin{equation} \psi(x,0) = (2\pi\sigma^{2})^{-1/4} \exp\!\left[-\frac{(x-x_0)^{2}}{4\sigma^{2}}\right] e^{\,ik_0(x-x_0)} . \label{eq:packet} \end{equation}

Its Fourier transform is obtained by completing the square in the exponent, one of the few transforms that can be done in a single line:

\begin{equation} A(k) = \frac{1}{\sqrt{2\pi}}\int_{-\infty}^{\infty} \psi(x,0)\,e^{-ik(x-x_0)}\,dx = \left(\frac{2\sigma^{2}}{\pi}\right)^{1/4} e^{-\sigma^{2}(k-k_0)^{2}} , \label{eq:Ak} \end{equation}

a Gaussian in \(k\) centred on \(k_0\) with width \(\Delta k = 1/2\sigma = 0.1\), narrow compared with \(k_0 = 1.6\).

Equation \(\eqref{eq:Ak}\) expresses the initial packet as a sum of plane waves. Each plane wave \(e^{ikx}\) in that sum is replaced by the exact scattering state \(\varphi_k(x)\) of \(\eqref{eq:ansatz}\), which is what that plane wave becomes in the presence of the barrier, and each evolves with its own phase \(e^{-iE_kt/\hbar}\). The exact solution of the original time-dependent problem is therefore

\begin{equation} \psi(x,t) = \frac{1}{\sqrt{2\pi}}\int_{-\infty}^{\infty} dk\; A(k)\,e^{-ikx_0}\,\varphi_k(x)\,e^{-i\hbar k^{2}t/2m} . \label{eq:evolution} \end{equation}

Every ingredient of \(\eqref{eq:evolution}\) is a closed-form expression. The integral over \(k\) is not: \(\varphi_k\) depends on \(k\) through \(\cosh\) and \(\sinh\) of \(\kappa(k)a\), and no elementary antiderivative exists. It is evaluated numerically, and it is the only step of this report that is. The distinction matters — no grid carries the differential equation, no time step is taken, and \(\psi\) at \(t = 55\) is obtained without first obtaining \(\psi\) at \(t = 54\). Each frame of the animation below is \(\eqref{eq:evolution}\) evaluated at one time, independently of every other.

6.1Predicting the result in advance

Because the plane waves evolve independently, each is transmitted with probability \(T(E_k)\), and the transmitted fraction of the whole packet is the average of \(T\) over the packet's own spectrum:

\begin{equation} T_{\text{packet}} = \int_{0}^{\infty} |A(k)|^{2}\, T\!\left(\frac{\hbar^{2}k^{2}}{2m}\right) dk = 0.2901 . \label{eq:Tpacket} \end{equation}

This is a prediction made before any evolution is carried out.

7Results

Evaluating \(\eqref{eq:evolution}\) on a grid of times and integrating \(|\psi|^{2}\) on each side of the barrier at \(t = 55\), when the reflected and transmitted parts have fully separated, gives

\begin{equation} R = 0.7099, \qquad T = 0.2901, \qquad R + T = 1.0000 . \label{eq:results} \end{equation}

The transmitted fraction agrees with the prediction \(\eqref{eq:Tpacket}\) to four decimal places, as it must: the two are the same integral evaluated in a different order.

Four instants of the collision, each computed directly from eq:evolution. The packet approaches; the ripples in the second panel are the incident and reflected waves superposed on one another in front of the barrier; the two parts then separate, leaving a large reflected packet moving left and a smaller transmitted packet moving right. All four panels share one vertical scale.
Figure 3. Four instants of the collision, each computed directly from \(\eqref{eq:evolution}\). The packet approaches; the ripples in the second panel are the incident and reflected waves superposed on one another in front of the barrier; the two parts then separate, leaving a large reflected packet moving left and a smaller transmitted packet moving right. All four panels share one vertical scale.
The collision in full. Above, the probability density |(x,t)|^2 with the barrier marked; below, the probability on each side of the barrier as it accumulates, so that R and T are watched forming rather than quoted at the end. Nothing crosses until the packet arrives near t=20, the transfer is complete by t35, and both curves are flat afterwards because the two parts have separated for good.
Figure 4. The collision in full. Above, the probability density \(|\psi(x,t)|^{2}\) with the barrier marked; below, the probability on each side of the barrier as it accumulates, so that \(R\) and \(T\) are watched forming rather than quoted at the end. Nothing crosses until the packet arrives near \(t=20\), the transfer is complete by \(t\approx35\), and both curves are flat afterwards because the two parts have separated for good.

7.1Why the packet does not transmit like a plane wave

Equation \(\eqref{eq:Tpacket}\) gives \(0.2901\), while \(\eqref{eq:T}\) evaluated at the packet's mean energy alone gives \(0.288\). The difference is not an error; it is the answer to part (d) of the problem.

The packet contains a spread of energies, and \(T(E)\) is not constant over that spread — it rises by a factor of \(11\) between the slowest and the fastest components present. Transmission therefore acts as a filter on the packet's spectrum,

\begin{equation} |A_{\text{transmitted}}(k)|^{2} \;\propto\; |A(k)|^{2}\,T(E_k), \label{eq:filter} \end{equation}

weighting it towards larger \(k\). Two consequences follow. First, \(\langle T(E)\rangle > T(\langle E\rangle)\), because \(T\) is convex across the packet: the fast half of the packet gains more than the slow half loses. Second, the transmitted packet is measurably faster than the incident one, its mean wavenumber rising from \(1.600\) to \(1.616\). Nothing accelerated; the barrier simply passed the fast components preferentially and reflected the slow ones.

Above: the transmission of each plane wave present in the packet, which rises by a factor of 11 across the packet's own width. Below: the packet's spectrum before and after the barrier, the second being the first multiplied by that curve and renormalised, with the shift in k marked.
Figure 5. Above: the transmission of each plane wave present in the packet, which rises by a factor of \(11\) across the packet's own width. Below: the packet's spectrum before and after the barrier, the second being the first multiplied by that curve and renormalised, with the shift in \(\langle k\rangle\) marked.

8Which course methods were used

9Conclusion

Matching exponential solutions across the two walls of a rectangular barrier gives the transmission coefficient in closed form, equation \(\eqref{eq:T}\), valid at every energy. From that single expression follow three results that require no further calculation: transmission is non-zero at every energy below the barrier, \(R+T=1\) holds identically rather than approximately, and the barrier becomes perfectly transparent above its top whenever its width is a whole number of half-wavelengths.

Expanding a Gaussian wave packet in these exact stationary states and evolving it by superposition reproduces the full collision, with \(T = 0.290\) and \(R = 0.710\) summing to unity, and the transmitted fraction is predicted exactly by averaging \(T(E)\) over the packet's spectrum before any evolution is performed. That spectral view also explains why a packet does not transmit like a plane wave of its mean energy: the barrier is a high-pass filter, and what emerges from it is faster than what arrived.

The method generalises directly. Any potential that is piecewise constant can be solved the same way, with one factor in the transfer-matrix product per layer, and a smoothly varying potential can be approximated to any accuracy by enough layers.

References

← All projects