Test case for charmm permissive read

wrappers/python/tests/TestCharmmFiles.py

@@ -142,6 +142,43 @@ class TestCharmmFiles(unittest.TestCase):
                 diff = norm(f1-f2)
                 self.assertTrue(diff < 0.1 or diff/norm(f1) < 1e-4)
 
+    def test_PermissiveRead(self):
+        """Compare permissive and strict reading of Charmm parameters"""
+
+        psf = CharmmPsfFile('systems/5dhfr_cube.psf')
+        pdb = PDBFile('systems/5dhfr_cube.pdb')
+
+        params_strict     = CharmmParameterSet('systems/par_all22_prot_with_mass.inp')
+        params_permissive = CharmmParameterSet('systems/par_all22_prot.inp', permissive=True)
+        # Box dimensions (found from bounding box)
+        psf.setBox(62.23*angstroms, 62.23*angstroms, 62.23*angstroms)
+
+        # Turn off charges so we only test the Lennard-Jones energies
+        for a in psf.atom_list:
+            a.charge = 0.0
+
+        # Now compute the full energy
+        plat = Platform.getPlatformByName('Reference')
+
+        system_strict     = psf.createSystem(params_strict    , nonbondedMethod=PME,
+                                  nonbondedCutoff=8*angstroms)
+        system_permissive = psf.createSystem(params_permissive, nonbondedMethod=PME,
+                                  nonbondedCutoff=8*angstroms)
+
+        con_strict     = Context(system_strict    , VerletIntegrator(2*femtoseconds), plat)
+        con_permissive = Context(system_permissive, VerletIntegrator(2*femtoseconds), plat)
+
+        con_strict.setPositions(pdb.positions)
+        con_permissive.setPositions(pdb.positions)
+
+        state_strict     = con_strict.getState(getEnergy=True, enforcePeriodicBox=True)
+        state_permissive = con_permissive.getState(getEnergy=True, enforcePeriodicBox=True)
+
+        ene_strict     = state_strict.getPotentialEnergy().value_in_unit(kilocalories_per_mole)
+        ene_permissive = state_permissive.getPotentialEnergy().value_in_unit(kilocalories_per_mole)
+        self.assertAlmostEqual(ene_strict, ene_permissive, delta=0.00001)
+
+
 if __name__ == '__main__':
     unittest.main()