RustOpenMM/openmm
1
0
Fork 0
mirror of https://github.com/openmm/openmm synced 2026-06-03 06:39:48 +09:00
Code Issues Packages Projects Releases Wiki Activity
Files
master
openmm/wrappers/python/tests/TestAPIUnits.py
Christian Clauss 82703dff73 Fix typos discovered by codespell in Python files (#5173)
2025-12-30 12:21:00 -08:00

1414 lines
65 KiB
Python
Raw Permalink Blame History

from __future__ import division
import unittest
from openmm import *
from openmm.app import *
from openmm.unit import *
import random
import math
class TestAPIUnits(unittest.TestCase):
"""Test the Simulation class"""
def assertAlmostEqualUnit(self, x1, x2):
self.assertAlmostEqual(x1._value, x2.value_in_unit(x1.unit))
def assertAlmostEqualUnitArray(self, a1, a2):
for x1, x2 in zip(a1._value, a2.value_in_unit(a1.unit)):
self.assertAlmostEqual(x1, x2)
def testHarmonicBondForce(self):
""" Tests HarmonicBondForce API features """
force = HarmonicBondForce()
force.addBond(0, 1, 1.0, 1.0)
force.addBond(2, 3, 1.0*angstroms,
1.0*kilocalories_per_mole/angstroms**2)
i, j, length, K = force.getBondParameters(0)
self.assertEqual(force.getNumBonds(), 2)
self.assertEqual(i, 0)
self.assertEqual(j, 1)
self.assertEqual(length, 1.0*nanometers)
self.assertEqual(K, 1*kilojoules_per_mole/nanometers**2)
self.assertIs(length.unit, nanometers)
self.assertIs(K.unit, kilojoules_per_mole/nanometers**2)
i, j, length, K = force.getBondParameters(1)
self.assertEqual(i, 2)
self.assertEqual(j, 3)
self.assertEqual(length, 1.0*angstroms)
self.assertAlmostEqualUnit(K, 1*kilocalories_per_mole/angstroms**2)
self.assertIs(length.unit, nanometers)
self.assertIs(K.unit, kilojoules_per_mole/nanometers**2)
def testHarmonicAngleForce(self):
""" Tests HarmonicAngleForce API features """
force = HarmonicAngleForce()
force.addAngle(0, 1, 2, 180*degrees,
1.0*kilocalories_per_mole/radians**2)
force.addAngle(1, 2, 3, math.pi/2, 1.0)
self.assertEqual(force.getNumAngles(), 2)
i, j, k, angle, K = force.getAngleParameters(0)
self.assertEqual(i, 0)
self.assertEqual(j, 1)
self.assertEqual(k, 2)
self.assertEqual(angle, 180*degrees)
self.assertEqual(K, 1.0*kilocalories_per_mole/radians**2)
self.assertIs(angle.unit, radians)
self.assertIs(K.unit, kilojoules_per_mole/radians**2)
i, j, k, angle, K = force.getAngleParameters(1)
self.assertEqual(i, 1)
self.assertEqual(j, 2)
self.assertEqual(k, 3)
self.assertEqual(angle, math.pi/2*radians)
self.assertEqual(K, 1.0*kilojoules_per_mole/radians**2)
self.assertIs(angle.unit, radians)
self.assertIs(K.unit, kilojoules_per_mole/radians**2)
def testPeriodicTorsionForce(self):
""" Tests PeriodicTorsionForce API features """
force = PeriodicTorsionForce()
force.addTorsion(0, 1, 2, 3, 1, math.pi, 1)
force.addTorsion(1, 2, 3, 4, 2, 180*degrees, 1*kilocalories_per_mole)
self.assertEqual(force.getNumTorsions(), 2)
i, j, k, l, per, phase, K = force.getTorsionParameters(0)
self.assertEqual(i, 0)
self.assertEqual(j, 1)
self.assertEqual(k, 2)
self.assertEqual(l, 3)
self.assertEqual(per, 1)
self.assertFalse(is_quantity(per))
self.assertEqual(phase, math.pi*radians)
self.assertEqual(K, 1*kilojoules_per_mole)
self.assertIs(phase.unit, radians)
self.assertIs(K.unit, kilojoules_per_mole)
i, j, k, l, per, phase, K = force.getTorsionParameters(1)
self.assertEqual(i, 1)
self.assertEqual(j, 2)
self.assertEqual(k, 3)
self.assertEqual(l, 4)
self.assertEqual(per, 2)
self.assertFalse(is_quantity(per))
self.assertEqual(phase, 180*degrees)
self.assertEqual(K, 1*kilocalories_per_mole)
self.assertIs(phase.unit, radians)
self.assertIs(K.unit, kilojoules_per_mole)
def testRBTorsionForce(self):
""" Tests the RBTorsionForce API features """
force = RBTorsionForce()
force.addTorsion(0, 1, 2, 3, 1, 2, 3, 4, 5, 6)
force.addTorsion(1, 2, 3, 4, 1*kilocalories_per_mole,
2*kilocalories_per_mole, 3*kilocalories_per_mole,
4*kilocalories_per_mole, 5*kilocalories_per_mole,
6*kilocalories_per_mole)
self.assertEqual(force.getNumTorsions(), 2)
i, j, k, l, c0, c1, c2, c3, c4, c5 = force.getTorsionParameters(0)
self.assertEqual(i, 0)
self.assertEqual(j, 1)
self.assertEqual(k, 2)
self.assertEqual(l, 3)
self.assertEqual(c0, 1*kilojoules_per_mole)
self.assertEqual(c1, 2*kilojoules_per_mole)
self.assertEqual(c2, 3*kilojoules_per_mole)
self.assertEqual(c3, 4*kilojoules_per_mole)
self.assertEqual(c4, 5*kilojoules_per_mole)
self.assertEqual(c5, 6*kilojoules_per_mole)
self.assertIs(c0.unit, kilojoules_per_mole)
self.assertIs(c1.unit, kilojoules_per_mole)
self.assertIs(c2.unit, kilojoules_per_mole)
self.assertIs(c3.unit, kilojoules_per_mole)
self.assertIs(c4.unit, kilojoules_per_mole)
self.assertIs(c5.unit, kilojoules_per_mole)
i, j, k, l, c0, c1, c2, c3, c4, c5 = force.getTorsionParameters(1)
self.assertEqual(i, 1)
self.assertEqual(j, 2)
self.assertEqual(k, 3)
self.assertEqual(l, 4)
self.assertAlmostEqualUnit(c0, 1*kilocalories_per_mole)
self.assertAlmostEqualUnit(c1, 2*kilocalories_per_mole)
self.assertAlmostEqualUnit(c2, 3*kilocalories_per_mole)
self.assertAlmostEqualUnit(c3, 4*kilocalories_per_mole)
self.assertAlmostEqualUnit(c4, 5*kilocalories_per_mole)
self.assertAlmostEqualUnit(c5, 6*kilocalories_per_mole)
self.assertIs(c0.unit, kilojoules_per_mole)
self.assertIs(c1.unit, kilojoules_per_mole)
self.assertIs(c2.unit, kilojoules_per_mole)
self.assertIs(c3.unit, kilojoules_per_mole)
self.assertIs(c4.unit, kilojoules_per_mole)
self.assertIs(c5.unit, kilojoules_per_mole)
def testNonbondedForce(self):
""" Tests the NonbondedForce API features """
force = NonbondedForce()
force.addParticle(1.0, 1.0, 1.0)
force.addParticle(1.0*coulombs, 1.0*angstroms,
1.0*kilocalories_per_mole)
self.assertEqual(force.getNumParticles(), 2)
charge, sigma, epsilon = force.getParticleParameters(0)
self.assertEqual(charge, 1.0*elementary_charge)
self.assertEqual(sigma, 1.0*nanometers)
self.assertEqual(epsilon, 1.0*kilojoules_per_mole)
self.assertIs(charge.unit, elementary_charge)
self.assertIs(sigma.unit, nanometers)
self.assertIs(epsilon.unit, kilojoules_per_mole)
charge, sigma, epsilon = force.getParticleParameters(1)
self.assertEqual(charge, 1.0*coulombs)
self.assertEqual(sigma, 1.0*angstroms)
self.assertEqual(epsilon, 1.0*kilocalories_per_mole)
self.assertIs(charge.unit, elementary_charge)
self.assertIs(sigma.unit, nanometers)
self.assertIs(epsilon.unit, kilojoules_per_mole)
force.setCutoffDistance(10*angstroms)
self.assertEqual(force.getCutoffDistance(), 1*nanometers)
force.setCutoffDistance(1)
self.assertEqual(force.getCutoffDistance(), 1*nanometers)
force.setSwitchingDistance(8*angstroms)
self.assertEqual(force.getSwitchingDistance(), 0.8*nanometers)
self.assertIs(force.getSwitchingDistance().unit, nanometer)
self.assertFalse(force.usesPeriodicBoundaryConditions())
force.setNonbondedMethod(NonbondedForce.NoCutoff)
self.assertFalse(force.usesPeriodicBoundaryConditions())
force.setNonbondedMethod(NonbondedForce.CutoffNonPeriodic)
self.assertFalse(force.usesPeriodicBoundaryConditions())
force.setNonbondedMethod(NonbondedForce.CutoffPeriodic)
self.assertTrue(force.usesPeriodicBoundaryConditions())
force.setNonbondedMethod(NonbondedForce.Ewald)
self.assertTrue(force.usesPeriodicBoundaryConditions())
force.setNonbondedMethod(NonbondedForce.PME)
self.assertTrue(force.usesPeriodicBoundaryConditions())
force.setNonbondedMethod(NonbondedForce.LJPME)
self.assertTrue(force.usesPeriodicBoundaryConditions())
def testConstantPotentialForce(self):
""" Tests the ConstantPotentialForce API features """
force = ConstantPotentialForce()
force.addParticle(1.0)
force.addParticle(1.0*coulomb)
self.assertEqual(force.getNumParticles(), 2)
charge = force.getParticleParameters(0)
self.assertAlmostEqualUnit(charge, 1.0*elementary_charge)
self.assertIs(charge.unit, elementary_charge)
charge = force.getParticleParameters(1)
self.assertAlmostEqualUnit(charge, 1.0*coulomb)
self.assertIs(charge.unit, elementary_charge)
force.setParticleParameters(0, 2.0*coulomb)
force.setParticleParameters(1, 2.0)
charge = force.getParticleParameters(0)
self.assertAlmostEqualUnit(charge, 2.0*coulomb)
charge = force.getParticleParameters(1)
self.assertAlmostEqualUnit(charge, 2.0*elementary_charge)
force.addException(0, 1, 1.0)
force.addException(1, 2, 1.0*coulomb*coulomb)
self.assertEqual(force.getNumExceptions(), 2)
index1, index2, chargeProd = force.getExceptionParameters(0)
self.assertEqual(index1, 0)
self.assertEqual(index2, 1)
self.assertAlmostEqualUnit(chargeProd, 1.0*elementary_charge*elementary_charge)
self.assertIs(chargeProd.unit, elementary_charge*elementary_charge)
index1, index2, chargeProd = force.getExceptionParameters(1)
self.assertEqual(index1, 1)
self.assertEqual(index2, 2)
self.assertAlmostEqualUnit(chargeProd, 1.0*coulomb*coulomb)
self.assertIs(chargeProd.unit, elementary_charge*elementary_charge)
force.setExceptionParameters(0, 3, 4, 2.0*coulomb*coulomb)
force.setExceptionParameters(1, 5, 6, 2.0)
index1, index2, chargeProd = force.getExceptionParameters(0)
self.assertEqual(index1, 3)
self.assertEqual(index2, 4)
self.assertAlmostEqualUnit(chargeProd, 2.0*coulomb*coulomb)
self.assertIs(chargeProd.unit, elementary_charge*elementary_charge)
index1, index2, chargeProd = force.getExceptionParameters(1)
self.assertEqual(index1, 5)
self.assertEqual(index2, 6)
self.assertAlmostEqualUnit(chargeProd, 2.0*elementary_charge*elementary_charge)
self.assertIs(chargeProd.unit, elementary_charge*elementary_charge)
force.addElectrode((0, 1, 2), 1.0, 1.0, 1.0)
force.addElectrode((3, 4, 5), 1.0*kilocalorie_per_mole/coulomb, 1.0*angstrom, 1.0/angstrom)
self.assertEqual(force.getNumElectrodes(), 2)
indices, potential, gaussianWidth, thomasFermiScale = force.getElectrodeParameters(0)
self.assertEqual(set(indices), {0, 1, 2})
self.assertAlmostEqualUnit(potential, 1.0*kilojoule_per_mole/elementary_charge)
self.assertAlmostEqualUnit(gaussianWidth, 1.0*nanometer)
self.assertAlmostEqualUnit(thomasFermiScale, 1.0/nanometer)
self.assertIs(potential.unit, kilojoule_per_mole/elementary_charge)
self.assertIs(gaussianWidth.unit, nanometer)
self.assertIs(thomasFermiScale.unit, nanometer**-1)
indices, potential, gaussianWidth, thomasFermiScale = force.getElectrodeParameters(1)
self.assertEqual(set(indices), {3, 4, 5})
self.assertAlmostEqualUnit(potential, 1.0*kilocalorie_per_mole/coulomb)
self.assertAlmostEqualUnit(gaussianWidth, 1.0*angstrom)
self.assertAlmostEqualUnit(thomasFermiScale, 1.0/angstrom)
self.assertIs(potential.unit, kilojoule_per_mole/elementary_charge)
self.assertIs(gaussianWidth.unit, nanometer)
self.assertIs(thomasFermiScale.unit, nanometer**-1)
force.setElectrodeParameters(0, (6, 7, 8), 2.0*kilocalorie_per_mole/coulomb, 2.0*angstrom, 2.0/angstrom)
force.setElectrodeParameters(1, (9, 10, 11), 2.0, 2.0, 2.0)
indices, potential, gaussianWidth, thomasFermiScale = force.getElectrodeParameters(0)
self.assertEqual(set(indices), {6, 7, 8})
self.assertAlmostEqualUnit(potential, 2.0*kilocalorie_per_mole/coulomb)
self.assertAlmostEqualUnit(gaussianWidth, 2.0*angstrom)
self.assertAlmostEqualUnit(thomasFermiScale, 2.0/angstrom)
indices, potential, gaussianWidth, thomasFermiScale = force.getElectrodeParameters(1)
self.assertEqual(set(indices), {9, 10, 11})
self.assertAlmostEqualUnit(potential, 2.0*kilojoule_per_mole/elementary_charge)
self.assertAlmostEqualUnit(gaussianWidth, 2.0*nanometer)
self.assertAlmostEqualUnit(thomasFermiScale, 2.0/nanometer)
force.setCutoffDistance(2.0)
cutoffDistance = force.getCutoffDistance()
self.assertAlmostEqualUnit(cutoffDistance, 2.0*nanometer)
self.assertIs(cutoffDistance.unit, nanometer)
force.setCutoffDistance(2.0*angstrom)
cutoffDistance = force.getCutoffDistance()
self.assertAlmostEqualUnit(cutoffDistance, 2.0*angstrom)
self.assertIs(cutoffDistance.unit, nanometer)
force.setEwaldErrorTolerance(1e-4)
self.assertEqual(force.getEwaldErrorTolerance(), 1e-4)
force.setEwaldErrorTolerance(2e-4)
self.assertEqual(force.getEwaldErrorTolerance(), 2e-4)
force.setPMEParameters(1.0, 10, 10, 10)
alpha, nx, ny, nz = force.getPMEParameters()
self.assertAlmostEqualUnit(alpha, 1.0/nanometer)
self.assertEqual(nx, 10)
self.assertEqual(ny, 10)
self.assertEqual(nz, 10)
self.assertIs(alpha.unit, nanometer**-1)
force.setPMEParameters(1.0/angstrom, 20, 20, 20)
alpha, nx, ny, nz = force.getPMEParameters()
self.assertAlmostEqualUnit(alpha, 1.0/angstrom)
self.assertEqual(nx, 20)
self.assertEqual(ny, 20)
self.assertEqual(nz, 20)
self.assertIs(alpha.unit, nanometer**-1)
self.assertTrue(force.usesPeriodicBoundaryConditions())
force.setExceptionsUsePeriodicBoundaryConditions(False)
self.assertFalse(force.getExceptionsUsePeriodicBoundaryConditions())
force.setExceptionsUsePeriodicBoundaryConditions(True)
self.assertTrue(force.getExceptionsUsePeriodicBoundaryConditions())
force.setConstantPotentialMethod(ConstantPotentialForce.CG)
self.assertEqual(force.getConstantPotentialMethod(), ConstantPotentialForce.CG)
force.setConstantPotentialMethod(ConstantPotentialForce.Matrix)
self.assertEqual(force.getConstantPotentialMethod(), ConstantPotentialForce.Matrix)
force.setUsePreconditioner(False)
self.assertFalse(force.getUsePreconditioner())
force.setUsePreconditioner(True)
self.assertTrue(force.getUsePreconditioner())
force.setUseChargeConstraint(False)
self.assertFalse(force.getUseChargeConstraint())
force.setUseChargeConstraint(True)
self.assertTrue(force.getUseChargeConstraint())
force.setChargeConstraintTarget(1.0)
target = force.getChargeConstraintTarget()
self.assertAlmostEqualUnit(target, 1.0*elementary_charge)
self.assertIs(target.unit, elementary_charge)
force.setChargeConstraintTarget(1.0*coulomb)
target = force.getChargeConstraintTarget()
self.assertAlmostEqualUnit(target, 1.0*coulomb)
self.assertIs(target.unit, elementary_charge)
force.setCGErrorTolerance(1.0)
target = force.getCGErrorTolerance()
self.assertAlmostEqualUnit(target, 1.0*kilojoule_per_mole/elementary_charge)
self.assertIs(target.unit, kilojoule_per_mole/elementary_charge)
force.setCGErrorTolerance(1.0*kilocalorie_per_mole/coulomb)
target = force.getCGErrorTolerance()
self.assertAlmostEqualUnit(target, 1.0*kilocalorie_per_mole/coulomb)
self.assertIs(target.unit, kilojoule_per_mole/elementary_charge)
force.setExternalField(Vec3(1.0, 2.0, 3.0))
field_x, field_y, field_z = force.getExternalField()
self.assertAlmostEqualUnit(field_x, 1.0*kilojoule_per_mole/(nanometer*elementary_charge))
self.assertAlmostEqualUnit(field_y, 2.0*kilojoule_per_mole/(nanometer*elementary_charge))
self.assertAlmostEqualUnit(field_z, 3.0*kilojoule_per_mole/(nanometer*elementary_charge))
self.assertIs(field_x.unit, kilojoule_per_mole/(nanometer*elementary_charge))
self.assertIs(field_y.unit, kilojoule_per_mole/(nanometer*elementary_charge))
self.assertIs(field_z.unit, kilojoule_per_mole/(nanometer*elementary_charge))
force.setExternalField(Vec3(1.0, 2.0, 3.0)*kilocalorie_per_mole/(angstrom*coulomb))
field_x, field_y, field_z = force.getExternalField()
self.assertAlmostEqualUnit(field_x, 1.0*kilocalorie_per_mole/(angstrom*coulomb))
self.assertAlmostEqualUnit(field_y, 2.0*kilocalorie_per_mole/(angstrom*coulomb))
self.assertAlmostEqualUnit(field_z, 3.0*kilocalorie_per_mole/(angstrom*coulomb))
self.assertIs(field_x.unit, kilojoule_per_mole/(nanometer*elementary_charge))
self.assertIs(field_y.unit, kilojoule_per_mole/(nanometer*elementary_charge))
self.assertIs(field_z.unit, kilojoule_per_mole/(nanometer*elementary_charge))
system = System()
system.addParticle(1.0)
system.addParticle(1.0)
force = ConstantPotentialForce()
force.addParticle(2.0*elementary_charge)
force.addParticle(0.0*elementary_charge)
force.addElectrode((1,), 0.0, 1.0, 0.0)
force.setUseChargeConstraint(True)
system.addForce(force)
context = Context(system, VerletIntegrator(0.001))
context.setPositions([[0.0, 0.0, 0.0], [1.0, 1.0, 1.0]])
charge1, charge2 = force.getCharges(context)
self.assertAlmostEqualUnit(charge1, 2.0*elementary_charge)
self.assertAlmostEqualUnit(charge2, -2.0*elementary_charge)
self.assertIs(charge1.unit, elementary_charge)
self.assertIs(charge2.unit, elementary_charge)
def testLCPOForce(self):
""" Tests the LCPOForce API features """
force = LCPOForce()
force.setSurfaceTension(10.0*kilocalorie_per_mole/bohr**2)
surfaceTension = force.getSurfaceTension()
self.assertAlmostEqualUnit(surfaceTension, 10.0*kilocalorie_per_mole/bohr**2)
self.assertIs(surfaceTension.unit, kilojoule_per_mole/nanometer**2)
force.setSurfaceTension(10.0)
surfaceTension = force.getSurfaceTension()
self.assertAlmostEqualUnit(surfaceTension, 10.0*kilojoule_per_mole/nanometer**2)
force.addParticle(1.0, 2.0, 3.0, 4.0, 5.0)
force.addParticle(6.0*bohr, 7.0, 8.0, 9.0, 10.0/bohr**2)
self.assertEqual(force.getNumParticles(), 2)
radius, p1, p2, p3, p4 = force.getParticleParameters(0)
self.assertAlmostEqualUnit(radius, 1.0*nanometer)
self.assertEqual(p1, 2.0)
self.assertEqual(p2, 3.0)
self.assertEqual(p3, 4.0)
self.assertAlmostEqualUnit(p4, 5.0/nanometer**2)
self.assertIs(radius.unit, nanometer)
self.assertIs(p4.unit, nanometer**-2)
radius, p1, p2, p3, p4 = force.getParticleParameters(1)
self.assertAlmostEqualUnit(radius, 6.0*bohr)
self.assertEqual(p1, 7.0)
self.assertEqual(p2, 8.0)
self.assertEqual(p3, 9.0)
self.assertAlmostEqualUnit(p4, 10.0/bohr**2)
self.assertIs(radius.unit, nanometer)
self.assertIs(p4.unit, nanometer**-2)
force.setParticleParameters(0, 11.0, 12.0, 13.0, 14.0, 15.0)
force.setParticleParameters(1, 16.0*bohr, 17.0, 18.0, 19.0, 20.0/bohr**2)
radius, p1, p2, p3, p4 = force.getParticleParameters(0)
self.assertAlmostEqualUnit(radius, 11.0*nanometer)
self.assertEqual(p1, 12.0)
self.assertEqual(p2, 13.0)
self.assertEqual(p3, 14.0)
self.assertAlmostEqualUnit(p4, 15.0/nanometer**2)
radius, p1, p2, p3, p4 = force.getParticleParameters(1)
self.assertAlmostEqualUnit(radius, 16.0*bohr)
self.assertEqual(p1, 17.0)
self.assertEqual(p2, 18.0)
self.assertEqual(p3, 19.0)
self.assertAlmostEqualUnit(p4, 20.0/bohr**2)
self.assertFalse(force.usesPeriodicBoundaryConditions())
force.setUsesPeriodicBoundaryConditions(True)
self.assertTrue(force.usesPeriodicBoundaryConditions())
force.setUsesPeriodicBoundaryConditions(False)
self.assertFalse(force.usesPeriodicBoundaryConditions())
def testCmapForce(self):
""" Tests the CMAPTorsionForce API features """
map1 = [random.random() for i in range(24*24)]
map2 = [random.random() for i in range(12*12)] * kilocalories_per_mole
force = CMAPTorsionForce()
force.addMap(24, map1)
force.addMap(12, map2)
force.addTorsion(0, 0, 1, 2, 3, 1, 2, 3, 4)
force.addTorsion(1, 5, 6, 7, 8, 6, 7, 8, 9)
force.addTorsion(0, 10, 11, 12, 13, 11, 12, 13, 14)
force.addTorsion(1, 15, 16, 17, 18, 16, 17, 18, 19)
self.assertEqual(force.getNumTorsions(), 4)
self.assertEqual(force.getNumMaps(), 2)
self.assertEqual(force.getMapParameters(0)[0], 24)
self.assertEqual(force.getMapParameters(1)[0], 12)
self.assertIs(force.getMapParameters(0)[1].unit, kilojoules_per_mole)
self.assertIs(force.getMapParameters(1)[1].unit, kilojoules_per_mole)
for x, y in zip(force.getMapParameters(0)[1], map1):
self.assertAlmostEqual(x.value_in_unit(kilojoules_per_mole), y)
for x, y in zip(force.getMapParameters(1)[1], map2):
self.assertAlmostEqualUnit(x, y)
def testCustomBondForce(self):
""" Tests the CustomBondForce API features """
force = CustomBondForce('1/2*k*(r-r0)^2')
force.addPerBondParameter('r0')
force.addBond(0, 1, [0.1])
force.addBond(1, 2, [1.0*angstroms])
self.assertEqual(force.getNumBonds(), 2)
i, j, req = force.getBondParameters(0)
self.assertEqual(i, 0)
self.assertEqual(j, 1)
self.assertEqual(req[0], 0.1)
i, j, req = force.getBondParameters(1)
self.assertEqual(i, 1)
self.assertEqual(j, 2)
self.assertEqual(req[0], 0.1)
def testCustomAngleForce(self):
""" Tests the CustomAngleForce API features """
force = CustomAngleForce('1/2*k*(theta-theta0)^2')
force.addPerAngleParameter('theta0')
force.addAngle(0, 1, 2, [math.pi / 2])
force.addAngle(3, 4, 5, [90*degrees])
self.assertEqual(force.getNumAngles(), 2)
i, j, k, theta = force.getAngleParameters(0)
self.assertEqual(i, 0)
self.assertEqual(j, 1)
self.assertEqual(k, 2)
self.assertEqual(theta[0], math.pi / 2)
i, j, k, theta = force.getAngleParameters(1)
self.assertEqual(i, 3)
self.assertEqual(j, 4)
self.assertEqual(k, 5)
self.assertEqual(theta[0], math.pi / 2)
def testCustomTorsionForce(self):
""" Tests the CustomTorsionForce API features """
force = CustomTorsionForce('1/2*k*(theta-theta0)^2')
force.addTorsion(0, 1, 2, 3, [math.pi])
force.addTorsion(4, 5, 6, 7, [180*degrees])
self.assertEqual(force.getNumTorsions(), 2)
i, j, k, l, theta = force.getTorsionParameters(0)
self.assertEqual(i, 0)
self.assertEqual(j, 1)
self.assertEqual(k, 2)
self.assertEqual(l, 3)
self.assertEqual(theta[0], math.pi)
i, j, k, l, theta = force.getTorsionParameters(1)
self.assertEqual(i, 4)
self.assertEqual(j, 5)
self.assertEqual(k, 6)
self.assertEqual(l, 7)
self.assertEqual(theta[0], math.pi)
def testCustomCompoundBondForce(self):
""" Tests the CustomCompoundBondForce API features """
force = CustomCompoundBondForce(4, 'kb*distance(p1, p2)*distance(p2, p3)*distance(p3, p4)+'
'ka*angle(p1, p2, p3)*angle(p2, p3, p4)+'
'kd*dihedral(p1, p2, p3, p4)')
force.addPerBondParameter('kb')
force.addPerBondParameter('ka')
force.addPerBondParameter('kd')
force.addBond([0, 1, 2, 3], [1.0, 2.0, 3.0])
force.addBond([4, 5, 6, 7], [1.0*kilocalories_per_mole/angstroms,
2.0*kilocalories_per_mole/radians,
3.0*kilocalories_per_mole]
)
self.assertEqual(force.getNumBonds(), 2)
(i, j, k, l), (kb, ka, kd) = force.getBondParameters(0)
self.assertEqual(i, 0)
self.assertEqual(j, 1)
self.assertEqual(k, 2)
self.assertEqual(l, 3)
self.assertEqual(kb, 1)
self.assertEqual(ka, 2)
self.assertEqual(kd, 3)
(i, j, k, l), (kb, ka, kd) = force.getBondParameters(1)
self.assertEqual(i, 4)
self.assertEqual(j, 5)
self.assertEqual(k, 6)
self.assertEqual(l, 7)
self.assertEqual(kb, 1*10*4.184)
self.assertEqual(ka, 2*4.184)
self.assertEqual(kd, 3*4.184)
def testCustomExternalForce(self):
""" Tests the CustomExternalForce API features """
force = CustomExternalForce('1/2*k*k2*((x-x0)^2+(y-y0)^2+(z-z0)^2)')
force.addGlobalParameter('k', 10*kilocalories_per_mole/angstroms**2)
force.addGlobalParameter('k2', 20)
force.addPerParticleParameter('x0')
force.addPerParticleParameter('y0')
force.addPerParticleParameter('z0')
force.addParticle(0, [1.0, 2.0, 3.0])
force.addParticle(1, [1.0*angstroms, 2.0*angstroms, 3.0*angstroms])
self.assertEqual(force.getNumParticles(), 2)
self.assertEqual(force.getNumGlobalParameters(), 2)
self.assertEqual(force.getGlobalParameterName(0), 'k')
self.assertEqual(force.getGlobalParameterName(1), 'k2')
self.assertEqual(force.getGlobalParameterDefaultValue(0), 1000*4.184)
self.assertEqual(force.getGlobalParameterDefaultValue(1), 20)
i, (x0, y0, z0) = force.getParticleParameters(0)
self.assertEqual(i, 0)
self.assertEqual(x0, 1)
self.assertEqual(y0, 2)
self.assertEqual(z0, 3)
i, (x0, y0, z0) = force.getParticleParameters(1)
self.assertEqual(i, 1)
self.assertAlmostEqual(x0, 1/10)
self.assertAlmostEqual(y0, 2/10)
self.assertAlmostEqual(z0, 3/10)
def testCustomGBForce(self):
""" Tests the CustomGBForce API features """
force = CustomGBForce()
force.addPerParticleParameter('q')
force.addPerParticleParameter('radius')
force.addPerParticleParameter('scale')
force.addGlobalParameter('extDiel', 78.5)
force.addGlobalParameter('intDiel', 1.0)
force.addComputedValue('I',
'step(r+sr2-or1)*0.5*(1/L-1/U+0.25*(1/U^2-1/L^2)*'
'(r-sr2^2/r)+0.5*log(L/U)/r+C);'
'U=r+sr2; C=2*(1/or1-1/L)*step(sr2-r-or1);'
'L=max(or1, D); D=abs(r-sr2); sr2=scale2*or2;'
'or1=radius1-0.009; or2=radius2-0.009',
CustomGBForce.ParticlePairNoExclusions
)
force.addComputedValue('B',
'1/(1/or-tanh(1*psi-0.8*psi^2+4.85*psi^3)/radius);'
'psi=I*or; or=radius-0.009',
CustomGBForce.SingleParticle
)
force.addEnergyTerm('-138.935456*(1/intDiel-1/extDiel)*q1*q2/f;'
'f=sqrt(r^2+B1*B2*exp(-r^2/(4*B1*B2)))',
CustomGBForce.ParticlePair
)
force.setNonbondedMethod(CustomGBForce.CutoffPeriodic)
force.addParticle([1.0, 0.1, 0.5])
force.addParticle([-1.0*coulombs, 1.0*angstroms, 0.5])
self.assertEqual(force.getNumParticles(), 2)
self.assertEqual(force.getNumComputedValues(), 2)
self.assertEqual(force.getNumEnergyTerms(), 1)
self.assertEqual(force.getNumPerParticleParameters(), 3)
self.assertEqual(force.getPerParticleParameterName(0), 'q')
self.assertEqual(force.getPerParticleParameterName(1), 'radius')
self.assertEqual(force.getPerParticleParameterName(2), 'scale')
self.assertTrue(force.usesPeriodicBoundaryConditions())
force.setNonbondedMethod(CustomGBForce.NoCutoff)
self.assertFalse(force.usesPeriodicBoundaryConditions())
force.setNonbondedMethod(CustomGBForce.CutoffNonPeriodic)
self.assertFalse(force.usesPeriodicBoundaryConditions())
q, rad, scale = force.getParticleParameters(0)
self.assertEqual(q, 1.0)
self.assertEqual(rad, 0.1)
self.assertEqual(scale, 0.5)
q, rad, scale = force.getParticleParameters(1)
self.assertEqual(q, -6.241509074460763e18) # very electronegative
self.assertEqual(rad, 0.1)
self.assertEqual(scale, 0.5)
force.setCutoffDistance(12*angstroms)
self.assertAlmostEqualUnit(force.getCutoffDistance(), 1.2*nanometers)
force.setCutoffDistance(1)
self.assertEqual(force.getCutoffDistance(), 1*nanometers)
def testCustomHBondForce(self):
""" Tests the CustomHbondForce API features """
force = CustomHbondForce('kd*(distance(a1,d1)-r0)^2 + '
'ka*(angle(a1,d1,d2)-theta0)^2')
force.addPerAcceptorParameter('r0')
force.addPerAcceptorParameter('ka')
force.addPerDonorParameter('theta0')
force.addPerDonorParameter('kd')
force.setCutoffDistance(10*angstroms)
self.assertEqual(force.getNumPerAcceptorParameters(), 2)
self.assertEqual(force.getNumPerDonorParameters(), 2)
self.assertEqual(force.getCutoffDistance(), 1*nanometers)
self.assertIs(force.getCutoffDistance().unit, nanometer)
force.addAcceptor(0, 1, 2, [0.2, 10.0])
force.addAcceptor(3, -1, -1, [4*angstroms,
20.0*kilocalories_per_mole/angstroms**2])
force.addDonor(4, 5, 6, [math.pi, 30])
force.addDonor(7, 8, -1, [180*degrees,
40*kilocalories_per_mole/radians**2])
self.assertEqual(force.getNumAcceptors(), 2)
self.assertEqual(force.getNumDonors(), 2)
i, j, k, (r0, ka) = force.getAcceptorParameters(0)
self.assertEqual(i, 0)
self.assertEqual(j, 1)
self.assertEqual(k, 2)
self.assertEqual(r0, 0.2)
self.assertEqual(ka, 10)
i, j, k, (r0, ka) = force.getAcceptorParameters(1)
self.assertEqual(i, 3)
self.assertEqual(j, -1)
self.assertEqual(k, -1)
self.assertEqual(r0, 0.4)
self.assertEqual(ka, 20*4.184*100)
i, j, k, (theta0, kd) = force.getDonorParameters(0)
self.assertEqual(i, 4)
self.assertEqual(j, 5)
self.assertEqual(k, 6)
self.assertEqual(theta0, math.pi)
self.assertEqual(kd, 30)
i, j, k, (theta0, kd) = force.getDonorParameters(1)
self.assertEqual(i, 7)
self.assertEqual(j, 8)
self.assertEqual(k, -1)
self.assertEqual(theta0, math.pi)
self.assertEqual(kd, 40*4.184)
self.assertFalse(force.usesPeriodicBoundaryConditions())
force.setNonbondedMethod(CustomHbondForce.NoCutoff)
self.assertFalse(force.usesPeriodicBoundaryConditions())
force.setNonbondedMethod(CustomHbondForce.CutoffNonPeriodic)
self.assertFalse(force.usesPeriodicBoundaryConditions())
force.setNonbondedMethod(CustomHbondForce.CutoffPeriodic)
self.assertTrue(force.usesPeriodicBoundaryConditions())
def testCustomNonbondedForce(self):
""" Tests the CustomNonbondedForce API features """
force = CustomNonbondedForce('4*diel*q1*q2/r+sqrt(p1*p2)/r^2-sqrt(m1*m2)/r^3')
force.addGlobalParameter('diel', 1.0)
force.addPerParticleParameter('q')
force.addPerParticleParameter('p')
force.addPerParticleParameter('m')
force.addParticle([1, 2, 3])
force.addParticle([1*coulombs, 2*kilocalories_per_mole*angstroms**2, 3*kilocalories_per_mole*angstroms**3])
force.addTabulatedFunction('f', Continuous1DFunction([1,2,3,4,5], 0.0, 2.0))
self.assertEqual(force.getNumParticles(), 2)
charge, sigma, epsilon = force.getParticleParameters(0)
self.assertEqual(charge, 1)
self.assertEqual(sigma, 2)
self.assertEqual(epsilon, 3)
charge, sigma, epsilon = force.getParticleParameters(1)
self.assertEqual(charge, 6.241509074460763e18) # very electronegative
self.assertEqual(sigma, 2*4.184/100)
self.assertAlmostEqual(epsilon, 3*4.184/1000)
force.setCutoffDistance(10*angstroms)
self.assertEqual(force.getCutoffDistance(), 1*nanometers)
force.setCutoffDistance(1)
self.assertEqual(force.getCutoffDistance(), 1*nanometers)
force.setSwitchingDistance(8*angstroms)
self.assertEqual(force.getSwitchingDistance(), 0.8*nanometers)
self.assertIs(force.getSwitchingDistance().unit, nanometer)
self.assertFalse(force.usesPeriodicBoundaryConditions())
force.setNonbondedMethod(CustomNonbondedForce.NoCutoff)
self.assertFalse(force.usesPeriodicBoundaryConditions())
force.setNonbondedMethod(CustomNonbondedForce.CutoffNonPeriodic)
self.assertFalse(force.usesPeriodicBoundaryConditions())
force.setNonbondedMethod(CustomNonbondedForce.CutoffPeriodic)
self.assertTrue(force.usesPeriodicBoundaryConditions())
self.assertIs(type(force.getTabulatedFunction(0)), Continuous1DFunction)
def testCustomManyParticleForce(self):
""" Tests the CustomManyParticleForce API features """
force = CustomManyParticleForce(3,
"C*(1+3*cos(theta1)*cos(theta2)*cos(theta3))/(r12*r13*r23)^3;"
"theta1=k1*angle(p1,p2,p3); theta2=k2*angle(p2,p3,p1); theta3=k3*angle(p3,p1,p2);"
"r12=distance(p1,p2); r13=distance(p1,p3); r23=distance(p2,p3)")
force.setPermutationMode(CustomManyParticleForce.SinglePermutation)
force.setTypeFilter(0, [0])
force.setTypeFilter(1, [1])
force.setTypeFilter(2, [2])
force.addGlobalParameter('C', 1.0*kilocalories_per_mole)
force.addPerParticleParameter('k')
self.assertEqual(force.getNumGlobalParameters(), 1)
self.assertEqual(force.getGlobalParameterName(0), 'C')
self.assertEqual(force.getGlobalParameterDefaultValue(0), 4.184)
self.assertEqual(force.getNumPerParticleParameters(), 1)
force.addParticle([10], 0)
force.addParticle([20], 1)
force.addParticle([30*kilocalories_per_mole], 2)
self.assertEqual(force.getNumParticles(), 3)
self.assertEqual(force.getParticleParameters(0)[0][0], 10)
self.assertEqual(force.getParticleParameters(1)[0][0], 20)
self.assertEqual(force.getParticleParameters(2)[0][0], 30*4.184)
def testATMForce(self):
"""Tests the ATMForce API features"""
force = ATMForce(0.1, 0.2, 0.3, 0.4, 0.5, 0.7, 0.6, 0.8, -1.0);
#particle 0: fixed displacements,
force.addParticle(Vec3(1, 2, 3), Vec3(4, 5, 6))
#particle 1: fixed displacements using a Transformation object
p = force.addParticle()
force.setParticleTransformation(p, FixedDisplacement(Vec3(7, 8, 9), Vec3(10, 11, 12)))
#particle 2: particle distance displacement
p = force.addParticle()
force.setParticleTransformation(p, ParticleOffsetDisplacement(1, 0))
#particle 3: stationary particle
force.addParticle()
self.assertEqual(0.1, force.getGlobalParameterDefaultValue(0))
self.assertEqual(0.2, force.getGlobalParameterDefaultValue(1))
self.assertEqual(0.3, force.getGlobalParameterDefaultValue(2))
self.assertEqual(0.4, force.getGlobalParameterDefaultValue(3))
self.assertEqual(0.5, force.getGlobalParameterDefaultValue(4))
self.assertEqual(0.7, force.getGlobalParameterDefaultValue(5))
self.assertEqual(0.6, force.getGlobalParameterDefaultValue(6))
self.assertEqual(0.8, force.getGlobalParameterDefaultValue(7))
self.assertEqual(-1.0, force.getGlobalParameterDefaultValue(8))
d1, d0 = force.getParticleParameters(0)
self.assertEqual(Vec3(1, 2, 3)*nanometers, d1)
self.assertEqual(Vec3(4, 5, 6)*nanometers, d0)
fixed_displacement_transformation = force.getParticleTransformation(1)
d1 = fixed_displacement_transformation.getFixedDisplacement1()
d0 = fixed_displacement_transformation.getFixedDisplacement0()
self.assertEqual(Vec3(7, 8, 9)*nanometers, d1)
self.assertEqual(Vec3(10, 11, 12)*nanometers, d0)
vectordistance_displacement_transformation = force.getParticleTransformation(2)
j1 = vectordistance_displacement_transformation.getDestinationParticle1()
i1 = vectordistance_displacement_transformation.getOriginParticle1()
j0 = vectordistance_displacement_transformation.getDestinationParticle0()
i0 = vectordistance_displacement_transformation.getOriginParticle0()
self.assertEqual( 1, j1)
self.assertEqual( 0, i1)
self.assertEqual(-1, j0)
self.assertEqual(-1, i0)
transformation = force.getParticleTransformation(3)
d1, d0 = force.getParticleParameters(3)
self.assertEqual(Vec3(0, 0, 0)*nanometers, d1)
self.assertEqual(Vec3(0, 0, 0)*nanometers, d0)
def testDrudeForce(self):
""" Tests the DrudeForce API features """
force = DrudeForce()
self.assertFalse(force.usesPeriodicBoundaryConditions())
force.addParticle(0, 1, -1, -1, -1, 1, 1, 0, 0)
force.addParticle(1, 2, 3, -1, -1, 1*elementary_charge, 1*angstrom**3, 0.5, 0)
force.addParticle(2, 3, 4, 5, 6, 1*elementary_charge, 10*angstrom**3, 0.5, 0.5)
force.addScreenedPair(0, 1, 0.5)
force.addScreenedPair(1, 2, 0.25)
force.addScreenedPair(0, 2, 0.125)
self.assertEqual(force.getNumParticles(), 3)
self.assertEqual(force.getNumScreenedPairs(), 3)
i, j, k, l, m, q, a, an12, an34 = force.getParticleParameters(0)
self.assertEqual(i, 0)
self.assertEqual(j, 1)
self.assertEqual(k, -1)
self.assertEqual(l, -1)
self.assertEqual(m, -1)
self.assertEqual(q, 1*elementary_charge)
self.assertEqual(a, 1*nanometer**3)
self.assertEqual(an12, 0)
self.assertEqual(an34, 0)
i, j, k, l, m, q, a, an12, an34 = force.getParticleParameters(1)
self.assertEqual(i, 1)
self.assertEqual(j, 2)
self.assertEqual(k, 3)
self.assertEqual(l, -1)
self.assertEqual(m, -1)
self.assertEqual(q, 1*elementary_charge)
self.assertAlmostEqualUnit(a, 1*angstrom**3)
self.assertEqual(an12, 0.5)
self.assertEqual(an34, 0)
i, j, k, l, m, q, a, an12, an34 = force.getParticleParameters(2)
self.assertEqual(i, 2)
self.assertEqual(j, 3)
self.assertEqual(k, 4)
self.assertEqual(l, 5)
self.assertEqual(m, 6)
self.assertEqual(q, 1*elementary_charge)
self.assertAlmostEqualUnit(a, 10*angstrom**3)
self.assertEqual(an12, 0.5)
self.assertEqual(an34, 0.5)
i, j, thole = force.getScreenedPairParameters(0)
self.assertEqual(i, 0)
self.assertEqual(j, 1)
self.assertEqual(thole, 0.5)
i, j, thole = force.getScreenedPairParameters(1)
self.assertEqual(i, 1)
self.assertEqual(j, 2)
self.assertEqual(thole, 0.25)
i, j, thole = force.getScreenedPairParameters(2)
self.assertEqual(i, 0)
self.assertEqual(j, 2)
self.assertEqual(thole, 0.125)
def testGeneralizedKirkwood(self):
""" Tests the AmoebaGeneralizedKirkwoodForce API features """
force = AmoebaGeneralizedKirkwoodForce()
self.assertEqual(force.getProbeRadius(), 0.14*nanometer) # default
force.setProbeRadius(0.16)
self.assertEqual(force.getProbeRadius(), 0.16*nanometer)
self.assertIs(force.getProbeRadius().unit, nanometer)
force.setProbeRadius(1.4*angstrom)
self.assertEqual(force.getProbeRadius(), 1.4*angstrom)
self.assertIs(force.getProbeRadius().unit, nanometer)
self.assertEqual(force.getSoluteDielectric(), 1.0) # default
force.setSoluteDielectric(2.0)
self.assertEqual(force.getSoluteDielectric(), 2.0)
self.assertEqual(force.getSolventDielectric(), 78.3) # default
force.setSolventDielectric(80)
self.assertEqual(force.getSolventDielectric(), 80)
self.assertAlmostEqualUnit(force.getSurfaceAreaFactor(),
-170.35173066268223*kilojoule_per_mole/nanometer**2) # default
force.setSurfaceAreaFactor(-1.0*kilocalorie_per_mole/angstrom**2)
self.assertAlmostEqualUnit(force.getSurfaceAreaFactor(),
-1.0*kilocalorie_per_mole/angstrom**2) # default
force.addParticle(1.0*coulomb, 1.0*angstroms, 0.5)
force.addParticle(1.0, 1.0, 0.4)
self.assertEqual(force.getNumParticles(), 2)
q, r, s, d, k = force.getParticleParameters(0)
self.assertAlmostEqualUnit(q, 1.0*coulomb)
self.assertIs(q.unit, elementary_charge)
self.assertEqual(r, 1.0*angstroms)
self.assertIs(r.unit, nanometer)
self.assertEqual(s, 0.5)
q, r, s, d, k = force.getParticleParameters(1)
self.assertAlmostEqualUnit(q, 1.0*elementary_charge)
self.assertIs(q.unit, elementary_charge)
self.assertEqual(r, 1.0*nanometer)
self.assertIs(r.unit, nanometer)
self.assertEqual(s, 0.4)
def testAmoebaTorsionTorsionForce(self):
""" Tests the AmoebaTorsionTorsionForce API features """
force = AmoebaTorsionTorsionForce()
grid1 = [[[i, j, random.random(), random.random(), random.random(), random.random()]
for j in range(-180, 180, 24)] for i in range(-180, 180, 24)]
kcal = kilocalories_per_mole
grid2 = [[[i*math.pi/180*radians, j*math.pi/180*radians, random.random()*kcal,
random.random()*kcal/radians, random.random()*kcal/radians,
random.random()*kcal/radians**2]
for j in range(-180, 0, 10)] for i in range(-180, 0, 10)]
grid3 = [[[i, j, random.random()*kcal]
for j in range(-180, 180, 24)] for i in range(-180, 180, 24)]
force.setTorsionTorsionGrid(0, grid1)
force.setTorsionTorsionGrid(1, grid2)
force.setTorsionTorsionGrid(2, grid3)
self.assertEqual(force.getNumTorsionTorsionGrids(), 3)
# Check the grids
g1 = force.getTorsionTorsionGrid(0)
for row1, row2 in zip(g1, grid1):
for column1, column2 in zip(row1, row2):
self.assertEqual(len(column1), len(column2))
for x1, x2 in zip(column1, column2):
self.assertEqual(x1, x2)
g2 = force.getTorsionTorsionGrid(1)
for row1, row2 in zip(g2, grid2):
for column1, column2 in zip(row1, row2):
self.assertEqual(column1[0], column2[0].value_in_unit(degree))
self.assertEqual(column1[1], column2[1].value_in_unit(degree))
self.assertEqual(column1[2], column2[2].value_in_unit(kilojoules_per_mole))
self.assertEqual(column1[3], column2[3].value_in_unit(kilojoules_per_mole/radian))
self.assertEqual(column1[4], column2[4].value_in_unit(kilojoules_per_mole/radian))
self.assertEqual(column1[5], column2[5].value_in_unit(kilojoules_per_mole/radian**2))
g3 = force.getTorsionTorsionGrid(2)
for row1, row2 in zip(g3, grid3):
for column1, column2 in zip(row1, row2):
self.assertEqual(len(column1), 6)
self.assertEqual(len(column2), 3)
self.assertEqual(column1[0], column2[0])
self.assertEqual(column1[1], column2[1])
self.assertEqual(column1[2], column2[2].value_in_unit(kilojoules_per_mole))
force.addTorsionTorsion(0, 1, 2, 3, 4, 5, 0)
force.addTorsionTorsion(1, 2, 3, 4, 5, 6, 1)
force.addTorsionTorsion(2, 3, 4, 5, 6, 7, 2)
self.assertEqual(force.getNumTorsionTorsions(), 3)
i, j, k, l, m, ch, g = force.getTorsionTorsionParameters(0)
self.assertEqual(i, 0)
self.assertEqual(j, 1)
self.assertEqual(k, 2)
self.assertEqual(l, 3)
self.assertEqual(m, 4)
self.assertEqual(ch, 5)
self.assertEqual(g, 0)
i, j, k, l, m, ch, g = force.getTorsionTorsionParameters(1)
self.assertEqual(i, 1)
self.assertEqual(j, 2)
self.assertEqual(k, 3)
self.assertEqual(l, 4)
self.assertEqual(m, 5)
self.assertEqual(ch, 6)
self.assertEqual(g, 1)
def testAmoebaVdwForce(self):
""" Tests the AmoebaVdwForce API features """
force = AmoebaVdwForce()
self.assertEqual(force.getSigmaCombiningRule(), 'CUBIC-MEAN')
force.setSigmaCombiningRule('ARITHMETIC')
self.assertEqual(force.getSigmaCombiningRule(), 'ARITHMETIC')
force.setSigmaCombiningRule('GEOMETRIC')
self.assertEqual(force.getSigmaCombiningRule(), 'GEOMETRIC')
self.assertEqual(force.getEpsilonCombiningRule(), 'HHG')
force.setEpsilonCombiningRule('HARMONIC')
self.assertEqual(force.getEpsilonCombiningRule(), 'HARMONIC')
force.setEpsilonCombiningRule('GEOMETRIC')
self.assertEqual(force.getEpsilonCombiningRule(), 'GEOMETRIC')
force.setEpsilonCombiningRule('ARITHMETIC')
self.assertEqual(force.getEpsilonCombiningRule(), 'ARITHMETIC')
self.assertTrue(force.getUseDispersionCorrection())
force.setUseDispersionCorrection(False)
self.assertFalse(force.getUseDispersionCorrection())
self.assertIs(force.getNonbondedMethod(), AmoebaVdwForce.NoCutoff)
self.assertFalse(force.usesPeriodicBoundaryConditions())
force.setNonbondedMethod(AmoebaVdwForce.CutoffPeriodic)
self.assertTrue(force.usesPeriodicBoundaryConditions())
self.assertIs(force.getNonbondedMethod(), AmoebaVdwForce.CutoffPeriodic)
force.setCutoff(10.0*angstroms)
self.assertEqual(force.getCutoff(), 10.0*angstroms)
self.assertIs(force.getCutoff().unit, nanometers)
force.addParticle(0, 0.1, 1.0, 1.0)
force.addParticle(1, 1.0*angstroms, 1.0*kilocalories_per_mole, 0.5)
force.addParticle(1, 0.8*angstroms, 2.0*kilocalories_per_mole, 0.25)
self.assertEqual(force.getNumParticles(), 3)
p, sig, eps, reduction, alchemical, type, scale = force.getParticleParameters(0)
self.assertEqual(p, 0)
self.assertEqual(sig, 0.1*nanometers)
self.assertIs(sig.unit, nanometers)
self.assertEqual(eps, 1.0*kilojoules_per_mole)
self.assertIs(eps.unit, kilojoules_per_mole)
self.assertEqual(reduction, 1.0)
self.assertEqual(type, -1)
self.assertEqual(scale, 1.0)
p, sig, eps, reduction, alchemical, type, scale = force.getParticleParameters(1)
self.assertEqual(p, 1)
self.assertEqual(sig, 1.0*angstroms)
self.assertIs(sig.unit, nanometers)
self.assertEqual(eps, 1.0*kilocalories_per_mole)
self.assertIs(eps.unit, kilojoules_per_mole)
self.assertEqual(reduction, 0.5)
self.assertEqual(type, -1)
self.assertEqual(scale, 1.0)
p, sig, eps, reduction, alchemical, type, scale = force.getParticleParameters(2)
self.assertEqual(p, 1)
self.assertAlmostEqualUnit(sig, 0.8*angstroms)
self.assertIs(sig.unit, nanometers)
self.assertEqual(eps, 2.0*kilocalories_per_mole)
self.assertIs(eps.unit, kilojoules_per_mole)
self.assertEqual(reduction, 0.25)
self.assertEqual(type, -1)
self.assertEqual(scale, 1.0)
def testAmoebaWcaDispersionForce(self):
""" Tests the AmoebaWcaDispersionForce API features """
force = AmoebaWcaDispersionForce()
self.assertEqual(force.getDispoff(), 0.1056*nanometer)
self.assertEqual(force.getAwater(), 33.428*nanometer**-3)
self.assertEqual(force.getEpsh(), 0.056484*kilojoule_per_mole)
self.assertEqual(force.getEpso(), 0.46024000000000004*kilojoule_per_mole)
self.assertEqual(force.getRminh(), 0.13275*nanometer)
self.assertEqual(force.getRmino(), 0.17025*nanometer)
self.assertEqual(force.getShctd(), 0.82)
self.assertEqual(force.getSlevy(), 1.0)
force.setDispoff(3*angstroms)
self.assertAlmostEqualUnit(force.getDispoff(), 3*angstroms)
self.assertIs(force.getDispoff().unit, nanometer)
force.setAwater(3*angstroms**-3)
self.assertAlmostEqualUnit(force.getAwater(), 3*angstroms**-3)
self.assertEqual(1*force.getAwater().unit, 1*nanometer**-3)
force.setEpsh(1*kilocalorie_per_mole)
self.assertEqual(force.getEpsh(), 1*kilocalorie_per_mole)
self.assertIs(force.getEpsh().unit, kilojoule_per_mole)
force.setEpso(1*kilocalorie_per_mole)
self.assertEqual(force.getEpso(), 1*kilocalorie_per_mole)
self.assertIs(force.getEpso().unit, kilojoule_per_mole)
force.setRminh(20*angstroms)
self.assertEqual(force.getRminh(), 20*angstroms)
self.assertIs(force.getRminh().unit, nanometer)
force.setRmino(30*angstroms)
self.assertEqual(force.getRmino(), 30*angstroms)
self.assertIs(force.getRmino().unit, nanometer)
force.setShctd(1)
self.assertEqual(force.getShctd(), 1)
force.setSlevy(2)
self.assertEqual(force.getSlevy(), 2)
force.addParticle(0.5, 1)
force.addParticle(3*angstroms, 1*kilocalorie_per_mole)
self.assertEqual(force.getNumParticles(), 2)
sig, eps = force.getParticleParameters(0)
self.assertEqual(sig, 0.5*nanometer)
self.assertIs(sig.unit, nanometer)
self.assertEqual(eps, 1*kilojoule_per_mole)
self.assertIs(eps.unit, kilojoule_per_mole)
sig, eps = force.getParticleParameters(1)
self.assertAlmostEqualUnit(sig, 3*angstrom)
self.assertIs(sig.unit, nanometer)
self.assertEqual(eps, 1*kilocalorie_per_mole)
self.assertIs(eps.unit, kilojoule_per_mole)
def testAmoebaMultipoleForce(self):
""" Tests the AmoebaMultipoleForce API features """
force = AmoebaMultipoleForce()
self.assertIs(force.getNonbondedMethod(), AmoebaMultipoleForce.NoCutoff)
self.assertFalse(force.usesPeriodicBoundaryConditions())
force.setNonbondedMethod(AmoebaMultipoleForce.PME)
self.assertIs(force.getNonbondedMethod(), AmoebaMultipoleForce.PME)
self.assertTrue(force.usesPeriodicBoundaryConditions())
self.assertEqual(force.getCutoffDistance(), 1*nanometer)
self.assertIs(force.getCutoffDistance().unit, nanometer)
force.setCutoffDistance(8*angstroms)
self.assertEqual(force.getCutoffDistance(), 8*angstrom)
self.assertIs(force.getCutoffDistance().unit, nanometer)
self.assertEqual(force.getPolarizationType(), AmoebaMultipoleForce.Mutual)
force.setPolarizationType(AmoebaMultipoleForce.Direct)
self.assertEqual(force.getPolarizationType(), AmoebaMultipoleForce.Direct)
force.setPolarizationType(AmoebaMultipoleForce.Mutual)
self.assertEqual(force.getPolarizationType(), AmoebaMultipoleForce.Mutual)
force.addMultipole(1.0, [0.5, 0, -0.5], list(range(9)),
AmoebaMultipoleForce.ZThenX, 1, 2, 3, 0.5, 0.5, 1.0)
force.addMultipole(1.0*elementary_charge, [0.5, 0, -0.5]*elementary_charge*angstrom,
list(range(9))*elementary_charge*angstrom**2,
AmoebaMultipoleForce.Bisector, 2, 3, 4, 0.5, 0.5, 1.0*angstrom**3)
self.assertEqual(force.getNumMultipoles(), 2)
q, mu, quad, ax, i, j, k, thole, damp, polarity = force.getMultipoleParameters(0)
self.assertEqual(q, 1.0*elementary_charge)
self.assertEqual(mu, (0.5, 0, -0.5)*elementary_charge*nanometer)
self.assertEqual(quad, tuple(range(9))*elementary_charge*nanometer**2)
self.assertEqual(ax, AmoebaMultipoleForce.ZThenX)
self.assertEqual(i, 1)
self.assertEqual(j, 2)
self.assertEqual(k, 3)
self.assertEqual(thole, 0.5)
self.assertEqual(damp, 0.5)
self.assertEqual(polarity, 1*nanometer**3)
q, mu, quad, ax, i, j, k, thole, damp, polarity = force.getMultipoleParameters(1)
self.assertEqual(q, 1.0*elementary_charge)
self.assertEqual(mu, (0.5, 0, -0.5)*elementary_charge*angstrom)
self.assertAlmostEqualUnitArray(quad, tuple(range(9))*elementary_charge*angstrom**2)
self.assertEqual(ax, AmoebaMultipoleForce.Bisector)
self.assertEqual(i, 2)
self.assertEqual(j, 3)
self.assertEqual(k, 4)
self.assertEqual(thole, 0.5)
self.assertEqual(damp, 0.5)
self.assertAlmostEqualUnit(polarity, 1*angstrom**3)
def testVerletIntegrator(self):
""" Tests the VerletIntegrator API features """
integrator = VerletIntegrator(1.0)
self.assertEqual(integrator.getStepSize(), 1.0*picosecond)
self.assertEqual(integrator.getConstraintTolerance(), 1e-5)
integrator.setConstraintTolerance(1e-6)
self.assertEqual(integrator.getConstraintTolerance(), 1e-6)
integrator = VerletIntegrator(1.0*femtoseconds)
self.assertEqual(integrator.getStepSize(), 1.0*femtoseconds)
def testVariableVerletIntegrator(self):
""" Tests the VariableVerletIntegrator API features """
integrator = VariableVerletIntegrator(0.1)
self.assertEqual(integrator.getErrorTolerance(), 0.1)
integrator.setErrorTolerance(0.01)
self.assertEqual(integrator.getErrorTolerance(), 0.01)
self.assertEqual(integrator.getConstraintTolerance(), 1e-5)
integrator.setConstraintTolerance(1e-6)
self.assertEqual(integrator.getConstraintTolerance(), 1e-6)
def testLangevinIntegrator(self):
""" Tests the LangevinIntegrator API features """
integrator = LangevinIntegrator(300, 0.1, 1)
self.assertEqual(integrator.getTemperature(), 300*kelvin)
self.assertEqual(integrator.getFriction(), 0.1/picosecond)
self.assertEqual(integrator.getStepSize(), 1*picosecond)
integrator = LangevinIntegrator(300*kelvin, 0.1/microsecond, 1*femtosecond)
self.assertEqual(integrator.getTemperature(), 300*kelvin)
self.assertAlmostEqualUnit(integrator.getFriction(), 0.1/microsecond)
self.assertEqual(integrator.getStepSize(), 1*femtosecond)
def testVariableLangevinIntegrator(self):
""" Tests the VariableLangevinIntegrator API features """
integrator = VariableLangevinIntegrator(300, 0.1, 0.1)
self.assertEqual(integrator.getTemperature(), 300*kelvin)
self.assertEqual(integrator.getFriction(), 0.1/picosecond)
self.assertEqual(integrator.getErrorTolerance(), 0.1)
integrator = VariableLangevinIntegrator(300*kelvin, 0.1/microsecond, 0.01)
self.assertEqual(integrator.getTemperature(), 300*kelvin)
self.assertAlmostEqualUnit(integrator.getFriction(), 0.1/microsecond)
self.assertEqual(integrator.getErrorTolerance(), 0.01)
def testBrownianIntegrator(self):
""" Tests the BrownianIntegrator API features """
integrator = BrownianIntegrator(300, 0.1, 1)
self.assertEqual(integrator.getTemperature(), 300*kelvin)
self.assertEqual(integrator.getFriction(), 0.1/picosecond)
self.assertEqual(integrator.getStepSize(), 1*picosecond)
integrator = BrownianIntegrator(300*kelvin, 0.1/microsecond, 1*femtosecond)
self.assertEqual(integrator.getTemperature(), 300*kelvin)
self.assertAlmostEqualUnit(integrator.getFriction(), 0.1/microsecond)
self.assertEqual(integrator.getStepSize(), 1*femtosecond)
def testNoseHooverIntegrator(self):
""" Tests the NoseHooverIntegrator API features """
integrator = NoseHooverIntegrator(300, 0.1, 1)
self.assertEqual(integrator.getTemperature(0), 300*kelvin)
self.assertEqual(integrator.getCollisionFrequency(), 0.1/picosecond)
self.assertEqual(integrator.getStepSize(), 1*picosecond)
integrator = NoseHooverIntegrator(300*kelvin, 0.1/microsecond, 1*femtosecond)
self.assertEqual(integrator.getTemperature(), 300*kelvin)
self.assertAlmostEqualUnit(integrator.getCollisionFrequency(), 0.1/microsecond)
self.assertEqual(integrator.getStepSize(), 1*femtosecond)
self.assertEqual(integrator.computeHeatBathEnergy(), 0.0*kilojoule_per_mole)
# Test setters
integrator.setTemperature(200*kelvin)
self.assertEqual(integrator.getTemperature(), 200*kelvin)
integrator.setCollisionFrequency(0.1/picosecond)
self.assertEqual(integrator.getCollisionFrequency(), 0.1/picosecond)
integrator.setRelativeTemperature(200*kelvin)
self.assertEqual(integrator.getRelativeTemperature(), 200*kelvin)
integrator.setRelativeCollisionFrequency(0.1/picosecond)
self.assertEqual(integrator.getRelativeCollisionFrequency(), 0.1/picosecond)
# Test bare constructor and addThermostat
integrator = NoseHooverIntegrator(1*femtosecond)
self.assertEqual(integrator.getStepSize(), 1*femtosecond)
integrator.addThermostat(300*kelvin, 0.1/microsecond, 3, 3, 3)
self.assertAlmostEqualUnit(integrator.getTemperature(), 300*kelvin)
self.assertAlmostEqualUnit(integrator.getCollisionFrequency(), 0.1/microsecond)
integrator = NoseHooverIntegrator(1*femtosecond)
integrator.addSubsystemThermostat([0], [], 300*kelvin, 0.1/microsecond, 1.0*kelvin, 1.0/microsecond, 3, 3, 3)
self.assertAlmostEqualUnit(integrator.getTemperature(), 300*kelvin)
self.assertAlmostEqualUnit(integrator.getCollisionFrequency(), 0.1/microsecond)
self.assertAlmostEqualUnit(integrator.getRelativeTemperature(), 1.0*kelvin)
self.assertAlmostEqualUnit(integrator.getRelativeCollisionFrequency(), 1.0/microsecond)
def testDrudeNoseHooverIntegrator(self):
""" Tests the DrudeNoseHooverIntegrator API features """
integrator = DrudeNoseHooverIntegrator(300, 0.1, 1.0, 1.0, 1)
self.assertEqual(integrator.getTemperature(0), 300*kelvin)
self.assertEqual(integrator.getCollisionFrequency(), 0.1/picosecond)
self.assertEqual(integrator.getRelativeTemperature(0), 1.0*kelvin)
self.assertEqual(integrator.getRelativeCollisionFrequency(), 1.0/picosecond)
self.assertEqual(integrator.getStepSize(), 1*picosecond)
integrator = DrudeNoseHooverIntegrator(300*kelvin, 0.1/microsecond, 5.0*kelvin, 0.01/microsecond, 1*femtosecond)
self.assertEqual(integrator.getTemperature(), 300*kelvin)
self.assertAlmostEqualUnit(integrator.getCollisionFrequency(), 0.1/microsecond)
self.assertEqual(integrator.getRelativeTemperature(), 5.0*kelvin)
self.assertAlmostEqualUnit(integrator.getRelativeCollisionFrequency(), 0.01/microsecond)
self.assertEqual(integrator.getStepSize(), 1*femtosecond)
def testAndersenThermostat(self):
""" Tests the AndersenThermostat API features """
force = AndersenThermostat(300*kelvin, 1/microsecond)
self.assertEqual(force.getDefaultTemperature(), 300*kelvin)
self.assertAlmostEqualUnit(force.getDefaultCollisionFrequency(), 1/microsecond)
force.setDefaultTemperature(298.15)
force.setDefaultCollisionFrequency(1)
self.assertEqual(force.getDefaultTemperature(), 298.15*kelvin)
self.assertAlmostEqualUnit(force.getDefaultCollisionFrequency(), 1/picosecond)
def testMonteCarloBarostat(self):
""" Tests the various Monte Carlo barostats API features """
def checkPressureUnits(force):
system = System()
system.addParticle(0.0)
system.addForce(force)
bonds = HarmonicBondForce()
bonds.setUsesPeriodicBoundaryConditions(True)
system.addForce(bonds)
context = Context(system, VerletIntegrator(0.001))
context.setPositions([Vec3(0, 0, 0)])
pressure = force.computeCurrentPressure(context)
self.assertTrue(is_quantity(pressure))
self.assertTrue(pressure.unit == bar)
force = MonteCarloBarostat(1.1*bar, 350*kelvin)
self.assertEqual(force.getDefaultPressure(), 1.1*bar)
self.assertEqual(force.getDefaultTemperature(), 350*kelvin)
checkPressureUnits(force)
force = MonteCarloFlexibleBarostat(1.1*bar, 350*kelvin)
self.assertEqual(force.getDefaultPressure(), 1.1*bar)
self.assertEqual(force.getDefaultTemperature(), 350*kelvin)
checkPressureUnits(force)
force = MonteCarloMembraneBarostat(1.0, 1.5, 300, MonteCarloMembraneBarostat.XYAnisotropic, MonteCarloMembraneBarostat.ZFixed, 25)
self.assertEqual(force.getDefaultPressure(), 1.0*bar)
self.assertEqual(force.getDefaultSurfaceTension(), 1.5*bar*nanometer)
self.assertEqual(force.getDefaultTemperature(), 300*kelvin)
self.assertEqual(force.getXYMode(), MonteCarloMembraneBarostat.XYAnisotropic)
self.assertEqual(force.getZMode(), MonteCarloMembraneBarostat.ZFixed)
self.assertEqual(force.getFrequency(), 25)
checkPressureUnits(force)
force = MonteCarloMembraneBarostat(1.1*bar, 2.0*bar*nanometer, 350*kelvin, MonteCarloMembraneBarostat.XYAnisotropic, MonteCarloMembraneBarostat.ZFixed, 25)
self.assertEqual(force.getDefaultPressure(), 1.1*bar)
self.assertEqual(force.getDefaultSurfaceTension(), 2.0*bar*nanometer)
self.assertEqual(force.getDefaultTemperature(), 350*kelvin)
force.setDefaultPressure(1.2*bar)
force.setDefaultSurfaceTension(2.5*bar*nanometer)
force.setDefaultTemperature(298.15)
self.assertEqual(force.getDefaultPressure(), 1.2*bar)
self.assertEqual(force.getDefaultSurfaceTension(), 2.5*bar*nanometer)
self.assertEqual(force.getDefaultTemperature(), 298.15*kelvin)
force = MonteCarloAnisotropicBarostat(Vec3(1.1, 2.2, 3.3)*bar, 350*kelvin)
self.assertEqual(force.getDefaultPressure(), Vec3(1.1, 2.2, 3.3)*bar)
self.assertEqual(force.getDefaultTemperature(), 350*kelvin)
checkPressureUnits(force)
def testDrudeSCFIntegrator(self):
""" Tests the DrudeSCFIntegrator API features """
integrator = DrudeSCFIntegrator(0.002)
self.assertEqual(integrator.getStepSize(), 0.002*picoseconds)
self.assertAlmostEqualUnit(integrator.getMinimizationErrorTolerance(),
1*kilojoule_per_mole/nanometer)
integrator.setMinimizationErrorTolerance(0.1*kilocalorie_per_mole/angstrom)
self.assertAlmostEqualUnit(integrator.getMinimizationErrorTolerance(),
0.1*kilocalorie_per_mole/angstrom)
def testDrudeLangevinIntegrator(self):
""" Tests the DrudeLangevinIntegrator API features """
integrator = DrudeLangevinIntegrator(300*kelvin, 1.0/microsecond,
100*kelvin, 1/picosecond, 0.1*femtosecond)
self.assertEqual(integrator.getTemperature(), 300*kelvin)
integrator.setTemperature(298.15)
self.assertEqual(integrator.getTemperature(), 298.15*kelvin)
self.assertAlmostEqualUnit(integrator.getFriction(), 1/microsecond)
integrator.setFriction(1)
self.assertEqual(integrator.getFriction(), 1/picosecond)
self.assertEqual(integrator.getDrudeTemperature(), 100*kelvin)
integrator.setDrudeTemperature(50)
self.assertEqual(integrator.getDrudeTemperature(), 50*kelvin)
self.assertEqual(integrator.getStepSize(), 0.1*femtosecond)
integrator.setStepSize(0.0005)
self.assertEqual(integrator.getStepSize(), 0.0005*picosecond)
self.assertEqual(integrator.getMaxDrudeDistance(), 0.02*nanometer)
integrator.setMaxDrudeDistance(0.05)
self.assertEqual(integrator.getMaxDrudeDistance(), 0.05*nanometer)
def testCustomIntegrator(self):
""" Tests the CustomIntegrator API features """
integrator = CustomIntegrator(1.0)
self.assertEqual(integrator.getStepSize(), 1.0*picosecond)
self.assertEqual(integrator.getConstraintTolerance(), 1e-5)
integrator.setConstraintTolerance(1e-6)
self.assertEqual(integrator.getConstraintTolerance(), 1e-6)
integrator = CustomIntegrator(1.0*femtoseconds)
self.assertEqual(integrator.getStepSize(), 1.0*femtoseconds)
def testQTBIntegrator(self):
""" Tests the LangevinIntegrator API features """
integrator = QTBIntegrator(300, 0.1, 0.001)
self.assertEqual(integrator.getTemperature(), 300*kelvin)
self.assertEqual(integrator.getFriction(), 0.1/picosecond)
self.assertEqual(integrator.getStepSize(), 0.001*picosecond)
integrator = QTBIntegrator(300*kelvin, 0.1/microsecond, 1*femtosecond)
self.assertEqual(integrator.getTemperature(), 300*kelvin)
self.assertAlmostEqualUnit(integrator.getFriction(), 0.1/microsecond)
self.assertEqual(integrator.getStepSize(), 1*femtosecond)
integrator.setSegmentLength(0.5)
self.assertEqual(integrator.getSegmentLength(), 0.5*picosecond)
integrator.setSegmentLength(100*femtosecond)
self.assertEqual(integrator.getSegmentLength(), 0.1*picosecond)
if __name__ == '__main__':
unittest.main()
Reference in New Issue View Git Blame Copy Permalink