MDALGO: Difference between revisions

From VASP Wiki
VisualWikitext
 
(82 intermediate revisions by 6 users not shown)
Line 1: Line 1:
{{TAGDEF|MDALGO|0 {{!}} 1 {{!}} 2 {{!}} 3 {{!}} 11 {{!}} 21 {{!}} 13|0}}
{{TAGDEF|MDALGO|0 {{!}} 1 {{!}} 2 {{!}} 3 {{!}} 4 {{!}} 5 {{!}} 11 {{!}} 21 {{!}} 13 |0}}


Description: {{TAG|MDALGO}} specifies the molecular dynamics simulation protocol (in case {{TAG|IBRION}}=0 and VASP was compiled with [[Precompiler_flags|-Dtbdyn]]).  
Description: {{TAG|MDALGO}} specifies the molecular-dynamics-simulation protocol (in case {{TAG|IBRION}}=0 and VASP was compiled with [[Precompiler options#-Dtbdyn|-Dtbdyn]]).  
----
----


== {{TAG|MDALGO}}=0: Standard molecular dynamics ==
== {{TAG|MDALGO}}=0: Standard molecular dynamics ==
Standard molecular dynamics ({{TAG|IBRION}}=0), the same behavior as if VASP were compiled without [[Precompiler_flags|-Dtbdyn]].
Should provide the same results as {{TAG|MDALGO}}=2 ({{TAG|Nose-Hoover thermostat}}). The difference is that it is a different implementation and works also without the precompiler flag [[Precompiler_flags|-Dtbdyn]].


== {{TAG|MDALGO}}=1: ''NVT''-ensemble with Andersen thermostat ==
=== NVE ensemble ===
<div id="Andersen"></div>
=== Andersen thermostat ===
''NVT''-simulation with Andersen thermostat. In the approach proposed by Andersen<ref name="Andersen80"/> the system is thermally coupled to a fictitious heat bath with the desired temperature. The coupling is represented by stochastic impulsive forces that act occasionally on randomly selected particles. The collision probability is defined as an average number of collisions per atom and time-step. This quantity can be controlled by the flag {{TAG|ANDERSEN_PROB}}. The total number of collisions with the heat-bath is written in the file {{FILE|REPORT}} for each MD step.


<div id="ConstrMD"></div>
*To perform a calculation in the NVE ensemble:
=== Constrained molecular dynamics ===
#Set the standard MD-related tags: {{TAG|IBRION}}=0, {{TAG|TEBEG}}, {{TAG|POTIM}}, and {{TAG|NSW}}.
{{TAG|Constrained molecular dynamics}}
#Set {{TAG|MDALGO}}=0 and set {{TAG|SMASS}}=-3.
 
== {{TAG|MDALGO}}=1: Andersen thermostat ==
 
For the description of the Andersen thermostat see: {{TAG|Andersen thermostat}}.
 
The Andersen thermostat is only available for the NVT ensemble.


=== Slow-growth approach ===
=== Standard molecular dynamics in ===
In general, constrained molecular dynamics generates biased statistical averages.
It can be shown that the correct average for a quantity <math>a(\xi)</math> can be obtained using the formula:
:<math>
a(\xi)=\frac{\langle |\mathbf{Z}|^{-1/2} a(\xi^*) \rangle_{\xi^*}}{\langle |\mathbf{Z}|^{-1/2}\rangle_{\xi^*}},
</math>
where <math>\langle ... \rangle_{\xi^*}</math> stands for the statistical average of the quantity enclosed in angular parentheses computed for a constrained ensemble and <math>Z</math> is a mass metric tensor defined as:
:<math>
Z_{\alpha,\beta}={\sum}_{i=1}^{3N} m_i^{-1} \nabla_i \xi_\alpha \cdot \nabla_i \xi_\beta, \, \alpha=1,...,r, \, \beta=1,...,r,
</math>
It can be shown that the free energy gradient can be computed using the equation:<ref name="Carter89"/><ref name="Otter00"/><ref name="Darve02"/><ref name="Fleurat05"/>
:<math>
\Bigl(\frac{\partial A}{\partial \xi_k}\Bigr)_{\xi^*}=\frac{1}{\langle|Z|^{-1/2}\rangle_{\xi^*}}\langle |Z|^{-1/2} [\lambda_k +\frac{k_B T}{2 |Z|} \sum_{j=1}^{r}(Z^{-1})_{kj} \sum_{i=1}^{3N} m_i^{-1}\nabla_i \xi_j \cdot \nabla_i |Z|]\rangle_{\xi^*},
</math>
where <math>\lambda_{\xi_k}</math> is the Lagrange multiplier associated with the parameter <math>{\xi_k}</math> used in the [[#SHAKE|SHAKE algorithm]].<ref name="Ryckaert77"/>


The free-energy difference between states (1) and (2) can be computed by integrating the free-energy gradients over a connecting path:
*For a standard molecular dynamics run with Anderson thermostat, one has to:
:<math>
#Set the standard MD-related tags: {{TAG|IBRION}}=0, {{TAG|TEBEG}}, {{TAG|POTIM}}, and {{TAG|NSW}}.
{\Delta}A_{1 \rightarrow 2} = \int_{{\xi(1)}}^{{\xi(2)}}\Bigl( \frac{\partial {A}} {\partial \xi} \Bigr)_{\xi^*} \cdot d{\xi}.
#Set {{TAG|MDALGO}}=1, and choose an appropriate setting for {{TAG|ANDERSEN_PROB}}
</math>
Note that as the free-energy is a state quantity, the choice of path connecting (1) with (2) is irrelevant.


=== Constrained molecular dynamics ===
For a description of constrained molecular dynamics see {{TAG|Constrained molecular dynamics}}.


* For a constrained molecular dynamics run with Andersen thermostat, one has to:
* For a constrained molecular dynamics run with Andersen thermostat, one has to:
#Set the standard MD-related tags: {{TAG|IBRION}}=0, {{TAG|TEBEG}}, {{TAG|POTIM}}, and {{TAG|NSW}}
#Set the standard MD-related tags: {{TAG|IBRION}}=0, {{TAG|TEBEG}}, {{TAG|POTIM}}, and {{TAG|NSW}}
#Set {{TAG|MDALGO}}=1, and choose an appropriate setting for {{TAG|ANDERSEN_PROB}}
#Set {{TAG|MDALGO}}=1, and choose an appropriate setting for {{TAG|ANDERSEN_PROB}}
#Define geometric constraints in the {{FILE|ICONST}}-file, and set the {{TAG|STATUS}} parameter for the constrained coordinates to 0
#Define geometric constraints in the {{FILE|ICONST}}-file, and set the STATUS parameter for the constrained coordinates to 0
#When the free-energy gradient is to be computed, set {{TAG|LBLUEOUT}}=.TRUE.
 
=== Slow-growth approach ===
For a description of slow-growth approach see {{TAG|Slow-growth approach}}.
* For a slow-growth simulation, one has to perform a calcualtion very similar to {{TAG|Constrained molecular dynamics}} but additionally the transformation velocity-related {{TAG|INCREM}}-tag for each geometric parameter with <tt>STATUS=0</tt> has to be specified:
#Set the standard MD-related tags: {{TAG|IBRION}}=0, {{TAG|TEBEG}}, {{TAG|POTIM}}, and {{TAG|NSW}}
#Set {{TAG|MDALGO}}=1, and choose an appropriate setting for {{TAG|ANDERSEN_PROB}}
#Define geometric constraints in the {{FILE|ICONST}}-file, and set the STATUS parameter for the constrained coordinates to 0
#When the free-energy gradient is to be computed, set {{TAG|LBLUEOUT}}=.TRUE.
#When the free-energy gradient is to be computed, set {{TAG|LBLUEOUT}}=.TRUE.


For a slow-growth simulation, one has to additionally:
<ol start="5">
<ol start="5">
<li>Specify the transformation velocity-related {{TAG|INCREM}}-tag for each geometric parameter with <tt>STATUS=0</tt></li>
<li>Specify the transformation velocity-related {{TAG|INCREM}}-tag for each geometric parameter with <tt>STATUS=0</tt>.</li>
</ol>
</ol>
VASP can handle multiple (even redundant) constraints. Note, however, that a too large number of constraints can cause problems with the stability of the [[#SHAKE|SHAKE algorithm]]. In problematic cases, it is recommended to use a looser convergence criterion (see {{TAG|SHAKETOL}}) and to allow a larger number of iterations (see {{TAG|SHAKEMAXITER}}) in the [[#SHAKE|SHAKE algorithm]]. Hard constraints may also be used in [[#Metadynamics|metadynamics simulations]] (see {{TAG|MDALGO}}=11 {{!}} 21). Information about the constraints is written onto the {{FILE|REPORT}}-file: check the lines following the string: <tt>Const_coord</tt>


=== Monitoring geometric parameters ===
=== Monitoring geometric parameters ===
Line 58: Line 52:
Sometimes it is desirable to terminate the simulation if all values of monitored parameters get larger that some predefined upper and/or lower limits. These limits can be set by the user by means of the {{TAG|VALUE_MAX}} and {{TAG|VALUE_MIN}}-tags.
Sometimes it is desirable to terminate the simulation if all values of monitored parameters get larger that some predefined upper and/or lower limits. These limits can be set by the user by means of the {{TAG|VALUE_MAX}} and {{TAG|VALUE_MIN}}-tags.


To monitor geometric parameters during an MD run:
*To monitor geometric parameters during an MD run:
#Set the standard MD-related tags: {{TAG|IBRION}}=0, {{TAG|TEBEG}}, {{TAG|POTIM}}, and {{TAG|NSW}}
#Set the standard MD-related tags: {{TAG|IBRION}}=0, {{TAG|TEBEG}}, {{TAG|POTIM}}, and {{TAG|NSW}}.
#Set {{TAG|MDALGO}}=1, and choose an appropriate setting for {{TAG|ANDERSEN_PROB}}
#Set {{TAG|MDALGO}}=1, and choose an appropriate setting for {{TAG|ANDERSEN_PROB}}.
#Define geometric constraints in the {{FILE|ICONST}}-file, and set the <tt>STATUS</tt> parameter for the constrained coordinates to 7
#Define geometric constraints in the {{FILE|ICONST}}-file, and set the <tt>STATUS</tt> parameter for the constrained coordinates to 7.
#Optionally, set the upper and/or lower limits for the coordinates, by means of the {{TAG|VALUE_MAX}} and {{TAG|VALUE_MIN}}-tags, respectively.
#Optionally, set the upper and/or lower limits for the coordinates, by means of the {{TAG|VALUE_MAX}} and {{TAG|VALUE_MIN}}-tags, respectively.


=== Metadynamics ===
For a description of metadynamics see {{TAG|Metadynamics}}.


== {{TAG|MDALGO}}=2: ''NVT''-ensemble with Nose-Hoover thermostat==
* For a metadynamics run with Andersen thermostat, one has to:
Nose-Hoover thermostat ({{TAG|SMASS}} needs to be specified in the {{FILE|INCAR}} file).
#Set the standard MD-related tags: {{TAG|IBRION}}=0, {{TAG|TEBEG}}, {{TAG|POTIM}}, and {{TAG|NSW}}.
#Set {{TAG|MDALGO}}=1 (or {{TAG|MDALGO}}=11 in VASP 5.x), and choose an appropriate setting for {{TAG|ANDERSEN_PROB}}.
#Set the parameters {{TAG|HILLS_H}}, {{TAG|HILLS_W}}, and {{TAG|HILLS_BIN}}.
#Define collective variables in the {{FILE|ICONST}}-file, and set the {{TAG|STATUS}} parameter for the collective variables to 5.
#If needed, define the bias potential in the {{FILE|PENALTYPOT}}-file.
 
The actual time-dependent bias potential is written to the {{FILE|HILLSPOT}}-file, which is updated after adding a new Gaussian. At the beginning of the simulation, VASP attempts to read the initial bias potential from the {{FILE|PENALTYPOT}}-file. For the continuation of a metadynamics run, copy {{FILE|HILLSPOT}} to {{FILE|PENALTYPOT}}. The values of all collective variables for each MD step are listed in {{FILE|REPORT}}-file, check the lines after the string <tt>Metadynamics</tt>.
 
=== Biased molecular dynamics ===
 
For a description of biased molecular dynamics see {{TAG|Biased molecular dynamics}}.
 
* For a biased molecular dynamics 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 (or {{TAG|MDALGO}}=11 in VASP 5.x), and choose an appropriate setting for {{TAG|ANDERSEN_PROB}}.
#In order 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 if necessary.
 
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>.
 
=== Special case: NVE ensemble ===
 
NVE ensemble calculations can be also run by selecting the Anderson thermostat and setting the update collision probability ({{TAG|ANDERSEN_PROB}}) to zero.
*To run an NVE ensemble:
#Set the standard MD-related tags: {{TAG|IBRION}}=0, {{TAG|TEBEG}}, {{TAG|POTIM}}, and {{TAG|NSW}}.
#Set {{TAG|MDALGO}}=1 and {{TAG|ANDERSEN_PROB}}=0.0.
 
== {{TAG|MDALGO}}=2: Nose-Hoover thermostat==
 
For the description of the Nose-Hoover thermostat see: {{TAG|Nose-Hoover thermostat}}.
 
The Nose-Hoover thermostat is currently only available for the NVT ensemble.
 
=== Standard molecular dynamics ===
 
* For a standard molecular dynamics 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, and choose an appropriate setting for {{TAG|SMASS}}.


=== Constrained molecular dynamics ===
=== Constrained molecular dynamics ===
{{TAG|Constrained molecular dynamics}}
For a description of constrained molecular dynamics see {{TAG|Constrained molecular dynamics}}.
 
* For a constrained molecular dynamics 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, and choose an appropriate setting for {{TAG|SMASS}}.
#Define geometric constraints in the {{FILE|ICONST}}-file, and set the STATUS parameter for the constrained coordinates to 0.
#When the free-energy gradient is to be computed, set {{TAG|LBLUEOUT}}=.TRUE.


=== Slow-growth approach ===
=== Slow-growth approach ===
For a short description of the slow-growth approach, see [[#Slowgro|its section]] under {{TAG|MDALGO}}=1.


* For a constrained molecular dynamics run with Nose-Hoover thermostat, one has to:
For a description of slow-growth approach see {{TAG|Slow-growth approach}}.
 
* For a slow-growth approach 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 the standard MD-related tags: {{TAG|IBRION}}=0, {{TAG|TEBEG}}, {{TAG|POTIM}}, and {{TAG|NSW}}
#Set {{TAG|MDALGO}}=2, and choose an appropriate setting for {{TAG|SMASS}}
#Set {{TAG|MDALGO}}=2, and choose an appropriate setting for {{TAG|SMASS}}
#Define geometric constraints in the {{FILE|ICONST}}-file, and set the <tt>STATUS</tt> parameter for the constrained coordinates to 0
#Define geometric constraints in the {{FILE|ICONST}}-file, and set the <tt>STATUS</tt> parameter for the constrained coordinates to 0
#When the free-energy gradient is to be computed, set {{TAG|LBLUEOUT}}=.TRUE.
#When the free-energy gradient is to be computed, set {{TAG|LBLUEOUT}}=.TRUE.
For a slow-growth simulation, one has to additionally:
<ol start="5">
<ol start="5">
<li>Specify the transformation velocity-related {{TAG|INCREM}}-tag for each geometric parameter with <tt>STATUS=0</tt></li>
<li>Specify the transformation velocity-related {{TAG|INCREM}}-tag for each geometric parameter with <tt>STATUS=0</tt></li>
</ol>
</ol>
VASP can handle multiple (even redundant) constraints. Note, however, that a too large number of constraints can cause problems with the stability of the [[#SHAKE|SHAKE algorithm]]. In problematic cases, it is recommended to use a looser convergence criterion (see {{TAG|SHAKETOL}}) and to allow a larger number of iterations (see {{TAG|SHAKEMAXITER}}) in the [[#SHAKE|SHAKE algorithm]]. Hard constraints may also be used in [[#Metadynamics|metadynamics simulations]] (see {{TAG|MDALGO}}=11 {{!}} 21). Information about the constraints is written onto the {{FILE|REPORT}}-file: check the lines following the string: <tt>Const_coord</tt>


=== Monitoring geometric parameters ===
=== Monitoring geometric parameters ===
Line 100: Line 137:
<div id="Langevin"></div>
<div id="Langevin"></div>


=== Metadynamics ===
For a description of metadynamics see {{TAG|Metadynamics}}.
* For a metadynamics 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 (or {{TAG|MDALGO}}=21 in VASP 5.x), and choose an appropriate setting for {{TAG|SMASS}}.
#Set the parameters {{TAG|HILLS_H}}, {{TAG|HILLS_W}}, and {{TAG|HILLS_BIN}}.
#Define collective variables in the {{FILE|ICONST}}-file, and set the <tt>STATUS</tt> parameter for the collective variables to 5.
#If needed, define the bias potential in the {{FILE|PENALTYPOT}}-file.
The actual time-dependent bias potential is written to the {{FILE|HILLSPOT}}-file, which is updated after adding a new Gaussian. At the beginning of the simulation, VASP attempts to read the initial bias potential from the {{FILE|PENALTYPOT}}-file. For the continuation of a metadynamics run, copy {{FILE|HILLSPOT}} to {{FILE|PENALTYPOT}}. The values of all collective variables for each MD step are listed in {{FILE|REPORT}}-file, check the lines after the string <tt>Metadynamics</tt>.
=== Biased molecular dynamics ===
For a description of biased molecular dynamics see {{TAG|Biased molecular dynamics}}.
* For a biased molecular dynamics 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 (or {{TAG|MDALGO}}=21 in VASP 5.x), and choose an appropriate setting for {{TAG|SMASS}}.
#In order 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.
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>.


== {{TAG|MDALGO}}=3: Langevin thermostat ==
== {{TAG|MDALGO}}=3: Langevin thermostat ==
'''Note:''' Geometric constraints and metadynamics are not available for Langevin dynamics.
 
For the description of the Langevin thermostat see: {{TAG|Langevin thermostat}}.
 
The Langevin thermostat is available for [[NVT ensemble|NVT]] and [[NpT ensemble|NpT]] ensembles.


=== ''NVT''-simulation with Langevin thermostat ===
=== ''NVT''-simulation with Langevin thermostat ===
<span id="LangevinEOM">
The Langevin thermostat<ref name="Allen91"/> maintains the temperature through a modification of Newton's equations of motion
:<math>
\dot{r_i} = p_i/m_i \qquad
\dot{p_i} = F_i - {\gamma}_i\,p_i + f_i,
</math>
</span>
where ''F<sub>i</sub>'' is the force acting on atom ''i'' due to the interaction potential, &gamma;<sub>i</sub> is a friction coefficient, and ''f<sub>i</sub>'' is a random force with dispersion &sigma;<sub>i</sub> related to &gamma;<sub>i</sub> through:
:<math>
\sigma_i^2 = 2\,m_i\,{\gamma}_i\,k_B\,T/{\Delta}t
</math>
with &Delta;''t'' being the time-step used in the MD to integrate the equations of motion. Obviously, Langevin dynamics is identical to the classical Hamiltonian in the limit of vanishing &gamma;.


*To run an ''NVT''-simulation with a Langevin thermostat, one has to:
*To run an ''NVT''-simulation with a Langevin thermostat, one has to:
#Set the standard MD-related tags: {{TAG|IBRION}}=0, {{TAG|TEBEG}}, {{TAG|POTIM}}, and {{TAG|NSW}}
{{Template:NVT_Langevin_thermostat_recipe}}
#Set {{TAG|ISIF}}=2
#Set {{TAG|MDALGO}}=3 to invoke the Langevin thermostat
#Specify friction coefficients for all species in the {{FILE|POSCAR}} file, by means of the {{TAG|LANGEVIN_GAMMA}}-tag.


<div id="ParrinelloRahman"></div>
=== ''NpT''-simulation with Langevin thermostat ===
=== ''NpT''-simulation with Langevin thermostat ===
The Langevin dynamics in the isobaric-isothermal ensemble is calculated by the method of Parrinello and Rahman{{cite|parrinello:prl:1980}}{{cite|parrinello:jap:1981}} (see {{TAG|NpT ensemble}} for more descriptions) combined with a {{TAG|Langevin thermostat}}.


In the method of Parrinello and Rahman<ref name="Parrinello80"/><ref name="Parrinello81"/>, the equations of motion for ionic and lattice degrees-of-freedom are derived from the following Lagrangian:
*To run an NpT-simulation (Parinello-Rahman dynamics) with a Langevin thermostat, one has to:
:<math>
{{Template:NpT_Langevin_thermostat_recipe}}
{\mathcal{L}}(s,h,\dot{s},\dot{h})=\frac{1}{2}\sum_i^N m_i \dot{s_i}^t\,G \dot{s_i}-V(s,h) +\frac{1}{2}W\,\mathrm{Tr}(\dot{h}^t \dot{h}) - p_\mathrm{ext}\Omega,
</math>
where ''s<sub>i</sub>'' is a position vector in fractional coordinates for atom ''i'', ''h'' is the matrix formed by lattice vectors, the tensor ''G'' is defined as ''G''=''h''<sup>t</sup>''h'', ''p''<sub>ext</sub> is the external pressure, &Omega; is the cell volume (&Omega;=det ''h''), and ''W'' is a constant with the dimensionality of mass. Integrating equations of motion based on the Parrinello-Rahman Lagrangian generates an ''NpH'' ensemble, where the enthalpy <math>H=E+p_\mathrm{ext}\Omega</math> is the constant of motion. The Parrinello-Rahman method can be combined with a Langevin thermostat<ref name="Allen91"/> to generate an ''NpT''-ensemble.
 
*To run an ''NpT''-simulation (Parinello-Rahman dynamics) with a Langevin thermostat, one has to:
#Set the standard MD-related tags: {{TAG|IBRION}}=0, {{TAG|TEBEG}}, {{TAG|POTIM}}, and {{TAG|NSW}}
#Set {{TAG|ISIF}}=3 to allow for relaxation of the cell volume and shape. At the moment, dynamics with ''fixed volume+variable shape'' ({{TAG|ISIF}}=4) or ''fixed shape+variable volume'' ({{TAG|ISIF}}=7) are not available.
#Set {{TAG|MDALGO}}=3 to invoke the Langevin thermostat
#Specify friction coefficients for all species in the {{FILE|POSCAR}} file, by means of the {{TAG|LANGEVIN_GAMMA}}-tag.
#Specify a separate set of friction coefficient for the lattice degrees-of-freedom, using the {{TAG|LANGEVIN_GAMMA_L}}-tag.
#Set a mass for the lattice degrees-of-freedom, using the {{TAG|PMASS}}-tag.
#Optionally, one may define an external pressure (in kB), by means of the {{TAG|PSTRESS}}-tag.


The temperatures listed in the {{FILE|OSZICAR}} are computed using the kinetic energy including contribution from both atomic and lattice degrees of freedom. The external pressure for a simulation can be computed as one third of the trace of the stress-tensor corrected for kinetic contributions, listed in the {{FILE|OUTCAR}} file. This can be achieved, ''e.g.'' using the following command:
Note that the advanced simulation protocols, such metadynamics, or geometric constraints and restraints are available within the {{TAG|MDALGO}}=3 setting.  
<code lang="text">
grep "Total+kin" OUTCAR| awk 'BEGIN {p=0.} {p+=($2+$3+$4)/3.} END {print "pressure (kB):",p}'
</code>
 
Important: In Parinello-Rahman<ref name="Parrinello80"/><ref name="Parrinello81"/> dynamics (''NpT''), the stress tensor is used to define forces on lattice degrees-of-freedom. In order to achieve a reasonable quality of sampling (and to avoid numerical problems), it is essential to eliminate [[Energy_vs_volume_Volume_relaxations_and Pulay_stress|Pulay stress]]. Unfortunately, this usually requires a rather large value of {{TAG|ENCUT}}. {{TAG|PREC}}=low, frequently used in ''NVT''-MD is not recommended for molecular dynamics with variable cell volume.


=== Stochastic boundary conditions ===
=== Stochastic boundary conditions ===
Line 158: Line 190:
#Set {{TAG|ISIF}}=2
#Set {{TAG|ISIF}}=2
#Set {{TAG|MDALGO}}=3 to invoke the Langevin thermostat
#Set {{TAG|MDALGO}}=3 to invoke the Langevin thermostat
#Prepare the {{FILE|POSCAR}} file in such a way that the Newtonian and Langevin atoms are treated as different species (even if they are chemically identical). In your {{FILE|POSCAR}}, use [[Selective Dynamics]] and the corresponding logical flags to define the frozen and moveable atoms.
#Prepare the {{FILE|POSCAR}} file in such a way that the Newtonian and Langevin atoms are treated as different species (even if they are chemically identical). In your {{FILE|POSCAR}}, use "selective dynamics" and the corresponding logical flags to define the frozen and moveable atoms.
#Specify friction coefficients &gamma;, for all species in the {{FILE|POSCAR}} file, by means of the {{TAG|LANGEVIN_GAMMA}}-tag: set the friction coefficients to 0 for all fixed and Newtonian atoms, and choose a proper &gamma; for the Langevin atoms.
#Specify friction coefficients &gamma;, for all species in the {{FILE|POSCAR}} file, by means of the {{TAG|LANGEVIN_GAMMA}}-tag: set the friction coefficients to 0 for all fixed and Newtonian atoms, and choose a proper &gamma; for the Langevin atoms.


Line 174: Line 206:
''i.e.'', for all non-Langevin atoms, &gamma; should be set to zero.
''i.e.'', for all non-Langevin atoms, &gamma; should be set to zero.


== {{TAG|MDALGO}}=11: Biased-MD and metadynamics with Andersen thermostat ==
<div id="multiAnderson"></div>
=== Andersen thermostat ===
For a short description of the Andersen thermostat see [[#Andersen|its section]] under {{TAG|MDALGO}}=1.
 
<div id="Metadynamics"></div>
=== Metadynamics ===
In metadynamics,<ref name="Laio02"/><ref name="Iannuzzi03"/> the bias potential
that acts on a selected number of geometric parameters (collective variables) &xi;={&xi;<sub>1</sub>, &xi;<sub>2</sub>, ...,&xi;<sub>''m''</sub>} is constructed on-the-fly during the simulation. The Hamiltonian for the metadynamics <math>\tilde{H}(q,p)</math> can be written as:
:<math>
\tilde{H}(q,p,t) = H(q,p) + \tilde{V}(t,\xi),
</math>
where <math>H(q,p)</math> is the Hamiltonian for the original (unbiased) system, and <math>\tilde{V}(t,\xi)</math> is the time-dependent bias potential. The latter term is usually defined as a sum of Gaussian hills with height ''h'' and width ''w'':
:<math>
\tilde{V}(t,\xi) = h \sum_{i=1}^{\lfloor t/t_G \rfloor} \exp{\left\{ -\frac{|\xi^{(t)}-\xi^{(i \cdot t_G)}|^2}{2
w^2} \right\}}.
</math>
In practice, <math>\tilde{V}(t,\xi)</math> is updated by adding a new Gaussian with a time increment ''t''<sub>G</sub>, which is typically one or two orders of magnitude greater than the time step used in the MD simulation.
 
In the limit of infinite simulation time, the bias potential is related to the free energy via:
:<math>
A(\xi) = - \lim_{t \to \infty} \tilde{V}(t,\xi) + const.
</math>
Practical hints as how to adjust the parameters used in metadynamics (''h'', ''w'', ''t''<sub>G</sub>) are given in Refs.<ref name="Ensing05"/><ref name="Laio05"/>.
 
The error estimation in free-energy calculations with metadynamics is discussed in Ref.<ref name="Laio05"/>.
 


* For a metadynamics run with Andersen thermostat, one has to:
== {{TAG|MDALGO}}=4: Nose Hoover chains thermostat ==
For a description of Nose Hoover chains theromastat, see {{TAG|NHC thermostat}}.
*To run an ''NVT''-simulation with a Nose Hoover chains thermostat, one has to:
#Set the standard MD-related tags: {{TAG|IBRION}}=0, {{TAG|TEBEG}}, {{TAG|POTIM}}, and {{TAG|NSW}}
#Set the standard MD-related tags: {{TAG|IBRION}}=0, {{TAG|TEBEG}}, {{TAG|POTIM}}, and {{TAG|NSW}}
#Set {{TAG|MDALGO}}=11, and choose an appropriate setting for {{TAG|ANDERSEN_PROB}}
#Set {{TAG|ISIF}}=2
#Set the parameters {{TAG|HILLS_H}}, {{TAG|HILLS_W}}, and {{TAG|HILLS_BIN}}
#Set {{TAG|MDALGO}}=4 to invoke the thermostat, select the number of thermostats in chain via {{TAG|NHC_NCHAINS}} and choose an appropriate setting for the thermostat parameter {{TAG|NHC_PERIOD}}.
#Define collective variables in the {{FILE|ICONST}}-file, and set the {{TAG|STATUS}} parameter for the collective variables to 5
Note that the advanced simulation protocols, such metadynamics, or geometric constraints and restraints are available within the {{TAG|MDALGO}}=4 setting.  
#If needed, define the bias potential in the {{FILE|PENALTYPOT}}-file
 
The actual time-dependent bias potential is written to the {{FILE|HILLSPOT}}-file, which is updated after adding a new Gaussian. At the beginning of the simulation, VASP attempts to read the initial bias potential from the {{FILE|PENALTYPOT}}-file. For the continuation of a metadynamics run, copy {{FILE|HILLSPOT}} to {{FILE|PENALTYPOT}}. The values of all collective variables for each MD step are listed in {{FILE|REPORT}}-file, check the lines after the string <tt>Metadynamics</tt>.
 
<div id="BiasedMD"></div>
 
=== Biased molecular dynamics ===
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 only on a selected internal parameter 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 such as umbrella sampling<ref name="Torrie77"/> use a bias potential to enhance sampling of &xi; in regions with low ''P''(&xi;<sub>''i''</sub>) such as transition regions of chemical
reactions.
The correct distributions are recovered afterwards using the equation for <math>\langle A \rangle_{H}</math> above.
 
A more detailed description of the method can be found in Ref.<ref name="FrenkelSmit"/>.
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.


== {{TAG|MDALGO}}=5: Canonical sampling through velocity-rescaling (CSVR thermostat) ==
{{NB|mind|This option is available as of VASP 6.4.3.}}


* For a biased molecular dynamics run with Andersen thermostat, one has to:
For a description of CSVR theromastat, see {{TAG|CSVR thermostat}}.
*To run an ''NVT''-simulation with a CSVR thermostat, one has to:
#Set the standard MD-related tags: {{TAG|IBRION}}=0, {{TAG|TEBEG}}, {{TAG|POTIM}}, and {{TAG|NSW}}
#Set the standard MD-related tags: {{TAG|IBRION}}=0, {{TAG|TEBEG}}, {{TAG|POTIM}}, and {{TAG|NSW}}
#Set {{TAG|MDALGO}}=11, and choose an appropriate setting for {{TAG|ANDERSEN_PROB}}
#Set {{TAG|ISIF}}=2
#In order to avoid updating of the bias potential, set {{TAG|HILLS_BIN}}={{TAG|NSW}}
#Set {{TAG|MDALGO}}=5 to invoke the thermostat and choose an appropriate setting for the thermostat parameter {{TAG|CSVR_PERIOD}}.
#Define collective variables in the {{FILE|ICONST}}-file, and set the <tt>STATUS</tt> parameter for the collective variables to 5
Note that the advanced simulation protocols, such metadynamics, or geometric constraints and restraints are available within the {{TAG|MDALGO}}=5 setting.
#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>.
== {{TAG|MDALGO}}=13: Multiple Andersen thermostats ==
 
== {{TAG|MDALGO}}=21: Biased-MD and metadynamics with Nose-Hoover Thermostat ==
Biased-molecular- or [[Metadynamics|meta-dynamics]] with Nose-Hoover Thermostat ({{TAG|SMASS}} needs to be specified in the {{FILE|INCAR}} file).
 
=== Metadynamics ===
For a short description of metadynamics see [[#Metadynamics|its section]] under {{TAG|MDALGO}}=11.
 
* For a metadynamics 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}}=21, and choose an appropriate setting for {{TAG|SMASS}}
#Set the parameters {{TAG|HILLS_H}}, {{TAG|HILLS_W}}, and {{TAG|HILLS_BIN}}
#Define collective variables in the {{FILE|ICONST}}-file, and set the <tt>STATUS</tt> parameter for the collective variables to 5
#If needed, define the bias potential in the {{FILE|PENALTYPOT}}-file
 
The actual time-dependent bias potential is written to the {{FILE|HILLSPOT}}-file, which is updated after adding a new Gaussian. At the beginning of the simulation, VASP attempts to read the initial bias potential from the {{FILE|PENALTYPOT}}-file. For the continuation of a metadynamics run, copy {{FILE|HILLSPOT}} to {{FILE|PENALTYPOT}}. The values of all collective variables for each MD step are listed in {{FILE|REPORT}}-file, check the lines after the string <tt>Metadynamics</tt>.
 
=== Biased molecular dynamics ===
For a short description of biased molecular dynamics see [[#BiasedMD|its section]] under {{TAG|MDALGO}}=11.
 
* For a biased molecular dynamics 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}}=21, and choose an appropriate setting for {{TAG|SMASS}}
#In order 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
 
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>.
 
<div id="multiAnderson"></div>
== {{TAG|MDALGO}}=13: Multiple Anderson thermostats ==
Up to three user-defined atomic subsystems may be coupled with independent Andersen thermostats<ref name="Andersen80"/> (see remarks under {{TAG|MDALGO}}=1 as well).
Up to three user-defined atomic subsystems may be coupled with independent Andersen thermostats<ref name="Andersen80"/> (see remarks under {{TAG|MDALGO}}=1 as well).
The {{FILE|POSCAR}} file must be organized such that the positions of atoms of subsystem ''i+1'' are defined after those for the subsystem ''i'', and the following flags must be set by the user:  
The {{FILE|POSCAR}} file must be organized such that the positions of atoms of subsystem ''i+1'' are defined after those for the subsystem ''i'', and the following flags must be set by the user:  
*{{TAG|NSUBSYS}}=[int array]  
*{{TAG|NSUBSYS}}=[int array]  
:Define the last atom for each subsystem (two or three values must be supplied). For instance, if total of 20 atoms is defined in the {{FILE|POSCAR}} file, and the initial 10 atoms belong to the subsystem 1, the next 7 atoms to the subsystem 2, and the last 3 atoms to the subsystem 3, {{TAG|NSUBSYS}} should be defined as follows:
:Define the last atom for each subsystem (two or three values must be supplied). For instance, if total of 20 atoms is defined in the {{FILE|POSCAR}} file, and the initial 10 atoms belong to the subsystem 1, the next 7 atoms to the subsystem 2, and the last 3 atoms to the subsystem 3, {{TAG|NSUBSYS}} should be defined as follows:
:<code lang="text">
 
:{{TAG|NSUBSYS}}= 10 17 20
:{{TAG|NSUBSYS}}= 10 17 20
:</code>
 
:Note that the last number in the previous example is actually redundant (clearly the last three atoms belong to the last subsystem) and does not have to be user-supplied.
:Note that the last number in the previous example is actually redundant (clearly the last three atoms belong to the last subsystem) and does not have to be user-supplied.


Line 336: Line 272:
<references>
<references>
<ref name="Andersen80">[http://dx.doi.org/10.1063/1.439486 H. C. Andersen, J. Chem. Phys. 72, 2384 (1980).]</ref>
<ref name="Andersen80">[http://dx.doi.org/10.1063/1.439486 H. C. Andersen, J. Chem. Phys. 72, 2384 (1980).]</ref>
<ref name="Ryckaert77">[http://dx.doi.org/10.1016/0021-9991(77)90098-5 J. P. Ryckaert, G. Ciccotti, and H. J. C. Berendsen, J. Comp. Phys. 23, 327 (1977).]</ref>
<ref name="Carter89">[http://dx.doi.org/10.1016/S0009-2614(89)87314-2 E. A. Carter, G. Ciccotti, J. T. Hynes, and R. Kapral, Chem. Phys. Lett. 156, 472 (1989).]</ref>
<ref name="Otter00">[http://dx.doi.org/10.1080/00268970009483348 W. K. Den Otter and W. J. Briels, Mol. Phys. 98, 773 (2000).]</ref>
<ref name="Darve02">[http://dx.doi.org/10.1080/08927020211975 E. Darve, M. A. Wilson, and A. Pohorille, Mol. Simul. 28, 113 (2002).]</ref>
<ref name="Fleurat05">[http://dx.doi.org/10.1063/1.1948367 P. Fleurat-Lessard and T. Ziegler, J. Chem. Phys. 123, 084101 (2005).]</ref>
<ref name="Allen91">M. P. Allen and D. J. Tildesley, ''Computer simulation of liquids'', Oxford university press: New York, 1991.</ref>
<ref name="Parrinello80">[http://dx.doi.org/10.1103/PhysRevLett.45.1196 M. Parrinello and A. Rahman, Phys. Rev. Lett. 45, 1196 (1980).]</ref>
<ref name="Parrinello81">[http://dx.doi.org/10.1063/1.328693 M. Parrinello and A. Rahman, J. Appl. Phys. 52, 7182 (1981).]</ref>
<ref name="Toton10">[http://dx.doi.org/10.1088/0953-8984/22/7/074205 D. Toton, C. D. Lorenz, N. Rompotis, N. Martsinovich, and L. Kantorovich, J. Phys.: Condens. Matter 22, 074205 (2010).]</ref>
<ref name="Toton10">[http://dx.doi.org/10.1088/0953-8984/22/7/074205 D. Toton, C. D. Lorenz, N. Rompotis, N. Martsinovich, and L. Kantorovich, J. Phys.: Condens. Matter 22, 074205 (2010).]</ref>
<ref name="Kantorovich08">[http://dx.doi.org/10.1103/PhysRevB.78.094305 L. Kantorovich and N. Rompotis, Phys. Rev. B 78, 094305 (2008).]</ref>
<ref name="Kantorovich08">[http://dx.doi.org/10.1103/PhysRevB.78.094305 L. Kantorovich and N. Rompotis, Phys. Rev. B 78, 094305 (2008).]</ref>
<ref name="Laio02">[http://dx.doi.org/10.1073/pnas.202427399 A. Laio and M. Parrinello, Proc. Natl. Acad, Sci. USA 99, 12562 (2002).]</ref>
<ref name="Iannuzzi03">[http://dx.doi.org/10.1103/PhysRevLett.90.238302 M. Iannuzzi, A. Laio, and M. Parrinello, Phys. Rev. Lett. 90, 238302 (2003).]</ref>
<ref name="Ensing05">[http://dx.doi.org/10.1021/jp045571i B. Ensing, A. Laio, M. Parrinello, and M. L. Klein, J. Phys. Chem. B 109, 6676 (2005).]</ref>
<ref name="Laio05">[http://dx.doi.org/10.1021/jp045424k A. Laio, A. Rodriguez-Fortea, F. L. Gervasio, M. Ceccarelli, and M. Parrinello, J. Phys. Chem. B 109, 6714 (2005).]</ref>
<ref name="Torrie77">[http://dx.doi.org/10.1016/0021-9991(77)90121-8 G. M. Torrie and J. P. Valleau, J. Comp. Phys. 23, 187 (1977).]</ref>
<ref name="FrenkelSmit">D. Frenkel and B. Smit, ''Understanding molecular simulations: from algorithms to applications'', Academic Press: San Diego, 2002.</ref>
</references>
</references>
----
----
[[The_VASP_Manual|Contents]]


[[Category:Molecular Dynamics]][[Category:Howto]]
[[Category:INCAR tag]][[Category:Molecular dynamics]][[Category:Advanced molecular-dynamics sampling]][[Category:Howto]]

Latest revision as of 09:27, 15 October 2024

MDALGO = 0 | 1 | 2 | 3 | 4 | 5 | 11 | 21 | 13
Default: MDALGO = 0 

Description: MDALGO specifies the molecular-dynamics-simulation protocol (in case IBRION=0 and VASP was compiled with -Dtbdyn).

MDALGO=0: Standard molecular dynamics

Should provide the same results as MDALGO=2 (Nose-Hoover thermostat). The difference is that it is a different implementation and works also without the precompiler flag -Dtbdyn.

NVE ensemble

  • To perform a calculation in the NVE ensemble:
  1. Set the standard MD-related tags: IBRION=0, TEBEG, POTIM, and NSW.
  2. Set MDALGO=0 and set SMASS=-3.

MDALGO=1: Andersen thermostat

For the description of the Andersen thermostat see: Andersen thermostat.

The Andersen thermostat is only available for the NVT ensemble.

Standard molecular dynamics in

  • For a standard molecular dynamics run with Anderson thermostat, one has to:
  1. Set the standard MD-related tags: IBRION=0, TEBEG, POTIM, and NSW.
  2. Set MDALGO=1, and choose an appropriate setting for ANDERSEN_PROB

Constrained molecular dynamics

For a description of constrained molecular dynamics see Constrained molecular dynamics.

  • For a constrained 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, and choose an appropriate setting for ANDERSEN_PROB
  3. Define geometric constraints in the ICONST-file, and set the STATUS parameter for the constrained coordinates to 0
  4. When the free-energy gradient is to be computed, set LBLUEOUT=.TRUE.

Slow-growth approach

For a description of slow-growth approach see Slow-growth approach.

  • For a slow-growth simulation, one has to perform a calcualtion very similar to Constrained molecular dynamics but additionally the transformation velocity-related INCREM-tag for each geometric parameter with STATUS=0 has to be specified:
  1. Set the standard MD-related tags: IBRION=0, TEBEG, POTIM, and NSW
  2. Set MDALGO=1, and choose an appropriate setting for ANDERSEN_PROB
  3. Define geometric constraints in the ICONST-file, and set the STATUS parameter for the constrained coordinates to 0
  4. When the free-energy gradient is to be computed, set LBLUEOUT=.TRUE.
  1. Specify the transformation velocity-related INCREM-tag for each geometric parameter with STATUS=0.

Monitoring geometric parameters

Geometric parameters with STATUS = 7 in the ICONST-file are monitored during the MD simulation. The corresponding values are written onto the REPORT-file, for each MD step, after the lines following the string Monit_coord.

Sometimes it is desirable to terminate the simulation if all values of monitored parameters get larger that some predefined upper and/or lower limits. These limits can be set by the user by means of the VALUE_MAX and VALUE_MIN-tags.

  • To monitor geometric parameters during an MD run:
  1. Set the standard MD-related tags: IBRION=0, TEBEG, POTIM, and NSW.
  2. Set MDALGO=1, and choose an appropriate setting for ANDERSEN_PROB.
  3. Define geometric constraints in the ICONST-file, and set the STATUS parameter for the constrained coordinates to 7.
  4. Optionally, set the upper and/or lower limits for the coordinates, by means of the VALUE_MAX and VALUE_MIN-tags, respectively.

Metadynamics

For a description of metadynamics see Metadynamics.

  • For a metadynamics run with Andersen thermostat, one has to:
  1. Set the standard MD-related tags: IBRION=0, TEBEG, POTIM, and NSW.
  2. Set MDALGO=1 (or MDALGO=11 in VASP 5.x), and choose an appropriate setting for ANDERSEN_PROB.
  3. Set the parameters HILLS_H, HILLS_W, and HILLS_BIN.
  4. Define collective variables in the ICONST-file, and set the STATUS parameter for the collective variables to 5.
  5. If needed, define the bias potential in the PENALTYPOT-file.

The actual time-dependent bias potential is written to the HILLSPOT-file, which is updated after adding a new Gaussian. At the beginning of the simulation, VASP attempts to read the initial bias potential from the PENALTYPOT-file. For the continuation of a metadynamics run, copy HILLSPOT to PENALTYPOT. The values of all collective variables for each MD step are listed in REPORT-file, check the lines after the string Metadynamics.

Biased molecular dynamics

For a description of biased molecular dynamics see Biased molecular dynamics.

  • 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 (or 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 if necessary.

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

Special case: NVE ensemble

NVE ensemble calculations can be also run by selecting the Anderson thermostat and setting the update collision probability (ANDERSEN_PROB) to zero.

  • To run an NVE ensemble:
  1. Set the standard MD-related tags: IBRION=0, TEBEG, POTIM, and NSW.
  2. Set MDALGO=1 and ANDERSEN_PROB=0.0.

MDALGO=2: Nose-Hoover thermostat

For the description of the Nose-Hoover thermostat see: Nose-Hoover thermostat.

The Nose-Hoover thermostat is currently only available for the NVT ensemble.

Standard molecular dynamics

  • For a standard 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, and choose an appropriate setting for SMASS.

Constrained molecular dynamics

For a description of constrained molecular dynamics see Constrained molecular dynamics.

  • For a constrained 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, and choose an appropriate setting for SMASS.
  3. Define geometric constraints in the ICONST-file, and set the STATUS parameter for the constrained coordinates to 0.
  4. When the free-energy gradient is to be computed, set LBLUEOUT=.TRUE.

Slow-growth approach

For a description of slow-growth approach see Slow-growth approach.

  • For a slow-growth approach 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, and choose an appropriate setting for SMASS
  3. Define geometric constraints in the ICONST-file, and set the STATUS parameter for the constrained coordinates to 0
  4. When the free-energy gradient is to be computed, set LBLUEOUT=.TRUE.
  1. Specify the transformation velocity-related INCREM-tag for each geometric parameter with STATUS=0

Monitoring geometric parameters

Geometric parameters with STATUS = 7 in the ICONST-file are monitored during the MD simulation. The corresponding values are written onto the REPORT-file, for each MD step, after the lines following the string Monit_coord.

Sometimes it is desirable to terminate the simulation if all values of monitored parameters get larger that some predefined upper and/or lower limits. These limits can be set by the user by means of the VALUE_MAX and VALUE_MIN-tags.

To monitor geometric parameters during an MD run:

  1. Set the standard MD-related tags: IBRION=0, TEBEG, POTIM, and NSW
  2. Set MDALGO=2, and choose an appropriate setting for SMASS
  3. Define geometric constraints in the ICONST-file, and set the STATUS parameter for the constrained coordinates to 7
  4. Optionally, set the upper and/or lower limits for the coordinates, by means of the VALUE_MAX and VALUE_MIN-tags, respectively.

Metadynamics

For a description of metadynamics see Metadynamics.

  • For a metadynamics 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 (or MDALGO=21 in VASP 5.x), and choose an appropriate setting for SMASS.
  3. Set the parameters HILLS_H, HILLS_W, and HILLS_BIN.
  4. Define collective variables in the ICONST-file, and set the STATUS parameter for the collective variables to 5.
  5. If needed, define the bias potential in the PENALTYPOT-file.

The actual time-dependent bias potential is written to the HILLSPOT-file, which is updated after adding a new Gaussian. At the beginning of the simulation, VASP attempts to read the initial bias potential from the PENALTYPOT-file. For the continuation of a metadynamics run, copy HILLSPOT to PENALTYPOT. The values of all collective variables for each MD step are listed in REPORT-file, check the lines after the string Metadynamics.

Biased molecular dynamics

For a description of biased molecular dynamics see Biased molecular dynamics.

  • 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 (or 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.

MDALGO=3: Langevin thermostat

For the description of the Langevin thermostat see: Langevin thermostat.

The Langevin thermostat is available for NVT and NpT ensembles.

NVT-simulation with Langevin thermostat

  • To run an NVT-simulation with a Langevin thermostat, one has to:
  1. Set the standard MD-related tags: IBRION=0, TEBEG, POTIM, and NSW.
  2. Set ISIF=2.
  3. Set MDALGO=3 to invoke the Langevin thermostat.
  4. Specify friction coefficients for all species in the POSCAR file, by means of the LANGEVIN_GAMMA-tag.

NpT-simulation with Langevin thermostat

The Langevin dynamics in the isobaric-isothermal ensemble is calculated by the method of Parrinello and Rahman[1][2] (see NpT ensemble for more descriptions) combined with a Langevin thermostat.

  • To run an NpT-simulation (Parinello-Rahman dynamics) with a Langevin thermostat, one has to:
  1. Set the standard MD-related tags: IBRION=0, TEBEG, POTIM, and NSW.
  2. Set ISIF=3 to allow for relaxation of the cell volume and shape. At the moment, dynamics with fixed volume+variable shape (ISIF=4) or fixed shape+variable volume (ISIF=7) are not available.
  3. Set MDALGO=3 to invoke the Langevin thermostat.
  4. Specify friction coefficients for all species in the POSCAR file, by means of the LANGEVIN_GAMMA-tag.
  5. Specify a separate set of friction coefficient for the lattice degrees-of-freedom, using the LANGEVIN_GAMMA_L-tag.
  6. Set a mass for the lattice degrees-of-freedom, using the PMASS-tag.
  7. Optionally, one may define an external pressure (in kB), by means of the PSTRESS-tag.

Note that the advanced simulation protocols, such metadynamics, or geometric constraints and restraints are available within the MDALGO=3 setting.

Stochastic boundary conditions

In some cases it is desirable to study approach of initially non-equilibrium system to equilibrium. Examples of such simulations include the impact problems when a particle with large kinetic energy hits a surface or calculation of friction force between two surfaces sliding with respect to each other. As shown by Toton et al.[3], this type of problems can be studied using the stochastic boundary conditions (SBC) derived from the generalized Langevin equation by Kantorovich and Rompotis.[4] In this approach, the system of interest is divided into three regions: (a) fixed atoms, (b) the internal (Newtonian) atoms moving according to Newtonian dynamics, and (c) a buffer region of Langevin atoms (i.e., atoms governed by Langevin equations of motion) located between (a) and (b).

The role of the Langevin atoms is to dissipate heat, while the fixed atoms are needed solely to create the correct potential well for the Langevin atoms to move in. The Newtonian region should include all atoms relevant to the process under study: in the case of the impact problem, for instance, the Newtonian region should contain atoms of the molecule hitting the surface and several uppermost layers of the material forming the surface. Performing molecular dynamics with such a setup guarantees that the system (possibly out of equilibrium initially) arrives at the appropriate canonical distribution.

  • To run a simulation with stochastic boundary conditions, one has to:
  1. Set the standard MD-related tags: IBRION=0, TEBEG, POTIM, and NSW
  2. Set ISIF=2
  3. Set MDALGO=3 to invoke the Langevin thermostat
  4. Prepare the POSCAR file in such a way that the Newtonian and Langevin atoms are treated as different species (even if they are chemically identical). In your POSCAR, use "selective dynamics" and the corresponding logical flags to define the frozen and moveable atoms.
  5. Specify friction coefficients γ, for all species in the POSCAR file, by means of the LANGEVIN_GAMMA-tag: set the friction coefficients to 0 for all fixed and Newtonian atoms, and choose a proper γ for the Langevin atoms.

Practical example

Consider a system consisting of 16 hydrogen and 48 silicon atoms. Suppose that eight silicon atoms are considered to be Langevin atoms and the remaining 32 Si atoms are either fixed or Newtonian atoms. The Langevin and Newtonian (or fixed) atoms should be considered as different species, i.e., the POSCAR-file should contain the line like this: 

Si H Si
40 16 8

As only the final eight Si atoms are considered to be Langevin atoms, the INCAR-file should contain the following line defining the friction coefficients: 

LANGEVIN_GAMMA = 0.0   0.0   10.0

i.e., for all non-Langevin atoms, γ should be set to zero.

MDALGO=4: Nose Hoover chains thermostat

For a description of Nose Hoover chains theromastat, see NHC thermostat.

  • To run an NVT-simulation with a Nose Hoover chains thermostat, one has to:
  1. Set the standard MD-related tags: IBRION=0, TEBEG, POTIM, and NSW
  2. Set ISIF=2
  3. Set MDALGO=4 to invoke the thermostat, select the number of thermostats in chain via NHC_NCHAINS and choose an appropriate setting for the thermostat parameter NHC_PERIOD.

Note that the advanced simulation protocols, such metadynamics, or geometric constraints and restraints are available within the MDALGO=4 setting.

MDALGO=5: Canonical sampling through velocity-rescaling (CSVR thermostat)

Mind: This option is available as of VASP 6.4.3.

For a description of CSVR theromastat, see CSVR thermostat.

  • To run an NVT-simulation with a CSVR thermostat, one has to:
  1. Set the standard MD-related tags: IBRION=0, TEBEG, POTIM, and NSW
  2. Set ISIF=2
  3. Set MDALGO=5 to invoke the thermostat and choose an appropriate setting for the thermostat parameter CSVR_PERIOD.

Note that the advanced simulation protocols, such metadynamics, or geometric constraints and restraints are available within the MDALGO=5 setting.

MDALGO=13: Multiple Andersen thermostats

Up to three user-defined atomic subsystems may be coupled with independent Andersen thermostats[5] (see remarks under MDALGO=1 as well). The POSCAR file must be organized such that the positions of atoms of subsystem i+1 are defined after those for the subsystem i, and the following flags must be set by the user:

Define the last atom for each subsystem (two or three values must be supplied). For instance, if total of 20 atoms is defined in the POSCAR file, and the initial 10 atoms belong to the subsystem 1, the next 7 atoms to the subsystem 2, and the last 3 atoms to the subsystem 3, NSUBSYS should be defined as follows:
NSUBSYS= 10 17 20
Note that the last number in the previous example is actually redundant (clearly the last three atoms belong to the last subsystem) and does not have to be user-supplied.
Simulation temperature for each subsystem
Collision probability for atoms in each subsystem. Only the values 0≤PSUBSYS≤1 are allowed.

Related Tags and Sections

IBRION, ISIF, SMASS, ANDERSEN_PROB, RANDOM_SEED, LBLUEOUT, SHAKETOL, SHAKEMAXITER, HILLS_H, HILLS_W, HILLS_BIN, INCREM, VALUE_MIN, VALUE_MAX, LANGEVIN_GAMMA, LANGEVIN_GAMMA_L, PMASS, NSUBSYS, TSUBSYS, PSUBSYS, ICONST, PENALTYPOT, HILLSPOT, REPORT

Examples that use this tag

References

  1. ↑ M. Parrinello and A. Rahman, Phys. Rev. Lett. 45, 1196 (1980).
  2. ↑ M. Parrinello and A. Rahman, J. Appl. Phys. 52, 7182 (1981).
  3. ↑ D. Toton, C. D. Lorenz, N. Rompotis, N. Martsinovich, and L. Kantorovich, J. Phys.: Condens. Matter 22, 074205 (2010).
  4. ↑ L. Kantorovich and N. Rompotis, Phys. Rev. B 78, 094305 (2008).
  5. ↑ H. C. Andersen, J. Chem. Phys. 72, 2384 (1980).
Retrieved from "https://www.vasp.at/wiki/index.php?title=MDALGO&oldid=26380"