Crystalline Si †
In this example, silicon in the diamond structure (space group number of 227) with the lattice constant of 10.30 Bohr, cutoff energy of 16 Ry (GMAX=4, GMAXP=8), and nonshifted 8 x 8 x 8 k-point, is used as an example of the SCF calculation of solid with the band gap.
SCF calculation †
- Input file (nfinp_1)
0 0 0 0 0 0 : I_CTRL(1:6) (DUMMY)
4.00 8.00 1 2 2 : GMAX, GMAXP, NTYP, NATM, NATM2
227 2 : NUM_SPACE_GROUP TYPE_BRAVIS_LATTICE
10.30 10.30 10.30 90.00 90.00 90.00 : A B C ALPHA BETA GAMMA
08 08 08 1 1 1 : N1 N2 N3 M1 M2 M3
0 0 : NCORD NINV : IWEI IMDTYP ITYP
0.00d0 0.00d0 0.00d0 1 1 1
0.25d0 0.25d0 0.25d0 1 1 1
14 0.50 28.09 6 1 0.2 : TYPE 1IATOMN,ALFA,AMION,ILOC,IVAN
0 0 0 0 0 : ICOND INIPOS INIVEL ININOS INIACC
0 1 : IPRE IPRI
20 20 0 84200.00 0 : NMD1 NMD2 LAST_ITER CPUMAX IFSTOP
6 1 : WAY_MIX MIX_WHAT
0 20 0.60 : ITER_START KBXMIX MIX_ALPHA
0.20 0.30 0.20 0.20 0.20 : DTIM1 DTIM2 DTIM3 DTIM4 DTIM
300.00 4 1 0.50D-09 : DTIO IMDALG IEXPL EDELTA
0.0002 0.50D+03 0 : WIDTH FORCCR ISTRESS
ggapbe 1 : XCTYPE KSPIN
2.00 : DESTM
102 : NBZTYP
4 4 4 : DUMMY
4 4 4 : DUMMY
8 : KEG
1 : NEXTST
0 : DUMMY
2 : IMSD
0 : EVALUATE_EKO_DIFF
2 : NPDOSAO
1
2
-15.00 5.00 0.20 501
2.400 0.20
0.20 14
0 0.00 : SM_N DOPPING
- Execution (interactive mode)
$ mpirun -np 6 ./STATE < nfinp_1 > nfout_1
- Execution (batch job)
$ qsub qsub.sh
- Output file (nfout_1)
Convergence of the total energy can be monitored by using grep as
$ grep ETOT: nfout_1
The result is:
ETOT: 1 -6.05513096 0.6055E+01 0.3203E-02
ETOT: 2 -7.84016187 0.1785E+01 0.5187E-02
ETOT: 3 -7.87270490 0.3254E-01 0.2825E-02
ETOT: 4 -7.87351715 0.8123E-03 0.6089E-03
ETOT: 5 -7.87355245 0.3530E-04 0.1887E-03
ETOT: 6 -7.87355822 0.5769E-05 0.2560E-04
ETOT: 7 -7.87355833 0.1069E-06 0.1066E-04
ETOT: 8 -7.87355833 0.4548E-08 0.1730E-05
ETOT: 9 -7.87355833 0.1916E-09 0.4194E-06
ETOT: 10 -7.87355833 0.2984E-10 0.1669E-07
ETOT: 11 -7.87355833 0.1951E-11 0.9524E-08
Converged total energy and its components
TOTAL ENERGY AND ITS COMPONENTS
TOTAL ENERGY = -7.87355833 A.U.
FREE ENERGY = -7.87355833 A.U.
KINETIC ENERGY = 3.01922457 A.U.
HARTREE ENERGY = 0.55014222 A.U.
XC ENERGY = -2.40098662 A.U.
LOCAL ENERGY = -0.84294976 A.U.
NONLOCAL ENERGY = 0.16885288 A.U.
EWALD ENERGY = -8.36784162 A.U.
PC ENERGY = 0.00000000 A.U.
ENTROPIC ENERGY = 0.00000000 A.U.
Band structure calculation †
To calculate the band structure, an SCF calculation is performed to obtain a converged charge density, and the non-SCF calculation along the symmetry lines is performed.
- Input file (nfinp_scf)
0 0 0 0 0 0 : INPUT_CTRL(1:6) (DUMMY)
6.00 12.00 1 2 2 : GMAX GMAXP NTYP NATM NATM2
227 2 : num_space_group, type
10.2347 10.2347 10.2347 90.0 90.0 90.0 : A B C ALPHA BETA GAMMA
4 4 4 2 2 2 : N1 N2 N3 M1 M2 M3
0 0 : NCORD, NINV
0.00 0.00 0.00 1 1 1 : CPS(1,1:3) IWEI IMDTYP ITYP
0.25 0.25 0.25 1 1 1 : CPS(2,1:3) IWEI IMDTYP ITYP
14 0.5000 28.09 6 1 0.0 : ATOMN ALFA AMION ILOC IVAN ZETA1
0 0 0 0 0 : ICOND INIPOS INIVEL ININOS INIACC
0 1 : IPRE IPRI
200 200 0 84200.00 0 : NMD1 NMD2 LAST_ITER CPUMAX IFSTOP
6 1 : WAY_MIX MIX_WHAT
0 20 0.60 : ITER_START KBXMIX MIX_ALPHA
0.20 0.30 0.20 0.20 0.20 : DTIM1 DTIM2 DTIM3 DTIM4 DTIM
300.00 4 1 0.50D-09 : DTIO IMDALG IEXPL EDELTA
-0.0020 0.50D+03 0 : WIDTH FORCCR ISTRESS
ggapbe 1 : xctype,kspin
2.00 : DETSTM
101 : NBZTYP
0 0 0 : NKX NKY NKZ (DUMMY)
0 0 0 : NKX2 NKY2 NKZ2 (DUMMY)
8 : KEG
1 : NEXTST
0 : (DUMMY)
2 : IMSD
0 : EVALUATE_EKO_DIFF
0 : NPDOSAO
0 0.00 : SM_N DOPING
- Input file (nfinp_band)
0 0 0 0 0 0 : INPUT_CTRL(1:6) (DUMMY)
6.00 12.00 1 2 2 : GMAX GMAXP NTYP NATM NATM2
227 2 : num_space_group, type
10.2347 10.2347 10.2347 90.0 90.0 90.0 : A B C ALPHA BETA GAMMA
4 4 4 2 2 2 : N1 N2 N3 M1 M2 M3
0 0 : NCORD, NINV
0.00 0.00 0.00 1 1 1 : CPS(1,1:3) IWEI IMDTYP ITYP
0.25 0.25 0.25 1 1 1 : CPS(2,1:3) IWEI IMDTYP ITYP
14 0.5000 28.09 6 1 0.0 : ATOMN ALFA AMION ILOC IVAN ZETA1
22 0 0 0 0 : ICOND INIPOS INIVEL ININOS INIACC
0 1 : IPRE IPRI
200 200 0 84200.00 0 : NMD1 NMD2 LAST_ITER CPUMAX IFSTOP
6 1 : WAY_MIX MIX_WHAT
0 20 0.60 : ITER_START KBXMIX MIX_ALPHA
0.20 0.30 0.20 0.20 0.20 : DTIM1 DTIM2 DTIM3 DTIM4 DTIM
300.00 4 1 0.50D-09 : DTIO IMDALG IEXPL EDELTA
-0.0100 0.50D+03 0 : WIDTH FORCCR ISTRESS
ggapbe 1 : xctype,kspin
2.00 : DETSTM
101 : NBZTYP
0 0 0 : NKX NKY NKZ (DUMMY)
0 0 0 : NKX2 NKY2 NKZ2 (DUMMY)
16 : KEG
1 : NEXTST
0 : (DUMMY)
2 : IMSD
0 : EVALUATE_EKO_DIFF
0 : NPDOSAO
0 0.00 : SM_N DOPING
&KPOINTS_BAND
NKSEG 6
KMESH 40 20 20 20 40 20
KPOINTS
0.000 0.000 0.000
0.000 0.500 0.500
0.250 0.500 0.750
0.500 0.500 0.500
0.000 0.000 0.000
0.375 0.375 0.750
0.500 0.500 1.000
&END
When the non-SCF calculation is performed, k-points and Kohn-Sham eigenvalues are printed in "energy.data." However, "energy.data" alone cannot be used to plot the band structure. For this purpose, "energy2band" included in the "util" directory is used.
Here a sample bash script, which generate the band structure is shown.
#!/bin/sh
#
# Fermi level obtained in the SCF calculation
ef=0.0790726449
# Number of k-points
nk=161
# Number of bands
nb=16
# Divisions for each segment (in this case we have 6 segments)
ik1=40
ik2=60
ik3=80
ik4=100
ik5=140
ik6=160
# Maximum energy
emin=-12.5
# Minimum energy
emax=7.5
#
nk1=`expr $ik1 + 2`
nk2=`expr $ik2 + 2`
nk3=`expr $ik3 + 2`
nk4=`expr $ik4 + 2`
nk5=`expr $ik5 + 2`
nk6=`expr $ik6 + 2`
#
e2b=~/STATE/util/bandutil/src/energy2band
#
${e2b} << EOF
${nb}
${nb}
${nk}
${ef}
EOF
#
k1=`head -${nk1} band.data | tail -1 | awk '{print $1}'`
k2=`head -${nk2} band.data | tail -1 | awk '{print $1}'`
k3=`head -${nk3} band.data | tail -1 | awk '{print $1}'`
k4=`head -${nk4} band.data | tail -1 | awk '{print $1}'`
k5=`head -${nk5} band.data | tail -1 | awk '{print $1}'`
k6=`head -${nk6} band.data | tail -1 | awk '{print $1}'`
#
xmin=0.0
xmax=${k6}
offset=-0.5
#
pos=`echo "scale=2; $emin + $offset" | bc`
#
gnuplot<<EOF
set xrange [$xmin:$xmax]
set yrange [$emin:$emax]
set arrow from ${k1},${emin} to ${k1},${emax} nohead lt 0
set arrow from ${k2},${emin} to ${k2},${emax} nohead lt 0
set arrow from ${k3},${emin} to ${k3},${emax} nohead lt 0
set arrow from ${k4},${emin} to ${k4},${emax} nohead lt 0
set arrow from ${k5},${emin} to ${k5},${emax} nohead lt 0
unset xtics
set mytics
set ylabel 'Energy (eV)'
set label '{/Symbol G}' at 0.0, ${pos}
set label 'X' at ${k1}, ${pos}
set label 'W' at ${k2}, ${pos}
set label 'L' at ${k3}, ${pos}
set label '{/Symbol G}' at ${k4}, ${pos}
set label 'K' at ${k5}, ${pos}
set label 'X' at ${k6}, ${pos}
plot 'band.data' with lines
EOF
Note that the energy origin is set to the valence band maximum, which is located at the Gamma point.
|