Add python interface tests

openmmapi/src/NoseHooverIntegrator.cpp

@@ -250,7 +250,9 @@ double NoseHooverIntegrator::computeKineticEnergy() {
 double NoseHooverIntegrator::computeHeatBathEnergy() {
     double energy = 0;
     for(auto &nhc : noseHooverChains) {
-        energy += nhcKernel.getAs<NoseHooverChainKernel>().computeHeatBathEnergy(*context, nhc);
+        if (context && (nhc.getDefaultNumDegreesOfFreedom() > 0)) {
+            energy += nhcKernel.getAs<NoseHooverChainKernel>().computeHeatBathEnergy(*context, nhc);
+        }
     }
     return energy;
 }

wrappers/python/src/swig_doxygen/swigInputConfig.py

@@ -475,7 +475,7 @@ UNITS = {
 ("MonteCarloMembraneBarostat", "MonteCarloMembraneBarostat") : (None, ("unit.bar", "unit.bar*unit.nanometer", "unit.kelvin", None, None, None)),
 ("MonteCarloMembraneBarostat", "getXYMode") : (None, ()),
 ("MonteCarloMembraneBarostat", "getZMode") : (None, ()),
-("DrudeLangevinIntegrator", "getDrudeFriction") : ("1/unit.picosecond", ()),
+("DrudeLangevinIntegrator", "getDrudeFriction") : ("unit.picosecond**-1", ()),
 ("DrudeSCFIntegrator", "getMinimizationErrorTolerance") : ("unit.kilojoules_per_mole/unit.nanometer", ()),
 ("RPMDIntegrator", "getContractions") : (None, ()),
 ("RPMDIntegrator", "getTotalEnergy") : ("unit.kilojoules_per_mole", ()),
@@ -485,14 +485,13 @@ UNITS = {
 #("NoseHooverChain", "getThermostatedPairs") : (None, ()),
 #("NoseHooverChain", "getThermostatedAtoms") : (None, ()),
 #("NoseHooverChain", "getDefaultYoshidaSuzukiWeights") : (None, ()),
-#("NoseHooverIntegrator", "NoseHooverIntegrator") : (None, ("unit.kelvin", "1/unit.picosecond", "unit,picosecond", None, None, None)),
 ("NoseHooverIntegrator", "setTemperature") : (None, ("unit.kelvin", None)),
 ("NoseHooverIntegrator", "setRelativeTemperature") : (None, ("unit.kelvin", None) ),
-("NoseHooverIntegrator", "setCollisionFrequency") : (None, ("1/unit.picosecond", None)),
-("NoseHooverIntegrator", "setRelativeCollisionFrequency") : (None, ("1/unit.picosecond", None)),
+("NoseHooverIntegrator", "setCollisionFrequency") : (None, ("unit.picosecond**-1", None)),
+("NoseHooverIntegrator", "setRelativeCollisionFrequency") : (None, ("unit.picosecond**-1", None)),
 ("NoseHooverIntegrator", "computeHeatBathEnergy") : ( "unit.kilojoules_per_mole", ()),
-("NoseHooverIntegrator", "addThermostat"): (None, ("unit.kelvin", "1/unit.picosecond", None, None, None)),
+("NoseHooverIntegrator", "addThermostat"): (None, ("unit.kelvin", "unit.picosecond**-1", None, None, None)),
 ("NoseHooverIntegrator", "addSubsystemThermostat"):
-    (None, (None, None, "unit.kelvin", "1/unit.picosecond", "unit.kelvin", "1/unit.picosecond", None, None, None))
+    (None, (None, None, "unit.kelvin", "unit.picosecond**-1", "unit.kelvin", "unit.picosecond**-1", None, None, None))
 }
 

wrappers/python/tests/TestAPIUnits.py

@@ -682,7 +682,7 @@ class TestAPIUnits(unittest.TestCase):
         force.setSolventDielectric(80)
         self.assertEqual(force.getSolventDielectric(), 80)
 
-        self.assertEqual(force.getSurfaceAreaFactor(),
+        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(),
@@ -1176,8 +1176,47 @@ class TestAPIUnits(unittest.TestCase):
         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 consructor 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):

wrappers/python/tests/TestIntegrators.py

 import unittest
+import warnings
 import tempfile
 from datetime import datetime, timedelta
 from simtk.openmm import *
@@ -120,5 +121,36 @@ class TestIntegrators(unittest.TestCase):
             context.setPositions(pdb.positions)
             integrator.step(10)
 
+    def testNoseHooverIntegrator(self):
+        """Test partial thermostating in the NoseHooverIntegrator (only API)"""
+        pdb = PDBFile('systems/alanine-dipeptide-explicit.pdb')
+        ff = ForceField('amber99sbildn.xml', 'tip3p.xml')
+        system = ff.createSystem(pdb.topology, cutoffMethod=PME)
+
+        integrator = NoseHooverIntegrator(1.0*femtosecond) 
+        integrator.addSubsystemThermostat(list(range(5)), [], 200*kelvin, 1/picosecond, 200*kelvin, 1/picosecond, 3,3,3)
+        con = Context(system, integrator)
+        con.setPositions(pdb.positions)
+
+        integrator.step(5)
+        self.assertNotEqual(integrator.computeHeatBathEnergy(), 0.0*kilojoule_per_mole)
+
+    def testDrudeNoseHooverIntegrator(self):
+        """Test the DrudeNoseHooverIntegrator"""
+        warnings.filterwarnings('ignore', category=CharmmPSFWarning)
+        psf = CharmmPsfFile('systems/ala3_solv_drude.psf')
+        crd = CharmmCrdFile('systems/ala3_solv_drude.crd')
+        params = CharmmParameterSet('systems/toppar_drude_master_protein_2013e.str')
+        # Box dimensions (cubic box)
+        psf.setBox(33.2*angstroms, 33.2*angstroms, 33.2*angstroms)
+
+        system = psf.createSystem(params, nonbondedMethod=PME, ewaldErrorTolerance=0.0005)
+        integrator = DrudeNoseHooverIntegrator(300*kelvin, 1.0/picosecond, 1*kelvin, 10/picosecond, 0.001*picoseconds)
+        con = Context(system, integrator)
+        con.setPositions(crd.positions)
+
+        integrator.step(5)
+        self.assertNotEqual(integrator.computeHeatBathEnergy(), 0.0*kilojoule_per_mole)
+
 if __name__ == '__main__':
     unittest.main()