openmmtools.mcmc.LangevinSplittingDynamicsMove

class openmmtools.mcmc.LangevinSplittingDynamicsMove(timestep=Quantity(value=1.0, unit=femtosecond), collision_rate=Quantity(value=10.0, unit=/picosecond), n_steps=1000, reassign_velocities=False, splitting='V R O R V', constraint_tolerance=1e-08, measure_shadow_work=False, measure_heat=False, **kwargs)[source]

Langevin dynamics segment with custom splitting of the operators and optional Metropolized Monte Carlo validation.

Besides all the normal properties of the LangevinDynamicsMove, this class implements the custom splitting sequence of the openmmtools.integrators.LangevinIntegrator. Additionally, the steps can be wrapped around a proper Generalized Hybrid Monte Carlo step to ensure that the exact distribution is generated.

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).

splittingstring, default: “V R O R V”

Sequence of “R”, “V”, “O” (and optionally “{”, “}”, “V0”, “V1”, …) substeps to be executed each timestep.

Forces are only used in V-step. Handle multiple force groups by appending the force group index to V-steps, e.g. “V0” will only use forces from force group 0. “V” will perform a step using all forces. “{” will cause metropolization, and must be followed later by a “}”.

constraint_tolerancefloat, default: 1.0e-8

Tolerance for constraint solver

measure_shadow_workboolean, default: False

Accumulate the shadow work performed by the symplectic substeps, in the global shadow_work

measure_heatboolean, default: False

Accumulate the heat exchanged with the bath in each step, in the global heat

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 = LangevinSplittingDynamicsMove()

or create a Langevin move with specified splitting.

>>> move = LangevinSplittingDynamicsMove(splitting="O { V R V } O")

Where this splitting is a 5 step symplectic integrator:

*. Ornstein-Uhlenbeck (O) interactions with the stochastic heat bath interactions *. Hybrid Metropolized step around the half-step velocity updates (V) with full position updates (R).

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.

splittingstr

Splitting applied to this integrator represented as a string.

constraint_tolerancefloat, default: 1.0e-8

Tolerance for constraint solver

measure_shadow_workboolean, default: False

Accumulate the shadow work performed by the symplectic substeps, in the global shadow_work

measure_heatboolean, default: False

Accumulate the heat exchanged with the bath in each step, in the global heat

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, splitting='V R O R V', constraint_tolerance=1e-08, measure_shadow_work=False, measure_heat=False, **kwargs)[source]

Methods

__init__([timestep, collision_rate, ...])

apply(thermodynamic_state, sampler_state[, ...])

Apply the Langevin dynamics MCMC move.

Attributes

context_cache