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 LinearAlgebra
using Plots
using PseudoPotentialData
using Unitful
using UnitfulAtomic
# 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, PseudoFamily("cp2k.nc.sr.pbe.v0_1.semicore.gth"))
atoms = [C, C]
# Run SCF
model = model_DFT(lattice, atoms, positions; functionals=PBE(), temperature)
basis = PlaneWaveBasis(model; Ecut, kgrid)
scfres = self_consistent_field(basis)
# Construct 2D path through Brillouin zone
kpath = irrfbz_path(model; dim=2, space_group_number=13) # graphene space group number
bands = compute_bands(scfres, kpath; kline_density=20)
plot_bandstructure(bands)
n Energy log10(ΔE) log10(Δρ) Diag Δtime --- --------------- --------- --------- ---- ------ 1 -11.15656705121 -0.60 5.9 191ms 2 -11.16020721823 -2.44 -1.30 1.1 127ms 3 -11.16039992622 -3.72 -2.33 2.4 147ms 4 -11.16041622339 -4.79 -3.21 2.3 138ms 5 -11.16041703063 -6.09 -3.47 2.7 148ms 6 -11.16041704404 -7.87 -3.63 1.1 113ms 7 -11.16041704963 -8.25 -3.98 1.1 111ms 8 -11.16041705117 -8.81 -4.51 1.9 123ms 9 -11.16041705142 -9.60 -5.01 1.9 125ms 10 -11.16041705145 -10.58 -5.38 1.7 451ms 11 -11.16041705145 -11.57 -5.87 1.7 118ms 12 -11.16041705145 -12.34 -6.23 2.0 127ms