QTBIntegrator

class openmm.openmm.QTBIntegrator(*args)

This integrator implements the Adaptive Quantum Thermal Bath (adQTB) algorithm as described in . This is a fast method for approximating nuclear quantum effects by applying a Langevin thermostat whose random force varies with frequency to match the expected energy of a quantum harmonic oscillator.

To compensate for zero point energy leakage, the spectrum of the random force is adjusted automatically. The trajectory is divided into short segments. At the end of each segment, the distribution of energy across frequencies is evaluated, and the friction coefficients used in generating the random force are adjusted to better match the target distribution. You can monitor this process by calling getAdaptedFriction(). It is important to equilibrate the simulation long enough for the friction coefficients to converge before beginning the production part of the simulation.

To make this process more robust, it is recommended to average the data over all particles that are expected to behave identically. To do this, you can optionally call setParticleType() to define a particle as having a particular type. The data for all particles of the same type is averaged and they are adjusted together. You also can set a different adaptation rate for each particle type. The more particles that are being averaged over, the higher the adaptation rate can reasonably be set, leading to faster adaptation.

Most properties of this integrator are fixed at Context creation time and can only be changed after that by reinitializing the Context. That includes the step size, friction coefficient, segment length, cutoff frequency, particle types, and adaptation rates. The only property of the integrator that can be freely changed in the middle of a simulation is the temperature, and even that requires some care. The new temperature will not take effect until the beginning of the next segment. Furthermore, changing the temperature is a potentially expensive operation, since it requires performing a calculation whose cost scales as the cube of the number of time steps in a segment.

One must be very careful when trying to compute velocity dependent thermodynamic quantities, such as the instantaneous temperature or instantaneous pressure. The standard calculations for these quantities assume the velocities follow a classical distribution. They do not produce correct results for an adQTB simulation, in which the velocities follow a quantum distribution.

__init__(self, temperature, frictionCoeff, stepSize) QTBIntegrator
__init__(self, other) QTBIntegrator

Create a QTBIntegrator.

Parameters:
  • temperature (float) – the temperature of the heat bath (in Kelvin)

  • frictionCoeff (float) – the friction coefficient which couples the system to the heat bath (in inverse picoseconds)

  • stepSize (float) – the step size with which to integrate the system (in picoseconds)

Methods

__init__(-> QTBIntegrator)

Create a QTBIntegrator.

getAdaptedFriction(self, particle)

Get the adapted friction coefficients for a particle.

getConstraintTolerance(self)

Get the distance tolerance within which constraints are maintained, as a fraction of the constrained distance.

getCutoffFrequency(self)

Get the cutoff frequency applied to the colored noise.

getDefaultAdaptationRate(self)

Get the default adaptation rate.

getFriction(self)

Get the friction coefficient which determines how strongly the system is coupled to the heat bath (in inverse ps).

getIntegrationForceGroups(self)

Get which force groups to use for integration.

getParticleTypes(self)

Get a map whose keys are particle indices and whose values are particle types.

getRandomNumberSeed(self)

Get the random number seed.

getSegmentLength(self)

Get the length of each segment used for adapting the noise spectrum.

getStepSize(self)

Get the size of each time step, in picoseconds.

getTemperature(self)

Get the temperature of the heat bath (in Kelvin).

getTypeAdaptationRates(self)

Get a map whose keys are particle types and whose values are adaptation rates.

setConstraintTolerance(self, tol)

Set the distance tolerance within which constraints are maintained, as a fraction of the constrained distance.

setCutoffFrequency(self, cutoff)

Set the cutoff frequency applied to the colored noise.

setDefaultAdaptationRate(self, rate)

Set the default adaptation rate.

setFriction(self, coeff)

Set the friction coefficient which determines how strongly the system is coupled to the heat bath (in inverse ps).

setIntegrationForceGroups(groups)

Set which force groups to use for integration.

setParticleType(self, index, type)

Set the type of a particle.

setRandomNumberSeed(self, seed)

Set the random number seed.

setSegmentLength(self, length)

Set the length of each segment used for adapting the noise spectrum.

setStepSize(self, size)

Set the size of each time step, in picoseconds.

setTemperature(self, temp)

Set the temperature of the heat bath (in Kelvin).

setTypeAdaptationRate(self, type, rate)

Set the adaptation rate to use for particles of a given type.

step(self, steps)

Advance a simulation through time by taking a series of time steps.

Attributes

thisown

The membership flag

property thisown

The membership flag

getTemperature(self) float

Get the temperature of the heat bath (in Kelvin).

Returns:

the temperature of the heat bath, measured in Kelvin

Return type:

float

setTemperature(self, temp)

Set the temperature of the heat bath (in Kelvin).

Parameters:

temp (float) – the temperature of the heat bath, measured in Kelvin

getFriction(self) float

Get the friction coefficient which determines how strongly the system is coupled to the heat bath (in inverse ps).

Returns:

the friction coefficient, measured in 1/ps

Return type:

float

setFriction(self, coeff)

Set the friction coefficient which determines how strongly the system is coupled to the heat bath (in inverse ps).

Parameters:

coeff (float) – the friction coefficient, measured in 1/ps

getSegmentLength(self) float

Get the length of each segment used for adapting the noise spectrum.

Returns:

the segment length, measured in ps

Return type:

float

setSegmentLength(self, length)

Set the length of each segment used for adapting the noise spectrum.

Parameters:

length (float) – the segment length, measured in ps

getCutoffFrequency(self) float

Get the cutoff frequency applied to the colored noise.

Returns:

the cutoff frequency, measured in 1/ps

Return type:

float

setCutoffFrequency(self, cutoff)

Set the cutoff frequency applied to the colored noise.

Parameters:

cutoff (float) – the cutoff frequency, measured in 1/ps

getParticleTypes(self) Mapping[int, int]

Get a map whose keys are particle indices and whose values are particle types. This contains only the particles that have been specifically set with setParticleType().

setParticleType(self, index, type)

Set the type of a particle. This is an arbitrary integer.

getDefaultAdaptationRate(self) float

Get the default adaptation rate. This is the rate used for particles whose rate has not been otherwise specified with setTypeAdaptationRate(). It determines how much the noise spectrum changes after each segment.

setDefaultAdaptationRate(self, rate)

Set the default adaptation rate. This is the rate used for particles whose rate has not been otherwise specified with setTypeAdaptationRate(). It determines how much the noise spectrum changes after each segment.

getTypeAdaptationRates(self) Mapping[int, float]

Get a map whose keys are particle types and whose values are adaptation rates. These are the rates used for particles whose rates have been specified with setTypeAdaptationRate(). The rate determines how much the noise spectrum changes after each segment.

setTypeAdaptationRate(self, type, rate)

Set the adaptation rate to use for particles of a given type. The rate determines how much the noise spectrum changes after each segment.

getAdaptedFriction(self, particle)

Get the adapted friction coefficients for a particle. The return value is a vector of length numFreq = (3*n+1)/2, where n is the number of time steps in a segment. Element i is the friction coefficient used in generating the random force with frequency i*pi/(numFreq*stepSize).

This method is guaranteed to return identical results for particles that have the same type.

Parameters:
  • particle (int) – the index of the particle for which to get the friction

  • friction (Sequence[float]) – the adapted friction coefficients used in generating the random force.

getRandomNumberSeed(self) int

Get the random number seed. See setRandomNumberSeed() for details.

setRandomNumberSeed(self, seed)

Set the random number seed. The precise meaning of this parameter is undefined, and is left up to each Platform to interpret in an appropriate way. It is guaranteed that if two simulations are run with different random number seeds, the sequence of random forces will be different. On the other hand, no guarantees are made about the behavior of simulations that use the same seed. In particular, Platforms are permitted to use non-deterministic algorithms which produce different results on successive runs, even if those runs were initialized identically.

If seed is set to 0 (which is the default value assigned), a unique seed is chosen when a Context is created from this Integrator. This is done to ensure that each Context receives unique random seeds without you needing to set them explicitly.

step(self, steps)

Advance a simulation through time by taking a series of time steps.

Parameters:

steps (int) – the number of time steps to take

getConstraintTolerance(self) float

Get the distance tolerance within which constraints are maintained, as a fraction of the constrained distance.

getIntegrationForceGroups(self) int

Get which force groups to use for integration. By default, all force groups are included. This is interpreted as a set of bit flags: the forces from group i will be included if (groups&(1<<i)) != 0.

getStepSize(self) float

Get the size of each time step, in picoseconds. If this integrator uses variable time steps, the size of the most recent step is returned.

Returns:

the step size, measured in ps

Return type:

float

setConstraintTolerance(self, tol)

Set the distance tolerance within which constraints are maintained, as a fraction of the constrained distance.

setIntegrationForceGroups(groups)

Set which force groups to use for integration. By default, all force groups are included.

Parameters:

groups (set or int) – a set of indices for which force groups to include when integrating the equations of motion. Alternatively, the groups can be passed as a single unsigned integer interpreted as a bitmask, in which case group i will be included if (groups&(1<<i)) != 0.

setStepSize(self, size)

Set the size of each time step, in picoseconds. If this integrator uses variable time steps, the effect of calling this method is undefined, and it may simply be ignored.

Parameters:

size (float) – the step size, measured in ps