Category:Biased molecular dynamics: Difference between revisions

From VASP Wiki
VisualWikitext
No edit summary
No edit summary
Line 1: Line 1:
== Theory ==
'''Biased molecular dynamics''' (MD) refers to advanced [[:Category:Molecular dynamics|MD-simulation methods]] that introduce a ''bias potential''.
*Biased molecular dynamics: {{TAG|Biased molecular dynamics}}.
 
For instance, umbrella sampling{{cite|torrie:jcp:1977}} uses a bias potential to enhance the sampling of a geometric parameter in regions with low probability density, e.g., transition regions of chemical reactions.
The correct distributions are recovered afterward using an exact relation between the biased and unbiased systems. D. Frenkel and B. Smit provide a more detailed description of the method in Ref. {{cite|frenkel:ap-book:2002}}.
'''Biased MD''' can also be used to introduce soft geometric constraints in which the controlled geometric parameter is not strictly constant. Instead, it oscillates in a narrow interval of values.
 
In VASP, the bias potentials are supported in both the [[NVT ensemble|NVT]] and [[NpT ensemble|NpT]] MD simulations, regardless of the particular [[:Category:Thermostats|thermostat]] and/or [[PSTRESS|barostat]] setting.
Different types of potentials can be freely combined. The values of all collective variables defined in the {{FILE|ICONST}} file for each MD step are listed in the {{FILE|REPORT}} file. Check the lines after the string <tt>Metadynamics</tt>.
[[File:Bias potentials.png|thumb|Graphical representation of (a) harmonic, (b) Fermi function-shaped, and (c) and Gauss function-shaped bias potentials.  ]]
 
== Theoretical background ==
 
The probability density for a geometric parameter &xi; of the system driven by a Hamiltonian:
:<math>
H(q,p) = T(p) + V(q), \;
</math>
with ''T''(''p''), and ''V''(''q'') being kinetic, and potential energies, respectively, can be written as:
:<math>
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}.
</math>
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 on one or multiple selected internal coordinates of the system &xi;=&xi;(''q''), the Hamiltonian takes a form:
:<math>
\tilde{H}(q,p) = H(q,p) + \tilde{V}(\xi),
</math>
and the probability density of &xi; in the biased ensemble is:
:<math>
\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}}.
</math>
It can be shown that the biased and unbiased averages are related via a simple formula:
:<math>
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}}}.
</math>
More generally, an observable <math>\langle A \rangle_{H}</math>:
:<math>
\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}
</math>
can be expressed in terms of thermal averages within the biased ensemble:
:<math>
\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}}}.
</math>
Simulation methods, e.g., umbrella sampling{{cite|torrie:jcp:1977}}, use a bias potential to enhance sampling of &xi; in regions with low ''P''(&xi;<sub>''i''</sub>), e.g., transition regions of chemical
reactions.
The correct distributions are recovered afterward using the equation for <math>\langle A \rangle_{H}</math> above. Ref.{{cite|frenkel:ap-book:2002}} provides a more detailed description of the method.
 
== Supported types of bias potentials ==
Presently, the following types of bias potential are supported:
 
=== Harmonic potentials ===
 
:A sum of Harmonic potentials (see curve (a) on the plot above)
::<math>
\tilde{V}(\xi_1,\dots,\xi_{M_8}) = \sum_{\mu=1}^{M}\frac{1}{2}\kappa_{\mu} (\xi_{\mu}(q)-\xi_{0\mu})^2 \;
</math>
:where the sum runs over all (<math>M_8</math>) coordinates the potential acts upon. The coordinates are defined in the {{FILE|ICONST}} file by setting the <code>status</code> to 8. The parameters of the potential are the force constant <math>\kappa_{\mu}</math> ({{TAG|SPRING_K}}) and the minimum or the potential <math>\xi_{0\mu}</math> {{TAG|SPRING_R0}}. These must be set in the {{FILE|INCAR}} file. Optionally, it is possible to change the value of <math>\xi_{0\mu}</math> every MD step at a constant rate defined via the {{FILE|INCAR}} tag {{TAG|SPRING_V0}}.
{{NB|mind|The number of items defined via {{TAG|SPRING_K}}, {{TAG|SPRING_R0}}, and {{TAG|SPRING_V0}} must be equal to <math>M_8</math>. Otherwise, the calculation terminates with an error message.|:}}
:This form of bias potential is employed in several simulation protocols, such as the umbrella sampling{{cite|torrie:jcp:1977}}, umbrella integration, or steered MD, and is useful also in cases where the <math>\xi_{\mu}</math> values need to be restrained.
 
=== Step function ===
 
:A sum of Fermi-like step functions (see curve (b) on the plot above)
::<math>
\tilde{V}(\xi_1,\dots,\xi_{M_4}) = \sum_{\mu=1}^{M_4}\frac{A_{\mu}}{1+\text{exp}\left [-D_{\mu}(\frac{\xi(q)}{\xi_{0\mu}} -1) \right ]} \;
</math>
:where the sum runs over all (<math>M_4</math>) coordinates the potential acts upon. The coordinates are defined in the {{FILE|ICONST}} file by setting the <code>status</code> to 4. The parameters of the potential are the height of the step (<math>A_{\mu}</math> set by {{TAG|FBIAS_A}}), the slope around the point <math>\xi_{0\mu}</math> (<math>D_{\mu}</math> set by {{TAG|FBIAS_D}}), and the position of the step (<math>\xi_{0\mu}</math> set by {{TAG|FBIAS_R0}}). These must be set in the {{FILE|INCAR}} file.
{{NB|mind|The number of items defined via {{TAG|FBIAS_A}}, {{TAG|FBIAS_D}}, and {{TAG|FBIAS_R0}} must be equal to <math>M_4</math>. Otherwise, the calculation terminates with an error message.|:}}
:This form of potential is suitable especially for imposing restrictions on the upper (or lower) limit of the value of <math>\xi</math>.
 
=== Gaussian potential ===
 
:A sum of Gauss functions (see curve (b) on the plot above)
::<math>
\tilde{V}(\xi_1,\dots,\xi_{M})  = \sum_{\nu=1}^{N_5}h_{\nu}\text{exp}\left [-\frac{\sum_{\mu=1}^{M_5}(\xi_{\mu}(q)-\xi_{0\nu,\mu})^2}{2w_{\nu}^2}  \right ], \;
</math>
:where <math>N_5</math> is the number of Gaussian functions and <math>M_5</math> is the number of coordinates the potential acts upon. The coordinates are defined in the {{FILE|ICONST}} file by setting the <code>status</code> to 5. The parameters of the potentials, <math>h_{\nu}</math>, <math>w_{\nu}</math>, and <math>\xi_{0\nu,\mu}</math> are defined in the {{FILE|PENALTYPOT}} file.
 
:This type of bias potential is primarily intended for use in metadynamics, but since Gaussians can be used as basis functions for more general shapes, they can also be used to prepare various atypically-shaped bias potentials. 
 
== Examples of usage ==
Let us consider the nucleophile substitution reaction of CH<math>_3</math>Cl with Cl<math>^-</math>. The reactant is a weak van-der-Waals complex. The corresponding {{FILE|POSCAR}} file reads 
 
vdW complex CH3Cl...Cl
1.00000000000000
12.0000000000000000    0.0000000000000000    0.0000000000000000
0.0000000000000000    12.0000000000000000    0.0000000000000000
0.0000000000000000    0.0000000000000000    12.0000000000000000
C H Cl
1 3 2
cart
5.91331371  7.11364924  5.78037960
5.81982231  8.15982106  5.46969017
4.92222130  6.65954232  5.88978969
6.47810398  7.03808479  6.71586385
4.32824726  8.75151396  7.80743202
6.84157897  6.18713289  4.46842049
 
Due to the weak interactions between CH<math>_3</math>Cl and Cl<math>^-</math>, the complex can collapse at high temperatures. This can be avoided by setting an upper bound for the length of the non-bonding Cl...C interactions.
This can be conveniently achieved by using a Fermi-like step-shaped bias potential. To this end, we need to define the Cl...C distance, i.e., the distance between the atoms 1 and 5, as a coordinate with status 4 in the {{FILE|ICONST}} file:
 
R 1 5 4 
 
Next, we need to specify the bias potential parameters {{TAG|FBIAS_A}}, {{TAG|FBIAS_D}}, and {{TAG|FBIAS_R0}} in the {{FILE|INCAR}} file. Since the bias potential acts only on one internal coordinate (<math>M_4=1</math>), we need to provide only one number for each of the tags. The setting
 
FBIAS_A  = 1
FBIAS_D  = 50
FBIAS_R0 = 3.5
 
ensures that repulsive bias forces steeply increase when the C...Cl distance is increased beyond about <math>3.2 \AA</math>. This causes a shortening of the distance in the next MD step. Notice that the bias force is essentially negligible for distances below <math>3 \AA</math>.
A careful adjustment of {{TAG|FBIAS_A}} and {{TAG|FBIAS_D}} is needed to ensure that (i) the bias force is large enough to effectively limit the value of <math>\xi</math>, and (ii) the interval of <math>\xi</math> values for which the bias forces are significant is broad enough to avoid overcoming via random fluctuations.     
A suitable setting can be found by noting that the maximal bias force of <math>\frac{D\,A}{4\xi_0}</math> is exerted on the system at the point <math>\xi = \xi_{0}</math>. This can be seen by inspecting the analytical expression for the potential.


== How to ==
== How to ==
*Biased molecular dynamics: {{TAG|Biased molecular dynamics}}.


=== Biased MD run using the Andersen thermostat and a Gaussian bias potential ===
:For a biased MD run with Andersen thermostat, one has to:
:#Set the standard MD-related tags: {{TAG|IBRION}}=0, {{TAG|TEBEG}}, {{TAG|POTIM}}, and {{TAG|NSW}}
:#Set {{TAG|MDALGO}}=1 ({{TAG|MDALGO}}=11 in VASP 5.x), and choose an appropriate setting for {{TAG|ANDERSEN_PROB}}
:#To avoid updating of the bias potential, set {{TAG|HILLS_BIN}}={{TAG|NSW}}
:#Define collective variables in the {{FILE|ICONST}} file, and set the <tt>STATUS</tt> parameter for the collective variables to 5
:#Define the bias potential in the {{FILE|PENALTYPOT}} file
=== Biased MD run using the Nose-Hoover thermostat and a Gaussian bias potential ===
:For a biased MD run with Nose-Hoover thermostat, one has to:
:#Set the standard MD-related tags: {{TAG|IBRION}}=0, {{TAG|TEBEG}}, {{TAG|POTIM}}, and {{TAG|NSW}}
:#Set {{TAG|MDALGO}}=2 ({{TAG|MDALGO}}=21 in VASP 5.x), and choose an appropriate setting for {{TAG|SMASS}}
:#To avoid updating the bias potential, set {{TAG|HILLS_BIN}}={{TAG|NSW}}
:#Define collective variables in the {{FILE|ICONST}} file, and set the <tt>STATUS</tt> parameter for the collective variables to 5
:#Define the bias potential in the {{FILE|PENALTYPOT}} file
The values of all collective variables for each MD step are listed in the {{FILE|REPORT}} file. Check the lines after the string <tt>Metadynamics</tt>.
== References ==
<references/>
----
----


[[Category:VASP|Constrained molecular dynamics]][[Category:Molecular dynamics]]
[[Category:VASP|Biased molecular dynamics]][[Category:Molecular dynamics]][[Category:Theory]][[Category:Howto]]

Revision as of 09:37, 24 April 2023

Biased molecular dynamics (MD) refers to advanced MD-simulation methods that introduce a bias potential.

For instance, umbrella sampling[1] uses a bias potential to enhance the sampling of a geometric parameter in regions with low probability density, e.g., transition regions of chemical reactions. The correct distributions are recovered afterward using an exact relation between the biased and unbiased systems. D. Frenkel and B. Smit provide a more detailed description of the method in Ref. [2]. Biased MD can also be used to introduce soft geometric constraints in which the controlled geometric parameter is not strictly constant. Instead, it oscillates in a narrow interval of values.

In VASP, the bias potentials are supported in both the NVT and NpT MD simulations, regardless of the particular thermostat and/or barostat setting. Different types of potentials can be freely combined. The values of all collective variables defined in the ICONST file for each MD step are listed in the REPORT file. Check the lines after the string Metadynamics.

Graphical representation of (a) harmonic, (b) Fermi function-shaped, and (c) and Gauss function-shaped bias potentials.

Theoretical background

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/":): {\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:

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, e.g., umbrella sampling[1], use a bias potential to enhance sampling of ξ in regions with low P(ξi), e.g., transition regions of chemical reactions. The correct distributions are recovered afterward 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. Ref.[2] provides a more detailed description of the method.

Supported types of bias potentials

Presently, the following types of bias potential are supported:

Harmonic potentials

A sum of Harmonic potentials (see curve (a) on the plot above)
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 {\tilde {V}}(\xi _{1},\dots ,\xi _{M_{8}})=\sum _{\mu =1}^{M}{\frac {1}{2}}\kappa _{\mu }(\xi _{\mu }(q)-\xi _{0\mu })^{2}\;}
where the sum runs over all (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 M_{8}} ) coordinates the potential acts upon. The coordinates are defined in the ICONST file by setting the status to 8. The parameters of the potential are the force constant 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 \kappa _{\mu }} (SPRING_K) and the minimum or the 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/":): {\displaystyle \xi_{0\mu}} SPRING_R0. These must be set in the INCAR file. Optionally, it is possible to change the value of 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\mu}} every MD step at a constant rate defined via the INCAR tag SPRING_V0.
Mind: The number of items defined via SPRING_K, SPRING_R0, and SPRING_V0 must be equal to 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 M_{8}} . Otherwise, the calculation terminates with an error message.
This form of bias potential is employed in several simulation protocols, such as the umbrella sampling[1], umbrella integration, or steered MD, and is useful also in cases where the 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_{\mu}} values need to be restrained.

Step function

A sum of Fermi-like step functions (see curve (b) on the plot above)
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 \tilde{V}(\xi_1,\dots,\xi_{M_4}) = \sum_{\mu=1}^{M_4}\frac{A_{\mu}}{1+\text{exp}\left [-D_{\mu}(\frac{\xi(q)}{\xi_{0\mu}} -1) \right ]} \; }
where the sum runs over all (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 M_4} ) coordinates the potential acts upon. The coordinates are defined in the ICONST file by setting the status to 4. The parameters of the potential are the height of the step (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 A_{\mu}} set by FBIAS_A), the slope around the point 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\mu}} (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 D_{\mu}} set by FBIAS_D), and the position of the step (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\mu}} set by FBIAS_R0). These must be set in the INCAR file.
Mind: The number of items defined via FBIAS_A, FBIAS_D, and FBIAS_R0 must be equal to 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 M_4} . Otherwise, the calculation terminates with an error message.
This form of potential is suitable especially for imposing restrictions on the upper (or lower) limit of the value of 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/":): \xi .

Gaussian potential

A sum of Gauss functions (see curve (b) on the plot above)
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 \tilde{V}(\xi_1,\dots,\xi_{M}) = \sum_{\nu=1}^{N_5}h_{\nu}\text{exp}\left [-\frac{\sum_{\mu=1}^{M_5}(\xi_{\mu}(q)-\xi_{0\nu,\mu})^2}{2w_{\nu}^2} \right ], \; }
where 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 N_5} is the number of Gaussian functions 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 M_5} is the number of coordinates the potential acts upon. The coordinates are defined in the ICONST file by setting the status to 5. The parameters of the potentials, 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 h_{\nu}} , 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 w_{\nu}} , 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\nu,\mu}} are defined in the PENALTYPOT file.
This type of bias potential is primarily intended for use in metadynamics, but since Gaussians can be used as basis functions for more general shapes, they can also be used to prepare various atypically-shaped bias potentials.

Examples of usage

Let us consider the nucleophile substitution reaction of CHFailed 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 _3} Cl with ClFailed 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/":): ^{-} . The reactant is a weak van-der-Waals complex. The corresponding POSCAR file reads 

vdW complex CH3Cl...Cl 
1.00000000000000
12.0000000000000000    0.0000000000000000    0.0000000000000000
0.0000000000000000    12.0000000000000000    0.0000000000000000
0.0000000000000000    0.0000000000000000    12.0000000000000000
C H Cl
1 3 2
cart
5.91331371  7.11364924  5.78037960
5.81982231  8.15982106  5.46969017
4.92222130  6.65954232  5.88978969
6.47810398  7.03808479  6.71586385
4.32824726  8.75151396  7.80743202
6.84157897  6.18713289  4.46842049

Due to the weak interactions between CHFailed 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 _3} Cl and ClFailed 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/":): ^{-} , the complex can collapse at high temperatures. This can be avoided by setting an upper bound for the length of the non-bonding Cl...C interactions. This can be conveniently achieved by using a Fermi-like step-shaped bias potential. To this end, we need to define the Cl...C distance, i.e., the distance between the atoms 1 and 5, as a coordinate with status 4 in the ICONST file: 

R 1 5 4

Next, we need to specify the bias potential parameters FBIAS_A, FBIAS_D, and FBIAS_R0 in the INCAR file. Since the bias potential acts only on one internal coordinate (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 M_4=1} ), we need to provide only one number for each of the tags. The setting 

FBIAS_A  = 1
FBIAS_D  = 50
FBIAS_R0 = 3.5

ensures that repulsive bias forces steeply increase when the C...Cl distance is increased beyond about 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 3.2 \AA} . This causes a shortening of the distance in the next MD step. Notice that the bias force is essentially negligible for distances 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 3 \AA} . A careful adjustment of FBIAS_A and FBIAS_D is needed to ensure that (i) the bias force is large enough to effectively limit the value of 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/":): \xi , and (ii) the interval of 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/":): \xi values for which the bias forces are significant is broad enough to avoid overcoming via random fluctuations. A suitable setting can be found by noting that the maximal bias force of 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 \frac{D\,A}{4\xi_0}} is exerted on the system at the point 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 = \xi_{0}} . This can be seen by inspecting the analytical expression for the potential.

How to

Biased MD run using the Andersen thermostat and a Gaussian bias potential

For a biased MD 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. 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

Biased MD run using the Nose-Hoover thermostat and a Gaussian bias potential

For a biased MD 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. To avoid updating 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. ↑ a b c G. M. Torrie and J. P. Valleau, J. Comp. Phys. 23, 187 (1977).
  2. ↑ a b D. Frenkel and B. Smit, Understanding molecular simulations: from algorithms to applications, Academic Press: San Diego, 2002.

Pages in category "Biased molecular dynamics"

The following 8 pages are in this category, out of 8 total.

F

I

P

S

Retrieved from "https://www.vasp.at/wiki/index.php?title=Category:Biased_molecular_dynamics&oldid=20284"