openmmtools.states.CompoundThermodynamicState

class openmmtools.states.CompoundThermodynamicState(thermodynamic_state, composable_states)[source]

Thermodynamic state composed by multiple states.

Allows to extend a ThermodynamicState through composition rather than inheritance.

The class dynamically inherits from the ThermodynamicState object given in the constructor, and it preserves direct access to all its methods and attributes. It is compatible also with subclasses of ThermodynamicState, but it does not support objects which make use of __slots__.

It is the user’s responsibility to check that IComposableStates are compatible to each other (i.e. that they do not depend on and/or modify the same properties of the system). If this is not the case, consider merging them into a single IComposableStates. If an IComposableState needs to access properties of ThermodynamicState (e.g. temperature, pressure) consider extending it through normal inheritance.

It is not necessary to explicitly inherit from IComposableState for compatibility as long as all abstract methods are implemented. All its attributes and methods will still be directly accessible unless they are masked by the main ThermodynamicState or by a IComposableState that appeared before in the constructor argument composable_states.

After construction, changing the original thermodynamic_state or any of the composable_states changes the state of the compound state.

Parameters:
thermodynamic_state : ThermodynamicState

The main ThermodynamicState which holds the OpenMM system.

composable_states : list of IComposableState

Each element represent a portion of the overall thermodynamic state.

Examples

Create an alchemically modified system.

>>> from openmmtools import testsystems, alchemy
>>> factory = alchemy.AbsoluteAlchemicalFactory(consistent_exceptions=False)
>>> alanine_vacuum = testsystems.AlanineDipeptideVacuum().system
>>> alchemical_region = alchemy.AlchemicalRegion(alchemical_atoms=range(22))
>>> alanine_alchemical_system = factory.create_alchemical_system(reference_system=alanine_vacuum,
...                                                              alchemical_regions=alchemical_region)
>>> alchemical_state = alchemy.AlchemicalState.from_system(alanine_alchemical_system)

AlchemicalState implement the IComposableState interface, so it can be used with CompoundThermodynamicState. All the alchemical parameters are accessible through the compound state.

>>> from simtk import openmm, unit
>>> thermodynamic_state = ThermodynamicState(system=alanine_alchemical_system,
...                                                 temperature=300*unit.kelvin)
>>> compound_state = CompoundThermodynamicState(thermodynamic_state=thermodynamic_state,
...                                                    composable_states=[alchemical_state])
>>> compound_state.lambda_sterics
1.0
>>> compound_state.lambda_electrostatics
1.0

You can control the parameters in the OpenMM Context in this state by setting the state attributes.

>>> compound_state.lambda_sterics = 0.5
>>> integrator = openmm.VerletIntegrator(1.0*unit.femtosecond)
>>> context = compound_state.create_context(integrator)
>>> context.getParameter('lambda_sterics')
0.5
>>> compound_state.lambda_sterics = 1.0
>>> compound_state.apply_to_context(context)
>>> context.getParameter('lambda_sterics')
1.0
Attributes:
barostat

The barostat associated to the system.

beta

Thermodynamic beta in units of mole/energy.

default_box_vectors

The default box vectors of the System (read-only).

is_periodic

True if the system is in a periodic box (read-only).

kT

Thermal energy per mole.

n_particles

Number of particles (read-only).

pressure

Constant pressure of the thermodynamic state.

system

The system in this thermodynamic state.

temperature

Constant temperature of the thermodynamic state.

volume

Constant volume of the thermodynamic state (read-only).

Methods

apply_to_context(context) Apply this compound thermodynamic state to the context.
create_context(integrator[, platform]) Create a context in this ThermodynamicState.
get_system(**kwargs) Manipulate and return the system.
get_volume([ignore_ensemble]) Volume of the periodic box (read-only).
is_context_compatible(context) Check compatibility of the given context.
is_state_compatible(thermodynamic_state) Check compatibility between ThermodynamicStates.
reduced_potential(context_state) Reduced potential in this thermodynamic state.
reduced_potential_at_states(context, …) Efficiently compute the reduced potential for a list of compatible states.
set_system(system[, fix_state]) Allow to set the system and fix its thermodynamic state.
__init__(thermodynamic_state, composable_states)[source]

Initialize self. See help(type(self)) for accurate signature.

Methods

__init__(thermodynamic_state, composable_states) Initialize self.
apply_to_context(context) Apply this compound thermodynamic state to the context.
create_context(integrator[, platform]) Create a context in this ThermodynamicState.
get_system(**kwargs) Manipulate and return the system.
get_volume([ignore_ensemble]) Volume of the periodic box (read-only).
is_context_compatible(context) Check compatibility of the given context.
is_state_compatible(thermodynamic_state) Check compatibility between ThermodynamicStates.
reduced_potential(context_state) Reduced potential in this thermodynamic state.
reduced_potential_at_states(context, …) Efficiently compute the reduced potential for a list of compatible states.
set_system(system[, fix_state]) Allow to set the system and fix its thermodynamic state.

Attributes

barostat The barostat associated to the system.
beta Thermodynamic beta in units of mole/energy.
default_box_vectors The default box vectors of the System (read-only).
is_periodic True if the system is in a periodic box (read-only).
kT Thermal energy per mole.
n_particles Number of particles (read-only).
pressure Constant pressure of the thermodynamic state.
system The system in this thermodynamic state.
temperature Constant temperature of the thermodynamic state.
volume Constant volume of the thermodynamic state (read-only).
apply_to_context(context)[source]

Apply this compound thermodynamic state to the context.

barostat

The barostat associated to the system.

Note that this is only a copy of the barostat, and you will need to set back the ThermodynamicState.barostat property for the changes to take place internally. If the pressure is allowed to fluctuate, this is None. Normally, you should only need to access the pressure and temperature properties, but this allows you to modify other parameters of the MonteCarloBarostat (e.g. frequency) after initialization. Setting this to None will place the system in an NVT ensemble.

beta

Thermodynamic beta in units of mole/energy.

create_context(integrator, platform=None)

Create a context in this ThermodynamicState.

The context contains a copy of the system. If the integrator is coupled to a heat bath (e.g. LangevinIntegrator), the system in the context will not have a thermostat, and vice versa if the integrator is not thermostated the system in the context will have a thermostat.

An integrator is considered thermostated if it exposes a method getTemperature(). A CompoundIntegrator is considered coupled to a heat bath if at least one of its integrators is. An exception is raised if the integrator is thermostated at a temperature different from the thermodynamic state’s.

Parameters:
integrator : simtk.openmm.Integrator

The integrator to use for Context creation. The eventual heat bath temperature must be consistent with the thermodynamic state.

platform : simtk.openmm.Platform, optional

Platform to use. If None, OpenMM tries to select the fastest available platform. Default is None.

Returns:
context : simtk.openmm.Context

The created OpenMM Context object.

Raises:
ThermodynamicsError

If the integrator has a temperature different from this ThermodynamicState.

Examples

When passing an integrator that does not expose getter and setter for the temperature, the context will be created with a thermostat.

>>> from simtk import openmm, unit
>>> from openmmtools import testsystems
>>> toluene = testsystems.TolueneVacuum()
>>> state = ThermodynamicState(toluene.system, 300*unit.kelvin)
>>> integrator = openmm.VerletIntegrator(1.0*unit.femtosecond)
>>> context = state.create_context(integrator)
>>> system = context.getSystem()
>>> [force.__class__.__name__ for force in system.getForces()
...  if 'Thermostat' in force.__class__.__name__]
['AndersenThermostat']

The thermostat is removed if we choose an integrator coupled to a heat bath.

>>> del context  # Delete previous context to free memory.
>>> integrator = openmm.LangevinIntegrator(300*unit.kelvin, 5.0/unit.picosecond,
...                                        2.0*unit.femtosecond)
>>> context = state.create_context(integrator)
>>> system = context.getSystem()
>>> [force.__class__.__name__ for force in system.getForces()
...  if 'Thermostat' in force.__class__.__name__]
[]
default_box_vectors

The default box vectors of the System (read-only).

get_system(**kwargs)[source]

Manipulate and return the system.

With default arguments, this is equivalent as the system property. By setting the arguments it is possible to obtain a modified copy of the system without the thermostat or the barostat.

Parameters:
remove_thermostat : bool

If True, the system thermostat is removed.

remove_barostat : bool

If True, the system barostat is removed.

Returns:
system : simtk.openmm.System

The system of this ThermodynamicState.

get_volume(ignore_ensemble=False)

Volume of the periodic box (read-only).

Parameters:
ignore_ensemble : bool, optional

If True, the volume of the periodic box vectors is returned even if the volume fluctuates.

Returns:
volume : simtk.unit.Quantity

The volume of the periodic box (units of length^3) or None if the system is not periodic or allowed to fluctuate.

is_context_compatible(context)[source]

Check compatibility of the given context.

Parameters:
context : simtk.openmm.Context

The OpenMM context to test.

Returns:
is_compatible : bool

True if this ThermodynamicState can be applied to context.

is_periodic

True if the system is in a periodic box (read-only).

is_state_compatible(thermodynamic_state)

Check compatibility between ThermodynamicStates.

The state is compatible if Contexts created by thermodynamic_state can be set to this ThermodynamicState through apply_to_context. The property is symmetric and transitive.

This is faster than checking compatibility of a Context object through is_context_compatible, and it should be preferred when possible.

Parameters:
thermodynamic_state : ThermodynamicState

The thermodynamic state to test.

Returns:
is_compatible : bool

True if the context created by thermodynamic_state can be converted to this state through apply_to_context().

Examples

States in the same ensemble (NVT or NPT) are compatible.

>>> from simtk import unit
>>> from openmmtools import testsystems
>>> alanine = testsystems.AlanineDipeptideExplicit()
>>> state1 = ThermodynamicState(alanine.system, 273*unit.kelvin)
>>> state2 = ThermodynamicState(alanine.system, 310*unit.kelvin)
>>> state1.is_state_compatible(state2)
True

States in different ensembles are not compatible.

>>> state1.pressure = 1.0*unit.atmosphere
>>> state1.is_state_compatible(state2)
False

States that store different systems (that differ by more than barostat and thermostat pressure and temperature) are also not compatible.

>>> alanine_implicit = testsystems.AlanineDipeptideImplicit().system
>>> state_implicit = ThermodynamicState(alanine_implicit, 310*unit.kelvin)
>>> state2.is_state_compatible(state_implicit)
False
kT

Thermal energy per mole.

n_particles

Number of particles (read-only).

pressure

Constant pressure of the thermodynamic state.

If the pressure is allowed to fluctuate, this is None. Setting this will automatically add/configure a barostat to the system. If it is set to None, the barostat will be removed.

reduced_potential(context_state)

Reduced potential in this thermodynamic state.

Parameters:
context_state : SamplerState or simtk.openmm.Context

Carry the configurational properties of the system.

Returns:
u : float

The unit-less reduced potential, which can be considered to have units of kT.

Notes

The reduced potential is defined as in Ref. [1]

u = eta [U(x) + p V(x) + mu N(x)]

where the thermodynamic parameters are

eta = 1/(kB T) is the inverse temperature p is the pressure mu is the chemical potential

and the configurational properties are

x the atomic positions U(x) is the potential energy V(x) is the instantaneous box volume N(x) the numbers of various particle species (e.g. protons of

titratible groups)

References

[1] Shirts MR and Chodera JD. Statistically optimal analysis of equilibrium states. J Chem Phys 129:124105, 2008.

Examples

Compute the reduced potential of a water box at 298 K and 1 atm.

>>> from openmmtools import testsystems
>>> waterbox = testsystems.WaterBox(box_edge=20.0*unit.angstroms)
>>> system, positions = waterbox.system, waterbox.positions
>>> state = ThermodynamicState(system=waterbox.system,
...                            temperature=298.0*unit.kelvin,
...                            pressure=1.0*unit.atmosphere)
>>> integrator = openmm.VerletIntegrator(1.0*unit.femtosecond)
>>> context = state.create_context(integrator)
>>> context.setPositions(waterbox.positions)
>>> sampler_state = SamplerState.from_context(context)
>>> u = state.reduced_potential(sampler_state)

If the sampler state is incompatible, an error is raised

>>> incompatible_sampler_state = sampler_state[:-1]
>>> state.reduced_potential(incompatible_sampler_state)
Traceback (most recent call last):
...
openmmtools.states.ThermodynamicsError: The sampler state has a different number of particles.

In case a cached SamplerState containing the potential energy and the volume of the context is not available, the method accepts a Context object and compute them with Context.getState().

>>> u = state.reduced_potential(context)
classmethod reduced_potential_at_states(context, thermodynamic_states)

Efficiently compute the reduced potential for a list of compatible states.

The user is responsible to ensure that the given context is compatible with the thermodynamic states.

Parameters:
context : openmm.Context

The OpenMM Context object with box vectors and positions set.

thermodynamic_states : list of ThermodynamicState

The list of thermodynamic states at which to compute the reduced potential.

Returns:
reduced_potentials : list of float

The unit-less reduced potentials, which can be considered to have units of kT.

Raises:
ValueError

If the thermodynamic states are not compatible to each other.

set_system(system, fix_state=False)[source]

Allow to set the system and fix its thermodynamic state.

With default arguments, this is equivalent to assign the system property, which raise an error if the system is in a different thermodynamic state.

Parameters:
system : simtk.openmm.System

The system to set.

fix_state : bool, optional

The thermodynamic state of the state will be fixed by all the composable states. Default is False.

system

The system in this thermodynamic state.

The returned system is a copy and can be modified without altering the internal state of ThermodynamicState. In order to ensure a consistent thermodynamic state, the system has a Thermostat force. You can use get_system() to obtain a copy of the system without the thermostat. The method create_context() then takes care of removing the thermostat when an integrator with a coupled heat bath is used (e.g. LangevinIntegrator).

It can be set only to a system which is consistent with the current thermodynamic state. Use set_system() if you want to correct the thermodynamic state of the system automatically before assignment.

temperature

Constant temperature of the thermodynamic state.

volume

Constant volume of the thermodynamic state (read-only).

If the volume is allowed to fluctuate, or if the system is not in a periodic box this is None.