AmoebaVdwForce

class openmm.openmm.AmoebaVdwForce(*args)

This class models van der Waals forces in the AMOEBA force field. It can use either buffered 14-7 potential or a Lennard-Jones 12-6 potential.

This class can operate in two different modes. In one mode, force field parameters are defined for each particle. When two particles interact, a combining rule is used to calculate the interaction parameters based on the parameters for the two particles. To use the class in this mode, call the version of addParticle() that takes sigma and epsilon values. It should be called once for each particle in the System.

In the other mode, each particle has a type index, and parameters are specified for each type rather than each individual particle. By default this mode also uses a combining rule, but you can override it by defining alternate parameters to use for specific pairs of particle types. To use the class in this mode, call the version of addParticle() that takes a type index. It should be called once for each particle in the System. You also must call addParticleType() once for each type. If you wish to override the combining for particular pairs of types, do so by calling addTypePair().

A unique feature of this class is that the interaction site for a particle does not need to be exactly at the particle’s location. Instead, it can be placed a fraction of the distance from that particle to another one. This is typically done for hydrogens to place the interaction site slightly closer to the parent atom. The fraction is known as the “reduction factor”, since it reduces the distance from the parent atom to the interaction site.

Support is also available for softcore interactions based on setting a per particle alchemical flag and setting the AmoebaVdwForce to use an “AlchemicalMethod” either Decouple or Annihilate. For Decouple, two alchemical atoms interact normally. For Annihilate, all interactions involving an alchemical atom are influenced. The softcore state is specified by setting a single Context parameter “AmoebaVdwLambda” between 0.0 and 1.0.

The softcore functional form can be modified by setting the softcore power (default of 5) and the softcore alpha (default of 0,7). For more information on the softcore functional form see Eq. 2 from: Jiao, D.; Golubkov, P. A.; Darden, T. A.; Ren, P., Calculation of protein-ligand binding free energy by using a polarizable potential. Proc. Natl. Acad. Sci. U.S.A. 2008, 105 (17), 6290-6295. .

__init__(self)AmoebaVdwForce
__init__(self, other)AmoebaVdwForce

Create an Amoeba VdwForce.

Methods

Lambda()

This is the name of the parameter which stores the current Amoeba vdW lambda value.

__init__(-> AmoebaVdwForce)

Create an Amoeba VdwForce.

addParticle(-> int)

Add the force field parameters for a vdw particle.

addParticleType(self, sigma, epsilon)

Add a particle type.

addTypePair(self, type1, type2, sigma, epsilon)

Add a type pair.

getAlchemicalMethod(self)

Get the method used for alchemical interactions.

getCutoff(self)

Get the cutoff distance.

getCutoffDistance(self)

Get the cutoff distance (in nm) being used for nonbonded interactions.

getEpsilonCombiningRule(self)

Get epsilon combining rule

getForceGroup(self)

Get the force group this Force belongs to.

getName(self)

Get the name of this Force.

getNonbondedMethod(self)

Get the method used for handling long range nonbonded interactions.

getNumParticleTypes(self)

Get the number of particle types.

getNumParticles(self)

Get the number of particles

getNumTypePairs(self)

Get the number of type pairs.

getParticleExclusions(self, particleIndex)

Get exclusions for specified particle

getParticleParameters(self, particleIndex)

Get the force field parameters for a vdw particle.

getParticleTypeParameters(self, typeIndex)

Get the force field parameters for a particle type.

getPotentialFunction(self)

Get the potential function to use.

getSigmaCombiningRule(self)

Get sigma combining rule

getSoftcoreAlpha(self)

Get the softcore alpha value.

getSoftcorePower(self)

Get the softcore power on lambda.

getTypePairParameters(self, pairIndex)

Get the force field parameters for a type pair.

getUseDispersionCorrection(self)

Get whether to add a contribution to the energy that approximately represents the effect of VdW interactions beyond the cutoff distance.

getUseParticleTypes(self)

Get whether parameters were specified by particle or by particle type.

setAlchemicalMethod(self, method)

Set the method used for handling long range nonbonded interactions.

setCutoff(self, cutoff)

Set the cutoff distance.

setCutoffDistance(self, distance)

Set the cutoff distance (in nm) being used for nonbonded interactions.

setEpsilonCombiningRule(self, …)

Set epsilon combining rule

setForceGroup(self, group)

Set the force group this Force belongs to.

setName(self, name)

Set the name of this Force.

setNonbondedMethod(self, method)

Set the method used for handling long range nonbonded interactions.

setParticleExclusions(self, particleIndex, …)

Set exclusions for specified particle

setParticleParameters(self, particleIndex, …)

Set the force field parameters for a vdw particle.

setParticleTypeParameters(self, typeIndex, …)

Set the force field parameters for a particle type.

setPotentialFunction(self, potential)

Set the potential function to use.

setSigmaCombiningRule(self, sigmaCombiningRule)

Set sigma combining rule

setSoftcoreAlpha(self, alpha)

Set the softcore alpha value (default = 0.7).

setSoftcorePower(self, n)

Set the softcore power on lambda (default = 5).

setTypePairParameters(self, pairIndex, …)

Set the force field parameters for a type pair.

setUseDispersionCorrection(self, useCorrection)

Set whether to add a contribution to the energy that approximately represents the effect of VdW interactions beyond the cutoff distance.

updateParametersInContext(self, context)

Update the per-particle parameters in a Context to match those stored in this Force object.

usesPeriodicBoundaryConditions(self)

Returns whether or not this force makes use of periodic boundary conditions.

Attributes

Annihilate

Buffered147

CutoffPeriodic

Decouple

LennardJones

NoCutoff

thisown

The membership flag

property thisown

The membership flag

static Lambda()std::string const &

This is the name of the parameter which stores the current Amoeba vdW lambda value.

getNumParticles(self)int

Get the number of particles

getNumParticleTypes(self)int

Get the number of particle types.

getNumTypePairs(self)int

Get the number of type pairs.

setParticleParameters(self, particleIndex, parentIndex, sigma, epsilon, reductionFactor, isAlchemical=False, typeIndex=- 1)

Set the force field parameters for a vdw particle.

Parameters
  • particleIndex (int) – the particle index

  • parentIndex (int) – the index of the parent particle

  • sigma (double) – vdw sigma

  • epsilon (double) – vdw epsilon

  • reductionFactor (double) – the fraction of the distance along the line from the parent particle to this particle at which the interaction site should be placed

  • isAlchemical (bool) – if true, this vdW particle is undergoing an alchemical change.

  • typeIndex (int) – the index of the particle type for this particle

getParticleParameters(self, particleIndex)

Get the force field parameters for a vdw particle.

Parameters

particleIndex (int) – the particle index

Returns

  • parentIndex (int) – the index of the parent particle

  • sigma (double) – vdw sigma

  • epsilon (double) – vdw epsilon

  • reductionFactor (double) – the fraction of the distance along the line from the parent particle to this particle at which the interaction site should be placed

  • isAlchemical (bool) – if true, this vdW particle is undergoing an alchemical change.

  • typeIndex (int) – the index of the particle type for this particle

addParticle(self, parentIndex, sigma, epsilon, reductionFactor, isAlchemical=False)int
addParticle(self, parentIndex, typeIndex, reductionFactor, isAlchemical=False)int

Add the force field parameters for a vdw particle. This version is used when parameters are defined by particle type.

Parameters
  • parentIndex (int) – the index of the parent particle

  • typeIndex (int) – the index of the particle type for this particle

  • reductionFactor (double) – the fraction of the distance along the line from the parent particle to this particle at which the interaction site should be placed

  • isAlchemical (bool) – if true, this vdW particle is undergoing an alchemical change.

Returns

index of added particle

Return type

int

addParticleType(self, sigma, epsilon)int

Add a particle type.

Parameters
  • sigma (double) – the sigma value for particles of this type

  • epsilon (double) – the epsilon value for particles of this type

Returns

the index of the particle type that was just added.

Return type

int

getParticleTypeParameters(self, typeIndex)

Get the force field parameters for a particle type.

Parameters

typeIndex (int) – the index of the particle type

Returns

  • sigma (double) – the sigma value for particles of this type

  • epsilon (double) – the epsilon value for particles of this type

setParticleTypeParameters(self, typeIndex, sigma, epsilon)

Set the force field parameters for a particle type.

Parameters
  • typeIndex (int) – the index of the particle type

  • sigma (double) – the sigma value for particles of this type

  • epsilon (double) – the epsilon value for particles of this type

addTypePair(self, type1, type2, sigma, epsilon)int

Add a type pair. This overrides the standard combining rule for interactions between particles of two particular types.

Parameters
  • type1 (int) – the index of the first particle type

  • type2 (int) – the index of the second particle type

  • sigma (double) – the sigma value for interactions between particles of these two types

  • epsilon (double) – the epsilon value for interactions between particles of these two types

Returns

the index of the type pair that was just added.

Return type

int

getTypePairParameters(self, pairIndex)

Get the force field parameters for a type pair. This overrides the standard combining rule for interactions between particles of two particular types.

Parameters

pairIndex (int) – the index of the type pair

Returns

  • type1 (int) – the index of the first particle type

  • type2 (int) – the index of the second particle type

  • sigma (double) – the sigma value for interactions between particles of these two types

  • epsilon (double) – the epsilon value for interactions between particles of these two types

setTypePairParameters(self, pairIndex, type1, type2, sigma, epsilon)

Set the force field parameters for a type pair. This overrides the standard combining rule for interactions between particles of two particular types.

Parameters
  • pairIndex (int) – the index of the type pair

  • type1 (int) – the index of the first particle type

  • type2 (int) – the index of the second particle type

  • sigma (double) – the sigma value for interactions between particles of these two types

  • epsilon (double) – the epsilon value for interactions between particles of these two types

setSigmaCombiningRule(self, sigmaCombiningRule)

Set sigma combining rule

Parameters

sigmaCombiningRule (string) – sigma combining rule: ‘ARITHMETIC’, ‘GEOMETRIC’. ‘CUBIC-MEAN’

getSigmaCombiningRule(self)std::string const &

Get sigma combining rule

Returns

sigmaCombiningRule sigma combining rule: ‘ARITHMETIC’, ‘GEOMETRIC’. ‘CUBIC-MEAN’

Return type

string

setEpsilonCombiningRule(self, epsilonCombiningRule)

Set epsilon combining rule

Parameters

epsilonCombiningRule (string) – epsilon combining rule: ‘ARITHMETIC’, ‘GEOMETRIC’. ‘HARMONIC’, ‘W-H’, ‘HHG’

getEpsilonCombiningRule(self)std::string const &

Get epsilon combining rule

Returns

epsilonCombiningRule epsilon combining rule: ‘ARITHMETIC’, ‘GEOMETRIC’. ‘HARMONIC’, ‘W-H’, ‘HHG’

Return type

string

getUseDispersionCorrection(self)bool

Get whether to add a contribution to the energy that approximately represents the effect of VdW interactions beyond the cutoff distance. The energy depends on the volume of the periodic box, and is only applicable when periodic boundary conditions are used. When running simulations at constant pressure, adding this contribution can improve the quality of results.

setUseDispersionCorrection(self, useCorrection)

Set whether to add a contribution to the energy that approximately represents the effect of VdW interactions beyond the cutoff distance. The energy depends on the volume of the periodic box, and is only applicable when periodic boundary conditions are used. When running simulations at constant pressure, adding this contribution can improve the quality of results.

getUseParticleTypes(self)bool

Get whether parameters were specified by particle or by particle type.

setParticleExclusions(self, particleIndex, exclusions)

Set exclusions for specified particle

Parameters
  • particleIndex (int) – particle index

  • exclusions (vector< int >) – vector of exclusions

getParticleExclusions(self, particleIndex)

Get exclusions for specified particle

Parameters

particleIndex (int) – particle index

Returns

exclusions – vector of exclusions

Return type

vector< int >

getCutoffDistance(self)double

Get the cutoff distance (in nm) being used for nonbonded interactions. If the NonbondedMethod in use is NoCutoff, this value will have no effect.

Returns

the cutoff distance, measured in nm

Return type

double

setCutoffDistance(self, distance)

Set the cutoff distance (in nm) being used for nonbonded interactions. If the NonbondedMethod in use is NoCutoff, this value will have no effect.

Parameters

distance (double) – the cutoff distance, measured in nm

setCutoff(self, cutoff)

Set the cutoff distance.

Deprecated

This method exists only for backward compatibility. Use setCutoffDistance() instead.

getCutoff(self)double

Get the cutoff distance.

Deprecated

This method exists only for backward compatibility. Use getCutoffDistance() instead.

getNonbondedMethod(self)OpenMM::AmoebaVdwForce::NonbondedMethod

Get the method used for handling long range nonbonded interactions.

setNonbondedMethod(self, method)

Set the method used for handling long range nonbonded interactions.

getPotentialFunction(self)OpenMM::AmoebaVdwForce::PotentialFunction

Get the potential function to use.

setPotentialFunction(self, potential)

Set the potential function to use.

setSoftcorePower(self, n)

Set the softcore power on lambda (default = 5).

getSoftcorePower(self)int

Get the softcore power on lambda.

setSoftcoreAlpha(self, alpha)

Set the softcore alpha value (default = 0.7).

getSoftcoreAlpha(self)double

Get the softcore alpha value.

getAlchemicalMethod(self)OpenMM::AmoebaVdwForce::AlchemicalMethod

Get the method used for alchemical interactions.

setAlchemicalMethod(self, method)

Set the method used for handling long range nonbonded interactions.

updateParametersInContext(self, context)

Update the per-particle parameters in a Context to match those stored in this Force object. This method provides an efficient method to update certain parameters in an existing Context without needing to reinitialize it. Simply call setParticleParameters() to modify this object’s parameters, then call updateParametersInContext() to copy them over to the Context.

The only information this method updates is the values of per-particle parameters. All other aspects of the Force (the nonbonded method, the cutoff distance, etc.) are unaffected and can only be changed by reinitializing the Context.

usesPeriodicBoundaryConditions(self)bool

Returns whether or not this force makes use of periodic boundary conditions.

Returns

true if nonbondedMethod uses PBC and false otherwise

Return type

bool

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.