openmmtools.integrators.AlchemicalNonequilibriumLangevinIntegrator¶

class openmmtools.integrators.AlchemicalNonequilibriumLangevinIntegrator(alchemical_functions=None, splitting='O { V R H R V } O', nsteps_neq=100, *args, **kwargs)[source]¶

Allows nonequilibrium switching based on force parameters specified in alchemical_functions. A variable named lambda is switched from 0 to 1 linearly throughout the nsteps of the protocol. The functions can use this to create more complex protocols for other global parameters.

Propagator is based on Langevin splitting, as described below.

One way to divide the Langevin system is into three parts which can each be solved “exactly:”

R: Linear “drift” / Constrained “drift”
Deterministic update of positions, using current velocities x <- x + v dt
V: Linear “kick” / Constrained “kick”
Deterministic update of velocities, using current forces v <- v + (f/m) dt

where f = force, m = mass
O: Ornstein-Uhlenbeck
Stochastic update of velocities, simulating interaction with a heat bath v <- av + b sqrt(kT/m) R

where a = e^(-gamma dt) b = sqrt(1 - e^(-2gamma dt)) R is i.i.d. standard normal
H: Hamiltonian update step

We can then construct integrators by solving each part for a certain timestep in sequence. (We can further split up the V step by force group, evaluating cheap but fast-fluctuating forces more frequently than expensive but slow-fluctuating forces. Since forces are only evaluated in the V step, we represent this by including in our “alphabet” V0, V1, …)

When the system contains holonomic constraints, these steps are confined to the constraint manifold.

References

[Nilmeier, et al. 2011] Nonequilibrium candidate Monte Carlo is an efficient tool for equilibrium simulation [Leimkuhler and Matthews, 2015] Molecular dynamics: with deterministic and stochastic numerical methods, Chapter 7

Examples

Create a nonequilibrium integrator to switch the center of a harmonic oscillator

>>> # Create harmonic oscillator testsystem
>>> from openmmtools import testsystems
>>> import openmm
>>> from openmm import unit
>>> testsystem = testsystems.HarmonicOscillator()
>>> # Create a nonequilibrium alchemical integrator
>>> alchemical_functions = { 'testsystems_HarmonicOscillator_x0' : 'lambda' }
>>> nsteps_neq = 100 # number of steps in the switching trajectory where lambda is switched from 0 to 1
>>> integrator = AlchemicalNonequilibriumLangevinIntegrator(temperature=300*unit.kelvin, collision_rate=1.0/unit.picoseconds, timestep=1.0*unit.femtoseconds,
...                                                         alchemical_functions=alchemical_functions, splitting="O { V R H R V } O", nsteps_neq=nsteps_neq,
...                                                         measure_shadow_work=True)
>>> # Create a Context
>>> context = openmm.Context(testsystem.system, integrator)
>>> # Run the whole switching trajectory
>>> context.setPositions(testsystem.positions)
>>> integrator.step(nsteps_neq)
>>> protocol_work = integrator.protocol_work # retrieve protocol work (excludes shadow work)
>>> total_work = integrator.total_work # retrieve total work (includes shadow worl)
>>> # Reset and run again
>>> context.setPositions(testsystem.positions)
>>> integrator.reset()
>>> integrator.step(nsteps_neq)

Attributes:

_kinetic_energystr: This is 0.5*m*v*v by default, and is the expression used for the kinetic energy

Methods

`addComputeGlobal`(self, variable, expression)	Add a step to the integration algorithm that computes a global value.
`addComputePerDof`(self, variable, expression)	Add a step to the integration algorithm that computes a per-DOF value.
`addComputeSum`(self, variable, expression)	Add a step to the integration algorithm that computes a sum over degrees of freedom.
`addComputeTemperatureDependentConstants`(...)	Wrap the ComputePerDof into an if-block executed only when kT changes.
`addConstrainPositions`(self)	Add a step to the integration algorithm that updates particle positions so all constraints are satisfied.
`addConstrainVelocities`(self)	Add a step to the integration algorithm that updates particle velocities so the net velocity along all constraints is 0.
`addGlobalVariable`(self, name, initialValue)	Define a new global variable.
`addPerDofVariable`(self, name, initialValue)	Define a new per-DOF variable.
`addTabulatedFunction`(self, name, function)	Add a tabulated function that may appear in expressions.
`addUpdateContextState`(self)	Add a step to the integration algorithm that allows Forces to update the context state.
`beginIfBlock`(self, condition)	Add a step which begins a new "if" block.
`beginWhileBlock`(self, condition)	Add a step which begins a new "while" block.
`deserialize_xml`(xml_serialization)	Shortcut to deserialize the XML representation and the restore interface.
`endBlock`(self)	Add a step which marks the end of the most recently begun "if" or "while" block.
`getComputationStep`(self, index)	Get the details of a computation step that has been added to the integration algorithm.
`getConstraintTolerance`(self)	Get the distance tolerance within which constraints are maintained, as a fraction of the constrained distance.
`getGlobalVariable`(self, index)	Get the current value of a global variable.
`getGlobalVariableByName`(self, name)	Get the current value of a global variable, specified by name.
`getGlobalVariableName`(self, index)	Get the name of a global variable.
`getIntegrationForceGroups`(self)	Get which force groups to use for integration.
`getKineticEnergyExpression`(self)	Get the expression to use for computing the kinetic energy.
`getNumComputations`(self)	Get the number of computation steps that have been added.
`getNumGlobalVariables`(self)	Get the number of global variables that have been defined.
`getNumPerDofVariables`(self)	Get the number of per-DOF variables that have been defined.
`getNumTabulatedFunctions`(self)	Get the number of tabulated functions that have been defined.
`getPerDofVariable`()
`getPerDofVariableByName`(self, name)	Get the value of a per-DOF variable, specified by name.
`getPerDofVariableName`(self, index)	Get the name of a per-DOF variable.
`getRandomNumberSeed`(self)	Get the random number seed.
`getStepSize`(self)	Get the size of each time step, in picoseconds.
`getTabulatedFunction`(-> TabulatedFunction)	Get a reference to a tabulated function that may appear in expressions.
`getTabulatedFunctionName`(self, index)	Get the name of a tabulated function that may appear in expressions.
`getTemperature`()	Return the temperature of the heat bath.
`get_acceptance_rate`()	Get acceptance rate for Metropolized integrators.
`get_heat`([dimensionless])	Get the current accumulated heat.
`get_protocol_work`([dimensionless])	Get the current accumulated protocol work.
`get_shadow_work`([dimensionless])	Get the current accumulated shadow work.
`get_total_work`([dimensionless])	Get the current accumulated total work.
`is_restorable`(openmm_object)	Check if the custom integrator or force has a restorable interface.
`is_thermostated`(integrator)	Return true if the integrator is a ThermostatedIntegrator.
`pretty_format`([as_list, step_types_to_highlight])	Generate a human-readable version of each integrator step.
`pretty_print`()	Pretty-print the computation steps of this integrator.
`reset`()	Reset all statistics, alchemical parameters, and work.
`reset_ghmc_statistics`()	Reset GHMC acceptance rate statistics.
`reset_heat`()	Reset heat.
`reset_protocol_work`()	Reset the protocol work.
`reset_shadow_work`()	Reset shadow work.
`restore_interface`(integrator)	Restore the original interface of a CustomIntegrator.
`setConstraintTolerance`(self, tol)	Set the distance tolerance within which constraints are maintained, as a fraction of the constrained distance.
`setGlobalVariable`(self, index, value)	Set the value of a global variable.
`setGlobalVariableByName`(self, name, value)	Set the value of a global variable, specified by name.
`setIntegrationForceGroups`(groups)	Set which force groups to use for integration.
`setKineticEnergyExpression`(self, expression)	Set the expression to use for computing the kinetic energy.
`setPerDofVariable`(self, index, values)	Set the value of a per-DOF variable.
`setPerDofVariableByName`(self, name, values)	Set the value of a per-DOF variable, specified by name.
`setRandomNumberSeed`(self, seed)	Set the random number seed.
`setStepSize`(self, size)	Set the size of each time step, in picoseconds.
`setTemperature`(temperature)	Set the temperature of the heat bath.
`step`(self, steps)	Advance a simulation through time by taking a series of time steps.

__init__(alchemical_functions=None, splitting='O { V R H R V } O', nsteps_neq=100, *args, **kwargs)[source]¶

Parameters:

alchemical_functionsdict of strings, optional, default=None

key: value pairs such as “global_parameter” : function_of_lambda where function_of_lambda is a Lepton-compatible string that depends on the variable “lambda” If not specified, no alchemical functions will be used.

splittingstring, default: “O { V R H R V } O”

Sequence of “R”, “V”, “O” (and optionally “{”, “}”, “V0”, “V1”, …) substeps to be executed each timestep. “H” increments the global parameter lambda by 1/nsteps_neq for each step and accumulates protocol work.

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

temperaturenp.unit.Quantity compatible with kelvin, default: 298.0*unit.kelvin

Fictitious “bath” temperature

collision_ratenp.unit.Quantity compatible with 1/picoseconds, default: 91.0/unit.picoseconds

Collision rate

timestepnp.unit.Quantity compatible with femtoseconds, default: 1.0*unit.femtoseconds

Integration timestep

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

nsteps_neqint, default: 100

Number of steps in nonequilibrium protocol. Default 100 This number cannot be changed without creating a new integrator.

Methods

`__init__`([alchemical_functions, splitting, ...])	Parameters:
`addComputeGlobal`(self, variable, expression)	Add a step to the integration algorithm that computes a global value.
`addComputePerDof`(self, variable, expression)	Add a step to the integration algorithm that computes a per-DOF value.
`addComputeSum`(self, variable, expression)	Add a step to the integration algorithm that computes a sum over degrees of freedom.
`addComputeTemperatureDependentConstants`(...)	Wrap the ComputePerDof into an if-block executed only when kT changes.
`addConstrainPositions`(self)	Add a step to the integration algorithm that updates particle positions so all constraints are satisfied.
`addConstrainVelocities`(self)	Add a step to the integration algorithm that updates particle velocities so the net velocity along all constraints is 0.
`addGlobalVariable`(self, name, initialValue)	Define a new global variable.
`addPerDofVariable`(self, name, initialValue)	Define a new per-DOF variable.
`addTabulatedFunction`(self, name, function)	Add a tabulated function that may appear in expressions.
`addUpdateContextState`(self)	Add a step to the integration algorithm that allows Forces to update the context state.
`beginIfBlock`(self, condition)	Add a step which begins a new "if" block.
`beginWhileBlock`(self, condition)	Add a step which begins a new "while" block.
`deserialize_xml`(xml_serialization)	Shortcut to deserialize the XML representation and the restore interface.
`endBlock`(self)	Add a step which marks the end of the most recently begun "if" or "while" block.
`getComputationStep`(self, index)	Get the details of a computation step that has been added to the integration algorithm.
`getConstraintTolerance`(self)	Get the distance tolerance within which constraints are maintained, as a fraction of the constrained distance.
`getGlobalVariable`(self, index)	Get the current value of a global variable.
`getGlobalVariableByName`(self, name)	Get the current value of a global variable, specified by name.
`getGlobalVariableName`(self, index)	Get the name of a global variable.
`getIntegrationForceGroups`(self)	Get which force groups to use for integration.
`getKineticEnergyExpression`(self)	Get the expression to use for computing the kinetic energy.
`getNumComputations`(self)	Get the number of computation steps that have been added.
`getNumGlobalVariables`(self)	Get the number of global variables that have been defined.
`getNumPerDofVariables`(self)	Get the number of per-DOF variables that have been defined.
`getNumTabulatedFunctions`(self)	Get the number of tabulated functions that have been defined.
`getPerDofVariable`()
`getPerDofVariableByName`(self, name)	Get the value of a per-DOF variable, specified by name.
`getPerDofVariableName`(self, index)	Get the name of a per-DOF variable.
`getRandomNumberSeed`(self)	Get the random number seed.
`getStepSize`(self)	Get the size of each time step, in picoseconds.
`getTabulatedFunction`(-> TabulatedFunction)	Get a reference to a tabulated function that may appear in expressions.
`getTabulatedFunctionName`(self, index)	Get the name of a tabulated function that may appear in expressions.
`getTemperature`()	Return the temperature of the heat bath.
`get_acceptance_rate`()	Get acceptance rate for Metropolized integrators.
`get_heat`([dimensionless])	Get the current accumulated heat.
`get_protocol_work`([dimensionless])	Get the current accumulated protocol work.
`get_shadow_work`([dimensionless])	Get the current accumulated shadow work.
`get_total_work`([dimensionless])	Get the current accumulated total work.
`is_restorable`(openmm_object)	Check if the custom integrator or force has a restorable interface.
`is_thermostated`(integrator)	Return true if the integrator is a ThermostatedIntegrator.
`pretty_format`([as_list, step_types_to_highlight])	Generate a human-readable version of each integrator step.
`pretty_print`()	Pretty-print the computation steps of this integrator.
`reset`()	Reset all statistics, alchemical parameters, and work.
`reset_ghmc_statistics`()	Reset GHMC acceptance rate statistics.
`reset_heat`()	Reset heat.
`reset_protocol_work`()	Reset the protocol work.
`reset_shadow_work`()	Reset shadow work.
`restore_interface`(integrator)	Restore the original interface of a CustomIntegrator.
`setConstraintTolerance`(self, tol)	Set the distance tolerance within which constraints are maintained, as a fraction of the constrained distance.
`setGlobalVariable`(self, index, value)	Set the value of a global variable.
`setGlobalVariableByName`(self, name, value)	Set the value of a global variable, specified by name.
`setIntegrationForceGroups`(groups)	Set which force groups to use for integration.
`setKineticEnergyExpression`(self, expression)	Set the expression to use for computing the kinetic energy.
`setPerDofVariable`(self, index, values)	Set the value of a per-DOF variable.
`setPerDofVariableByName`(self, name, values)	Set the value of a per-DOF variable, specified by name.
`setRandomNumberSeed`(self, seed)	Set the random number seed.
`setStepSize`(self, size)	Set the size of each time step, in picoseconds.
`setTemperature`(temperature)	Set the temperature of the heat bath.
`step`(self, steps)	Advance a simulation through time by taking a series of time steps.

Attributes

`BlockEnd`
`ComputeGlobal`
`ComputePerDof`
`ComputeSum`
`ConstrainPositions`
`ConstrainVelocities`
`IfBlockStart`
`UpdateContextState`
`WhileBlockStart`
`acceptance_rate`	Get acceptance rate for Metropolized integrators.
`global_variable_names`	The set of global variable names defined for this integrator.
`heat`
`is_metropolized`	Return True if this integrator is Metropolized, False otherwise.
`kT`	The thermal energy in openmm.Quantity
`protocol_work`	Total protocol work in energy units.
`shadow_work`
`thisown`	The membership flag
`total_work`	Total work (protocol work plus shadow work) in energy units.

openmmtools.integrators.AlchemicalNonequilibriumLangevinIntegrator¶

openmmtools

Navigation

Related Topics