box stabilisation

bin/openmmMD.py

@@ -111,7 +111,7 @@ def setup_system(modeller, forcefield, solvmol: str, no_restraints: bool):
     return system
 
 def setup_simulation(modeller, system):
-    integrator = mm.LangevinMiddleIntegrator(300*kelvin, 1/picosecond, 0.0001*picoseconds)
+    integrator = mm.LangevinMiddleIntegrator(300*kelvin, 1/picosecond, 0.001*picoseconds)
     simulation = app.Simulation(modeller.topology, system, integrator)
     simulation.context.setPositions(modeller.positions)
     return simulation
@@ -134,6 +134,7 @@ def energy_minimization(modeller):
 def md_nvt(simulation, csvname: str, totalsteps: int, reprate: int, pdbname):
     init_state = simulation.context.getState(getEnergy=True, getPositions=True)
     init_pe = init_state.getPotentialEnergy().value_in_unit(kilocalories_per_mole)
+    simulation.context.setVelocitiesToTemperature(300*unit.kelvin)
     simulation.reporters.append(app.PDBReporter(pdbname, reprate))
     simulation.reporters.append(app.StateDataReporter(stdout, reprate, step=True, potentialEnergy=True, temperature=True, volume=True))
     prepdf = {'Step':[], 'Potential Energy (kJ/mole)':[], 'Temperature (K)':[], 'Box Volume (nm^3)':[]}