h2_quantum_espresso.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "#Shell command libraries\n",
    "import shlex, subprocess\n",
    "import re"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "#A useful function to run shell commands\n",
    "\n",
    "def run_shell_command(cmdline, output=True):\n",
    "    print(\"run shell command: \", cmdline, \"\\n\")\n",
    "    args = shlex.split(cmdline)\n",
    "    # let's use the new run command\n",
    "    p = subprocess.run(args, stdout=subprocess.PIPE, stderr=subprocess.STDOUT)\n",
    "    if output:\n",
    "        print(p.stdout.decode('utf-8'))\n",
    "        return p.stderr"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "run shell command:  ls \n",
      "\n",
      "1d_band_structure.ipynb\n",
      "2d_band_structure.ipynb\n",
      "ebt617e_env.yml\n",
      "ebt617e_ex1.ipynb\n",
      "h2\n",
      "h2_quantum_espresso.ipynb\n",
      "hydrogen-wave-function.ipynb\n",
      "LICENSE\n",
      "README.md\n",
      "spec-file_ebt617e_env.txt\n",
      "\n"
     ]
    }
   ],
   "source": [
    "#Testing with a simple command\n",
    "p = run_shell_command('ls', output=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Quantum espresso input file with a variable \"d\"\n",
    "def create_input(d=1, inputfile=\"h2/h2.in\"):\n",
    "    filedata = f\"\"\"\n",
    "    &control\n",
    "    calculation = 'scf'\n",
    "    prefix = 'h2'\n",
    "    pseudo_dir = './'\n",
    "    /\n",
    "    &system\n",
    "        ibrav=1, celldm(1) =20, \n",
    "        nat=2, ntyp=1,\n",
    "        ecutwfc=10.0\n",
    "    /\n",
    "    &electrons\n",
    "        conv_thr=1e-6\n",
    "    /\n",
    "    &ions\n",
    "        ion_dynamics='bfgs'\n",
    "    /\n",
    "    ATOMIC_SPECIES\n",
    "     H  1.00784  H.pbe-kjpaw_psl.1.0.0.UPF\n",
    "    ATOMIC_POSITIONS (angstrom)\n",
    "     H 0.00 0.00 0.00\n",
    "     H {d:.2} 0.00 0.00\n",
    "    K_POINTS (automatic)\n",
    "      1 1 1 0 0 0 \n",
    "    \"\"\"\n",
    "    \n",
    "    # Write out the file to be compiledb\n",
    "    with open(inputfile, 'w') as file:\n",
    "        file.write(filedata)\n",
    "    \n",
    "    return"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Test the create input function\n",
    "create_input(d=1.1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "    &control\n",
      "    calculation = 'scf'\n",
      "    prefix = 'h2'\n",
      "    pseudo_dir = './'\n",
      "    /\n",
      "    &system\n",
      "        ibrav=1, celldm(1) =20, \n",
      "        nat=2, ntyp=1,\n",
      "        ecutwfc=10.0\n",
      "    /\n",
      "    &electrons\n",
      "        conv_thr=1e-6\n",
      "    /\n",
      "    &ions\n",
      "        ion_dynamics='bfgs'\n",
      "    /\n",
      "    ATOMIC_SPECIES\n",
      "     H  1.00784  H.pbe-kjpaw_psl.1.0.0.UPF\n",
      "    ATOMIC_POSITIONS (angstrom)\n",
      "     H 0.00 0.00 0.00\n",
      "     H 1.1 0.00 0.00\n",
      "    K_POINTS (automatic)\n",
      "      1 1 1 0 0 0 \n",
      "    "
     ]
    }
   ],
   "source": [
    "!cat h2/h2.in"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Run Quantum Espresso\n",
    "def run_quantum_espresso(cmdline):\n",
    "    print(\"run command: \", cmdline, \"\\n\")\n",
    "    args = shlex.split(cmdline)\n",
    "    # let's use the new run command\n",
    "    p = subprocess.run(args, stdout=subprocess.PIPE, stderr=subprocess.STDOUT, cwd=\"h2/\")\n",
    "    return p"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [],
   "source": [
    "infile = 'h2.in'\n",
    "cmdline = 'pw.x -inp '+ infile"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "'pw.x -inp h2.in'"
      ]
     },
     "execution_count": 20,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "cmdline"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "run command:  pw.x -inp h2.in \n",
      "\n"
     ]
    }
   ],
   "source": [
    "out = run_quantum_espresso(cmdline)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {},
   "outputs": [],
   "source": [
    "outstr = out.stdout.decode()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "     Program PWSCF v.6.7GPU starts on  1Apr2021 at 10: 0:54 \n",
      "\n",
      "     This program is part of the open-source Quantum ESPRESSO suite\n",
      "     for quantum simulation of materials; please cite\n",
      "         \"P. Giannozzi et al., J. Phys.:Condens. Matter 21 395502 (2009);\n",
      "         \"P. Giannozzi et al., J. Phys.:Condens. Matter 29 465901 (2017);\n",
      "         \"P. Giannozzi et al., J. Chem. Phys. 152 154105 (2020);\n",
      "          URL http://www.quantum-espresso.org\", \n",
      "     in publications or presentations arising from this work. More details at\n",
      "     http://www.quantum-espresso.org/quote\n",
      "\n",
      "     Parallel version (MPI), running on     1 processors\n",
      "\n",
      "     MPI processes distributed on     1 nodes\n",
      "     Reading input from h2.in\n",
      "\n",
      "     Current dimensions of program PWSCF are:\n",
      "     Max number of different atomic species (ntypx) = 10\n",
      "     Max number of k-points (npk) =  40000\n",
      "     Max angular momentum in pseudopotentials (lmaxx) =  3\n",
      " end of file reached, closing tag not found\n",
      "     file H.pbe-kjpaw_psl.1.0.0.UPF: wavefunction(s)  1S renormalized\n",
      "\n",
      "     Subspace diagonalization in iterative solution of the eigenvalue problem:\n",
      "     a serial algorithm will be used\n",
      "\n",
      " \n",
      "     G-vector sticks info\n",
      "     --------------------\n",
      "     sticks:   dense  smooth     PW     G-vecs:    dense   smooth      PW\n",
      "     Sum        1281    1281    325                34265    34265    4337\n",
      " \n",
      "     Using Slab Decomposition\n",
      " \n",
      "\n",
      "\n",
      "     bravais-lattice index     =            1\n",
      "     lattice parameter (alat)  =      20.0000  a.u.\n",
      "     unit-cell volume          =    8000.0000 (a.u.)^3\n",
      "     number of atoms/cell      =            2\n",
      "     number of atomic types    =            1\n",
      "     number of electrons       =         2.00\n",
      "     number of Kohn-Sham states=            1\n",
      "     kinetic-energy cutoff     =      10.0000  Ry\n",
      "     charge density cutoff     =      40.0000  Ry\n",
      "     scf convergence threshold =      1.0E-06\n",
      "     mixing beta               =       0.7000\n",
      "     number of iterations used =            8  plain     mixing\n",
      "     Exchange-correlation= SLA  PW   PBX  PBC\n",
      "                           (   1   4   3   4   0   0   0)\n",
      "\n",
      "     celldm(1)=  20.000000  celldm(2)=   0.000000  celldm(3)=   0.000000\n",
      "     celldm(4)=   0.000000  celldm(5)=   0.000000  celldm(6)=   0.000000\n",
      "\n",
      "     crystal axes: (cart. coord. in units of alat)\n",
      "               a(1) = (   1.000000   0.000000   0.000000 )  \n",
      "               a(2) = (   0.000000   1.000000   0.000000 )  \n",
      "               a(3) = (   0.000000   0.000000   1.000000 )  \n",
      "\n",
      "     reciprocal axes: (cart. coord. in units 2 pi/alat)\n",
      "               b(1) = (  1.000000  0.000000  0.000000 )  \n",
      "               b(2) = (  0.000000  1.000000  0.000000 )  \n",
      "               b(3) = (  0.000000  0.000000  1.000000 )  \n",
      "\n",
      "\n",
      "     PseudoPot. # 1 for H  read from file:\n",
      "     ./H.pbe-kjpaw_psl.1.0.0.UPF\n",
      "     MD5 check sum: d3667862644f9dbf5ba3e2108d7283f9\n",
      "     Pseudo is Projector augmented-wave, Zval =  1.0\n",
      "     Generated using &quot;atomic&quot; code by A. Dal Corso  v.6.3MaX\n",
      "     Shape of augmentation charge: PSQ\n",
      "     Using radial grid of  929 points,  2 beta functions with: \n",
      "                l(1) =   0\n",
      "                l(2) =   0\n",
      "     Q(r) pseudized with 0 coefficients \n",
      "\n",
      "\n",
      "     atomic species   valence    mass     pseudopotential\n",
      "        H              1.00     1.00784     H ( 1.00)\n",
      "\n",
      "      8 Sym. Ops. (no inversion) found\n",
      "\n",
      "\n",
      "\n",
      "   Cartesian axes\n",
      "\n",
      "     site n.     atom                  positions (alat units)\n",
      "         1           H   tau(   1) = (   0.0000000   0.0000000   0.0000000  )\n",
      "         2           H   tau(   2) = (   0.1039349   0.0000000   0.0000000  )\n",
      "\n",
      "     number of k points=     1\n",
      "                       cart. coord. in units 2pi/alat\n",
      "        k(    1) = (   0.0000000   0.0000000   0.0000000), wk =   2.0000000\n",
      "\n",
      "     Dense  grid:    34265 G-vectors     FFT dimensions: (  45,  45,  45)\n",
      "\n",
      "     Estimated max dynamical RAM per process >      25.03 MB\n",
      "\n",
      "     Initial potential from superposition of free atoms\n",
      "\n",
      "     starting charge    1.99998, renormalised to    2.00000\n",
      "\n",
      "     negative rho (up, down):  1.158E-01 0.000E+00\n",
      "     Starting wfcs are    2 randomized atomic wfcs\n",
      "     Checking if some PAW data can be deallocated... \n",
      "\n",
      "     total cpu time spent up to now is        0.4 secs\n",
      "\n",
      "     Self-consistent Calculation\n",
      "\n",
      "     iteration #  1     ecut=    10.00 Ry     beta= 0.70\n",
      "     Davidson diagonalization with overlap\n",
      "     ethr =  1.00E-02,  avg # of iterations =  2.0\n",
      "\n",
      "     negative rho (up, down):  7.557E-02 0.000E+00\n",
      "\n",
      "     total cpu time spent up to now is        0.5 secs\n",
      "\n",
      "     total energy              =      -2.16379671 Ry\n",
      "     estimated scf accuracy    <       0.03194867 Ry\n",
      "\n",
      "     iteration #  2     ecut=    10.00 Ry     beta= 0.70\n",
      "     Davidson diagonalization with overlap\n",
      "     ethr =  1.60E-03,  avg # of iterations =  1.0\n",
      "\n",
      "     negative rho (up, down):  5.916E-02 0.000E+00\n",
      "\n",
      "     total cpu time spent up to now is        0.6 secs\n",
      "\n",
      "     total energy              =      -2.16505268 Ry\n",
      "     estimated scf accuracy    <       0.00116559 Ry\n",
      "\n",
      "     iteration #  3     ecut=    10.00 Ry     beta= 0.70\n",
      "     Davidson diagonalization with overlap\n",
      "     ethr =  5.83E-05,  avg # of iterations =  2.0\n",
      "\n",
      "     negative rho (up, down):  5.523E-02 0.000E+00\n",
      "\n",
      "     total cpu time spent up to now is        0.6 secs\n",
      "\n",
      "     total energy              =      -2.16521795 Ry\n",
      "     estimated scf accuracy    <       0.00003104 Ry\n",
      "\n",
      "     iteration #  4     ecut=    10.00 Ry     beta= 0.70\n",
      "     Davidson diagonalization with overlap\n",
      "     ethr =  1.55E-06,  avg # of iterations =  2.0\n",
      "\n",
      "     negative rho (up, down):  5.537E-02 0.000E+00\n",
      "\n",
      "     total cpu time spent up to now is        0.6 secs\n",
      "\n",
      "     total energy              =      -2.16521814 Ry\n",
      "     estimated scf accuracy    <       0.00000592 Ry\n",
      "\n",
      "     iteration #  5     ecut=    10.00 Ry     beta= 0.70\n",
      "     Davidson diagonalization with overlap\n",
      "     ethr =  2.96E-07,  avg # of iterations =  1.0\n",
      "\n",
      "     negative rho (up, down):  5.536E-02 0.000E+00\n",
      "\n",
      "     total cpu time spent up to now is        0.7 secs\n",
      "\n",
      "     End of self-consistent calculation\n",
      "\n",
      "          k = 0.0000 0.0000 0.0000 (  4337 PWs)   bands (ev):\n",
      "\n",
      "    -8.7708\n",
      "\n",
      "     highest occupied level (ev):    -8.7708\n",
      "\n",
      "!    total energy              =      -2.16521868 Ry\n",
      "\n",
      "     total all-electron energy =        -2.165210 Ry\n",
      "     estimated scf accuracy    <          1.6E-09 Ry\n",
      "\n",
      "     The total energy is the sum of the following terms:\n",
      "     one-electron contribution =      -3.00761466 Ry\n",
      "     hartree contribution      =       1.60426958 Ry\n",
      "     xc contribution           =      -1.15863407 Ry\n",
      "     ewald contribution        =       0.39697969 Ry\n",
      "     one-center paw contrib.   =      -0.00021922 Ry\n",
      "\n",
      "     convergence has been achieved in   5 iterations\n",
      "\n",
      "     Writing output data file ./h2.save/\n",
      " \n",
      "     init_run     :      0.09s CPU      0.27s WALL (       1 calls)\n",
      "     electrons    :      0.22s CPU      0.25s WALL (       1 calls)\n",
      "\n",
      "     Called by init_run:\n",
      "     wfcinit      :      0.01s CPU      0.06s WALL (       1 calls)\n",
      "     potinit      :      0.04s CPU      0.08s WALL (       1 calls)\n",
      "     hinit0       :      0.01s CPU      0.06s WALL (       1 calls)\n",
      "\n",
      "     Called by electrons:\n",
      "     c_bands      :      0.03s CPU      0.04s WALL (       5 calls)\n",
      "     sum_band     :      0.04s CPU      0.04s WALL (       5 calls)\n",
      "     v_of_rho     :      0.16s CPU      0.17s WALL (       6 calls)\n",
      "     newd         :      0.01s CPU      0.02s WALL (       6 calls)\n",
      "     PAW_pot      :      0.00s CPU      0.02s WALL (       6 calls)\n",
      "     mix_rho      :      0.01s CPU      0.02s WALL (       5 calls)\n",
      "\n",
      "     Called by c_bands:\n",
      "     init_us_2    :      0.00s CPU      0.00s WALL (      11 calls)\n",
      "     cegterg      :      0.03s CPU      0.04s WALL (       5 calls)\n",
      "\n",
      "     Called by *egterg:\n",
      "     cdiaghg      :      0.00s CPU      0.04s WALL (      13 calls)\n",
      "     h_psi        :      0.03s CPU      0.04s WALL (      14 calls)\n",
      "     s_psi        :      0.00s CPU      0.00s WALL (      14 calls)\n",
      "     g_psi        :      0.00s CPU      0.00s WALL (       8 calls)\n",
      "\n",
      "     Called by h_psi:\n",
      "     h_psi:calbec :      0.00s CPU      0.00s WALL (      14 calls)\n",
      "     vloc_psi     :      0.03s CPU      0.03s WALL (      14 calls)\n",
      "     add_vuspsi   :      0.00s CPU      0.00s WALL (      14 calls)\n",
      "\n",
      "     General routines\n",
      "     calbec       :      0.00s CPU      0.00s WALL (      19 calls)\n",
      "     fft          :      0.04s CPU      0.09s WALL (      65 calls)\n",
      "     ffts         :      0.00s CPU      0.00s WALL (       5 calls)\n",
      "     fftw         :      0.03s CPU      0.03s WALL (      35 calls)\n",
      " \n",
      "     Parallel routines\n",
      " \n",
      "     PWSCF        :      0.37s CPU      0.68s WALL\n",
      "\n",
      " \n",
      "   This run was terminated on:  10: 0:55   1Apr2021            \n",
      "\n",
      "=------------------------------------------------------------------------------=\n",
      "   JOB DONE.\n",
      "=------------------------------------------------------------------------------=\n",
      "\n"
     ]
    }
   ],
   "source": [
    "print( outstr )"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "-2.16521868"
      ]
     },
     "execution_count": 24,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "float(re.search(\"!    total energy              =      (.*)Ry\", outstr).group(0).split()[-2])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 34,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "run command:  pw.x -inp h2.in \n",
      "\n",
      "run command:  pw.x -inp h2.in \n",
      "\n",
      "run command:  pw.x -inp h2.in \n",
      "\n",
      "run command:  pw.x -inp h2.in \n",
      "\n",
      "run command:  pw.x -inp h2.in \n",
      "\n",
      "run command:  pw.x -inp h2.in \n",
      "\n",
      "run command:  pw.x -inp h2.in \n",
      "\n",
      "run command:  pw.x -inp h2.in \n",
      "\n",
      "run command:  pw.x -inp h2.in \n",
      "\n",
      "run command:  pw.x -inp h2.in \n",
      "\n",
      "run command:  pw.x -inp h2.in \n",
      "\n"
     ]
    }
   ],
   "source": [
    "dlist = np.linspace(0.4,1.4,11)\n",
    "\n",
    "infile = 'h2.in'\n",
    "\n",
    "cmdline = 'pw.x -inp '+ infile\n",
    "\n",
    "elist=[]\n",
    "for d in dlist:\n",
    "    create_input(d=d)\n",
    "    out = run_quantum_espresso(cmdline)\n",
    "    outstr = out.stdout.decode()\n",
    "    energy = float(re.search(\"!    total energy              =      (.*)Ry\", outstr).group(0).split()[-2])\n",
    "    elist.append(energy)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 35,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[<matplotlib.lines.Line2D at 0x7fdbadf20af0>]"
      ]
     },
     "execution_count": 35,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXwAAAD4CAYAAADvsV2wAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8vihELAAAACXBIWXMAAAsTAAALEwEAmpwYAAAh7ElEQVR4nO3deXTU533v8fdXG5JAICGBhARiGRaBMYuRMYg4doIX7DgGu3E2J81yEx/ftLlJk7Z24uTmnDbJqZulaU/am7pObrr4OnEbszjewI1jOwbZZhGY3SADloRAEkIsArTM9/4hQWQiIcGM5qeZ+bzO0YGZeTS/7zPGH/30/J7n+Zm7IyIiiS8l6AJERCQ2FPgiIklCgS8ikiQU+CIiSUKBLyKSJNKCLuBSCgoKfNKkSUGXISISNzZt2tTo7mN6e21IB/6kSZPYuHFj0GWIiMQNMzvY12sa0hERSRIKfBGRJKHAFxFJEgp8EZEkEVHgm9k9ZrbDzMJmVt5HmxlmVtXj64SZfTmS44qIyOWLdJbOduBu4J/7auDue4B5AGaWCtQCKyM8roiIXKaIAt/ddwGY2UC/ZSmw3937nDYUDZsONlNZ3cSiKfksmJg3mIcSEYkbsZ6H/1Hg8Us1MLP7gPsASktLL/sAldWNfOLR1wm7k5GWwmOfW6TQFxFhAGP4ZvaCmW3v5Wv55RzIzDKAO4H/vFQ7d3/E3cvdvXzMmF4Xi13SxgPNdISdsEN7R5jK6qbLfg8RkUTU7xm+u98UpWPdBmx29yNRer9eLQ4VkGJ7CTukp6WwaEr+YB5ORCRuxHJa5sfoZzgnGhZMzOPz108B4G8/NEfDOSIi3SKdlnmXmdUAi4Gnzez57ueLzeyZHu2ygZuBJyM53kB9bGHX2H/LmY5YHE5EJC5EFPjuvtLdx7v7MHcvdPdbu5+vc/fbe7Rrdfd8d2+JtOCBmJifTfGoTDbsb4zF4URE4kJCrrQ1MxaHCtiwv4lwWDdpFxGBBA18gIpQPs2t7eyuPxl0KSIiQ0LiBv7Urtk56zWsIyICJHDgjxuVxZSC4azfr3n4IiKQwIEPsDiUz2vVTbR3hoMuRUQkcAkd+BWhAk63dfJmbUwmB4mIDGkJHfiLpowGYP0+jeOLiCR04OePGEZZUY7G8UVESPDAB1gytYCNB5s5294ZdCkiIoFK+MCvCOXT1hFm86HmoEsREQlUwgf+wsmjSU0xNmhYR0SSXMIHfk5mOleXjOJVXbgVkSSX8IEPsGRqPltrWjh1TrtnikjySorArwgV0Bl23nj7WNCliIgEJikCf8HEPDJSU7SvjogktaQI/Mz0VK6ZmKv5+CKS1JIi8KFrWGfn4RM0n24LuhQRkUAkTeAvmZqPO1RW6yxfRJJT0gT+nPG5ZGekalhHRJJW0gR+emoKCyeP1oVbEUlaSRP40LXNwv6G09S3nA26FBGRmEuywC8AYEO1zvJFJPkkVeDPGjeSUVnprN+ncXwRST5JFfgpKcbiKfms39+EuwddjohITCVV4ANUTM2n9vgZ3jl2JuhSRERiKvkCP5QPwKuarSMiSSbpAj80ZgRjc4ZpPr6IJJ2kC3wzoyKUz4b9jRrHF5GkElHgm9k9ZrbDzMJmVn6Jdn/W3W67mT1uZpmRHDdSFaECGk+18dbRU0GWISISU5Ge4W8H7gZe7quBmZUA/wsod/fZQCrw0QiPG5HF3eP463UXLBFJIhEFvrvvcvc9A2iaBmSZWRqQDdRFctxITRidzYTRWbyqcXwRSSKDPobv7rXA94FDwGGgxd3X9tXezO4zs41mtrGhoWHQ6loSKqCyuonOsMbxRSQ59Bv4ZvZC99j7xV/LB3IAM8sDlgOTgWJguJl9oq/27v6Iu5e7e/mYMWMG2o/LtjiUz8mzHeyoaxm0Y4iIDCVp/TVw95siPMZNwNvu3gBgZk8CFcB/RPi+Ebkwjr+/iTnjc4MsRUQkJmIxLfMQsMjMss3MgKXArhgc95LG5mQybewIXtWFWxFJEpFOy7zLzGqAxcDTZvZ89/PFZvYMgLu/BvwXsBl4s/uYj0RUdZRUhPJ548Ax2jrCQZciIjLoIp2ls9Ldx7v7MHcvdPdbu5+vc/fbe7T7lruXuftsd/+ku5+LtPBoqJhawNn2MFXvHA+6FBGRQZd0K217WjQ5HzN0FywRSQpJHfijstOZXTxK++qISFJI6sCHrnH8LYeaaW3rCLoUEZFBpcCfWkB7p7PxQHPQpYiIDKqkD/xrJ+WRlmIa1hGRhJf0gZ+dkcb80lw26MKtiCS4pA98gMWhAt6sbaHlTHvQpYiIDBoFPl0XbsMOr1VrWEdEEpcCH5hfmktmeorG8UUkoSnwgWFpqVw7aTQbFPgiksAU+N0Wh/LZc+QkDSeHxK4PIiJRp8DvVhEqAGCDxvFFJEEp8LvNLh5JTmaapmeKSMJS4HdLS03husn5unArIglLgd9DRSifg02t1DS3Bl2KiEjUKfB7qJjaddtDzdYRkUSkwO9h+tgc8odnaFhHRBKSAr+HlBRjcSif9fsbcfegyxERiSoF/kUqQgUcOXGO6sbTQZciIhJVCvyLVIS6xvE1rCMiiUaBf5GJ+dkUj8rUfHwRSTgK/IuYGYtDBWzY30Q4rHF8EUkcCvxeLJmaT3NrO7vqTwRdiohI1Cjwe7E4pPn4IpJ4FPi9GDcqiykFw3XhVkQSigK/D4tD+bxW3UR7ZzjoUkREokKB34clUws43dbJtpqWoEsREYkKBX4fFk05P46v6ZkikhgiCnwzu8fMdphZ2MzKL9HuS2a2vbvtlyM5ZqyMHp7BzHEjNY4vIgkj0jP87cDdwMt9NTCz2cDngYXAXOAOM5sW4XFjoiKUz8aDzZxt7wy6FBGRiEUU+O6+y9339NNsJlDp7q3u3gG8BNwVyXFjpSKUT1tHmM0Hm4MuRUQkYrEYw98OvNfM8s0sG7gdmNBXYzO7z8w2mtnGhoaGGJTXt4WTR5OaYhrWEZGEkNZfAzN7ASjq5aWH3H11f9/v7rvM7GFgHXAK2Ap0XKL9I8AjAOXl5YHubZCTmc6c8aNYv78RmBFkKSIiEes38N39pkgP4u4/BX4KYGbfBWoifc9YqQjl85OXqjl1roMRw/r9uEREhqyYTMs0s7Hdf5bSdZH38VgcNxoqQgV0hp3X39awjojEt0inZd5lZjXAYuBpM3u++/liM3umR9NfmdlO4CngT9w9bq6CLpiYR0ZaCuv3KfBFJL5FNEbh7iuBlb08X0fXxdnzj6+P5DhBykxPZUFpni7cikjc00rbAagI5bPz8AmaT7cFXYqIyBVT4A9AxdSubRYqq3WWLyLxS4E/AHPG5zI8I5VXta+OiMQxBf4ApKemsHDyaI3ji0hcU+APUEWogOqG09S3nA26FBGRK6LAH6ALtz2s1rCOiMQnBf4AzRo3klFZ6ZqPLyJxS4E/QCkpxuIp+azf34R7oFv8iIhcEQX+ZVgyNZ/a42c4dKw16FJERC6bAv8yLA4VAGi2jojEJQX+ZQiNGc7YnGEKfBGJSwr8y2BmVITy2bC/UeP4IhJ3FPiXqSJUQOOpNvYeORV0KSIil0WBf5nO76uzXtssiEicUeBfpvF52ZSOztY4vojEHQX+FagI5VNZ3URnWOP4IhI/FPhXYHEon5NnO9he2xJ0KSIiA6bAvwIVmo8vInFIgX8FxuQMY3rhCF24FZG4osC/QhWhAt44cIy2jnDQpYiIDIgC/wotDuVztj1M1TvHgy5FRGRAFPhXaNHkfFIMXt2nYR0RiQ8K/Cs0Kjud2SWj2KALtyISJxT4EVgcymfLO820tnUEXYqISL8U+BGoCBXQ3ulsPNAcdCkiIv1S4Efg2kl5pKUYr2p6pojEAQV+BLIz0phfmqtxfBGJCwr8CFWECthe20JLa3vQpYiIXFJEgW9m3zOz3Wa2zcxWmlluH+2WmdkeM9tnZg9GcsyhpiKUT9jhtbd1li8iQ1ukZ/jrgNnuPgfYC3zt4gZmlgr8I3AbMAv4mJnNivC4Q8a80lwy01O0r46IDHkRBb67r3X383MSK4HxvTRbCOxz92p3bwN+ASyP5LhDybC0VK6dNFr76ojIkBfNMfzPAs/28nwJ8E6PxzXdzyWMilABe4+couHkuaBLERHpU7+Bb2YvmNn2Xr6W92jzENABPNbbW/TyXJ93DjGz+8xso5ltbGhoGEgfAlcR6rrt4YZqDeuIyNCV1l8Dd7/pUq+b2aeAO4Cl7t5bkNcAE3o8Hg/UXeJ4jwCPAJSXl8fFLaWuKh5JTmYaG/Y3cufc4qDLERHpVaSzdJYBDwB3untrH83eAKaZ2WQzywA+CqyJ5LhDTVpqCtdNzteFWxEZ0iIdw/8xkAOsM7MqM/sJgJkVm9kzAN0Xdf8UeB7YBTzh7jsiPO6QUxHK52BTKzXNff3cExEJVr9DOpfi7lP7eL4OuL3H42eAZyI51lC3ZOrvb3v44fLsgKsREflDWmkbJdMLRzAyM42fv3qATQe1mZqIDD0K/CjZfOg4p851sPPwCe59tFKhLyJDjgI/Siqrmzg/R+lce5hKTdEUkSFGgR8li6bkMyy96+N0YM74UcEWJCJyEQV+lCyYmMdjn1vEJxdNBNAZvogMORHN0pF3WzAxjwUT8zh5tp1HX3mbe6+bSHFuVtBliYgAOsMfFH9+6wwc+P7ze4IuRUTkAgX+IBifl81nl0zmyS21bK9tCbocERFAgT9ovvC+EHnZ6Xz76Z30vsWQiEhsKfAHycjMdL5803Qqq4/xm91Hgy5HRESBP5g+fl0pUwqG891ndtHRGQ66HBFJcgr8QZSemsKDt5Wxv+E0j7/xTv/fICIyiBT4g+zmWYUsnDyaH63by8mz7UGXIyJJTIE/yMyMb3xgJk2n2/jJS/uDLkdEkpgCPwbmjM9lxbxiHn3lbeqOnwm6HBFJUgr8GLmwGGutFmOJSDAU+DFyYTHWZi3GEpFgKPBjSIuxRCRICvwY0mIsEQmSAj/GtBhLRIKiwI+xnouxfqHFWCISQwr8AJxfjPV3WowlIjGkwA+AFmOJSF82HWzmH/77LTYdbI76eyvwA6LFWCLSU0trO3/73C7u+cl6frhuL/c+Whn10FfgB0iLsUSSW2tbB2u21vG5f93Itd95gX/6bTXh7hnb7R3hqN8bW/e0DdD5xVj//PJ+PrtkMrNLRgVdkogMsraOMC/vbWDN1jpe2HWE1rZOCkcO45OLJzK9cATfWrOD9o4w6WkpLJqSH9VjK/AD9oX3hfjlG4f4ztO7+H+fvw4zC7okEYmyzrDzWnUTa7bW8ez2elrOtJObnc7yeSXcObeYhZNHk5rS9f/+1LE5VFY3sWhKPgsm5kW1DgV+wM4vxvrWmh38ZvdRls4sDLokEYkCd6fqneOs2VrH09sOc/TkObIzUrllViF3zivmPVPHkJH2h6PqCybmRT3oz1PgDwEfv66Uf11/gO8+s4sbpo8hLVWXVkTi1e76E6ypquOpbXW8c+wMGakp3DhjDHfOK2ZpWSFZGamB1RZR4JvZ94APAm3AfuAz7n68l3Y/A+4Ajrr77EiOmYjOL8a679838Ys33uETiyYGXZKIXIZDTa2s2VrLmq117D1yihSDJVML+OL7p3HrVUWMykoPukQg8jP8dcDX3L3DzB4GvgY80Eu7nwM/Bv4twuMlrAt3xnphL8vnFZOTOTT+gYhI746eOMtT2w6zZmsdW985DkD5xDz+avlV3H71OApGDAu2wF5EFPjuvrbHw0rgQ320e9nMJkVyrERnZjx0+0yW/+Or/PNL1fz5rTOCLklELnK8tY1nt9ezpqqOyrebcIdZ40by4G1l3DFnHOPzsoMu8ZKiOYb/WeCXkb6Jmd0H3AdQWloa6dvFlbkTclk+r5h/eaWaj19XSnFuVtAliSStTQebqaxuYt6EXBpPnWNNVR0vv9VAe6czuWA4X3z/NO6cW8zUsSOCLnXA+g18M3sBKOrlpYfcfXV3m4eADuCxSAty90eARwDKy8uTbtP4v7h1Bs9ur+f7a/fwww/PC7ockaRUWd3IJ3/6Ou2dv4+gopGZfLpiEnfOLWF2yci4nELdb+C7+02Xet3MPkXXBdmlrrt6REyLsUSCEQ47rx84xuqqWp7cXHsh7A34yLUT+O5dV5OSEn8h31Oks3SW0XWR9gZ3b41OSaLFWCKxs+vwCVZV1fJUVR11LWfJSk/lusmjqaxuojPspKelcE/5hLgPe4h8DP/HwDBgXXcoVbr7/WZWDDzq7rcDmNnjwI1AgZnVAN9y959GeOyE1XMx1ot7jvL+Mi3GEommmuZW1mytY/WWOvYcOUlqivHeaQU8cFsZN88qJDsj7cIY/mCseA2KDeVRmPLyct+4cWPQZQSivTPMrX/3MikpxnNful6LsUQidLy1jaffPMzqLXW8fuAYANeU5rJifgkfuHoc+UNwGuWVMLNN7l7e22taaTtEaTGWSOTOtHXywq4jrK6q46W9R2nvdEJjhvPVm6ezfF4JpflDexpltCnwhzAtxhK5fB2dYdbvb2JVVS3Pb6/ndPdulJ+umMTyeSVcVRyfM2yiQYE/hGkxlsjAuDvbalq6Lr5uPUzjqXPkZKbxgTnjWDGvhOum5F/YjTKZKfCHuJ6Lse5dVMq4UVqMJXLegcbTrKqqZXVVHW83niYjNYX3l41lxfxibpwxlsz04DYqG4oU+HHg/GKs7z2vxVgiDSfP8ettdayq6trDxgwWTc7n/humsGz2uCGzUdlQpMCPA1qMJcnq/NTIuRNyOdJyllVVtby6r5Fw9x42X7+9jA/OLdZvvgOkwI8TWowlyeb1t5u499HX3rW9wfi8LL5w41SWzytmWmFOgNXFJwV+nNBiLEkG7s7WmhZWbq7hiY0179re4OPXlfLtFbN1shMBBX4c+f2dsXbz3mm6M5YkjkNNrayqqmXVllqqG0+TkZZC+cQ83jhwjHD39gZ3XzNeYR8hBX4c0WIsSSTHW9v49bbDrNpSy8aDzT0uvoZYdnURIzPTE3J7gyAp8OOMFmNJPDvX0cmLu4/y5OZaXtzTtfJ1euEIHlhWxvJ5xX9wD4jBvKF3MlLgxxktxpJ4Ew47mw418+TmWp7eVseJsx2MyRnGpxZPYsX85F75GmsK/DikxVgSD/Y3nGLVllpWbqmlpvkMWempLJtdxF3zS6gI5esaVAAU+HHqwp2xnt/LDz48N+hyRABoPHWOp7bWsXJLLdtqWkgxeM+0MXz1luncMquI4cMUOUHSpx+nei7G+sySSVqMJYE509bJul1HWLm5hpffaqQz7FxVPJJvfGAmd84tZuzIzKBLlG4K/DimxVgSlM6wU1ndxJOba3lu+2FOt3VSPCqT+947hbvmlzBdi6KGJAV+HNNiLImV89Mjx43KZE/9SVZX1VF/4iw5w9K4Y04xK+aXcN3k0QlxG8BEpsCPc+cXY31z9Q521J2gIlSgaWwSVWt31POFxzbTEe5a9Zpq8L6ysXzzjlksnakdKeOJAj/OpaemcE/5eB5+bg8/WLuXf0zbx2OfX6TQl4i0nGnn2TcPs6qqlsrqYxeeN+D+G0P8xa1lwRUnV0yBnwDC7hjgwNmOME9urlHgy2U72961KGpVVS0v7m6grTPMlILhfKR8AquqaunoDJOelqKhwzimwE8Ai6YUMCx9H20dYdzh8dcPMTYnkz95X0hzneWSzl98XbWllue213PyXNeiqE8smsiK+cVcXTIKM+PD107QFgcJwNy9/1YBKS8v940bNwZdRlw4f1FtzvhR/GpTDauq6lgwMY8ffWQeE0Yn142a5dLcne21J1hdVctT2+o4cuIcI4alsWx2ESvmlbA4pNsBxjMz2+Tu5b2+psBPTKuravnGyu048FfLr+Ku+SWatpnkDjadZnVVHauqaqluOE16qnHjjLGsmFeii68J5FKBryGdBLV8XgnXlObxlSeq+MoTW3lxTwPfXjFbt39LMo2nzvH0tq6Lr1sOHQfgusmj+fz1U7htdhG52RnBFigxpTP8BNcZdn7y0n7+bt1exuYM44cfmceiKflBlyWD6PS5DtburGfVljp+t69r5WtZUQ4r5pdw59w/3JFSEouGdISt7xzny7+s4kDTae6/IcSf3TSdjDRd0E0U7Z1hXnmrgVVb6li38whn2jspyc1i+bxils8rYUaRVr4mCw3pCHMn5PLrL76Hv/71Tv7Pb/fzu7ca+dFH5xEaMyLo0uQKuTubDjazqqqWp7cdprm1ndzsdO6+poQV80tYUJqnla/yLjrDT0LPba/nwSe3ca49zDfvmMXHFk7QBd04cH4m1vjcLPYe7dreoKb5DJnpKdw8q4gV84q5ftoY/eaW5AbtDN/Mvgd8EGgD9gOfcffjF7WZAPwbUASEgUfc/e8jOa5EZtnsIuaX5vLVJ7by9ZVv8uKeo/zN3VeTP2JY0KVJL9ydlVtq+cv/2nZhewMDrp8+hq/cPJ1bripihLYdlgGI6AzfzG4BfuPuHWb2MIC7P3BRm3HAOHffbGY5wCZghbvv7O/9dYY/uMJh52evvs3fPreHUdnpfP+eudwwfUzQZQnQ0Rlm08Fm1u48wrqdRzh0rPXCawb8zxtD/OUybW8gf2jQzvDdfW2Ph5XAh3ppcxg43P33k2a2CygB+g18GVwpKcbnrp9CRaiAL/9yC5/62et8ZskkHlhWpjnZATjT1snLbzWwbucR/nvXEZpb28lITWHJ1Hxuv7qI//vqgQvbGyydqe0N5PJFbQzfzJ4Cfunu/3GJNpOAl4HZ7n6ijzb3AfcBlJaWLjh48GBU6pNLO9veyd88u5ufrz/AjMIc/v5j8ygrGhl0WQmv6dQ5/nv3UdbuOMLv9jVwtj3MyMw0ls4s5OZZhbx3+pgLwzXnx/C1vYFcSkTTMs3sBbrG3y/2kLuv7m7zEFAO3O19vKGZjQBeAr7j7k8OpHAN6cTei3uO8hf/uY0TZ9t5cFkZn66YpJkeUXag8TTruodqNh48RtiheFQmt1xVxM2zClk4eTTp2gNJrtCgzsM3s08B9wNL3b21jzbpwK+B5939hwN9bwV+MBpPnePBX23jhV1HuX5aAT+4Z65uUxcBd2dbTQvrdh5h7c569h45BcDMcSO5eVYht8wq5KrikZopJVExaIFvZsuAHwI3uHtDH20M+FfgmLt/+XLeX4EfHHfnsdcO8e2nd5KVnsrDfzSHW67q7Rc96U1bR5jK6ibW7qznhZ1HqT9xltQUY+Gk0dw8q2u4RpvayWAYzMDfBwwDmrqfqnT3+82sGHjU3W83s/cArwBv0jUtE+Dr7v5Mf++vwA/evqOn+NIvtrCj7gQfW1jKN++YSXaGpgD25sTZdn67p+ui6293H+XkuQ6y0lO5YfoYbp5VyPvLxpI3XHvXyODS1goSkbaOMD9Yt4dHXq5mcv5wfvTRecwZnxt0WYHpefG0JDeLdbuOsHZHPZXVTbR3OgUjMlhaVsgtVxWyZGqBZjxJTCnwJSrW72/kq09speHkOf7s5uncf0Mo6fZNX7+vkU///A3aO8JgcP5/n8kFw7mle6hmfmle0n0uMnQo8CVqWlrb+frKN3n6zcMsnDyazy6ZxP6G0wk3VTAcdmqaz7Cr/gS7D59kd/0J9tSfpLrx9LvaXT+tgG99cBahMSN00VWGBAW+RJW786vNtXxj5Zuc7QhjQFqq8fAfzWHpzMK423O/5Uw7e+q7Qn1Xd7jvrT/J6bZOAMxg4uhsyopGMjIrjZVbagmHnfS0FB77nG4YL0OLdsuUqDIzPrRgPHvqT/Avr7yNA+2dzlee2ApAzrA0inOzKMnLojg3k5LcbIpzMxmfl0VxbhZjczIDGfLo6AzzduNpdtWfZPfhE+zu/rOu5eyFNqOy0ikryuGe8gmUFeVQNm4k0wtHvOtC9UeuLdUCKIlLCny5Ystmj+PfKw/S3hEmNTWFP33fVDLTU6g7fpaa5jPUHT/DpoPNtJxpf9f3paUYRaMyKcnN6vrq/kFQkvv7P7MyrvxCp7vTcOpc11n74ZMXhmX2HT1FW2f4Qg2hMSO4dvJoyopGUjYuh5lFIykcOazfoZkFE/MU9BKXFPhyxRZMzOOxzy3q92z31LkO6o6fofb4GWq7fxDUHu/6s7K6ifoTZwlfNLI4enhG9w+AP/wNoSQ3i9HDM9h86Di/e6uB4lFZYHSdsXeHe9PptgvvVThyGDOKRnL9tALKxuVQVjSS0JgR2kZYko7G8CVwHZ1hjpw8964fBhf/cGjtHk8/LyPVaO90ev7rzUxPYUZhV6DPKMq5EO6jNfddkojG8GVIS0tNuTC80xt3p+VM+7t+CDz95mHeONAMQIrBZ5dM5mu3z9R0SJFLUODLkGdm5GZnkJudwVXFowC4enwu9z5aSXtH13bBt109TmEv0g8FvsSlgV4/EJHfU+BL3NJsGZHLo2kKIiJJQoEvIpIkFPgiIklCgS8ikiQU+CIiSUKBLyKSJIb01gpm1gAcvMJvLwAao1hOPFCfE1+y9RfU58s10d3H9PbCkA78SJjZxr72k0hU6nPiS7b+gvocTRrSERFJEgp8EZEkkciB/0jQBQRAfU58ydZfUJ+jJmHH8EVE5N0S+QxfRER6UOCLiCSJuA58M1tmZnvMbJ+ZPXiJdteaWaeZfSiW9Q2GgfTZzG40syoz22FmL8W6xmjrr89mNsrMnjKzrd19/kwQdUaTmf3MzI6a2fY+Xjcz+4fuz2SbmV0T6xqjaQD9vbe7n9vMbL2ZzY11jdHWX597tItefrl7XH4BqcB+YAqQAWwFZvXR7jfAM8CHgq57sPsM5AI7gdLux2ODrjsGff468HD338cAx4CMoGuPsN/vBa4Btvfx+u3As4ABi4DXgq55kPtbAeR1//22eO/vQPrc3Saq+RXPZ/gLgX3uXu3ubcAvgOW9tPsi8CvgaCyLGyQD6fPHgSfd/RCAu8d7vwfSZwdyzMyAEXQFfkdsy4wud3+Zrn70ZTnwb96lEsg1s3GxqS76+uuvu6939+buh5XA+JgUNogG8N8Yopxf8Rz4JcA7PR7XdD93gZmVAHcBP4lhXYOp3z4D04E8M/utmW0ysz+OWXWDYyB9/jEwE6gD3gS+5O7h2JQXmIF8Lonqf9D1201CG4z8iudbHPZ2x+qL55j+CHjA3Tu7Tv7i3kD6nAYsAJYCWcAGM6t0972DXdwgGUifbwWqgPcDIWCdmb3i7icGubYgDeRzSThm9j66Av89QdcSAz8iyvkVz4FfA0zo8Xg8XWd4PZUDv+j+sAqA282sw91XxaTC6BtIn2uARnc/DZw2s5eBuUC8Bv5A+vwZ4G+8a9Bzn5m9DZQBr8emxEAM5HNJKGY2B3gUuM3dm4KuJwainl/xPKTzBjDNzCabWQbwUWBNzwbuPtndJ7n7JOC/gC/EcdjDAPoMrAauN7M0M8sGrgN2xbjOaBpInw/R9RsNZlYIzACqY1pl7K0B/rh7ts4ioMXdDwdd1GAxs1LgSeCTcfzb6mUZjPyK2zN8d+8wsz8FnqfrSvbP3H2Hmd3f/XqijNtfMJA+u/suM3sO2AaEgUfd/ZLTvoayAf53/mvg52b2Jl1DHQ+4e1xvp2tmjwM3AgVmVgN8C0iHC31+hq6ZOvuAVrp+y4lbA+jv/wbygX/qPuPt8DjfQXMAfY7+Mbun/oiISIKL5yEdERG5DAp8EZEkocAXEUkSCnwRkSShwBcRSRIKfBGRJKHAFxFJEv8fR1aVQZFzpJMAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.plot(dlist,elist,'.-')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "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.8.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}