用awk和gnuplot可视化VASP计算材料的能带结构（以Pt金属为例）

用awk和gnuplot可视化VASP计算材料的能带结构（以Pt金属为例）

作者：     文章来源：    更新时间：2007-9-16

VASP是固体物理程序中非常稳定，性能优秀的软件之一，被广泛的用于固体性质的计算。在这里，我们以Pt金属为例，介绍如何用VASP计算并且可视化VASP的能带计算结果。

进行能带计算首先要进行正常的自洽的计算。主要目的是得到能带计算所需要的基态的电荷密度文件。所以必须保证LCHARG=.TRUE.。可以用如下的INCAR文件和KPOINTS文件。

INCAR

Pt-bulk

PREC = High

NSW    = 50        number of steps for

IBRION =2         conjugate-gradient

ISMEAR = 1

SIGMA = 0.2

LCHARG = .TRUE.



KPOINTS

kpoints

0

Monkhorst

12     12 12

0     0     0

进行完自恰计算后，用以下的INCAR文件进行能带结构的计算。关键是要保持自洽的电荷密度不变，即ICHARG=11。

Pt-bulk

PREC = High

NWRITE = 1

Ionic Relaxation

NSW    = 0        number of steps for IOM

IBRION = -1         conjugate-gradient algorithm

DOS related values:

ISMEAR = 0

ICHARG = 11



KPOINTS文件如下，我们选取了5个高对称方向，每个方向上21个K点来构建能带图。

kpoints

21      ! 21 intersections

Line-mode

cart

0 0 0   ! gamma

1 0 0   ! X


1 0 0   ! X

1 0.5 0   ! W


1 0.5 0 ! W

0.5 0.5 0.5   ! L



0.5 0.5 0.5   ! L

0 0 0   ! gamma



0 0 0   ! gamma

1 1 0   ! K



得到的OUTCAR文件里就包含了各个指定K点所对应的能级，也就是我们常说的能带结构。

E-fermi :   7.4279     XC(G=0): -13.2605     alpha+bet :-12.7584

   add alpha+bet to get absolut eigen values



k-point   1 :       0.0000    0.0000    0.0000

band No. band energies     occupation

      1      -2.4770      2.00000

      2       3.7825      2.00000

      3       3.7825      2.00000

      4       3.7825      2.00000

      5       5.7731      2.00000

      6       5.7731      2.00000

      7      22.8395      0.00000

      8      27.5743      0.00000

      9      27.5744      0.00000



k-point   2 :       0.0000    0.0250    0.0250

band No. band energies     occupation

      1      -2.4444      2.00000

      2       3.7614      2.00000

      3       3.7990      2.00000

      4       3.7990      2.00000

      5       5.7500      2.00000

      6       5.7816      2.00000

      7      22.8879      0.00000

      8      27.4306      0.00000

      9      27.5713      0.00000

…....

中间省略部分

…….

k-point 105 :       0.5000    0.5000    1.0000

band No. band energies     occupation

      1       0.6467      2.00000

      2       1.1067      2.00000

      3       7.3042      1.61813

      4       7.8021      0.00814

      5       7.8021      0.00814

      6       9.2022      0.00000

      7      15.7878      0.00000

      8      18.8190      0.00000

      9      18.8190      0.00000

可见，OUTCAR文件里包含了我们在KPOINTS中指定的105个K点。

将VASP的OUTCAR里的能带信息可视化通常是个比较头疼的事，如果要将这些信息粘贴到Origin之类的软件里进行可视化，需要大量的数据格式处理。这里我们教大家用awk文本处理程序和gnuplot画图程序简单快捷地画VASP的能带结构图。我们在这里将不详细介绍awk或者gnuplot的语法和功能，而仅仅给出一个用它们画VASP能带图的解决方案。

首先我们将OUTCAR里的从k-point 105到k-point 105的部分存成另外一个文件，不妨称之为Band.raw。

ExtractBand.awk

#!/usr/bin/awk –f

/k-point/ {

    printf "%d ", $2

    getline

    for ( i=1; i<=9; i++ ) {

        getline

        printf "%8.3f ", $2

    }

    printf "\n"

}



然后我们运行如下的awk script。./ExtractBand.awk Band.raw > band.dat
这个awk script的功能是将OUTCAR里的输出信息进行格式化，以便于gnuplot画图
OUTCAR里的原始格式

……

k-point 98 :       0.3250    0.3250    0.6500

band No. band energies     occupation

      1       1.8846      2.00000

      2       2.8157      2.00000

      3       3.7368      2.00000

      4       5.4560      2.00000

      5       6.6622      2.00000

      6      11.3258      0.00000

      7      16.2724      0.00000

      8      16.7466      0.00000

      9      24.2900      0.00000



k-point 99 :       0.3500    0.3500    0.7000

band No. band energies     occupation

      1       1.6640      2.00000

      2       2.5342      2.00000

      3       4.1600      2.00000

      4       5.7507      2.00000

      5       6.9403      1.99944

      6      12.1407      0.00000

      7      15.0639      0.00000

      8      16.8554      0.00000

      9      23.0580      0.00000

……

经过ExtractBand.awk处理后的band.dat的格式

……

98    1.885    2.816    3.737    5.456    6.662   11.326   16.272   16.747   24.290

99    1.664    2.534    4.160    5.751    6.940   12.141   15.064   16.855   23.058

……



然后我们运行以下的gnuplot脚本文件，来得到能带图。运行plotBand.sh，就得到band.eps，如图所示。

plotBand.sh

#!/bin/bash

x1=21

x2=41

x3=61

x4=81

x5=101

ymin=-10

ymax=20

ylabel=-10.6

ef_v=7.3089   #pls use E-fermi in the self consistent calculation

gnuplot << EOF

set term post landscape enhanced color dashed defaultplex "Helvetica" 14

set output 'band.ps'

set key outside

set title "Pt Band Structure (VASP)"

set ylabel "Energy (eV)"

set noxtics

set xrange [0:$x5]

set yrange [$ymin : $ymax]

set label "{/Symbol G}" at 0, $ylabel

set label "{/Symbol D}" at 0.5*$x1, $ylabel

set label "X" at $x1, $ylabel

set label "Z" at 0.5*($x1+$x2), $ylabel

set label "W" at $x2, $ylabel

set label "Q" at 0.5*($x2+$x3), $ylabel

set label "L" at $x3, $ylabel

set label "{/Symbol L}" at 0.5*($x3+$x4), $ylabel

set label "{/Symbol G}" at $x4, $ylabel

set label "{/Symbol S}" at $x4+0.4*($x5-$x4), $ylabel

set label "{/Symbol K}" at $x4+($x5-$x4)/2**0.5, $ylabel

set label "X" at $x5, $ylabel

set label "E_f" at -5.0, 0

set arrow from 0, 0 to $x5, 0 nohead

set arrow from $x1, $ymin to $x1, $ymax nohead

set arrow from $x2, $ymin to $x2, $ymax nohead

set arrow from $x3, $ymin to $x3, $ymax nohead

set arrow from $x4, $ymin to $x4, $ymax nohead

set arrow from $x4+1+($x5-$x4)/2**0.5,$ymin to $x4+1+($x5-$x4)/2**0.5, $ymax nohead

plot \

"band.dat" u 1:(\$2-$ef_v) w lines lt 1 t "", "band.dat" u 1:(\$3-$ef_v) w lines lt 1 t "",\

"band.dat" u 1:(\$4-$ef_v) w lines lt 1 t "", "band.dat" u 1:(\$5-$ef_v) w lines lt 1 t "",\

"band.dat" u 1:(\$6-$ef_v) w lines lt 1 t "", "band.dat" u 1:(\$7-$ef_v) w lines lt 1 t "",\

"band.dat" u 1:(\$8-$ef_v) w lines lt 1 t "", "band.dat" u 1:(\$9-$ef_v) w lines lt 1 t ""

EOF

最终得到的Pt的能带图如下