IBRION=5 and IBRION=6

IBRION=5, is only supported starting from VASP.4.5. IBRION=6, is only supported starting from VASP.5.1. Both flags allow to determine the Hessian matrix (matrix of the second derivatives of the energy with respect to the atomic positions) and the vibrational frequencies of a system. Only zone centered ($\Gamma$-point) frequencies are calculated automatically and printed after

 Eigenvectors and eigenvalues of the dynamical matrix

To calculate the Hessian matrix, finite differences are used, i.e. each ion is displaced in the direction of each Cartesian coordinate, and from the forces the Hessian matrix is determined. The two modes differ in the way symmetry is considered. For IBRION=5, all atoms are displaced in all three Cartesian directions, resulting in a significant computational effort even for moderately sized high symmetry systems. For IBRION=6, however only symmetry inequivalent displacements are considered, and the reminder of the Hessian matrix is filled using symmetry considerations.

Selective dynamics are presently only supported for IBRION=5; in this case, only those components of the Hessian matrix are calculated for which the selective dynamics tags are set to .TRUE. Contrary to the conventional behavior, the selective dynamics tags now refer to the Cartesian components of the Hessian matrix. For the following POSCAR file, for instance,

Cubic BN
   3.57
 0.0 0.5 0.5
 0.5 0.0 0.5
 0.5 0.5 0.0
   1 1
selective
Direct
 0.00 0.00 0.00  F F F
 0.25 0.25 0.25  T F F

atom 2 is displaced in the $\hat x$-direction only, and only the $\hat x$ component of the second atom of the Hessian matrix is calculated.

Three parameters influence the determination of the Hessian matrix. The parameter NFREE determines how many displacements are used for each direction and ion, and POTIM determines the step size. The step size is defaulted to 0.015 Å, if too large values are supplied in the input file. Expertise shows that this is a very reasonable compromise. NFREE=2 uses central difference, i.e. each ion is displaced in each direction by a small positive and negative displacement

$\pm$ POTIM $\times$ $\hat x$ , $\pm$ POTIM $\times$ $\hat y$, $\pm$ POTIM $\times$ $\hat z$

For NFREE=4, four displacement are used

$\pm$ POTIM $\times$ $\hat x$ and $\pm$ 2 POTIM $\times$ $\hat x$
$\pm$ POTIM $\times$ $\hat y$ and $\pm$ 2 POTIM $\times$ $\hat y$
...

For NFREE=1, only a single displacement is applied (it is strongly recommend to avoid NFREE=1).

Finally, IBRION=6 and ISIF$\ge$3 allows to calculate the elastic constants. The elastic tensor is determined by performing six finite distortions of the lattice and deriving the elastic constants from the strain-stress relationship [4]. The elastic tensor is calculated both, for rigid ions, as well, as allowing for relaxation of the ions. The elastic moduli for rigid ions are written after the line

SYMMETRIZED ELASTIC MODULI (kBar)

The ionic contributions are determined by inverting the ionic Hessian matrix and multiplying with the internal strain tensor [5], and the corresponding contributions are written after the lines:

 ELASTIC MODULI CONTR FROM IONIC RELAXATION (kBar)

The final elastic moduli including both, the contributions for distortions with rigid ions and the contributions from the ionic relaxations, are summarized at the very end.

 TOTAL ELASTIC MODULI (kBar)

There are a few caveats to this approach: most notably the plane wave cutoff needs to be sufficiently large to converge the stress tensor. This is usually only achieved if the default cutoff is increased by roughly 30 %, but it is strongly recommended to increase the cutoff systematically (e.g. in steps of 15 %), until full convergence is achieved.

Mind: In some older versions, NSW (number of ionic steps) must be set to 1 in the INCAR file, since NSW=0 resets the IBRION tag to $-1$ regardless of the value supplied in the INCAR file.

A final problem concerns the symmetry treatment in VASP.4.6. VASP determines the symmetry for the displaced configurations correctly, but unfortunately VASP does not change the set of k-points automatically (often the lower symmetry of configurations with displaced ions would require one to use more $k-$points). Hence, for accurate calculations, the symmetry must be switched off, or a $k-$point set which has not been reduced using symmetry considerations must be applied. VASP.5.1 changes the k-point set on the fly and the previous restriction does not apply.