openmmtools.integrators.ThermostatedIntegrator

class openmmtools.integrators.ThermostatedIntegrator(temperature, *args, **kwargs)[source]

Add temperature functions to a CustomIntegrator.

This class is intended to be inherited by integrators that maintain the stationary distribution at a given temperature. The constructor adds a global variable named “kT” defining the thermal energy at the given temperature. This global variable is updated through the temperature setter and getter.

It also provide a utility function to handle per-DOF constants that must be computed only when the temperature changes.

Notice that the CustomIntegrator internally stored by a Context object will loose setter and getter and any extra function you define. The same happens when you copy your integrator. You can restore the methods with the static method ThermostatedIntegrator.restore_interface().

Parameters:
temperatureunit.Quantity

The temperature of the integrator heat bath (temperature units).

timestepunit.Quantity

The timestep to pass to the CustomIntegrator constructor (time units).

Examples

We can inherit from ThermostatedIntegrator to automatically define setters and getters for the temperature and to add a per-DOF constant “sigma” that we need to update only when the temperature is changed.

>>> import openmm
>>> from openmm import unit
>>> class TestIntegrator(ThermostatedIntegrator):
...     def __init__(self, temperature=298.0*unit.kelvin, timestep=1.0*unit.femtoseconds):
...         super(TestIntegrator, self).__init__(temperature, timestep)
...         self.addPerDofVariable("sigma", 0)  # velocity standard deviation
...         self.addComputeTemperatureDependentConstants({"sigma": "sqrt(kT/m)"})
...

We instantiate the integrator normally.

>>> integrator = TestIntegrator(temperature=350*unit.kelvin)
>>> integrator.getTemperature()
Quantity(value=350.0, unit=kelvin)
>>> integrator.setTemperature(380.0*unit.kelvin)
>>> integrator.getTemperature()
Quantity(value=380.0, unit=kelvin)
>>> integrator.getGlobalVariableByName('kT')
3.1594995390636815

Notice that a CustomIntegrator loses any extra method after a serialization cycle.

>>> integrator_serialization = openmm.XmlSerializer.serialize(integrator)
>>> deserialized_integrator = openmm.XmlSerializer.deserialize(integrator_serialization)
>>> deserialized_integrator.getTemperature()
Traceback (most recent call last):
...
AttributeError: type object 'object' has no attribute '__getattr__'

We can restore the original interface with a class method

>>> ThermostatedIntegrator.restore_interface(integrator)
True
>>> integrator.getTemperature()
Quantity(value=380.0, unit=kelvin)
>>> integrator.setTemperature(400.0*unit.kelvin)
>>> isinstance(integrator, TestIntegrator)
True
Attributes:
global_variable_names

The set of global variable names defined for this integrator.

kT

The thermal energy in openmm.Quantity

thisown

The membership flag

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.

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.

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__(temperature, *args, **kwargs)[source]

Methods

__init__(temperature, *args, **kwargs)

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.

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.

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

global_variable_names

The set of global variable names defined for this integrator.

kT

The thermal energy in openmm.Quantity

thisown

The membership flag