v0.6.3/examples/gaas_surface.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "source": [
    "# Modelling a gallium arsenide surface\n",
    "\n",
    "This example shows how to use the\n",
    "[atomistic simulation environment](https://wiki.fysik.dtu.dk/ase/index.html),\n",
    "or ASE for short,\n",
    "to set up and run a particular calculation of a gallium arsenide surface.\n",
    "ASE is a Python package to simplify the process of setting up,\n",
    "running and analysing results from atomistic simulations across different simulation codes.\n",
    "By means of [ASEconvert](https://github.com/mfherbst/ASEconvert.jl) it is seamlessly\n",
    "integrated with the AtomsBase ecosystem and thus available to DFTK via our own\n",
    "AtomsBase integration.\n",
    "\n",
    "In this example we will consider modelling the (1, 1, 0) GaAs surface separated by vacuum."
   ],
   "metadata": {}
  },
  {
   "cell_type": "markdown",
   "source": [
    "Parameters of the calculation. Since this surface is far from easy to converge,\n",
    "we made the problem simpler by choosing a smaller `Ecut` and smaller values\n",
    "for `n_GaAs` and `n_vacuum`.\n",
    "More interesting settings are `Ecut = 15` and `n_GaAs = n_vacuum = 20`."
   ],
   "metadata": {}
  },
  {
   "outputs": [],
   "cell_type": "code",
   "source": [
    "miller = (1, 1, 0)   # Surface Miller indices\n",
    "n_GaAs = 2           # Number of GaAs layers\n",
    "n_vacuum = 4         # Number of vacuum layers\n",
    "Ecut = 5             # Hartree\n",
    "kgrid = (4, 4, 1);   # Monkhorst-Pack mesh"
   ],
   "metadata": {},
   "execution_count": 1
  },
  {
   "cell_type": "markdown",
   "source": [
    "Use ASE to build the structure:"
   ],
   "metadata": {}
  },
  {
   "outputs": [],
   "cell_type": "code",
   "source": [
    "using ASEconvert\n",
    "\n",
    "a = 5.6537  # GaAs lattice parameter in Ångström (because ASE uses Å as length unit)\n",
    "gaas = ase.build.bulk(\"GaAs\", \"zincblende\"; a)\n",
    "surface = ase.build.surface(gaas, miller, n_GaAs, 0, periodic=true);"
   ],
   "metadata": {},
   "execution_count": 2
  },
  {
   "cell_type": "markdown",
   "source": [
    "Get the amount of vacuum in Ångström we need to add"
   ],
   "metadata": {}
  },
  {
   "outputs": [],
   "cell_type": "code",
   "source": [
    "d_vacuum = maximum(maximum, surface.cell) / n_GaAs * n_vacuum\n",
    "surface = ase.build.surface(gaas, miller, n_GaAs, d_vacuum, periodic=true);"
   ],
   "metadata": {},
   "execution_count": 3
  },
  {
   "cell_type": "markdown",
   "source": [
    "Write an image of the surface and embed it as a nice illustration:"
   ],
   "metadata": {}
  },
  {
   "outputs": [
    {
     "output_type": "execute_result",
     "data": {
      "text/plain": "\u001b[0m\u001b[1mPython None\u001b[22m"
     },
     "metadata": {},
     "execution_count": 4
    }
   ],
   "cell_type": "code",
   "source": [
    "ase.io.write(\"surface.png\", surface * pytuple((3, 3, 1)), rotation=\"-90x, 30y, -75z\")"
   ],
   "metadata": {},
   "execution_count": 4
  },
  {
   "cell_type": "markdown",
   "source": [
    "<img src=\"https://docs.dftk.org/stable/surface.png\" width=500 height=500 />"
   ],
   "metadata": {}
  },
  {
   "cell_type": "markdown",
   "source": [
    "Use the `pyconvert` function from `PythonCall` to convert to an AtomsBase-compatible system.\n",
    "These two functions not only support importing ASE atoms into DFTK,\n",
    "but a few more third-party data structures as well.\n",
    "Typically the imported `atoms` use a bare Coulomb potential,\n",
    "such that appropriate pseudopotentials need to be attached in a post-step:"
   ],
   "metadata": {}
  },
  {
   "outputs": [
    {
     "output_type": "execute_result",
     "data": {
      "text/plain": "FlexibleSystem(As₂Ga₂, periodic = TTT):\n    bounding_box      : [ 3.99777        0        0;\n                                0  3.99777        0;\n                                0        0  21.2014]u\"Å\"\n\n    Atom(Ga, [       0,        0,  8.48055]u\"Å\")\n    Atom(As, [ 3.99777,  1.99888,  9.89397]u\"Å\")\n    Atom(Ga, [ 1.99888,  1.99888,  11.3074]u\"Å\")\n    Atom(As, [ 1.99888, 6.96512e-16,  12.7208]u\"Å\")\n"
     },
     "metadata": {},
     "execution_count": 5
    }
   ],
   "cell_type": "code",
   "source": [
    "using DFTK\n",
    "system = attach_psp(pyconvert(AbstractSystem, surface);\n",
    "                    Ga=\"hgh/pbe/ga-q3.hgh\",\n",
    "                    As=\"hgh/pbe/as-q5.hgh\")"
   ],
   "metadata": {},
   "execution_count": 5
  },
  {
   "cell_type": "markdown",
   "source": [
    "We model this surface with (quite large a) temperature of 0.01 Hartree\n",
    "to ease convergence. Try lowering the SCF convergence tolerance (`tol`)\n",
    "or the `temperature` or try `mixing=KerkerMixing()`\n",
    "to see the full challenge of this system."
   ],
   "metadata": {}
  },
  {
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "n     Energy            log10(ΔE)   log10(Δρ)   Diag   Δtime\n",
      "---   ---------------   ---------   ---------   ----   ------\n",
      "  1   -16.58798704964                   -0.58    5.0         \n",
      "  2   -16.72530584784       -0.86       -1.01    1.0    519ms\n",
      "  3   -16.73059540617       -2.28       -1.57    2.2    662ms\n",
      "  4   -16.73127301563       -3.17       -2.16    1.9    629ms\n",
      "  5   -16.73133330732       -4.22       -2.61    2.1    672ms\n",
      "  6   -16.73133601714       -5.57       -2.96    1.6    549ms\n",
      "  7   -16.73109187482   +   -3.61       -2.61    2.2    667ms\n",
      "  8   -16.73130566476       -3.67       -2.97    2.6    674ms\n",
      "  9   -16.73131992822       -4.85       -3.12    2.4    649ms\n",
      " 10   -16.73134018248       -4.69       -4.21    1.9    606ms\n"
     ]
    }
   ],
   "cell_type": "code",
   "source": [
    "model = model_DFT(system, [:gga_x_pbe, :gga_c_pbe],\n",
    "                  temperature=0.001, smearing=DFTK.Smearing.Gaussian())\n",
    "basis = PlaneWaveBasis(model; Ecut, kgrid)\n",
    "\n",
    "scfres = self_consistent_field(basis, tol=1e-4, mixing=LdosMixing());"
   ],
   "metadata": {},
   "execution_count": 6
  },
  {
   "outputs": [
    {
     "output_type": "execute_result",
     "data": {
      "text/plain": "Energy breakdown (in Ha):\n    Kinetic             5.8593901 \n    AtomicLocal         -105.6097400\n    AtomicNonlocal      2.3494860 \n    Ewald               35.5044300\n    PspCorrection       0.2016043 \n    Hartree             49.5611544\n    Xc                  -4.5976614\n    Entropy             -0.0000035\n\n    total               -16.731340182481"
     },
     "metadata": {},
     "execution_count": 7
    }
   ],
   "cell_type": "code",
   "source": [
    "scfres.energies"
   ],
   "metadata": {},
   "execution_count": 7
  }
 ],
 "nbformat_minor": 3,
 "metadata": {
  "language_info": {
   "file_extension": ".jl",
   "mimetype": "application/julia",
   "name": "julia",
   "version": "1.8.5"
  },
  "kernelspec": {
   "name": "julia-1.8",
   "display_name": "Julia 1.8.5",
   "language": "julia"
  }
 },
 "nbformat": 4
}