Requests for technical support from the VASP group should be posted in the VASP-forum.

# Biased molecular dynamics

Jump to navigation Jump to search

The probability density for a geometric parameter ξ of the system driven by a Hamiltonian:

${\displaystyle H(q,p)=T(p)+V(q),\;}$

with T(p), and V(q) being kinetic, and potential energies, respectively, can be written as:

${\displaystyle P(\xi _{i})={\frac {\int \delta {\Big (}\xi (q)-\xi _{i}{\Big )}\exp \left\{-H(q,p)/k_{B}\,T\right\}dq\,dp}{\int \exp \left\{-H(q,p)/k_{B}\,T\right\}dq\,dp}}=\langle \delta {\Big (}\xi (q)-\xi _{i}{\Big )}\rangle _{H}.}$

The term ${\displaystyle \langle X\rangle _{H}}$ stands for a thermal average of quantity X evaluated for the system driven by the Hamiltonian H.

If the system is modified by adding a bias potential ${\displaystyle {\tilde {V}}(\xi )}$ acting only on a selected internal parameter of the system ξ=ξ(q), the Hamiltonian takes a form:

${\displaystyle {\tilde {H}}(q,p)=H(q,p)+{\tilde {V}}(\xi ),}$

and the probability density of ξ in the biased ensemble is:

${\displaystyle {\tilde {P}}(\xi _{i})={\frac {\int \delta {\Big (}\xi (q)-\xi _{i}{\Big )}\exp \left\{-{\tilde {H}}(q,p)/k_{B}\,T\right\}dq\,dp}{\int \exp \left\{-{\tilde {H}}(q,p)/k_{B}\,T\right\}dq\,dp}}=\langle \delta {\Big (}\xi (q)-\xi _{i}{\Big )}\rangle _{\tilde {H}}}$

It can be shown that the biased and unbiased averages are related via a simple formula:

${\displaystyle P(\xi _{i})={\tilde {P}}(\xi _{i}){\frac {\exp \left\{{\tilde {V}}(\xi )/k_{B}\,T\right\}}{\langle \exp \left\{{\tilde {V}}(\xi )/k_{B}\,T\right\}\rangle _{\tilde {H}}}}.}$

More generally, an observable ${\displaystyle \langle A\rangle _{H}}$:

${\displaystyle \langle A\rangle _{H}={\frac {\int A(q)\exp \left\{-H(q,p)/k_{B}\,T\right\}dq\,dp}{\int \exp \left\{-H(q,p)/k_{B}\,T\right\}dq\,dp}}}$

can be expressed in terms of thermal averages within the biased ensemble:

${\displaystyle \langle A\rangle _{H}={\frac {\langle A(q)\,\exp \left\{{\tilde {V}}(\xi )/k_{B}\,T\right\}\rangle _{\tilde {H}}}{\langle \exp \left\{{\tilde {V}}(\xi )/k_{B}\,T\right\}\rangle _{\tilde {H}}}}.}$

Simulation methods such as umbrella sampling[1] use a bias potential to enhance sampling of ξ in regions with low Pi) such as transition regions of chemical reactions. The correct distributions are recovered afterwards using the equation for ${\displaystyle \langle A\rangle _{H}}$ above.

A more detailed description of the method can be found in Ref.[2]. Biased molecular dynamics can be used also to introduce soft geometric constraints in which the controlled geometric parameter is not strictly constant, instead it oscillates in a narrow interval of values.

## Andersen thermostat

• For a biased molecular dynamics run with Andersen thermostat, one has to:
1. Set the standard MD-related tags: IBRION=0, TEBEG, POTIM, and NSW
2. Set MDALGO=1 (MDALGO=11 in VASP 5.x), and choose an appropriate setting for ANDERSEN_PROB
3. In order to avoid updating of the bias potential, set HILLS_BIN=NSW
4. Define collective variables in the ICONST-file, and set the STATUS parameter for the collective variables to 5
5. Define the bias potential in the PENALTYPOT-file

## Nose-Hoover thermostat

• For a biased molecular dynamics run with Nose-Hoover thermostat, one has to:
1. Set the standard MD-related tags: IBRION=0, TEBEG, POTIM, and NSW
2. Set MDALGO=2 (MDALGO=21 in VASP 5.x), and choose an appropriate setting for SMASS
3. In order to avoid updating of the bias potential, set HILLS_BIN=NSW
4. Define collective variables in the ICONST-file, and set the STATUS parameter for the collective variables to 5
5. Define the bias potential in the PENALTYPOT-file

The values of all collective variables for each MD step are listed in the REPORT-file, check the lines after the string Metadynamics.

## References

1. G. M. Torrie and J. P. Valleau, J. Comp. Phys. 23, 187 (1977).
2. D. Frenkel and B. Smit, Understanding molecular simulations: from algorithms to applications, Academic Press: San Diego, 2002.