Geometry Optimization

For geometry optimizations, all settings from the Description of Input Files page can be used. Besides them, some others must be added.

Unfortunately, a geometry optimization with VASP can be challenging, especially if large systems periodic with many smooth degrees of freedom (e.g., long alkyl chains adsorbed on a meta surface) are treated. Therefore, a separate Tricks section has been added after the Settings session.

Settings

The following list gives a keywords that should be added to an INCAR file besides those in the Description of Input Files for a geometry optimization.

• IBRION = 1 or 2 If 1 is set, the quasi-Newton optimzation algirithm is chosen. It is usually the fastest one to converge, but requires a resonable starting structure. Else, explosions might occur. If IBRION = 2 is set, the simpler conjugate-gradient algorithm is chosen. It can be used for bad starting structures. If the system get stucked and stagnates for one of the algorithms, simply try the other if it improves the situation.

• ISIF = [number] How the atomic and lattice degrees of freedom shall be treated. For a system in a rigid unit cell, where only atomic positions shall be optimized, ISIF = 0 is the choice. If only the cell shall be optimized, choose either ISIF = 6 (shape and volume changed) or ISIF = 7 (volume changed, shape kept). If the cell shape/volume and the atomic positions shall be optimized, use ISIF = 3. The new option ISIF = 8 keeps the cell shape constant and allows only the volume to relax.

• EDIFFG = [value] The criterion for geometric convergence. A value often used in publications is 0.02 (EDIFFG = -0.02) or 0.01 eV/Å (EDIFFG = -0.01) for the maximum gradient component, such that all components of the actual gradient for the atoms activated for optimization (T T T with selective dynamcs) must be below this value (with a minus sign!). If the calculation does not converge and stucks, e.g., at 0.028, criterion of 0.03 eV/Å os also acceptable and reported in literature. Alternatively, the change in total energy can be used as convergence criterion (now without a minus sign). Then, convergence is reached if the total energy between two geometry optimization steps is lower than this value. The value should not be larger than 1E-4 (EDIFFG = 1E-4) but rather be 1E-5 or 1E-6. The energy criterion should only be used if the force criterion is impossible to fulfill for reasonable values (e.g., for huge systems of several 100 atoms)

• NSW = [number] Number of ionic steps, i.e., geometry optimization steps. Usually, quite a large number is needed (even thousands of steps). If you have enough computation time at hand, it is useful to directly set this value to NSW = 500 or NSW = 1000.

• POTIM = [value] The size of a geometry optimization step. Values around POTIM = 0.5 perform reasonable in many cases. If, however, convergence is tricky, the step might be set smaller. Too large steps will lead to explosions!

• NFREE = 15 Number of remembered ionic steps for approximate buildup of the Hessian matrix within the quasi-Newton algorithm. 10-20 is the recommended value for rather large systems in the VASP wiki.

• ISYM = [-1 or 2] For systems without clear symmetry (most large molecular or adsorbate systems), the symmetry should be turned off (ISYM = -1). For highly symmetric systems such as clean metal surfaces, the parameter should be set to (ISYM = 2).

Tricks

Some tricks on the trade might help here, if you :

• Use selective dynamics and set as many degrees as freedom as possible to frozen (F F F), but not the atoms that are crucial for the properties of interest. Usually, lower layers of a metal surface are frozen if an adsorption shall be calculated. The script build_adsorbates.py has an option included to freeze all subtrate atoms that are far enough away from any of the adsorbate atoms.

• Use a large number of ionic steps (see below). Even if the progress of a calculation seems to be nonexistent, it often happens that after 100-200 steps of stagnation, the deciding last movement occurs. The total number of steps might then exceed 1000 or more. If you do your calculations on a cluster with limited walltime, the script opt_long.sh might be ideal for you.

• Raise the plane wave cutoff (ENCUT) of the calculation, e.g., to 600 instead of 400 eV for a systems containing hydrocarbons.

• Raise the NELMIN number to 8 or more. It might happen that the system is trapped and the SCF always converges after 2-3 cycles but nothing goes forward. This change can sometimes release the system.