The nudged elastic band (NEB) is a method for finding saddle points and minimum energy paths between known reactants and products. The method works by optimizing a number of intermediate images along the reaction path. Each image finds the lowest energy possible while maintaining equal spacing to neighboring images. This constrained optimization is done by adding spring forces along the band between images and by projecting out the component of the force due to the potential perpendicular to the band. Difference from the implementation in VASPThere are a few improvements to the NEB method which are not yet included in the current version of vasp. A new tangent definition and a climbing image method combine to allow for the more accurate finding of saddle points using the NEB with fewer images than the original method. The setup and operation of this implementation can be identical to what is described in the vasp manual under the elastic band section. The new tangent is implemented by default, and the climbing image method can be turned by setting LCLIMB=.TRUE. in the INCAR file. Climbing Image The climbing image is a small modification to the NEB method in which the highest energy image is driven up to the saddle point. This image does not feel the spring forces along the band. Instead, the true force at this image along the tangent is inverted. In this way, the image tries to maximize it's energy along the band, and minimize in all other directions. When this image converges, it will be at the exact saddle point. Because the highest image is moved to the saddle point and it does not feel the spring forces, the spacing of images on either side of this image will be different. It can be important to do some minimization with the regular NEB method before this flag is turned on, both to have a good estimate of the reaction co-ordinate around the saddle point, and so that the highest image is close to the saddle point. If the maximum image is initially very far from the saddle point, and the climbing image was used from the outset, the path would develop very different spacing on either side of the saddle point. To use the climbing image, set LCLIMB=.True. Example The graph on the right shows an NEB calculation (blue) and a climbing image cNEB calculation* (red). The system is an Al adatom on an Al(100) surface. The process is an exchange between the adatom and a substraight atom, leading to adatom diffusion. Notice how the climbing image calculation has shifted the position of the images (by compressing the images on the left) so that one image sits right at the saddle point. *The cNEB energies have been shifted by 0.05 eV so that the two curves are distinct. Standard Parameters Parameter | Default Value | Description | ICHAIN | 0 | Indicates which method to run. NEB (ICHAIN=0) is the default | IMAGES | | Number of images in the band, excluding endpoints | SPRING | -5.0 eV/A2 | Spring force between images; negative value turns on nudging | LCLIMB | .TRUE. | Flag to turn on the climbing image algorithm | LTANGENTOLD | .FALSE. | Flag to turn on the old central difference tangent | LDNEB | .FALSE. | Flag to turn on modified doubble nudging | LNEBCELL | .FALSE. | Flag to turn on SS-NEB. Used with ISIF=3 and IOPT=3. | JACOBIAN | (Ω/N) 1/3N1/2 | Controls weight of lattice to atomic motion. Ω is volume and N is the number of atoms. |
FORCE BASED OPTIMIZERS
| | The NEB and min-mode following (dimer/Lanczos) saddle point finding methods use a force projection in order to direct the optimizers towards minimum energy paths and saddle points. These modifications to the force mean that the energy is no longer consistent with the force being optimized. Because of this, we can only use optimizers that are solely based upon the force (and not the energy). The quasi-Newton and quick-min (IBRION=1 and 3 respectively) optimizers that are built into VASP are both force-based, but the conjugate-gradient method (IBRION=2) is not. Here, we present a set of optimizers that are all force-based so they can be used with the NEB and min-mode following methods. To use them, the INCAR must set IBRION=3 and POTIM=0, to disable the built in optimizers. Then, the IOPT parameter will select one of the following methods. This version of quick-min is essential the same as what has been implemented in vasp. The conjugate-gradient method is different in that is uses a Newton's line optimizer, and are entirely force based. The LBFGS is also different in that the NEB can be optimized globally, instead of image-by-image. The FIRE optimizer is an interesting new optimizer which has similarities to quick-min, but tends to be faster. The steepest descent method is provided primarily for testing. We recommend using CG or LBFGS when accurate forces are available. This is essential for evaluating curvatures. For high forces (far from the minimum) or inaccurate forces (close to the minimum) the quick-min or FIRE methods are recommended. These two methods do not rely on curvatures, and tend to be less aggressive, better behaved, but also less efficient than CG/LBFGS. Here is a recent paper discussing the performance of these different optimizers with the NEB. The min-mode following methods can only be used with one of these methods. If IOPT is not set in a dimer run, the job with die with a corresponding error message. Optimizer input parameters The following parameters are read from the INCAR file. (IOPT = 0) Use VASP optimizers specified from IBRION (default) (IOPT = 1) LBFGS = Limited-memory Broyden-Fletcher-Goldfarb-Shanno (IOPT = 2) CG = Conjugate Gradient (IOPT = 3) QM = Quick-Min (IOPT = 4) SD = Steepest Descent (IOPT = 5) BFGS = Full Broyden-Fletcher-Goldfarb-Shanno (IOPT = 6) MD = Molecular dynamics based on Verlet-Velocity (IOPT = 7) FIRE = Fast Inertial Relaxation Engine Required ParametersParameter | Default | Description | IOPT | 0 | Default is to use the VASP optimizers | IBRION | 3 | Specify that VASP do molecular dynamics (with zero time step) | POTIM | 0 | Zero time step so that VASP does not move the ions | LBFGS Parameters (IOPT = 1)Parameter | Default | Description | MAXMOVE | 0.2 | Maximum allowed step size for translation | ILBFGSMEM | 20 | Number of steps saved when building the inverse Hessian matrix | LGLOBAL | .TRUE. | Optimize the NEB globally instead of image-by-image | INVCURV | 0.01 | Initial inverse curvature, used to construct the inverse Hessian matrix | LLINEOPT | .FALSE. | Use a force based line minimizer for translation | FDSTEP | 5E-3 | Finite difference step size for line optimizer | CG Parameters (IOPT = 2)Parameter | Default | Description | MAXMOVE | 0.2 | Maximum allowed step size for translation | FDSTEP | 5E-3 | Finite difference step size to calculate curvature | QM Parameters (IOPT = 3)Parameter | Default | Description | MAXMOVE | 0.2 | Maximum allowed step size for translation | TIMESTEP | 0.1 | Dynamical time step | SD Parameters (IOPT = 4)Parameter | Default | Description | MAXMOVE | 0.2 | Maximum allowed step size for translation | SDALPHA | 0.1 | Ratio between force and step size | BFGS Parameters (IOPT = 5)Parameter | Default | Description | MAXMOVE | 0.2 | Maximum allowed step size for translation | FDSTEP | 5E-3 | Finite difference step size to calculate curvature | LINEMIN | .TRUE. | Use force based line minimizer for direction | BFGSINVCURV | 0.0001 | Scales Hessian for initial move | BFGSDFP | .FALSE. | Use BFGS updating formula to update the Hessian | MD Parameters (IOPT = 6)Parameter | Default | Description | MAXMOVE | 0.2 | Maximum allowed step size for translation | TIMESTEP | 0.1 | Dynamical time step | MDTHERMO | 0 | Thermostat used for MD 0=Anderson, 1=Langevin | MDTEMP | 300 | Temperature in Kelvin for MD | MDGAMMA | 0.1 | Friction for Langevin dynamics | MDALPHA | 0.5 | Controls how aggressively velocity is reset | MDCOLLISION | dt | For Anderson thermostat, the probability of collision is 1-exp(-dt/MDCOLLISION) | FIRE Parameters (IOPT = 7)Parameter | Default | Description | MAXMOVE | 0.2 | Maximum allowed step size for translation | TIMESTEP | 0.1 | Initial dynamical time step | FTIMEMAX | 1.0 | Maximum dynamical time step allowed | FTIMEDEC | 0.5 | Factor to decrease dt | FTIMEINC | 1.1 | Factor to increase dt | FALPHA | 0.1 | Parameter that controls velocity damping | FALPHADEC | 0.99 | Factor to decrease alpha | FNMIN | 5 | Minimum number of iterations before adjusting alpha and dt |
|
计算过渡态先要摆正心态,不急于下手。步骤如下: (1)做模型,初态IS和终态FS,分别结构优化到基态; (2)线形插入images: nebmake.pl POSCAR.IS POSCAR.FS N N为image个数。 (3)nebmovie.pl,生成movie.xyz。用Xcrysden --xyz movie.xyz 反复观看动画,仔细检查过程的合理性。这里要提醒,POSCAR.IS 和POSCAR.FS中原子坐标列表的顺序必须对应。 (4)写INCAR,选IOPT。注意,最好忘记vasp自带的NEB,而全部改用包含vtstool的vasp. IBRION=3,POTIM=0关闭vasp自带的NEB功能。 (5)过渡态计算第一个离子步最耗时,也最容易出问题,也是模型设计合理性检验的首要环节。所以可以选小一些的ENCUT,可以不用考虑自旋(ISPIN=1),也不用考虑DFT+U。而且用最快最粗糙的算法(IOPT=3,其他默认)。 (6)带vtstool的vasp-ClNEB(NEB)过渡态计算ICHAIN=0作为入口,这个也是默认的。LCLIMB=TRUE也是默认的。如果不要climb image,可以设置LCLIMB = False. (7)收敛判据EDIFFG<0。过渡态计算要以力为收敛判据,而不是能量。一般EDIFFG=-0.05就可以接受,-0.02或者-0.01更好。但是作为开始的过渡态计算,可以设置很宽的收敛条件,如EDIFFG=-1. (8)初步过渡态收敛后,修改INCAR中的优化器(IOPT),并修改相应参数(参考vtstool官方论坛),EDIFFG改小(如-0.05),然后运行vfin.pl,这个脚本自动帮你准备在原来的基础上继续运行新的过渡态计算(完成cp CONTCAR POSCAR, 保留电荷密度和波函数的操作)。 (9)过渡态如何验算虚频呢? | | | FORCES: max atom, RMS |
|