From 9b11039b608aa94b7e6021970e8427515892513a Mon Sep 17 00:00:00 2001 From: Lester Hedges Date: Mon, 23 Aug 2021 11:56:48 +0100 Subject: [PATCH] Fixed somd-freenrg OpenMM random number generator seeding bug. --- CHANGELOG | 5 ++++- wrapper/Tools/OpenMMMD.py | 46 +++++++++++++++++++-------------------- 2 files changed, 26 insertions(+), 25 deletions(-) diff --git a/CHANGELOG b/CHANGELOG index 46946296d..5627b5580 100644 --- a/CHANGELOG +++ b/CHANGELOG @@ -3,7 +3,7 @@ SIRE changelog devel branch: - [2021.1.0] Feb 2021 - Added support for multiple combining rules in SOMD + [2021.1.0] Aug 2021 - Added support for multiple combining rules in SOMD (@SofiaBariami). Added support for triclinic simulation boxes. Convert Ryckaert-Bellememans form dihedral functions from GROMACS to Fourier series to allow conversion to AMBER format. @@ -24,6 +24,9 @@ devel branch: conda version of libcpuid. Added support for generating PDB CONECT records from a Sire.Mol.Connectivity object. Fixed issue with PMEMD skipping torsions with zero periodicity. + Fixed random number seeding bug in somd-freenrg, which + resulted in the OpenMM generator being seeded with the same + seed for each cycle of the simulation. [2020.1.0] July 2020 - Fixed bug in WaterView program to ensure that a molecule is extracted from the returned list. Stable sorting diff --git a/wrapper/Tools/OpenMMMD.py b/wrapper/Tools/OpenMMMD.py index e6c796091..428479567 100644 --- a/wrapper/Tools/OpenMMMD.py +++ b/wrapper/Tools/OpenMMMD.py @@ -1,4 +1,3 @@ - #################################################################################################### # # # RUN SCRIPT to perform an MD simulation in Sire with OpenMM # @@ -87,8 +86,10 @@ nmoves = Parameter("nmoves", 1000, """Number of Molecular Dynamics moves to be performed during the simulation.""") -random_seed = Parameter("random seed", None, """Random number seed. Set this if you - want to have reproducible simulations.""") +debug_seed = Parameter("debug seed", 0, """Debugging seed number seed. Set this if you + want to reproduce a single cycle. Don't use this seed for production simulations + since the same seed will be used for all cycles! A value of zero means that a unique + seed will be generated for each cycle.""") ncycles = Parameter("ncycles", 1, """The number of MD cycles. The total elapsed time will be nmoves*ncycles*timestep""") @@ -455,7 +456,7 @@ def setupForcefields(system, space): return system -def setupMoves(system, random_seed, GPUS): +def setupMoves(system, debug_seed, GPUS): print("Setting up moves...") @@ -507,11 +508,11 @@ def setupMoves(system, random_seed, GPUS): moves = WeightedMoves() moves.add(mdmove, 1) - if (not random_seed): - random_seed = RanGenerator().randInt(100000, 1000000) - print("Generated random seed number %d " % random_seed) + if debug_seed != 0: + debug_seed = RanGenerator().randInt(100000, 1000000) + print("Generated debugging seed number %d " % debug_seed) - moves.setGenerator(RanGenerator(random_seed)) + moves.setGenerator(RanGenerator(debug_seed)) return moves @@ -1154,7 +1155,7 @@ def setupForceFieldsFreeEnergy(system, space): return system -def setupMovesFreeEnergy(system, random_seed, GPUS, lam_val): +def setupMovesFreeEnergy(system, debug_seed, GPUS, lam_val): print ("Setting up moves...") @@ -1165,7 +1166,7 @@ def setupMovesFreeEnergy(system, random_seed, GPUS, lam_val): solute_fromdummy = system[MGName("solute_ref_fromdummy")] Integrator_OpenMM = OpenMMFrEnergyST(molecules, solute, solute_hard, solute_todummy, solute_fromdummy) - Integrator_OpenMM.setRandomSeed(random_seed) + Integrator_OpenMM.setRandomSeed(debug_seed) Integrator_OpenMM.setIntegrator(integrator_type.val) Integrator_OpenMM.setFriction(inverse_friction.val) # Only meaningful for Langevin/Brownian integrators Integrator_OpenMM.setPlatform(platform.val) @@ -1205,7 +1206,7 @@ def setupMovesFreeEnergy(system, random_seed, GPUS, lam_val): #This calls the OpenMMFrEnergyST initialise function Integrator_OpenMM.initialise() velocity_generator = MaxwellBoltzmann(temperature.val) - velocity_generator.setGenerator(RanGenerator(random_seed)) + velocity_generator.setGenerator(RanGenerator(debug_seed)) mdmove = MolecularDynamics(molecules, Integrator_OpenMM, timestep.val, {"velocity generator":velocity_generator}) @@ -1215,12 +1216,11 @@ def setupMovesFreeEnergy(system, random_seed, GPUS, lam_val): moves = WeightedMoves() moves.add(mdmove, 1) - if (not random_seed): - random_seed = RanGenerator().randInt(100000, 1000000) - - print("Generated random seed number %d " % random_seed) + if debug_seed != 0: + debug_seed = RanGenerator().randInt(100000, 1000000) + print("Generated debugging seed number %d " % debug_seed) - moves.setGenerator(RanGenerator(random_seed)) + moves.setGenerator(RanGenerator(debug_seed)) return moves @@ -1429,12 +1429,11 @@ def run(): system = setupForcefields(system, space) - if random_seed.val: - ranseed = random_seed.val + if debug_seed.val == 0: + ranseed = debug_seed.val else: ranseed = RanGenerator().randInt(100000, 1000000) - - print("Setting up the simulation with random seed %s" % ranseed) + print("Setting up the simulation with debugging seed %s" % ranseed) moves = setupMoves(system, ranseed, gpu.val) @@ -1587,12 +1586,11 @@ def runFreeNrg(): system = setupForceFieldsFreeEnergy(system, space) - if random_seed.val: - ranseed = random_seed.val + if debug_seed.val == 0: + ranseed = debug_seed.val else: ranseed = RanGenerator().randInt(100000, 1000000) - - print("Setting up the simulation with random seed %s" % ranseed) + print("Setting up the simulation with debugging seed %s" % ranseed) moves = setupMovesFreeEnergy(system, ranseed, gpu.val, lambda_val.val)