Biased molecular dynamics: Difference between revisions

From VASP Wiki
VisualWikitext
No edit summary
No edit summary
Line 11: Line 11:
The term <math>\langle X \rangle_H</math> stands for a thermal average of quantity ''X'' evaluated for the system driven by the Hamiltonian ''H''.
The term <math>\langle X \rangle_H</math> 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 <math>\tilde{V}(\xi)</math> acting only on a selected internal parameter of the system &xi;=&xi;(''q''), the Hamiltonian takes a form:
If the system is modified by adding a bias potential <math>\tilde{V}(\xi)</math> acting on one or multiple selected internal coordinates of the system &xi;=&xi;(''q''), the Hamiltonian takes a form:
:<math>
:<math>
\tilde{H}(q,p) = H(q,p) + \tilde{V}(\xi),
\tilde{H}(q,p) = H(q,p) + \tilde{V}(\xi),
Line 42: Line 42:
Presently, the following types of bias potential are supported:
Presently, the following types of bias potential are supported:


*Gauss function
 
*harmonic potential (see curve (a) on the plot below)
:<math>
:<math>
V_{bias}(\xi) = h\,\text{exp}\left [-\frac{(\xi(q)-\xi_0)^2}{2w^2}  \right ], \;
V_{bias}(\xi) = \sum_{\mu}^{M}\frac{1}{2}\kappa_{\mu} (\xi(q)-\xi_{0,\mu})^2 \;
</math>
</math>
where the sum runs all coordinates the potential acts upon, which are defined in the {{FILE|ICONST}}-file by setting the status to 8.
The parameters <math>\kappa</math> (force constant) and <math>\xi_0</math> (minimum of potential)
The potential acts on coordinates with status 8 defined in the {{FILE|ICONST}}-file. 


*harmonic potential
*Gauss function (see curve (b) on the plot below)
:<math>
:<math>
V_{bias}(\xi) = \frac{1}{2}\kappa (\xi(q)-\xi_0)^2 \;
V_{bias}(\xi) = h\,\text{exp}\left [-\frac{(\xi(q)-\xi_0)^2}{2w^2}  \right ], \;
</math>
</math>




*Fermi function
 
*Fermi function (see curve (b) on the plot below)
:<math>
:<math>
V_{bias}(\xi) = \frac{A}{1+\text{exp}\left [-D\frac{\xi(q)}{\xi_0} -1 \right ]} \;
V_{bias}(\xi) = \frac{A}{1+\text{exp}\left [-D\frac{\xi(q)}{\xi_0} -1 \right ]} \;

Revision as of 09:51, 6 April 2023

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

Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://www.vasp.at/wiki/restbase/vasp.at/v1/":): H(q,p)=T(p)+V(q),\;

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

Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://www.vasp.at/wiki/restbase/vasp.at/v1/":): 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 Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://www.vasp.at/wiki/restbase/vasp.at/v1/":): \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 Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://www.vasp.at/wiki/restbase/vasp.at/v1/":): {\tilde {V}}(\xi ) acting on one or multiple selected internal coordinates of the system ξ=ξ(q), the Hamiltonian takes a form:

Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://www.vasp.at/wiki/restbase/vasp.at/v1/":): {\tilde {H}}(q,p)=H(q,p)+{\tilde {V}}(\xi ),

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

Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://www.vasp.at/wiki/restbase/vasp.at/v1/":): {\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:

Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://www.vasp.at/wiki/restbase/vasp.at/v1/":): 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 Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://www.vasp.at/wiki/restbase/vasp.at/v1/":): \langle A\rangle _{{H}} :

Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://www.vasp.at/wiki/restbase/vasp.at/v1/":): \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:

Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://www.vasp.at/wiki/restbase/vasp.at/v1/":): \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 P(ξi) such as transition regions of chemical reactions. The correct distributions are recovered afterwards using the equation for Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://www.vasp.at/wiki/restbase/vasp.at/v1/":): \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.

Supported types of bias potentials

Presently, the following types of bias potential are supported:


  • harmonic potential (see curve (a) on the plot below)
Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://www.vasp.at/wiki/restbase/vasp.at/v1/":): {\displaystyle V_{bias}(\xi) = \sum_{\mu}^{M}\frac{1}{2}\kappa_{\mu} (\xi(q)-\xi_{0,\mu})^2 \; }

where the sum runs all coordinates the potential acts upon, which are defined in the ICONST-file by setting the status to 8. The parameters Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://www.vasp.at/wiki/restbase/vasp.at/v1/":): \kappa (force constant) and Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://www.vasp.at/wiki/restbase/vasp.at/v1/":): {\displaystyle \xi_0} (minimum of potential) The potential acts on coordinates with status 8 defined in the ICONST-file.

  • Gauss function (see curve (b) on the plot below)
Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://www.vasp.at/wiki/restbase/vasp.at/v1/":): {\displaystyle V_{bias}(\xi) = h\,\text{exp}\left [-\frac{(\xi(q)-\xi_0)^2}{2w^2} \right ], \; }


  • Fermi function (see curve (b) on the plot below)
Failed to parse (SVG (MathML can be enabled via browser plugin): Invalid response ("Math extension cannot connect to Restbase.") from server "https://www.vasp.at/wiki/restbase/vasp.at/v1/":): {\displaystyle V_{bias}(\xi) = \frac{A}{1+\text{exp}\left [-D\frac{\xi(q)}{\xi_0} -1 \right ]} \; }

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.
Retrieved from "https://www.vasp.at/wiki/index.php?title=Biased_molecular_dynamics&oldid=20085"