||
关注:
1) INCAR解读
2) NEB原理及使用方法:http://sobereva.com/44
3) 新手入门
一些认识摘录:
NEB 是 chain-of-states 搜索方法的一种,经过长时间的发展形成现在的准确性好,稳定性好的 CINEB 方法。做法就是在反应物到产物之间插入一系列结构,共插入 P-1 个,反应物编号为 0,产编号物为 P。不同的是优化不是对每个点孤立地优化,而是优化一个函数,每一步所有点一起运动。
CI-NEB方法只需要很少的点,比如包含初、末态总共5个甚至3个点就能准确定位过渡态,是最有效率的寻找过渡态的方法之一。
用VTST中CI-NEB做示范,(首先要确保编译VASP的时候安装了VTST:IBRION = 3, POTIM = 0,这是VTST识别并启动VTST优化算法的标志。
IBRION = 3 # Disbale VASP default optimizers which are basesd on Energy
POTIM = 0 # Disbale VASP default optimizers which are basesd on Energy
【有利于爬坡,寻找能量最高点的过渡态?】
LCLIMB = .TRUE. 爬坡即CI-NEB
注意用的核的个数必须整除插点的数量,(为保证最大的并行效率节点数最好等于差点的个数)。
计算过程中,可以用nebef.pl
或者nebefs.pl
命令(不带参数)查看计算收敛情况:
~/yexq/work/neb2020-al333-he-ois $ nebefs.pl
Image Force Stress Volume Magnet Rel Energy
0 0.00000000 0.00000000 1797.06 0.00000000 0.00000000
1 1.42943000 0.36451400 1797.61 0.00000000 0.23100000
2 7.22470900 4.55774200 1797.80 0.00000000 0.55940000
3 40.22819000 0.84127900 1821.86 0.00000000 3.74250000
4 38.93652200 8.65748600 1843.88 0.00000000 61.48530000
5 657.29611100 8.26608500 1825.69 0.00000000 266.41220000
6 33.40365500 0.85736200 1821.87 0.00000000 4.43150000
7 7.22617100 4.55904700 1797.80 0.00000000 0.55930000
8 1.42975700 0.36507000 1797.61 0.00000000 0.23110000
9 0.00000000 0.00000000 1797.06 0.00000000 0.00000000
做NEB的正确方法是:
1.优化反应物与产物结构;
2.用线性插值方法,在反应物和产物间插入若干结构(image)(通常不会只有一个)
3.按NEB方法设置几个参数计算就可以了。计算完成后OUTCAR里有收敛提示
4.那些脚本用不用都可以,0x文件夹下的每个OUTCAR里都有前后两个image间的距离,可以用来做reaction path
5.得到的是不是过渡态,还需要用频率来检验
NEB: distance to prev, next image, angle between1.501237 1.501237 那几个images的能量值在哪查看?
5. 频率的验证是不是进行Gamma点 vibrational frequency的计算?
我的做法是把neb算好的文件夹中的CONTCAR取出来改名为POSCAR
然后POTCAR和KPOINTS不变
INCAR中将NSW改为0(之前计算NEB时设的是2000)
IBRION改为5(在手册上看到这个是计算频率的)
NSW=1
IBRION=5
NFREE=2
POTIM=0.05
如果你使用的vasp为5.2版本,改IBRION为6,vasp会为你寻找体系对称性,可以极大减少计算量
INCAR解读-from Jinfan
SYSTEM = XO2
#starting parameter
ISTART = 0
ICHARG = 2
INIWAV = 1
#Default: INIWAV = 1 Fill wavefunction arrays with random numbers. Use whenever possible
LCHARG=.FALSE.
LWAVE=.FALSE.
#PREC = Accurate
GGA = PE
#VOSKOWN = 1
# If VOSKOWN is set to 1 the interpolation formula according to Vosko, Wilk and Nusair[53] is used. This usually enhances the magnetic moments and the magnetic energies.
IVWD = 12
#spin
#ISPIN = 2
#MAGMOM = 45*1 45*0
#electronic
ENCUT = 500
EDIFF=5E-5
NELM = 200 #default, if not converged then test the parameters
NELMIN = 3
#ALGO = Normal #eletronic minimization algorithum
#IALGO = 38
#mixing paramters
#ionic
EDIFFG=-0.05
NSW=2000
IBRION=3
POTIM=0
ISIF = 2
#ISYM = 2
#dos
ISMEAR = 1
SIGMA = 0.2
#others
LREAL = Auto
#ALGO = Fast
NPAR = 8
#VOSKOWN = 1
LDAU=.TRUE.
LDAUTYPE=2
LDAUL=3 -1 -1
LDAUU=4.5 0 0
LDAUJ=0.5 0 0
#LDAUPRINT=2
#LWAVE = .FALSE.
#LCHARG = .TRUE.
#LREAL = Auto
LMAXMIX = 6
#NPAR = 4
#FOR NEB
IOPT=3
ICHAIN=0
LCLIMB=.TRUE.
SPRING=-5
IMAGES=7
LTANGENTOLD = .FALSE.
LDNEB = .FALSE.
TIMESTEP = 0.1
MAXMOVE = 0.2
参数解释:
IBRION=3
For difficult relaxation problems it is recommended to use the conjugate gradient algorithm (IBRION=2), which presently possesses the most reliable backup routines. Damped molecular dynamics (IBRION=3) are often useful when starting from very bad initial guesses.
https://www.researchgate.net/post/NEB_Calculation_using_VASP
Dear Berna Akgenc,
I will try my best to answer your question.
You need to have your initial and final structures well relaxed. In this case, it means POSCAR00 and POSCAR06 . Then you use the VTST script nebmake.pl to create the number of images. I understand you want to have 5 images thus you do:
nebmake.pl POSCAR00 POSCAR06 5.
thus you will have 7 folder (00,01,02,03,04,05,06). For further analysis you need to put the OUTCARs for POSCAR00 and POSCAR06 in the folders 00 and 06.
Then you need to make sure that your vasp executable is complied with VTST (http://theory.cm.utexas.edu/vtsttools/neb.html#neb). unless you can not make a CI-NEB.
Then in your parent directory where you have all the (00,01,02,03,04,05,06) directories, you put your INCAR, KPOINTS and POTCAR file. But this time, your INCAR file will be a little be different as you need to include the NEB tags on it. And me personally use the spring constant be -5.
A sample of INCAR file is here in the following:
-------------------------------------------------------
SYSTEM = NEB-Calc
ISYM = 0 # switch on symmetry
#Start parameter for this Run:
ISTART = 0
ICHARG = 2
#Electronic Relaxation 1
ENCUT = 500 eV
IALGO = 38
EDIFF = 1.0E-6
NELMIN = 5
#Ionic Relaxation
EDIFFG = -1.0E-3
NSW = 500
#DOS related values:
ISMEAR = 0
SIGMA = 0.05
#Cell Volume + Cell shape
ISIF = 3
#Other Calculation parameters
PREC = ACCURATE # Advised to do so
LREAL = .FALSE.
LWAVE = .FALSE.
LCHARG = .FALSE
#NEB Related
IMAGES = 5 # Number of NEB images between the fixed endpoints
SPRING = -5.0 # The spring constant, in eV/Ang^2 between the images; negative value turns on nudging
ICHAIN = 0 # Indicates which method to run. NEB (ICHAIN=0) is the default
LCLIMB = .TRUE. # Flag to turn on the climbing image algorithm
LNEBCELL = .TRUE. # Flag to turn on SS-NEB. Used with ISIF=3 and IOPT=3
# Must set IOPT = 3 or 7 when using LNEBCELL=.TRUE.
IOPT = 3 # QM (Quick-Min) force based optimizers is used and good for high force far from minmum
# O(default vasp),1(LBFGS),2(CG),4(SD),7(FIRE)
IBRION = 3 # Disbale VASP default optimizers which are basesd on Energy
POTIM = 0 # Disbale VASP default optimizers which are basesd on Energy
【有利于爬坡,寻找能量最高点的过渡态?】
#MAXMOVE = 0.2 # Maximum allowed step size for translation(Default comes with IOPT=3)
TIMESTEP = 0.1 # Dynamical time step (Default comes with IOPT=3)
LSCALAPACK = .FALSE.
NPAR = 4
-------------------------------------------------------
This INCAR is for calculation allowing a variable CELL with Climbing NEB.
Hope it helps.
Regards,
Abel
新手入门
http://bbs.keinsci.com/thread-7489-1-1.html
最好使用vtst tool版本的vasp,里面有很多好用的算法,比原版vasp好用。编译方法见官网:Transition State Tools for VASP
以下NEB计算的方法只介绍vtst tool版本的vasp
NEB计算步骤如下:
1. 准备好初末构型的vasp文件: is.vasp,fs.vasp,以及脚本:nebmake.pl(来源:Transition State Tools for VASP),然后输入命令:
> nebmake.pl is.vasp fs.vasp 8
复制代码
最后一个8代表生成8个image,如果运行成功,会在当前文件夹生成编号为00-09的十个文件夹,每个文件夹包含一个POSCAR。编号00和09的文件夹对应is.vasp和fs.vasp。
之后把配置文件INCAR,KPOINTS,POTCAR放到当前文件夹。INCAR文件需要有以下字段:
IMAGES=8
LCLIMB = .TRUE.
复制代码
IMAGES代表有8个image,LCLIMB是是否开启climbing image方法,详见官网。
nebmake.pl会用线性插值的方法在is.vasp和fs.vasp之间插入若干个image,这会导致其中某些image的构型极度不合理,比如某两个原子重叠在了一起,需要调整。可以自己把01-08文件夹中的POSCAR文件下载下来,调整后再上传,但比较麻烦。推荐使用我写的脚本nebinfo(来源:Mabinogiysk/VASP-script),输入image的数量之后,会生成一个压缩包NEBINFO.tar.gz. 里面包含有每个image的vasp格式文件,gjf格式文件以及xyz格式的动画文件(可用VMD观看)。
调整好image之后在当前文件夹执行vasp命令就可以开始NEB计算。
刚开始算的时候,由于构型不合理,因此选用最速下降法会比较快,也不容易出错。此时选择IOPT=7,跑个二三十步之后再换用IOPT=1,当然也可以从头到尾只用IOPT=1.
对于复杂的过渡态,NEB计算需要非常多的image才能收敛。因此我一般只跑个五六十步,根据经验来看,后面基本不会有太大变化了。
之后可以采用vtst tool中的脚本来输出最后一步NEB的信息,如能量、切线力等,用法请到官网查找。也可以使用我脚本仓库中的nebvtst.pl脚本:
> nebvtst.pl 8
复制代码
该脚本将输出NEB计算每一步每个image的信息:
steps: 60
images: 1 -314.70152 -0.20180 0.07899
images: 2 -314.66766 -0.10107 0.08997
images: 3 -314.21920 -0.03095 0.06242
images: 4 -314.82366 0.36034 0.08089
images: 5 -315.02908 0.17184 0.10215
images: 6 -315.14834 0.10812 0.06171
images: 7 -315.21699 0.05799 0.05365
images: 8 -315.25374 0.02400 0.12696
复制代码
三列数字从左到右分别是能量,切线力和受力最大的原子所受的力(都小于EDIFFG时NEB计算就收敛了)。
其中最重要的是切线力,切线力为零的点就是过渡态。正常的NEB计算切线力应该是从负值到正值或者从正值到负值,并且能量应该先升高再下降。不满足这两个条件时意味着计算很可能出错了。
对于上面这个计算,可以看到image3的切线力很接近0了,也就意味着这个image跟过渡态很接近了,可以把它单独拿出来做个优化(IBRION=1),有可能就收敛到过渡态了。
如果所有的image的切线力都比0差很多,如下面这样:
steps: 42
images: 1 -314.71090 -0.19732 0.67034
images: 2 -314.68128 -0.10245 0.18674
images: 3 -314.22095 -0.23067 0.47518
images: 4 -314.80703 0.24911 0.40551
images: 5 -315.01767 0.16162 0.27819
images: 6 -315.13853 0.10266 0.38136
images: 7 -315.21270 0.04747 0.28035
images: 8 -315.25231 0.02491 0.24664
复制代码
就只能把零两边的image,即image3和image4拿出来重新用nebmake.pl插值再跑NEB,直到找到切线力接近零的点。这是原始版本的vasp做过渡态的方法,如果编译了vtst tool工具,采用climbing image方法,即设置LCLIMB = .TRUE. 基本上跑一遍NEB之后都会有一个image的切线力接近零,避免了重复插值。因此用vasp做过渡态一定要编译上vtst tool.
另外经常遇到的情况是切线力接近零的image拿出来做优化之后,发现力总是降不下去。此时NEB方法就无能为力了,需要另一个强大的过渡态搜索方法Dimer上场了。这就留做下一篇帖子来说了。
http://blog.sina.com.cn/s/blog_b364ab230102vghk.html
先来原理:
Eyring equation(势垒与温度决定反应速率)
NEB方法(撒点,弹簧力;详情参考下面给出的文献)
CL-NEB方法(爬坡的NEB,使得其中一个点跑到势能面鞍点,即过渡态)
VASP-vtst提供的NEB优化算法(IOPT参数)
vtst脚本包
一些用perl语言写的脚本,不错的处理neb计算的工具箱。详情请进:http://theory.cm.utexas.edu/vtsttools/scripts.html
后面所有后缀名为pl的脚本都来自这里。
步骤一:初末态结构优化
建立两个文件夹ini(初态的意思),fin(末态),每个文件夹放入vasp计算必备的四个文件(INCAR,POSCAR,KPOINTS,POTCAR),其中的两个POSCAR对应未优化的初末态。确保两个文件夹里面除POSCAR外,其他文件完全一样。
对于两个POSCAR,注意每一行的原子一一对应。若有固定位置的原子(比方做表面计算,要固定底部1-2层原子),最好让两个POSCAR里面的固定原子位置完全相同。
下图是我的两个优化结构
ini fin
用dist.pl检查两个优化后结构(两个CONTCAR)的相似程度(每个对应原子的初末态距离的平方和,再开根号)。
具体命令为: dist.pl ini/CONTCAR fin/CONTCAR
若返回值小于5A,一般可以进行下一步。我的情况:
步骤二:插点
插点的数目取决于前面dist.pl的返回值,一般插点数目可取(dist.pl返回值/0.8)。
具体插点的命令: nebmake.pl ini/CONTCAR fin/CONTCAR N
最后的N表示插点数目。插入点的算法为线性插值,详情请进前面给的vtst脚本链接。我这里为简单插一个点,执行完命令后出现文件夹00,01,02.
00表示初态,里面放的是ini/CONTCAR,02表示末态,放的是fin/CONTCAR, 01是插入的点。三个文件夹里面的文件名称都是POSCAR。
输出中提示讲初末态对应的OUTCAR复制到对应的文件夹中,以便后续数据分析。
检查插入点的合理性:
输入命令 nebmovie.pl 0
参数0表示用POSCAR生成xyz文件;还可取1,为用CONTCAR生成。
可以看到每个文件夹内多了一个后缀名为xyz的文件。将这几个xyz文件合成一个movie.xyz文件即可在VMD中动画演示。我更倾向于使用jmol来查看xyz文件,不过VMD可能大家更熟悉。
用nebmovie.pl,没有直接生成movie.xyz是因为从官方主页下载的脚本默认当使用nebmovie.pl后面不带参数或者参数为0的时候(即使用POSCAR产生xyz),不输出movie.xyz。 你想让脚本自动生成movie.xyz而不是自己去合并各个文件夹里面的xyz文件,需要自行修改nebmovie.pl文件。很简单,把倒数第二个if语句整个用#注释掉或者直接删掉即可。
也可以直接用VESTA查看01的POSCAR。
进一步可用命令nebavoid.pl 1,确保中间插入的点每一个原子间距都大于1A。该命令的参数1表示最小允许间距,可取小数。
上面输出表示没有原子间距小于1A。
原子间距太小说明结构有问题。要么是因为初末态POSCAR里的原子没有对应好,要么是因为nebmake.pl线性插值的方法不适用于估计反应路径路径。如果是前者,调整好对应原子,也就是CONTCAR里面相应行进行调换;如果是后者,用Materials studio之类的软件调整好结构,再放入对应文件夹。
步骤三:其它文件的准备
在当前目录下面放入KPOINTS,POTCAR及INCAR文件。要求KPOINTS和POTCAR与ini、fin文件夹中的一样;INCAR中的基本参数也与初末态计算的INCAR一样,另外再加入NEB计算的相关参数。
我的结构优化的INCAR:
我的NEB计算的INCAR:
INCAR再多说两句。优化算法若想使用vtst的,需要如上设置IBRION=3,POTIM=0,IOPT设置见本文第四个图。若使用vasp自带的优化器,可以使用IBRION=1或者3,不要取2;POTIM取合适的值(0.01-0.5范围内去尝试)。
这样就可以提交NEB任务了。
步骤四:后续处理
计算过程中,可以用nebef.pl命令(不带参数)查看计算收敛情况。
输出中,第二列即为最大受力(force of images in the neb),第三列为相应结构的能量。
前面INCAR中EDIFFG设置为-0.02(一般结构优化可取-0.02,NEB计算取-0.03),表示当所有结构最大力小于0.02eV/A时,结束计算。上图中的力都小于0.02,此时已经完全收敛,计算结束。
另一个可以用来观察收敛情况的命令是nebbarrier.pl(同样不带参数)。该命令没有输出到屏幕的内容,而是生成neb.dat文件。
neb.dat文件第二列表示距离(即临近两结构的dist.pl的计算结果),第三列表示能量(以初态能量为参考值),第四列为力(forces along the neb)。
EDIFFG参数对应的力是nebef.pl输出中的force of images in the neb.
计算完成后,使用命令nebresult.pl.
nebresult.pl做的事情如其所输出表明的,完成了nebbarrier.pl, nebspline.pl, nebef.pl, nebmovie.pl, nebjmovie.pl, nebconverge.pl还有对各文件夹中的OUTCAR打包压缩。生成了很多文件。其中mep.eps是以dist.pl距离为横坐标,能量为纵坐标画出的能势垒图,如下:
该文件一般查看图片的软件打不开,推荐使用EPS/PS viewer.
生成的spline.dat文件是对上面几个点的拟合曲线数据,已经有了上图,这些数据就意义不大了。
生成的vaspgr文件夹内是各个插点结构的收敛图(同样可用EPS/PS viewer打开),如下(我的只插入一个结构,所以只有一个图):
结束语:
NEB计算,收敛是个核心问题,可以参见我的涉及NEB计算的其他博文。
http://blog.wangruixing.cn/2019/08/19/cineb/
过渡态(transition state)搜索,优化到一阶鞍点。是MEP路径上能量最高的点(一阶鞍点)。在过渡态位置上,仅在沿着反应路径上是能量极大点,而在于之正交的其他所有方向上都是能量极小点。
energy pathway),矛盾得到解决。
movie.xyz 直接拖入Jmol可以看所有结构。
计算过程中,可以用nebef.pl
或者nebefs.pl
命令(不带参数)查看计算收敛情况:
nebefs.pl
:(有多少个点就会显示几行,这是插了4个点的算例)Image Force Stress Volume Magnet Rel Energy
0 0.00000000 0.00000000 1797.06 0.00000000 0.00000000
1 1.42943000 0.36451400 1797.61 0.00000000 0.23100000
2 7.22470900 4.55774200 1797.80 0.00000000 0.55940000
3 40.22819000 0.84127900 1821.86 0.00000000 3.74250000
4 38.93652200 8.65748600 1843.88 0.00000000 61.48530000
5 657.29611100 8.26608500 1825.69 0.00000000 266.41220000
6 33.40365500 0.85736200 1821.87 0.00000000 4.43150000
7 7.22617100 4.55904700 1797.80 0.00000000 0.55930000
8 1.42975700 0.36507000 1797.61 0.00000000 0.23110000
9 0.00000000 0.00000000 1797.06 0.00000000 0.00000000
步骤十:频率计算验证虚频(下节课)
IBRION = 5 |
nebresults.pl也是VTST脚本,用来总结neb结果
nebspline.pl,
nebef.pl,
nebmovie.pl,
nebjmovie.pl,
nebconverge.pl还有对各文件夹中的OUTCAR打包压缩。生成了很多文件。
其中mep.eps是以dist.pl距离为横坐标,能量为纵坐标画出的能势垒图,如下:
Archiver|手机版|科学网 ( 京ICP备07017567号-12 )
GMT+8, 2023-9-29 15:33
Powered by ScienceNet.cn
Copyright © 2007- 中国科学报社