ICODE-ADC/lammps
4
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

updates to hyper doc pages

Browse Source
This commit is contained in:
Steve Plimpton
2018-10-11 11:25:19 -06:00
parent 56598fcd0b
commit 4a6f088c0b
2 changed files with 192 additions and 103 deletions
Download Patch File Download Diff File Expand all files Collapse all files

154
doc/src/fix_hyper_global.txt
View File

@ -33,9 +33,9 @@ dynamics of the system in a manner that effectively accelerates time.
This is in contrast to the "fix hyper/local"_fix_hyper_local.html This is in contrast to the "fix hyper/local"_fix_hyper_local.html
command, which can be user to perform a local hyperdynamics (LHD) command, which can be user to perform a local hyperdynamics (LHD)
simulation, by adding a local bias potential to multiple pairs of simulation, by adding a local bias potential to multiple pairs of
atoms at each timestep. GHD can speed-up a small simulation with up atoms at each timestep. GHD can time accelerate a small simulation
to a few 100 atoms. For larger systems, LHD is needed to achieve good with up to a few 100 atoms. For larger systems, LHD is needed to
speed-ups. achieve good time acceleration.
For a system that undergoes rare transition events, where one or more For a system that undergoes rare transition events, where one or more
atoms move over an energy barrier to a new potential energy basin, the atoms move over an energy barrier to a new potential energy basin, the
@ -43,23 +43,24 @@ effect of the bias potential is to induce more rapid transitions.
This can lead to a dramatic speed-up in the rate at which events This can lead to a dramatic speed-up in the rate at which events
occurs, without altering their relative frequencies, thus leading to occurs, without altering their relative frequencies, thus leading to
an overall increase in the elapsed real time of the simulation as an overall increase in the elapsed real time of the simulation as
compared to simulating for the same number of timesteps with normal compared to running for the same number of timesteps with normal MD.
MD. See the "hyper"_hyper.html doc page for a more general discussion See the "hyper"_hyper.html doc page for a more general discussion of
of hyperdynamics and citations that explain both GHD and LHD. hyperdynamics and citations that explain both GHD and LHD.
The equations and logic used by this fix to perform GHD are as The equations and logic used by this fix to perform GHD follow the
follows. This follows the description of GHD given in description given in "(Voter2013)"_#Voter2013ghd. The bond-boost form
"(Voter2013)"_#Voter2013ghd. The bond-boost form of a bias potential of a bias potential for HD is due to Miron and Fichthorn as described
for GHD is due to Miron and Fichthorn as described in in "(Miron)"_#Mironghd. In LAMMPS we use a simplified version of
"(Miron)"_#Mironghd. Here we use a simplified version of bond-boost bond-boost GHD where a single bond in the system is biased at any one
GHD where a single bond is biased at any one timestep. timestep.
Bonds are defined between every pair of I,J atoms whose R0ij distance Bonds are defined between each pair of I,J atoms whose R0ij distance
is less than {cutbond}, when the system is in a quenched state is less than {cutbond}, when the system is in a quenched state
(minimum) energy. Note that these are not "bonds" in a covalent (minimum) energy. Note that these are not "bonds" in a covalent
sense. A bond is simply any pair of atoms that meet the distance sense. A bond is simply any pair of atoms that meet the distance
criterion. {Cutbond} is an argument to this fix; it is discussed more criterion. {Cutbond} is an argument to this fix; it is discussed
below. below. A bond is only formed if one or both of the I.J atoms are in
the specified group.
The current strain of bond IJ (when running dynamics) is defined as The current strain of bond IJ (when running dynamics) is defined as
@ -74,16 +75,16 @@ Vij = Vmax * (1 - (Eij/q)^2) for abs(Eij) < qfactor
= 0 otherwise :pre = 0 otherwise :pre
where the prefactor {Vmax} and the cutoff {qfactor} are arguments to where the prefactor {Vmax} and the cutoff {qfactor} are arguments to
this fix; they are discussed more below. This functional form is an this fix; they are discussed below. This functional form is an
inverse parabola centered at 0.0 with height Vmax and which goes to inverse parabola centered at 0.0 with height Vmax and which goes to
0.0 at +/- qfactor. 0.0 at +/- qfactor.
Let Emax = the maximum of abs(Eij) for all IJ bonds in the system on a Let Emax = the maximum of abs(Eij) for all IJ bonds in the system on a
given timestep. On that step, Vij is added as a bias potential to given timestep. On that step, Vij is added as a bias potential to
only the two atoms IJ in the bond with strain Emax, call it Vij(max). only the single bond with strain Emax, call it Vij(max). Note that
Note that Vij(max) will be 0.0 if Emax >= qfactor on that timestep. Vij(max) will be 0.0 if Emax >= qfactor on that timestep. Also note
Also note that Vij(max) is added to the normal interatomic potential that Vij(max) is added to the normal interatomic potential that is
that is computed between all atoms in the system. computed between all atoms in the system at every step.
The derivative of Vij(max) with respect to the position of each atom The derivative of Vij(max) with respect to the position of each atom
in the Emax bond gives a force Fij(max) acting on the bond as in the Emax bond gives a force Fij(max) acting on the bond as
@ -91,19 +92,21 @@ in the Emax bond gives a force Fij(max) acting on the bond as
Fij(max) = - dVij(max)/dEij = 2 Vmax Eij / qfactor^2 for abs(Eij) < qfactor Fij(max) = - dVij(max)/dEij = 2 Vmax Eij / qfactor^2 for abs(Eij) < qfactor
= 0 otherwise :pre = 0 otherwise :pre
which can be decomposed into an equal and opposite force acting (only) which can be decomposed into an equal and opposite force acting on
on the two I,J atoms in the Emax bond. only the two I,J atoms in the Emax bond.
The boost factor in time for the system is given each timestep I by The time boost factor for the system is given each timestep I by
Bi = exp(beta * Vij(max)) :pre Bi = exp(beta * Vij(max)) :pre
where beta = 1/kTequil, and {Tequil} is the temperature of the system where beta = 1/kTequil, and {Tequil} is the temperature of the system
and an arguement to this fix. Note that Bi >= 1 at every step. and an arguement to this fix. Note that Bi >= 1 at every step.
NOTE: Must use langevin set to Tequil also. LAMMPS does not NOTE: To run GHD, the input script must also use the "fix
check that you do this. langevin"_fix_langevin.html command to thermostat the atoms at the
as maintained by a thermostat for constant NVT dynamics. same {Tequil} as specified by this fix, so that the system is running
constant-temperature (NVT) dynamics. LAMMPS does not check that this
is done.
The elapsed system time t_hyper for a GHD simulation running for {N} The elapsed system time t_hyper for a GHD simulation running for {N}
timesteps is simply timesteps is simply
@ -114,34 +117,59 @@ where dt is the timestep size defined by the "timestep"_timestep.html
command. The effective time acceleration due to GHD is thus t_hyper / command. The effective time acceleration due to GHD is thus t_hyper /
N*dt, where the denominator is elapsed time for a normal MD run. N*dt, where the denominator is elapsed time for a normal MD run.
Note that in global hyperdynamics, the boost factor varies from Note that in GHD, the boost factor varies from timestep to timestep.
timestep to timestep. Likewise, which bond has Emax strain and thus Likewise, which bond has Emax strain and thus which pair of atoms the
which pair of atoms the bias potential is added to, will also vary bias potential is added to, will also vary from timestep to timestep.
from timestep to timestep. This in contrast to local hyperdynamics (LHD) where the boost factor
is an input parameter; see the "fix hyper/local"_fix_hyper_local.html
doc page for details.
:line :line
Here is additional information on the input parameters for GHD.
The {cutbond} argument is the cutoff distance for defining bonds The {cutbond} argument is the cutoff distance for defining bonds
between pairs of nearby atoms. A pair of atoms in their equilibrium, between pairs of nearby atoms. A pair of I,J atoms in their
minimum-energy config, which are a distance Rij < cutbond, are defined equilibrium, minimum-energy configuration, which are separated by a
as a bonded pair. This should normally be roughly ~25% larger than distance Rij < {cutbond}, are flagged as a bonded pair. Setting
the nearest-neighbor distance in a crystalline lattice. {cubond} to be ~25% larger than the nearest-neighbor distance in a
crystalline lattice is a typical choice for solids, so that bonds
exist only between nearest neighbor pairs.
The {qfactor} argument is the limiting strain at which the bias The {qfactor} argument is the limiting strain at which the bias
potential goes to 0.0. potential goes to 0.0. It is dimensionless, so a value of 0.3 means a
bond distance can be up to 30% larger or 30% smaller than the
equilibrium (quenched) R0ij distance and the two atoms in the bond
could still experience a non-zero bias force.
If {qfactor} is set too large, then transitions from one energy basin If {qfactor} is set too large, then transitions from one energy basin
to another are affected because the bias potneial is non-zero at the to another are affected because the bias potneial is non-zero at the
transitions. If {qfactor} is set too small than little boost is transition state (e.g. saddle point). If {qfactor} is set too small
achieved because some the Eij strain of some bond is (nearly) always than little boost is achieved because the Eij strain of some bond in
exceeding {qfactor). Note the system will (nearly) always exceed {qfactor). A value of 0.3 for
{qfactor} is typically a reasonable value.
The {Vmax} argument is the prefactor on the bias potential. The {Vmax} argument is the prefactor on the bias potential. It should
It should be set to a value smaller than be set to a value less than the smallest barrier height for an event
to occur. Otherwise the applied bias potential may be large enough
(when added to the interatomic potential) to produce a local energy
basin with a maxima in the center. This can produce artificial energy
minima in the same basin that trap an atom. Or if {Vmax} is even
larger, it may induce an atom(s) to rapidly transition to another
energy basin. Both cases are "bad dynamics" which violate the
assumptions of GHD that guarantee an accelerated time-accurate
trajectory of the system.
The {Tequil} argument is part of the beta term in the exponential The {Tequil} argument is the temperature at which the system is
factor that determines how much boost is achieved as a function of the simulated; see the comment above about the "fix
bias potential. langevin"_fix_langevin.html thermostatting. It is also part of the
beta term in the exponential factor that determines how much boost is
achieved as a function of the bias potential.
In general, the lower the value of {Tequil} and the higher the value
of {Vmax}, the more boost will be achieved by the GHD algorithm.
:line
[Restart, fix_modify, output, run start/stop, minimize info:] [Restart, fix_modify, output, run start/stop, minimize info:]
@ -152,7 +180,7 @@ The "fix_modify"_fix_modify.html {energy} option is supported by this
fix to add the energy of the bias potential to the the system's fix to add the energy of the bias potential to the the system's
potential energy as part of "thermodynamic output"_thermo_style.html. potential energy as part of "thermodynamic output"_thermo_style.html.
This fix computes a global scalar and global vector of length 10, This fix computes a global scalar and global vector of length 11,
which can be accessed by various "output which can be accessed by various "output
commands"_Section_howto.html#howto_15. The scalar is the magnitude of commands"_Section_howto.html#howto_15. The scalar is the magnitude of
the bias potential (energy units) applied on the current timestep. the bias potential (energy units) applied on the current timestep.
@ -164,30 +192,30 @@ The vector stores the following quantities:
4 = ID of second atom in the max-strain bond 4 = ID of second atom in the max-strain bond
5 = average # of bonds/atom on this step :ul 5 = average # of bonds/atom on this step :ul
6 = fraction of steps with bias = 0.0 during this run 6 = fraction of timesteps with bias = 0.0 during this run
7 = max drift distance of any atom during this run (distance units) 7 = max drift distance of any atom during this run (distance units)
8 = max bond length during this run (distance units) :ul 8 = max bond length during this run (distance units) :ul
9 = cummulative hyper time since fix created (time units) 9 = cummulative hyper time since fix was defined (time units)
10 = cummulative count of event timesteps since fix created 10 = cummulative count of event timesteps since fix was defined
11 = cummulative count of atoms in events since fix created :ul 11 = cummulative count of atoms in events since fix was defined :ul
The first 5 quantities are for the current timestep. The quantities The first 5 quantities are for the current timestep. Quantities 6-8
6-8 are for the current run. The quantities 9-11 are cummulative are for the current hyper run. Quantities 9-11 are cummulative across
across multiple runs (since the fix was defined in the input script). multiple runs (since the fix was defined in the input script).
For value 10, events are checked for by the "hyper"_hyper.html command For value 10, events are checked for by the "hyper"_hyper.html command
once every {Nevent} timesteps. This value is the count of the number once every {Nevent} timesteps. This value is the count of those
of timesteps on which one (or more) events was detected. It is NOT timesteps on which one (or more) events was detected. It is NOT the
the number of distinct events, since more than one physical event may number of distinct events, since more than one physical event may
occur in the same {Nevent} time window. occur in the same {Nevent} time window.
For value 11, each time the "hyper"_hyper.html command checks for an For value 11, each time the "hyper"_hyper.html command checks for an
event, the event compute it uses will flag zero or more atoms as event, it invokes a compute to flag zero or more atoms as
participating in an event. E.g. atoms that have displaced more than participating in an event. E.g. atoms that have displaced more than
some distance from the previous quench state. This value is the count some distance from the previous quench state. Value 11 is the
of the number of atoms participating in any of the events that were cummulative count of the number of atoms participating in any of the
found. events that were found.
The scalar and vector values calculated by this fix are all The scalar and vector values calculated by this fix are all
"intensive". "intensive".
@ -199,8 +227,8 @@ minimization"_minimize.html.
[Restrictions:] [Restrictions:]
This command can only be used if LAMMPS was built with the REPLICA This command can only be used if LAMMPS was built with the REPLICA
package. See the "Build package"_Build_package.html doc package. See the "Build package"_Build_package.html doc page for more
page for more info. info.
[Related commands:] [Related commands:]
@ -208,9 +236,11 @@ page for more info.
[Default:] None [Default:] None
:link(Voter2013) :line
:link(Voter2013ghd)
[(Voter2013)] S. Y. Kim, D. Perez, A. F. Voter, J Chem Phys, 139, [(Voter2013)] S. Y. Kim, D. Perez, A. F. Voter, J Chem Phys, 139,
144110 (2013). 144110 (2013).
:link(Miron) :link(Mironghd)
[(Miron)] R. A. Miron and K. A. Fichthorn, J Chem Phys, 119, 6210 (2003). [(Miron)] R. A. Miron and K. A. Fichthorn, J Chem Phys, 119, 6210 (2003).

141
doc/src/fix_hyper_local.txt
View File

@ -18,7 +18,7 @@ cutbond = max distance at which a pair of atoms is considered bonded (distance u
qfactor = max strain at which bias potential goes to 0.0 (unitless) :l qfactor = max strain at which bias potential goes to 0.0 (unitless) :l
Vmax = estimated height of bias potential (energy units) :l Vmax = estimated height of bias potential (energy units) :l
Tequil = equilibration temperature (temperature units) :l Tequil = equilibration temperature (temperature units) :l
Dcut = min distance between boosted bonds (distance units) :l Dcut = minimum distance between boosted bonds (distance units) :l
alpha = boostostat relaxation time (time units) :l alpha = boostostat relaxation time (time units) :l
boost = desired time boost factor (unitless) :l boost = desired time boost factor (unitless) :l
zero or more keyword/value pairs may be appended :l zero or more keyword/value pairs may be appended :l
@ -40,66 +40,115 @@ fix 1 all hyper/local 1.0 0.3 0.8 300.0 :pre
[Description:] [Description:]
This fix is meant to be used with the "hyper"_hyper.html command to This fix is meant to be used with the "hyper"_hyper.html command to
perform a bond-boost hyperdynamics simulation. The role of this fix perform a bond-boost local hyperdynamics (LHD) simulation. The role
is a select a multiple pairs of atoms within the system at each of this fix is to a select multiple pairs of atoms in the system at
timestep to add a local bias potential to, which will alter their each timestep to add a local bias potential to, which will alter the
dynamics. This is in contrast to the "fix dynamics of the system in a manner that effectively accelerates time.
hyper/global"_fix_hyper_global.html command, which adds a global bias This is in contrast to the "fix hyper/global"_fix_hyper_global.html
potential to a single pair of atoms at each timestep. command, which can be user to perform a global hyperdynamics (GHD)
simulation, by adding a global bias potential to a single pair of
atoms at each timestep. GHD can time accelerate a small simulation
with up to a few 100 atoms. For larger systems, LHD is needed to
achieve good time acceleration.
For a system that undergoes rare transition events, where one or more For a system that undergoes rare transition events, where one or more
atoms move across an energy barrier to a new potential energy basin, atoms move over an energy barrier to a new potential energy basin, the
the effect of the bias potential is to induce more rapid transitions. effect of the bias potential is to induce more rapid transitions.
This can lead to a dramatic speed-up in the rate at which events This can lead to a dramatic speed-up in the rate at which events
occurs, without altering their relative frequencies, thus leading to occurs, without altering their relative frequencies, thus leading to
an overall dramatic speed-up in the effective wall-clock time of the an overall increase in the elapsed real time of the simulation as
simulation. compared to running for the same number of timesteps with normal MD.
See the "hyper"_hyper.html doc page for a more general discussion of
hyperdynamics and citations that explain both GHD and LHD.
Cite various papers. The equations and logic used by this fix to perform LHD follow the
description given in "(Voter2013)"_#Voter2013lhd. The bond-boost form
of a bias potential for HD is due to Miron and Fichthorn as described
in "(Miron)"_#Mironlhd.
The current strain of a bond IJ is defined as To understand this description, you should first read the description
of the GHD algorithm on the "fix hyper/global"_fix_hyper/global.html
doc page. This description of LHD builds on the GHD description.
Eij = (Rij - R0ij) / R0ij The definition of bonds, Eij, and Emax are the same for GHD and LHD.
The formula for Vij(max) and Fij(max) are also the same except for a
pre-factor Cij, explained below.
Emax = is the max of the absolute value of Eij for all IJ bonds. The bias energy Vij applied to a bond IJ with maximum strain is
dVij = Vmax * (1 - (Eij/q)^2) for abs(Eij) < q Vij(max) = Cij * Vmax * (1 - (Eij/q)^2) for abs(Eij) < qfactor
= 0 otherwise = 0 otherwise :pre
Delta Vbias = minimum of dVij for all IJ bonds The derivative of Vij(max) with respect to the position of each atom
in the IJ bond gives a force Fij(max) acting on the bond as
The boost factor B = exp(beta * Delta Vbias) Fij(max) = - dVij(max)/dEij = 2 Cij Vmax Eij / qfactor^2 for abs(Eij) < qfactor
for a single timestep. = 0 otherwise :pre
Thus the accumulated hypertime is simply which can be decomposed into an equal and opposite force acting on
only the two I,J atoms in the IJ bond.
The key difference is that in GHD this bias energy and force is added
(on a particular timestep) to only one bond (pair of atoms) in the
system, the bond which has maximum strain Emax.
In LHD, the bias energy and force
t_hyper = Sum (i = 1 to Nsteps) Bi * dt
NOTE: Add eqs for boostostat and boost coeff. NOTE: Add eqs for boostostat and boost coeff.
Explain how many bonds are boosted simultaneously Explain how many bonds are boosted simultaneously
and how to choose boost factor and initial Vmax. and how to choose boost factor and initial Vmax.
:line
Here is additional information on the input parameters for GHD.
The {cutbond} argument is the cutoff distance for defining bonds The {cutbond} argument is the cutoff distance for defining bonds
between pairs of nearby atoms. A pair of atoms in their equilibrium, between pairs of nearby atoms. A pair of I,J atoms in their
minimum-energy config, which are a distance Rij < cutbond, are equilibrium, minimum-energy configuration, which are separated by a
defined as a bonded pair. distance Rij < {cutbond}, are flagged as a bonded pair. Setting
{cubond} to be ~25% larger than the nearest-neighbor distance in a
crystalline lattice is a typical choice for solids, so that bonds
exist only between nearest neighbor pairs.
The {qfactor} argument is the limiting strain at which The {qfactor} argument is the limiting strain at which the bias
the Vbias (the bias potential) goes to 0.0. potential goes to 0.0. It is dimensionless, so a value of 0.3 means a
bond distance can be up to 30% larger or 30% smaller than the
equilibrium (quenched) R0ij distance and the two atoms in the bond
could still experience a non-zero bias force.
If qfactor is too big, then transitions are affected b/c If {qfactor} is set too large, then transitions from one energy basin
the bias energy is non-zero at the transitions. If it is to another are affected because the bias potneial is non-zero at the
too small than not must boost is achieved for a large system transition state (e.g. saddle point). If {qfactor} is set too small
with many bonds (some bonds Eij always exceeds qfactor). than little boost can be achieved because the Eij strain of some bond in
the system will (nearly) always exceed {qfactor). A value of 0.3 for
{qfactor} is typically a reasonable value.
The {Vmax} argument is the initial prefactor on the bias potential. What about requested boost
Should be chosen as estimate of final. Will be adjusted by
boost cooeficient
The {Tequil} argument is part of the beta term in the
exponential factor that determines how much boost is The {Vmax} argument is the prefactor on the bias potential. It should
be set to a value less than the smallest barrier height for an event
to occur. Otherwise the applied bias potential may be large enough
(when added to the interatomic potential) to produce a local energy
basin with a maxima in the center. This can produce artificial energy
minima in the same basin that trap an atom. Or if {Vmax} is even
larger, it may induce an atom(s) to rapidly transition to another
energy basin. Both cases are "bad dynamics" which violate the
assumptions of GHD that guarantee an accelerated time-accurate
trajectory of the system.
The {Tequil} argument is the temperature at which the system is
simulated; see the comment above about the "fix
langevin"_fix_langevin.html thermostatting. It is also part of the
beta term in the exponential factor that determines how much boost is
achieved as a function of the bias potential. achieved as a function of the bias potential.
In general, the lower the value of {Tequil} and the higher the value
of {Vmax}, the more boost will be achieved by the GHD algorithm.
The {Dcut} argument is the distance required between two bonds for The {Dcut} argument is the distance required between two bonds for
them to be selected as both being boosted. The distance is between them to be selected as both being boosted. The distance is between
the center points of each bond. Actually between any pair of atoms in the center points of each bond. Actually between any pair of atoms in
@ -122,6 +171,8 @@ bonded atoms may become unphysically large, leading to bad dynamics.
If chosen too small, the hyperdynamics run may be inefficient in the If chosen too small, the hyperdynamics run may be inefficient in the
sense that events take a long time to occur. sense that events take a long time to occur.
:line
[Restart, fix_modify, output, run start/stop, minimize info:] [Restart, fix_modify, output, run start/stop, minimize info:]
No information about this fix is written to "binary restart No information about this fix is written to "binary restart
@ -165,11 +216,10 @@ is run (distance units)
22 = cummulative count of atoms in events since fix created 22 = cummulative count of atoms in events since fix created
23 = cummulative # of new bonds since fix created :ul 23 = cummulative # of new bonds since fix created :ul
The first quantities (1-5) are for the current timestep. The The first quantities (1-5) are for the current timestep. Quantities
quantities 6-19 are for the current hyper run. They are reset each 6-19 are for the current hyper run. They are reset each time a new
time a new hyper run is performed. The quantities 20-23 are hyper run is performed. Quantities 20-23 are cummulative across
cummulative across multiple hyper runs, They are only set to initial multiple runs (since the fix was defined in the input script).
values once, when this fix is defined in the input script.
For value 6, the numerator is a count of all biased bonds on every For value 6, the numerator is a count of all biased bonds on every
timestep whose bias value = 0.0. The denominator is the count of all timestep whose bias value = 0.0. The denominator is the count of all
@ -198,3 +248,12 @@ LAMMPS"_Section_start.html#start_3 section for more info.
[Default:] None [Default:] None
:line
:link(Voter2013lhd)
[(Voter2013)] S. Y. Kim, D. Perez, A. F. Voter, J Chem Phys, 139,
144110 (2013).
:link(Mironlhd)
[(Miron)] R. A. Miron and K. A. Fichthorn, J Chem Phys, 119, 6210 (2003).