Crystalline Al †
In this example, how to determine the equilibrium lattice constant is described.
Since the stress is not yet implemented in the official version of STATE, we need to optimize the lattice parameters manually. To do so, we calculate the total energy as a function of lattice constant and obtain the one which gives the lowest total energy
Cell optimization †
- Job script (qsub_seq.sh)
Here, instead of creating input files one by one, let us use a simple job script, which perform the job (semi-)automatically
#$ -S /bin/sh
#$ -cwd
#$ -pe fillup 6
#$ -N Al
#disable OPENMP parallelism
export OMP_NUM_THREADS=1
# execuable of the STATE code
ln -fs ${HOME}/STATE/src/STATE .
# pseudopotential data
ln -fs ${HOME}/STATE/gncpp/pot.Al_pbe1 fort.37
# launch STATE
for a in $(seq 7.50 0.10 7.90)
do
INPUT_FILE=nfinp_a${a}
OUTPUT_FILE=nfout_a${a}
cat > ${INPUT_FILE} << EOF
0 0 0 0 0 0 : I_CTRL(1:6) (DUMMY)
4.00 8.00 1 1 1 : GMAX GMAXP NTYP NATM NATM2
221 2 : NUM_SPACE_GROUP TYPE_BRAVIS_LATTICE
${a} ${a} ${a} 90.00 90.00 90.00 : A B C ALPHA BETA GAMMA
12 12 12 1 1 1 : K_MESH
0 0 : NCORD NINV : IWEI IMDTYP ITYP
0.00 0.00 0.00 1 0 1 : POS(1:3) IWEI IMDTYP ITYP
13 0.50 26.98 6 1 0.20 : IATOMN ALFA AMION ILOC IVAN
0 0 0 0 0 : ICOND INIPOS INIVEL ININOS INIACC
0 1 : IPRE IPRI
30 30 0 84200.00 0 : NMD1 NMD2 LAST_ITER CPUMAXjIFSTOP
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 : DESTM
101 : NBZTYP
4 4 4 : DUMMY
4 4 4 : DUMMY
6 : KEG
1 : NEXTST
0 : DUMMY
2 : IMSD
0 : EVALUATE_EKO_DIFF
0 : NPDOSAO
0 0.00 : SM_N DOPPING
EOF
mpirun -np $NSLOTS ./STATE < ${INPUT_FILE} > ${OUTPUT_FILE}
done
After the calculations, let us use "getetot.sh" placed in the working directory to get
7.5 -2.07149508
7.6 -2.07210987
7.7 -2.07223992
7.8 -2.07193773
7.9 -2.07125357
Use gnuplot or whatever software to fit the total energy to any function.
Band strucure calculation †
To perform the band structure calculation, we need to a SCF calculation, followed by a non-SCF calculation along the high symmetry line.
- Input file for SCF (nfinp_dav)
0 0 0 0 0 0 : I_CTRL(1:6) (DUMMY)
4.00 8.00 1 1 1 : GMAX GMAXP NTYP NATM NATM2
221 2 : NUM_SPACE_GROUP TYPE
7.60 7.60 7.60 90.00 90.00 90.00 : A B C ALPHA BETA GAMMA
06 06 06 1 1 1 : N1 N2 N3 M1 M2 M3
0 0 : NCORD, NINV
0.00 0.00 0.00 1 0 1 : CPS(1,1:3) IWEI IMDTYP ITYP
13 0.50 26.98 6 1 0.2 : IATOMN ALFA AMION ILOC IVAN
0 0 0 0 0 : ICOND INIPOS INIVEL ININOS INIACC
0 1 : IPRE IPRI
30 30 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.002 0 0.50D+03 0 : WIDTH FORCCR ISTRESS
ggapbe 1 : XCTYPE KSPIN
2.00 : DESTM
101 : NBZTYPE
4 4 4 : NKX NKY NKZ (DUMMY)
4 4 4 : NKX2 NKY2 NKZ2 (DUMMY)
6 : KEG
1 : NEXTST
0 : (DUMMY)
2 : IMSD
0 : EVALUATE_EKO_DIFF
0 : NPDOSAO
0 0.000 : SM_N DOPPING (DUMMY)
- Input file for the band structure calculation (nfinp_band)
0 0 0 0 0 0 : I_CTRL(1:6) (DUMMY)
4.00 8.00 1 1 1 : GMAX GMAXP NTYP NATM NATM2
221 2 : NUM_SPACE_GROUP TYPE
7.60 7.60 7.60 90.00 90.00 90.00 : A B C ALPHA BETA GAMMA
06 06 06 1 1 1 : N1 N2 N3 M1 M2 M3
0 0 : NCORD, NINV
0.00 0.00 0.00 1 0 1 : CPS(1,1:3) IWEI IMDTYP ITYP
13 0.50 26.98 6 1 0.2 : IATOMN ALFA AMION ILOC IVAN
22 0 0 0 0 : ICOND INIPOS INIVEL ININOS INIACC
0 1 : IPRE IPRI
30 30 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 : DESTM
101 : NBZTYPE
4 4 4 : NKX NKY NKZ (DUMMY)
4 4 4 : NKX2 NKY2 NKZ2 (DUMMY)
6 : KEG
1 : NEXTST
0 : (DUMMY)
2 : IMSD
0 : EVALUATE_EKO_DIFF
0 : NPDOSAO
0 0.000 : SM_N DOPPING (DUMMY)
&KPOINTS_BAND
NKSEG 4
KMESH 40 20 20 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
&END
In the band structure calculation, we set ICOND=22.
In addtion, high symmetry points and lines should be specified.
This can be done by adding the section &KPOINTS_BAND ... &END.
Here NKSEG is the number of k-point segment (and the number of k-points which define the k-point segments is NKSEG+1), KMESH is the k-point meshes for each segment, and NKSEG+1 k-points in the unit of the reciprocal lattice vector are specified after the keyword "KPOINTS."
When the band structure calculation converges, "energy.data" can be found in the working directory, which contains the k-points and Kohn-Sham eigenvalues.
To plot the band structure with gnuplot or xmgrace, I use "energy2band" contained in util/bandutil. By executing the "energy2band," you are prompted to input the number of bands, the number of bands to be plotted, the number of k-points, and the Fermi level.
The number of bands and k-points are extracted in the band structure calculation, while the Fermi level is taken from the previous SCF calculation.
Below, a gnuplot script to plot the band structure is shown (band.gp)
G=0.00000000
X=0.82673491
W=1.24010236
L=1.82469222
G2=2.54066566
emin=-12.25
emax=4.25
xmin=0.0
xmax=G2
set xrange [xmin:xmax]
set yrange [emin:emax]
set arrow from G,emin to G,emax nohead lt 0
set arrow from X,emin to X,emax nohead lt 0
set arrow from W,emin to W,emax nohead lt 0
set arrow from L,emin to L,emax nohead lt 0
set arrow from G2,emin to G2,emax nohead lt 0
set xzeroaxis lt 0
unset xtics
set ylabel 'Energy (eV)'
plot 'band.data' with points pt 7 ps 2
pause -1
Type
$ gnuplot band.gp
to view the band structure.
![[PukiWiki]](image/logo_state-senri3_80x80.png "[PukiWiki]"){kind=logo}
