HippoNonbondedForce

class openmm.openmm.HippoNonbondedForce(*args)

This class implements all nonbonded interactions in the HIPPO force field: electrostatics, induction, charge transfer, dispersion, and repulsion. Although some of these are conceptually distinct, they share parameters in common and are most efficiently computed together. For example, the same multipole definitions are used for both electrostatics and Pauli repulsion. Therefore, all of them are computed by a single Force object.

To use it, create a HippoNonbondedForce object, then call addParticle() once for each particle. After an entry has been added, you can modify its force field parameters by calling setParticleParameters(). This will have no effect on Contexts that already exist unless you call updateParametersInContext().

You also can specify “exceptions”, particular pairs of particles whose interactions should be reduced or completely omitted. Call addException() to define exceptions.

__init__(self)HippoNonbondedForce
__init__(self, other)HippoNonbondedForce

Create a HippoNonbondedForce.

Methods

__init__(-> HippoNonbondedForce)

Create a HippoNonbondedForce.

addException(self, particle1, particle2, …)

Add an interaction to the list of exceptions that should be calculated differently from other interactions.

addParticle(self, charge, dipole, …)

Add the nonbonded force parameters for a particle.

getCutoffDistance(self)

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

getDPMEParameters(self)

Get the parameters to use for dispersion PME calculations.

getDPMEParametersInContext(self, context)

Get the parameters being used for dispersion PME in a particular Context.

getEwaldErrorTolerance(self)

Get the error tolerance for Ewald summation.

getExceptionParameters(self, index)

Get the scale factors for an interaction that should be calculated differently from others.

getExtrapolationCoefficients(self)

Get the coefficients for the mu_0, mu_1, mu_2, …, mu_n terms in the extrapolation algorithm for induced dipoles.

getForceGroup(self)

Get the force group this Force belongs to.

getInducedDipoles(self, context)

Get the induced dipole moments of all particles.

getLabFramePermanentDipoles(self, context)

Get the fixed dipole moments of all particles in the global reference frame.

getName(self)

Get the name of this Force.

getNonbondedMethod(self)

Get the method used for handling long-range nonbonded interactions.

getNumExceptions(self)

Get the number of exceptions.

getNumParticles(self)

Get the number of particles in the potential function.

getPMEParameters(self)

Get the parameters to use for PME calculations.

getPMEParametersInContext(self, context)

Get the parameters being used for PME in a particular Context.

getParticleParameters(self, index)

Get the nonbonded force parameters for a particle.

getSwitchingDistance(self)

Get the distance at which the switching function begins to reduce the repulsion and charge transfer interactions.

setCutoffDistance(self, distance)

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

setDPMEParameters(self, alpha, nx, ny, nz)

Set the parameters to use for dispersion PME calculations.

setEwaldErrorTolerance(self, tol)

Get the error tolerance for Ewald summation.

setExceptionParameters(self, index, …)

Set the scale factors for an interaction that should be calculated differently from others.

setExtrapolationCoefficients(self, coefficients)

Set the coefficients for the mu_0, mu_1, mu_2, …, mu_n terms in the extrapolation algorithm for induced dipoles.

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.

setPMEParameters(self, alpha, nx, ny, nz)

Set the parameters to use for PME calculations.

setParticleParameters(self, index, charge, …)

Set the nonbonded force parameters for a particle.

setSwitchingDistance(self, distance)

Set the distance at which the switching function begins to reduce the repulsion and charge transfer interactions.

updateParametersInContext(self, context)

Update the particle and exception 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

Bisector

NoAxisType

NoCutoff

PME

ThreeFold

ZBisect

ZOnly

ZThenX

thisown

The membership flag

property thisown

The membership flag

getNumParticles(self)int

Get the number of particles in the potential function.

getNumExceptions(self)int

Get the number of exceptions.

getNonbondedMethod(self)OpenMM::HippoNonbondedForce::NonbondedMethod

Get the method used for handling long-range nonbonded interactions.

setNonbondedMethod(self, method)

Set the method used for handling long-range nonbonded interactions.

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

getSwitchingDistance(self)double

Get the distance at which the switching function begins to reduce the repulsion and charge transfer interactions. This must be less than the cutoff distance.

setSwitchingDistance(self, distance)

Set the distance at which the switching function begins to reduce the repulsion and charge transfer interactions. This must be less than the cutoff distance.

getExtrapolationCoefficients(self)vectord

Get the coefficients for the mu_0, mu_1, mu_2, …, mu_n terms in the extrapolation algorithm for induced dipoles.

setExtrapolationCoefficients(self, coefficients)

Set the coefficients for the mu_0, mu_1, mu_2, …, mu_n terms in the extrapolation algorithm for induced dipoles.

Parameters

coefficients (vector< double >) – a vector whose mth entry specifies the coefficient for mu_m. The length of this vector determines how many iterations are performed.

getPMEParameters(self)

Get the parameters to use for PME calculations. If alpha is 0 (the default), these parameters are ignored and instead their values are chosen based on the Ewald error tolerance.

Returns

  • alpha (double) – the separation parameter

  • nx (int) – the number of grid points along the X axis

  • ny (int) – the number of grid points along the Y axis

  • nz (int) – the number of grid points along the Z axis

getDPMEParameters(self)

Get the parameters to use for dispersion PME calculations. If alpha is 0 (the default), these parameters are ignored and instead their values are chosen based on the Ewald error tolerance.

Returns

  • alpha (double) – the separation parameter

  • nx (int) – the number of dispersion grid points along the X axis

  • ny (int) – the number of dispersion grid points along the Y axis

  • nz (int) – the number of dispersion grid points along the Z axis

setPMEParameters(self, alpha, nx, ny, nz)

Set the parameters to use for PME calculations. If alpha is 0 (the default), these parameters are ignored and instead their values are chosen based on the Ewald error tolerance.

Parameters
  • alpha (double) – the separation parameter

  • nx (int) – the number of grid points along the X axis

  • ny (int) – the number of grid points along the Y axis

  • nz (int) – the number of grid points along the Z axis

setDPMEParameters(self, alpha, nx, ny, nz)

Set the parameters to use for dispersion PME calculations. If alpha is 0 (the default), these parameters are ignored and instead their values are chosen based on the Ewald error tolerance.

Parameters
  • alpha (double) – the separation parameter

  • nx (int) – the number of grid points along the X axis

  • ny (int) – the number of grid points along the Y axis

  • nz (int) – the number of grid points along the Z axis

getPMEParametersInContext(self, context)

Get the parameters being used for PME in a particular Context. Because some platforms have restrictions on the allowed grid sizes, the values that are actually used may be slightly different from those specified with setPmeGridDimensions(), or the standard values calculated based on the Ewald error tolerance. See the manual for details.

Parameters

context (Context) – the Context for which to get the parameters

Returns

  • alpha (double) – the separation parameter

  • nx (int) – the number of grid points along the X axis

  • ny (int) – the number of grid points along the Y axis

  • nz (int) – the number of grid points along the Z axis

getDPMEParametersInContext(self, context)

Get the parameters being used for dispersion PME in a particular Context. Because some platforms have restrictions on the allowed grid sizes, the values that are actually used may be slightly different from those specified with setPMEParameters(), or the standard values calculated based on the Ewald error tolerance. See the manual for details.

Parameters

context (Context) – the Context for which to get the parameters

Returns

  • alpha (double) – the separation parameter

  • nx (int) – the number of grid points along the X axis

  • ny (int) – the number of grid points along the Y axis

  • nz (int) – the number of grid points along the Z axis

addParticle(self, charge, dipole, quadrupole, coreCharge, alpha, epsilon, damping, c6, pauliK, pauliQ, pauliAlpha, polarizability, axisType, multipoleAtomZ, multipoleAtomX, multipoleAtomY)int

Add the nonbonded force parameters for a particle. This should be called once for each particle in the System. When it is called for the i’th time, it specifies the parameters for the i’th particle.

Parameters
  • charge (double) – the particle’s charge

  • dipole (vector< double >) – the particle’s molecular dipole (vector of size 3)

  • quadrupole (vector< double >) – the particle’s molecular quadrupole (vector of size 9)

  • coreCharge (double) – the charge of the atomic core

  • alpha (double) – controls the width of the particle’s electron density

  • epsilon (double) – sets the magnitude of charge transfer

  • damping (double) – sets the length scale for charge transfer

  • c6 (double) – the coefficient of the dispersion interaction

  • pauliK (double) – the coefficient of the Pauli repulsion interaction

  • pauliQ (double) – the charge used in computing the Pauli repulsion interaction

  • pauliAlpha (double) – the width of the particle’s electron density for computing the Pauli repulsion interaction

  • polarizability (double) – atomic polarizability

  • axisType (int) – the particle’s axis type

  • multipoleAtomZ (int) – index of first atom used in defining the local coordinate system for multipoles

  • multipoleAtomX (int) – index of second atom used in defining the local coordinate system for multipoles

  • multipoleAtomY (int) – index of third atom used in defining the local coordinate system for multipoles

Returns

the index of the particle that was added

Return type

int

getParticleParameters(self, index)

Get the nonbonded force parameters for a particle.

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

  • charge (double) – the particle’s charge

  • dipole (vector< double >) – the particle’s molecular dipole (vector of size 3)

  • quadrupole (vector< double >) – the particle’s molecular quadrupole (vector of size 9)

  • coreCharge (double) – the charge of the atomic core

  • alpha (double) – controls the width of the particle’s electron density

  • epsilon (double) – sets the magnitude of charge transfer

  • damping (double) – sets the length scale for charge transfer

  • c6 (double) – the coefficient of the dispersion interaction

  • pauliK (double) – the coefficient of the Pauli repulsion interaction

  • pauliQ (double) – the charge used in computing the Pauli repulsion interaction

  • pauliAlpha (double) – the width of the particle’s electron density for computing the Pauli repulsion interaction

  • polarizability (double) – atomic polarizability

  • axisType (int) – the particle’s axis type

  • multipoleAtomZ (int) – index of first atom used in defining the local coordinate system for multipoles

  • multipoleAtomX (int) – index of second atom used in defining the local coordinate system for multipoles

  • multipoleAtomY (int) – index of third atom used in defining the local coordinate system for multipoles

setParticleParameters(self, index, charge, dipole, quadrupole, coreCharge, alpha, epsilon, damping, c6, pauliK, pauliQ, pauliAlpha, polarizability, axisType, multipoleAtomZ, multipoleAtomX, multipoleAtomY)

Set the nonbonded force parameters for a particle.

Parameters
  • index (int) – the index of the particle for which to set parameters

  • charge (double) – the particle’s charge

  • dipole (vector< double >) – the particle’s molecular dipole (vector of size 3)

  • quadrupole (vector< double >) – the particle’s molecular quadrupole (vector of size 9)

  • coreCharge (double) – the charge of the atomic core

  • alpha (double) – controls the width of the particle’s electron density

  • epsilon (double) – sets the magnitude of charge transfer

  • damping (double) – sets the length scale for charge transfer

  • c6 (double) – the coefficient of the dispersion interaction

  • pauliK (double) – the coefficient of the Pauli repulsion interaction

  • pauliQ (double) – the charge used in computing the Pauli repulsion interaction

  • pauliAlpha (double) – the width of the particle’s electron density for computing the Pauli repulsion interaction

  • polarizability (double) – atomic polarizability

  • axisType (int) – the particle’s axis type

  • multipoleAtomZ (int) – index of first atom used in defining the local coordinate system for multipoles

  • multipoleAtomX (int) – index of second atom used in defining the local coordinate system for multipoles

  • multipoleAtomY (int) – index of third atom used in defining the local coordinate system for multipoles

addException(self, particle1, particle2, multipoleMultipoleScale, dipoleMultipoleScale, dipoleDipoleScale, dispersionScale, repulsionScale, chargeTransferScale, replace=False)int

Add an interaction to the list of exceptions that should be calculated differently from other interactions. If all scale factors are set to 0, this will cause the interaction to be completely omitted from force and energy calculations.

Parameters
  • particle1 (int) – the index of the first particle involved in the interaction

  • particle2 (int) – the index of the second particle involved in the interaction

  • multipoleMultipoleScale (double) – the factor by which to scale the Coulomb interaction between fixed multipoles

  • dipoleMultipoleScale (double) – the factor by which to scale the Coulomb interaction between an induced dipole and a fixed multipole

  • dipoleDipoleScale (double) – the factor by which to scale the Coulomb interaction between induced dipoles

  • dispersionScale (double) – the factor by which to scale the dispersion interaction

  • repulsionScale (double) – the factor by which to scale the Pauli repulsion

  • chargeTransferScale (double) – the factor by which to scale the charge transfer interaction

  • replace (bool) – determines the behavior if there is already an exception for the same two particles. If true, the existing one is replaced. If false, an exception is thrown.

Returns

the index of the exception that was added

Return type

int

getExceptionParameters(self, index)

Get the scale factors for an interaction that should be calculated differently from others.

Parameters
  • index (int) – the index of the interaction for which to get parameters

  • particle1 (int) – the index of the first particle involved in the interaction

  • particle2 (int) – the index of the second particle involved in the interaction

  • multipoleMultipoleScale (double) – the factor by which to scale the Coulomb interaction between fixed multipoles

  • dipoleMultipoleScale (double) – the factor by which to scale the Coulomb interaction between an induced dipole and a fixed multipole

  • dipoleDipoleScale (double) – the factor by which to scale the Coulomb interaction between induced dipoles

  • dispersionScale (double) – the factor by which to scale the dispersion interaction

  • repulsionScale (double) – the factor by which to scale the Pauli repulsion

  • chargeTransferScale (double) – the factor by which to scale the charge transfer interaction

setExceptionParameters(self, index, particle1, particle2, multipoleMultipoleScale, dipoleMultipoleScale, dipoleDipoleScale, dispersionScale, repulsionScale, chargeTransferScale)

Set the scale factors for an interaction that should be calculated differently from others. If all scale factors are set to 0, this will cause the interaction to be completely omitted from force and energy calculations.

Parameters
  • index (int) – the index of the interaction for which to set parameters

  • particle1 (int) – the index of the first particle involved in the interaction

  • particle2 (int) – the index of the second particle involved in the interaction

  • multipoleMultipoleScale (double) – the factor by which to scale the Coulomb interaction between fixed multipoles

  • dipoleMultipoleScale (double) – the factor by which to scale the Coulomb interaction between an induced dipole and a fixed multipole

  • dipoleDipoleScale (double) – the factor by which to scale the Coulomb interaction between induced dipoles

  • dispersionScale (double) – the factor by which to scale the dispersion interaction

  • repulsionScale (double) – the factor by which to scale the Pauli repulsion

  • chargeTransferScale (double) – the factor by which to scale the charge transfer interaction

getEwaldErrorTolerance(self)double

Get the error tolerance for Ewald summation. This corresponds to the fractional error in the forces which is acceptable. This value is used to select the grid dimensions and separation (alpha) parameter so that the average error level will be less than the tolerance. There is not a rigorous guarantee that all forces on all atoms will be less than the tolerance, however.

This can be overridden by explicitly setting an alpha parameter and grid dimensions to use.

setEwaldErrorTolerance(self, tol)

Get the error tolerance for Ewald summation. This corresponds to the fractional error in the forces which is acceptable. This value is used to select the grid dimensions and separation (alpha) parameter so that the average error level will be less than the tolerance. There is not a rigorous guarantee that all forces on all atoms will be less than the tolerance, however.

This can be overridden by explicitly setting an alpha parameter and grid dimensions to use.

getLabFramePermanentDipoles(self, context)

Get the fixed dipole moments of all particles in the global reference frame.

Parameters

context (Context) – the Context for which to get the fixed dipoles

Returns

dipoles – the fixed dipole moment of particle i is stored into the i’th element

Return type

vector< Vec3 >

getInducedDipoles(self, context)

Get the induced dipole moments of all particles.

Parameters

context (Context) – the Context for which to get the induced dipoles

Returns

dipoles – the induced dipole moment of particle i is stored into the i’th element

Return type

vector< Vec3 >

updateParametersInContext(self, context)

Update the particle and exception 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.

This method has several limitations. The only information it updates is the parameters of particles and exceptions. All other aspects of the Force (the nonbonded method, the cutoff distance, etc.) are unaffected and can only be changed by reinitializing the Context. Furthermore, only the scale factors for an exception can be changed; the pair of particles involved in the exception cannot change. Finally, this method cannot be used to add new particles or exceptions, only to change the parameters of existing ones.

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.