Class CustomManyParticleForce
Be aware that the cost of evaluating an N-particle interaction increases very rapidly with N. Values larger than N=3 are rarely used.
We refer to a set of particles for which the energy is being evaluated as p1, p2, p3, etc. The energy expression may depend on the following variables and functions:
- x1, y1, z1, x2, y2, z2, etc.: The x, y, and z coordinates of the particle positions. For example, x1 is the x coordinate of particle p1, and y3 is the y coordinate of particle p3.
- distance(p1, p2): the distance between particles p1 and p2 (where "p1" and "p2" may be replaced by the names of whichever particles you want to calculate the distance between).
- angle(p1, p2, p3): the angle formed by the three specified particles.
- dihedral(p1, p2, p3, p4): the dihedral angle formed by the four specified particles.
- arbitrary global and per-particle parameters that you define.
To use this class, create a CustomManyParticleForce object, passing an algebraic expression to the constructor that defines the interaction energy of each set of particles. Then call addPerParticleParameter() to define per-particle parameters, and addGlobalParameter() to define global parameters. The values of per-particle parameters are specified as part of the system definition, while values of global parameters may be modified during a simulation by calling Context::setParameter().
Next, call addParticle() once for each particle in the System to set the values of its per-particle parameters. The number of particles for which you set parameters must be exactly equal to the number of particles in the System, or else an exception will be thrown when you try to create a Context. After a particle has been added, you can modify its parameters by calling setParticleParameters(). This will have no effect on Contexts that already exist unless you call updateParametersInContext().
Multi-particle interactions can be very expensive to evaluate, so they are usually used with a cutoff distance. The exact interpretation of the cutoff depends on the permutation mode, as discussed below.
CustomManyParticleForce also lets you specify "exclusions", particular pairs of particles whose interactions should be omitted from force and energy calculations. This is most often used for particles that are bonded to each other. If you specify a pair of particles as an exclusion, all sets that include those two particles will be omitted.
As an example, the following code creates a CustomManyParticleForce that implements an Axilrod-Teller potential. This is an interaction between three particles that depends on all three distances and angles formed by the particles.
CustomManyParticleForce force = new CustomManyParticleForce(3,
"C*(1+3*cos(theta1)*cos(theta2)*cos(theta3))/(r12*r13*r23)^3;" +
"theta1=angle(p1,p2,p3); theta2=angle(p2,p3,p1); theta3=angle(p3,p1,p2);" +
"r12=distance(p1,p2); r13=distance(p1,p3); r23=distance(p2,p3)");
force.setPermutationMode(CustomManyParticleForce.SinglePermutation);
This force depends on one parameter, C. The following code defines it as a global parameter:
force.addGlobalParameter("C", 1.0);
Notice that the expression is symmetric with respect to the particles. It only depends on the products cos(theta1)*cos(theta2)*cos(theta3) and r12*r13*r23, both of which are unchanged if the labels p1, p2, and p3 are permuted. This is required because we specified SinglePermutation as the permutation mode. (This is the default, so we did not really need to set it, but doing so makes the example clearer.) In this mode, the expression is only evaluated once for each set of particles. No guarantee is made about which particle will be identified as p1, p2, etc. Therefore, the energy must be symmetric with respect to exchange of particles. Otherwise, the results would be undefined because permuting the labels would change the energy.
-
Field Summary
-
Constructor Summary
ConstructorsConstructorDescriptionCustomManyParticleForce
(int particlesPerSet, String energy) Create a CustomManyParticleForce. -
Method Summary
Modifier and TypeMethodDescriptionint
addExclusion
(int particle1, int particle2) Add an exclusion for a pair of particles.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
(com.sun.jna.ptr.PointerByReference parameters, int type) Add a particle to the Force.int
addPerParticleParameter
(com.sun.jna.Pointer name) Add a new per-particle parameter that the interaction may depend on.int
Add a new per-particle parameter that the interaction may depend on.int
addTabulatedFunction
(com.sun.jna.Pointer name, com.sun.jna.ptr.PointerByReference function) Add a tabulated function that may appear in the energy expression.int
addTabulatedFunction
(String name, com.sun.jna.ptr.PointerByReference function) Add a tabulated function that may appear in the energy expression.void
createExclusionsFromBonds
(com.sun.jna.ptr.PointerByReference bonds, int bondCutoff) Identify exclusions based on the molecular topology.void
destroy()
Destroy the force.double
Get the cutoff distance (in nm) being used for interactions.Get the algebraic expression that gives the interaction energy for each set of particles.void
getExclusionParticles
(int index, com.sun.jna.ptr.IntByReference particle1, com.sun.jna.ptr.IntByReference particle2) Get the particles in an exclusion.void
getExclusionParticles
(int index, IntBuffer particle1, IntBuffer particle2) Get the particles in an exclusion.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 method used for handling long range nonbonded interactions.int
Get the number of exclusions.int
Get the number of global parameters that the interaction depends on.int
Get the number of particles for which force field parameters have been defined.int
Get the number of particles involved in each interaction.int
Get the number of per-particle parameters that the interaction depends on.int
Get the number of tabulated functions that have been defined.void
getParticleParameters
(int index, com.sun.jna.ptr.PointerByReference parameters, com.sun.jna.ptr.IntByReference type) Get the parameters for a particle.void
getParticleParameters
(int index, com.sun.jna.ptr.PointerByReference parameters, IntBuffer type) Get the parameters for a particle.int
Get the permutation mode.getPerParticleParameterName
(int index) Get the name of a per-particle parameter.com.sun.jna.ptr.PointerByReference
getTabulatedFunction
(int index) Get a reference to a tabulated function that may appear in the energy expression.getTabulatedFunctionName
(int index) Get the name of a tabulated function that may appear in the energy expression.void
getTypeFilter
(int index, com.sun.jna.ptr.PointerByReference types) Get the type filter for the specified type.void
setCutoffDistance
(double distance) Set the cutoff distance (in nm) being used for interactions.void
setEnergyFunction
(com.sun.jna.Pointer energy) Set the algebraic expression that gives the interaction energy for each set of particles.void
setEnergyFunction
(String energy) Set the algebraic expression that gives the interaction energy for each set of particles.void
setExclusionParticles
(int index, int particle1, int particle2) Set the particles in an exclusion.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
setNonbondedMethod
(int method) Set the method used for handling long range nonbonded interactions.void
setParticleParameters
(int index, com.sun.jna.ptr.PointerByReference parameters, int type) Set the parameters for a particle.void
setPermutationMode
(int mode) Set the permutation mode.void
setPerParticleParameterName
(int index, com.sun.jna.Pointer name) Set the name of a per-particle parameter.void
setPerParticleParameterName
(int index, String name) Set the name of a per-particle parameter.void
setTypeFilter
(int index, com.sun.jna.ptr.PointerByReference types) Set the type filter for the specified type.void
updateParametersInContext
(Context context) Update the per-particle parameters and tabulated functions in a Context to match those stored in this Force object.boolean
Check if the force uses periodic boundary conditions.Methods inherited from class ffx.openmm.Force
getForceGroup, getForceIndex, getName, getPointer, setForceGroup, setForceIndex, setName
-
Constructor Details
-
CustomManyParticleForce
Create a CustomManyParticleForce.- Parameters:
particlesPerSet
- The number of particles involved in each interaction.energy
- The algebraic expression that gives the interaction energy for each set of particles.
-
-
Method Details
-
addExclusion
public int addExclusion(int particle1, int particle2) Add an exclusion for a pair of particles.- Parameters:
particle1
- The index of the first particle.particle2
- The index of the second particle.- Returns:
- The index of the exclusion that was added.
-
addGlobalParameter
Add a new global parameter that the interaction may depend on.- Parameters:
name
- The name of the parameter.defaultValue
- 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.- Parameters:
name
- The name of the parameter.defaultValue
- The default value of the parameter.- Returns:
- The index of the parameter that was added.
-
addParticle
public int addParticle(com.sun.jna.ptr.PointerByReference parameters, int type) Add a particle to the Force.- Parameters:
parameters
- The list of parameters for the new particle.type
- The particle type.- Returns:
- The index of the particle that was added.
-
addPerParticleParameter
Add a new per-particle parameter that the interaction may depend on.- Parameters:
name
- The name of the parameter.- Returns:
- The index of the parameter that was added.
-
addPerParticleParameter
public int addPerParticleParameter(com.sun.jna.Pointer name) Add a new per-particle parameter that the interaction may depend on.- Parameters:
name
- The name of the parameter.- Returns:
- The index of the parameter that was added.
-
addTabulatedFunction
Add a tabulated function that may appear in the energy expression.- Parameters:
name
- The name of the function as it appears in expressions.function
- A TabulatedFunction object defining the function.- Returns:
- The index of the function that was added.
-
addTabulatedFunction
public int addTabulatedFunction(com.sun.jna.Pointer name, com.sun.jna.ptr.PointerByReference function) Add a tabulated function that may appear in the energy expression.- Parameters:
name
- The name of the function as it appears in expressions.function
- A TabulatedFunction object defining the function.- Returns:
- The index of the function that was added.
-
createExclusionsFromBonds
public void createExclusionsFromBonds(com.sun.jna.ptr.PointerByReference bonds, int bondCutoff) Identify exclusions based on the molecular topology.- Parameters:
bonds
- The set of bonds based on which to construct exclusions.bondCutoff
- Pairs of particles that are separated by this many bonds or fewer are added as exclusions.
-
destroy
public void destroy()Destroy the force. -
getCutoffDistance
public double getCutoffDistance()Get the cutoff distance (in nm) being used for interactions.- Returns:
- The cutoff distance, measured in nm.
-
getEnergyFunction
Get the algebraic expression that gives the interaction energy for each set of particles.- Returns:
- The energy expression.
-
getExclusionParticles
public void getExclusionParticles(int index, com.sun.jna.ptr.IntByReference particle1, com.sun.jna.ptr.IntByReference particle2) Get the particles in an exclusion.- Parameters:
index
- The index of the exclusion for which to get particles.particle1
- The index of the first particle in the exclusion (output).particle2
- The index of the second particle in the exclusion (output).
-
getExclusionParticles
Get the particles in an exclusion.- Parameters:
index
- The index of the exclusion for which to get particles.particle1
- The index of the first particle in the exclusion (output).particle2
- The index of the second particle in the exclusion (output).
-
getGlobalParameterDefaultValue
public double getGlobalParameterDefaultValue(int index) Get the default value of a global parameter.- Parameters:
index
- The index of the parameter for which to get the default value.- Returns:
- The parameter default value.
-
getGlobalParameterName
Get the name of a global parameter.- Parameters:
index
- The index of the parameter for which to get the name.- Returns:
- The parameter name.
-
getNonbondedMethod
public int getNonbondedMethod()Get the method used for handling long range nonbonded interactions.- Returns:
- The nonbonded method.
-
getNumExclusions
public int getNumExclusions()Get the number of exclusions.- Returns:
- The number of exclusions.
-
getNumGlobalParameters
public int getNumGlobalParameters()Get the number of global parameters that the interaction depends on.- Returns:
- The number of parameters.
-
getNumParticles
public int getNumParticles()Get the number of particles for which force field parameters have been defined.- Returns:
- The number of particles.
-
getNumParticlesPerSet
public int getNumParticlesPerSet()Get the number of particles involved in each interaction.- Returns:
- The number of particles per set.
-
getNumPerParticleParameters
public int getNumPerParticleParameters()Get the number of per-particle parameters that the interaction depends on.- Returns:
- The number of parameters.
-
getNumTabulatedFunctions
public int getNumTabulatedFunctions()Get the number of tabulated functions that have been defined.- Returns:
- The number of functions.
-
getParticleParameters
public void getParticleParameters(int index, com.sun.jna.ptr.PointerByReference parameters, com.sun.jna.ptr.IntByReference type) Get the parameters for a particle.- Parameters:
index
- The index of the particle for which to get parameters.parameters
- The list of parameters (output).type
- The particle type (output).
-
getParticleParameters
public void getParticleParameters(int index, com.sun.jna.ptr.PointerByReference parameters, IntBuffer type) Get the parameters for a particle.- Parameters:
index
- The index of the particle for which to get parameters.parameters
- The list of parameters (output).type
- The particle type (output).
-
getPerParticleParameterName
Get the name of a per-particle parameter.- Parameters:
index
- The index of the parameter for which to get the name.- Returns:
- The parameter name.
-
getPermutationMode
public int getPermutationMode()Get the permutation mode.- Returns:
- The permutation mode.
-
getTabulatedFunction
public com.sun.jna.ptr.PointerByReference getTabulatedFunction(int index) Get a reference to a tabulated function that may appear in the energy expression.- Parameters:
index
- The index of the function to get.- Returns:
- The TabulatedFunction object defining the function.
-
getTabulatedFunctionName
Get the name of a tabulated function that may appear in the energy expression.- Parameters:
index
- The index of the function to get.- Returns:
- The name of the function as it appears in expressions.
-
getTypeFilter
public void getTypeFilter(int index, com.sun.jna.ptr.PointerByReference types) Get the type filter for the specified type.- Parameters:
index
- The particle type index.types
- The allowed types for interactions (output).
-
setCutoffDistance
public void setCutoffDistance(double distance) Set the cutoff distance (in nm) being used for interactions.- Parameters:
distance
- The cutoff distance, measured in nm.
-
setEnergyFunction
Set the algebraic expression that gives the interaction energy for each set of particles.- Parameters:
energy
- The energy expression.
-
setEnergyFunction
public void setEnergyFunction(com.sun.jna.Pointer energy) Set the algebraic expression that gives the interaction energy for each set of particles.- Parameters:
energy
- The energy expression.
-
setExclusionParticles
public void setExclusionParticles(int index, int particle1, int particle2) Set the particles in an exclusion.- Parameters:
index
- The index of the exclusion for which to set particles.particle1
- The index of the first particle in the exclusion.particle2
- The index of the second particle in the exclusion.
-
setGlobalParameterDefaultValue
public void setGlobalParameterDefaultValue(int index, double defaultValue) Set the default value of a global parameter.- Parameters:
index
- The index of the parameter for which to set the default value.defaultValue
- The default value of the parameter.
-
setGlobalParameterName
Set the name of a global parameter.- Parameters:
index
- The index of the parameter for which to set the name.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 for which to set the name.name
- The name of the parameter.
-
setNonbondedMethod
public void setNonbondedMethod(int method) Set the method used for handling long range nonbonded interactions.- Parameters:
method
- The nonbonded method.
-
setParticleParameters
public void setParticleParameters(int index, com.sun.jna.ptr.PointerByReference parameters, int type) Set the parameters for a particle.- Parameters:
index
- The index of the particle for which to set parameters.parameters
- The list of parameters for the particle.type
- The particle type.
-
setPerParticleParameterName
Set the name of a per-particle parameter.- Parameters:
index
- The index of the parameter for which to set the name.name
- The name of the parameter.
-
setPerParticleParameterName
public void setPerParticleParameterName(int index, com.sun.jna.Pointer name) Set the name of a per-particle parameter.- Parameters:
index
- The index of the parameter for which to set the name.name
- The name of the parameter.
-
setPermutationMode
public void setPermutationMode(int mode) Set the permutation mode.- Parameters:
mode
- The permutation mode.
-
setTypeFilter
public void setTypeFilter(int index, com.sun.jna.ptr.PointerByReference types) Set the type filter for the specified type.- Parameters:
index
- The particle type index.types
- The allowed types for interactions.
-
updateParametersInContext
Update the per-particle parameters and tabulated functions in a Context to match those stored in this Force object.- Parameters:
context
- The Context in which to update the parameters.
-
usesPeriodicBoundaryConditions
public boolean usesPeriodicBoundaryConditions()Check if the force uses periodic boundary conditions.- Overrides:
usesPeriodicBoundaryConditions
in classForce
- Returns:
- True if the force uses periodic boundary conditions.
-