vasp计算吸附要多久完成(vasp计算声子谱教程)
vasp计算吸附要多久完成(vasp计算声子谱教程)2.2 phonopy code编译计算所用的软件为VASP与phonopy code。直接法,或称frozen-phonon方法,是通过在优化后的平衡结构中引入原子位移,计算作用在原子上的Hellmann-Feynman力,进而由动力学矩阵算出声子色散曲线。直接法的缺陷在于它要求声子波矢与原胞边界supersize正交,或者原胞足够大使得Hellmann-Feynman力在原胞外可以忽略不计。这使得对于复杂系统,如对称性高的晶体、合金、超晶格等材料需要采用超原胞。超原胞的采用使计算量急剧增加,极大的限制了该方法的使用。1987年,Baroni、Giannozzi和Testa提出了DFPT方法。DFPT通过计算系统能量对外场微扰的响应来求出晶格动力学性质。该方法最大的优势在于它不限定微扰的波矢与原胞边界正交,不需要超原胞也可以对任意波矢求解。Castep、Vasp等采用的是一种linear
1、前言
在第一性原理计算过程中考察体系稳定性是经常遇到的一个问题,也是大部分审稿人容易提问的地方;通常,使用声子谱研究体系的动力学稳定性,使用分子动力学研究体系的热稳定性;另外,可以借助于一些其它方法进一步说明体系的稳定性,例如研究钙钛矿体系时通过离子半径判断稳定性,研究异质结时通过计算结合能判断稳定性等。
2、声子谱计算方法
声子谱的计算主要有两种方法,一种是直接法,另一种是微扰密度泛函方法(DFPT)。
直接法,或称frozen-phonon方法,是通过在优化后的平衡结构中引入原子位移,计算作用在原子上的Hellmann-Feynman力,进而由动力学矩阵算出声子色散曲线。直接法的缺陷在于它要求声子波矢与原胞边界supersize正交,或者原胞足够大使得Hellmann-Feynman力在原胞外可以忽略不计。这使得对于复杂系统,如对称性高的晶体、合金、超晶格等材料需要采用超原胞。超原胞的采用使计算量急剧增加,极大的限制了该方法的使用。
1987年,Baroni、Giannozzi和Testa提出了DFPT方法。DFPT通过计算系统能量对外场微扰的响应来求出晶格动力学性质。该方法最大的优势在于它不限定微扰的波矢与原胞边界正交,不需要超原胞也可以对任意波矢求解。Castep、Vasp等采用的是一种linear response theory 的方法(或者称为 density perturbation functional theory,DFPT),直接计算出原子的移动而导致的势场变化,再进一步构造出动力学矩阵。进而计算出声子谱。
2.1 计算环境的搭建
计算所用的软件为VASP与phonopy code。
2.2 phonopy code编译
(1) 搭建Anaconda环境,下载链接(https://www.anaconda.com)。
(2) 安装Anaconda
bash Anaconda3-2019.07-Linux-x86_64.sh
(3) 写入环境变量
echo "export PATH=/……/anaconda3/bin:$PATH" >> ~/.bashrc
source ~/.bashrc
Note: /……/表示自己服务器路径
(3) 下载phonopy code,下载链接(https://pypi.org/project/phonopy/)
(4) 编译phonopy code
tar zxvf phonopy-2.2.0.tar.gz
cd phonopy-2.2.0
python3 setup.py install --user
(5) 写入phonopy的环境变量
echo "export PATH=/……/phonopy-2.2.0/build/scripts-3.7:$PATH" >> ~/.bashrc
source ~/.bashrc
Note:/……/表示自己服务器路径
(6) 检查phonopy code是否安装完成
$ phonopy
2.3 声子谱计算步骤
(1) 结构优化
声子谱的计算需要对原胞结构做高精度的充分优化,否则很容易出现虚频;下面以二维InSe为例,详细说明声子谱的计算步骤。
INCAR
SYSTEM = 2D_InSe
ISTART = 0
NWRITE = 2
PREC = Accurate
ENCUT = 500
GGA = PE
NSW = 200
ISIF = 3
ISYM = 2
NBLOCK = 1
KBLOCK = 1
IBRION = 2
NELM = 80
EDIFF = 1E-08
EDIFFG = -0.001
ALGO = Normal
LDIAG = .TRUE.
LREAL = .FALSE.
ISMEAR = 0
SIGMA = 0.02
ICHARG = 2
LPLANE = .TRUE.
NPAR = 4
LSCALU = .FALSE.
NSIM = 4
LWAVE = .FALSE.
LCHARG = .FALSE.
ICORELEVEL = 1
POSCAR
2D_InSe
1.0
4.0836000443 0.0000000000 0.0000000000
-2.0418000221 3.5365013772 0.0000000000
0.0000000000 0.0000000000 25.3775005341
In Se
2 2
Direct
0.666670026 0.333330010 0.589979992
0.666670026 0.333330010 0.478789968
0.333329978 0.666669998 0.428439986
0.333329978 0.666669998 0.640340007
KPOINTS
Monkhorst Pack of Gamma centered
0
Gamma
13 13 1
0.0 0.0 0.0
POTCAR
cat In_d/POTCAR Se/POTCAR > POTCAR
OPTCELL
100
110
000
不懂OPTCELL的设置请看刘博博客:
http://blog.wangruixing.cn/2019/05/05/constr/
(2) 通过phonopy code建立超胞
a) 将优化后的结构文件拷贝为POSCAR
cp CONTCAR POSCAR
b) 建立超胞
$ phonopy -d --dim="4 4 1"
$ ls
phonopy_disp.yaml POSCAR-002 POSCAR-005 POSCAR-008 POSCAR-011 POSCAR-014 SPOSCAR
POSCAR POSCAR-003 POSCAR-006 POSCAR-009 POSCAR-012 POSCAR-015
POSCAR-001 POSCAR-004 POSCAR-007 POSCAR-010 POSCAR-013 POSCAR-016
Note: 建立超胞之前,一定要在MS中找好对称性,否则会增加大量的计算量
c) POSCAR和SPOSCAR的重命名
cp POSCAR POSCAR-unitcell
cp SPOSCAR POSCAR
(3) 计算力学Hessian矩阵
力学Hessian矩阵可由VASP计算得到,该矩阵保存于VASP的输出文件vasprun.xml中;计算过程中KPOINTS文件根据服务器计算能力,可适当增加,POTCAR文件无需改变,INCAR文件设置如下:
INCAR
SYSTEM = 2D_InSe
ISTART = 0
NWRITE = 2
IBRION = 8
NSW = 1
IALGO = 38
NELM = 200
EDIFF = 1E-07
EDIFFG = -0.001
ISMEAR = 0
SIGMA = 0.02
ENCUT = 500
PREC = Accurate
LREAL = .FALSE.
LWAVE = .FALSE.
LCHARG = .FALSE.
ADDGRID = .TRUE
(4) 数据处理
a) 根据vasprun.xml文件生成力学文件FORCE_CONSTRAINS
phonopy --fc vasprun.xml
b) 编辑band.conf文件,该文件给出了高对称点路径的信息
ATOM_NAME = In Se
DIM = 4 4 1
BAND = 0.5 0.0 0.0 0.0 0.0 0.0 0.333333 0.333333 0.0 0.5 0.0 0.0
FORCE_CONSTANTS = READ
Note:
line 1:元素名称,顺序与POSCAR保持一致
line 2:建立超胞的大小
line 3:高对称点路径,每组高对称点之间用两个空格隔开
line 4:表示读取力常数文件FORCE_CONSTANTS
c) 生成band.yaml文件
phonopy --dim="4 4 1" -c POSCAR-unitcell band.conf
Note:4 4 1为建立超胞的大小
d) 得到声子谱数据文件PBAND.dat文件
phonopy-bandplot --gnuplot> PBAND.dat
Note: 高对称点标注说明:phonopy软件默认在两个高对称点之间打点51个,且在PBAND.dat中每组高对称点数据以一个空行分割,据此即可得到完整的声子谱图像
e) 计算结果
文献结果:
3、后记
如果计算完的声子谱其他地方都好,就是在Γ点有一点点虚频,那么这个材料很可能是稳定的,只是你优化做得不够好,进一步提高优化的精度可消除这一点点虚频。
对于二维材料,如果在Γ点出现很小的虚频,基本可以认为这个材料是稳定的,大部分二维材料都会有此现象;尤其是VASP结合phonopy code计算二维材料的声子谱在Γ点更是容易出现虚频;使用Quantum-ESPRESSO的PWSCF和PH模块计算声子谱对于内存的需求较小一些,且对于二维材料的声子谱计算更友好一些,后期会详细介绍Quantum-ESPRESSO计算声子谱的具体步骤。
请大家关注贺勇的个人博客网站:
https://yh-phys.github.io