In this tutorial, the converge the Fermi level with respect to the number of k-points is discussed.
First of all, cell relaxation is performed by using "vc-relax" using the following input file
&control calculation = 'vc-relax' restart_mode = 'from_scratch', pseudo_dir = '/home/ikutaro/QE/pseudo', outdir = './tmp' prefix = 'gr' tstress = .true. tprnfor = .true. / &system ibrav = 4 A = 2.46 C = 10.00 nat = 2 ntyp = 1 ecutwfc = 40.0 ecutrho = 400.0 occupations = 'smearing' smearing = 'mp' degauss = 0.02 / &electrons diagonalization = 'david' conv_thr = 1.0e-12 mixing_beta = 0.5 / &ions / &cell press = 0.0d0 press_conv_thr = 0.5d0 cell_dofree = 2Dxy / ATOMIC_SPECIES C 12.01070 c_pbe_v1.2.uspp.F.UPF ATOMIC_POSITIONS (crystal) C 0.0000000000000000 0.0000000000000000 0.0000000000000000 C 0.3333333333333333 0.6666666666666667 0.0000000000000000 K_POINTS (automatic) 12 12 01 0 0 0
In this example, GBRV ultrasoft pseudopotential, cutoff energies of 40 and 400 Ry, and nonshifted 12 x 12 uniform k-point mesh is used. The vacuum of 10 Angstrom is inserted to separate periodic images. For the cell relaxation in two dimension, 'cell_dofree=2Dxy' is used.
Having relaxed the cell parameters, let us start the band structure calculation. In this tutorial, band structure calculations are performed by varying the k-point mesh in the SCF calculation to learn how the number of k-points or k-point density is crucial to determine the position of Fermi level and thereby the location of the Dirac point at the K point.
&control calculation = 'scf' restart_mode = 'from_scratch', pseudo_dir = '/home/ikutaro/QE/pseudo', outdir = './tmp' prefix = 'gr' tstress = .true. tprnfor = .true. / &system ibrav = 4 A = 2.46575870 C = 10.00 nat = 2 ntyp = 1 ecutwfc = 40.0 ecutrho = 400.0 occupations = 'smearing' smearing = 'mp' degauss = 0.02 nbnd = 8 / &electrons diagonalization = 'david' conv_thr = 1.0e-12 mixing_beta = 0.5 / &ions / &cell press = 0.0d0 press_conv_thr = 0.5d0 cell_dofree = 2Dxy / ATOMIC_SPECIES C 12.01070 c_pbe_v1.2.uspp.F.UPF ATOMIC_POSITIONS (crystal) C 0.0000000000000000 0.0000000000000000 0.0000000000000000 C 0.3333333333333333 0.6666666666666667 0.0000000000000000 K_POINTS (automatic) 12 12 01 0 0 0
&control calculation = 'bands' restart_mode = 'from_scratch', pseudo_dir = '/home/ikutaro/QE/pseudo', outdir = './tmp' prefix = 'gr' tstress = .true. tprnfor = .true. / &system ibrav = 4 A = 2.46575870 C = 10.00 nat = 2 ntyp = 1 ecutwfc = 40.0 ecutrho = 400.0 occupations = 'smearing' smearing = 'mp' degauss = 0.02 nbnd = 8 / &electrons diagonalization = 'david' conv_thr = 1.0e-12 mixing_beta = 0.5 / &ions / &cell press = 0.0d0 press_conv_thr = 0.5d0 cell_dofree = 2Dxy / ATOMIC_SPECIES C 12.01070 c_pbe_v1.2.uspp.F.UPF ATOMIC_POSITIONS (crystal) C 0.0000000000000000 0.0000000000000000 0.0000000000000000 C 0.3333333333333333 0.6666666666666667 0.0000000000000000 K_POINTS (crystal_b) 4 0.00000000 0.00000000 0.00000000 20.0 0.66666667 -0.33333333 0.00000000 20.0 0.50000000 0.00000000 0.00000000 20.0 0.00000000 0.00000000 0.00000000 0.0In this example, the energy band is draw along the Gamma-K-M-Gamma line. The symmetry points are given in the unit of basic reciprocal lattice vectors.
&bands outdir = './tmp' prefix = 'gr' filband = 'gr.bands' lsym = .true /Note that the same number of processors as the previous band structure calculation should be used. With this input, "gr.bands.gnu" and "gr.bands.rap" are generated. For use with gnuplot, "gr.bands.gnu" can be used. One can use the following script ("band_plot.sh") to generate the band plot with gnuplot
#!/bin/sh ef=`grep Fermi gr.scf.out | awk '{print $5}'` echo $ef gnuplot <<EOF ef=${ef} xmin=0.0 xmax=1.5784 ymin=-21.0 ymax=4.0 g1=0.0 k=0.6671 m=1.0007 g2=1.5784 offset=0.5 set arrow from k,ymin to k,ymax nohead lt 0 set arrow from m,ymin to m,ymax nohead lt 0 unset xtics set ylabel 'E - E_F (eV)' set label '{/Symbol G}' at g1,ymin-offset center set label 'K' at k,ymin-offset center set label 'M' at m,ymin-offset center set label '{/Symbol G}' at g2,ymin-offset center set xrange [xmin:xmax] set yrange [ymin:ymax] set xzeroaxis lt 0 plot 'gr.bands.gnu' using (\$1):(\$2-ef) pause -1 EOSand type
sh band_plot.shHere the output file for the corresponding SCF calculation is assumed to be "gr.scf.out." If you want to save the band structure change the terminal in gnuplot as
#!/bin/sh ef=`grep Fermi gr.scf.out | awk '{print $5}'` echo $ef gnuplot <<EOF ef=${ef} xmin=0.0 xmax=1.5784 ymin=-21.0 ymax=4.0 g1=0.0 k=0.6671 m=1.0007 g2=1.5784 offset=0.5 set arrow from k,ymin to k,ymax nohead lt 0 set arrow from m,ymin to m,ymax nohead lt 0 unset xtics set ylabel 'E - E_F (eV)' set label '{/Symbol G}' at g1,ymin-offset center set label 'K' at k,ymin-offset center set label 'M' at m,ymin-offset center set label '{/Symbol G}' at g2,ymin-offset center set xrange [xmin:xmax] set yrange [ymin:ymax] set xzeroaxis lt 0 set terminal postscript eps set output 'gr.bands.eps' plot 'gr.bands.gnu' using (\$1):(\$2-ef) pause -1 EOSThe figure for the band structure is written to 'gr.bands.eps.'
The position of the Dirac cone at the K-point should be located at the Fermi level. However, the position of the Fermi level is very sensitive to the number of k-point. Examine how many k-points (how dense the k-point mesh) starting from, for e.g., 6x6x1, 8x8x1, ... 12x12x1.