*Crystalline Al [#jb02f7c9] 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 [#w76c80cd] - 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 [#rb127cad] 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 $ gnuplot band.gp to view the band structure.