Transition State search using CI-NEB by VASP-VTST
- Xijun Wang
- Jun 1, 2023
- 5 min read
The transition state (TS) of a chemical reaction is a particular configuration along the reaction coordinate. It is defined as the state corresponding to the highest potential energy along this reaction coordinate. Left figure is an illustration of the potential energy surface (PES) of a system, where the two local minimas (reactant and product) are connected by a red-marked reaction coordinate. The maximum along such a path is called a TS, which is a first order saddle point, a maximum in one coordinate and minima in all others.
Climbing Image Nudged Elastic Band (CI-NEB)
The nudged elastic band (NEB) is a widely-used method for searching TS on the PES. The method works by optimizing a number of intermediate images along the reaction path connecting reactants and products. 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.
The climbing image is an advanced version of the NEB method in which the highest energy image is driven up to the saddle point. In stead of adding spring forces, the inverted true forces along the tangent is added to this image. Such that the image's energy can be maximized along the band, and minimized in all other directions. With this approach, the exact saddle point (TS) can be found when this image converges.
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.
Installation of VTST codes into VASP I would like to skip this step because: 1. VTST is so popular that it can be easily loaded on most clusters by some simple module commands like: module load vasp/5.3.3.vtst. 2. Very detialed step-by-step tutorials has been provied by the official site of VTST: https://theory.cm.utexas.edu/vtsttools/installation.html
Example: NH3 activation on the surface of V2O5. Using CI-NEB implemented by VASP-VTST, the target of this example is to search the involved TS and calculate the energy required for the NH3 activation (*NH3 → *NH2 + *H), in which a proton transfers from the N site on NH3 to an oxygen site on the surface of V2O5.
Step 1. Prepartion Download all the wide-used scripts from the office-site of VTST: https://theory.cm.utexas.edu/vtsttools/scripts.html
tar -xvf vtstscripts.tgz
cp * ~/bin
chmod +x ~/bin/*Step 2. IS and FS Optimization Optimize the initial structure (IS, *NH3) and final structure (FS, *NH2 + *H) using VASP respectively prior to TS search. Copy the generated CONTCAR of IS and FS to a new directory as IS.vasp and FS.vasp, respectively.
cp IS/CONTCAR /your_directory/IS.vasp
cp FS/CONTCAR /your_directory/FS.vaspStep 3. Check the serial number of atoms in IS.vasp and FS.vasp It is required that the serial number of atoms in IS.vasp and FS.vasp are the same. A script ~/bin/dist.pl can be used to check how similar the IS.vasp and the FS.vasp are, which generates the average distance of the same atom in IS.vasp and FS.vasp. In most cases the output value is smaller than 5 Å . In addition, it is required to replace the first line of IS.vasp with the atomic type in the order of POTCAR as shown below.
[xwang99@login04 test]$ dist.pl IS.vasp FS.vasp
2.24955063206987
[xwang99@login04 test]$ head IS.vasp
V O N H # The first line of IS.vasp must be the atomic type in the order of POTCAR
1.00000000000000
11.5241003035999992 0.0000000000000000 0.0000000000000000
0.0000000000000000 12.6374998092999995 0.0000000000000000
0.0000000000000000 0.0000000000000000 20.0000000000000000
V O N H
18 45 1 3
Selective dynamics
Direct
0.1495900019999965 0.0364399999999989 0.0001299999999986 F F FStep 4. Interpolate images between IS and FS This step can be easily achieved by using ~/bin/nebmake.pl IS.vasp FS.vasp N, where N is the number of intermediate images as shown below. The generated sub-directories 00 contains the POSCAR of the initial geometry and 0(N+1) contains the POSCAR of the final geometry. The sub-directories 01 - 0N contains the POSCAR of intermediate images. After generating the sub-directories 01/, 02/, 03/, ... check the validity of these intermediate images (POSCAR) by visualizing them in VESTA.
[xwang99@login04 test]$ nebmake.pl IS.vasp FS.vasp 3 # Here 3 is the number of images, which can be replaced with any number you want.
filetype1: vasp5
filetype2: vasp5
OK, ALL SETUP HERE
FOR LATER ANALYSIS, PUT OUTCARs IN FOLDERS 00 and 04 !!!Step 5. Prepare INCAR, POTCAR, KPOINTS In INCAR, it is required to add/modify the related keywords of CI-NEB, including but not limitted to: ------------------------------------------------------------------------------------------------------------------------- #NEB LCLIMB=.True. ICHAIN=0 # Indicates which method to run. NEB (ICHAIN=0) is the default IOPT=0 # = 1, 2, 3, 4, 7 ... for different optimizers IBRION=1 # If you want to use IOPT = 1, 2, 3 ..., set IBRION=3 POTIM=0.2 # Change the step size based on your systems IMAGES = 3 # The number of images SPRING = -5 # The spring constant, in eV/Ang^2 between the images; negative value turns on nudging EDIFFG=-0.02 # Convergence criterion of force NSW=500 # The maximum steps allowed ISIF=2 # Fix the lattice constants or not ------------------------------------------------------------------------------------------------------------------------- POTCAR and KPOINTS are the same as structural optimization of IS and FS. Also, a script is required to submit the job. After preparing all the following files, submit the job.
[xwang99@login04 test]$ ls
00 01 02 03 04 FS.vasp INCAR IS.vasp KPOINTS POTCAR script-vasp5.4-NEBStep 6. Check the output files After the job ends normally, copy the IS/OUTCAR and FS/OUTCAR into /your_directory/00/ and 04/, respectively, which could make it easier to obtain the energy of IS, TS, and FS. Subsequently, double-check the validity of the obtained images (CONTCAR) in 01/, 02/, and 03/ by visualizing them in VESTA. If everything looks good, use nebresults.pl to generate all the information about this job. As an example, the highlighted information are the energies of the IS.vasp, interpolated images. and FS.vasp, in which the one with the highest energy is the TS.
[xwang99@login04 test]$ nebresults.pl
Unziping the OUTCARs ... done
Do nebbarrier.pl ; nebspline.pl
Do nebef.pl
Do nebmovie.pl
Do nebjmovie.pl
Do nebconverge.pl
Warning: empty x range [0:0], adjusting to [-1:1]
Warning: empty y range [0.035996:0.035996], adjusting to [0.035636:0.036356]
Warning: empty y2 range [-463.462:-463.462], adjusting to [-468.097:-458.828]
Warning: empty x range [0:0], adjusting to [-1:1]
Warning: empty y range [0.049215:0.049215], adjusting to [0.0487229:0.0497071]
Warning: empty y2 range [-462.495:-462.495], adjusting to [-467.12:-457.87]
Warning: empty x range [0:0], adjusting to [-1:1]
Warning: empty y range [0.040008:0.040008], adjusting to [0.0396079:0.0404081]
Warning: empty y2 range [-462.893:-462.893], adjusting to [-467.522:-458.264]
Zipping the OUTCARs again ... done
Forces and Energy:
0 0.000000 -463.541800 0.000000
1 0.015996 -463.462200 0.079600
2 0.019215 -462.495000 1.046800
3 0.010008 -462.893300 0.648500
4 0.000000 -462.999100 0.542700
Extremum 1 found at image 1.999395 with energy: 1.046749
Extremum 2 found at image 3.000000 with energy: 0.648451
Extremum 3 found at image 3.666667 with energy: 159015.437564Some Empirical Tips for TS searching
In most cases, the automatically interpolated images are not as good as the ones you built based on your chemical intuition. Feel free to adjust the interpolated images.
I would suggest to try IOPT=0 IBRION=1 first. If it doesn't work properly, try IOPT=3 IBRION=3 and then try IOPT=1, 2, 4, 7
In some cases, it would be efficient to use IOPT=3 IBRION=3 to find a rough TS configuration, and then switch to IOPT=0 IBRION=1 to find the exact TS. Remenber to cp CONTCAR POSCAR in each interpolated image directory 01/, 02/, ...
Keep monitoring the status of the jobs.
If you found this tutorial helpful, consider joining our community by making a donation. Your support will enable me to continue creating valuable content while also supporting my baby's care. Together, we can make a meaningful impact. Thank you!
Comments