Class ATMForce
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) https://doi.org/10.1021/acs.jcim.1c01129
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(new OpenMM_Vec3(0, 0, 0)); atm.addParticle(new OpenMM_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(new OpenMM_Vec3(1, 0, 0), new OpenMM_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.
-
Field Summary
-
Constructor Summary
Constructors -
Method Summary
Modifier and TypeMethodDescriptionvoid
addEnergyParameterDerivative
(com.sun.jna.Pointer name) Request that this Force compute the derivative of its energy with respect to a global parameter.void
Request that this Force compute the derivative of its energy with respect to a global parameter.int
Add a Force whose energy will be computed by the ATMForce.int
addGlobalParameter
(com.sun.jna.Pointer name, double defaultValue) Add a new global parameter that the interaction may depend on.int
addGlobalParameter
(String name, double defaultValue) Add a new global parameter that the interaction may depend on.int
addParticle
(edu.uiowa.jopenmm.OpenMM_Vec3 displacement1, edu.uiowa.jopenmm.OpenMM_Vec3 displacement0) Add a particle to the force.void
destroy()
Destroy the force.Get the energy function.getEnergyParameterDerivativeName
(int index) Get the name of an energy parameter derivative.com.sun.jna.ptr.PointerByReference
getForce
(int index) Get a force by index.double
getGlobalParameterDefaultValue
(int index) Get the default value of a global parameter.getGlobalParameterName
(int index) Get the name of a global parameter.int
Get the number of energy parameter derivatives.int
Get the number of forces.int
Get the number of global parameters.int
Get the number of particles.void
getParticleParameters
(int index, edu.uiowa.jopenmm.OpenMM_Vec3 displacement1, edu.uiowa.jopenmm.OpenMM_Vec3 displacement0) Get the parameters for a particlevoid
getPerturbationEnergy
(Context context, com.sun.jna.ptr.DoubleByReference u0, com.sun.jna.ptr.DoubleByReference u1, com.sun.jna.ptr.DoubleByReference energy) Returns the current perturbation energy.void
getPerturbationEnergy
(Context context, DoubleBuffer u0, DoubleBuffer u1, DoubleBuffer energy) Returns the current perturbation energy.void
setEnergyFunction
(com.sun.jna.Pointer energy) Set the energy function.void
setEnergyFunction
(String energy) Set the energy function.void
setGlobalParameterDefaultValue
(int index, double defaultValue) Set the default value of a global parameter.void
setGlobalParameterName
(int index, com.sun.jna.Pointer name) Set the name of a global parameter.void
setGlobalParameterName
(int index, String name) Set the name of a global parameter.void
setParticleParameters
(int index, edu.uiowa.jopenmm.OpenMM_Vec3 displacement1, edu.uiowa.jopenmm.OpenMM_Vec3 displacement0) Set the parameters for a particlevoid
updateParametersInContext
(Context context) Update the per-particle parameters in a Context to match those stored in this Force object.boolean
Returns whether or not this force makes use of periodic boundary conditions.Methods inherited from class ffx.openmm.Force
getForceGroup, getForceIndex, getName, getPointer, setForceGroup, setForceIndex, setName
-
Constructor Details
-
ATMForce
Create an ATMForce object.- Parameters:
energy
- an algebraic expression giving the energy of the system as a function of u0 and u1, the energies before and after displacement
-
ATMForce
public ATMForce(double lambda1, double lambda2, double alpha, double uh, double w0, double umax, double ubcore, double acore, double direction) 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
- the default value of the Lambda1 parameter (dimensionless). This should be a number between 0 and 1.lambda2
- the default value of the Lambda2 parameter (dimensionless). This should be a number between 0 and 1.alpha
- the default value of the Alpha parameter (kJ/mol)^-1uh
- the default value of the Uh parameter (kJ/mol)w0
- the default value of the W0 parameter (kJ/mol)umax
- the default value of the Umax parameter (kJ/mol)ubcore
- the default value of the Ubcore parameter (kJ/mol)acore
- the default value of the Acore parameter dimensionless)direction
- the default value of the Direction parameter (dimensionless). This should be either 1 for the forward transfer, or -1 for the backward transfer.
-
-
Method Details
-
addEnergyParameterDerivative
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
- the name of the parameter
-
addEnergyParameterDerivative
public void addEnergyParameterDerivative(com.sun.jna.Pointer 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
- the name of the parameter
-
addForce
Add a Force whose energy will be computed by the ATMForce.- Parameters:
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
-
addGlobalParameter
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
- the name of the parameterdefaultValue
- the default value of the parameter- Returns:
- the index of the parameter that was added
-
addGlobalParameter
public int addGlobalParameter(com.sun.jna.Pointer name, double defaultValue) 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
- the name of the parameterdefaultValue
- the default value of the parameter- Returns:
- the index of the parameter that was added
-
addParticle
public int addParticle(edu.uiowa.jopenmm.OpenMM_Vec3 displacement1, edu.uiowa.jopenmm.OpenMM_Vec3 displacement0) 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.
- Parameters:
displacement1
- the displacement of the particle for the target state in nmdisplacement0
- the displacement of the particle for the initial state in nm- Returns:
- the index of the particle that was added
-
destroy
public void destroy()Destroy the force. -
getEnergyFunction
Get the energy function.- Returns:
- The energy function.
-
getEnergyParameterDerivativeName
Get the name of an energy parameter derivative.- Parameters:
index
- The index of the parameter derivative.- Returns:
- The name of the parameter derivative.
-
getForce
public com.sun.jna.ptr.PointerByReference getForce(int index) Get a force by index.- Parameters:
index
- The index of the force.- Returns:
- The force at the specified index.
-
getGlobalParameterDefaultValue
public double getGlobalParameterDefaultValue(int index) Get the default value of a global parameter.- Parameters:
index
- The index of the parameter.- Returns:
- The default value of the parameter.
-
getGlobalParameterName
Get the name of a global parameter.- Parameters:
index
- The index of the parameter.- Returns:
- The name of the parameter.
-
getNumEnergyParameterDerivatives
public int getNumEnergyParameterDerivatives()Get the number of energy parameter derivatives.- Returns:
- The number of energy parameter derivatives.
-
getNumForces
public int getNumForces()Get the number of forces.- Returns:
- The number of forces.
-
getNumGlobalParameters
public int getNumGlobalParameters()Get the number of global parameters.- Returns:
- The number of global parameters.
-
getNumParticles
public int getNumParticles()Get the number of particles.- Returns:
- The number of particles.
-
getParticleParameters
public void getParticleParameters(int index, edu.uiowa.jopenmm.OpenMM_Vec3 displacement1, edu.uiowa.jopenmm.OpenMM_Vec3 displacement0) Get the parameters for a particle- Parameters:
index
- the index in the force for the particle for which to get parametersdisplacement1
- the displacement of the particle for the target state in nmdisplacement0
- the displacement of the particle for the initial state in nm
-
getPerturbationEnergy
public void getPerturbationEnergy(Context context, com.sun.jna.ptr.DoubleByReference u0, com.sun.jna.ptr.DoubleByReference u1, com.sun.jna.ptr.DoubleByReference energy) Returns the current perturbation energy.- Parameters:
context
- the Context for which to return the energyu0
- on exit, the energy of the non-displaced stateu1
- on exit, the energy of the displaced stateenergy
- on exit, the value of this force's energy function
-
getPerturbationEnergy
public void getPerturbationEnergy(Context context, DoubleBuffer u0, DoubleBuffer u1, DoubleBuffer energy) Returns the current perturbation energy.- Parameters:
context
- the Context for which to return the energyu0
- on exit, the energy of the non-displaced stateu1
- on exit, the energy of the displaced stateenergy
- on exit, the value of this force's energy function
-
setEnergyFunction
Set the energy function.- Parameters:
energy
- The energy function.
-
setEnergyFunction
public void setEnergyFunction(com.sun.jna.Pointer energy) Set the energy function.- Parameters:
energy
- The energy function.
-
setGlobalParameterDefaultValue
public void setGlobalParameterDefaultValue(int index, double defaultValue) Set the default value of a global parameter.- Parameters:
index
- The index of the parameter.defaultValue
- The default value of the parameter.
-
setGlobalParameterName
Set the name of a global parameter.- Parameters:
index
- The index of the parameter.name
- The name of the parameter.
-
setGlobalParameterName
public void setGlobalParameterName(int index, com.sun.jna.Pointer name) Set the name of a global parameter.- Parameters:
index
- The index of the parameter.name
- The name of the parameter.
-
setParticleParameters
public void setParticleParameters(int index, edu.uiowa.jopenmm.OpenMM_Vec3 displacement1, edu.uiowa.jopenmm.OpenMM_Vec3 displacement0) Set the parameters for a particle- Parameters:
index
- the index in the force of the particle for which to set parametersdisplacement1
- the displacement of the particle for the target state in nmdisplacement0
- the displacement of the particle for the initial state in nm
-
updateParametersInContext
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
public boolean usesPeriodicBoundaryConditions()Returns whether or not this force makes use of periodic boundary conditions.- Overrides:
usesPeriodicBoundaryConditions
in classForce
- Returns:
- True if the force uses periodic boundary conditions.
-