|
The dynamical matrix code in VASP allows for the calculation of
harmonic frequencies of and the prefactor of a reaction.
The formalism was developed by Vineyard, and the Arrhenius rate
includes the Vineyard prefactor. The basic idea is that
the potential is assumed to be harmonic for both the initial
state and the transition state. With this assumption, one can
formulate the migration prefactor (the product of the attempt
frequency and the exponential of the migration entropy over
k) as the ratio of the products of the frequencies of the
normal modes of the initial state over the transition state.
Another way of saying this is that the entropy of the initial
and transition state is evaluated using a harmonic approximation.
To be able to do a dynamical matrix calculation for a prefactor, one needs
to know beforehand both the initial state and the transition
state. The nudged elastic band (NEB)
method and the dimer method
are two ways of locating transition states.
The dynamical matrix calculation may also be used to calculate the harmonic frequencies
and modes of a given molecule whose geometry has already been optimized in VASP. This dynamical
matrix calculation requires only beforehand knowledge of the equilibrium optimized geometry of the
initial state.
In implementing the dynamical matrix calculation in VASP, several
support programs have been written, both to help set up the
calculation as well as process the output. A page describing these
programs, most of which are Perl scripts, can be found
here.
Setting up a dynamical matrix calculation
This is a flow chart of the steps that need to be done in order to do
a dynamical matrix calculation.
CALCULATING THE FIRST DYNAMICAL MATRIX
- Compile a version of VASP with the vtst source.
- if you want to run VASP, set up the directory structure in a
multiple-geometries-at-once mode in the same way that the nudged
elastic band is run, then create the number of directories you want
to run, from 01 to NN, where NN is the number of directories you
want. To simply calculate vibrational frequencies of a molecule, this step
is not necessary and may be ignored.
- Select the degrees of freedom to calculate in a DISPLACECAR file. This file is a
3 by N list of displacements where N is the number of atoms. A typical displacement and entry
in this list is on the order of 0.001. This file may also be generated with
either the dymselsph.pl script (to select atoms in a radius
about central atom(s) when calculating harmonic frequencies of a molecule)
or the dymseldsp.pl script (to select the atoms that are displaced the most between two POSCARs).
- Set up the INCAR file. You need to have the following variables set:
ICHAIN=1, POTIM=0.0, IMAGES=NN (only if running in multiple geometries at once mode from step 2),
NSW=(Degrees of freedom)/N+1. The reason for the +1 is because each directory will calculate the forces on
the initial structure before doing any displacements.
For instance, in the calculation of ethane
modes, the NSW tag would be 25 (3 degrees of freedom for each of 8 atoms, plus 1). For ethylene,
the NSW tag would be 19. No IMAGES tag would be added unless running in parallel as specified
in step 2.
| Parameter |
Value |
Description |
| ICHAIN | 1 | Run the dynamical matrix code |
| IMAGES | NN | Number of parallel images to run, if desired as in step 2 above. Otherwise, do not add. |
| NSW | (DOF)/N+1 | Number of ionic steps |
| IBRION | 3 | Tell VASP to run dynamics |
| POTIM | 0.0 | with a timestep of zero (ie, do nothing) |
Note: Choosing the number of time steps.
The dynamical matrix calculation may be set up run in parallel, both
with multiple processors working on one degree of freedom, and with
different degrees of freedom in different directories. This
parallelization is identical to that of the nudged elastic band. If
there are N degrees of freedom and directories in which they are
being calculated, then each directory needs to do N/M time steps.
However, each directory does one time step without a displacement, one
time step on the initial structure, and so each directory needs to do
one more time step. Thus, the total number of time steps needed
is: N/M + 1. Also note that M = IMAGES of the INCAR file.
- Run VASP after verifying that POTCAR, POSCAR, INCAR, DISPLACECAR, and KPOINTS files all
exist and are properly configured. To ensure accuracy of calculated modes, the geometry in POSCAR
should have first been optimized in VASP. All INCAR settings such as ISMEAR, SIGMA, PREC, etc used
in the geometry optimization should be present in the dynamical matrix INCAR as well.
- Extract the force constants from the output. You can use the
dymmatrix.pl script to do this. Please see the documentation
for this script for its arguments. The output of this script is
the dynamical matrix, its eigenvalues, the normal mode frequencies
and the normal mode eigenvectors.
- If all calculated frequencies are all zero, check the EDIFFG tag in dynamical matrix INCAR.
AFTER geometry minimization when running the dynamical matrix code, this INCAR tag may be set
arbitrarily small to prevent premature termination of the job (as if minimum structural accuracy
had been obtained). The number of force calls made should be equal to the number specified in the
NSW tag.
- Check for convergence of frequencies with respect to all parameters.
- Calculate the prefactor (not for simple frequency analysis). The script dymprefactor.pl does this.
It takes as input the eigenvalue files of both the minimum and the
transition state. The output is the prefactor in inverse cm.
CALCULATING ADDITIONAL DEGREES OF FREEDOM FOR THE SAME DISPLACEMENT
This completes one prefactor calculation. If you want to do another
calculation with the same system, but including more degrees of freedom,
then you need to:
- Save the OUTCAR files from the previous calculation.
- Create the new DISPLACECAR as above.
- Generate the difference between the new and old DISPLACECAR with
the dymcmpdisp.pl program. You will run the dynamical matrix
code with the output of this program.
- Modify NSW in your INCAR file as before.
- Run Vasp.
- Extract the forces from the OUTCAR files. Use dymmatrix.pl
to do this. To do this, you will need the DISPLACECAR from the
original calculation, the DISPLACECAR from this calculation, and
all of the corresponding OUTCARs. Again, the second and
subsequent DISPLACECARs should only contain the new degrees of
freedom not calculated before.
|