Sampling phonon spectra from molecular-dynamics simulations

Finite temperature phonon renormalization with velocity correlation functions

To account for finite temperature phonon renormalization, finite temperature molecular dynamics (MD) trajectories of the system under consideration are required. Additionally, harmonic phonon eigenvectors are necessary, as the velocities obtained from the MD simulations must be projected onto these eigenvectors. The molecular dynamics simulations should be conducted in the NVE ensemble, which can be set up using the following example INCAR-file.

INCAR

#INCAR molecular-dynamics tags NVE ensemble 
IBRION = 0                   # choose molecular-dynamics 
MDALGO = 1                   # using Andersen thermostat
ISIF = 2                     # compute stress tensor but do not change box volume/shape 
TEBEG = 300                  # set temperature 
NSW = 10000                  # number of time steps 
POTIM = 1.0                  # time step in femto seconds 
ANDERSEN_PROB = 0.0          # setting Andersen collision probability to zero to get NVE enseble

The forces during the molecular dynamics run can be either obtained by DFT or by machine learning force fields. The trajectories of the molecular dynamics have to be started from structures which were equilibrated for the conditions of interest. The harmonic phonon eigenvectors can be computed with VASP by performing a 3 step procedure

Step 1: Compute the force constants

Using finite differences with IBRION = 5, 6.
Using DFPT with IBRION = 7, 8.

Step 2: Provide q-points file for which the projected velocity correlation functions should be computed.
Step 3: Compute the eigenvectors by setting LPHON_DISPERSION = True and PHON_NWRITE = -3.

Further information can be found on the following page. External tools as for example phonopy may also be considered. To compute the power spectra of the Fourier transformed projected velocity autocorrelations

[math]\displaystyle{ |G_{\nu}(\mathbf{q},\omega)|^{2}=\sum_{I,\alpha}\sum_{J,\beta} \int\left( \varepsilon_{I\nu}^{\beta}(\mathbf{q}) \sqrt{M_{I}}v_{I}^{\alpha}(t') \right )\left( \varepsilon_{J\nu}^{\beta}(\mathbf{q}) \sqrt{M_{J}}v_{J}^{\beta}(t'') \right )e^{i\mathbf{q} \cdot (\mathbf{R}_{I}(t')-\mathbf{R}_{J}(t''))}e^{-i\omega (t'-t'')}d(t'-t'') }[/math]

external tools are required. The following table summarizes a small list of codes which can compute projected velocity correlation functions from VASP output.

code	publication
DSLEAP	Lahnsteiner et.al.
phq	Zhang et.al.
DynaPhoPy	Carreras et.al.

Finite temperature phonon renormalization with velocity correlation functions

INCAR

Related tags and articles