Merge pull request #466 from choderalab/update_sigma

Improve alchemical path for atoms being created/destroyed

perses/annihilation/new_relative.py

@@ -1,3 +1,4 @@
+import sys
 import simtk.openmm as openmm
 import simtk.openmm.app as app
 import simtk.unit as unit
@@ -1243,7 +1244,7 @@ def handle_nonbonded(self):
                 [charge, sigma, epsilon] = old_system_nonbonded_force.getParticleParameters(old_index)
 
                 #add the particle to the hybrid custom sterics and electrostatics.
-                self._hybrid_system_forces['core_sterics_force'].addParticle([sigma, epsilon, 1.0, 0.0])
+                self._hybrid_system_forces['core_sterics_force'].addParticle([sigma, epsilon, sigma, 0.0])
                 self._hybrid_system_forces['core_electrostatics_force'].addParticle([charge, 0.0])
 
                 #Add the particle to the regular nonbonded force as required, but zero out interaction
@@ -1256,7 +1257,7 @@ def handle_nonbonded(self):
                 [charge, sigma, epsilon] = new_system_nonbonded_force.getParticleParameters(new_index)
 
                 #add the particle to the hybrid custom sterics and electrostatics
-                self._hybrid_system_forces['core_sterics_force'].addParticle([1.0, 0.0, sigma, epsilon])
+                self._hybrid_system_forces['core_sterics_force'].addParticle([sigma, 0.0, sigma, epsilon])
                 self._hybrid_system_forces['core_electrostatics_force'].addParticle([0.0, charge])
 
                 #Add the particle to the regular nonbonded force as required, but zero out interaction

perses/dispersed/relative_setup.py

@@ -142,7 +142,7 @@ def __init__(self, ligand_file, old_ligand_index, new_ligand_index, forcefield_f
                 barostat = openmm.MonteCarloBarostat(self._pressure, self._temperature, self._barostat_period)
             else:
                 barostat = None
-            self._system_generator = SystemGenerator(forcefield_files, barostat=barostat, forcefield_kwargs={'nonbondedMethod' : self._nonbonded_method})
+            self._system_generator = SystemGenerator(forcefield_files, barostat=barostat, forcefield_kwargs={'nonbondedMethod' : self._nonbonded_method, 'constraints': app.HBonds})
         else:
             self._system_generator = SystemGenerator(forcefield_files)