Brief Description
I found that if I add a log reporter to unbiased alanine-dipeptide simulation, the histogram logger would record one time more xi values than the time steps.
Environment
I followed the set up in installation note in midway3. I ran my job in sinteractive mode in gpu partition.
Code to reproduce
(this code is adapted from example for spectral abf)
#!/usr/bin/env python3
"""
SpectralABF simulation of Alanine Dipeptide in vacuum with OpenMM and PySAGES.
Example command to run the simulation `python3 alanine-dipeptide.py --time-steps 1000`
For other supported commandline parameters, check `python3 alanine-dipeptide.py --help`
"""
# %%
import argparse
import sys
import time
import matplotlib.pyplot as plt
import numpy
import pysages
from pysages.approxfun import compute_mesh
from pysages.colvars import DihedralAngle
from pysages.methods import Unbiased, HistogramLogger
from pysages.utils import try_import
openmm = try_import("openmm", "simtk.openmm")
unit = try_import("openmm.unit", "simtk.unit")
app = try_import("openmm.app", "simtk.openmm.app")
# %%
pi = numpy.pi
kB = unit.BOLTZMANN_CONSTANT_kB * unit.AVOGADRO_CONSTANT_NA
kB = kB.value_in_unit(unit.kilojoules_per_mole / unit.kelvin)
T = 298.15 * unit.kelvin
dt = 2.0 * unit.femtoseconds
adp_pdb = "adp-vacuum.pdb"
# %%
def generate_simulation(pdb_filename=adp_pdb, T=T, dt=dt):
pdb = app.PDBFile(pdb_filename)
ff = app.ForceField("amber99sb.xml")
cutoff_distance = 1.0 * unit.nanometer
topology = pdb.topology
system = ff.createSystem(
topology, constraints=app.HBonds, nonbondedMethod=app.PME, nonbondedCutoff=cutoff_distance
)
# Set dispersion correction use.
forces = {}
for i in range(system.getNumForces()):
force = system.getForce(i)
forces[force.__class__.__name__] = force
forces["NonbondedForce"].setUseDispersionCorrection(True)
forces["NonbondedForce"].setEwaldErrorTolerance(1.0e-5)
positions = pdb.getPositions(asNumpy=True)
integrator = openmm.LangevinIntegrator(T, 1 / unit.picosecond, dt)
integrator.setRandomNumberSeed(42)
# platform = openmm.Platform.getPlatformByName(platform)
# simulation = app.Simulation(topology, system, integrator, platform)
simulation = app.Simulation(topology, system, integrator)
simulation.context.setPositions(positions)
simulation.minimizeEnergy()
simulation.reporters.append(
app.StateDataReporter('log',
1,
step=True,
time=True,
potentialEnergy=True,
kineticEnergy=True,
totalEnergy=True,
temperature=True))
return simulation
# %%
def get_args(argv):
available_args = [
("time-steps", "t", int, 5e5, "Number of simulation steps"),
]
parser = argparse.ArgumentParser(description="Example script to run Spectral ABF")
for (name, short, T, val, doc) in available_args:
parser.add_argument("--" + name, "-" + short, type=T, default=T(val), help=doc)
return parser.parse_args(argv)
# %%
def main(argv=[]):
args = get_args(argv)
cvs = [DihedralAngle([4, 6, 8, 14])]
# Method
method = Unbiased(cvs)
tic = time.perf_counter()
callback = HistogramLogger(1)
run_result = pysages.run(method, generate_simulation, args.time_steps, callback)
toc = time.perf_counter()
print(f"Completed the simulation in {toc - tic:0.4f} seconds.")
np.savetxt("dihedral.txt", run_result.callbacks[0].data)
# %%
if __name__ == "__main__":
main(sys.argv[1:])I changed phi-psi to one dihedral, changed abf to unbiased, and added a histogramLogger and a openmm log reporter.
And I execute
python alanine-dipeptide.py --time-steps 10
this would lead to a dihedral.txt file with 20 lines, and a log file with 11 lines (the first line is the head).
dihedral.txt:
-1.356931805610656738e+00
-1.358228564262390137e+00
-1.358228564262390137e+00
-1.360592961311340332e+00
-1.360592961311340332e+00
-1.361593604087829590e+00
-1.361593604087829590e+00
-1.363622307777404785e+00
-1.363622307777404785e+00
-1.365397930145263672e+00
-1.365397930145263672e+00
-1.365929484367370605e+00
-1.365929484367370605e+00
-1.368099331855773926e+00
-1.368099331855773926e+00
-1.369502067565917969e+00
-1.369502067565917969e+00
-1.370507478713989258e+00
-1.370507478713989258e+00
-1.371272802352905273e+00
log:
#"Step","Time (ps)","Potential Energy (kJ/mole)","Kinetic Energy (kJ/mole)","Total Energy (kJ/mole)","Temperature (K)"
1,0.002,-90.7515703562076,0.14926064088649582,-90.6023097153211,0.7039972076315866
2,0.004,-90.66084159644197,0.31577375944470987,-90.34506783699726,1.489366812121981
3,0.006,-90.52912772925447,0.4592449535266496,-90.06988277572782,2.166057729495585
4,0.008,-90.40165580542634,0.5241025778814219,-89.87755322754492,2.4719627970887195
5,0.01,-90.16859306128572,0.8659084612736478,-89.30268460001207,4.0841117603452854
6,0.012,-89.96650565894197,1.1798341972753406,-88.78667146166663,5.56476225357841
7,0.014,-89.90150321753572,1.4788492678198963,-88.42265394971582,6.975085654663139
8,0.016,-89.8797441843326,1.8018517950549722,-88.07789238927762,8.498547405067642
9,0.018000000000000002,-89.8254228952701,2.0839733593165874,-87.7414495359535,9.829191520443416
10,0.020000000000000004,-89.55366386206697,1.9868355998769403,-87.56682826219003,9.371035163918704
you can see that many consecutive lines are repeated in the dihedral.txt. If I print out the counter at the end:
print(run_result.callbacks[0].counter), this would print 20, but the total nsteps is only 10.
If I comment out the openmm log reporter, everything would behave normally.
dihedral.txt:
-1.356931805610656738e+00
-1.358228564262390137e+00
-1.360592961311340332e+00
-1.361593604087829590e+00
-1.363622307777404785e+00
-1.365397930145263672e+00
-1.365929484367370605e+00
-1.368099331855773926e+00
-1.369502067565917969e+00
-1.370507478713989258e+00
log:
#"Step","Time (ps)","Potential Energy (kJ/mole)","Kinetic Energy (kJ/mole)","Total Energy (kJ/mole)","Temperature (K)"
1,0.002,-90.7515703562076,0.14926064088649582,-90.6023097153211,0.7039972076315866
2,0.004,-90.66084159644197,0.31577375944470987,-90.34506783699726,1.489366812121981
3,0.006,-90.52912772925447,0.4592449535266496,-90.06988277572782,2.166057729495585
4,0.008,-90.40165580542634,0.5241025778814219,-89.87755322754492,2.4719627970887195
5,0.01,-90.16859306128572,0.8659084612736478,-89.30268460001207,4.0841117603452854
6,0.012,-89.96650565894197,1.1798341972753406,-88.78667146166663,5.56476225357841
7,0.014,-89.90150321753572,1.4788492678198963,-88.42265394971582,6.975085654663139
8,0.016,-89.8797441843326,1.8018517950549722,-88.07789238927762,8.498547405067642
9,0.018000000000000002,-89.8254228952701,2.0839733593165874,-87.7414495359535,9.829191520443416
10,0.020000000000000004,-89.55366386206697,1.9868355998769403,-87.56682826219003,9.371035163918704
![pre-commit-ci[bot] avatar](https://avatars.githubusercontent.com/in/68672?v=4 "pre-commit-ci[bot]")
React
Typescript
Django
Laravel
Facebook
Microsoft
Google
Alibaba
D3
Tencent