openmmtools.mcmc.LangevinDynamicsMove¶
- class openmmtools.mcmc.LangevinDynamicsMove(timestep=Quantity(value=1.0, unit=femtosecond), collision_rate=Quantity(value=10.0, unit=/picosecond), n_steps=1000, reassign_velocities=False, constraint_tolerance=1e-08, **kwargs)[source]¶
Langevin dynamics segment as a (pseudo) Monte Carlo move.
This move assigns a velocity from the Maxwell-Boltzmann distribution and executes a number of Maxwell-Boltzmann steps to propagate dynamics. This is not a true Monte Carlo move, in that the generation of the correct distribution is only exact in the limit of infinitely small timestep; in other words, the discretization error is assumed to be negligible. Use HybridMonteCarloMove instead to ensure the exact distribution is generated.
The OpenMM LangevinMiddleIntegrator, based on BAOAB [1], is used.
Warning
The LangevinMiddleIntegrator generates velocities that are half a timestep lagged behind the positions.
Warning
No Metropolization is used to ensure the correct phase space distribution is sampled. This means that timestep-dependent errors will remain uncorrected, and are amplified with larger timesteps. Use this move at your own risk!
- Parameters:
- timestepopenmm.unit.Quantity, optional
The timestep to use for Langevin integration (time units, default is 1*openmm.unit.femtosecond).
- collision_rateopenmm.unit.Quantity, optional
The collision rate with fictitious bath particles (1/time units, default is 10/openmm.unit.picoseconds).
- n_stepsint, optional
The number of integration timesteps to take each time the move is applied (default is 1000).
- reassign_velocitiesbool, optional
If True, the velocities will be reassigned from the Maxwell-Boltzmann distribution at the beginning of the move (default is False).
- constraint_tolerancefloat, optional
Fraction of the constrained distance within which constraints are maintained for the integrator (default is 1e-8).
References
[1] Leimkuhler B and Matthews C. Robust and efficient configurational molecular sampling via Langevin dynamics. https://doi.org/10.1063/1.4802990 [2] Leimkuhler B and Matthews C. Efficient molecular dynamics using geodesic integration and solvent–solute splitting. https://doi.org/10.1098/rspa.2016.0138
Examples
First we need to create the thermodynamic state and the sampler state to propagate. Here we create an alanine dipeptide system in vacuum.
>>> from openmm import unit >>> from openmmtools import testsystems >>> from openmmtools.states import SamplerState, ThermodynamicState >>> test = testsystems.AlanineDipeptideVacuum() >>> sampler_state = SamplerState(positions=test.positions) >>> thermodynamic_state = ThermodynamicState(system=test.system, temperature=298*unit.kelvin)
Create a Langevin move with default parameters
>>> move = LangevinDynamicsMove()
or create a Langevin move with specified parameters.
>>> move = LangevinDynamicsMove(timestep=0.5*unit.femtoseconds, ... collision_rate=20.0/unit.picoseconds, n_steps=10)
Perform one update of the sampler state. The sampler state is updated with the new state.
>>> move.apply(thermodynamic_state,sampler_state,context_cache=context_cache) >>> np.allclose(sampler_state.positions, test.positions) False
The same move can be applied to a different state, here an ideal gas.
>>> test = testsystems.IdealGas() >>> sampler_state = SamplerState(positions=test.positions) >>> thermodynamic_state = ThermodynamicState(system=test.system, ... temperature=298*unit.kelvin) >>> move.apply(thermodynamic_state,sampler_state,context_cache=context_cache) >>> np.allclose(sampler_state.positions, test.positions) False
- Attributes:
- timestepopenmm.unit.Quantity
The timestep to use for Langevin integration (time units).
- collision_rateopenmm.unit.Quantity
The collision rate with fictitious bath particles (1/time units).
- n_stepsint
The number of integration timesteps to take each time the move is applied.
- reassign_velocitiesbool
If True, the velocities will be reassigned from the Maxwell-Boltzmann distribution at the beginning of the move.
- constraint_tolerancefloat
Fraction of the constrained distance within which constraints are maintained for the integrator.
Methods
apply
(thermodynamic_state, sampler_state[, ...])Apply the Langevin dynamics MCMC move.
- __init__(timestep=Quantity(value=1.0, unit=femtosecond), collision_rate=Quantity(value=10.0, unit=/picosecond), n_steps=1000, reassign_velocities=False, constraint_tolerance=1e-08, **kwargs)[source]¶
Methods
__init__
([timestep, collision_rate, ...])apply
(thermodynamic_state, sampler_state[, ...])Apply the Langevin dynamics MCMC move.
Attributes
context_cache