6. 2D materials - lwoncheol/Kioupakis_Group_Wiki GitHub Wiki
ℹ️ These instructions cover how to run basic DFT calculations for 2D materials
- For 2D material calculations, we need to include enough empty space (enough vacuum thickness) between each periodic image
- Additionally, if a material is polar along z-axis, so that it exhibits finite dipole along z-axis, we need to correct this dipole to cancel out artifacts of the periodic boundary conditions.
For this example, we use monolayer hexagonal boron nitride (mBN). mBN is polar along in-plane direction, but doesn’t show polarity along z-axis.
- First, we run scf calculation
&control
prefix = 'BN'
calculation = 'scf'
restart_mode = 'from_scratch'
wf_collect = .true.
tstress = .true.
tprnfor = .true.
outdir = './'
wfcdir = './'
pseudo_dir = './'
verbosity = 'high'
/
&system
ibrav = 0
celldm(1) = 4.731873875
nat = 2
ntyp = 2
nbnd = 16
ecutwfc = 70.0
/
&electrons
electron_maxstep = 500
conv_thr = 1.0d-10
mixing_mode = 'plain'
mixing_beta = 0.7
mixing_ndim = 8
diagonalization = 'david'
diago_david_ndim = 4
diago_full_acc = .true.
/
CELL_PARAMETERS
1.0000000000 0.0000000000 0.0000000000
-0.5000000000 0.8660254038 0.0000000000
0.0000000000 0.0000000000 8.0000000000
&ions
/
&cell
/
ATOMIC_SPECIES
B 10.811 b.cpi.UPF
N 14.007 n.cpi.UPF
ATOMIC_POSITIONS (crystal)
B 0.33333333333333 0.66666666666667 0.500000000
N 0.66666666666667 0.33333333333333 0.500000000
K_POINTS automatic
24 24 1 1 1 0- Note that we introduced thick vacuum region (about 20 angstroms) to prevent the interaction between periodic supercell images
- Second, we run post-processing calculation
&INPUTPP
prefix='BN',
plot_num=0,
filplot='Charge'
/
&INPUTPP
/
&plot
nfile=1
filepp(1)='Charge'
iflag=3,
fileout="Charge.cube"
output_format=6
/- Next, we run
pp.x < pp_charge > pp.outto calculate the charge distribution.
- Third, we prepare a script (
surface.inp) to analyze the charge density. In most cases, you don’t have to modify this script, so just copy and paste the following code.
header
inputfilename Charge.cube
inputfileformat cube
outputfilename Charge.pov
outputfileformat povray
isovalue 0.99
sign positive
power 1
algorithm cube
smooth T
box F
basis F
uc F
uco
0 0 0
ucv
0.0 0.0 0.0
0.0 0.0 0.0
0.0 0.0 0.0
ucu latvec
sc T
sco
0 0 0
scv
1.0 0.0 0.0
0.0 1.0 0.0
0.0 0.0 1.0
scu latvec
sct T- Finally, we run
./03-calculate_lattice.shto calculate the optimal c-axis lattice constant so that the half-size of the supercell contains 99% of the total charge density.
#!/bin/bash -l
module load berkeleygw
surface.x surface.inp &> lattice.out
# get x,y,z lengths from lattice.out
export dx2=`grep dx lattice.out | awk '{print $9;}'`
export dy2=`grep dy lattice.out | awk '{print $9;}'`
export dz2=`grep dz lattice.out | awk '{print $9*2;}'`
echo "How large does the c-axis have to be so that c/2 contains 99% of the charge density?"
echo $dz2 "a.u."For this example, we use monolayer BiTeI, which contains finite dipole along z-axis, due to the lack of inversion symmetry along z-axis. In this case, the artificial dipole correction in the vacuum region is necessary to cancel out artifacts of the periodic boundary conditions.
&control
prefix = 'BiTeI'
calculation = 'scf'
restart_mode = 'from_scratch'
wf_collect = .true.
tstress = .true.
tprnfor = .true.
outdir = './'
wfcdir = './'
pseudo_dir = './'
verbosity = 'high'
tefield = .true.
dipfield = .true.
/
&system
ibrav = 0
celldm(1) = 8.1998990096
nat = 3
ntyp = 3
nbnd = 12
ecutwfc = 50.0
edir = 3
emaxpos = 0.96667
eopreg = 0.01
eamp = 0.000
/
&electrons
electron_maxstep = 500
conv_thr = 1.0d-10
mixing_mode = 'plain'
mixing_beta = 0.7
mixing_ndim = 8
diagonalization = 'david'
diago_david_ndim = 4
diago_full_acc = .true.
/
&ions
/
&cell
/
ATOMIC_SPECIES
Bi 208.98 Bi.cpi.UPF
Te 127.6 Te.cpi.UPF
I 126.90447 I.cpi.UPF
CELL_PARAMETERS (alat= 8.1998990096)
1.0000000000 0.0000000000 0.0000000000
-0.5000000000 0.8660254038 0.0000000000
0.0000000000 0.0000000000 5.0000000000
ATOMIC_POSITIONS (crystal)
Bi 0.000000000 0.000000000 0.500000000
Te 0.666666667 0.333333333 0.420920294
I 0.333333333 0.666666667 0.597040655
K_POINTS automatic
6 6 3 1 1 1- Note that here we included 4 additional flags for the dipole correction
- edir = 3
- emaxpos = 0.96667
- eopreg = 0.01
- eamp = 0.000
- More details about the dipole correction scheme can be found at https://doi.org/10.1021/acs.nanolett.7b03003