ATMForce¶
- class openmm.openmm.ATMForce(*args)¶
The ATMForce class implements the Alchemical Transfer Method (ATM) for OpenMM. ATM is used to compute the binding free energies of molecular complexes and of other equilibrium processes. ATM and its implementation are described in the open access article:
Solmaz Azimi, Sheenam Khuttan, Joe Z. Wu, Rajat K. Pal, and Emilio Gallicchio. Relative Binding Free Energy Calculations for Ligands with Diverse Scaffolds with the Alchemical Transfer Method. J. Chem. Inf. Model. 62, 309 (2022)
Refer to the publication above for a detailed description of the ATM method and the parameters used in this API and please cite it to support our work if you use this software in your research.
The ATMForce implements an arbitrary potential energy function that depends on the potential energies (u0 and u1) of the system before and after a set of atoms are displaced by a specified amount. For example, you might displace a molecule from the solvent bulk to a receptor binding site to simulate a binding process. The potential energy function typically also depends on one or more parameters that are dialed to implement alchemical transformations.
To use this class, create an ATMForce object, passing an algebraic expression to the constructor that defines the potential energy. This expression can be any combination of the variables u0 and u1. Then call addGlobalParameter() to define the parameters on which the potential energy expression depends. The values of global parameters may be modified during a simulation by calling Context::setParameter(). Next, call addForce() to add Force objects that define the terms of the potential energy function that change upon displacement. Finally, call addParticle() to specify the displacement applied to each particle. Displacements can be changed by calling setParticleParameters(). As any per-particle parameters, changes in displacements take effect only after calling updateParametersInContext().
As an example, the following code creates an ATMForce based on the change in energy of two particles when the second particle is displaced by 1 nm in the x direction. The energy change is dialed using an alchemical parameter Lambda, which in this case is set to 1/2:
ATMForce *atmforce = new ATMForce("u0 + Lambda*(u1 - u0)"); atm->addGlobalParameter("Lambda", 0.5); atm->addParticle(Vec3(0, 0, 0)); atm->addParticle(Vec3(1, 0, 0)); CustomBondForce* force = new CustomBondForce("0.5*r^2"); atm->addForce(force);
Expressions may involve the operators + (add), - (subtract), * (multiply), / (divide), and ^ (power), and the following functions: sqrt, exp, log, sin, cos, sec, csc, tan, cot, asin, acos, atan, atan2, sinh, cosh, tanh, erf, erfc, min, max, abs, floor, ceil, step, delta, select. All trigonometric functions are defined in radians, and log is the natural logarithm. step(x) = 0 if x is less than 0, 1 otherwise. delta(x) = 1 if x is 0, 0 otherwise. select(x,y,z) = z if x = 0, y otherwise.
If instead of the energy expression the ATMForce constructor specifies the values of a series of parameters, the default energy expression is used:
select(step(Direction), u0, u1) + ((Lambda2-Lambda1)/Alpha)*log(1+exp(-Alpha*(usc-Uh))) + Lambda2*usc + W0; usc = select(step(u-Ubcore), (Umax-Ubcore)*fsc+Ubcore, u), u); fsc = (z^Acore-1)/(z^Acore+1); z = 1 + 2*(y/Acore) + 2*(y/Acore)^2; y = (u-Ubcore)/(Umax-Ubcore); u = select(step(Direction), 1, -1)*(u1-u0)
which is the same as the soft-core softplus alchemical potential energy function in the Azimi et al. paper above.
The ATMForce is then added to the System as any other Force
system.addForce(atmforce);
after which it will be used for energy/force evaluations for molecular dynamics and energy optimization. You can call getPerturbationEnergy() to query the values of u0 and u1, which are needed for computing free energies.
In most cases, particles are only displaced in one of the two states evaluated by this force. It computes the change in energy between the current particle coordinates (as stored in the Context) and the displaced coordinates. In some cases, it is useful to apply displacements to both states. You can do this by providing two displacement vectors to addParticle():
atm->addParticle(Vec3(1, 0, 0), Vec3(-1, 0, 0));
In this case, u1 will be computed after displacing the particle in the positive x direction, and u0 will be computed after displacing it in the negative x direction.
This class also has the ability to compute derivatives of the potential energy with respect to global parameters. Call addEnergyParameterDerivative() to request that the derivative with respect to a particular parameter be computed. You can then query its value in a Context by calling getState() on it.
- __init__(self, energy) → ATMForce¶
- __init__(self, lambda1, lambda2, alpha, uh, w0, umax, ubcore, acore, direction) → ATMForce
- __init__(self, other) → ATMForce
Create an ATMForce object with the default softplus energy expression. The values passed to this constructor are the default values of the global parameters for newly created Contexts. Their values can be changed by calling setParameter() on the Context using the parameter names defined by the Lambda1(), Lambda2(), etc. methods below.
- Parameters
lambda1 (double) – the default value of the Lambda1 parameter (dimensionless). This should be a number between 0 and 1.
lambda2 (double) – the default value of the Lambda2 parameter (dimensionless). This should be a number between 0 and 1.
alpha (double) – the default value of the Alpha parameter (kJ/mol)^-1
uh (double) – the default value of the Uh parameter (kJ/mol)
w0 (double) – the default value of the W0 parameter (kJ/mol)
umax (double) – the default value of the Umax parameter (kJ/mol)
ubcore (double) – the default value of the Ubcore parameter (kJ/mol)
acore (double) – the default value of the Acore parameter dimensionless)
direction (double) – the default value of the Direction parameter (dimensionless). This should be either 1 for the forward transfer, or -1 for the backward transfer.
Methods
Acore
()Returns the name of the global parameter corresponding to acore.
Alpha
()Returns the name of the global parameter corresponding to alpha.
Returns the name of the global parameter corresponding to direction.
Lambda1
()Returns the name of the global parameter corresponding to lambda1.
Lambda2
()Returns the name of the global parameter corresponding to lambda2.
Ubcore
()Returns the name of the global parameter corresponding to ubcore.
Uh
()Returns the name of the global parameter corresponding to uh.
Umax
()Returns the name of the global parameter corresponding to umax.
W0
()Returns the name of the global parameter corresponding to w0.
__init__
(-> ATMForce -> ATMForce)Create an ATMForce object with the default softplus energy expression.
addEnergyParameterDerivative
(self, name)Request that this Force compute the derivative of its energy with respect to a global parameter.
addForce
(self, force)Add a Force whose energy will be computed by the ATMForce.
addGlobalParameter
(self, name, defaultValue)Add a new global parameter that the interaction may depend on.
addParticle
(self, displacement1[, displacement0])Add a particle to the force.
getEnergyFunction
(self)Get the algebraic expression that gives the energy of the system
getEnergyParameterDerivativeName
(self, index)Get the name of a global parameter with respect to which this Force should compute the derivative of the energy.
getForce
(self, index)return the force from index
getForceGroup
(self)Get the force group this Force belongs to.
getGlobalParameterDefaultValue
(self, index)Get the default value of a global parameter.
getGlobalParameterName
(self, index)Get the name of a global parameter.
getName
(self)Get the name of this Force.
Get the number of global parameters with respect to which the derivative of the energy should be computed.
getNumForces
(self)Get the number of Forces included in the ATMForce.
getNumGlobalParameters
(self)Get the number of global parameters that the interaction depends on.
getNumParticles
(self)Get the number of particles managed by ATMForce.
getParticleParameters
(self, index)Get the parameters for a particle
getPerturbationEnergy
(self, context)Returns the current perturbation energy.
setEnergyFunction
(self, energy)Set the algebraic expression that gives the energy of the system
setForceGroup
(self, group)Set the force group this Force belongs to.
setGlobalParameterDefaultValue
(self, index, …)Set the default value of a global parameter.
setGlobalParameterName
(self, index, name)Set the name of a global parameter.
setName
(self, name)Set the name of this Force.
setParticleParameters
(self, index, displacement1)Set the parameters for a particle
updateParametersInContext
(self, context)Update the per-particle parameters in a Context to match those stored in this Force object.
Returns whether or not this force makes use of periodic boundary conditions.
Attributes
The membership flag
- property thisown¶
The membership flag
- static Lambda1() → std::string const &¶
Returns the name of the global parameter corresponding to lambda1. The value assigned to this parameter should be a number between 0 and 1.
- static Lambda2() → std::string const &¶
Returns the name of the global parameter corresponding to lambda2. The value assigned to this parameter should be a number between 0 and 1.
- static Alpha() → std::string const &¶
Returns the name of the global parameter corresponding to alpha. The value assigned to this parameter should be in units of (kJ/mol)^-1.
- static Uh() → std::string const &¶
Returns the name of the global parameter corresponding to uh. The value assigned to this parameter should be in units of (kJ/mol).
- static W0() → std::string const &¶
Returns the name of the global parameter corresponding to w0. The value assigned to this parameter should be in units of (kJ/mol).
- static Umax() → std::string const &¶
Returns the name of the global parameter corresponding to umax. The value assigned to this parameter should be in units of (kJ/mol).
- static Ubcore() → std::string const &¶
Returns the name of the global parameter corresponding to ubcore. The value assigned to this parameter should be in units of (kJ/mol).
- static Acore() → std::string const &¶
Returns the name of the global parameter corresponding to acore.
- static Direction() → std::string const &¶
Returns the name of the global parameter corresponding to direction. The value assigned to this parameter should be either 1 for the forward transfer, or -1 for the backward transfer.
- getNumParticles(self) → int¶
Get the number of particles managed by ATMForce.
This should be the same number of particles as the System
- getNumForces(self) → int¶
Get the number of Forces included in the ATMForce.
- getNumGlobalParameters(self) → int¶
Get the number of global parameters that the interaction depends on.
- getNumEnergyParameterDerivatives(self) → int¶
Get the number of global parameters with respect to which the derivative of the energy should be computed.
- getEnergyFunction(self) → std::string const &¶
Get the algebraic expression that gives the energy of the system
- setEnergyFunction(self, energy)¶
Set the algebraic expression that gives the energy of the system
- addForce(self, force) → int¶
Add a Force whose energy will be computed by the ATMForce.
- Parameters
force (Force *) – the Force to the be added, which should have been created on the heap with the “new” operator. The ATMForce takes over ownership of it, and deletes the Force when the ATMForce itself is deleted.
- Returns
The index within ATMForce of the force that was added
- Return type
int
- addParticle(self, displacement1, displacement0=Vec3()) → int¶
Add a particle to the force.
All of the particles in the System must be added to the ATMForce in the same order as they appear in the System.
- getParticleParameters(self, index)¶
Get the parameters for a particle
- setParticleParameters(self, index, displacement1, displacement0=Vec3())¶
Set the parameters for a particle
- addGlobalParameter(self, name, defaultValue) → int¶
Add a new global parameter that the interaction may depend on. The default value provided to this method is the initial value of the parameter in newly created Contexts. You can change the value at any time by calling setParameter() on the Context.
- Parameters
name (string) – the name of the parameter
defaultValue (double) – the default value of the parameter
- Returns
the index of the parameter that was added
- Return type
int
- getGlobalParameterName(self, index) → std::string const &¶
Get the name of a global parameter.
- Parameters
index (int) – the index of the parameter for which to get the name
- Returns
the parameter name
- Return type
string
- setGlobalParameterName(self, index, name)¶
Set the name of a global parameter.
- Parameters
index (int) – the index of the parameter for which to set the name
name (string) – the name of the parameter
- getGlobalParameterDefaultValue(self, index) → double¶
Get the default value of a global parameter.
- Parameters
index (int) – the index of the parameter for which to get the default value
- Returns
the parameter default value
- Return type
double
- setGlobalParameterDefaultValue(self, index, defaultValue)¶
Set the default value of a global parameter.
- Parameters
index (int) – the index of the parameter for which to set the default value
defaultValue (double) – the default value of the parameter
- addEnergyParameterDerivative(self, name)¶
Request that this Force compute the derivative of its energy with respect to a global parameter. The parameter must have already been added with addGlobalParameter().
- Parameters
name (string) – the name of the parameter
- getEnergyParameterDerivativeName(self, index) → std::string const &¶
Get the name of a global parameter with respect to which this Force should compute the derivative of the energy.
- Parameters
index (int) – the index of the parameter derivative, between 0 and getNumEnergyParameterDerivatives()
- Returns
the parameter name
- Return type
string
- updateParametersInContext(self, context)¶
Update the per-particle parameters in a Context to match those stored in this Force object. This method should be called after updating parameters with setParticleParameters() to copy them over to the Context. The only information this method updates is the values of per-particle parameters. The number of particles cannot be changed.
- usesPeriodicBoundaryConditions(self) → bool¶
Returns whether or not this force makes use of periodic boundary conditions.
- getPerturbationEnergy(self, context)¶
Returns the current perturbation energy.
- Parameters
context (Context) – the Context for which to return the energy
u1 (double) – on exit, the energy of the displaced state
u0 (double) – on exit, the energy of the non-displaced state
energy (double) – on exit, the value of this force’s energy function
- getForceGroup(self) → int¶
Get the force group this Force belongs to.
- getName(self) → std::string const &¶
Get the name of this Force. This is an arbitrary, user modifiable identifier. By default it equals the class name, but you can change it to anything useful.
- setForceGroup(self, group)¶
Set the force group this Force belongs to.
- Parameters
group (int) – the group index. Legal values are between 0 and 31 (inclusive).
- setName(self, name)¶
Set the name of this Force. This is an arbitrary, user modifiable identifier. By default it equals the class name, but you can change it to anything useful.