PBE and HSE06 band structures of GaAs

Version: 2016.4

In this tutorial you will learn the basics of how to use the FHI-aims calculator in ATK. FHI-aims is an all-electron electronic structure simulations package, so it does not employ any pseudopotentials. Moreover, it provides access to a wide range of beyond-GGA methods, including hybrid density functionals, the random-phase approximation for total energies, and the GW approximation for excited states. See the FHI-aims website for more details.



We shall consider the band structure of GaAs using the PBE and HSE06 density functionals, the latter of which employs a screened variant of exact exchange. The HSE06 can be computationally intense, but is often more accurate than standard GGA methods. See also Refs. [HSE03], [PMH+06], and [KVIS06].


The video below considers silicon as an example semiconductor, while the written below it gives instructions using GaAs instead. Do not let this confuse you; the difference in workflow is very small.

ATK 2017: The video was made before the default ATK datafile format was changed to HDF5. Therefore, ATK 2017 users should just use the ”.hdf5” extension wherever the tutorial text says ”.nc”.

Many more video tutorials are available on the QuantumWise YouTube channel.


Your ATK license file must include the ATKFHIaims feature in order to use FHI-aims with ATK. Please contact us for any inquiries.


This tutorial uses rather low accuracy settings in order to keep the computational time to a minimum, for instructional purposes. Please remember that settings related to k-point grids, basis sets, and integration grids should always be carefully checked for convergence in publication-quality scientific work.

PBE band structure

  • Open the builder_icon Builder and set up the GaAs crystal (zincblende structure).
  • Send the structure to the script_generator_icon Scripter and add the calculator_icon New Calculator and analysis_icon Analysis ‣ Bandstructure script blocks.
  • Open the New Calculator block and select the FHI-aims calculator.
  • Select the PBE density functional and set up a 12x12x12 k-point grid.
  • In Control Parameters, use the system default settings for “Insulator”.

Save the script as pbe.py (but do not close the Script Generator window) and run the script using the job_manager_icon Job Manager. Alternatively, you can run the FHI-aims calculation from the command line:

atkpython pbe.py > pbe.log

The calculation should take only a few minutes, and the GaAs band structure will pop up on the VNL LabFloor. You can visualize it using the Bandstructure Analyzer and zoom in around the Fermi level:


HSE06 band structure

You will now compute the same band structure using HSE06 and compare it to the PBE results. HSE06 can be computationally heavy for periodic systems, so we will here use a coarser 6x6x6 k-point grid.


In the following, the PBE and HSE06 calculations are conveniently done after each other in a single Python script. This illustrates one of the advantages of running FHI-aims calculations through ATK Python scripts: Composite workflows are easily created using the flexibility of Python scripting, see more at Python in ATK.

  • PBE calculation:
    • Go back to the script_generator_icon Script Generator window used above. Set the k-point grid to 6x6x6, and change the band structure route to L, G, X.
  • HSE calculation:
    • Add one more pair of calculator_icon New Calculator and analysis_icon Bandstructure blocks.
    • Adjust the settings to make them identical to those for PBE, but chooce HSE06 as xc functional.

Finally, change the default output file name to hse.nc, and save the script as hse.py.


Run the calculations, preferably in parallel using e.g. 4 MPI processes:


You can of course obtain the same in a terminal:

QWDIR/libexec/mpiexec.hydra -n 4 atkpython hse.py > hse.log

where QWDIR is the VNL-ATK installation directory.

The calculations should take roughly 30 minutes in total if parallelized over 4 CPU cores on a modern laptop.

When finished, you can mark both band structure objects on the LabFloor and use the Compare Data plugin to visualize both band structures in one plot:

../../_images/pic_31.png ../../_images/pic_30.png

The two band structures are clearly different; the HSE06 conduction bands (green, object ID gID003) are pushed up relative to the PBE ones, leading to a predicted band gap of 1.36 eV, not too far from the experimental value of 1.52 eV.


Please be aware that the results above were obtained using Medium accuracy (“Light” integration grids and “Tier 1” basis sets) and a somewhat coarse k-point grid. Calculations with Tight settings, more basis functions, and more k-points may give different results. See more below.

Also note that band structure calculations should usually be preceeded by a first principles geometry optimization of the structure.

Basis sets and integration grids


Fig. 103 GaAs band gap calculated with different FHI-aims settings.

As seen in the figure above, FHI-aims results do in general depend on the computational settings, for example basis sets and integration grids. “Tight” means more refined integrations grids than “Light”, while “Tier 2” has more basis functions than “Tier 1”. Note that increased computational quality requires in general also more computing resources (CPU hours and/or memory).

The calculations were performed using the script GaAs.py, which is reproduced below. The first 18 script lines (highlighted) set up a Python loop over different combinations of basis sets and levels of integration accuracy. SCF and band structure calculations for all 6 combinations are hereby executed in a single ATK job.

# ====================================================
basis_sets  = [Light, Tight]
tiers       = [0, 1, 2]
# ====================================================

for basis in basis_sets:
    for tier in tiers:
        out = 'GaAs_pbe_%s_tier-%i.nc' % (basis.speciesDefaultAsString(), tier)

        # Basis Set
        basis_set = FHIaimsBasisSet()
        species_bases = [
            FHIaimsBasisSetSpecies(Arsenic, basis, tier),
            FHIaimsBasisSetSpecies(Gallium, basis, tier),

        # -------------------------------------------------------------
        # Bulk Configuration
        # -------------------------------------------------------------
        # Set up lattice
        lattice = FaceCenteredCubic(5.6537*Angstrom)
        # Define elements
        elements = [Gallium, Arsenic]
        # Define coordinates
        fractional_coordinates = [[ 0.  ,  0.  ,  0.  ],
                                  [ 0.25,  0.25,  0.25]]
        # Set up configuration
        bulk_configuration = BulkConfiguration(

        # -------------------------------------------------------------
        # Calculator
        # -------------------------------------------------------------
        k_point_sampling = MonkhorstPackGrid(9,9,9)
        numerical_accuracy_parameters = NumericalAccuracyParameters(
        iteration_control_parameters = IterationControlParameters(

        # FHI-aims Parameters
        fhiaims_parameters = FHIaimsParameters()
        fhiaims_parameters.addControlKeywordValuePair("occupation_type","gaussian 0.000001")

        # Working Directory
        working_directory = '.'

        # FHI-aims calculator
        calculator = FHIaimsCalculator(


        # -------------------------------------------------------------
        # Bandstructure
        # -------------------------------------------------------------
        bandstructure = Bandstructure(
        nlsave(out, bandstructure)

The script plot.py may then be used to read the band structures, extract the band gap, plot the results, and save the plot as a PNG file:

basis_sets  = [Light, Tight]
tiers       = [0, 1, 2]
gaps   = []
labels = []
for tier in tiers:
    for basis in basis_sets:
        out = 'GaAs_pbe_%s_tier-%i.nc' % (basis.speciesDefaultAsString(), tier)
        bandstructure = nlread(out, Bandstructure)[0]
        gap = bandstructure._directBandGap().inUnitsOf(eV)
        label = "%s\n" % basis.speciesDefaultAsString().capitalize()
        label += "Tier %i" % tier

import pylab
x = range(len(labels))
pylab.plot(x, gaps, 'bo-')
pylab.ylabel('GaAs direct band gap (eV)')
pylab.title('XC=PBE using 9x9x9 k-points')
ax = pylab.gca()

You can try to run the scripts yourself – the calculations should take ca. 30 mins in total if run in serial on a modern laptop.

Germanium band structure

For completeness, you can also calculate the band structure of Ge using HSE06: Download the script germanium.py and run the calculation. It takes roughly 25 minutes when parallelized on 16 cores, but can of course be run on fewer cores. The result is illustrated below.



[HSE03]Jochen Heyd, Gustavo E. Scuseria, and Matthias Ernzerhof. Hybrid functionals based on a screened coulomb potential. J. Chem. Phys., 118(18):8207–8215, 2003. doi:10.1063/1.1564060.
[KVIS06]Aliaksandr V. Krukau, Oleg A. Vydrov, Artur F. Izmaylov, and Gustavo E. Scuseria. Influence of the exchange screening parameter on the performance of screened hybrid functionals. J. Chem. Phys., 2006. doi:10.1063/1.2404663.
[PMH+06]J. Paier, M. Marsman, K. Hummer, G. Kresse, I. C. Gerber, and J. G. Angyan. Screened hybrid density functionals applied to solids. J. Chem. Phys., 2006. doi:10.1063/1.2187006.