{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "\n# Band Structures of Bulk Structures\n\nHere, we calculate the band structure of relaxed bulk crystal structures.\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "import matplotlib.pyplot as plt\nfrom gpaw import GPAW, PW\n\nfrom ase.build import bulk\nfrom ase.dft.dos import DOS" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Setting up bulk structures\n\nFor more details regarding the set-up and relaxation of bulk structures,\ncheck out the `bulk` tutorial.\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "atoms = bulk('Ag')" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Bulk DFT calculation\n\nFor periodic DFT calculations we should generally use a number of\nk-points which properly samples the Brillouin zone.\nMany calculators including GPAW and Aims\naccept the ``kpts`` keyword which can be a tuple such as\n``(4, 4, 4)``. In GPAW, the planewave mode\nis very well suited for smaller periodic systems.\nUsing the planewave mode, we should also set a planewave cutoff (in eV):\n\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "calc = GPAW(\n mode=PW(350), kpts=[8, 8, 8], txt='gpaw.bulk_Ag.txt', setups={'Ag': '11'}\n)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Here we have used the ``setups`` keyword to specify that we want the\n11-electron PAW dataset instead of the default which has 17 electrons,\nmaking the calculation faster.\n\n(In principle, we should be sure to converge both kpoint sampling\nand planewave cutoff -- I.e., write a loop and try different samplings\nso we know both are good enough to accurately describe the quantity\nwe want.)\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "atoms.calc = calc\nprint(\n 'Bulk {0} potential energy = {1:.3f}eV'.format(\n atoms.get_chemical_formula(), atoms.get_potential_energy()\n )\n)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We can save the ground-state into a file\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "ground_state_file = 'bulk_Ag_groundstate.gpw'\ncalc.write(ground_state_file)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Density of states\n\nHaving saved the ground-state, we can reload it for ASE to extract\nthe density of states:\n\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "calc = GPAW(ground_state_file)\ndos = DOS(calc, npts=800, width=0)\nenergies = dos.get_energies()\nweights = dos.get_dos()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Calling the DOS class with ``width=0`` means ASE calculates the DOS using\nthe linear tetrahedron interpolation method, which takes time but gives a\nnicer representation. If the width is nonzero (e.g., 0.1 (eV)) ASE uses a\nsimple Gaussian smearing with that width, but we would need more k-points\nto get a plot of the same quality.\n\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "fig, ax = plt.subplots()\nax.axvline(0.0, linestyle='--', color='black', alpha=0.5)\nax.plot(energies, weights)\nax.set_xlabel('Energy - Fermi Energy (eV)')\nax.set_ylabel('Density of States (1/eV)')\nfig.tight_layout()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Time for analysis: Which parts of the spectrum do you think originate\n(mostly) from s electrons? And which parts (mostly) from d electrons?\n\nAs we probably know, the d-orbitals in a transition metal atom are\nlocalized close to the nucleus while the s-electron is much more\ndelocalized.\n\nIn bulk systems, the s-states overlap a lot and therefore split into a\nvery broad band over a wide energy range. d-states overlap much less\nand therefore also split less: They form a narrow band with a\nvery high DOS. Very high indeed because there are 10 times as\nmany d electrons as there are s electrons.\n\nSo to answer the question, the d-band accounts for most of the states\nforming the big, narrow chunk between -6.2 eV to -2.6 eV. Anything outside\nthat interval is due to the much broader s band.\n\nThe DOS above the Fermi level may not be correct, since the SCF\nconvergence criterion (in this calculation)\nonly tracks the convergenece of occupied states.\nHence, the energies over the Fermi level 0 are probably wrong.\n\n\nWhat characterizes the noble metals Cu, Ag, and Au, is that the d-band\nis fully occupied. I.e.: The whole d-band lies below the Fermi level\n(energy=0).\nIf we had calculated any other transition metal, the Fermi level would\nlie somewhere within the d-band.\n\n

Note

We could calculate the s, p, and d-projected DOS to see more\n conclusively which states have what character.\n In that case we should look up the GPAW documentation, or other\n calculator-specific documentation. So let's not do that now.

\n\n\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Band structure\n\nLet's calculate the band structure of silver.\n\nFirst we need to set up a band path. Our favourite image search\nengine can show us some reference graphs. We might find band\nstructures from both Exciting and GPAW with Brillouin-zone path:\n\n$\\mathrm{W L \\Gamma X W K}$.\n\nLuckily ASE knows these letters and can also help us\nvisualize the reciprocal cell:\n\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "lat = atoms.cell.get_bravais_lattice()\nprint(lat.description())" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "lat.plot_bz(show=True)\nplt.show()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "In general, the :mod:`ase.lattice` module provides\n:class:`~ase.lattice.BravaisLattice` classes used to represent each\nof the 14 + 5 Bravais lattices in 3D and 2D, respectively.\nThese classes know about the high-symmetry k-points\nand standard Brillouin-zone paths\n(using the [AFlow](http://aflowlib.org/) conventions).\nYou can visualize the band path object in bash via the command:\n\n $ ase reciprocal path.json\n\nYou can ``print()`` the band path object to see some basic information\nabout it, or use its :meth:`~ase.dft.kpoints.BandPath.write` method\nto save the band path to a json file such as :file:`path.json`.\n\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "path = atoms.cell.bandpath('WLGXWK', density=10)\npath.write('path.json')\nprint(path)" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "path.plot()\nplt.show()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Once we are sure we have a good path with a reasonable number of k-points,\nwe can run the band structure calculation.\nHow to trigger a band structure calculation depends\non which calculator we are using, so we would typically consult\nthe documentation for that calculator (ASE will one day provide\nshortcuts to make this easier with common calculators):\n\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "calc = GPAW(ground_state_file)\ncalc = calc.fixed_density(kpts=path, symmetry='off')" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We have here told GPAW to use our bandpath for k-points, not to\nperform symmetry-reduction of the k-points, and to fix the electron\ndensity.\n\nThen we trigger a new calculation, which will be non-selfconsistent,\nand extract and save the band structure:\n\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "bs = calc.band_structure()\nbs.write('bs.json')\nprint(bs)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We can plot the band structure in python, or in the terminal\nfrom a file using:\n\n $ ase band-structure bs.json\n\nThe plot will show the Fermi level as a dotted line\n(but does not define it as zero like the DOS plot before).\nLooking at the band structure, we see the complex tangle of what must\nbe mostly d-states from before, as well as the few states with lower energy\n(at the $\\Gamma$ point) and higher energy (crossing the Fermi level)\nattributed to s.\n\n\n" ] }, { "cell_type": "code", "execution_count": null, "metadata": { "collapsed": false }, "outputs": [], "source": [ "ax = bs.plot()\nax.set_ylim(-2.0, 30.0)\nplt.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.13.7" } }, "nbformat": 4, "nbformat_minor": 0 }