notebook/quantum-mechanics/1quantumwell.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# **Numerical Solution of the Schrödinger Equation for a One-Dimensional Quantum Well**\n",
    "\n",
    "**Authors:** Dou Du, Taylor James Baird and Giovanni Pizzi\n",
    "\n",
    "<i class=\"fa fa-home fa-2x\"></i><a href=\"../index.ipynb\" style=\"font-size: 20px\"> Go back to index</a>\n",
    "\n",
    "**Source code:** https://github.com/osscar-org/quantum-mechanics/blob/master/notebook/quantum-mechanics/1quantumwell.ipynb\n",
    "\n",
    "This notebook solves numerically the quantum-mechanical problem of a single rectangular one-dimensional quantum well, and displays interactively the eigenfunctions (plotted at the height of the corresponding eigenvalues).\n",
    "\n",
    "<hr style=\"height:1px;border:none;color:#cccccc;background-color:#cccccc;\" />"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## **Goals**\n",
    "* Investigate quantum confinement by manipulating the depth and width of a finite one-dimensional quantum well.\n",
    "* Learn how to solve the Schrödinger equation by matrix diagonalization.\n",
    "* Understand quantum tunnelling in regions where the energy is lower than the potential."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## **Background theory** \n",
    "\n",
    "[More on the background theory.](./theory/theory_1quantumwell.ipynb)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## **Tasks and exercises**\n",
    "\n",
    "1. Keep the depth at the default value (-0.2) and move the slider of the width from the\n",
    "smallest value (0.1) to the largest value (2.0). Do you see any change for the number\n",
    "of the states? How does the eigenvalue of the ground state (the lowest eigenstate) change \n",
    "with the increasing of the width of the quantum well and why?\n",
    "\n",
    "    <details>\n",
    "    <summary style=\"color: red\">Solution</summary>\n",
    "    <div style=\"border:blue; border-width:3px; border-style:outset;\">  \n",
    "        The wider the well potential is, the more (bound state) wavefunctions you will obtain. The reason for this can be traced back to the relationship between the eigenenergies of the bound states and the width of the well: $E_n \\propto \\frac{1}{L^2}$. Therefore, increasing the width of the well decreases the values of the various eigenergies. Note that the \"zero\" of the wavefunction is shifted, for every wavefunction, at the energy of the corresponding eigenvalue.\n",
    "    </div>\n",
    "    </details>\n",
    "    \n",
    "2. Keep the width at the default value (1.2) and move the slider of the depth from\n",
    "the smallest value (-1.0) to the largest value (0.0). Do you see any change for the \n",
    "number of the states? How does the eigenvalue of the ground state (the lowest eigenstate)\n",
    "change with the increasing of the depth of the quantum well and why?\n",
    "\n",
    "    <details>\n",
    "    <summary style=\"color: red\">Solution</summary>\n",
    "    <div style=\"border:blue; border-width:3px; border-style:outset;\">\n",
    "        The deeper the well potential is, the more wavefunctions you will obtain. Mathematically, the explanation for this lies in the relationship between the values of the eigenenergies and the potential depth. This is a transcendental relationship and has to be solved numerically or graphically but one may note that by increasing the well depth, the values of the various energy eigenvalues become progressively more negative. Some important remarks are that, regardless of the well depth, there is always at least one bound solution in one dimension if the well depth is greater than zero, and that as the depth is made progressively larger, the energy levels shall approach those of the infinite square well. For more information see \n",
    "                <a href=\"https://en.wikipedia.org/wiki/Finite_potential_well\">Wikipedia</a>.\n",
    "        Note that the \"zero\" of the wavefunction is shifted, for every wavefunction, \n",
    "        at the energy of the corresponding eigenvalue.\n",
    "    </div>\n",
    "    </details>\n",
    "\n",
    "3. Investigate the role played by quantum confinement in this 1D potential well system by varying the width of the well via the \"width\" slider and observing the spacing of the energy levels of the system.\n",
    "\n",
    "    <details>\n",
    "    <summary style=\"color: red\">Solution</summary>\n",
    "    <div style=\"border:blue; border-width:3px; border-style:outset;\">\n",
    "        When the dimension of the confining potential becomes comparable to the de Broglie wavelength of the confined particle then we begin to observe so-called confinement effects. These manifest as a discretization of the bound state energy levels (when the de Broglie wavelength and well width are not comparable the particle behaves essentially like a free particle, with a continuous energy spectrum). As mentioned in task 1, the relationship between the energy levels and the width of the well is $E_n \\propto \\frac{1}{L^2}$; consequently, as we decrease the width of the well, the spacing of the energy levels becomes larger and larger. As an application of this effect we note that in a nanostructure the effect of increased confinement (reducing the size of the structure) is to increase the band gap of that structure. This can be leveraged to manipulate the opto-electronic properties of the nanostructure. \n",
    "        Furthermore, we note that in a classical potential, the lowest energy state would be at the bottom of the well. In quantum mechanics, the lowest energy state has a higher energy \n",
    "        than the bottom of the quantum well, as you can see in the figure below. \n",
    "        The reason for this energy difference is what is typically called the quantum \n",
    "        confinement effect. Observe also how the lowest energy state changes with \n",
    "        the width of the potential well.\n",
    "    </div>\n",
    "    </details>\n",
    "    \n",
    "4. Keeping the sliders for the width and the height of the well at their default values and selecting the \"Probability density\" option, investigate the probability of a confined particle to observed outside the potential well.\n",
    "  <details>\n",
    "    <summary style=\"color: red\">Solution</summary>\n",
    "    <div style=\"border:blue; border-width:3px; border-style:outset;\">\n",
    "       One can see that even although the potential outside the well is greater than the energy of the particle, there exists a non-zero probability of finding said particle at a position outside the well. This possibility of the particle \"tunnelling\" through the potential barrier in this way is a purely quantum mechanical effect. In the classical picture of this system one would never expect to observe the particle outside the confines of the potential well.\n",
    "    </div>\n",
    "    </details>\n",
    "  \n",
    "5. Please read the background theory section and try to understand the core of the\n",
    "numerical algorithm. Why is the diagonalization of the Hamiltonian matrix key to \n",
    "solve the Schrödinger equation?\n",
    "\n",
    "    <details>\n",
    "    <summary style=\"color: red\">Solution</summary>\n",
    "    <div style=\"border:blue; border-width:3px; border-style:outset;\">\n",
    "        As shown in the background theory section, the Schrödinger equation is an \n",
    "        eigenvalue equation. In mathematics, diagonalization of the matrix allows us\n",
    "        to obtain its eigenvalues and eigenfunctions. You can read more \n",
    "        about the theory on \n",
    "        <a href=\"https://en.wikipedia.org/wiki/Diagonalizable_matrix\">Wikipedia</a>.\n",
    "    </div>\n",
    "    </details>\n",
    "\n",
    "<hr style=\"height:1px;border:none;color:#cccccc;background-color:#cccccc;\" />"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## **Interactive visualization**\n",
    "(be patient, it might take a few seconds to load)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%matplotlib widget\n",
    "\n",
    "from numpy import linspace, sqrt, ones, arange, diag, argsort, zeros\n",
    "from scipy.linalg import eigh_tridiagonal\n",
    "import matplotlib.pyplot as plt\n",
    "from ipywidgets import FloatSlider, jslink, VBox, HBox, Button, Label, RadioButtons"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "colors = ['#a6cee3','#1f78b4','#b2df8a','#33a02c','#fb9a99','#e31a1c','#fdbf6f','#ff7f00','#cab2d6','#6a3d9a','#ffff99','#b15928']\n",
    "ixx = 0\n",
    "\n",
    "def singlewell_potential(x, width, depth):\n",
    "    x1 = zeros(len(x))\n",
    "    for i in range(len(x)):\n",
    "        if x[i] > - width/2.0 and x[i] < width/2.0:\n",
    "            x1[i] = depth\n",
    "    return x1\n",
    "    \n",
    "\n",
    "def diagonalization(hbar, L, N, pot=singlewell_potential, width = 0.1, depth = -0.1):\n",
    "    \"\"\"Calculate sorted eigenvalues and eigenfunctions by diagonalization of the Hamiltonian matrix. \n",
    "\n",
    "       Input:\n",
    "         hbar: reduced Planck constant\n",
    "         L: set discretized interval as [-L,L] \n",
    "         N: number of grid points, i.e., size of the H matrix \n",
    "         pot: python function returning the value of the potential energy\n",
    "         x0: center of the quantum well\n",
    "         width: the width of the quantum well\n",
    "         depth: the depth of the quantum well\n",
    "       Ouput:\n",
    "         ew: sorted eigenvalues (array of length N)\n",
    "         ef: sorted eigenfunctions, ef[:,i] (size N*N)\n",
    "         x:  grid points (array of length N)\n",
    "         dx: grid spacing (called Delta in the theory above)\n",
    "         V:  values of the potential discretized on the grid x (array of length N)\n",
    "    \"\"\"\n",
    "    x = linspace(-L, L, N+2)[1:N+1]                 # grid points\n",
    "    dx = x[1] - x[0]                                # grid spacing\n",
    "    V = pot(x, width, depth)\n",
    "    z = hbar**2 /2.0/dx**2                         # coefficient\n",
    "\n",
    "    ew, ef = eigh_tridiagonal(V+2.0*z, -z*ones(N-1))\n",
    "    ew = ew.real                                    # real part of the eigenvalues\n",
    "    ind = argsort(ew)                               # Indices f. sort. Array\n",
    "    ew = ew[ind]                                    # Sort the ew by ind\n",
    "    ef = ef[:, ind]                                 # Sort the columns \n",
    "    ef = ef/sqrt(dx)                                # Correct standardization \n",
    "    return ew, ef, x, dx, V\n",
    "\n",
    "\n",
    "def plot_eigenfunctions(ax, ew, ef, x, V, width=1, Emax=0.05, fak= 5.0):\n",
    "    \"\"\"Create the full plot.\n",
    "    \n",
    "    Plot the lowest squared eigenfunctions 'ef' at the level of the eigenvalues\n",
    "    'ew' in the plot area 'ax', and the potential 'V(x)'.\n",
    "    \"\"\"\n",
    "    if psi_x.value == \"Wavefunction\":\n",
    "        y_text = r'$V(x)$ and $\\psi(x)$'\n",
    "\n",
    "        fig.suptitle('Numerical Solution ($\\psi$) Of The One Dimensional Schrödinger Equation', fontsize = 13)\n",
    "    else:\n",
    "        y_text = r'$V(x)$ and $|\\psi(x)|^2$'\n",
    "        fig.suptitle('Numerical Solution ($\\psi^2$) Of The One Dimensional Schrödinger Equation', fontsize = 13)\n",
    "        \n",
    "    fak = fak/100.0\n",
    "    \n",
    "    ax[0].axhspan(0.0, Emax, facecolor='lightgrey')\n",
    "    \n",
    "    ax[0].set_xlim([min(x), max(x)])\n",
    "    ax[0].set_ylim([min(V)-0.05, Emax])\n",
    "    \n",
    "    ax[0].set_xlabel(r'$x/a$', fontsize = 10)\n",
    "    ax[0].set_ylabel(y_text, fontsize = 10)\n",
    "    \n",
    "    ax[1].set_xlim([min(x), max(x)])\n",
    "    ax[1].set_ylim([min(V)-0.05, Emax])\n",
    "    \n",
    "    ax[1].yaxis.set_label_position(\"right\")\n",
    "    ax[1].yaxis.tick_right()\n",
    "    \n",
    "    ax[1].get_xaxis().set_visible(False)\n",
    "    ax[1].set_ylabel(r'$\\rm{\\ Eigenvalues}$', fontsize = 10)\n",
    "    \n",
    "    indmax = sum(ew<=0.0)                       \n",
    "    if not hasattr(width, \"__iter__\"):           \n",
    "        width = width*ones(indmax)               \n",
    "    for i in arange(indmax):                     \n",
    "        if psi_x.value == \"Wavefunction\":\n",
    "            ax[0].plot(x, fak*(ef[:, i])+ew[i], linewidth=width[i]+.1, color=colors[i%len(colors)])\n",
    "        else:\n",
    "            ax[0].plot(x, fak*abs(ef[:, i])**2+ew[i], linewidth=width[i]+.1, color=colors[i%len(colors)])\n",
    "\n",
    "        ax[1].plot(x, x*0.0+ew[i], linewidth=width[i]+2.5, color=colors[i%len(colors)])\n",
    "        \n",
    "    ax[0].plot(x, V, c='k', linewidth=1.6)\n",
    "    "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "mu = 0.06                                            # Potential parameter\n",
    "L = 1.5                                              # x range [-L,L]\n",
    "N = 200                                              # Number of grid points\n",
    "hbar = 0.06                                          # Reduced Planck constant\n",
    "sigma_x = 0.1                                        # Width of the Gaussian function\n",
    "zeiten = linspace(0.0, 10.0, 400)                    # time\n",
    "\n",
    "\n",
    "swidth = FloatSlider(value = 1.2, min = 0.1, max = 2.0, description = 'Width: ')\n",
    "sdepth = FloatSlider(value = -0.2, min = -1.0, max = -0.1, step = 0.05, description = 'Depth: ')\n",
    "sfak = FloatSlider(value = 5, min = 1.0, max = 10.0, step = 1.0, description = r'Zoom factor: ')\n",
    "update = Button(description=\"Show all\")\n",
    "psi_x = RadioButtons(options=[\"Wavefunction\", \"Probability density\"], value=\"Wavefunction\", description=\"Show:\")\n",
    "\n",
    "width = 1.2\n",
    "depth = -0.2\n",
    "fak = 5.0\n",
    "\n",
    "ew, ef, x, dx, V = diagonalization(hbar, L, N, width = width, depth = depth)\n",
    "    \n",
    "fig, ax = plt.subplots(1, 2, figsize=(7,5), gridspec_kw={'width_ratios': [10, 1]})\n",
    "fig.canvas.header_visible = False\n",
    "fig.canvas.layout.width = \"750px\"\n",
    "\n",
    "fig.suptitle(r'Numerical Solution ($\\psi^2$) of the Schrödinger Equation for a 1D Quantum Well', fontsize = 12)\n",
    "plot_eigenfunctions(ax, ew, ef, x, V)\n",
    "\n",
    "plt.show()\n",
    "\n",
    "def on_update_click(b):\n",
    "    for i in ax[0].lines:\n",
    "        i.set_alpha(1.0)\n",
    "    for i in ax[1].lines:\n",
    "        i.set_alpha(1.0)\n",
    "    try:\n",
    "        ann.remove()\n",
    "        ann1.remove()\n",
    "    except:\n",
    "        pass\n",
    "\n",
    "def on_width_change(change):\n",
    "    global ew, ef, x, dx, V\n",
    "    ax[0].lines.clear()\n",
    "    ax[1].lines.clear()\n",
    "    \n",
    "    try:\n",
    "        ann.remove()\n",
    "        ann1.remove()\n",
    "    except:\n",
    "        pass\n",
    "\n",
    "    ew, ef, x, dx, V = diagonalization(hbar, L, N, width = swidth.value, depth = sdepth.value)\n",
    "    plot_eigenfunctions(ax, ew, ef, x, V, fak = sfak.value)\n",
    "\n",
    "def on_depth_change(change):\n",
    "    global ew, ef, x, dx, V\n",
    "    ax[0].lines.clear()\n",
    "    ax[1].lines.clear()\n",
    "    \n",
    "    try:\n",
    "        ann.remove()\n",
    "        ann1.remove()\n",
    "    except:\n",
    "        pass\n",
    "\n",
    "    ew, ef, x, dx, V = diagonalization(hbar, L, N, width = swidth.value, depth = sdepth.value)\n",
    "    plot_eigenfunctions(ax, ew, ef, x, V, fak = sfak.value)\n",
    "    \n",
    "def on_xfak_change(change):\n",
    "    ax[0].lines.clear()\n",
    "    ax[1].lines.clear()\n",
    "    \n",
    "    try:\n",
    "        ann.remove()\n",
    "        ann1.remove()\n",
    "    except:\n",
    "        pass\n",
    "\n",
    "    plot_eigenfunctions(ax, ew, ef, x, V, fak = sfak.value)\n",
    "\n",
    "def on_press(event):\n",
    "    global ann, ann1, ixx\n",
    "    \n",
    "    ixx = min(enumerate(ew), key = lambda x: abs(x[1]-event.ydata))[0]\n",
    "    \n",
    "    for i in range(len(ax[1].lines)):\n",
    "        ax[0].lines[i].set_alpha(0.1)\n",
    "        ax[1].lines[i].set_alpha(0.1)\n",
    "        ax[0].lines[i].set_linewidth(1.1)\n",
    "        \n",
    "    ax[0].lines[ixx].set_alpha(1.0)\n",
    "    ax[1].lines[ixx].set_alpha(1.0)\n",
    "    ax[0].lines[ixx].set_linewidth(2.0)\n",
    "    \n",
    "    try:\n",
    "        ann.remove()\n",
    "        ann1.remove()\n",
    "    except:\n",
    "        pass\n",
    "    \n",
    "    ann = ax[0].annotate(s = 'n = ' + str(ixx+1), xy = (0, ew[ixx]), xytext = (-0.15, ew[ixx]), xycoords = 'data', color='k', size=15)\n",
    "    ann1 = ax[1].annotate(s = str(\"{:.3f}\".format(ew[ixx])), xy = (0, ew[ixx]), xytext = (-1.2, ew[ixx]+0.005), xycoords = 'data', color='k', size=9)\n",
    "\n",
    "cid = fig.canvas.mpl_connect('button_press_event', on_press)\n",
    "\n",
    "swidth.observe(on_width_change, names = 'value')\n",
    "sdepth.observe(on_depth_change, names = 'value')\n",
    "sfak.observe(on_xfak_change, names = 'value')\n",
    "update.on_click(on_update_click)\n",
    "psi_x.observe(on_width_change, names = 'value')\n",
    "\n",
    "label1 = Label(value=\"(click on a state to select it)\");\n",
    "\n",
    "display(HBox([swidth, sdepth]), sfak, HBox([psi_x, update, label1]))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* **Width:** the width of the quantum well.\n",
    "* **Depth:** the depth of the quantum well.\n",
    "* **Zoom factor:** the zoom factor for the magnitude of the eigenfunctions."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<hr style=\"height:1px;border:none;color:#cccccc;background-color:#cccccc;\" />\n",
    "\n",
    "## **Legend**\n",
    "\n",
    "(How to use the interactive visualization)\n",
    "\n",
    "### Controls\n",
    "There are three sliders you can tune interactively. Two of them adjust \n",
    "the width and depth of the potential well. The zoom factor slider is used \n",
    "to adjust the normalization factor of the wavefunctions. This factor is used for visualization \n",
    "purposes only, to avoid overlap of the wavefunctions whilst simultaneously preventing the wavefunctions from being so narrow that you cannot appreciate their shape.\n",
    "\n",
    "### Interactive figure\n",
    "In the main visualization part, there are two subplots. The wide figure on \n",
    "the left shows the well potential and the square moduli of the eigenfunctions \n",
    "$|\\psi|^2$ of different states, plotted at the energy of the corresponding \n",
    "eigenvalues. The narrow figure on the right shows the eigenvalues.\n",
    "You can select either an eigenfunction or an eigenvalue by clicking on \n",
    "it in one of the two subplots. The eigenfunction and corresponding eigenvalue \n",
    "will be highlighted, and the corresponding eigenenergy and value of the index \n",
    "$n$ will be shown. All other eigenstates will fade out. By clicking the \"Show all\" \n",
    "button, you can reset the plot and show all states again."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.9.12"
  },
  "voila": {
   "authors": "Dou Du and Giovanni Pizzi"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}