simulations always blow up

[juanmit]

Hi,

I am trying to do an openmm test internally on FAH, and the NVT simulation could run successfully on my workstation.
However, it keeps blowing-up when I grab it from FAH internally. I think there should be something wrong with the *.xml files. Could anyone help me out? Thanks.

--Hongbin

from simtk.openmm.app import *
from simtk.openmm import *
from simtk.unit import *

pdb = PDBFile('alanine.pdb')

forcefield = ForceField('amber99sb.xml', 'tip3p.xml')
system = forcefield.createSystem(pdb.topology, nonbondedMethod=PME, nonbondedCutoff=1*nanometer, constraints=HBonds)
integrator = LangevinIntegrator(300*kelvin, 1/picosecond, 0.002*picoseconds)

simulation = Simulation(topology, system, integrator)
simulation.context.setPositions(pdb.positions)
simulation.context.setVelocitiesToTemperature(300*kelvin)

with open('system.xml', 'w') as f:
    system_xml = XmlSerializer.serialize(system)
    f.write(system_xml)

with open('integrator.xml', 'w') as f:
    integrator_xml = XmlSerializer.serialize(integrator)
    f.write(integrator_xml)

state = simulation.context.getState(getPositions=True, getVelocities=True, getForces=True, getEnergy=True, getParameters=True, enforcePeriodicBox=True)

with open('state.xml', 'w') as f:
    f.write(XmlSerializer.serialize(state))
    print('saved state.xml')

[peastman]

Can you post all the files needed to reproduce this? Also, can you give more details about exactly what you're seeing? How quickly does it blow up? Does it happen on the very first step, or does it take a while? And does it always take the same amount of time, or is it variable? Is the energy of the initial conformation reasonable? What platform are you using? What hardware? What version of OpenMM?

[juanmit]

Hi Peter, Thanks for replying.

Here is all info on my workstation openmm=7.0.1 gpu_drives=375.26 GTX980 ubuntu=15.04

And attached includes all codes&results of the system. The starting structure was pre-equilibrated with gromacs. Please let me know if you need anything else.

--Hongbin

result.tar.gz

[peastman]

Could you answer the other questions from my previous post too? I need to know what to look for.

[juanmit]

Hi Peter,

1. Can you post all the files needed to reproduce this? I have put all the files in the tar.gz.(attached yesterday) There will be two folders(naming differently) consists of codes/results from FAH internal tests and codes/results from local runs.

2. can you give more details about exactly what you're seeing? The water diffuses out of box. See the attached figure.

3. How quickly does it blow up? Does it happen on the very first step, or does it take a while? And does it always take the same amount of time, or is it variable? Is the energy of the initial conformation reasonable? For the runs directly on the workstation, everything works well. The internal tests on FAH can run but the water molecules start to diffuse to very far at the very beginning. The initial energy = ~-149566.167185 kJ/mol, which I think it is a good starting point, and I have equilibrated everything using Gromacs beforehand.

4.What platform are you using? What hardware? What version of OpenMM? openmm=7.0.1 gpu_drives=375.26 GTX980 ubuntu=15.04

I am new to openmm21 on FAH, and I think this error comes from the system setup(*xml files). Hopefully you can help me out. Thanks a lot.

--Hongbin

[peastman]

Wait, what do you mean when you say it "blows up"? Do you just mean the water diffuses away from the origin? That may just mean nothing more than that you need to reimage the molecules to put them into a single periodic box (for example, with Trajectory.image_molecules() in MDTraj).

[juanmit]

top = GromacsTopFile('topol.top', periodicBoxVectors=gro.getPeriodicBoxVectors())

I thought this line of codes defines the PBC in the simulation. And also, when I run the simulation directly on my workstation, there were no pbc problems with same codes(all the molecules stay in the box). Are there codes I am missing in the project.xml that defines PBC? Thanks for helping.

[peastman]

Period boundary conditions refer to how forces are computed, not to where particles can go. When the forces are computed, it considers each particle to represent an infinite grid of identical particles repeating through all of space. You can still place a particle anywhere, not just within any one arbitrarily chosen box. But as far as the force is concerned, it doesn't matter. A particle at position x and one at x+box_width are both considered to represent exactly the same infinite grid of particles, and the forces on them are the same. But for other purposes it can make a big difference. For example, if you're measuring the diffusion rate of water molecules, you need to know which of those positions it's now out, since it had to travel different distances to get to them.

[juanmit]

Thanks. I know what you are trying to say about PBC. What I am confused is that with the same setup, the outputs from FAH-internal and LOCAL are different. If that is the system setup issue, I hope to find out where I went wrong and edit my codes.

[peastman]

The DCDReporter is wrapping all particle coordinates into a single box when it saves them. If you add the argument enforcePeriodicBox=False to its constructor, you'll get the same result as on FAH.

[cxhernandez]

The DCDReporter is wrapping all particle coordinates into a single box when it saves them. If you add the argument enforcePeriodicBox=False to its constructor, you'll get the same result as on FAH.

Not sure if this is a question better suited for @jcoffland, but does the FAH core itself have an option for this?

[jcoffland]

@cxhernandez , no all simulation options are set on the WS in the project configuration.

[juanmit]

@peastman @cxhernandez I think everything works well now. I really appreciate your help on it. Thanks!!!! --Hongbin