加载中…NEB计算反应路径力、角、距离和坐标的意义(附图及vststool脚本说明)
(2012-05-17 10:48:19)
标签：
杂谈 |
分类： Vasp |
相信做渡态计算的虫子大部分读过Sobereva的《过渡态、反应路径的计算方法及相关问题》。我反复读过，还是觉得困难重重。因为，这篇文章并不提供详细的数学矢量推导，初学者很难获得具体细节性的理解。我有如下几个问题希望路过的虫子给予评论，看我理解的对不对：
（1）什么是反应坐标（reaction coordinates）？
从包含N个原子的反应物到产物的过程记为：A->B
我们只知道A和B的所有原子坐标，插入若干中间点后得到一条反应路径。这条路径上，若不考虑原子的固定，其内坐标变量为3N-6.我们常常以一条线来表示反应路径，这条线是不是3N-6维空间中的线呢？如果再加上能量维度，这条反应路径是不是3N-5维空间中的曲线？
（2）在(C)NEB计算中，各image是怎样移动的？
我的初浅理解，是各个image在垂直反应路径方向受势能力而移动，这个与结构有过程类似：小小的移动一个离子步，然后求能量对坐标的变化率，得到image在垂直反应路径的势能力。而在反应路径方向，各个image之间受弹簧力作用。
一般三个紧邻的image之间成一定的（0~180）夹角，中间的image受到前后两个弹簧的合力。优化过程为把image受到的弹簧力分解，保留image在反应路径上切线方向的力。
所以，每个image有两个运动：受势能力驱动的，垂直于反应路径上的运动；和受弹簧力驱动的，沿着反应路径切线方向的运动。
优化最终的结果，应该是每个image受到的势能力和弹簧力小于设定标准。那么，具体到vasp的NEB计算nebbarrier.pl的输出第四列，说是沿着反应路径方向（从反应物指向产物的有向曲线）的受力。这个力为正，表示沿着反应路径想产物移动能量增加，否则相反。
这里提到的力，包括势能力(A1)和弹簧力(B1)，加上二者在反应路径切线平行方向(A11,B11)和垂直方向(A12,B12)的投影，力的收敛标准（EDIFFG）到底是以什么为对象呢？
（3）image之间的夹角
三个相邻的image成一个角，这个角的角度如何定义？如果三个image依次命名为i-1,i和i+1，对应坐标为R(i-1),R(i)和R(i+1)，角度是不是矢量到与矢量的方向一致所需转动的最小角度（0-180度）？
（4）image之间的距离
三维空间两点距离是对应坐标项差的平方和的开方。对于原子个数众多且变化的体系，为了距离的可比较，似乎求的是二阶矩：对应原子坐标项的差的平方和的开方对原子个数的平均。
（5)image的受力和收敛力判据到底指什么？
下面9楼和10楼贴了两个跟力相关的帖子。nebef.pl的屏幕输出得到的力是每个image中的FORCES: max atom.难道收敛判据是每个image中原子最大受力？为什么不是每个image所受的弹簧力沿切线方向分量和势能力沿垂直切线方向分量？或者后两者的合力？
我觉得，反应路径上image的受力和每个原子的受力是有差别的。最小能量路径上每个image的原子受力有可能就是较大的，因为image不是真的能量最低点，所以原子受力有可能是较大的（比如大于0.05eV/A）.只要反应路径平行方向的弹簧力和反应路径垂直方向的势能力小于设定标准，就应该是收敛了。
如果一定要解释用image的原子最大受力做收敛标准的合理性，即EDIFFG=-0.05表示$(grep "FORCES: max atom"\) <= 0.05的话，那可能是因为每个原子的受力小到一定程度，自然对应的弹簧力和势能力投影就小了。但是弹簧力是虚拟的力，似乎并不会体现在image的原子受力上。原子受力应该仅仅是实际的势能力。
I have check the code nebbarrier.pl, where I can find the $force='grep "NEB: projections" $directories/OUTCAR | tail -1', and then
@force=split(/\s+/,$force);
$force=$force;
In one of my OUTCAR, the line "NEB: projections" read like : NEB: projections on to tangent (spring, REAL) -0.425069 0.137351
So, the last $force=0.137351.
As results, the nebef.pl and nebbarrier.pl provide different forces. The former corresponds to "FORCES: max atom", and the latter to REAL force projected to the tangent of the reacton path, excluding the spring forces.
I am confused! And, I do not know what kind of forces, or what components of forces push the images and atoms in the images down to the MEP.
http://edu.emuch.net/attachment/b7/f6/537452_1303822905_728.jpg
图一 OUTCAR中的两个力项，红色划线部分分别对应nebef.pl和nebbarrier.pl的力。
nebef.pl ------------> FORCES: max atom; nebbarrier.pl ----------> NEB: projections
这里似乎是说，nebbarrier.pl抽取的力是REAL force在反应路径切线方向分量，不包括一丁点儿spring force. 是这样的吗？
http://edu.emuch.net/attachment/6f/dc/537452_1303821850_642.jpg
图二 nebef.pl 脚本中关于力的部分。 见$force.
http://edu.emuch.net/attachment/01/dc/537452_1303821860_502.jpg
图三nebbarrier.pl 脚本中关于力的部分。 见$force.
http://edu.emuch.net/attachment/2b/58/537452_1303833508_382.jpg
nebef.pl是弹簧力沿反应路径切线方向分量+势能力垂直切线方向分量；
nebbarrier.pl是势能力沿反应路径切线方向分量。
，小木虫的字体有点问题。
这两篇文章就是讲NEB具体怎么做的。引用都高得离谱。
Last edited by cenwanglai on 2011-4-26 at 09:11 ]
nebbarrier.pl是势能力沿反应路径切线方向分量。
http://edu.emuch.net/attachment/2b/58/537452_1303833508_382.jpg
Last edited by cenwanglai on 2011-4-27 at 00:04 ]
tangent的确定依赖于前后两个image与当前image的矢量（R(i+1)-R(i)）和(R(i)-R(i-1))。如果image比较少，即image之间的距离与势能面的起伏尺度相比，前者较大的话，所得切不准确，由此会带来计算收敛的困难吧？
因此，如果反应路径如果收敛困难，某些image之间距离较大，且计算过程中受力起伏变化较大的话，是不是可以通过适当增加image来加速收敛呢？
Graeme他们似乎推荐image之间的距离不超过1A. 我的计算资源有限，给了四个image,image之间的距离有2~3A的。这样如果收敛得到一条路径和能垒，这个能垒是否可信呢？
如果你要准确的描述一个单峰势垒过程，算头尾至少要5个点（由傅里叶变换得到的描述一个单峰周期的基本精度要求）。多峰过程依次类推。另外NEB方法处理多峰过程个人感觉有问题，因为NEB方法会义无反顾的跳进稳态的谷中，个人感觉这个和实际的物理过程不相符。不过还好，不影响势垒计算。
至于你说的image数目的选取，至少要保证上面数目。如果条件允许，最好能够做到平均每个原子的移动小于一个原子的大小即0.1-0.2A，不过一般做不到。
另外NEB有两点要注意，初始路径要仔细的设计。保证是你想要的初猜路径。还有不要一开始就是用CINEB否则初猜路径会被破坏，当路径稳定之后，再使用。切忌，初猜路径非常重要，是决定你能否收敛和收敛之后结果是否可信的关键！
另外如果出现N多的峰或者收敛趋势很不明显时，我建议你仔细检查一下初始路径！
（1）什么是反应坐标（reaction coordinates）？
从包含N个原子的反应物到产物的过程记为：A->B
我们只知道A和B的所有原子坐标，插入若干中间点后得到一条反应路径。这条路径上，若不考虑原子的固定，其内坐标变量为3N-6.我们常常以一条线来表示反应路径，这条线是不是3N-6维空间中的线呢？如果再加上能量维度，这条反应路径是不是3N-5维空间中的曲线？
（2）在(C)NEB计算中，各image是怎样移动的？
我的初浅理解，是各个image在垂直反应路径方向受势能力而移动，这个与结构有过程类似：小小的移动一个离子步，然后求能量对坐标的变化率，得到image在垂直反应路径的势能力。而在反应路径方向，各个image之间受弹簧力作用。
一般三个紧邻的image之间成一定的（0~180）夹角，中间的image受到前后两个弹簧的合力。优化过程为把image受到的弹簧力分解，保留image在反应路径上切线方向的力。
所以，每个image有两个运动：受势能力驱动的，垂直于反应路径上的运动；和受弹簧力驱动的，沿着反应路径切线方向的运动。
优化最终的结果，应该是每个image受到的势能力和弹簧力小于设定标准。那么，具体到vasp的NEB计算nebbarrier.pl的输出第四列，说是沿着反应路径方向（从反应物指向产物的有向曲线）的受力。这个力为正，表示沿着反应路径想产物移动能量增加，否则相反。
这里提到的力，包括势能力(A1)和弹簧力(B1)，加上二者在反应路径切线平行方向(A11,B11)和垂直方向(A12,B12)的投影，力的收敛标准（EDIFFG）到底是以什么为对象呢？
（3）image之间的夹角
三个相邻的image成一个角，这个角的角度如何定义？如果三个image依次命名为i-1,i和i+1，对应坐标为R(i-1),R(i)和R(i+1)，角度是不是矢量到与矢量的方向一致所需转动的最小角度（0-180度）？
（4）image之间的距离
三维空间两点距离是对应坐标项差的平方和的开方。对于原子个数众多且变化的体系，为了距离的可比较，似乎求的是二阶矩：对应原子坐标项的差的平方和的开方对原子个数的平均。
（5)image的受力和收敛力判据到底指什么？
下面9楼和10楼贴了两个跟力相关的帖子。nebef.pl的屏幕输出得到的力是每个image中的FORCES: max atom.难道收敛判据是每个image中原子最大受力？为什么不是每个image所受的弹簧力沿切线方向分量和势能力沿垂直切线方向分量？或者后两者的合力？
我觉得，反应路径上image的受力和每个原子的受力是有差别的。最小能量路径上每个image的原子受力有可能就是较大的，因为image不是真的能量最低点，所以原子受力有可能是较大的（比如大于0.05eV/A）.只要反应路径平行方向的弹簧力和反应路径垂直方向的势能力小于设定标准，就应该是收敛了。
如果一定要解释用image的原子最大受力做收敛标准的合理性，即EDIFFG=-0.05表示$(grep "FORCES: max atom"\) <= 0.05的话，那可能是因为每个原子的受力小到一定程度，自然对应的弹簧力和势能力投影就小了。但是弹簧力是虚拟的力，似乎并不会体现在image的原子受力上。原子受力应该仅仅是实际的势能力。
I have check the code nebbarrier.pl, where I can find the $force='grep "NEB: projections" $directories/OUTCAR | tail -1', and then
@force=split(/\s+/,$force);
$force=$force;
In one of my OUTCAR, the line "NEB: projections" read like : NEB: projections on to tangent (spring, REAL) -0.425069 0.137351
So, the last $force=0.137351.
As results, the nebef.pl and nebbarrier.pl provide different forces. The former corresponds to "FORCES: max atom", and the latter to REAL force projected to the tangent of the reacton path, excluding the spring forces.
I am confused! And, I do not know what kind of forces, or what components of forces push the images and atoms in the images down to the MEP.
http://edu.emuch.net/attachment/b7/f6/537452_1303822905_728.jpg
图一 OUTCAR中的两个力项，红色划线部分分别对应nebef.pl和nebbarrier.pl的力。
nebef.pl ------------> FORCES: max atom; nebbarrier.pl ----------> NEB: projections
这里似乎是说，nebbarrier.pl抽取的力是REAL force在反应路径切线方向分量，不包括一丁点儿spring force. 是这样的吗？
http://edu.emuch.net/attachment/6f/dc/537452_1303821850_642.jpg
图二 nebef.pl 脚本中关于力的部分。 见$force.
http://edu.emuch.net/attachment/01/dc/537452_1303821860_502.jpg
图三nebbarrier.pl 脚本中关于力的部分。 见$force.
http://edu.emuch.net/attachment/2b/58/537452_1303833508_382.jpg
nebef.pl是弹簧力沿反应路径切线方向分量+势能力垂直切线方向分量；
nebbarrier.pl是势能力沿反应路径切线方向分量。
G. Henkelman and H. Jo ´nsson, J. Chem. Phys. 113, 9978
(2000)
G. Henkelman, B. P. Uberuaga and H. Jónsson, J. Chem. Phys. 113, 9901 (2000)
做NEB怎么能不看这两篇开创性的文章呢！
G. Henkelman, B. P. Uberuaga and H. Jónsson, J. Chem. Phys. 113, 9901 (2000)
做NEB怎么能不看这两篇开创性的文章呢！
Originally posted by dxcharlary at 2011-04-25 2101:
G. Henkelman and H. Jo ´nsson, J. Chem. Phys. 113, 9978 (2000)
G. Henkelman, B. P. Uberuaga and H. Jónsson, J. Chem. Phys. 113, 9901 (2000)
做NEB怎么能不看这两篇开创性的文章呢！
第一篇文章是G. Henkelman and H. Jónsson,J. Chem. Phys. 113, 9978
(2000) G. Henkelman and H. Jo ´nsson, J. Chem. Phys. 113, 9978 (2000)
G. Henkelman, B. P. Uberuaga and H. Jónsson, J. Chem. Phys. 113, 9901 (2000)
做NEB怎么能不看这两篇开创性的文章呢！
这两篇文章就是讲NEB具体怎么做的。引用都高得离谱。
Last edited by cenwanglai on 2011-4-26 at 09:11 ]
Originally posted by dxcharlary at 2011-04-25 2107:
第一篇文章是G. Henkelman and H. Jónsson,J. Chem. Phys. 113, 9978 (2000) ，小木虫的字体有点问题。
这两篇文章就是讲NEB具体怎么做的。引用都高得离谱。
Last edited by cenwanglai on 2011-4 ...
请教：ó这种用什么字体输入呢？第一篇文章是G. Henkelman and H. Jónsson,J. Chem. Phys. 113, 9978 (2000)
这两篇文章就是讲NEB具体怎么做的。引用都高得离谱。
Last edited by cenwanglai on 2011-4 ...
Originally posted by cenwanglai at 2011-04-26 0902:
请教：ó这种用什么字体输入呢？
我也不知道，直接粘的。请教：ó这种用什么字体输入呢？
Originally posted by dxcharlary at 2011-04-26 1326:
我也不知道，直接粘的。
提供的文献确实讲述的很清楚。早上一口气读完两篇。爽快！我也不知道，直接粘的。
在vasp官网上发现这个帖子。他使用vasp自带的NEB来算的。这个OUTCAR似乎与编译了VSTSTCODE的CINEB计算输出OUTCAR在力上的表述有很大差别。
路过的有可以提供解答的不？
how to find force acting on an image
Dear all,
I use regular NEB method to investigate atom diffusion process. I encounter a convergence problem. I read some papers on NEB methed and knew that the total force acting on an image is the sum of the spring force along the local tangent and the true fore perpendicular to the local tangent. I also note that the OUTCAR file gives
"energy of chain is (eV)",
"tangential force (eV/A)',
"CHAIN-FORCE (eV/Angst)",
"CHAIN+TOTAL (eV/Angst)",
"g(Force)" and "g(Stress)".
What is total force acting on an image (How to find this total force in OUTCAR)? If I set EDIFFG=-0.05 for NEB calculation, does it means that the chain convergence criteria is the fore actiong on an image are less than 0.05 eV/A?
PS: I use regular NEB implanted in VASP 4.6.
Last edited by cenwanglai on 2011-4-26 at 15:00 ]
路过的有可以提供解答的不？
how to find force acting on an image
Dear all,
I use regular NEB method to investigate atom diffusion process. I encounter a convergence problem. I read some papers on NEB methed and knew that the total force acting on an image is the sum of the spring force along the local tangent and the true fore perpendicular to the local tangent. I also note that the OUTCAR file gives
"energy of chain is (eV)",
"tangential force (eV/A)',
"CHAIN-FORCE (eV/Angst)",
"CHAIN+TOTAL (eV/Angst)",
"g(Force)" and "g(Stress)".
What is total force acting on an image (How to find this total force in OUTCAR)? If I set EDIFFG=-0.05 for NEB calculation, does it means that the chain convergence criteria is the fore actiong on an image are less than 0.05 eV/A?
PS: I use regular NEB implanted in VASP 4.6.
Last edited by cenwanglai on 2011-4-26 at 15:00 ]
http://edu.emuch.net/attachment/d4/36/537452_1303801343_182.jpg
http://edu.emuch.net/attachment/4a/9c/537452_1303801350_304.jpg
上面两个图是编译了vststcode计算CNEB的OUTCAR中的截图。
I have read G. Henkelman's two papers on CNEB publised in 2000, and several posts in this forum. I am still confusing on what the forces convergent criterion standing for.
In the screen output when running the nebef.pl, the second column seem to be the forces from /grep "max atom" OUTCAR/.
But, in my opinion, the forces shoud be the spring force parallel to th the tangent and the real force perpendicular to the tangent, or the sum of the two parts, which is added to the image, rather than any single atom.
Can any body give me some advice?
Last edited by cenwanglai on 2011-4-26 at 15:09 ]
http://edu.emuch.net/attachment/4a/9c/537452_1303801350_304.jpg
上面两个图是编译了vststcode计算CNEB的OUTCAR中的截图。
I have read G. Henkelman's two papers on CNEB publised in 2000, and several posts in this forum. I am still confusing on what the forces convergent criterion standing for.
In the screen output when running the nebef.pl, the second column seem to be the forces from /grep "max atom" OUTCAR/.
But, in my opinion, the forces shoud be the spring force parallel to th the tangent and the real force perpendicular to the tangent, or the sum of the two parts, which is added to the image, rather than any single atom.
Can any body give me some advice?
Last edited by cenwanglai on 2011-4-26 at 15:09 ]
麻烦路过的近来讨论一下吧！！！！
vasp有自带的NEB吗？没用过，不大清楚。
作用在单个原子上面和作用在整个image上面没什么太大的区别。因为所有原子坐标组成一个hilbert空间。所以整个image的受力可以看成是每个原子的受力的直和。那个原子的受力应该是计算了弹簧力之后的受力吧。这个以前没太注意过，你可以单独算一下那个image，然后比较一下受力。
收敛标准一定应该是合力收敛性，如果只是势场力收敛的话，那不是都跑到稳态附近或是鞍点附近去了！
作用在单个原子上面和作用在整个image上面没什么太大的区别。因为所有原子坐标组成一个hilbert空间。所以整个image的受力可以看成是每个原子的受力的直和。那个原子的受力应该是计算了弹簧力之后的受力吧。这个以前没太注意过，你可以单独算一下那个image，然后比较一下受力。
收敛标准一定应该是合力收敛性，如果只是势场力收敛的话，那不是都跑到稳态附近或是鞍点附近去了！
http://theory.cm.utexa ... ools/scripts/
usage: nebef.pl
output: energy and force of images in the neb
usage: nebbarrier.pl
output: energy, distance, and forces along the neb
注意这里面的措辞。
nebef给出的应该是势场力，而nebbarrier.pl给出的是合力，是判断收敛的依据。
至于是原子的还是image整体的，意义不大。
usage: nebef.pl
output: energy and force of images in the neb
usage: nebbarrier.pl
output: energy, distance, and forces along the neb
注意这里面的措辞。
nebef给出的应该是势场力，而nebbarrier.pl给出的是合力，是判断收敛的依据。
至于是原子的还是image整体的，意义不大。
令外我现在用的程序是我和我师兄自己写的。austin的那个程序我也只用过那么一两次，你问得太具体的我也不清楚。
Originally posted by dxcharlary at 2011-04-26 2237:
http://theory.cm.utexa ... ools/scripts/
usage: nebef.pl
output: energy and force of images in the neb
usage: nebbarrier.pl
output: energy, distance, and forces along the neb ...
nebef.pl是弹簧力沿反应路径切线方向分量+势能力垂直切线方向分量；http://theory.cm.utexa ... ools/scripts/
usage: nebef.pl
output: energy and force of images in the neb
usage: nebbarrier.pl
output: energy, distance, and forces along the neb ...
nebbarrier.pl是势能力沿反应路径切线方向分量。
http://edu.emuch.net/attachment/2b/58/537452_1303833508_382.jpg
Last edited by cenwanglai on 2011-4-27 at 00:04 ]
Originally posted by dxcharlary at 2011-04-26 2311:
令外我现在用的程序是我和我师兄自己写的。austin的那个程序我也只用过那么一两次，你问得太具体的我也不清楚。
dxcharlary,
CINEB来计算反应路径的方法中，对local 令外我现在用的程序是我和我师兄自己写的。austin的那个程序我也只用过那么一两次，你问得太具体的我也不清楚。
因此，如果反应路径如果收敛困难，某些image之间距离较大，且计算过程中受力起伏变化较大的话，是不是可以通过适当增加image来加速收敛呢？
Graeme他们似乎推荐image之间的距离不超过1A. 我的计算资源有限，给了四个image,image之间的距离有2~3A的。这样如果收敛得到一条路径和能垒，这个能垒是否可信呢？
Originally posted by cenwanglai at 2011-04-27 1711:
dxcharlary, CINEB来计算反应路径的方法中，对local tangent的确定依赖于前后两个image与当前image的矢量（R(i+1)-R(i)）和(R(i)-R(i-1))。如果image比较少，即image之间的距离与势能面的起伏尺度相比，前
...
image之间的距离与势能面的起伏尺度相比，前者较大的话，如果你运气好的话，经过很多步之后，有可能收敛，我遇到过这样的情况，但是结果根本就不对。dxcharlary, CINEB来计算反应路径的方法中，对local
如果你要准确的描述一个单峰势垒过程，算头尾至少要5个点（由傅里叶变换得到的描述一个单峰周期的基本精度要求）。多峰过程依次类推。另外NEB方法处理多峰过程个人感觉有问题，因为NEB方法会义无反顾的跳进稳态的谷中，个人感觉这个和实际的物理过程不相符。不过还好，不影响势垒计算。
至于你说的image数目的选取，至少要保证上面数目。如果条件允许，最好能够做到平均每个原子的移动小于一个原子的大小即0.1-0.2A，不过一般做不到。
另外NEB有两点要注意，初始路径要仔细的设计。保证是你想要的初猜路径。还有不要一开始就是用CINEB否则初猜路径会被破坏，当路径稳定之后，再使用。切忌，初猜路径非常重要，是决定你能否收敛和收敛之后结果是否可信的关键！
另外如果出现N多的峰或者收敛趋势很不明显时，我建议你仔细检查一下初始路径！
16楼: Originally posted by cenwanglai at
2011-04-27 0032:
nebef.pl是弹簧力沿反应路径切线方向分量+势能力垂直切线方向分量；
nebbarrier.pl是势能力沿反应路径切线方向分量。
http://edu.emuch.net/attachment/2b/58/537452_1303833508_382.jpg
Last edited by cenwanglai on 20 ...
请教一下，由nebspline.pl得到的spline.dat中第四列数据是什么含义？怎么拟合的?nebef.pl是弹簧力沿反应路径切线方向分量+势能力垂直切线方向分量；
nebbarrier.pl是势能力沿反应路径切线方向分量。
http://edu.emuch.net/attachment/2b/58/537452_1303833508_382.jpg
Last edited by cenwanglai on 20 ...
喜欢
0
赠金笔
前一篇：分三段