【更新1.1版本】二维材料的二维能带图绘制以及有效质量计算软件band2d.py

ggdh

听说2D材料最近要火？我博后期间做这个写了个程序。关键是能够产生2D能带图，1维材料只能产生1维能带图，2维材料可以产生2维能带和1维能带。看2维能带能够直观掌握整个布里渊区的信息。而3维材料的3维能带图估计只有用VR跑到布里渊区里面才能看个contour版的3D能带图。。。现在把这个脚本分享出来，亲生的脚本希望能产生更多价值。

需求：
vasp
python 3.8.3
mayavi 4.7.3

1.1版新特性
1.1版主要增加了一些风格的控制
1. --shift 选项位移能带，--shift homo/fermi/0.04 分别是把vbm位移到0，fermi位移到0，以及位移0.04eV
2. --style3d 选项可以设置2维能带图的一些风格，包括是否显示边框，坐标轴，填色风格，尺寸大小等等（具体用法见帮助）
3. --rangeC 选项设置contour图的色彩范围（具体用法见帮助）
4. --contour 选项设置二维contour填色图等高线的数量和位置（具体用法见帮助）
5. -j pc 配合 --range3d 可以显示某个能带的contour填色图

主要功能：
对于二维材料:
a.计算1维能带结构，并绘制如下所示的图(左边是能带图，右边是K路径和布里渊区示意图）：


b.计算2维能带结构，并绘制如下所示的图（这里只绘制了top of valence band）：



c.计算有效质量，并在能带平面图上标出方向（分别是空穴和电子的有效质量）：



安装说明:
先安装anaconda作为python环境
https://www.continuum.io/downloads
这里装python3.8.3版本的，可以装windows或者linux版本
安装好之后再linux Termial 或者是windows cmd 里面输入：
pip install mayavi
pip install PyQt5

使用说明：
1. 准备工作
用vasp优化好你的二维结构，或者是拿优化好的结构算个单点，例子输入文件如下，注意你的slab必须在xy平面，而z方向是真空层，同时晶体的ab轴在xy平面内，c轴和z平行：

INCAR文件如下：
PREC=Normal
ISMEAR=0
SIGMA=0.05
LREAL=AUTO
NSW=100
IBRION=2
ISIF=4

KPOINT文件如下：
0
Gamma
4 4 1
计算完成后产生的OUTCAR IBZKPT文件会被本程序读取。


2.检查当前情况以及设置条件
将band2d.py 拷贝到刚才vasp计算的文件夹下，运行
python band2d.py -j c
或者在linux 下把band2d.py 放到PATH路径下，直接运行
band2d.py -j c
对当前目录进行检查，并产生输出如下（红色是我的说明文字，不是输出的部分）
10 k points found in OUTCAR  #在OUTCAR中找到10个K点
no mesh_band data found #没有找到网格数据，网格数据是用来产生2D能带的
no mass_points data found #没有找到质量数据，质量数据是用来计算有效质量的
no points detected for peripheral band #没有找到外周能带数据点，外周能带就是布里渊区的边
no points detected for radical band #没有找到径向能带数据点，径向能带就是从Gamma点到布里渊区边缘的能带
0 k points found in 2D data #2维数据里面找到了0个K点，2维数据就是上面的mesh,mass,peripheral,radical点之和，因为这一步还没有计算，所以没有找到任何点。
Using eigen values from OUTCAR to determine the following infor #从OUTCAR里面的K点数据来计算下表中的数据
Gap       HOMO      LUMO      H_width   L_width   H_kpt               L_kpt
1.1188    -4.9184   -3.7996   0.1497    0.2720    0.0,0.0             0.0,0.0 #分别是能带，HOMO，LUMO，价带宽度，导带宽度，价带顶点坐标，导带底坐标（这两个坐标程序会自动选取价带顶和导带底，也可以手动选。有可能2d能带算了以后会找到一个更高的价带顶，那么就需要重新算有效质量）
step size for h, e are 0.060, 0.035 with delte E 0.020000 #计算有效质量的步长，因为价带宽度小一些（平滑一些），所以改变相同的能量（0.02）所需要的步长就大一些，用-s选项设置能量变化
angle between a b is 59.9979
length of a, b is 0.4273, 0.4273
hex reciprocal lattice detected with precision 1 #根据ab和夹角检测到是一个正六变形的晶胞, 这里对称性高的话，下面所需要算的k点数量也少，可以用-p选项提高precision 来强制降低对称性
number of kpts along a is 8 with step size (-G) of 0.050000 #根据间隔是0.05计算出沿着a轴要算8个k点， 这部分k数量点包括了peripheral 和 radical，决定了1维能带是否平滑，用-G选项改变间隔
number of grid pt in x is 4 with step size (-g) of 0.050000 #根据间隔是0.05 x和y方向需要各算4个k点， 这部分k数量点包括了mesh，决定了2维能带是否平滑，用-g选项改变间隔
number of grid pt in y is 4 with step size (-g) of 0.050000

另外还可以产生一个图片如下，其中的小黑点表明了mesh k点的密集度
，

3.产生输入文件
python band2d.py -j i
然后程序会创建新的文件夹，复制WAVECAR和CHGCAR，产生INCAR，KPOINTS
完成后会发现当面目录下多了一串 2D开头的文件夹



4.计算
cd到每个2D开头的目录下，运行vasp计算。这里我是用脚本批量提交到服务器上算。输入vasps 2D* 就行。但是每个服务器的作业调度系统不一样。所以这个脚本要根据具体的服务器来编写。

5.产生结果文件
算完以后将所有的文件拷贝到本地机器上，一般我是把CHGCAR和WAVECAR 删了再拷。
for i in 2D*
do
rm $i/WAVECAR
rm $i/CHGCAR
done
然后在当前目录（当前目录包括了所有2D开头的文件夹以及之前单点/优化计算的结果）下运行
python band2d.py -j p
在终端会产生下面的输出：
10 k points found in OUTCAR
Warning! Number of points (12) in mesh line 0 not equal to 0.
Change the grid_np_x to 12 now
12 X 4 points detected for mesh grid
25 X 2 points detected for eff. mass grid
7 points detected for peripheral band
11 13 points detected for radical band
129 k points found in 2D data #可以看到这里找到了很多新算的k点
Using 2D data to determine the following infor #采用2Dk点来决定下面的属性，可以看到导带和价带都变宽了，这是布里渊区取样更全导致。
Gap       HOMO      LUMO      H_width   L_width   H_kpt               L_kpt
1.1188    -4.9262   -3.8074   0.2324    0.4061    0.0,0.0             0.0,0.0
step size for h, e are 0.040, 0.025 with delte E 0.020000
angle between a b is 59.9979
length of a, b is 0.4273, 0.4273
hex reciprocal lattice detected with precision 1
number of kpts along a is 8 with step size (-G) of 0.050000
number of grid pt in x is 12 with step size (-g) of 0.050000
number of grid pt in y is 4 with step size (-g) of 0.050000
12 X 4 points detected for mesh grid
25 X 2 points detected for eff. mass grid
7 points detected for peripheral band
11 13 points detected for radical band
然后会依次产生开头出现的那四个图。调整大小后保存就可以了。另外会产生如下的数据文件
DATA_BAND1D.out  绘制1维能带的数据文件（包含所有的能带）
DATA_BAND2D.out  绘制2维能带的数据文件（只包含当前绘制的能带，在--range3d里面设置）
DATA_MASS_K.out  用来计算有效质量的K点能量数据
DATA_MESH_K.out  用来产生2维能带的K点能量数据
DATA_PERI_K.out  用来产生1维外周能带的K点能量数据（外周能带是指那些不经过Γ点的K路径）
DATA_RADI_K.out 用来产生1维径向能带的K点能量数据（径向能带是指那些经过Γ点的K路径）
注意：
使用band2d.py -h 查看可用的选项和帮助
下面举一个使用例子
band2d.py -j p -g 0.03,0.03 -G 0.04 --range2d h5,l2 --range3d h2,h=-5,0  -p  0120340
说明：-j p 运行处理数据的任务，p for process
-g 0.03,0.03 对于2维能带，在x,y方向上的取样间隔是0.03
-G 0.04 对于1维能带，取样间隔是 0.04 （单位都是angstrom^-1)
--range2d h5,l2对于1维能带，绘制HOMO-5 到 LUMO+5 的能带
--range3d h2,h=-5,0 对于2维能带，绘制HOMO-2 和HOMO的能带，然后把z坐标的范围改成-5到0eV
-m h1,h,e 计算HOMO-1，HOMO，LUMO 三条能带的有效质量
-p 0120340 对于1维能带，k点的路径是从Γ点（0）出发，到外周的第一个点（1），然后逆时针转到外周第二个点（2），然后回到Γ，然后到外周的第三个点（3）。。。

遇到bug，或者有使用问题，或者有新功能建议在本帖留言或者加本人QQ：32589927。

不取俗名

很感兴趣，楼主介不介意放一两个例子，有相关的2D DOS的文章那就更好了，多谢

helpme

这个估计会用到，先备份一个。

lindlar

谢谢，学习啦

luckstar18

good, vary good!
Thanks for sharing!

lindlar

先收藏了，非常好的应用，谢谢啦！

solarman

很好的东西，赞一个，不过应该称之为二维布里渊区的三维能带图更合适些？

追枫少年

您好，我是编程的小白，在尝试您的脚本的时候提示如下错误：
  File "band2d.py", line 731
    for i in range(len(self.radi_frac))}.copy()
      ^
SyntaxError: invalid syntax
望解答是什么原因？先拜谢了

ggdh

追枫少年 发表于 2017-5-10 09:55
您好，我是编程的小白，在尝试您的脚本的时候提示如下错误：
  File "band2d.py", line 731
    for i in ...

你装了anaconda没有
具体哪一步运行出错的？用 python band2d.py 的方式运行试试看

追枫少年

ggdh 发表于 2017-5-10 17:52
你装了anaconda没有
具体哪一步运行出错的？用 python band2d.py 的方式运行试试看

我安装的是anaconda2  版本对这个脚本有影响吗？
我尝试了 Python band2d.py的命令   还是提示相同的错误

ggdh

追枫少年 发表于 2017-5-10 19:44
我安装的是anaconda2  版本对这个脚本有影响吗？
我尝试了 Python band2d.py的命令   还是提示相同的错 ...

anaconda2 是没有问题的
你试试输入
python --version 看看版本是不是2.7.12
另外 使用dos2unix band2d.py
另外 重新下载band2d.py
最后可以加我qq 32589927 ，有问题可以实时跟我反馈

追枫少年

ggdh 发表于 2017-5-11 23:23
anaconda2 是没有问题的
你试试输入
python --version 看看版本是不是2.7.12

非常感谢得到您的指点

shuihuntian

谢谢楼主分享，请问能处理hse06的结果（计算有效质量）吗？

ggdh

shuihuntian 发表于 2017-5-29 09:09
谢谢楼主分享，请问能处理hse06的结果（计算有效质量）吗？

可以

luobo123456

谢谢楼主分享 请问可以把图片的数据导出来吗？