Quantum Mechanics

In the last lecture, we have modeled electromagnetic waves not by solving the wave equation, but by taking the solutions of wave equations like a plane wave or a spherical wave. Today we will solve a wave equation, but not for electromagnetic waves, but for matter waves. We will solve the stationary Schrödinger equation with the implicit solution scheme, which we have already applied for the diffusion equation. With the help of that we will tackle the particle in a box, the harmonic oscillator and the periodic potential. All of these problems have also analytical solutions, thus we do not need the numerical solution in principle. But it will give us some practice on how to tackle such problems. As not all of you might be familiar with the physical description of quantum mechanics, we will give a short introduction into this field first.

[4]:
import numpy as np
from scipy.sparse import diags
from scipy.sparse.linalg import eigsh
import matplotlib.pyplot as plt

%matplotlib inline
%config InlineBackend.figure_format = 'retina'

# default values for plotting
plt.rcParams.update({'font.size': 12,
                     'axes.titlesize': 18,
                     'axes.labelsize': 16,
                     'axes.labelpad': 14,
                     'lines.linewidth': 1,
                     'lines.markersize': 10,
                     'xtick.labelsize' : 16,
                     'ytick.labelsize' : 16,
                     'xtick.top' : True,
                     'xtick.direction' : 'in',
                     'ytick.right' : True,
                     'ytick.direction' : 'in',})

Quantum Mechanics in a Nutshell

Quantum Mechanics assumes that all particles propagate as waves. They are described by a wavefunction \(\Psi(x,t)\). A quantum mechanical object thus posseses an amplitude and a phase which propagate in space and time. One could see the wavefunction in analogy to the electric field \(\vec{E}(x)\) of an electromagnetic wave. As the square of the electric field describes the propagation of energy of a wave, the square magnitude of the wavefunction, i.e. \(|\Psi(x,t)|^2\), describes the propagation of probability density of the quantum mechanical wave. The wavefunction itself is thus just the probability amplitude.

Time dependent Schrödinger equation

The dynamics of a quantum mechanical wave is described,for example, by the time dependent Schrödinger equation

\begin{equation} -i\hbar\frac{\partial \Psi(x,t)}{\partial t} = \left ( \frac{-\hbar^2 }{2m}\frac{\partial^2}{\partial x^2}+V(x,t) \right ) \Psi(x,t) \end{equation}

whis is written here for one dimension only.

The bracket on the right side of the above equation contains the so-called Hamilton operator \(\hat{H}\). The Hamilton operator \(\hat{H}\) contains the energy operators for the kinetic and potential energies and represents the total energy of the system.

\begin{equation} \hat{H}=\left ( \frac{-\hbar^2 }{2m}\frac{\partial^2}{\partial x^2}+V(x,t) \right ) \end{equation}

Stationary Schrödinger equation

Our first problems will be stationary problems. We will not ask for the temporal development of the quantum object. We will rather ask, what solutions without time dependence are possible. In general this is much like the question asking what kind of standing waves are possible on a string or in an optical resonator. In quantum mechanics the boundaries, which define the standing waves are formed by the potential energy \(V(x)\).

We therefore also need the stationary Schrödinger equation, where the left side of the time dependent Schrödinger equation does not depend on time, hence is constant in time. This stationary (time-independent) Schrödinger equation is

\begin{equation} \hat{H}\Psi(x)=E\Psi(x) \end{equation}

The Hamilton operator \(\hat{H}\) gives a recipe how to calculate the energies for a given wavefunction \(\Psi(x)\) in terms of derivates or multiplications by functions. If this recipe reduces to a multiplication of the wave function with a number \(E\), then these wavefunctions are eigenfunction of the Hamilton operator and the values of \(E\) are the eigenvalues of the problem, i.e. the time-independent solutions of this differential equation.

Recap: Implicit Solution

According to our above description, the Hamilton operator \(\hat{H}\) contains two parts, a second derivative in the position, which represents the kinetic energy and the potential energy operator \(V(x)\), which is in the simplest case just a function of \(x\).

\begin{equation} \left ( \frac{-\hbar^2 }{2m}\frac{\partial^2}{\partial x^2}+V(x) \right ) \Psi(x)=E\Psi(x) \end{equation}

Since we want to apply our implicit solution scheme (Cranck Nicolson), we want to represent both parts as matrices.

Kinetic energy

We remember that we can write the second derivative of our wavefunction \(\Psi(x)\) in the finite difference approximation as

\begin{equation} \Psi^{''}(x)=\frac{\Psi(x+\delta x)-2\Psi(x)+\Psi(x-\delta x)}{\delta x^{2}} \end{equation}

If we want to evaluate the wavefunction at certain positions \(x_{i}\), then this second derivative translates into

\(T\Psi=\frac{d^2}{dx^2}\Psi=\frac{1}{\delta x^2} \begin{bmatrix} -2 & 1 & 0 & 0 & 0 & 0\\ 1 & -2 & 1 & 0 & 0 & 0\\ 0 & 1 & -2 & 1 & 0 & 0\\ 0 & 0 & 1 & -2 & 1 & 0\\ 0 & 0 & 0 & 1 & -2 & 1\\ 0 & 0 & 0 & 0 & 1 & -2\\ \end{bmatrix} \begin{bmatrix} \Psi(x_{1})\\ \Psi(x_{2})\\ \Psi(x_{3})\\ \Psi(x_{4})\\ \Psi(x_{5})\\ \Psi(x_{6}) \end{bmatrix}\)

if we just use 6 positions. Please remember, that in the version above, we have imposed already boundary conditions in the first and the last row, which are \(\Psi(x_{0})=0\) and \(\Psi(x_{7})=0\). If we multiply this matrix by \(-\hbar^{2}/2m\), we obtain the kinetic energy for an object of mass \(m\).

Potential energy

The potential energy values are just values at the diagonal of the matrix

\(V\Psi= \begin{bmatrix} V(x_{1}) & 0 & 0 & 0 & 0 & 0\\ 0 & V(x_{2}) & 0 & 0 & 0 & 0\\ 0 & 0 & V(x_{3}) & 0 & 0 & 0\\ 0 & 0 & 0 & V(x_{4}) & 0 & 0\\ 0 & 0 & 0 & 0 & V(x_{5}) & 0\\ 0 & 0 & 0 & 0 & 0 & V(x_{6})\\ \end{bmatrix} \begin{bmatrix} \Psi(x_{1})\\ \Psi(x_{2})\\ \Psi(x_{3})\\ \Psi(x_{4})\\ \Psi(x_{5})\\ \Psi(x_{6}) \end{bmatrix}\)

an you may insert the specific potential energy values for your particular problem here.

Our final problem \(\hat{H}\Psi=E\Psi\) will thus have the following shape

\begin{equation} \begin{bmatrix} -2+V(x_{1}) & 1 & 0 & 0 & 0 & 0\\ 1 & -2+V(x_{2}) & 1 & 0 & 0 & 0\\ 0 & 1 & -2+V(x_{3}) & 1 & 0 & 0 \\ 0 &0 & 1 & -2+V(x_{4}) & 1 & 0 \\ 0 & 0 & 0 & 1 & -2+V(x_{5}) & 1 \\ 0 & 0 & 0 & 0 & 1 & -2+V(x_{6}) \\ \end{bmatrix} \begin{bmatrix} \Psi(x_{1})\\ \Psi(x_{2})\\ \Psi(x_{3})\\ \Psi(x_{4})\\ \Psi(x_{5})\\ \Psi(x_{6}) \end{bmatrix}=E \begin{bmatrix} \Psi(x_{1})\\ \Psi(x_{2})\\ \Psi(x_{3})\\ \Psi(x_{4})\\ \Psi(x_{5})\\ \Psi(x_{6}) \end{bmatrix} \end{equation}

where I skipped the prefactor \(-\hbar^2/2m\), to fit the matrices on one line. Yet I did not succeed. This is the final system of coupled equations which we can supply to any matrix solver. We will use a solver from the scipy.linalg module. In case we have special boundary conditions, we need to take them into account and replace the first and the last line for example with the particular boundary conditions.

[ ]: