This example plots the band structure of graphene, a 2D material. 2D band structures are not supported natively (yet), so we manually build a custom path in reciprocal space.
using DFTK
using Unitful
using UnitfulAtomic
using LinearAlgebra
# Define the convergence parameters (these should be increased in production)
L = 20 # height of the simulation box
kgrid = [6, 6, 1]
Ecut = 15
temperature = 1e-3
# Define the geometry and pseudopotential
a = 4.66 # lattice constant
a1 = a*[1/2,-sqrt(3)/2, 0]
a2 = a*[1/2, sqrt(3)/2, 0]
a3 = L*[0 , 0 , 1]
lattice = [a1 a2 a3]
C1 = [1/3,-1/3,0.0] # in reduced coordinates
C2 = -C1
positions = [C1, C2]
C = ElementPsp(:C, psp=load_psp("hgh/pbe/c-q4"))
atoms = [C, C]
# Run SCF
model = model_PBE(lattice, atoms, positions; temperature)
basis = PlaneWaveBasis(model; Ecut, kgrid)
scfres = self_consistent_field(basis)
# Construct 2D path through Brillouin zone
sgnum = 13 # Graphene space group number
kpath = irrfbz_path(model; dim=2, sgnum)
plot_bandstructure(scfres, kpath; kline_density=20)
n Energy log10(ΔE) log10(Δρ) Diag Δtime --- --------------- --------- --------- ---- ------ 1 -11.15665868980 -0.60 6.0 2 -11.16018340015 -2.45 -1.30 1.0 162ms 3 -11.16039703433 -3.67 -2.33 2.1 209ms 4 -11.16041671291 -4.71 -3.27 2.6 202ms 5 -11.16041704730 -6.48 -3.44 3.1 240ms 6 -11.16041704968 -8.62 -3.60 1.0 149ms 7 -11.16041705061 -9.03 -3.85 1.3 157ms 8 -11.16041705110 -9.31 -4.25 1.9 174ms 9 -11.16041705132 -9.66 -4.60 2.0 185ms 10 -11.16041705139 -10.15 -4.88 2.1 191ms 11 -11.16041705144 -10.36 -5.24 2.6 202ms 12 -11.16041705145 -10.86 -5.83 3.0 223ms 13 -11.16041705145 -12.09 -6.38 2.6 205ms Computing bands along kpath: Γ -> M -> K -> Γ Diagonalising Hamiltonian kblocks: 100%|████████████████| Time: 0:00:03