Band Unfolding in Siesta


Sara G Mayo (Universidad Autónoma de Madrid)


For background, see the slide presentation by Sara Mayo.


In this session we present the Siesta band unfolding algorithm, a post-processing tool that allows to visualize energy bands of non-periodic systems in their reference pure system primi tive Brillouin zone (PBZ). This method consists on two steps: a full unfolding of the bands from the supercell Brillouin zone (SBZ) to all reciprocal space and a refolding into the PBZ. Both fully-unfolded bands and refolded bands can provide useful information of the di erent materials under study.


The unfolding algorithm is a post-processing tool. It is run independently after a Siesta calculation, employing its output files. The source code can be found in Siesta folder /Util/Unfolding along with documentation and some model examples, some of which we will cover in this session.

Set fdf file

Your fdf file is common for Siesta and Unfold and must include the following blocks:

SaveHS true # makes siesta write the .HSX output file
LatticeConstant 5.430 Ang # used by siesta and unfold
BandLinesScale pi/a # used in UnfoldedBandLines blocks

%block LatticeVectors # simulation (super)cell
1.000 0.000 0.000 # (this is an 8-atom Si supercell)
0.000 1.000 0.000
0.000 0.000 1.000
%endblock LatticeVectors

%block UnfoldedBandLines
320 -20.0 20.0 eV # energy: nE, Emin, Emax, units
1 3.0 3.0 3.0 # mustBeOne, first_qVector
100 0.0 0.0 0.0 \Gamma # nq_line, endLine_q, qLabel(optional)
120 6.0 0.0 0.0
1 0.0 0.0 0.0 \Gamma # paths begin with a single point
100 4.0 4.0 0.0
%endblock UnfoldedBandLines

Up to this point, the algorithm will only fully-unfold the bands. If you recquire both fully-unfolded and refolded bands, the fdf file must include as well the refolding blocks:

# cutoff for (the square of) refolding G vectors
RefoldingGcutoff 25 Ry

%block RefoldingLatticeVectors
0.000 0.500 0.500 # primitive unit cell, units of LatticeConstant
0.500 0.000 0.500
0.500 0.500 0.000
%endblock RefoldingLatticeVectors

This example unfolds the bands of a 8-atom Si supercell into the fcc Si cell. A defect can be introduced by removing one atom of the supercell.

Run Siesta and unfold

Perform a Siesta calculation as usual, then check for the .ion, .EIG and .HSX outputs (and the .psf files) before running unfold:

unfold < sysLabel.fdf > unfold.out

You may run unfold in parallel for long calculations. Outputs are sysLabel.unfoldedBands and sysLabel.refoldedBands (if RefoldingLatticeVectors block is present). If more than one path is specified in UnfoldingBandLines, the files are suffixed with .pathN. Some Matlab plotting routines are available (see this section).

CPU time

The CPU time used by unfold is proportional to the number of q vectors specified in Unfolding BandLines, but it does not increase with the energy resolution. Since the parallelization over q vectors is almost perfect, execution wall-clock time is inversely proportional to the number of MPI nodes, provided they are less than the number of q points in UnfoldedBandLines. Bands can be refolded to any Brillouin zone, even to one not commensurate with the simulation cell, but this requires a new diagonalization for each q+G vector, increasing CPU time.

Proposed examples


Most of the provided model systems involve long Siesta calculations. You can run them yourself or use the provided output files to save this time. The links to those files are given below as appropriate.

  • Si bulk and vacancy in an 8-atom supercell

    A fast and simple example. Run Siesta unfold for a periodic system and plot the refolded bands with plotdos.m. Conventional bands are recovered. As unfold yields a (local) density of states, additional information about degeneracies appears in the plot. This will be your reference system for the defective Si crystal. Next, run unfold for the defective system and compare both. You can plot the DOS difference with plotdiff.m.

  • Amorphous Si

    Amorphous Si is a completely non-periodic system, yet its bands can be fully unfolded and refolded into the PBZ of the crystal to recover some information, such as the band gap width. Run unfold using the provided Siesta output files:

    Plot the fully unfolded bands and try to obtain the effective mass by fitting to a free-electron dispersion relation.

  • Divacancy in graphene

    We propose a defect in a graphene layer as an example of a 2D system. The divacancy breaks the symmetry and thus the six K points split in two kinds. Use the provided Siesta output files to refold the bands around one of the K0 points and observe the effect on the Dirac cone. Periodicity in the z-direction is artificially induced by the DFT calculation, therefore a full unfolding in this dimension does not make sense. We must refold in the z dimension to correct this. To obtain the fully-unfolded bands of graphene, we set small x and y refolding lattice vectors while the z ones must be equal. In this way, the x and y refolding G vectors in reciprocal space are excluded by the cutoff.

    Extra: LDOS surfaces can be visualized by plotmesh.m. To generate a LDOS surface, create the file meshgen.dat with format:

    ne emin emax
    nx x0 xend
    ny y0 yend

    then run the routine meshgen and include the output in your .fdf file. Here, nq = nx*ny.


SG Mayo, F Yndurain, JM Soler 2020, J. Phys.: Condens. Matter 32 205902.

More theoretical background and the reference plots for the proposed exercises are available in the slide presentation referenced above.


If the utility meshgen is not available as an executable, you can get it by:

cc -o meshgen meshgen.c

Guide to unfolding plotting scripts

Directory scripts contains some plotting scripts useful for this tutorial:

  • The original Matlab scripts by Jose Soler and Sara Mayo:

    • plotmesh.m

    • plotdos.m

    • plotdiff.m

    • plotUandR.m

  • A python+matplotlib version of plotUandR.m, to plot unfolded/refolded bands, called`, by Nils Wittemeier:


    If you are in the School VM, you need to do first:

    workon siesta_school

    to access the right python modules.


    ./ si8.refoldedBands --kpoints 0 20 45 75 --klabel 'L' '$\Gamma$' 'X' '$\Gamma$'

    Control of saturation and contrast can be easily done by changing the color map:

    ./ system.refoldedBands -cmap 'ColorMapName'

    A full list of color maps can be found in the matplotlib documentation here. The ‘Greys’ settings gives a reasonable contrast.

  • A version of plotdos.m for octave, by Vladimir Dikan, called plotdos-octave.m.

    Usage example:

    octave --persist plotdos-octave.m

    (One has to edit the top of the file for the proper system label).

More contributions are welcome!