Tutorials/Exercise/TI_01/Exercise_TI-01.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "a28357c6",
   "metadata": {},
   "outputs": [],
   "source": [
    "import sisl\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "%matplotlib inline\n",
    "import functools"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "d4f313be",
   "metadata": {},
   "source": [
    "# Quantum Hall Effect\n",
    "In this exercise, we will build on **TB_07** and **A_03** to simulate the quantum Hall effect. In **A_03**, we learned how to created a Hall bar device and in **TB_07**, how to deal with magnetic fields.\n",
    "\n",
    "Here, we will extract the Hall resistance from the transmissions calculated with `TBtrans` using the Landauer-Büttiker formalism.\n",
    "\n",
    "## Exercise Overview:\n",
    "1. Create a Hall bar (see **A_03**)\n",
    "2. Construct Hamiltonians and add magnetic fields (see **TB_07**)\n",
    "3. Calculate the transmission with `TBtrans`\n",
    "4. Extract the Hall resistance $(R_H)$.\n",
    "4. Extract the logitudinal resistance $(R_L)$.\n",
    "\n",
    "## Exercise"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "287ce57a",
   "metadata": {},
   "source": [
    "### 1. Create a Hall bar\n",
    "In order to be able to observe the quantum Hall effect, the size of the Hall bar needs to be big enough. For a 4(6) lead device reasonable dimensions are: \n",
    "1. 4-lead device (square lattice):\n",
    "   - Width of electrodes (perpendicular to the semi-infinite axis): 30 atoms\n",
    "   - Offset of the electrodes 2(3) from the corner of the central part: >= 5 atoms\n",
    "\n",
    "\n",
    "2. 6-lead device (square lattice): \n",
    "   - Width of electrodes (perpendicular to the semi-infinite axis): 30 atoms\n",
    "   - Spacing between electrodes on the same side: 50 atoms\n",
    "   - Offset of the electrodes 2(3,5,6) from the corner of the central part: > 25 atoms\n",
    " \n",
    "\n",
    "3. 6-lead Graphene Hall bar:\n",
    "   - Width of electrodes (perpendicular to the semi-infinite axis): ~30 Å\n",
    "   - Spacing between electrodes on the same side: ~50 Å\n",
    "   - Offset of the electrodes 2(3,5,6) from the corner of the central part: ~ 15 Å\n",
    "   \n",
    "|4-lead device |6-lead device |\n",
    "|:-------------|:-------------|\n",
    "|<img src=\"img/Bi2D-BHex-c.png\" alt=\"BHex-c\" style=\"width: 400px;\"/>|<img src=\"img/Bi2D-BHex-a.png\" alt=\"BHex-a\" style=\"width: 400px;\"/>|"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "15cf46df",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Create a Hall bar"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ddd90a12",
   "metadata": {},
   "source": [
    "### 2. Construct Hamiltonian and add magnetic fields\n",
    "\n",
    "The required field strengths may vary depending on the size of the Hall bar. We should start with a corse grid, and create a finer grid once we have identified the correct range. A good starting point might by `B = 1 / np.arange(10,31)`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "889e23fd",
   "metadata": {},
   "outputs": [],
   "source": [
    "def peierls(self, ia, atoms, atoms_xyz=None, B=None):\n",
    "    idx = self.geometry.close(ia, R=[0.1, 1.01], atoms=atoms, atoms_xyz=atoms_xyz)\n",
    "    # Onsite\n",
    "    self[ia, idx[0]] = 4\n",
    "    \n",
    "    # Hopping\n",
    "    if B == 0:\n",
    "        self[ia, idx[1]] = -1\n",
    "    else:\n",
    "        xyz = self.geometry.xyz[ia]\n",
    "        dxyz = self.geometry[idx[1]]\n",
    "        self[ia, idx[1]] = - np.exp(-0.5j * B * (dxyz[:, 0] - xyz[0])*(dxyz[:,1] + xyz[1]))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "4c615c6d",
   "metadata": {},
   "outputs": [],
   "source": [
    "H0 = sisl.Hamiltonian(geom, dtype=np.float64)\n",
    "H0.construct(functools.partial(peierls, B=0.))\n",
    "\n",
    "HB = sisl.Hamiltonian(geom, dtype=np.complex128)\n",
    "\n",
    "rec_phis = np.arange(10,31)\n",
    "for rec_phi in rec_phis:\n",
    "    HB.construct(functools.partial(peierls, B=1/rec_phi))\n",
    "    dH = ..."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "2a5971f1",
   "metadata": {},
   "source": [
    "### 3. Calculate the transmission with TBtrans\n",
    "\n",
    "The folder of this exercise contains the skeleton of an input file for a 4-lead (RUN-4.fdf) and 6-lead device (RUN-6.fdf), as well as a script to run TBtrans for all values of the magnetic field (run.sh).\n",
    "\n",
    "*Depending on the size of the Hall bar, this step might require a considerable amount of time.*\n",
    "\n",
    "### 4. Extract the Hall resistance $(R_H)$\n",
    "\n",
    "In Hall effect we measure the build up of a potential difference between the measurement electrodes as response to an electric current. The Hall resistance ($R_H$) in a 4 lead Hall bar like the one shown above is given by\n",
    "\n",
    "$$R_H = \\left.\\frac{V_2-V_3}{I_1}\\right|_{I_2 = I_3 = 0}.$$\n",
    "\n",
    "In `TBtrans` the chemical potentials of all electrodes can be specified and the currents are calculated as a response to the biases. Rather than trying many combinations of chemical potentials to find the Hall resistance, we use transmission curves to calculate $R_H$.\n",
    "\n",
    "To start, we express the lead currents $I_i$ in terms of applied biases $V_i$ and the transmissions $T_{ij}$ between leads $i$ and $j$\n",
    "$$I_i = \\sum_j G_{ij} (V_i - V_j)\\quad\\text{where}\\quad G_{ij} = \\frac{2e^2}{h} T_{ij}.$$\n",
    "We rewrite this relation as\n",
    "\n",
    "$$\\mathbf{I} = \\mathcal{G} \\mathbf{V}\n",
    "\\quad\\text{, where}\\quad\n",
    "\\mathcal{G}_{ii} = \\sum_{i\\neq j} G_{ij}\n",
    "\\quad\\text{and}\\quad\n",
    "\\mathcal{G}_{ij} = - G_{ij}.$$\n",
    "\n",
    "Since the currents only depend on bias differences, we can set one of them to zero without loss of generality (here $V_4 = 0$ ). Further, Kirchhoff's current law allow us to eliminate one of current (here $I_4=-I_1-I_2-I_3$). This leaves us with an invertible $3\\times3$ matrix equation. \n",
    "\n",
    "Using the inverse $\\mathbf{R}$ of $\\mathcal{G}$, we can express $V_2$ and $V_3$ in terms of the lead currents $I_i$ and calculate the Hall conductance:\n",
    "\n",
    "$$\n",
    "\\mathbf{V} = \\mathcal{G}^{-1}\\mathbf{I} = \\mathbf{R} \\mathbf{I}\\quad\\Rightarrow\\quad\n",
    "V_i = R_{i1} I_1 + R_{i2} I_2 + R_{i3} I_3,\n",
    "$$\n",
    "\n",
    "and finally, we find the Hall resistance:\n",
    "\n",
    "\\begin{align}\n",
    "R_H &= \\left.\\frac{R_{21}I_1+R_{22}I_2 + R_{23}I_3-(R_{31}I_1+R_{32}I_2 + R_{33}I_3)}{I_1}\\right|_{I_2 = I_3 = 0} \\\\\n",
    "&= R_{21}-R_{31}\n",
    "\\end{align}\n",
    "\n",
    "The derivation for the 6-lead device is analogous and yields:\n",
    "$$\n",
    "R_H = R_{21}-R_{61}\n",
    "$$\n",
    "\n",
    "<br>\n",
    "<br>\n",
    "<br>\n",
    "\n",
    "*If everything is set up correctly, the quantization of the Hall resistance should be visible.*\n",
    "\\begin{align}\n",
    "R_H &= \\frac{h}{2ne^2}  &&\\text{for square lattice Hall bar} & n\\in\\mathbb{N}\\\\\n",
    "R_H &= \\frac{h}{2(2n-1)e^2} &&\\text{for graphene Hall bar}   & n\\in\\mathbb{N}\n",
    "\\end{align}\n",
    "\n",
    "### 5. Extract longitudinal resistance $R_L$\n",
    "The longitudinal resistance can be extracted using the same approach used for the Hall resistance. With a 6-lead Hall bar we can replace the $V_3$ with $V_6$ and get it immediately.\n",
    "\n",
    "$$\n",
    "R_L = R_{21} - R_{31}\n",
    "$$\n",
    "\n",
    "With a 4-lead Hall bar we need to create a new device with electrodes 2 and 3 on the same side of the Hall bar. \n",
    "\n",
    "*If the energy mesh in TBTrans and the mesh for magnetic field strength are fine engough, spikes in $R_L$ should be observable at each step of $R_H$.*"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "68eb4a61",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Create short-hand function to open files\n",
    "gs = sisl.get_sile\n",
    "# No magnetic field\n",
    "tbt0 = gs('M_0/siesta.TBT.nc')\n",
    "# All magnetic fields in increasing order\n",
    "tbts = [gs('M_{}/siesta.TBT.nc'.format(rec_phi)) for rec_phi in rec_phis]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "158d601b",
   "metadata": {},
   "outputs": [],
   "source": [
    "def calc_G(tbtsile, n):\n",
    "    # Construct G    \n",
    "    G = np.zeros((tbtsile.nE, n, n))\n",
    "    for i in range(n):\n",
    "        for j in range(n):\n",
    "            if i == j: continue\n",
    "            Gij = tbtsile.transmission(i,j)\n",
    "            G[:,i,j] = -Gij\n",
    "            G[:,i,i] += Gij\n",
    "    return G\n",
    "\n",
    "# Number of electrodes used in the Hall bar\n",
    "nEl = 6\n",
    "\n",
    "# Calculate Conductance matrix G\n",
    "G = np.asarray([calc_G(tbt, nEl) for tbt in tbts])\n",
    "\n",
    "# Remove one row and column from the matrix\n",
    "referenceElectrode = 3\n",
    "G = np.delete(np.delete(G, referenceElectrode, axis=2), referenceElectrode, axis=3)\n",
    "\n",
    "# Invert G to get resistance matrix R\n",
    "# Hint: If the matrix inversion might fail for some values of rec_phi, because the magentic field is too strong.\n",
    "R = np.linalg.inv(G)\n",
    "\n",
    "if nEl == 6:\n",
    "    # For 6-lead device \n",
    "    RH = R[:,:,1,0] - R[:,:,4,0]   # R_21 - R_61\n",
    "    RL = R[:,:,1,0] - R[:,:,2,0]   # R_21 - R_31\n",
    "elif nEl == 4:\n",
    "    # For 4-lead device\n",
    "    RH = R[:,:,1,0] - R[:,:,2,0]   # R_21 - R_31"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "f7f21a11",
   "metadata": {},
   "outputs": [],
   "source": [
    "# Plot the Hall resistance\n",
    "#  - as a function of energy for a fixed magnetic field strength\n",
    "#  - as a function of magnetic field strength for a fixed energy\n",
    "\n",
    "# E_idx = ...\n",
    "# phi_idx = ...\n",
    "\n",
    "fig, axs = plt.subplots(1,2, sharey=True, figsize=(12,4))\n",
    "\n",
    "axs[0].set_title('$\\phi = {:.5f}$'.format(1/rec_phis[phi_idx]))\n",
    "axs[0].set_ylabel(r'$R \\left[h/2e^2\\right]$')\n",
    "axs[0].set_xlabel('$E/t$')\n",
    "axs[0].axvline(x=E[E_idx], ls='--', c='k')\n",
    "# axs[0].plot(...)\n",
    "\n",
    "axs[1].set_title(f'E = {E[E_idx]:.5f}')\n",
    "axs[1].axvline(x=1/rec_phis[phi_idx], ls='--', c='k')\n",
    "axs[1].set_xlabel('$\\phi/\\phi_0$')\n",
    "# axs[1].plot(...)\n",
    "\n",
    "axs[1].legend(loc='center right', bbox_to_anchor=(1.5,0.5), ncol=1)\n",
    "\n",
    "plt.show()"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "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.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}