forked from lijiext/lammps
301 lines
14 KiB
Plaintext
301 lines
14 KiB
Plaintext
"LAMMPS WWW Site"_lws - "LAMMPS Documentation"_ld - "LAMMPS Commands"_lc :c
|
||
|
||
:link(lws,http://lammps.sandia.gov)
|
||
:link(ld,Manual.html)
|
||
:link(lc,Section_commands.html#comm)
|
||
|
||
:line
|
||
|
||
tad command :h3
|
||
|
||
[Syntax:]
|
||
|
||
tad N t_event T_lo T_hi delta tmax compute-ID \
|
||
seed keyword value ... :pre
|
||
|
||
N = # of timesteps to run (not including dephasing/quenching) :ulb,l
|
||
t_event = timestep interval between event checks :l
|
||
T_lo = temperature at which event times are desired :l
|
||
T_hi = temperature at which MD simulation is performed :l
|
||
delta = desired confidence level for stopping criterion :l
|
||
tmax = reciprocal of lowest expected preexponential factor (time units) :l
|
||
compute-ID = ID of the compute used for event detection :l
|
||
zero or more keyword/value pairs may be appended :l
|
||
keyword = {min} or {neb} or {min_style} or {neb_style} or {neb_log} :l
|
||
{min} values = etol ftol maxiter maxeval
|
||
etol = stopping tolerance for energy (energy units)
|
||
ftol = stopping tolerance for force (force units)
|
||
maxiter = max iterations of minimize
|
||
maxeval = max number of force/energy evaluations
|
||
{neb} values = ftol N1 N2 Nevery
|
||
etol = stopping tolerance for energy (energy units)
|
||
ftol = stopping tolerance for force (force units)
|
||
N1 = max # of iterations (timesteps) to run initial NEB
|
||
N2 = max # of iterations (timesteps) to run barrier-climbing NEB
|
||
Nevery = print NEB statistics every this many timesteps
|
||
{min_style} value = {cg} or {hftn} or {sd} or {quickmin} or {fire}
|
||
{neb_style} value = {quickmin} or {fire}
|
||
{neb_log} value = file where NEB statistics are printed :pre
|
||
|
||
:ule
|
||
|
||
[Examples:]
|
||
|
||
tad 2000 50 1800 2300 0.01 0.01 event 54985
|
||
tad 2000 50 1800 2300 0.01 0.01 event 54985 &
|
||
min 1e-05 1e-05 100 100 &
|
||
neb 0.0 0.01 200 200 20 &
|
||
min_style cg &
|
||
neb_style fire &
|
||
neb_log log.neb :pre
|
||
|
||
[Description:]
|
||
|
||
Run a temperature accelerated dynamics (TAD) simulation. This method
|
||
requires two or more partitions to perform NEB transition state
|
||
searches.
|
||
|
||
TAD is described in "this paper"_#Voter by Art Voter. It is a method
|
||
that uses accelerated dynamics at an elevated temperature to generate
|
||
results at a specified lower temperature. A good overview of
|
||
accelerated dynamics methods for such systems is given in "this review
|
||
paper"_#Voter2 from the same group. In general, these methods assume
|
||
that the long-time dynamics is dominated by infrequent events i.e. the
|
||
system is is confined to low energy basins for long periods,
|
||
punctuated by brief, randomly-occurring transitions to adjacent
|
||
basins. TAD is suitable for infrequent-event systems, where in
|
||
addition, the transition kinetics are well-approximated by harmonic
|
||
transition state theory (hTST). In hTST, the temperature dependence of
|
||
transition rates follows the Arrhenius relation. As a consequence a
|
||
set of event times generated in a high-temperature simulation can be
|
||
mapped to a set of much longer estimated times in the low-temperature
|
||
system. However, because this mapping involves the energy barrier of
|
||
the transition event, which is different for each event, the first
|
||
event at the high temperature may not be the earliest event at the low
|
||
temperature. TAD handles this by first generating a set of possible
|
||
events from the current basin. After each event, the simulation is
|
||
reflected backwards into the current basin. This is repeated until
|
||
the stopping criterion is satisfied, at which point the event with the
|
||
earliest low-temperature occurrence time is selected. The stopping
|
||
criterion is that the confidence measure be greater than
|
||
1-{delta}. The confidence measure is the probability that no earlier
|
||
low-temperature event will occur at some later time in the
|
||
high-temperature simulation. hTST provides an lower bound for this
|
||
probability, based on the user-specified minimum pre-exponential
|
||
factor (reciprocal of {tmax}).
|
||
|
||
In order to estimate the energy barrier for each event, the TAD method
|
||
invokes the "NEB"_neb.html method. Each NEB replica runs on a
|
||
partition of processors. The current NEB implementation in LAMMPS
|
||
restricts you to having exactly one processor per replica. For more
|
||
information, see the documentation for the "neb"_neb.html command. In
|
||
the current LAMMPS implementation of TAD, all the non-NEB TAD
|
||
operations are performed on the first partition, while the other
|
||
partitions remain idle. See "this section"_Section_howto.html#howto_5
|
||
of the manual for further discussion of multi-replica simulations.
|
||
|
||
|
||
A TAD run has several stages, which are repeated each time an event is
|
||
performed. The logic for a TAD run is as follows:
|
||
|
||
while (time remains):
|
||
while (time < tstop):
|
||
until (event occurs):
|
||
run dynamics for t_event steps
|
||
quench
|
||
run neb calculation using all replicas
|
||
compute tlo from energy barrier
|
||
update earliest event
|
||
update tstop
|
||
reflect back into current basin
|
||
execute earliest event :pre
|
||
|
||
Before this outer loop begins, the initial potential energy basin is
|
||
identified by quenching (an energy minimization, see below) the
|
||
initial state and storing the resulting coordinates for reference.
|
||
|
||
Inside the inner loop, dynamics is run continuously according to
|
||
whatever integrator has been specified by the user, stopping every
|
||
{t_event} steps to check if a transition event has occurred. This
|
||
check is performed by quenching the system and comparing the resulting
|
||
atom coordinates to the coordinates from the previous basin.
|
||
|
||
A quench is an energy minimization and is performed by whichever
|
||
algorithm has been defined by the {min} and {min_style} keywords or
|
||
their default values. Note that typically, you do not need to perform
|
||
a highly-converged minimization to detect a transition event.
|
||
|
||
The event check is performed by a compute with the specified
|
||
{compute-ID}. Currently there is only one compute that works with the
|
||
TAD commmand, which is the "compute
|
||
event/displace"_compute_event_displace.html command. Other
|
||
event-checking computes may be added. "Compute
|
||
event/displace"_compute_event_displace.html checks whether any atom in
|
||
the compute group has moved further than a specified threshold
|
||
distance. If so, an "event" has occurred.
|
||
|
||
The neb calculation is similar to that invoked by the "neb"_neb.html
|
||
command, except that the final state is generated internally, instead
|
||
of being read in from a file. The TAD implementation provides default
|
||
values for the NEB settings, which can be overridden using the {neb}
|
||
and {neb_style} keywords.
|
||
|
||
:line
|
||
|
||
A key aspect of the TAD method is setting the stopping criterion
|
||
appropriately. If this criterion is too conservative, then many
|
||
events must be generated before one is finally executed. Conversely,
|
||
if this criterion is too aggressive, high-entropy high-barrier events
|
||
will be over-sampled, while low-entropy low-barrier events will be
|
||
under-sampled. If the lowest pre-exponential factor is known fairly
|
||
accurately, then it can be used to estimate {tmax}, and the value of
|
||
{delta} can be set to the desired confidence level e.g. {delta} = 0.05
|
||
corresponds to 95% confidence. However, for systems where the dynamics
|
||
are not well characterized (the most common case), it will be
|
||
necessary to experiment with the values of {delta} and {tmax} to get a
|
||
good trade-off between accuracy and performance.
|
||
|
||
A second key aspect is the choice of {t_hi}. A larger value greatly
|
||
increases the rate at which new events are generated. However, too
|
||
large a value introduces errors due to anharmonicity (not accounted
|
||
for within hTST). Once again, for any given system, experimentation is
|
||
necessary to determine the best value of {t_hi}.
|
||
|
||
:line
|
||
|
||
Five kinds of output can be generated during a TAD run: event
|
||
statistics, NEB statistics, thermodynamic output by each replica, dump
|
||
files, and restart files.
|
||
|
||
Event statistics are printed to the screen and master log.lammps file
|
||
each time an event is executed. The quantities are the timestep, CPU
|
||
time, global event number {N}, local event number {M}, event status,
|
||
energy barrier, time margin, {t_lo} and {delt_lo}. The timestep is
|
||
the usual LAMMPS timestep, which corresponds to the high-temperature
|
||
time at which the event was detected, in units of timestep. The CPU
|
||
time is the total processor time since the start of the TAD run. The
|
||
global event number {N} is a counter that increments with each
|
||
executed event. The local event number {M} is a counter that resets to
|
||
zero upon entering each new basin. The event status is {E} when an
|
||
event is executed, and is {D} for an event that is detected, while
|
||
{DF} is for a detected event that is also the earliest (first) event
|
||
at the low temperature.
|
||
|
||
The time margin is the ratio of the high temperature time in the
|
||
current basin to the stopping time. This last number can be used to
|
||
judge whether the stopping time is too short or too long (see above).
|
||
|
||
{t_lo} is the low-temperature event time when the current basin was
|
||
entered, in units of timestep. del{t_lo} is the time of each detected
|
||
event, measured relative to {t_lo}. {delt_lo} is equal to the
|
||
high-temperature time since entering the current basin, scaled by an
|
||
exponential factor that depends on the hi/lo temperature ratio and the
|
||
energy barrier for that event.
|
||
|
||
On lines for executed events, with status {E}, the global event number
|
||
is incremented by one, and the timestep, local event number, energy
|
||
barrier, {t_lo}, and {delt_lo} match the last event with status {DF}
|
||
in the immediately preceding block of detected events.
|
||
|
||
The NEB statistics are written to the file specified by the {neb_log}
|
||
keyword. If the keyword value is "none", then no NEB statistics are
|
||
printed out. The statistics are written every {Nevery} timesteps. See
|
||
the "neb"_neb.html command for a full description of the NEB
|
||
statistics. When invoked from TAD, NEB statistics are never printed to
|
||
the screen.
|
||
|
||
Because the NEB calculation must run on multiple partitions, LAMMPS
|
||
produces additional screen and log files for each partition,
|
||
e.g. log.lammps.0, log.lammps.1, etc. For the TAD command, these
|
||
contain the thermodynamic output of each NEB replica. In addition, the
|
||
log file for the first partition, log.lammps.0, will contain
|
||
thermodynamic output from short runs and minimizations corresponding
|
||
to the dynamics and quench operations, as well as a line for each new
|
||
detected event, as described above.
|
||
|
||
After the TAD command completes, timing statistics for the TAD run are
|
||
printed in each replica's log file, giving a breakdown of how much CPU
|
||
time was spent in each stage (NEB, dynamics, quenching, etc).
|
||
|
||
Any "dump files"_dump.html defined in the input script will be written
|
||
to during a TAD run at timesteps when an event is executed. This
|
||
means the the requested dump frequency in the "dump"_dump.html command
|
||
is ignored. There will be one dump file (per dump command) created
|
||
for all partitions. The atom coordinates of the dump snapshot are
|
||
those of the minimum energy configuration resulting from quenching
|
||
following the executed event. The timesteps written into the dump
|
||
files correspond to the timestep at which the event occurred and NOT
|
||
the clock. A dump snapshot corresponding to the initial minimum state
|
||
used for event detection is written to the dump file at the beginning
|
||
of each TAD run.
|
||
|
||
If the "restart"_restart.html command is used, a single restart file
|
||
for all the partitions is generated, which allows a TAD run to be
|
||
continued by a new input script in the usual manner. The restart file
|
||
is generated after an event is executed. The restart file contains a
|
||
snapshot of the system in the new quenched state, including the event
|
||
number and the low-temperature time. The restart frequency specified
|
||
in the "restart"_restart.html command is interpreted differently when
|
||
performing a TAD run. It does not mean the timestep interval between
|
||
restart files. Instead it means an event interval for executed
|
||
events. Thus a frequency of 1 means write a restart file every time
|
||
an event is executed. A frequency of 10 means write a restart file
|
||
every 10th executed event. When an input script reads a restart file
|
||
from a previous TAD run, the new script can be run on a different
|
||
number of replicas or processors.
|
||
|
||
Note that within a single state, the dynamics will typically
|
||
temporarily continue beyond the event that is ultimately chosen, until
|
||
the stopping criterionis satisfied. When the event is eventually
|
||
executed, the timestep counter is reset to the value when the event
|
||
was detected. Similarly, after each quench and NEB minimization, the
|
||
timestep counter is reset to the value at the start of the
|
||
minimization. This means that the timesteps listed in the replica log
|
||
files do not always increase monotonically. However, the timestep
|
||
values printed to the master log file, dump files, and restart files
|
||
are always monotonically increasing.
|
||
|
||
:line
|
||
|
||
[Restrictions:]
|
||
|
||
This command can only be used if LAMMPS was built with the REPLICA
|
||
package. See the "Making LAMMPS"_Section_start.html#start_3 section
|
||
for more info on packages.
|
||
|
||
{N} setting must be integer multiple of {t_event}.
|
||
|
||
Runs restarted from restart files written during a TAD run will only
|
||
produce identical results if the user-specified integrator supports
|
||
exact restarts. So "fix nvt"_fix_nh.html will produce an exact
|
||
restart, but "fix langevin"_fix_langevin.html will not.
|
||
|
||
This command cannot be used when any fixes are defined that keep track
|
||
of elapsed time to perform time-dependent operations. Examples
|
||
include the "ave" fixes such as "fix
|
||
ave/spatial"_fix_ave_spatial.html. Also "fix
|
||
dt/reset"_fix_dt_reset.html and "fix deposit"_fix_deposit.html.
|
||
|
||
[Related commands:]
|
||
|
||
"compute event/displace"_compute_event_displace.html,
|
||
"min_modify"_min_modify.html, "min_style"_min_style.html,
|
||
"run_style"_run_style.html, "minimize"_minimize.html,
|
||
"temper"_temper.html, "neb"_neb.html,
|
||
"prd"_prd.html
|
||
|
||
[Default:]
|
||
|
||
The option defaults are {min} = 0.1 0.1 40 50, {neb} = 0.01 100 100
|
||
10, {min_style} = {cg}, {neb_style} = {quickmin}, and {neb_log} =
|
||
"none"
|
||
|
||
:line
|
||
|
||
:link(Voter)
|
||
[(Voter)] S<>rensen and Voter, J Chem Phys, 112, 9599 (2000)
|
||
|
||
:link(Voter2)
|
||
[(Voter2)] Voter, Montalenti, Germann, Annual Review of Materials
|
||
Research 32, 321 (2002).
|